Formation and Initiation of Erupting Flux Rope and Embedded Filament Driven by Photospheric Converging Motion

The MHD simulation presented here is performed in a rectangular Cartesian geometry, ignoring the curvature of the solar surface. The evolution of the system can be adequately described by a 2.5D thermodynamic MHD model including gravity, heat conduction, radiative cooling, background heating, viscosity, and resistivity (Priest 1982):

$$\begin{eqnarray}&&\displaystyle \frac{\partial \rho }{\partial t}+\triangledown \cdot (\rho {\boldsymbol{v}})=0,\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial (\rho {\boldsymbol{v}})}{\partial t}+\triangledown \cdot \left[\rho {\boldsymbol{vv}}+\left(p+\displaystyle \frac{{{\boldsymbol{B}}}^{2}}{8\pi }\right){\boldsymbol{I}}-\displaystyle \frac{{\boldsymbol{BB}}}{4\pi }\right]=\rho {\boldsymbol{g}}+\triangledown \cdot {\boldsymbol{\tau }},\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial e}{\partial t}+\triangledown \cdot \left[\left(e-\displaystyle \frac{{{\boldsymbol{B}}}^{2}}{8\pi }+p\right){\boldsymbol{v}}+\displaystyle \frac{c{\boldsymbol{E}}\times {\boldsymbol{B}}}{4\pi }\right]\\ &&\quad =\rho {\boldsymbol{g}}\cdot {\boldsymbol{v}}+\triangledown \cdot ({\boldsymbol{\kappa }}\cdot \triangledown T)-Q+H+\triangledown \cdot ({\boldsymbol{\tau }}\cdot {\boldsymbol{v}}),\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t}+\triangledown \cdot ({\boldsymbol{vB}}-{\boldsymbol{Bv}})=-c\triangledown \times (\eta {\boldsymbol{J}}).\end{eqnarray} \tag{ 4 }$$

Here $e=[(3/2)p+(1/2)\rho {{\boldsymbol{v}}}^{2}+{{\boldsymbol{B}}}^{2}/8\pi ]$ is the total plasma energy per unit volume, ${\boldsymbol{E}}=[\eta {\boldsymbol{J}}-({\boldsymbol{v}}\times {\boldsymbol{B}})/c]$ is the electric field, η is the resistivity, Q is the radiative loss rate, c is the light speed in vacuum, and H is the background heating rate. The viscous stress tensor ${\boldsymbol{\tau }}$ has the components ${\tau }_{{ij}}=\mu [(1/2)({\partial }_{j}{v}_{i}\,+{\partial }_{i}{v}_{j})-(1/3)({\delta }_{{ij}}\triangledown \cdot {\boldsymbol{v}})]$, where μ is the dynamic viscosity. All the equations are written in Gaussian units in this paper if not specified otherwise.

As per Xia et al. (2012), we assume a perfect gas law and consider fully ionized plasma with a 10:1 abundance of hydrogen and helium. The mass ratio of hydrogen and helium atoms is 1:4, thus the mass density ratio of hydrogen and helium is 10:4. The mass density ρ of our model atmosphere is

$$\begin{eqnarray}&&\rho =1.4{m}_{{\rm{p}}}{n}_{{\rm{H}}},\end{eqnarray} \tag{ 5 }$$

where m_p is the proton mass and n_H is the number density of hydrogen. Similarly, due to the charge neutrality, the number density ratio of ions and electrons is 11:12, so the gas pressure p is

$$\begin{eqnarray}&&p=2.3{n}_{{\rm{H}}}{k}_{{\rm{B}}}T=\displaystyle \frac{\rho {k}_{{\rm{B}}}T}{{\mu }_{w}{m}_{{\rm{H}}}},\end{eqnarray} \tag{ 6 }$$

where k_B is the Boltzmann constant and ${\mu }_{w}=1.4/2.3$ is the mean molecular weight.

It should be noted that we will solve Equations (1)–(6) in a nondimensionalized fashion. To nondimensionalize the equations, each variable is divided by its normalizing unit. The normalizing units of variables are given in Table 1, and the units of other derived variables are listed as well.

Table 1.  Normalization Units

Symbol	Quantity	Unit	Value
$x,y,z$	Length	L0	${10}^{9}\,\mathrm{cm}$
T	Temperature	T0	$1.0\times {10}^{6}\,{\rm{K}}$
ρ	Density	${\rho }_{0}$	$2.3417\times {10}^{-15}\,{\rm{g}}\,{\mathrm{cm}}^{-3}$
p	Pressure	${p}_{0}=({\rho }_{0}{k}_{{\rm{B}}}{T}_{0})/({\mu }_{w}{m}_{{\rm{H}}})$	$0.3175\,\mathrm{erg}\,{\mathrm{cm}}^{-3}$
e	Energy density	p0	$0.3175\,\mathrm{erg}\,{\mathrm{cm}}^{-3}$
${\boldsymbol{B}}$	Magnetic induction	${B}_{0}=\sqrt{4\pi {p}_{0}}$	$1.9976\,\ \mathrm{Gauss}$
${\boldsymbol{v}}$	Velocity	${v}_{0}={B}_{0}/\sqrt{4\pi {\rho }_{0}}$	$1.1645\times {10}^{7}\,\mathrm{cm}\,{{\rm{s}}}^{-1}$
t	Time	${t}_{0}={L}_{0}/{v}_{0}$	$85.8746\,{\rm{s}}$
η	Resistivity	$(4\pi {L}_{0}^{2})/({c}^{2}{t}_{0})$	$1.6282\times {10}^{-5}\,{\rm{s}}$
${\boldsymbol{E}}$	Electric field	$({B}_{0}{L}_{0})/({t}_{0}c)$	$7.759\times {10}^{-4}\,\mathrm{statvolt}\,{\mathrm{cm}}^{-1}$
${\boldsymbol{J}}$	Current density	$({B}_{0}c)/(4\pi {L}_{0})$	$4.7689\,\mathrm{statamp}\,{\mathrm{cm}}^{-2}$

Download table as:  ASCIITypeset image

We set the gravitational acceleration as ${\boldsymbol{g}}=-{g}_{0}\,{R}_{\odot }^{2}/{({R}_{\odot }+y)}^{2}{{\boldsymbol{e}}}_{y}$, with the solar radius as ${R}_{\odot }=6.961\times {10}^{10}\,\mathrm{cm}$ and the solar surface gravitational acceleration as ${g}_{0}=2.74\times {10}^{4}\,\mathrm{cm}\,{{\rm{s}}}^{-2}$. The anisotropic thermal conduction along the magnetic field lines and the radiative cooling terms are treated in the same way as Xia & Keppens (2016). We use the radiative cooling coefficient calculated by Colgan et al. (2008) with a lower cutoff temperature of $2\times {10}^{4}\,{\rm{K}}$. To obtain a self-consistent thermally structured corona, we augment this setup with a background heating rate decaying exponentially with height (Fang et al. 2013),

$$\begin{eqnarray}&&H={c}_{0}\exp \left(-\displaystyle \frac{y}{{L}_{\mathrm{Heat}}}\right),\end{eqnarray} \tag{ 7 }$$

with ${c}_{0}=1.0\times {10}^{-4}\,\mathrm{erg}\,{\mathrm{cm}}^{-3}\,{{\rm{s}}}^{-1}$ and ${L}_{\mathrm{Heat}}=50\,\mathrm{Mm}$. The coefficient c₀ is taken as the same value as that in Fang et al. (2013, 2015). Xia et al. (2011) assumed that the heating rate decreased exponentially with the distance away from the nearest footpoint in their one-dimensional magnetic loop, taking the heating rate at the footpoints as $3.0\,\times {10}^{-4}\,\mathrm{erg}\,{\mathrm{cm}}^{-3}\,{{\rm{s}}}^{-1}$. Karpen et al. (2003) assumed a small and spatially uniform background heating rate of $1.5\,\times {10}^{-4}\,\mathrm{erg}\,{\mathrm{cm}}^{-3}\,{{\rm{s}}}^{-1}$. Our setup for the background heating rate is thus similar to these previous works. We stress here that the background heating term compensates for the radiative loss and thermal conduction, maintaining the typical temperature and density of the solar atmosphere. Besides the background heating term, Xia et al. (2012) and Fang et al. (2013) used a localized extra heating term, of the order of ${10}^{-2}\,\mathrm{erg}\,{\mathrm{cm}}^{-3}\,{{\rm{s}}}^{-1}$, to make the chromospheric material evaporate into the corona. However, we did not utilize the localized extra heating term in our simulation and the only background heating term does not lead to any evaporation-condensation process in the corona in our simulation.

In our model, we set a relatively high resistivity in the chromosphere due to the collision effect of particles in this region, which will lead to an effective reconnection in the chromosphere to ensure that the reconnected magnetic flux is sufficient for the FR formation in the lower atmosphere. The dimensionless resistivity in the chromosphere is assumed to be 0.002, corresponding to a real resistivity of $2930.7\,{\rm{\Omega }}\,{\rm{m}}$, while it is set as 0.0002 in the corona, corresponding to a real resistivity of $293.1\,{\rm{\Omega }}\,{\rm{m}}$. The two resistivities are connected by a thin transition layer with a thickness of w_tr. We have

