2D MHD AND 1D HD MODELS OF A SOLAR FLARE—A COMPREHENSIVE COMPARISON OF THE RESULTS

As a 1D numerical code we adopt the modified NRL code (Mariska et al. 1982, 1989), in which we implement the new radiative losses and heating functions, with a mesh of the radiative loss rates calculated for the given range of temperatures (10⁴–10⁸ K) and densities (10⁸–10¹⁴ cm⁻³) with the use of the CHIANTI code (Dere et al. 1997; Landi et al. 2006). The original NRL code is written for the isothermal chromosphere, which provides an insufficient amount of plasma available at the feet of the loop in the course of the massive chromospheric evaporation process. Thus, energy fluxes consistent with the medium or large solar flares result in a massive evaporation of this isothermal chromosphere and in a nonphysical downward motion of the feet even if the chromosphere was set unrealistically thick. In order to solve the problem, likewise in the 2D model (see Section 3.2), we implement the realistic temperature model C7 of the solar atmosphere (Avrett & Loeser 2008), which is smoothly extended downward by the BP04 model of the solar interior (Bahcall & Pinsonneault 2004; Figure 3). This realistic model ensures a big enough storage of plasma within the loop. All other aspects of the adopted NRL model are unchanged; among others, this code does not include viscosity and diffusion, but it takes into account radiative losses of the optically thin plasma, which can be switched off if necessary, and the whole plasma is assumed to be fully ionized. Energy losses due to emission of the optically thick plasma are not taken into account, but Milligan et al. (2014) have already shown that the energy losses due to emission of the flaring loop in the optically thick chromospheric lines range only from 3.5% to 15% of the total energy deposited to the flaring loop (and even less during the impulsive phase solely). The modified code was thoroughly verified by comparison of the obtained results with those obtained with the original code (Mariska et al. 1982, 1989).

Zoom In Zoom Out Reset image size

We developed a set of models, where in the framework of a particular model the energy deposition region has a limited length and volume and is located at a selected, fixed height. The diversity of the positions of the energy deposition volumes reflects various properties of the NTEs' beams, particularly various hardness of the NTEs' energy spectra. The actual position of the energy deposition volume along the loop in the 2D models (see Section 3.2) is selected in the vicinity of the position of the energy deposition region that was estimated with the use of the 1D HD model of the flare (Falewicz et al. 2011). In the case of 2D models, we applied Gaussian distributions of energy deposition in both horizontal and vertical directions. In the horizontal direction a Gaussian distribution of energy deposition is qualitatively similar to the observed spatial distribution of the HXR emission at the feet in real flaring loops. Usually the HXR emission is highest in the centers of these emission sources and decreases toward the edges of the same sources, but an exact distribution of the emission is event dependent. In the vertical direction a Gaussian distribution of energy deposition rates was adopted in both 1D and 2D models as a simplest model of an effective distribution of the NTEs' precipitation depths. The implemented temporal evolution of a delivered energy flux mimics typical variations of energy fluxes observed in solar flares, with gradual growth, and impulsive and gradual decay phases. Such simplification is acceptable as we do not aim to build an exact model of any specific solar flare, but instead our goal is to compare 1D and 2D models of plasma dynamics and evolution of a flaring loop. Guided by this purpose, we construct an initial, quasi-stationary preflare model of the flaring loop, using basic geometrical and thermodynamic parameters estimated from RHESSI and GOES data: half-length, cross section, initial pressure at the base of the transition region (TR), temperature, emission measure, mean electron density, and GOES class (Table 1). The initial pressure at the base is set in order to obtain a static plasma having a temperature at the top of the loop that is close to the observed one. The models are calculated for the periods lasting from the very beginning of the impulsive phase to far beyond the maximum of the SXR emissions, lasting usually for about 300 s.

Two-dimensional numerical models of the solar flare are calculated using the FLASH code (Fryxell et al. 2000; Lee & Deane 2009). In this code third-order unsplit Godunov and Riemann solvers, various slope limiters, and an adaptive mesh refinement are implemented. It was verified by numerical experiments that numerically induced flows remain at physically acceptable levels for the minmod slope limiter and the Roe Riemann solver (Toro 2006), and for the refinement strategy based on controlling the gradients in mass density.

The simulation box is set as (−10 Mm, 10 Mm) × (−10 Mm, 10 Mm), which spans 20 Mm along both the horizontal and vertical directions. In the Cartesian coordinate system, we represent the horizontal axis as "x" and the vertical axis as "y," and all plasma quantities remain invariant along the third direction (z), i.e., ∂/∂z = 0. We impose boundary conditions by fixing in time all plasma quantities at all four boundary layers to their equilibrium values. The only exception is the bottom boundary, where we replace the equilibrium plasma quantities by implementing plasma heating, which was modeled by gas pressure variations defined as

