NUMERICAL STUDY ON IN SITU PROMINENCE FORMATION BY RADIATIVE CONDENSATION IN THE SOLAR CORONA

The corona model is set up in a rectangular box, in which the x- and y-axis are horizontal and vertical, respectively, and the z-axis is orthogonal to the x–y plane. The corona is stratified under uniform temperature (${T}_{\mathrm{cor}}=1\;\mathrm{MK}$) and gravity (${g}_{\mathrm{cor}}=2.7\times {10}^{4}\;\mathrm{cm}\;{\rm{s}}^{-2}$) conditions,

$$\begin{eqnarray}&&\rho (y)={\rho }_{\mathrm{cor}}\mathrm{exp}\left[-{\displaystyle \frac{{{mg}}_{\mathrm{cor}}}{{k}_{\rm{B}}{T}_{\mathrm{cor}}}}y\right],\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&p(y)={\displaystyle \frac{{k}_{\rm{B}}}{m}}\rho (y){T}_{\mathrm{cor}},\end{eqnarray} \tag{ 2 }$$

where ${\rho }_{\mathrm{cor}}=3.2\times {10}^{-15}\;\rm{g}\;{\mathrm{cm}}^{-3}$, k_B is Boltzmann's constant, and m is the mean molecular mass. We set $m={m}_{p}$ with m_p being the proton mass. The force-free arcade field is described as

$$\begin{eqnarray}&&{B}_{x}=-\left({\displaystyle \frac{2{L}_{a}}{\pi a}}\right){B}_{a}\mathrm{cos}\left({\displaystyle \frac{\pi }{2{L}_{a}}}x\right)\mathrm{exp}\left[-{\displaystyle \frac{y}{a}}\right],\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&{B}_{y}={B}_{a}\mathrm{sin}\left({\displaystyle \frac{\pi }{2{L}_{a}}}x\right)\mathrm{exp}\left[-{\displaystyle \frac{y}{a}}\right],\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{B}_{z}=-\sqrt{1-{\left({\displaystyle \frac{2{L}_{a}}{\pi a}}\right)}^{2}}{B}_{a}\mathrm{cos}\left({\displaystyle \frac{\pi }{2{L}_{a}}}x\right)\mathrm{exp}\left[-{\displaystyle \frac{y}{a}}\right],\end{eqnarray} \tag{ 5 }$$

where ${B}_{a}=3\;\rm{G}$ is the field strength at the footpoint, ${L}_{a}=12\;\mathrm{Mm}$ is the width, and $a=31\;\mathrm{Mm}$ is the magnetic scale height of the arcade field. Initially, the system exists in mechanical equilibrium.

The left and right boundaries are subjected to symmetric (for $\rho ,p,{v}_{y},{B}_{y}$) or anti-symmetric (for ${v}_{x},{v}_{z},{B}_{x},{B}_{z}$) boundary conditions. A free boundary condition is applied to the top. In the region below y = 0, the converging and shearing motions are introduced. We test three types of the footpoint motions. One is the converging motion without shearing; the others are converging motion with shearing that increases the magnetic shear of the arcade field or with anti-shearing that decreases the shear of the arcade field. In all the cases, the velocity components v_x and v_y within this region are set as follows,

$$\begin{eqnarray}&&{v}_{x}=-{v}_{0}(t){\displaystyle \frac{x}{{L}_{a}/4}},\;{v}_{y}=0,\;\;(0\leqslant x\lt {L}_{a}/4,\;y\lt 0)\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}\begin{array}{ccc}{v}_{x} & = & -{v}_{0}(t){\displaystyle \frac{({L}_{a}-x)}{3{L}_{a}/4}},\\ & & {v}_{y}=0,\;\;({L}_{a}/4\leqslant x\leqslant {L}_{a},\;y\lt 0),\end{array}\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&{v}_{0}(t)={v}_{00},(0\lt t\lt {t}_{1})\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&{v}_{0}(t)={v}_{00}({t}_{2}-t)/({t}_{2}-{t}_{1}),({t}_{1}\leqslant t\lt {t}_{2})\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&{v}_{0}(t)=0,(t\geqslant {t}_{2})\end{eqnarray} \tag{ 10 }$$