$$\begin{eqnarray}&&\eta ={\eta }_{\mathrm{chro}}+\displaystyle \frac{({\eta }_{\mathrm{cor}}-{\eta }_{\mathrm{chro}})}{2}\left[\tanh \left(\displaystyle \frac{y-{h}_{\mathrm{tr}}}{{w}_{\mathrm{tr}}}\right)+1\right],\end{eqnarray} \tag{ 8 }$$

where ${\eta }_{\mathrm{chro}}=0.002$, ${\eta }_{\mathrm{cor}}=0.0002$, ${h}_{\mathrm{tr}}=0.23$, and ${w}_{\mathrm{tr}}\,=0.05$, with the units listed in Table 1. Taking the local Alfvén speed, the resistivity, and the typical length scale of $L={10}^{9}\,\mathrm{cm}$, we estimate the Lundquist number ${S}_{L}=4\pi {{Lv}}_{A}/{c}^{2}\eta$ as ${10}^{4}\sim {10}^{5}$ in the corona and 10³ in the chromosphere.

Note that the photosphere and the chromosphere are made up of low-temperature weakly ionized plasmas. In a more self-consistent consideration, the deep photosphere should be treated by dynamo theory, and the generalized Ohm's law should be taken into account for the reconnection in the chromosphere (Somov 2013). Our one-fluid MHD model is a relatively simplified model. There is some artificiality in the choice of the form of the resistivity in this simulation, but the choice of the resistivity keeps the feature predicted by the Spitzer resistivity model: the resistivity in the chromosphere is higher than that in the corona. The Spitzer resistivity model is a more realistic model; however, it cannot be directly used in our simulation since the value predicted by the Spitzer resistivity model is much smaller than the numerical resistivity. The dimensionless viscosity μ is set as a constant at 0.0001 over the entire simulation box.

The computational box is located in the x–y Cartesian plane. The x-axis is parallel to the solar surface, while the y-axis is perpendicular to the photosphere. The computational domain is x_min ≤ x ≤ x_max and y_min ≤ y ≤ y_max, where x_min = −10.0, x_max = 0.0, and y_min = 0, y_max = 25.0, with the length unit at 10 Mm. Using a mirroring boundary condition at x = 0, we obtain a doubled box with an effective domain of −10 ≤ x ≤ 10 and 0 ≤ y ≤ 25.

Following Fang et al. (2013), the initial magnetic topology now adopts a linear force-free magnetic field, which represents a sheared magnetic arcade characterized by a constant angle of ${\theta }_{0}=30^\circ$ and a typical magnetic field strength of ${B}_{0}=20\,{\rm{G}}$ as follows:

$$\begin{eqnarray}&&{B}_{x}=-{B}_{0}\cos \left(\displaystyle \frac{\pi x}{{l}_{{\rm{B}}}}\right)\sin {\theta }_{0}\exp \left(-\displaystyle \frac{\pi y\sin {\theta }_{0}}{{l}_{{\rm{B}}}}\right),\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&{B}_{y}={B}_{0}\sin \left(\displaystyle \frac{\pi x}{{l}_{{\rm{B}}}}\right)\exp \left(-\displaystyle \frac{\pi y\sin {\theta }_{0}}{{l}_{{\rm{B}}}}\right),\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&{B}_{z}=-{B}_{0}\cos \left(\displaystyle \frac{\pi x}{{l}_{{\rm{B}}}}\right)\cos {\theta }_{0}\exp \left(-\displaystyle \frac{\pi y\sin {\theta }_{0}}{{l}_{{\rm{B}}}}\right),\end{eqnarray} \tag{ 11 }$$

where ${l}_{{\rm{B}}}=50\,\mathrm{Mm}$. The velocity is set as ${\boldsymbol{v}}=0$ at the beginning. For the initial thermal structure, we set a chromosphere with a temperature of ${10}^{4}\,{\rm{K}}$ below a height of 2.3 Mm and choose a vertically stratified temperature profile with a constant vertical thermal conduction flux (i.e., $\kappa \partial T/\partial y=2\times {10}^{5}\,\mathrm{erg}\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$, where $\kappa ={10}^{-6}{T}^{5/2}\,\mathrm{erg}\,{\mathrm{cm}}^{-1}\,{{\rm{s}}}^{-1}\,{{\rm{K}}}^{-3.5}$) above the height, following Xia et al. (2012). The initial density is then determined by assuming a hydrostatic atmosphere with a number density of $2.552\times {10}^{14}\,{\mathrm{cm}}^{-3}$ at the bottom of the boundary layer, i.e., solving the hydrostatic equation,

$$\begin{eqnarray}&&\displaystyle \frac{{dp}}{{dy}}=\rho g(y),\end{eqnarray} \tag{ 12 }$$

where $g{(y)=-{g}_{0}{R}_{\odot }^{2}/({R}_{\odot }+y)}^{2}$. The relation between pressure and density is given by the nondimensionalized ideal gas law,

$$\begin{eqnarray}&&\displaystyle \frac{p}{T}=\rho .\end{eqnarray} \tag{ 13 }$$

In order to solve the equations numerically, we introduce a new variable c_g:

$$\begin{eqnarray}&&{c}_{g}=\displaystyle \frac{{\rm{\Delta }}y({g}_{i,j-1}+{g}_{i,j}\,)}{4}=\displaystyle \frac{{\rm{\Delta }}y}{2}{g}_{i,j-\tfrac{1}{2}}.\end{eqnarray} \tag{ 14 }$$

Then the following finite-difference equation is adopted to obtain the pressure:

$$\begin{eqnarray}&&{p}_{i,j}=\displaystyle \frac{{p}_{i,j-1}+{c}_{g}\cdot {\rho }_{i,j-1}}{1-\tfrac{{c}_{g}}{{T}_{j}}}.\end{eqnarray} \tag{ 15 }$$

2.3.1. Top Boundary

Suppose that j_max is the vertical index of the cells in the uppermost layer in the physical region and that i, j represent the indices of a ghost cell in the upper ghost layer where $j\gt {j}_{\max }$ and the antisymmetric relation ${{\boldsymbol{v}}}_{i,j}=-{{\boldsymbol{v}}}_{i,2{j}_{\max }-j}$ is used for velocity at the top boundary, which ensures that no mass flux passes through the top boundary. The magnetic field at the top boundary is fixed as the same linear force-free field in the physical region initially. Continuous condition, which means that the gradient is kept at zero by copying the variable values from the edge of the mesh into the ghost cells, is applied for the temperature ${T}_{i,j}={T}_{i,{j}_{\max }}$ at the top boundary. By hydrostatic assumption, we can thus derive the exponential formula of pressure. The integral in the exponential index becomes a summation:

$$\begin{eqnarray}&&{p}_{i,j}={p}_{i,{j}_{\max }}\exp \left[\displaystyle \frac{-1}{{T}_{i,{j}_{\max }}}({g}_{i,{j}_{\max }}+{{\rm{\Sigma }}}_{k={j}_{\max }+1}^{k=j}{g}_{i,k}){\rm{\Delta }}y\right].\end{eqnarray} \tag{ 16 }$$

In the same way, the density at the top boundary is prescribed by

$$\begin{eqnarray}&&{\rho }_{i,j}={\rho }_{i,{j}_{\max }}\exp \left[\displaystyle \frac{-1}{{T}_{i,{j}_{\max }}}({g}_{i,{j}_{\max }}+{{\rm{\Sigma }}}_{k={j}_{\max }+1}^{k=j}{g}_{i,k}){\rm{\Delta }}y\right].\end{eqnarray} \tag{ 17 }$$

2.3.2. Bottom Boundary

The bottom boundary conditions are arranged as follows. Our boundary motion is chosen to converge toward the inversion line x = 0 at the following velocity:

$$\begin{eqnarray}&&{v}_{x}=\left\{\begin{array}{c}{-v}_{\max }\sin (2\pi x/{L}_{a}),| x| \leqslant {L}_{d}\\ 0,| x| \gt {L}_{d}\end{array}\right.,\end{eqnarray} \tag{ 18 }$$

where ${L}_{d}={L}_{a}/2$, ${L}_{a}=50\,\mathrm{Mm}$, ${v}_{\max }=7\,\mathrm{km}\,{{\rm{s}}}^{-1}$ is the amplitude of the driving velocity, and ${v}_{y}={v}_{z}=0$. The driving velocity mimics the photospheric converging motion. The driving photospheric velocity is switched off initially and the system first evolves through a relaxation procedure with the background heating, radiation, resistivity, viscosity, gravity, and thermal conduction until the system has reached a mechanical and thermodynamic equilibrium state, which is characterized by the temporal invariability, except for small fluctuations, of the total energy $E=\int [(1/2)\rho$ $({v}_{x}^{2}+{v}_{y}^{2}+{v}_{z}^{2})+(3/2)p+(1/8\pi )$ $({B}_{x}^{2}+{B}_{y}^{2}+{B}_{z}^{2})]{dxdy}$. The height of the chromosphere-corona transition region rises from 2.3 Mm to about 3 Mm after relaxation. After the relaxation procedure, the driving velocity is switched on. The instant when the driving velocity is switched on is counted as 0 s in this simulation. Then the amplitude of the driving velocity v_max increases linearly with time from 0 to ${v}_{\max }=7\,\mathrm{km}\,{{\rm{s}}}^{-1}$ within 500 s and keeps 7 km s⁻¹ from then on.The peak flow speed is about an order of magnitude higher than that observed in the photosphere, which is less than 1 km s⁻¹ (Kosovichev & Zharkova 2001). The value of the peak flow speed is reasonable since it is less than the sound speed of 130 km s⁻¹ and the Alfvén speed of 80 km s⁻¹ in the chromosphere. The boundary magnetic field is obtained by the use of the symmetric relation ${{\boldsymbol{B}}}_{i,j}={{\boldsymbol{B}}}_{i,2{j}_{\min }-j}$, where j_min represents the lowermost layer of the physical region and j represents the lower ghost layer, $j\lt {j}_{\min }$. Then we apply the divergence-free condition at the bottom boundary:

$$\begin{eqnarray}&&{B}_{y;i,j}={B}_{y;i,j+2}+\displaystyle \frac{{\rm{\Delta }}y}{{\rm{\Delta }}x}({B}_{x;i+1,j+1}-{B}_{x;i-1,j+1}).\end{eqnarray} \tag{ 19 }$$

The density and pressure are determined by assuming a hydrostatic atmosphere with a number density of $2.552\times {10}^{14}\,{\mathrm{cm}}^{-3}$ at the bottom of the boundary layer, which is close to the chromospheric density.

The governing equations are solved with a total variation diminishing (TVD) Lax-Friedrichs scheme and a three-step Runge–Kutta method in time integration. The grid refinement/coarsening criteria are based on the Lohner's error estimation (Lohner 1987), which evaluates the variations (gradients) of specified quantities by calculating their discrete second derivatives. The implementation of the AMR technique in the code is described in Keppens et al. (2012). In this simulation, the AMR criterion is based on density nonsmoothness. An effective resolution of 8192 × 20480 is attained by using 7 AMR levels, with an equivalent spatial resolution of 24.414 km in both directions. The effective numerical resistivity can be estimated as $v{\rm{\Delta }}{x}^{2}/L$, where v, Δx, and L are the typical values of velocity, grid spacing, and length (Arber et al. 2007). In dimensionless units, the velocity value is taken as the typical value of the Alfvén speed in the corona, about 10, and the length and grid spacing are taken as 1 and 0.0024, respectively. Then the effective numerical resistivity is estimated as 10⁻⁵, smaller than our physical resistivity. The discrete divergence of the magnetic field $\triangledown \cdot {\boldsymbol{B}}$ does not exactly vanish and may increase with time, leading to the wrong magnetic topologies (Dedner et al. 2002). One needs to handle the solenoidality constraint. As suggested by Powell et al. (1999) and Keppens et al. (2003), error-related source terms are added when writing

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t}+\triangledown \cdot ({\boldsymbol{vB}}-{\boldsymbol{Bv}})=-c\triangledown \times (\eta {\boldsymbol{J}})\\ &&\qquad -(\triangledown \cdot {\boldsymbol{B}}){\boldsymbol{v}}+\triangledown [{C}_{d}{({\rm{\Delta }}l)}^{2}\triangledown \cdot {\boldsymbol{B}}],\end{eqnarray} \tag{ 20 }$$

where ${({\rm{\Delta }}l)}^{2}={({\rm{\Delta }}x{\rm{\Delta }}y)}^{2}/[{({\rm{\Delta }}x)}^{2}+{({\rm{\Delta }}y)}^{2}]$ and ${C}_{d}=0.5/{t}_{0}$. To make the energy equation consistent with this treatment, two terms,  $-(\triangledown \cdot {\boldsymbol{B}}){\boldsymbol{B}}\cdot {\boldsymbol{v}}$ and ${\boldsymbol{B}}\cdot \triangledown [{C}_{d}{({\rm{\Delta }}l)}^{2}\triangledown \cdot {\boldsymbol{B}}]$, are added to the right of Equation (3). These ensure a controlled diffusion and advection of monopole errors in a manner that is easily ensured on an AMR grid.

The snapshots in Figure 1 capture the evolution of the density (left half) and the temperature (right half) with magnetic field lines (solid black lines) overlaid at eight instants of the evolution. Note that only a centered region of the simulation box corresponding to [−5, 5] × [0, 10] with a length unit of 10 Mm where the formation and eruption of the FR occurs is presented in Figure 1. Figure 2 illustrates the mass density distribution in the region [−0.6, 0.6] × [0, 5] with a length unit of 10 Mm at 16 instants depicting the prominence formation process and the current sheet evolution.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As shown in Figure 1(a), the temperature and the density distribution are stratified. From the chromosphere to the corona, the temperature increases with the height while the density decreases. The rapid increase of the temperature from 10⁴ to 1 MK occurs in the height interval between 2.5 and 12.5 Mm. The temperature increases from 1 to 2 MK as the height rises from 12.5 to 70 Mm. At a height of 100 Mm, the temperature reaches 2.4 MK. We now impose a converging motion to the footpoints of the magnetic arcades as described in Section 2.3.2. The mimicked photospheric converging motion drives the chromospheric material to flow to the magnetic reversal line. Owing to magnetic flux freezing, the magnetic fields of opposite polarities are carried to the magnetic polarity inversion line by the converging flow. As shown in Figure 1, the opposite-polarity magnetic fields collide in the polarity inversion line, and subsequently reconnection is driven in the chromosphere, leading to the formation of an FR at 1992.3 s. The newly formed FR is buried in the chromosphere, forming a concave-downward bulge in the chromosphere around the polarity inversion line. More magnetic fields of opposite polarity are sustainably advected to the magnetic polarity inversion line, and reconnection continues. The center of the FR rises from 2.8 Mm at 1992.3 s to 12.5 Mm at 2988.4 s, and the velocity increases from 0 to 0.01 Mm s⁻¹ at the same instants. Along with a sheared arcade, an FR is one of the magnetic structures capable of supporting a filament. During the process of the FR rising from the chromosphere into the corona, some cool and dense chromospheric plasma embedded in the FR is lifted into the corona, giving rise to the prominence formation process, which is discussed in Section 3.2. As shown in Figure 1, the lifted plasma has a comparatively lower temperature and higher density than the corona.

Afterward, the newly formed FR rises upward, evolving into the corresponding structure of a CME, namely the three-component morphology. At 2988.4 s, the FR detaches itself from the chromosphere, then a CME bubble is gradually formed and eruption is initiated. As the FR commences to rise, the magnetic field lines are severely stretched and a ray-like current sheet develops below the FR. It is widely accepted that the reconnection process in the current sheet accounts for the instantaneous energy release in the solar eruptive events, during which a large amount of free magnetic energy is converted into thermal and bulk energy of plasma. In our numerical experiment, since 3984.6 s the temperature of the thin current sheet is higher or at least equal to the ambient coronal plasma, much higher than in the chromospheric material.

As the FR rises, the vertical current sheet following the FR grows in length. Resistive MHD theory predicts that the current sheet becomes unstable once the diffusion region gets long enough (Furth et al. 1963). As shown in Figure 2, at ~4122 s a series of magnetic islands of different sizes appears on the current sheet, which marks the beginning of resistive unsteady reconnection. As indicated in Figure 2(o), the magnetic island above and connected with the flare loops may be identified as the cusp-like structure in observations, which was first discovered by Yohkoh (Masuda et al. 1994).

The reconnection outflows both moving upward and downward are heated to high temperatures, carrying kinetic energy and heated plasma out of the current sheet. The reconnection outflows in this simulation are always Alfvénic, which is consistent with the prediction by the Sweet–Parker model and the Petschek model that the reconnection outflows are of the same order of magnitude as that of the Alfvén speed in the inflow region. The upward-moving outflow hits the bottom of the FR and eventually enters into the FR bubble, heating the plasma there. Since in most of the domain the magnetic flux is frozen into the plasma and the thermal conduction perpendicular to the magnetic field lines is always neglected, both the heat and the mass entering into the FR would only move around the border of the FR rather than diffuse inward to the center of the FR or outward to the ambient coronal plasma. The mass and heat rapidly diffuse along the magnetic field lines around the border of the FR, forming a hot outer layer whose temperature reaches 14 MK at 4505.0 s, inconsistent with the observation. The overestimated temperature of the outer layer is attributed to our 2.5D setup, where there is no thermal conduction along the z-direction since the gradient of the temperature along the z-direction is zero. However, heat can be conducted along the axis of an FR on the Sun, transferring energy to the cooler chromosphere. The hot outer layer around the FR also appears in the two-dimensional simulation by Mei et al. (2012). The temperature and density of the cavity of the FR are 0.85 ~ 1.2 MK and $2.8\sim 3.5\times {10}^{-16}\,{\rm{g}}\,{\mathrm{cm}}^{-3}$, respectively, at 4505.0 s. Meanwhile, the momentum carried by the upward-moving reconnection outflow is transferred to the FR, accelerating the rising FR. This mechanical impulse distorts the prominence, seen in the bottom of the prominence at 3984.6 s. A bundle of flare loops below the current sheet structure appears at 4336.7 s. Since we mainly focus on the CME initiation and the beginning of CME acceleration, which are discussed in Section 3.3, the simulation stops at 4505.0 s. It should be noted that in this simulation the highest ratio of the speed of the CME front to the local Alfvén speed is ~0.15 at 4505.0 s. Thus the CME in this simulation cannot be identified as a fast CME but it may evolve to a fast CME if the computation continues.

The global evolution of the system can be clearly divided into three distinct stages. The system goes through an almost ideal evolution in the first stage that starts at 0 s and ends at 1992.3 s, the instant the reconnection is launched and the FR starts to form. The second stage lasts from 1992.3 to ~4122 s, during which the reconnection is driven by the converging flow and the FR grows in size. After ~4122 s, the ignition of the resistive instabilities marks the beginning of the third stage.