$$\begin{eqnarray}p(x,y,t) & = & p(x)\left(1+{A}_{{\rm{p}}}\displaystyle \sum _{i=1}^{2}\mathrm{exp}\left(\displaystyle \frac{-{(x-{x}_{i})}^{2}}{{w}_{x}^{2}}\right)\right.\\ & & \times \mathrm{exp}\left.\left(\displaystyle \frac{-{(y-{y}_{0})}^{2}}{{w}_{y}^{2}}\right)f(t)\right).\end{eqnarray} \tag{ 1 }$$

Here A_p is the relative amplitude of the pressure signals, x₁ and x₂ are the horizontal positions of these two signals, y₀ is their vertical position, and w_x and w_y are the horizontal and vertical widths of the signals. We set A_p = 60, x₁ = −x₂ = 6 Mm, w_x = 0.5 Mm, and w_y = 1 Mm and hold them fixed but allow y₀ to vary.

Henceforth, time is expressed in seconds and counted from the moment the energy is deposited at the feet of the flaring loop, that is, at t = 0 s. The symbol f(t) denotes the temporal profile of the applied gas pressure signal, which is given by the following formula:

$$\begin{eqnarray*}f(t)=\left\{\begin{array}{lr}1-\mathrm{exp}(-t/\tau ) & :t\lt {\tau }_{\mathrm{max}}\\ \mathrm{exp}(-(t-{\tau }_{\mathrm{max}})/\tau ) & :t\geqslant {\tau }_{\mathrm{max}},\end{array}\right.\end{eqnarray*}$$

where τ = 60 s is the growth and decay time of the time profile. τ_max = 30 s is the time at which f(t) reaches a maximum. The function f(t) is shown in Figure 4.

Zoom In Zoom Out Reset image size

The FLASH code solves a set of nonideal and nonadiabatic MHD equations, which are written in the nonconservative form:

$$\begin{eqnarray}&&\displaystyle \frac{\partial \varrho }{\partial t}+{\rm{\nabla }}\cdot (\varrho {\boldsymbol{V}})=0,\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}\varrho \displaystyle \frac{\partial {\boldsymbol{V}}}{\partial t}+\varrho \left({\boldsymbol{V}}\cdot {\rm{\nabla }}\right){\boldsymbol{V}} & = & -{\rm{\nabla }}p+\displaystyle \frac{1}{\mu }({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}\\ & & +\varrho {\boldsymbol{g}}\ +{{\boldsymbol{F}}}_{{\upsilon }_{{\rm{v}}}},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial p}{\partial t}+{\rm{\nabla }}\cdot (p{\boldsymbol{V}})=(1-\gamma )p{\rm{\nabla }}\cdot {\boldsymbol{V}}+{{\boldsymbol{F}}}_{\kappa },p=\displaystyle \frac{{k}_{{\rm{B}}}}{m}\varrho T,\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t}=\eta {{\rm{\nabla }}}^{2}{\boldsymbol{B}}+{\rm{\nabla }}\times ({\boldsymbol{V}}\times {\boldsymbol{B}}),{\rm{\nabla }}\cdot {\boldsymbol{B}}=0,\end{eqnarray} \tag{ 5 }$$

where ϱ is the mass density, ${\boldsymbol{V}}$ the flow velocity, ${\boldsymbol{B}}$ the magnetic field, p the gas pressure, T temperature, γ = 5/3 the adiabatic index, κ and υ_v the coefficients of thermal conduction and plasma viscosity, ${{\boldsymbol{F}}}_{{\upsilon }_{{\rm{v}}}}$ the viscous term and ${{\boldsymbol{F}}}_{\kappa }$ the thermal conduction term, ${\boldsymbol{g}}$ = (0, −g, 0) gravitational acceleration with its value g = 274 m s⁻², m the mean particle mass, and k_B Boltzmann's constant.

For reasons of simplicity, the medium is assumed to be invariant along the z-direction (∂/∂z = 0), and z-components of the preset background plasma velocity and magnetic field are equal to zero, i.e., V_z = B_z = 0. As a result of this assumption, Alfvén waves are removed from the system, which is able to guide the magnetoacoustic-gravity waves only (e.g., Nakariakov & Verwichte 2005).

We assume that the solar atmosphere in its background state is in static equilibrium (${\boldsymbol{V}}$ = 0) with a force-free magnetic field:

$$\begin{eqnarray}&&({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}=0,\end{eqnarray} \tag{ 6 }$$

and the pressure gradient is balanced by the gravity force,

$$\begin{eqnarray}&&-{\rm{\nabla }}p+\varrho {\boldsymbol{g}}=0.\end{eqnarray} \tag{ 7 }$$

Note that the temperature profile of the initially static plasma determines directly the initial distribution of mass density and gas pressure (see Figure 3). This includes an abrupt decrease of the mass density and pressure in a vicinity of the TR. Using the ideal gas law and the y-component of the hydrostatic balance condition indicated by Equation (7), the equilibrium gas pressure and mass density are given as follows:

$$\begin{eqnarray}&&p(y)={p}_{{\rm{0}}}\;{\rm{exp}}\left[-{\displaystyle \int }_{{y}_{{\rm{r}}}}^{y}\displaystyle \frac{{{dy}}^{\prime }}{{\rm{\Lambda }}({y}^{\prime })}\right],\varrho (y)=\displaystyle \frac{p(y)}{g{\rm{\Lambda }}(y)},\end{eqnarray} \tag{ 8 }$$

where

$$\begin{eqnarray}&&{\rm{\Lambda }}(y)=\displaystyle \frac{{k}_{{\rm{B}}}T(y)}{{mg}}\end{eqnarray} \tag{ 9 }$$

is the pressure scale height, p₀ denotes the gas pressure at the reference level that we choose in the solar corona at y_r = 10 Mm, and y = 0 corresponds to the base of the photosphere.

The force-free condition of Equation (6) is satisfied by the current-free magnetic field, i.e., ∇ × ${\boldsymbol{B}}$ = 0. The solenoidal condition ∇ · ${\boldsymbol{B}}$ = 0 is automatically fulfilled if we introduce the following magnetic flux function ${\boldsymbol{A}}$ = A z:

$$\begin{eqnarray}&&{\boldsymbol{B}}={\rm{\nabla }}\times {\boldsymbol{A}}\left[\displaystyle \frac{\partial A}{\partial y},-\displaystyle \frac{\partial A}{\partial x},0\right].\end{eqnarray} \tag{ 10 }$$

Here we choose the z-component of ${\boldsymbol{A}}$ as

$$\begin{eqnarray}&&A(x,y)={B}_{{\rm{0}}}{\lambda }_{{\rm{B}}}\ {\rm{exp}}\left(-\displaystyle \frac{y}{{\lambda }_{{\rm{B}}}}\right)\mathrm{cos}\left(\displaystyle \frac{x}{{\lambda }_{{\rm{B}}}}\right),\end{eqnarray} \tag{ 11 }$$

where λ_B = L/π is the magnetic scale height, with L being half of the arcade width. We choose and hold fixed L = 19 Mm and specify the magnetic field B₀ by requiring that at the reference point (x = 0 Mm, y = 10 Mm) the Alfvén speed, ${c}_{{\rm{A}}}(x,y)=| {\boldsymbol{B}}(x,y)| /\sqrt{\mu \varrho (y)},$ is 10 times larger than the speed of sound, ${c}_{{\rm{s}}}(y)=\sqrt{\gamma p(y)/\varrho (y)}.$ For these settings ${B}_{{\rm{0}}}\cong 11.4$ G and the magnetic field is predominantly horizontal around x = 0 Mm, while around (x = 11, y = 1.75) Mm it is oblique with a π/4 angle to the solar surface. As a result, the topology of the magnetic field mimics the closed and bent field lines of the active region magnetic arcade. The plasma parameter $\beta =P/({B}^{2}/2\pi )$ is equal to 0.016 at y = 3 Mm and 0.028 at y = 5 Mm. As β < 1, the trajectories of the evaporating plasma are fully controlled by the magnetic field lines. The topology of the modeled flaring loop, i.e., separation of the feet, half-length, and diameter, is not predefined but results directly from the actual position and size of the heating regions (Figure 5).

Zoom In Zoom Out Reset image size

A complex, multiscale structure of the adaptive mesh is applied, which has the finest spatial size in the vicinity of steepest profiles of the mass density and in the region of interest that is chosen below the altitude y = 11 Mm (Figure 6). As every block is divided into 8 × 8 numerical cells, on the first level of grid refinement a spatial resolution is equal to 2.5 Mm and is refined gradually into 2 × 2 twice smaller numerical cells, forming a nonuniform numerical grid. The resulting grid is coarse in the solar corona but very fine in the low corona and the chromosphere, with the finest spatial resolution of 0.08 Mm (Figure 6).

Zoom In Zoom Out Reset image size

A transient energy of the solar flare is delivered to the flaring loop by the beams of high-energy NTEs that deposit their energy in the dense plasma, usually in the vicinity of the TR, in the chromosphere, and sometimes even in the photosphere. The local depositions of energy are defined and modeled by a variable-in-time, volumetric increase of the plasma pressure, launched in a finite volume of the magnetic structure at the bottom boundary of the numerical box (see Equation (1) and cf. Figure 6). The pressure signal is of Gaussian shape in both horizontal and vertical directions, and its amplitude, A_p, is 60 times higher than the local ambient pressure. Temporal variations of the pressure mimic a typical time profile of the X-ray emission of a solar flares (Figure 4). However, the amount of the delivered energy depends on actual physical parameters of the plasma. It is equal to 10,900 erg cm⁻³ at the height of h = 1.24 Mm above the temperature minimum of the model atmosphere and is equal to only 2685 erg cm⁻³ at h = 1.74 Mm. While the actual deposited energy depends both on the parameters of the NTE beams and on the local plasma characteristics, we investigate a set of models, where the altitudes of the energy deposition regions vary between h = 1.2 Mm and h = 1.8 Mm above the temperature minimum with the spatial step of Δh = 0.1 Mm. This is in accordance with our former results of our previous investigation of the same flare (Falewicz et al. 2011), in which the height of the heating region was defined as a position of the maximum of the heating caused by NTEs and estimated using Fisher's approximation (Fisher 1989).

In this part of the paper we consider a flare with an energy deposition region located in the upper chromosphere, h = 1.24 Mm above the temperature minimum. Owing to the initial Gaussian spatial distribution of the local pressure impulse, which mimics the local energy deposition of the NTEs during the early phase of the flare, the evaporating plasma rises toward the loop top in all the models, forming in the 2D models the heterogeneous columns of hot and dense cores and relatively cool and less dense ambient layers. However, an exclusion of the thermal conduction in the IMHD models exaggerates even more local differences in temperature perpendicular to the loop's axis. When the chromospheric evaporation is fully developed, i.e., about 20 s after the onset of the noticeable X-ray emissions of the loop in the 1–8 Å band, the highest temperature of the evaporating plasma at the representative height of 1.5 Mm above the TR is obtained in the NRL model (about T_NRL = 2.5 MK), while the temperature in the axial part of the loop is of the order of T_IMHD = 1.3 MK in the case of the IMHD model and ${T}_{\mathrm{VDC}}\cong 0.1$ MK in the case of the VDC model. Representative temperatures of the ambient layer of the loops in both 2D models are of the order of chromospheric temperatures only, while these layers of the loop are formed by chromospheric-type plasma pushed from below without any significant heating. Mass density of the hot plasma is ${\varrho }_{\mathrm{NRL}}\cong 6\times {10}^{-14}\;{\rm{g}}\;{\mathrm{cm}}^{-3}$ in the NRL and in the axial part of the 2D VDC model, but it is 4 times lower in the axial part of the 2D IMHD model. The latter may be due to neglected viscosity in this model. The raising pillars of plasma caused compression of the plasma contained in the upper part of the loop, gradually forming hot and dense loop-top sources of X-rays. Twenty seconds after the onset of the X-ray emission, their temperature is still moderate and equal to T_IMHD = 1.3 MK, T_VDC = 1.5 MK, and T_NRL = 2 MK in the IMHD, VDC, and NRL models, respectively.

The X-ray 1–8 Å emission attains its maximum at t_maxIMHD = 130 s in the IMHD model, at t_maxVDC = 110 s in the VDC model, and at t_maxNRL = 150 s in the NRL model. The overall properties of still-ascending pillars of the evaporating plasma at that time remain roughly similar to the plasma properties during the previous stage of evolution. The highest temperature of the evaporating plasma at the representative height of h = 1.5 Mm above the TR is again noted in the NRL model (about 3.7 MK), while the temperatures in the axial part of the loop are of the order of T_IMHD = 1.1 MK in the case of the IMHD model and T_VDC = 0.1 MK only in the case of the VDC model. At this evolution stage, the highest temperatures of the plasma occur in the well-formed but still compressed loop-top regions of the loops, co-spatial with the loop-top sources of the X-ray emissions, where the temperatures reach T_IMHD = 3.5 MK in the IMHD model, T_VDC = 1.5 MK in the case of the VDC model, and peaked at T_NRL = 5.2 MK in the NRL model. Plasma densities are ϱ_IMHD = 4 × 10⁻¹⁴ g cm⁻³, ϱ_VDC = 2 × 10⁻¹⁴ g cm⁻³, and ϱ_NRL = 1.4 × 10⁻¹³ g cm⁻³. For both 2D models, the calculated temperatures of the loop-top plasma are in a good agreement with the temperature derived from the calculated synthetic X-ray fluxes emitted by the modeled loops, i.e., GOES-like temperature, being equal to T_IMHD-GOES = 3.5 MK and T_VDC-GOES = 3.4 MK for the IMHD and the VDC models, respectively. The GOES-like temperature describes an average or representative temperature of the flaring plasma. According to the algorithm of White (2005), the GOES-like temperature is based on the ratio of the integral X-ray fluxes recorded in the 0.5–4 Å and 1–8 Å bands by the GOES satellite. In the NRL model the same temperature peaks 40 s before the maximum of the emitted 1–8 Å flux, reaching T_NRL-GOES = 5.2 MK at T_NRL = 110 s. A nearly equal maximum of the GOES-like temperature occurs at ${t}_{\mathrm{NRL}}\cong 150\;{\rm{s}},$ during the maximum of the emission in the 1–8 Å band.

During the late, gradual phase of the evolution, at t = 200 s in the case of the IMHD and VDC models and t = 250 s in the case of the NRL model, the loop-top regions expanded, and the plasma began to flow down along the legs of the loops in the NRL model, as well as in the central parts of the legs in the 2D models, with the velocity of 40–60 km s⁻¹ for the IMHD model, 10–30 km s⁻¹ for the VDC model, and 20–40 km s⁻¹ for the NRL model. However, in the case of 2D models the plasma in the ambient layer of the loop still moves toward the top of the loop, with the velocity of the order of 10 km s⁻¹. The temperature of the expanding plasma in the loop-top region gradually decreases, reaching at that moment a value of T_IMHD = 1.5 MK in the IMHD model, T_VDC = 1.8 MK in the case of the VDC model, and T_NRL = 3.2 MK in the NRL model, which is still the highest amount of the various models we discussed. The X-ray emission of the loop-top source vanished at that moment, and the X-ray emission becomes dominant by the emission of the plasma in the legs of the loops. At the height of h = 1.5 Mm above the TR, the plasma temperature in the central part of the loop reaches ${T}_{\mathrm{IMHD}}\cong 1.2$ MK in the case of the IMHD model and T_VDC = 0.2 MK only in the case of the VDC model. The temperature of the plasma at the same height in the NRL model is equal to T_NRL = 1.3 MK. Plasma densities are equal to ϱ_IMHD = 2 × 10⁻¹⁴ g cm⁻³, ϱ_VDC = 5 × 10⁻¹⁴ g cm⁻³, and ϱ_NRL = 1.8 × 10⁻¹³ g cm⁻³. The representative temperature of the ambient plasma layer of the loops in both 2D models is of the order of the chromospheric temperatures only. As a result of an isotropic thermal conduction applied in the FLASH code, the total width of the hot loop gradually grows in time when modeled in 2D, reaching a value of about 2500 km after 200 s of its evolution.

In this part of the paper we consider the flare with an energy deposition region located in the lower corona, h = 1.74 Mm above the temperature minimum (see Figure 3). It is heated in the same way as in the model described in Section 4.1. Spatial distributions of the temperature and mass density are presented in Figures 10 and 11, and representative values of the main physical parameters of the plasma are presented in Table 2. Resulting X-ray fluxes in the 0.5–4 Å and 1–8 Å bands and derived GOES-like temperatures and densities are presented in Figure 13. The beginning of the noticeable X-ray emissions of the loop in the 1–8 Å band occurs in all models (IMHD, VDC, and NRL) on much faster timescales than in the case of the models with the heating region located at the height of h = 1.24 Mm above the TR, namely, at ${t}_{\mathrm{IMHD}}=9\;{\rm{s}}$ for the IMHD model, t_VDC = 8 s for the VDC model, and t_NRL = 5 s for the NRL model (see Figure 13, right panels).

Zoom In Zoom Out Reset image size

The global evolution of the main physical parameters of the plasma confined in the flaring loop occurs in the same way in all models with the heating region at h = 1.74 Mm above the temperature minimum. During the early phase of the flare, the evaporating plasma rises toward the loop top, forming in 2D models the columns of very hot plasma in cores and much cooler plasma in the outside layers. When the chromospheric evaporation is fully developed (again, about 20 s after the onset of the noticeable X-ray emissions of the loop in the 1–8 Å band), the highest temperature of the evaporating plasma at the representative height of h = 1.5 Mm above the TR is noted for the IMHD model, being equal to T_IMHD = 10 MK, while for the VDC model, the temperature is T_VDC = 2.1 MK. In the case of the NRL model, the temperature of the plasma at h = 1.5 Mm is equal to the intermediate value of T_NRL = 5.6 MK.

Significant differences between the temperatures calculated in various models result from application or neglecting various physical mechanisms of energy redistribution such as viscosity and thermal conduction. For instance, exclusion of thermal conduction results in the strong local overheating of the loop core in the IMHD model. The representative temperatures of the outside layer of the loops in the IMHD and VDC models are T_IMHD = 0.1 MK and T_VDC = 0.5 MK, respectively. The very hot plasma in the core of the loop of the IMHD model at h = 1.5 Mm is at the same time very rarefied, with ϱ_IMHD = 5.6 × 10⁻¹⁶ g cm⁻³. The cold and dense envelopes at the edge of the loop resulted from the applied model of plasma heating, with Gaussian-like spatial distribution of the impulse, where plasma on the edges is less heated and expands less than in the central part of the loop. However, an exclusion of thermal conduction in the IMHD model prevents any energy transfers perpendicular to the loop's axis and exaggerates even more local differences in temperature. Much colder plasma in the core of the loop of the VDC model has ϱ_VDC = 2.5 × 10⁻¹⁵ g cm⁻³, but the densest plasma occurs in the NRL model at the same stage of the evolution, reaching ϱ_NRL = 7.6 × 10⁻¹⁴ g cm⁻³. Note that for the IMHD and VDC models at the same stage of the evolution massive upflows of the evaporating plasma take place, while in the NRL model the plasma is still static at the reference level of h = 1.5 Mm above the TR, contemporarily moving up in the lower parts of the loop. Conversely to the models described above, for the 2D models the rising pillars of the plasma do not cause any important compression of the plasma present in the upper parts of the loop, where loop tops remain moderately cold, with T_IMHD = T_VDC = 0.5 MK. As a result of plasma substantial compression, the top of the loop in the NRL model is much hotter, with T_NRL = 10.6 MK. This is consistent with the observational Soft X-ray Telescope (SXT/Yohkoh) data, which show the temperature of the loop-top sources calculated from the ratio of images in the two thickest filters to be typically ≲15 MK (McTiernan et al. 1993; Doschek 1999).

The maxima of the X-ray 1–8 Å emissions occur at t_maxIMHD = 90 s in the IMHD model, at t_maxVDC = 70 s in the VDC model, and at t_maxNRL = 150 s in the NRL model. The substantial difference between relatively fast occurrence of the maxima in the IMHD and VDC models and delayed maximum in the NRL model results mostly from the very extended in time but rather slow increase of the 1–8 Å emissions in the NRL model, which still lasts between t = 90 s and t = 150 s. However, during the early phases of the evolution of plasma in the framework of all three models, the calculated 1–8 Å fluxes rise in a very similar way and up to roughly t = 70 s. The difference is caused by differences in the dynamics of the plasma motions and densities along the loops in 1D and 2D models, with the mass densities of the plasma in the loop top being equal to ϱ_IMHD = 8 × 10⁻¹⁵ g cm⁻³, ϱ_VDC = 1 × 10⁻¹⁴ g cm⁻³, and ϱ_NRL = 2.6 × 10⁻¹³ g cm⁻³. In the case of the IMHD model, the plasma flows downward along the central part of the loop, but it still rises in the ambient layer toward the top of the loop. In the VDC model, the plasma flows down through the whole legs of the loop, while in the NRL model it moves toward the loop top.

For both 2D models, the calculated plasma temperatures at the loop top are equal to T_IMHD = 1.8 MK and T_VDC = 2.5 MK, and they are much lower than the temperatures derived from the calculated synthetic X-ray fluxes emitted by the modeled loops. For example, the GOES-like temperatures are equal to T_IMHD-GOES = 13 MK and T_VDC-GOES = 7 MK for the IMHD and VDC models, respectively. For the NRL model, the temperature of the loop-top plasma is equal to T_NRL = 11.5 MK at t = 150 s, and it almost fits the GOES-like temperature of T_NRL-GOES = 10.5 MK. The GOES-like temperature peaks in the NRL model about 100 s before the maximum of the emitted 1–8 Å flux, reaching T_NRL-GOES = 15 MK at t = 55 s.

During the late, gradual decay phase of the evolution, at t = 200 s in the case of the IMHD and VDC models and t = 250 s in the case of the NRL model, the loop-top regions expand, and the plasma flows down along the legs of the loops in all the models. The temperature of the expanding plasma remains on the constant level of T_IMHD = 1.8 MK in the IMHD model, or slowly decreases, reaching at that moment T_VDC = 1.4 MK and T_VDC = 7.3 MK for the VDC and NRL models, respectively. Again, the difference in the temporal variations of the loop-top temperatures results from differences in the dynamics of the plasma motions and mass densities along the loops in the 1D and 2D models. At the representative height of h = 1.5 Mm above the TR, the plasma temperature in the axial part of the loops is of the order of T_IMHD = 10 MK in the case of the IMHD model and only of T_VDC = 4.5 MK in the case of the VDC model. For the NRL model the temperature at the same height is equal to T_NRL = 4.1 MK. Plasma densities are ϱ_IMHD = 1 × 10⁻¹⁵ g cm⁻³, ϱ_VDC = 1.6 × 10⁻¹⁵ g cm⁻³, and ϱ_NRL = 4.8 × 10⁻¹³ g cm⁻³. Besides, the representative temperatures of the ambient plasma in all the 2D models remain relatively high, and they are equal to T_IMHD = 3.1 MK and T_VDC = 1.2 MK, respectively.

A comparative analysis of the IMHD and VDC models reveals an overall similarity of temporal variations of the most basic flare parameters, such as GOES light curves, derived temperatures, and emission measures. In the case of models that are heated in the upper chromosphere, the IMHD model leads to high-temperature apex volume and relatively cold upper parts and hot lower parts of the legs. The temperature of the loop-top for the VDC model is similar to that for the IMHD model, but the temperatures along whole legs are essentially homogeneous, owing to effective thermal conduction. Temperature gradients for the IMHD model are inevitably very large, while these gradients for the VDC model are smoothed as a result of the thermal conduction. Numerous small-scale variations of the plasma temperature observed for the IMHD model at the legs of the flaring loop are not seen for the VDC model, being smoothed out efficiently by viscosity and thermal conduction. Comparing both models at the same evolution time, the volume of the apex X-ray source for the VDC model is smaller than for the IMHD model. Local mass density in the loop-top source is relatively low but with abrupt spatial variations for the IMHD model, while for the VDC model the mass density is higher, and all small-scale variations are smoothed by viscosity. These differences are particularly well seen for time between t = 60 s and t = 90 s.

In the case of the models with heating above the TR temperature, gradients in the IMHD model are inevitably very strong, while these gradients in the VDC model are smoothed as a result of thermal conduction action. Similar differences could be noticed also for mass densities, which show abrupt local variations for the IMHD model but have a relatively smooth distribution for the VDC model. The solely strong mass density variations for the VDC model during the late phase of the evolution occur on the forehead edge of the expanding apex kernel, where the matter is moving down. In both models the loop-top sources exhibit similar volumes.

The applied model of local plasma heating mimics positions and limited volumes of plasma in the legs of the flaring loop, where the NTEs effectively deposit their energy in real flares. However, in our models energy is not deposited in any other part of the flaring loop, including also the loop-top region. Temporal variations of the deposited energy flux are modeled by variations of the local pressure amplitude in rough accordance and, at least qualitatively, with a typical time profile of the HXR emission observed for a typical solar flare.

Despite local variations of the main physical parameters of the plasma, such as temperature, actual local mass density, etc., discussed above, and plasma parameters derived from simulated X-ray emissions of the flare, i.e., GOES-like temperatures and emission measures, the most important difference between modeled and real solar flares is the higher brightness of the loop legs in X-ray than the brightness of the loop-top sources—quite the opposite of the real observations (see Figure 15). The reason is the absence of a heating component of the loop-top source, originated from NTEs, which was not included in the simulation and low plasma mass density despite its relatively high temperature. This striking difference could be attributed to an absence of heating of the plasma deposited along the whole loop, especially in the loop-top region. The modeled densities of the apex volumes are also relatively low, being, despite relatively high temperatures, a not very efficient source of X-rays. However, the model also appropriately fits to the compact loop flares, where the flare is triggered near the loop's feet during low atmospheric reconnection in the TR or inner corona, and heats the plasma to excite its evolution and multiwavelength emissions (EUV, SXR, etc.; e.g., Kumar et al. 2011).

Zoom In Zoom Out Reset image size

Specific spatial distribution of the applied heating model, being a reasonably good model of the spatial distribution of the energy flux inside a beam of the NTEs and automatically mimicking a multithread structure of the real flaring loops, causes the core parts of the modeled loops to be usually hotter and denser than ambient regions. In the case of other specific spatial distributions of the heating energy, the resulting courses of the modeled flares could be substantially different.

The heating of the spatially limited volume of the plasma and during finite periods leads to the evolution of the heating pulses, which we presented here as Gaussian. The implementation of such pulses in the full MHD model results in their propagation in the stratified atmosphere, and under suitable ambient conditions these thermal perturbations generate associated slow waves that can be dissipated by thermal conduction (Ofman & Wang 2002). It should be noted that the release of plasma heating results in perturbed gas pressure or temperature and mass density. The evolution of Alfvén waves is likely not allowed in the considered model. The combination of gas pressure and magnetic perturbations during the flare process may trigger fast magnetoacoustic-gravity waves. Within the given timescale of a few tens of seconds, these waves leave the simulation domain, thus essentially not contributing to energy transfer on the given timescale of a few hundred seconds. In the low β plasma slow magnetoacoustic waves are guided along magnetic field lines, and they may contribute to the partial energy fulfillment of the localized plasma of the medium, where they evolve (see theoretical papers by De Moortel & Hood 2003; Dwivedi & Srivastava 2008). For instance, in the present flaring loop system, where the plasma pressure is perturbed by flare energy release, the evolution of slow magnetoacustic waves is likely along the magnetic field lines, and the presence of thermal conduction and viscosity results in their dissipation (Dwivedi & Srivastava 2008), leading to the plasma heating. Under suitable conditions, some of these waves may steepen in the form of shocks, which propagate along the magnetic field lines (see Srivastava & Murawski 2012). Therefore, in conclusion, slow magnetoacoustic waves may dissipate and heat the plasma.

Numerous previous hydrodynamic simulations of flaring loops failed to reproduce high temperatures and mass densities typical for the gradual phase of the solar flares. In particular, during the gradual phase of the flare a much slower decay of the SXR flux is usually observed than could be expected from analytical estimation of energy losses of the flaring plasma by radiation and by conduction (Cargill et al. 1995) or from HD numerical modeling of the flares (e.g., Reale et al. 1997), proving that during this phase of the flares significant heating is also present. Similar conclusions were already drawn by McTiernan et al. (1993) and Jiang et al. (2006). In order to solve this problem, Warren (2006) applied multithread, time-dependent 1D hydrodynamic simulations of solar flares as a multitude of tiny, parallel loops for investigations of the so-called Masuda flare, which took place on 1992 January 13. He showed that modeling a flare as a sequence of independently heated threads instead of as a single loop may resolve the discrepancies between the simulations and observations, but at the cost of some ad hoc parameters (number of loops, energy share by particular loops, etc.). Additionally, Falewicz et al. (2011) and Falewicz (2014) showed that energy delivered by NTEs was fully sufficient to fulfill the energy budgets of the plasma during the preheating and impulsive phases of analyzed flares, as well as during the decay phase. They concluded that in the case of the investigated flares there was no need to use any additional heating mechanisms other than heating by NTEs. These two results allow one to model a flaring loop as a macroscopic loop having a multithread internal structure and heated by a variable in time energy flux transported to the feet of the loop by NTEs only. However, the multithread internal structure of the flaring loop is poorly resolved by contemporary instruments, and there are no means to investigate individual heating episodes and thermodynamic evolution of the individual threads. These observational limitations cause fundamental problems in numerical modeling of multithread flaring loops, including a proper selection of a relevant structure and number of treads, i.e., filling factor of the observed macroscopic loop, as well as a proper model of spatial and temporal energy deposition in the threads.

In this work we analyzed thoroughly a broad set of 1D HD and 2D MHD models of a solar flare, derived in order to compare and understand the energy transfer, plasma redistribution, and physical conditions of the various models of the flaring loops. Initial physical and geometrical parameters of the flaring loop, applied as a working model of the typical solar flare, are set according to the observations of an M1.8 flare recorded in the AR 10126 active region on 2002 September 20 between 09:21 UT and 09:50 UT. The modeled loop was embedded into the realistic temperature model of the solar atmosphere, smoothly extended downward below the photosphere. The nonideal 1D models include thermal conduction and radiative losses of the optically thin plasma as energy-loss mechanisms, while the nonideal 2D models include viscosity and thermal conduction as energy-loss mechanisms only. The global evolution of the main physical properties and associated parameters of the plasma confined in the flaring loops occurred in the same way in all the models. The thorough comparison of the main physical characteristics of the 1D HD and 2D MHD models of the flare shows that the basic properties of the flaring plasma are acceptably well reproduced using both models, ensuring that results obtained already using a simplified 1D approach remain valuable for understanding the chief physical properties of the solar flares.

The IMHD models, often applied in modeling of the solar flares, apply a very simplified approximation of the real processes occurring during the flares, as well as of local and macroscopic properties of flaring plasma, if compared with much more realistic VDC models. Some nonrealistic properties of plasma modeled in the IMHD models are thus the direct effect of basic assumptions and limitation of the IMHD methodology. However, IMHD models can be applied, with great success, for investigations of some problems, such as interactions of the waves in the loops having a nonstiff magnetic field or plasma–wave interactions.

We show that simplified 2D models, having a continuous distribution of the physical parameters of the plasma across the loop and powered by a heating flux variable in time as well as along and across the loop, are an extreme borderline case of a multithread internal structure of the flaring loop with a filling factor equal to 1. Therefore, such models might mimic to some extent the subtle structure and variations of the plasma parameters better than their 1D counterparts, revealing processes that are inherently absent in 1D models. However, the whole complexity in time and space of the processes occurring in the real multithread flaring loops is far beyond their scope.

In our approach used in the present work the 2D energy deposition region located inside the limited volume of the magnetic field automatically defines a multithread loop (or multithread arcade) in a preexisting magnetic field, filled in various ways by evaporating plasma. In this approach, the only parameters to be selected/determined are the size of the heating area (and spatial distribution of the energy flux in a frame of the source, if necessary) and the amount of delivered energy.