where ${v}_{00}=12\;\mathrm{km}\;{\rm{s}}^{-1}$, ${t}_{1}=1000\;\rm{s}$, and ${t}_{2}=1500\;\rm{s}$. In the case of no shearing motion, ${v}_{z}=0$. The shearing and the anti-shearing motions are set as follows,

$$\begin{eqnarray}&&{v}_{z}=\pm {v}_{0}(t){\displaystyle \frac{x}{{L}_{a}/4}},\;\;(0\leqslant x\lt {L}_{a}/4,\;y\lt 0)\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&{v}_{z}=\pm {v}_{0}(t){\displaystyle \frac{({L}_{a}-x)}{3{L}_{a}/4}},\;\;({L}_{a}/4\leqslant x\leqslant {L}_{a},\;y\lt 0)\end{eqnarray} \tag{ 12 }$$

where the plus and minus signs represent the shearing and anti-shearing motions, respectively. The magnetic fields are computed with the induction equation by coupling with the given converging and shearing motions. The magnetic fields at the bottom boundary are fixed to the initial state. The gas pressure and density are fixed by assuming hydrostatic equilibrium at a constant temperature of ${T}_{\mathrm{cor}}={10}^{6}\;\rm{K}$.

The basic equations are

$$\begin{eqnarray}&&{\displaystyle \frac{\partial \rho }{\partial t}}+\boldsymbol{v}\cdot \bf{\nabla }\rho =-\rho \bf{\nabla }\cdot \boldsymbol{v},\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}\begin{array}{ccc}{\displaystyle \frac{\partial {e}_{\mathrm{th}}}{\partial t}}+\boldsymbol{v}\cdot \bf{\nabla }{e}_{\mathrm{th}} & = & -({e}_{\mathrm{th}}+p)\bf{\nabla }\\ & & \cdot \boldsymbol{v}+\bf{\nabla }\cdot \left(\kappa {T}^{5/2}\boldsymbol{b}\boldsymbol{b}\cdot \bf{\nabla }T\right)-L,\end{array}\end{eqnarray} \tag{ 14 }$$

$$\begin{eqnarray}&&{e}_{\mathrm{th}}={\displaystyle \frac{p}{\gamma -1}},\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}&&{\displaystyle \frac{\partial \boldsymbol{v}}{\partial t}}+\boldsymbol{v}\cdot \bf{\nabla }v=-{\displaystyle \frac{1}{\rho }}\bf{\nabla }p+{\displaystyle \frac{1}{4\pi \rho }}(\bf{\nabla }\times \boldsymbol{B})\times \boldsymbol{B}+\boldsymbol{g},\end{eqnarray} \tag{ 16 }$$

$$\begin{eqnarray}&&{\displaystyle \frac{\partial \boldsymbol{B}}{\partial t}}=-c\bf{\nabla }\times \boldsymbol{E},\end{eqnarray} \tag{ 17 }$$

$$\begin{eqnarray}&&\boldsymbol{E}=-{\displaystyle \frac{1}{c}}\boldsymbol{v}\times \boldsymbol{B}+{\displaystyle \frac{4\pi \eta }{{c}^{2}}}\boldsymbol{J},\end{eqnarray} \tag{ 18 }$$

$$\begin{eqnarray}&&\boldsymbol{J}={\displaystyle \frac{c}{4\pi }}\bf{\nabla }\times \boldsymbol{B},\end{eqnarray} \tag{ 19 }$$

where $\kappa =2\times {10}^{-6}\;\mathrm{erg}\;{\mathrm{cm}}^{-1}\;{\rm{s}}^{-1}\;{\rm{K}}^{-7/2}$ is the coefficient of thermal conduction, $\boldsymbol{b}$ is a unit vector along the magnetic field, and $\boldsymbol{g}=(0,-{g}_{\mathrm{cor}},0)$ is the gravitational acceleration and η is the magnetic diffusion rate. The temperature is computed by the following equation of state,