Prominences (or historically called filaments when observed on the solar disk) are one of the most intriguing structures in the solar atmosphere. Observations suggest that prominence plasmas are two orders of magnitude denser and cooler than their ambient coronal plasmas and have properties similar to those of chromospheric material. It was ultimately realized that the formation of a prominence requires more material than is available in the ambient corona and that the prominence plasma should originate from the chromosphere (Pikel'Ner 1971; Saito & Tandberg-Hanssen 1973). The physical processes of prominence formation and evolution still remain under active investigation. Four leading models have been proposed to account for prominence formation: the injection model (Zirker et al. 1998; Wang 1999; Chae 2003), the magnetothermal convection model (Berger et al. 2011), the evaporation-condensation model (Serio et al. 1981; Dahlburg et al. 1998; Karpen et al. 2003, 2006; Xia et al. 2012), and the levitation model (Rust & Kumar 1994; Deng et al. 2000).

As illustrated in Figure 2, our model shows how cool and dense chromospheric material is directly lifted upward into the solar corona by the newly formed FR through magnetic forces and is suspended by magnetic tension, which leads to the formation of a prominence. This scenario points to a model where cool prominence resides in the upward concave dip of the FR and is lifted up with rising magnetic fields (Rust & Kumar 1994), consistent with the levitation model, which does not rely on in situ condensation and formation. Similar formation processes have been reported from several observations (Lites 2005; Yelles Chaouche et al. 2012).

The opposite-polarity magnetic fields and the mass are advected to the polarity inversion line by the converging flow. At 2713.6 s, a bulk of mass has been raised into the corona by the upward-retracting field. Because the choice of a closed bottom boundary prohibits flux submergence and the model is 2.5D, the mass advected by the converging flow to the polarity inversion line can neither move downward through the bottom boundary nor escape horizontally along the polarity inversion line. It is likely that the mass collected at the polarity inversion line and raised by the upward-retracting field is overestimated in this work. The top region of the prominence is immersed into the rarefied corona while the lower region is still trapped in the chromosphere. From ~2885.4 s, our prominence begins to detach from the chromosphere and then slowly rises up at an almost constant speed ~13 km s⁻¹, followed by a long ray-like current sheet with a relatively higher density than the corona. We did not observe any uplifted chromospheric material draining along the axial magnetic fields back to the solar surface due to the translational invariance along the z-direction in the 2.5D setup, and the prominence mass in this simulation may be overestimated. After ~3700 s, the prominence and the FR rise in an accelerated fashion and the FR expands rapidly in size while the area of the prominence changes little. As shown in Figures 2(f)–(i), when the upward-moving reconnection outflow hits the bottom of the FR and enters into it, the lower part of the prominence is perturbed and deformed by the energetic reconnection jet. The flowing plasma spreads horizontally at a relatively fast speed, forming a low-density bubble-like structure at the lower part of the prominence (the tail of the prominence). The bubble-like structure is obvious in Figures 2(i)–(k) and then it is flattened gradually and fully merges into the prominence at ~4328.1 s, perturbing the prominence. After ~4190.7 s, the bottom of the prominence appears to be detached from the bottom of the FR as the FR expands rapidly. The reconnection jet continues to hit the bottom of the FR, forming a rather nonsteady and complex flow pattern.

Note that in this paper the FR is defined as a region where the value of the magnetic flux function A, which determines the magnetic field in the x–y plane by $\triangledown \times (A{{\boldsymbol{e}}}_{z})$, is larger than that of the uppermost X-point in the current sheet. So the bottom of the FR always connects the uppermost X-point in the current sheet. As shown in Figure 2(i), the bottom of the FR is elongated and connected to the X-point in the chromosphere, forming a tail-like structure. The tail-like structure is considered to be a part of the current sheet since it carries a much higher electric current than its ambient plasma and is directly connected to the X-point. The current sheet fragmentation, which starts from 4122 s, cuts off the lower part of the FR that is squeezed into the current sheet with other parts of the FR, creating magnetic islands. This shows that some parts of the FR can be fragmented into magnetic islands and move into the current sheet, carrying mass, magnetic flux, and electric current from the FR to the current sheet. Similar separation processes are also observed in the later evolution. We note that the density of the current sheet is higher than its ambient corona before 4122 s, the instant that resistive instabilities are ignited. The denser ray-like current sheet is fragmented into multiple plasmoids. These islands moving upward or downward carry the mass away from the current sheet. In comparing Figures 2(m) with (e), we note that the density of the current sheet, at 4414.0 s, is less than 3400.6 s. There are always interactions between the FR and the current sheet. On the one hand, the reconnection outflows, together with the magnetic islands produced in the current sheet, hit the FR, taking mass and magnetic flux from the current sheet to the FR. On the other hand, the lower part of the FR that extends into the current sheet can be separated from the FR by reconnection, forming magnetic islands with mass and magnetic flux from the FR. The mass and magnetic flux of the FR still change with time after the FR detaches from the chromosphere. It should be noted that the mass carried away from the FR by the magnetic islands is small compared to the mass of the FR while the electric current carried away from the FR by the magnetic islands is large, which can lead to the sudden decrease of the electric current of the FR, as is discussed in Section 3.3.

Prominence plasma is dynamic rather than static, and we note an interesting morphology change of the prominence during the rising process. As shown in Figures 2(a)–(d), the prominence expands first to fill the entire concave-upward part of the FR. At ~3400.6 s, the upper part of the prominence contracts. In Figure 2(g), we see that the width of the prominence varies with height. The width of the prominence changes from 5 Mm at 2885.4 s to 10 Mm at 3400.6 s, 8 Mm at 4328.1 s, and finally to 5 Mm at 4505.0 s. The length/width aspect ratio of the prominence changes from 0.18 at 2885.4 s to 1.1 at 3400.6 s, 2.375 at 4328.1 s, and eventually to 3.8 at 4505.0 s. Here we define the width and length of the prominence as the maximum horizontal and vertical distance ranges within the prominence region. The prominence presents a sheet-like structure at last. This shape change is driven by gravity and Lorentz force, which is discussed in Section 3.5.

Figure 3 shows the temporal evolution of various physical properties of the FR and the prominence. The vertical dashed lines in every panel indicate the instant 4122 s, when the first magnetic islands appear in the current sheet underneath, and the instant 4190.7 s, respectively.

Zoom In Zoom Out Reset image size

The top left panel shows the temporal evolution of the area of the cross-section of the FR and the prominence. The area of the FR increases from 0.95 Mm² at ~2000 s to $2250\,{\mathrm{Mm}}^{2}$ at ~4500 s, indicating the growth and the expansion of the FR. Compared to the great change of the area of the FR, the area of the prominence remains almost the same beginning at 3050 s, the instant when the material lifted up by the FR is detached from the chromosphere and the prominence forms. The area of the FR increases smoothly before 4122 s. After 4190.7 s, a drastic increase in the area of the FR begins, with some erratic fluctuations that are caused by the merging and separating of the magnetic islands and the FR.

The temporal evolution of the mass of the FR and the prominence are plotted in the top right panel of Figure 3. The mass is calculated as $m=\int \rho (x,y,z){dxdy}$, where the integral is over the region of the FR and the prominence. The region of the prominence is selected as the region with a temperature less than 25,000 K, a density higher than $2.3417\times {10}^{-14}\,{\rm{g}}\,{\mathrm{cm}}^{-3}$, a height above 6 Mm, and an electric current density less than $1.271\times {10}^{-3}\,{\rm{A}}\,{{\rm{m}}}^{-2}$. The last condition excludes the current sheet. The region of the FR is selected as the region where the value of the magnetic flux function A, which determines the magnetic field in the x–y plane by $\triangledown \times (A{{\boldsymbol{e}}}_{z})$, is larger than that of the uppermost X-point in the current sheet. The mass calculated here has the unit g cm⁻¹, representing the mass per unit length along the z-direction along which the variation of the physical variables is neglected in our 2.5D simulation. The mass of the prominence takes about 97% of the mass of the FR, suggesting that the mass of the FR is dominated by the mass of the prominence. The FR mass curve rises almost linearly before 3050 s but gradually flattens after 3050 s, the instant when the prominence detaches from the chromosphere and the mass exchange between the FR and the chromosphere is reduced. We determine the mass of the prominence per unit length as being about 1.6 × 10⁵ g. If the typical length of the prominence in the z-direction is taken as 2.0 × 10⁹ cm, then the mass of the prominence formed in our simulation will be estimated as 3.2 × 10¹⁴ g, in the same range as is given by Gilbert et al. (2005).

The middle left panel of Figure 3 presents the magnetic flux through the cross-sections in the x–y plane while the middle right panel of Figure 3 shows the electric current integrated over the cross-sections of the FR and the prominence in the x–y plane. There are obvious changes of the curves at 4122 s and 4190.7 s related to the current sheet dynamics. The magnetic flux of the FR (blue line in the middle left panel) increases linearly and smoothly from 2700 to 4122 s. After 4190.7 s, the FR magnetic flux curve rises with a bigger slope, accompanied by the fluctuations generated by the merging and separating of the magnetic islands with the FR. The magnetic flux through the prominence decreases after 4190.7 s. The decrease is due to the change of the relative position between the prominence and the FR. As pointed out previously, the bottom of the prominence leaves the bottom of the FR as the FR expands rapidly after 4190.7 s. Within the FR, the prominence moves inward and is detached from the border of the FR. The relatively weak z-component of the magnetic field in the inner region of the FR explains the decrease.

The curves of the electric current in the FR and the prominence are also changed at 4122 and 4190.7 s. The electric current of the FR drops to a lower value at 4122 s while the curve of the current in the prominence transits more smoothly. As shown in Figure 2(i), the lower part of the FR is elongated and squeezed into the current sheet, forming a tail-like structure that carries a much higher electric current density than other parts of the FR. The tail-like structure is fragmented into pieces and is disconnected with the FR at 4122 s, leading to an immediate loss of the current of the FR. From 4190.7 s, the curve of the electric current of the FR goes up with spiky structures since the magnetic islands generated by the unsteady bursty reconnection in the current sheet are ejected and ultimately merge into the FR, carrying electric current. The electric current in the prominence goes down for the same reason as the decrease of the magnetic flux. The prominence moves inward to the inner region of the FR where the current density is lower.

The bottom panels of Figure 3 show the kinematic evolution of the FR and the prominence. The bottom left panel shows the height of the mass center of the FR, the mass center of the prominence, and the FR center (O-point) plotted versus time. The plot shows that the curve of the FR mass center and the curve of the prominence mass center are almost overlaid, indicating once more that the mass of the prominence takes a large proportion of the mass of the FR. We note that the FR grows in size but does not rise before 2700 s. As the FR grows in size, the height of the FR center (O-point) and of the mass center increase. After 2700 s, the FR rises and eventually detaches itself from the chromosphere. The three curves in the bottom left panel of Figure 3 rise linearly from 2700 to 3700 s. The slopes of the three curves in this plot increase after about 3800 s, and after 4190.7 s the slopes increase more rapidly. The FR center (O-point) profile is similar to the profile of the mass center. Figure 3, bottom right panel, shows the velocity of the mass centers of the FR and the prominence, and their similarity verifies that the FR kinematics are dominated by the prominence. Before 2700 s, the velocity increase indicates the increasingly rapid expansion and size growth of the FR. The velocity of the FR decreases from 15.8 km s⁻¹ at 2700 s to 13.8 km s⁻¹ at 3700 s. The FR undergoes a slow rise with an almost constant but decreasing velocity of about 13–16 km s⁻¹ with fluctuations from 2700 to 3700 s and starts to accelerate gradually after 3700 s. An impulsive acceleration occurs after 4190.7 s with an acceleration value of 100 m s⁻². As noted by Zhang et al. (2001, 2004), the kinematic evolution of an impulsive CME contains a first phase (initiation phase) of slow rise (<80 km s⁻¹ ), a second phase (acceleration phase) with an acceleration of $100\mbox{--}500\,{\rm{m}}\,{{\rm{s}}}^{-2}$, and a third phase (propagation phase) with a constant or declining speed. Our simulation covers the first phase and the transition to the second phase. The transition from the first phase to the second phase starts at 4122 s, exactly the same instant when resistive instabilities occur in the current sheet. Note that many earlier CME models did not account for the embedded filament matter, and our model shows that their joint evolution is a unified ingredient in the overall dynamics.

Magnetic reconnection is an essential part of understanding solar flares. The Petschek model predicts fast reconnection rates comparable to observations, but it is not a self-consistent model within the limit of a small η (Biskamp 1997). A locally enhanced anomalous resistivity is often used to eliminate this inconsistency in the diffusion region in numerical simulations (Chen & Shibata 2000). However, the physical origin of the anomalous resistivity is still under debate. As described in Section 2.1, the anomalous resistivity is not adopted in our model.

As illustrated in Figures 1 and 2, the current sheet does not appear until ~2988.4 s, when the FR center moves with a velocity of ~10 km s⁻¹. It should be noted that we start our simulation with a 4 AMR level, and the grids are gradually refined to higher AMR levels as time goes on. Once the current sheet is formed, the grids are refined to the highest resolution. The AMR level of the grids covering the chromosphere and the transition region is at least 4, making the transitions of the density, the temperature, the plasma beta, and the resistivity smooth between the chromosphere and the corona. Note that if there is some nonsmoothness of density, the grids will be quickly refined to higher AMR levels. The current sheet evolution is a dynamical process in our simulation. The current sheet extends in length when the FR moves upward as shown in Figures 1 and 2, developing in a Sweet–Parker style before ~4122 s when the resistive instabilities are invoked.

3.4.1. Principal X-point

At every stage of the current sheet evolution since it appears, X-points always exist while O-points (in the center of each plasmoid) appear at the same time with the appearance of plasmoids. There exists a PX that has the maximal reconnection rate among all the X-points. The so-called absolute reconnection rate is measured as the reconnecting electric field (E_rec) at the X-point for a two-dimensional system (Priest & Forbes 2000), i.e.,

$$\begin{eqnarray}&&{E}_{\mathrm{rec}}=\eta {j}_{z}+\displaystyle \frac{{v}_{y}{B}_{x}}{c}-\displaystyle \frac{{v}_{x}{B}_{y}}{c},\end{eqnarray} \tag{ 21 }$$

with j_z representing the electric current density perpendicular outward to the x–y plane. The electric current density j_z may also be directed inward to the x–y plane but is identified as the reversed current, which we discuss in Section 3.4.4. We use this to locate the PX as the point with the maximal E_rec in the current sheet. The PX is indicated as a large green "X" in Figure 4.

Zoom In Zoom Out Reset image size

3.4.2. Thickness and Aspect Ratio

As pointed out by Biskamp (1997), the profile $j={j}_{\max }/{\cosh }^{2}(x/\delta )$ is found to fit the electric current density fairly well in the diffusion layer, where j_max is the maximal electric current density across the current sheet horizontally and δ is the half thickness of the current sheet. This profile is adopted to fit the current sheet along the line crossing the PX. As shown in Figure 4, the PX is marked as a big green "X," and the fitting result at 3984.6 s is illustrated in Figure 5(a). As indicated in Figure 5(b), the thickness of the current sheet becomes thinner and thinner before 4190.7 s, and the fluctuation becomes obvious after 4190.7 s. Not until 4190.7 s does the current sheet thin down to two to three cells. The thinnest thickness of the current sheet observed in this simulation is 68.01 km, larger than our grid spacing of 24.414 km. The current sheet does not form until 2500 s. Before 2500 s, there is a single X-point with the relatively concentrated electric current. The fitted thickness of the current distribution is larger than 1000 km before 2500 s. As the magnetic flux is continuously carried to the magnetic reversal line by the converging flow, the current distribution becomes more and more condensed around the single X-point. The single X-point eventually collapses at about 2500 s, leading to the formation of the current sheet. The thickness of the current sheet has its real significance after 2500 s in Figure 5(b).

Zoom In Zoom Out Reset image size

As the reconnection proceeds, the aspect ratio of the current sheet increases. Not until the aspect ratio reaches about 76.9 at ~4122 s do the first magnetic islands appear. The linear theory of resistive instability predicts that the tearing mode is invoked when the current sheet aspect ratio exceeds 2π (Furth et al. 1963; Priest & Forbes 2000). Nonlinear theories and simulations show that the aspect ratio could even exceed 100 before the tearing mode occurs (Loureiro et al. 2007; Shen et al. 2011). Biskamp (1997) gave a relation between the Lundquist number of the current sheet (S₀) and the aspect ratio (A):

$$\begin{eqnarray}&&A={S}_{0}^{1/2}{\left(1+\displaystyle \frac{\mu }{\rho \eta }\right)}^{-1/4}\approx {S}_{0}^{1/2},\end{eqnarray} \tag{ 22 }$$

where μ is the dynamic viscosity and η is the resistivity. Substituting A = 76.9 into the equation, we estimate the Lundquist number of the current sheet as 5911 when the resistive instabilities are invoked. The Lundquist number of the current sheet can also be estimated by ${S}_{0}={{lv}}_{A}/\eta$, where all the quantities are nondimensionalized by the units listed in Table 1. Taking ${v}_{A}=1.3$, the local Alfvén speed at a position a half thickness away from the center of the current sheet, the length of the current sheet l = 1.04, and η = 0.0002, we have ${S}_{0}=6760$, which is comparable to 5911.

3.4.3. Reconnection Rate

The reconnection rate near the PX is measured in two forms in our simulation: the absolute reconnection rate and the relative reconnection rate. As mentioned previously, the absolute reconnection rate is measured by the reconnecting electric field E_rec. The relative reconnection rate is characterized by the external Alfvén Mach number M_Ae, which is defined as the ratio between the external inflow speed v_ie and the local Alfvén speed ${v}_{\mathrm{Ae}}={B}_{e}/\sqrt{4\pi \rho }$, where B_e is the magnetic field strength in the external inflow region. In the classical steady-state reconnection theory, the electric field is dominated by the convective electric field $({v}_{{ie}}{B}_{e})/c$ outside the diffusion layer while the resistive electric field $\eta {j}_{z}$ is dominated inside the diffusion layer, and they should be matched at the boundary of the diffusion layer (Priest & Forbes 2000), so we have

$$\begin{eqnarray}&&{E}_{\mathrm{rec}}=\eta {j}_{z}=\displaystyle \frac{{v}_{{ie}}{B}_{e}}{c}=\displaystyle \frac{{v}_{{ie}}{B}_{e}}{c}\displaystyle \frac{{v}_{\mathrm{Ae}}{B}_{e}}{{v}_{\mathrm{Ae}}{B}_{e}}=\displaystyle \frac{{v}_{\mathrm{Ae}}{B}_{e}}{c}{M}_{\mathrm{Ae}}.\end{eqnarray} \tag{ 23 }$$

Equation (23) gives the relation between the absolute reconnection rate E_rec and the relative reconnection rate M_Ae,

$$\begin{eqnarray}&&{M}_{\mathrm{Ae}}=\displaystyle \frac{{{cE}}_{\mathrm{rec}}}{{v}_{\mathrm{Ae}}{B}_{e}}.\end{eqnarray} \tag{ 24 }$$

To calculate the relative reconnection rate, i.e., the Alfvén Mach number M_Ae, a small rectangle 0.1 Mm in width by 0.05 Mm in height near the PX, as indicated by the white box in Figure 4, is selected to evaluate the average velocity and Alfvén velocity in the external inflow region. The selected rectangle is outside the diffusion region, on the left of the half thickness of the current sheet. We refer to the region within the rectangle selected as "the external inflow region." We compute the ratio of the average fluid velocity and the Alfvén speed within this region and plot it as the blue curve in the middle right panel of Figure 6. The relative reconnection rate M_Ae increases and reaches its maximum of 0.093 at 2800 s, then decreases to 0.076 at 3700 s and increases again before the resistive instabilities occur. After the resistive instabilities occur, the relative reconnection rate presents a spiky structure, indicating the unsteady and bursty feature.

Zoom In Zoom Out Reset image size

The reconnecting electric field E_rec at the PX is evaluated on average over a small square containing the PX at its center with its width equal to the thickness of the current sheet, which is indicated as a yellow box in Figure 7(b). We refer to the region within the box selected as "the internal region around the PX." We plot the reconnecting electric field ${E}_{\mathrm{rec}}=({v}_{y}{B}_{x}-{v}_{x}{B}_{y})/c+\eta {j}_{z}$, the convective term $({v}_{y}{B}_{x}-{v}_{x}{B}_{y})/c$, and the resistive term $\eta {j}_{z}$ versus time in the middle left panel of Figure 6, denoted by the red, green, and blue solid curves, respectively. The reconnecting electric field is dominated by the resistive field before 4190.7 s while the convective field dominates after 4190.7 s as shown in the middle left panel of Figure 6. The reconnecting electric field E_rec increases before 3200 s and decreases from 8.9 V m⁻¹ at 3300 s to 7.8 V m⁻¹ at 4122 s. From 4190.7 s on, the E_rec increases drastically and fluctuates violently. The bottom left panel of Figure 6 shows the variations of the current density j_z directed outward to the x–y plane and the resistivity η in the internal region around the PX with time as the PX moves its position. The current density increases before 4190.7 s while the resistivity increases between 2000 and 2200 s, keeps constant between 2200 and 2600 s, and decreases beginning at 2600 s. The variations of the resistivity and the current density clearly explain the changes of the resistive term $\eta {j}_{z}$ and the reconnecting electric field E_rec before 4190.7 s. After 4190.7 s, the resistivity becomes weak as the PX moves to a higher position and the convective field dominates. The upper right panel of Figure 6 shows the variations of inflow speed and the Alfvén speed in the external inflow region. The drastic increase of the convective field indicated in the middle left panel of Figure 6 suggests a stronger inflow after 4190.7 s, which is also indicated by the increase of the inflow speed after 4190.7 s in the upper right panel of Figure 6. The temporal behavior and the order of magnitude of the absolute reconnection rate from 4122 s on resemble the reconnecting electric field derived in an M1.0 flare on 2000 September 12 by Qiu et al. (2004). The impulsive increase of E_rec from 4190.7 s onward may correlate with the microwave light curve and the hard X-ray impulse during the flare process (Qiu et al. 2004).

Zoom In Zoom Out Reset image size

As indicated in the upper left panel of Figure 6 as well as in Figure 4, the PX is always under the chromosphere before 4190.7 s and jumps from the chromosphere to the corona at 4190.7 s, which means that the resistive instabilities have sufficiently developed so that the unsteady reconnection due to the resistive instabilities in the corona is faster than the steady reconnection in the chromosphere. There is a time lag of about 69 s between the instant 4122 s, at which the first magnetic islands appear on the current sheet, and 4190.7 s, at which the PX jumps from the chromosphere to the corona. Simple estimation shows that the Alfvén transit time of the current sheet ${t}_{A}=l/{v}_{A}=67{\rm{s}}$ is comparable to the 69 s time lag, where l = 10 Mm is the length of the current sheet and ${v}_{A}=150\,\mathrm{km}\,{{\rm{s}}}^{-1}$ is the local Alfvén speed at a position a half thickness away from the center of the current sheet. The local Alfvén speed v_A in the current sheet is slower than the ambient corona speed due to the weaker magnetic field and the higher density as pointed out in Section 3.2. When the fragmentation happens in the corona, the chromosphere will respond 67s later, based on this estimation. So we consider the beginning of the fragmentation and the jump of the PX as a single event, and we do not distinguish them anymore. The PX is close to the polarity inversion line in the chromosphere before 4190.7 s, and the reconnection there may relate to transient brightenings (Harrison et al. 1985; Harrison 1986; Chifor et al. 2007), which have been long associated with the initiation of CMEs. Chifor et al. (2007) suggested that the slow reconnection prior to a full eruption is associated with X-ray transient brightenings, which are observable precursors in the early stages of eruptions. The localized X-ray brightenings are observed between two and 50 minutes by Chifor et al. (2007) before the impulsive phase of the flare and the FR acceleration, which is consistent with the time duration of the second stage of our system evolution from 1992.3 to ~4122 s.

To compare the absolute reconnection rate with the relative reconnection rate, we normalize E_rec dividing by $({v}_{\mathrm{Ae}}{B}_{e})/c$. The external Alfvén velocity v_Ae and the external magnetic field B_e are evaluated as the average values in the external inflow region. The time profile of the normalized E_rec is plotted as the red curve in the middle right panel of Figure 6, which matches the time profile of the Alfvén Mach number perfectly well.

Shown in the upper right panel of Figure 6, the Alfvén speed increases much more than the inflow speed in the external inflow region after 4190.7 s, which explains the decrease in the Alfvén Mach number after 4190.7 s. The bottom right panel of Figure 6 shows the variation of the magnetic field strength evaluated in the external inflow region, which decreases after 4190.7 s while the Alfvén speed (shown in the upper right panel of Figure 6) in the external inflow region increases, indicating that the decrease of the density, which is proportional to the square of the ratio of the magnetic field strength and the Alfvén speed, is the determining factor for the increase in the Alfvén speed. We conclude that the decrease in the Alfvén Mach number after 4190.7 s is mainly due to the decrease of mass density. The Alfvén Mach numbers in the chromosphere and in the corona are not comparable owing to the density difference between the two atmospheric layers.

After resistive instabilities are invoked, the long current sheet is fragmented into shorter pieces. The shortest growth time of resistive instabilities is given by ${({\tau }_{d}{\tau }_{A})}^{0.5}$ (Priest & Forbes 2000; Aschwanden 2005), where ${\tau }_{d}={l}^{2}/\eta$ is the diffusion timescale and ${\tau }_{A}=l/{v}_{A}$ is the Alfvén transit time. The shortest length of the current sheet that is resolvable in the simulation should be at least 2 times the grid spacing, about 24 km in our simulation. In the units listed in Table 1, taking l = 0.0048, η = 0.0002, and v_A = 1.3, the local Alfvén speed at a position a half thickness away from the center of the current sheet, we estimate the shortest growth time resolvable in our simulation as 0.0206, corresponding to a physical value of 1.77 s, comparable to the time duration of the spiked time structure of X-ray pulses, about 0.2–2 s (Cheng et al. 2012; Aschwanden 2005), which means that our simulation can resolve a timescale comparable to the spiked time structure of the X-ray pulses. Most of the spiked time structures of the X-ray pulses occur during the rise phase and peak times of the flares (Cheng et al. 2012). Particle acceleration is closely related to the dynamics of magnetic islands. Particles can be temporarily trapped inside the magnetic islands, which prolongs the acceleration time (Aschwanden 2005). Particles can gain kinetic energy by reflecting from the ends of the contracting magnetic islands (Drake et al. 2006; Guidoni et al. 2016). The nonthermal particles are associated with the radio and X-ray emission during the flare process. Detailed calculations and more numerical experiments are required to study the particle acceleration and radiation in the future.

Since ~4122 s, plasmoids of different scales are produced. The first magnetic island moves upward at a speed of about 180 km s⁻¹. Along with the tearing mode, coalescence, i.e., the merging of the neighboring magnetic islands into a single one, is another essential process occurring in the dynamical evolution of the current sheet. The first coalescence occurs at 4150 s.

3.4.4. Reversed Currents

3.4.5. Formation of Visible Arcade Structure

Observations have shown that the bulk of the acceleration of a CME occurs near the Sun, below $2\mbox{--}4\,{R}_{\odot }$ (MacQueen & Fisher 1983; Vršnak 2001; Chen & Krall 2003). Flare-associated CMEs attain faster speeds than solitary prominence-associated CMEs (St. Cyr et al. 1999; Chen & Krall 2003). The acceleration of a flare-associated CME is characterized by its "impulsiveness" (MacQueen & Fisher 1983). Why is the CME accelerated impulsively? Why does the kinematic evolution of the FR change from the initiation phase to the acceleration phase as pointed out in Section 3.3? Where does the impulsiveness of the accelerations of CMEs originate from? To answer these questions, we examine the forces acting upon the FR.

3.5.1. Global Momentum Balance

The total force ${{\boldsymbol{F}}}_{\mathrm{Total}}$ acting upon the FR, which we treat as a moving fluid parcel, is equal to the time rate of change of its momentum, which can be described by Lagrangian derivatives, i.e.,

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Total}}=\displaystyle \frac{d}{{dt}}{\int }_{\mathrm{FR}}\rho {\boldsymbol{v}}{dV}.\end{eqnarray} \tag{ 25 }$$