$$\begin{eqnarray}&&T={\displaystyle \frac{m}{{k}_{\rm{B}}}}{\displaystyle \frac{p}{\rho }}.\end{eqnarray} \tag{ 20 }$$

For fast reconnection, we adopt the following form of the anomalous resistivity (e.g., Yokoyama & Shibata 1994),

$$\begin{eqnarray}&&\eta =0,\quad \left(J\lt {J}_{c}\right)\end{eqnarray} \tag{ 21 }$$

$$\begin{eqnarray}&&\eta ={\eta }_{0}{\left(J/{J}_{c}-1\right)}^{2},\left(J\geqslant {J}_{c}\right)\end{eqnarray} \tag{ 22 }$$

where ${\eta }_{0}=3.6\times {10}^{13}\;{\mathrm{cm}}^{2}\;{\rm{s}}^{-1}$ and ${J}_{c}=25\;{\mathrm{dyn}}^{1/2}\;\mathrm{cm}\;{\rm{s}}^{-1}$. We restrict η to ${\eta }_{\mathrm{max}}=18.0\times {10}^{13}\;{\mathrm{cm}}^{2}\;{\rm{s}}^{-1}$. L in Equation (14) represents the net cooling term, which is expressed as

$$\begin{eqnarray}&&L={n}^{2}\rm{\Lambda\ }(T)-H-\eta {J}^{2},\end{eqnarray} \tag{ 23 }$$

where $n=\rho /m$ is the number density, $\rm{\Lambda\ }(T)$ is the radiative loss function of optically thin plasma, H is the background heating rate, and $\eta {J}^{2}$ represents ohmic heating. We adopt a simplified radiative loss function in Figure 2(a) (Hildner 1974) and simulate two different coronal heating models. In one model, the heating rate H depends on the local magnetic energy density (magnetic pressure) at temperatures above the cutoff temperature ${T}_{c}$ and is balanced out by the cooling rate when it falls below ${T}_{c}$,

$$\begin{eqnarray}&&H={\alpha }_{m}{P}_{m},\left(T\geqslant {T}_{c}\right)\end{eqnarray} \tag{ 24 }$$

$$\begin{eqnarray}&&H={n}^{2}\rm{\Lambda\ }(T),\left(T\lt {T}_{c}\right).\end{eqnarray} \tag{ 25 }$$

Initial thermal equilibrium required that

$$\begin{eqnarray}&&{\alpha }_{m}={\displaystyle \frac{{n}_{\mathrm{cor}}^{2}\rm{\Lambda\ }({T}_{\mathrm{cor}})}{{B}_{a}^{2}/8\pi }}\mathrm{exp}\left[-2\left({\displaystyle \frac{{{mg}}_{\mathrm{cor}}}{{k}_{\rm{B}}{T}_{\mathrm{cor}}}}-{\displaystyle \frac{1}{a}}\right)y\right],\end{eqnarray} \tag{ 26 }$$

where ${n}_{\mathrm{cor}}={\rho }_{\mathrm{cor}}/m=2.0\times {10}^{9}\;{\mathrm{cm}}^{-3}$. Selecting $a={k}_{\rm{B}}{T}_{\mathrm{cor}}/({{mg}}_{\mathrm{cor}})$, we obtain ${\alpha }_{m}=6.2\times {10}^{-4}\;{\mathrm{cm}}^{3}\;{\rm{s}}^{-1}$ (i.e., a constant). Equation (25) prevents the temperature from decreasing below ${T}_{c}$. In the other model, the heating rate H depends on the local density and height at temperatures above ${T}_{c}$ and balances with the cooling rate below the cutoff temperature ${T}_{c}$.

$$\begin{eqnarray}&&H={\alpha }_{d}(y)n,\left(T\geqslant {T}_{c}\right)\end{eqnarray} \tag{ 27 }$$

$$\begin{eqnarray}&&H={n}^{2}\rm{\Lambda\ }(T),\left(T\lt {T}_{c}\right).\end{eqnarray} \tag{ 28 }$$

For the initial thermal equilibrium, we have ${\alpha }_{d}(y)={n}_{\mathrm{cor}}\rm{\Lambda\ }({T}_{\mathrm{cor}})\mathrm{exp}[-{{mg}}_{\mathrm{cor}}y/({k}_{\rm{B}}{T}_{\mathrm{cor}})]$ where ${n}_{\mathrm{cor}}\rm{\Lambda\ }({T}_{\mathrm{cor}})=1.1\times {10}^{-13}\;\mathrm{erg}\;{\rm{s}}^{-1}$. Note that ohmic heating is always much smaller than background heating in these settings. ${T}_{c}$ is one of the model parameters. For each heating and shearing model, we vary ${T}_{c}$ as 0.3, 0.4, and 0.5 MK. The investigated cases are summarized in Table 1.

Zoom In Zoom Out Reset image size

Table 1.  The Presence of Radiative Condensation in Each Case

Case	Heating	Shearing	Tc (MK)	Condensation
M1	...	−	0.3	Yes
M2	...	−	0.4	Yes
M3	$H={\alpha }_{m}{P}_{m}$	−	0.5	Yes
M4	...	0	0.5	No
M5	...	+	0.5	No
D1	...	−	0.3	Yes
D2	...	−	0.4	Yes
D3	$H={\alpha }_{d}\rho$	−	0.5	Yes
D4	...	0	0.5	Yes
D5	...	+	0.5	Yes

Note. The second column shows the heating model. The third column shows the shearing model: plus sign (+), minus sign (−), and 0 represent the shearing, anti-shearing, and no shearing cases, respectively. The fourth column shows the cutoff temperature. The fifth column shows the presence of condensation.

Download table as:  ASCIITypeset image

The numerical scheme is based on the CIP-MOCCT method (Kudoh et al. 1998). To prevent numerical oscillations in the advection term (left-hand side of Equations (13)–(16)), we apply the rational CIP scheme (Xiao et al. 1996), in which the rational function replaces the cubic function as an interpolation function of CIP. The thermal conduction, radiative cooling, and heating terms are included as source terms in non-advection phases of CIP method. Thermal conduction is explicitly solved by a second-order accurate slope-limiting method (Sharma & Hammett 2007).

The required grid size for preventing numerical thermal instability can be estimated by the formula presented in Koyama & Inutsuka (2004). When applying this formula to our simulation, we assume that the density of the condensations is 10 times the coronal density,

$$\begin{eqnarray}&&\rm{\Delta\ }\lt {\displaystyle \frac{{\lambda }_{F}}{3}}={\displaystyle \frac{1}{3}}\sqrt{{\displaystyle \frac{\kappa {T}_{c}^{7/2}}{{n}_{c}^{2}\rm{\Lambda\ }({T}_{c})}}},\end{eqnarray} \tag{ 29 }$$

where Δ is the grid size, ${\lambda }_{F}$ is the Field length, ${T}_{c}$ is the cutoff temperature (also assume to be the temperature of condensations), and ${n}_{c}=10{n}_{\mathrm{cor}}$ is the assumed density of condensations. Figure 2(b) shows the required grid size against temperature. Under the constraint of Equation (29), the grid size is determined to be $\rm{\Delta\ }\;=$ 120, 60, and 30 km in the cases of ${T}_{c}\;=\;$0.5, 0.4, and 0.3 MK, respectively. The horizontal grid size $\rm{\Delta\ }x$ is uniform. The vertical grid size $\rm{\Delta\ }y$ is uniform up to $y\lt 24\;\mathrm{Mm}$ and increases to $\rm{\Delta\ }y=1.8\times {10}^{3}\;\mathrm{km}$ at $y\gt 24\;\mathrm{Mm}$ to ensure a sufficiently distant upper boundary.