Reynolds's Transport Theorem (Wesseling 2001) gives

$$\begin{eqnarray}&&\displaystyle \frac{d}{{dt}}{\int }_{\mathrm{FR}}\rho {\boldsymbol{v}}{dV}={\int }_{\mathrm{FR}}\left[\displaystyle \frac{\partial }{\partial t}(\rho {\boldsymbol{v}})+\triangledown \cdot (\rho {\boldsymbol{vv}})\right]{dV}.\end{eqnarray} \tag{ 26 }$$

The two terms on the right-hand side of the Equation (26) have clear physical significances. The first term represents the local rate of change of the momentum. The second term can be rewritten as ${\int }_{\mathrm{FR}}\triangledown \cdot (\rho {\boldsymbol{vv}}){dV}={\oint }_{\mathrm{FR}}\rho {\boldsymbol{v}}{v}_{n}{dS}$, where v_n is the velocity component perpendicular outward to the surface of the FR, dS is the surface element of the FR, and the integral is taken over the surface of the FR. The second term represents the momentum transported outward of the FR per unit time—that is, the force exerted on the surrounding fluid by the FR, which represents the interaction between the FR and the surroundings through the surface.

Inserting Equation (1), the mass continuity equation, Equation (25) can be rewritten as

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Total}}=\displaystyle \frac{d}{{dt}}{\int }_{\mathrm{FR}}\rho {\boldsymbol{v}}{dV}={\int }_{\mathrm{FR}}\rho \displaystyle \frac{d{\boldsymbol{v}}}{{dt}}{dV}.\end{eqnarray} \tag{ 27 }$$