In M-cases where the heating rate is proportional to the magnetic pressure, the anti-shearing and shearing motion causes a larger discrepancy between the cooling and heating rates, establishing the cooling- and heating-dominant thermal imbalances, respectively. Figure 7 shows the cooling and heating rates along the y-axis at time = $500\;\rm{s}$ in case M3–M5. In case M3 (anti-shearing case), the cooling rate is the highest of the three cases, whereas the heating rate is the lowest. Case M5 (shearing case) has the opposite tendency, where the cooling rate is the lowest and the heating rate is the highest. We discuss how these tendencies are established from a perspective of the evolution of the magnetic field and its influence on the cooling and heating rates.

Zoom In Zoom Out Reset image size

In the anti-shearing case, the flux rope is pushed upward by the reconnection and confined by the downward magnetic tension of the overlying arcade field. Figure 8 shows the vertical tension force along the y-axis at time = $500\;\rm{s}$ in case M3 (anti-shearing), M4 (no shearing), and M5 (shearing). The downward magnetic tension force in the anti-shearing case is larger than those in the shearing and no shearing cases. The reason the anti-shearing case has the larger downward tension force is as follows: the shear angle is related to the height of the arcade field as

$$\begin{eqnarray}&&a={\displaystyle \frac{{L}_{a}}{\pi \mathrm{cos}\theta }},\end{eqnarray} \tag{ 30 }$$

where a is the arcade height, $\theta =\mathrm{arctan}({B}_{z}/{B}_{x})$ is the shear angle against the positive x-axis in the xz-plane, L_a is the arcade width, and we assume a linear force-free arcade field for simplicity. Equation (30) indicates that the reduction of magnetic shear by the anti-shearing motion makes the arcade field shorter. Due to the shrink of the arcade field, the flux rope in the anti-shearing case experiences a larger downward magnetic tension force.

Zoom In Zoom Out Reset image size

The flux rope in the anti-shearing case is pinched by the overlying arcade field and pushed by the upward magnetic tension of the reconnected loops, leading to a density enhancement by compression. Figure 9 shows the time evolution of area and average density inside the dashed–dotted line in Figure 6. The average density is computed by dividing the total mass inside the dashed–dotted line by the corresponding area. The total mass is almost the same between these three cases and is conserved after the flux rope is detached from the bottom boundary. In Figure 9(b), the anti-shearing case shows the enhancement of the average density. This result indicates that the density enhancement inside the flux rope in the anti-shearing case is due to the compression by the magnetic tension as well as the converging motion.

Zoom In Zoom Out Reset image size

The anti-shearing motion led to the higher cooling rate as well as the lower heating rate. The anti-shearing motion reduces the magnetic component parallel to the PIL, corresponding to B_z in our simulations, and decreases the magnetic pressure, leading to the lowest heating rate. Thus, the interior of the flux rope in case M3 suffers from unstoppable cooling most efficiently. In case M5 (shearing case), on the other hand, the heating rate increases as the magnetic pressure increases, and the enhancement of the cooling rate is too small to overwhelm the heating rate. Eventually, the shearing motion is likely to establish the heating-dominant thermal imbalance, resulting only in the formation of hot flux rope. These results suggest that when the coronal heating depends on the magnetic energy density, the continuous shearing motion cannot trigger the radiative condensation, whereas the temporary anti-shearing motion can indeed trigger it, leading to the prominence formation.

From the results of the D-cases in which the heating is proportional to the density, it is found that the cooling-dominant thermal imbalance can be achieved only by the difference between the heating and cooling rates in the dependence on the density. Figure 10 shows the cooling and the heating rates in the cases D3–D5. The cooling rate in case D3 (anti-shearing case) is highest because the density enhancement is largest in the anti-shearing case, as discussed in the previous subsection. Unlike case M3, however, the heating rate in case D3 is also highest. The combination of the higher cooling rate and the higher heating rate can be a disadvantage to enhancing the thermal imbalance. What actually sets the cooling-dominant thermal imbalance in D-cases is the different density dependence between the heating ($H\propto \rho$) and cooling terms ($R\propto {\rho }^{2}$). The increment of the heating rate against the density increase is always smaller than that of the cooling rate. Consequently, the cooling overwhelms the heating, and the radiative condensation is triggered in all D-cases. These results suggest that if the coronal heating is only proportional to the density, the flux ropes in the corona are likely to cool down. The exact modeling of the coronal heating is required to confirm this mechanism.

Zoom In Zoom Out Reset image size

In evaporation–condensation models (Antiochos et al. 1999; Karpen et al. 2006), radiative condensation is triggered by the transport of dense plasma from lower altitudes. Our model is consistent with this mechanism, but it adopts a different approach to suppress the thermal conduction effects. In the evaporation–condensation model, thermal conduction effects are suppressed by preparing long magnetic field lines, whereas in our model, they are suppressed by closed loop formations. In other words, our model altered the magnetic field topology. Both our study and that of Choe & Lee (1992) assume that photospheric motion triggers the prominence formation. However, the direction of motion differs between the studies. In Choe & Lee (1992), radiative condensation arises by shearing motion along the PIL and subsequent expansion of the arcade field. In our study, converging motion toward the PIL and subsequent formation of closed magnetic fields (flux rope) are essential precursors of radiative condensation. Moreover, whereas their models reproduce normal-polarity prominences, our models form inverse-polarity prominences. The formation of the inverse-polarity prominence is found both in the results of our model and that of Linker et al. (2001). The difference is that our model does not require chromospheric plasma injection to trigger the radiative condensation; hence, our model reproduce the in situ condensation. Xia et al. (2014a) performed three-dimensional simulation of the in situ condensation inside the flux rope. The flux rope system in Xia et al. (2014a) was created by the converging and shearing motion in the isothermal simulation of Xia et al. (2014b). Because the strategy used to establish cooling-dominant thermal imbalance in Xia et al. (2014a) was the parameterized heating, it was still unclear why the thermal imbalance was created inside the flux rope. We ensure that the cooling-dominant thermal imbalance can be created inside the flux rope by two different mechanisms: the anti-shearing motion or the different dependence on the local density between the heating and cooling rates. We find that the direction of the shear motion strongly affected the thermal processes when the heating depends on the magnetic energy density, and it controls the presence of in situ condensation. When the heating is proportional to the density, on the other hand, the cooling-dominant thermal imbalance is always established by the difference between the cooling and heating terms in the dependence on density, leading to radiative condensation.

The formation mechanism through the anti-shearing motion smoothly connects to some of the eruptive models. These include instability models such as kink and torus instability (Kliem et al. 2004; Török & Kliem 2005; Kliem & Török 2006; Fan & Gibson 2007). The kink instability is likely to be triggered because the anti-shearing motion increases the twist of the flux rope. The torus instability also fit because the reduction of shear can alleviate the critical height of torus instability. The other model that fit with our mechanism is the reversed shear model of Kusano et al. (2004). We stopped the anti-shearing motion before the shear reversal. If we continued to impose the anti-shearing motion, the flux rope formation, radiative condensation, and eruption would have successively occurred even without the converging motion. Both the shearing and anti-shearing motions have been considered as the possible process of eruptions (Kusano 2002; Amari et al. 2003). Our results suggest that continuous shearing motion results in the eruption of the hot flux rope corresponding to flare eruption, whereas the anti-shearing motion causes the radiative condensation of the prominence formation. The eruptive mechanism coupling with flux emergence is also the subject to investigate (Chen & Shibata 2000; Shiota et al. 2005; Kaneko & Yokoyama 2014). The results of our simulation with the good mechanical balance can be directly adopted as the initial condition for these eruptive studies. The model including the dense materials of the prominence may reveal the more detailed structure of coronal mass ejection such as the dense core and its evolution.