The MHD momentum equation reads

$$\begin{eqnarray}&&\rho \displaystyle \frac{d{\boldsymbol{v}}}{{dt}}=-\triangledown p+\triangledown \cdot {\boldsymbol{\tau }}+\displaystyle \frac{1}{c}{\boldsymbol{J}}\times {\boldsymbol{B}}+\rho {\boldsymbol{g}}.\end{eqnarray} \tag{ 28 }$$

Substituting the MHD momentum equation into Equation (27), we obtain

$$\begin{eqnarray}{{\boldsymbol{F}}}_{\mathrm{Total}} & = & \displaystyle \int \rho \displaystyle \frac{d{\boldsymbol{v}}}{{dt}}{dV}=\displaystyle \int (-\triangledown p+\triangledown \cdot {\boldsymbol{\tau }}\\ & & +\left.\displaystyle \frac{1}{c}{\boldsymbol{J}}\times {\boldsymbol{B}}+\rho {\boldsymbol{g}}\right){dV}.\end{eqnarray} \tag{ 29 }$$

Integrating Equation (29) over the region of interest (i.e., the FR region), we obtain the forces acting upon the FR as follows:

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Pressure}}=-\int \triangledown {pdxdy},\end{eqnarray} \tag{ 30 }$$

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Lorentz}}=\int \displaystyle \frac{1}{c}{\boldsymbol{J}}\times {\boldsymbol{B}}{dxdy},\end{eqnarray} \tag{ 31 }$$

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Viscous}}=\int \triangledown \cdot {\boldsymbol{\tau }}{dxdy},\end{eqnarray} \tag{ 32 }$$

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Gravity}}=\int \rho {\boldsymbol{g}}{dxdy},\end{eqnarray} \tag{ 33 }$$