We compare the simulated process of in situ condensation with the observation by Berger et al. (2012) by synthetic emission through the extreme ultraviolet (EUV) filters of the Solar Dynamics Observatory Atmospheric Imaging Assembly (SDO/AIA). The emission of a certain wavelength channel i is expressed as

$$\begin{eqnarray}&&{D}_{i}=\int {n}_{e}^{2}{K}_{i}(T){dl},\end{eqnarray} \tag{ 31 }$$

where D_i, n_e, K_i(T), and l represent the pixel response to the photon flux, electron number density, temperature response function of the AIA filters, and the distance along the line of sight (LOS), respectively. The temperature response functions are obtained from aia_get_response.pro in the Solar Software Library. To compute the synthetic emission, we choose the z-axis as the LOS. In our 2.5 dimensional simulation, the physical variables are assumed to be constant along the z-axis; hence, the integration along the LOS is altered to the multiplication of an arbitrary constant l_z as follows,

$$\begin{eqnarray}&&{D}_{i}={n}_{e}^{2}{K}_{i}(T){l}_{z}.\end{eqnarray} \tag{ 32 }$$

Figure 11 is the time evolution of the synthetic emission of each EUV filter. Note that Figure 11 covers the time only before the temperature reaches the cutoff. Figure 11 clearly shows that the areas of strong emission progressively shift from the 193 ${\AA}$ channel, through 171 ${\AA}$ channel, to 131 ${\AA}$ channel both in time and altitude. The emission shift was claimed to be evidence of in situ condensation in Berger et al. (2012), and our simulation has validated their argument. The shift from the 211 ${\AA}$ channel to the 193 ${\AA}$ channel is not so clear because the initial coronal temperature is set to 1 MK. The shift in altitude in Figure 11 results from the descent of the flux rope and high density plasma. It can also be realized by some other mechanisms such as magnetic reconnection in the weighted dip (Liu et al. 2012; Keppens & Xia 2014), or cross-field slippage (Low et al. 2012a, 2012b). These mechanisms are not realized in our simulations because we set the cutoff temperature such that the simulated prominence density cannot be sufficiently large.

Zoom In Zoom Out Reset image size

Figure 12 shows the relationship between the temperature and the maximum density of the condensations obtained from our results. The lower the condensation temperature, the higher the density of the condensates. The resultant density is independent of the heating model. The temperature and density of the prominence are related by the following scaling formula,

$$\begin{eqnarray}&&\left({\displaystyle \frac{n}{2\times {10}^{9}\;{\mathrm{cm}}^{-3}}}\right)=85{\left({\displaystyle \frac{T}{{10}^{4}\;\rm{K}}}\right)}^{-1}.\end{eqnarray} \tag{ 33 }$$

The scaling relationship (Equation (33)) is plotted as a solid line in Figure 12. This formula indicates that the dense condensation results from re-stratification along the magnetic field lines with temperature ${T}_{c}$. Because the mass is conserved along the field lines and confined to the length scale of the scale height ${k}_{\rm{B}}{T}_{c}/({{mg}}_{\mathrm{cor}})$, the condensate density is inversely proportional to temperature. From Equation (33), we concluded that our model can reproduce the observed prominence density, which is 10–100 times larger than the surrounding coronal density at typical temperatures (below ${10}^{4}\;\rm{K}$).

Zoom In Zoom Out Reset image size

The direct simulation without the temperature cutoff is necessary to verify the realistic prominence density and compare DEM with observation. Including the chromosphere in our model is also reserved for future work. Due to the thermal conduction from the corona to chromosphere, the background heating rate must be larger than the cooling rate. When the initial coronal density is lower, the magnitude of the initial thermal imbalance is lower, and the present model without the initial thermal imbalance will be validated. It is the next step for our model to be extended to three dimensions. Self-consistent three dimensional MHD simulation including coronal thermal conduction is still a challenging issue. The comparison with observations including the LOS effect (Luna et al. 2012; Xia et al. 2014a) must contribute to the further understanding of the prominence formation.