$$\begin{eqnarray}&&{{\boldsymbol{F}}}_{\mathrm{Total}}={{\boldsymbol{F}}}_{\mathrm{Pressure}}+{{\boldsymbol{F}}}_{\mathrm{Lorentz}}+{{\boldsymbol{F}}}_{\mathrm{Viscous}}+{{\boldsymbol{F}}}_{\mathrm{Gravity}}.\end{eqnarray} \tag{ 34 }$$

The integral is only taken over the FR region on the x–y plane. The forces calculated represent the forces exerted on the body of interest per unit length with the unit (N m⁻¹) in the mks system of units. The force balance of the FR and the prominence along the y-axis is illustrated in the left and right panels of Figure 8, respectively. The left panel of Figure 8 shows that the total force acting on the FR is directed downward before 3700 s, which means that the acceleration of the FR is directed downward. Then the question arises as to why the FR goes up. As pointed out in Section 3.3, the rise of the center (O-point) of the FR and its mass center is due to the FR size growing and expanding before 2700 s. The converging flow continuously pushes the magnetic field to the magnetic reversal line, and reconnection leads to the growth in size of the FR. Although the total force is directed downward before 2700 s, the center (O-point) of the FR and its mass center rise as the FR expands in size. A direct consequence of the downward-directed total force is that the mass within the FR is concentrated in the lower part of the FR. After 2700 s, the total force acting on the FR oscillates with a small amplitude and is directed downward before 3700 s as is shown in the left panel of Figure 8. The oscillating amplitude of the total force acting on the FR is small. The bottom right panel of Figure 3 shows that the velocity curves are almost flattened from 2700 to 3700 s, indicating a near-zero acceleration, which is consistent with the small total force during this period. From 2700 to 3700 s the FR as a whole, which always stays close to the force-equilibrium point, can be considered to undergo a series of quasi-equilibrium states as the magnetic field is driven to evolve by the converging motion. We thus conclude that the rise of the FR from 2700 to 3700 s is due to the rise of the force-equilibrium point. At 2700 s, the FR reaches a force-equilibrium position as it continuously grows in size. From 2700 to 3700 s, the equilibrium point moves upward and our FR moves with the equilibrium point. A helical FR may suddenly jump from one equilibrium point to another as a result of catastrophe (Forbes & Isenberg 1991; Lin & Forbes 2000). The FR here jumps to a higher position of equilibrium. Then the process repeats—the FR loses equilibrium and jumps to a higher position as the equilibrium point moves up.

Zoom In Zoom Out Reset image size

During the process, the current sheet grows in length and reconnection cuts the magnetic field lines passing over the FR and tied to the photosphere. The Lorentz force eventually dominates the gravity after 3700 s, accelerating the FR. As discussed previously, the phase transition starts after 4122 s. After 4122 s, the Lorentz force and the pressure force with spiky structures increase drastically as shown in Figure 8. The spiky structure of the total force after 4190.7 s comes from the pressure force. The viscous force is negligible. Unlike with the FR, the pressure force acting upon the prominence is negligible throughout the time. The Lorentz force increases with time and lifts the prominence up into the corona.

We need to emphasize that the significant acceleration of the FR begins almost simultaneously with the beginning of the unsteady bursty regime of reconnection at 4190.7 s. The change of the reconnection mechanism from the Sweet–Parker to the unsteady bursty regime is the underlying physical reason for the phase transition of the FR kinematic evolution. Our simulation complements the major limitation on the understanding of the origin of CMEs: that the initiation and the early acceleration of CMEs are largely not observed, which leads to confusion on the relationship between CMEs and surface phenomena such as flares and filament eruptions (Zhang et al. 2004). Karpen et al. (2012) noticed that the phase transition of the numerical breakout model is attributed to the change in the form of the reconnection from the preliminary tearing-like form characterized by the formation of small, nearly stationary plasmoids (magnetic islands) to a strongly dynamic form with Alfvénic jets and multiple islands ejected from the current sheet. We stress that in the initiation phase of our model there are no magnetic islands, and the phase transition occurs after the first magnetic islands form, different from what is proposed by Karpen et al. (2012).

3.5.2. Local Force Distribution

As mentioned above and indicated in the upper panel of Figure 8, the spiky structure of the total force after 4122 s originates from the pressure force. The question arises as to where the spiky structure of the pressure force originates from. To answer this question, local force distribution at different instants should be closely inspected. Figure 9 illustrates the distribution of the y-component of the pressure force in the x–y plane. As can be seen from Figure 9, the colors in the figure are relatively faint before 4190.7 s, indicating that the pressure force is comparatively weak before this instant, which is consistent with Figure 8 and what is mentioned in Section 3.5.1. The pressure force increases drastically after 4190.7 s, and the force distribution presents a circle-layer pattern. As illustrated in Figure 9(f), from the outermost layer of the FR to the center of the FR the pressure force changes direction several times. The forces directed upward and downward are distributed in the neighboring circle-layers in a staggered form. If the force is directed upward in the upper part of a circle-layer, it will be directed downward in the corresponding lower part and vice versa. The force in the central region of the FR is weaker than that in the border of the FR. The number of circle-layers increases as time goes on, especially after 4122 s. From 4242.2 to 4493.0 s, the number of circle-layers increases a good deal, which indicates the drastic change within the FR after 4190.7 s. The pressure force exerted on the FR is dynamic rather than static, which drives the local evolution of the FR. The upward and downward forces may cancel each other at one instant but not another, so the total pressure force integrated over the FR region may vanish at some instant but not another, giving rise to the spiky structure of the pressure force illustrated in the upper panel of Figure 8. We note that the pressure force in the prominence region is weak compared to other parts of the FR, which explains the almost vanished pressure force exerted on the prominence as illustrated in the lower panel of Figure 8. Indicated in the bottom panel of Figure 8, the total force acting on the prominence is dominated by the Lorentz force and gravity. As we mentioned in Section 3.2, the shape change of the prominence is driven by gravity and the Lorentz force.

Zoom In Zoom Out Reset image size

The investigation of the evolution of the magnetic energy, the kinetic energy, and the internal energy of the system provides a deeper insight into the energy accumulation and conversion of the system. In this section, we discuss the magnetic energy accumulation, dissipation, and conversion, and we inspect the conservation of magnetic energy on the region [−5, 5] × [0, 10] (length unit ${L}_{0}=10\,\mathrm{Mm}$). The conservation law of the magnetic energy reads

$$\begin{eqnarray}&&\int \displaystyle \frac{\partial }{\partial t}\left(\displaystyle \frac{{{\boldsymbol{B}}}^{2}}{8\pi }\right){dVdt}=\int {P}_{m}{dVdt},\end{eqnarray} \tag{ 35 }$$

$$\begin{eqnarray}{P}_{m} & = & -\,{\boldsymbol{E}}\cdot {\boldsymbol{J}}-\triangledown \cdot {\boldsymbol{S}}-(\triangledown \cdot {\boldsymbol{B}}){\boldsymbol{B}}\cdot {\boldsymbol{v}}\\ & & +{\boldsymbol{B}}\cdot \triangledown [{C}_{d}{({\rm{\Delta }}l)}^{2}\triangledown \cdot {\boldsymbol{B}}],\end{eqnarray} \tag{ 36 }$$

where P_m is the magnetic power density, ${\boldsymbol{S}}=(c/4\pi ){\boldsymbol{E}}\times {\boldsymbol{B}}$ is the Poynting vector, and ${\boldsymbol{E}}\cdot {\boldsymbol{J}}=\eta {{\boldsymbol{J}}}^{2}+(1/c)({\boldsymbol{J}}\times {\boldsymbol{B}})\cdot {\boldsymbol{v}}$ is the electrical power density. Now we inspect the magnetic energy conservation over the region [−5, 5] × [0, 10]. The left panel of Figure 10 indicates the temporal evolution of the relative magnetic energy $\left(\int ({\boldsymbol{B}}\cdot {\boldsymbol{B}})/(8\pi ){dxdy}-{E}_{m0}\right)$, represented by the solid black line; the relative kinetic energy $(\int (1/2)(\rho {\boldsymbol{v}}\cdot {\boldsymbol{v}}){dxdy}-{E}_{k0})$, represented by the dashed purple line; and the relative internal energy $\left(\int (3p/2){dxdy}-{E}_{i0}\right)$, represented by the dashed orange line, where E_m0, E_k0, and E_i0 are the magnetic energy, the kinetic energy, and the internal energy at the instant 0 s and all the integrals are taken over the region [−5, 5] × [0, 10] with the length unit 10 Mm. The net magnetic energy gain $\int {P}_{m}{dxdydt}$ (red solid line), the energy transported by Poynting flux $-\int \triangledown \cdot {\boldsymbol{S}}{dxdydt}$ (green solid line), the dissipation by Joule heating $\int \eta {{\boldsymbol{J}}}^{2}{dxdydt}$ (yellow solid line), and the magnetic energy loss due to Lorentz work $-\int (1/c)({\boldsymbol{J}}\times {\boldsymbol{B}})\cdot {\boldsymbol{v}}{dxdydt}$ (blue solid line) are also plotted in the left panel of Figure 10. Theoretically, the net magnetic energy gain $\int {P}_{m}{dxdydt}$ and the relative magnetic energy $\int ({\boldsymbol{B}}\cdot {\boldsymbol{B}})/(8\pi ){dxdy}-{E}_{m0}$ should be equal. Obviously, in the left panel of Figure 10, the solid black line representing the relative magnetic energy and the solid red line representing the net magnetic energy gain are not overlaid on each other. We calculate the relative difference between the two, i.e.,

$$\begin{eqnarray}&&{d}_{r}=\left(\int {P}_{m}{dxdydt}+{E}_{m0}-\int \displaystyle \frac{{{\boldsymbol{B}}}^{2}}{8\pi }{dxdy}\right)/\int \displaystyle \frac{{{\boldsymbol{B}}}^{2}}{8\pi }{dxdy}.\end{eqnarray} \tag{ 37 }$$

The relative difference d_r is plotted over time in the bottom right panel of Figure 10. Initially the relative difference d_r is of the order of 0.001, then it increases from 0.001 to 0.01 as time goes on, and finally it keeps the order of 0.01, which means that the numerical code preserves the conservation of magnetic energy at an accuracy of 0.01. The increase of the relative difference d_r from about 1000 to 3000 s may be attributed to the gradual refining of the computational grids during this period, and more experiments and investigations should be conducted in the future.

Zoom In Zoom Out Reset image size

As can be seen from the left panel of Figure 10, once the converging driving starts, the relative magnetic energy decreases and reaches its minimum at 850 s. The relative magnetic energy increases almost linearly from 850 to 3100 s, indicating the energy accumulation driven by the converging motion. The first O-point, i.e., the FR center, appears at 1992.3 s, indicating the ignition of reconnection. Even if the reconnection starts, the relative magnetic energy continues to increase. The decrease and increase of the relative magnetic energy from 0 to 3100 s are dominated by Poynting flux (green solid line), which means that the energy is transported outside the region by Poynting flux from 0 to 850 s and then the energy is continuously transported into the region by Poynting flux from 850 to 3100 s, indicating the magnetic energy accumulation process.

The relative magnetic energy declines after 3100 s because the magnetic energy is converted to Joule heating and Lorentz work, indicating the magnetic energy dissipation and conversion process. It should be noted that the dissipation of Joule heating is overestimated in this simulation due to the artificiality of the resistivity, which is described in Section 2.1. The kinetic energy increases from 0 to ~4122 s while the internal energy increases from 0 to ~2500 s and decreases from ~2500 to ~4122 s. Besides Joule heating, the internal energy is influenced by radiative loss and thermal conduction, so the increasing Joule heating does not directly lead to the increase of the internal energy. The decrease in the magnetic energy becomes even steeper after ~4122 s while the internal and the kinetic energy increase drastically after this instant, suggesting the phase transition at the instant.

We have included the two terms added to the energy equation after the divergence cleaning: the Linde term $\int {\boldsymbol{B}}\cdot \triangledown [{C}_{d}{({\rm{\Delta }}l)}^{2}\triangledown \cdot {\boldsymbol{B}}]{dxdy}$ and the Powell term $\int -(\triangledown \,\cdot {\boldsymbol{B}}){\boldsymbol{B}}\cdot {\boldsymbol{v}}{dxdy}$, when calculating the net magnetic energy gain $\int {P}_{m}{dxdydt}$ (red solid curve in the left panel of Figure 10). The order of magnitude of the Powell term is ${10}^{14}\,{\mathrm{erg}}^{-1}\,{\mathrm{cm}}^{-1}\,{{\rm{s}}}^{-1}$, two orders of magnitude higher than the Linde term and two orders of magnitude lower than the other terms in Equation (36), which means that the Linde and the Powell terms have minor influences on the energy calculation. The Linde term and the Powell term are not plotted in the left panel of Figure 10 since they are too small compared to other quantities.

Aly (1991) and Sturrock (1991) conjectured that among a series of force-free magnetic fields with a simple topology, i.e., all the field lines are unknotted and connected to the boundary on which the normal component of the magnetic field is fixed, the so-called open field with its field lines connecting the boundary to infinity, which is current-free everywhere except on a set of current sheets in equilibrium, maximizes the magnetic energy over the upper half space. We point out that the Aly–Sturrock conjecture, which is valid only when the magnetic field is simply connected and force-free, is not applicable to our case since the field lines of the initial simple magnetic topology become knotted and disconnected to the boundary due to the reconnection process in our resistive MHD simulation. Figure 9 illustrates clearly that the Lorentz force is nonvanished at least in some parts of the simulated region, indicating that the magnetic field deviates from the initial force-free state and becomes non-force-free as the system evolves. We should also emphasize that the normal component, i.e., the y-component, of the magnetic field at the bottom boundary, which is described in Section 2.3.2, varies with time so that the upper and lower bounds of the magnetic energy of the simply connected force-free fields predicted by the Aly–Sturrock conjecture change with time.

The final snapshot, i.e., the snapshot at 4505.0 s, of our 2.5D simulation is visualized in three dimensions. An iso-density surface at $\rho =468.34\times {10}^{-13}\,{\rm{g}}\,{\mathrm{cm}}^{-3}$, represented by the color yellow, is plotted in Figure 11. This iso-density surface gives a clear sketch of the prominence and chromosphere. Several bundles of the magnetic field lines (solid black lines) are plotted that clearly show the twisted magnetic FR, the sheared flare arcade, and a plasmoid. In 2.5D simulation, there are no differences between an FR and a magnetic island. In 2.5D simulation, the flare loops are an infinite arcade along the polarity inversion line, and only a small section is plotted here.

Zoom In Zoom Out Reset image size