Extreme-ultraviolet and X-Ray Emission of Turbulent Solar Flare Loops

Turbulence has been observed in flare loops and is believed to be crucial for the acceleration of particles and in the emission of X-ray photons in flares, but how the turbulence is produced is still an open question. A scenario proposed by Fang et al. suggests that fast evaporation flows from flare loop footpoints can produce turbulence in the looptop via the Kelvin–Helmholtz instability (KHI). We revisit and improve on this scenario and study how the KHI turbulence influences extreme-ultraviolet (EUV) and X-ray emission. A 2.5D numerical simulation is performed in which we incorporate the penetration of high-energy electrons as a spatio-temporal dependent trigger for chromospheric evaporation flows. EUV, soft X-ray (SXR), and hard X-ray (HXR) emission are synthesized based on the evolving plasma parameters and given energetic electron spectra. KHI turbulence leads to clear brightness fluctuations in the EUV, SXR, and HXR emission, with the SXR light curve demonstrating a clear quasi-periodic pulsation (QPP) with period of 26 s. This QPP derives from a locally trapped, fast standing wave that resonates in between KHI vortices. The spectral profile of the Fe xxi 1354 line is also synthesized and found to be broadened due to the turbulent motion of plasma. HXR tends to mimic the variation of SXR flux and the footpoint HXR spectrum is flatter than the looptop HXR spectrum.


Introduction
Looptop hard X-ray (HXR) sources in flare loops are frequently reported in observations of solar flares (e.g., Masuda et al. 1994). It is widely believed that the HXR photons from the looptop source mainly come from bremsstrahlung (e.g., Kontar et al. 2011). Energetic electrons play a crucial role in this mechanism, as most of the energy of the HXR photons comes from electron bremsstrahlung. Stochastic acceleration by turbulence is suggested to be an efficient way to produce high-energy electrons in the environment of solar flares (e.g., Miller et al. 1996). In observations, turbulence has been observed in multiple flare events (e.g., Doschek et al. 2014;Kontar et al. 2017). Furthermore, turbulent magnetic fields can trap fast electrons into the loop apex, which helps increase the flux of HXR emission.
Although turbulence influences the looptop HXR source, how this turbulence is generated is still an open question. Fang et al. (2016) provided a possible answer to this question: the observed turbulence can be produced by Kelvin-Helmholtz instability (KHI). The KHI is triggered by fast and dense chromospheric evaporation flows. In typical coronal loops, we find low plasma density and low plasma beta conditions (much lower than 1), and we will initiate our model with such conditions. It has been noted that KHI is then likely triggered in the plane perpendicular to the axial magnetic field (Terradas et al. 2010;Antolin et al. 2014). Different from normal coronal loops, flare loops become filled with high-density plasma (10 10 -10 12 cm −3 ), which leads to a beta higher than 1 and a low Alfvén speed of several hundreds km s −1 . Because evaporation flows with upflow speed up to 500 km s −1 are not rare in solar flares (Nitta et al. 2012), it is then also possible to get KHI triggered by evaporation flows in flare loops. Ofman & Sui (2006) also suggested that KHI can be triggered in flare loops, although in their 2.5D scenario, KHI is triggered by reconnection outflows rather than by evaporation flows. Furthermore, KHI triggered in a solar blowout jet has been recently observed in Li et al. (2018) in an active region, where a loop like structure is rolling up at its edges. Ruan et al. (2018) surveyed the evaporation flow based trigger of KHI parametrically, and found that KHI turbulence is very likely to appear in short-duration, impulsive flares.
In this Letter, we study how KHI turbulence influences the emission of flare loops. The production of KHI turbulence triggered by evaporation flows in a flare loop is simulated and then extreme-ultraviolet (EUV), soft X-ray (SXR), and HXR views on the flare loop are studied. We improve on Ruan et al. (2018) and Fang et al. (2016) in a variety of ways: instead of a parametrized heating of the chromospheric layers, we now extend established 1D models (Allred et al. 2005) for energetic electron fluxes to our 2.5D study. We then obtain a fully realistic, time-and space-dependent heating source that triggers the evaporation flows in the flare loop. How the evaporation is produced, and how the HXR emission is synthesized, is introduced in Section 2. The main results are analyzed in Section 3 and summarized in Section 4.

Model Improvements
We performed our simulation with the MPI-parallelized Adaptive Mesh Refinement Versatile Advection Code (MPI-AMRVAC; Keppens et al. 2012;Porth et al. 2014;Xia et al. 2018). Magnetohydrodynamic (MHD) equations are solved in a 2.5D simulation box with a domain of  -30 Mm  x 30 Mm and 0 Mmy40 Mm. Adaptive mesh refinement is used and the smallest grid cells have a size of about 20 km×20 km. We use the same governing equations as in Ruan et al. (2018), in which anisotropic thermal conduction, gravity, a parametric, constant background heating, and optically thin radiative cooling are included. The expression of the initial force-free magnetic field is the same as in Ruan et al. (2018), but the horizontal size of the simulation box is reduced to L 0 =60 Mm and the angle between the apex of the magnetic loops and the neutral line is changed to θ 0 =60°to get similar loop configuration in x-y plane. Photosphere, chromosphere, transition region, and corona are included.

Flare Energy Deposition
A magnetic flux loop is selected to deposit flare energy from the overlying reconnection site, which is deliberately not (yet) included in our model. The footpoints of this flare loop are located at 21 Mmx22 Mm and −22 Mmx−21 Mm at the lower boundary. Flare energy is deposited into this loop based on the (changing) density profile along magnetic field lines. This reproduces how chromospheric footpoints are heated due to penetration of high-energy electrons. The energy flux carried by downward fast electrons at a height of h=10 Mm (where the loop width is 2 Mm) is where t 0 is the time when the energy flux reaches a maximum value, and we set τ 1 =30 s and τ 2 =180 s. The energy flux amplitude F 0 is set to 5.5×10 10 erg cm −2 s −1 for the left foot, and 6.5×10 10 erg cm −2 s −1 for the right foot, to achieve asymmetric heating. The change of downward energy flux due to the changing loop width has been taken into account. Multiple magnetic field lines through both footpoint regions are selected and their position is calculated and recorded. This is used to quantify the footpoint heating added as a spatiotemporal heat source in the simulation. The heating rate below h=10 Mm is calculated using the method of Allred et al. (2005). We assume that downward fast electrons have a double power-law energy distribution and an initial pitch angle of 60°. The energy distribution has a spectral index δ l =3 above the cutoff energy E c =20 keV and below the break energy E b =100 keV. The spectral index above E b is set to δ u =4. According to Allred et al. (2005), the energy distribution of input fast electrons is then and the heating rate is given by , m 0 is the cosine of the initial pitch angle that we fix to 0.5, The (minute) motion of the field lines in the bottom loop parts during the simulation is neglected for this quantification. The plasma is assumed to be fully ionized and consist of protons and electrons to calculate the heating rate. Figure 1 displays instantaneous profiles of number density and heating rate along the left loop axis at t=0 s and t=22 s: the heating rate adapts to the change of plasma density.

HXR Synthesis
HXR emission from few keV to a few hundred keV in flares is determined by electron-ion bremsstrahlung (Kontar et al. 2011). The HXR spectrum j(ò) (photons cm −3 s −1 keV −1 ) produced by electron-ion bremsstrahlung is calculated with where n i is the number density of target ions, ò is photon energy, and E is electron energy. The method to calculate the local fast electron spectrum F l (E) parallel to field line direction and corresponding cosine pitch angle μ(E) will be introduced below. The Kramers bremsstrahlung cross section is adopted (Kontar et al. 2011), given by where Z is ion charge and σ 0 =7.9×10 −25 cm −2 keV. The looptop has a density of ∼10 10 cm −3 . Plasma in this region can be regarded as thin target for electrons above E c =20 keV, thus the change of spectral profile can be neglected. We assume that the fast electrons spectrum F l (E) in regions where T>5 MK but Q e =0 (looptop) always has the double power-law distribution from Equation (2). The local energy flux carried by the electrons is assumed to be proportional to the local magnetic field magnitude: (1) is set to 6×10 10 erg cm −2 s −1 . This assumption ensures that the number of fast electrons is conserved when the electrons move along the field lines. The value of μ in Equation (4) is set to 0.5 and the change of pitch angle in this upper loop region has not been taken into account.
In contrast, regions where the heating rate Q e >0 (footpoints) are treated as thick target for fast electrons. Fast electron spectra F l (E) have a distribution of Equation (2) at h=10 Mm and become flatter and flatter owing to collisions when the electrons move downward. Here we estimate the local electron spectrum with an equation that describes the evolution of their average energy in high-density plasma (Emslie 1978). We use an equation stating how initial particle energy E 0 and initial pitch angle μ 0 (set 0.5, other parameters refer to Equation (3)) vary with column depth through In this thick target interaction, pitch angle changes must also be accounted for (Emslie 1978), and the cosine of the pitch angle in Equation (4) of a particle with initial energy E 0 at column depth N is given by Starting with the given distribution at h=10 Mm, we then obtain the evolution of the spectrum for varying column depths N. We select discrete energy values E i (0)=20+i keV (for i=0, 1, 2, K 280) and obtain energies E i (N) according to Equation (6) by conserving the integral under the energy spectrum: the average spectral density In this way, we can get the local fast electron spectra in the thick target footpoint regions. The local HXR spectra can then be estimated with the information of F l and μ. The flare loop is assumed to have a depth of 5 Mm in the out-of-plane direction.

Results
When flare energy is deposited into the footpoints, evaporation flows with high density, high temperature, and high speed are produced. The plasma beta inside the loop increases from its initial ∼0.001 to order unity and above, because of the evaporation of plasma. Evaporation flows interact with each other and produce KHI turbulence in this plasma environment. In contrast, Ofman & Sui (2006) set the initial beta to 4 in their simulation of the triggering of KHI by reconnection outflows in flare loops. The KHI we find here is of a completely different origin, and is fully consistent with expectations from linear theory and nonlinear simulations Henri et al. 2013). The evolution of the flare loop is similar to that in case 1 of Ruan et al. (2018), and is shown in X-rays in Figure 2.  field regions develop and suppress the effect of thermal conduction, leading to an uneven distribution of temperature. SXR images determined by thermal bremsstrahlung tend to be brighter in regions where the magnetic field is weaker, as plasma is hotter there. In contrast, HXR images determined by nonthermal bremsstrahlung are brighter in regions where magnetic field is stronger (e.g., at the footpoints), as converging magnetic field lines increase the flux of fast electrons and the density of plasma.

SXR and HXR Images and Light Curves
We find a clear quasi-periodic oscillation with period of about 26 s in the temporal profile of the SXR flux (Figure 3, left panel). The oscillations are obvious when t>4 minutes, and last for multiple periods. In the timeframe indicated between the red dotted lines, we compute the normalized correlation coefficient between the detrended total SXR flux and the local SXR flux all throughout the arcade. This shows us clearly where the quasi-periodic oscillation resides (Figure 4(a)): regions where this normalized correlation coefficient is higher than 0.7 are mainly near the apex center. The average Alfvén speed of the regions is about 60 km s −1 and the average acoustic speed is about 500 km s −1 , making the loop apex a region of high plasma beta. When a square area (1 Mm× 1 Mm) in this region is selected (red box in Figure 4(a)), the 26 s quasi-periodic oscillation is also obvious in the profile of average density and temperature (Figure 4(b)). Furthermore, the average density and temperature vary in phase with the average SXR flux.
Oscillations with periods shorter than 60 s in flare loops are normally interpreted as (global) standing fast sausage modes (Nakariakov et al. 2003;Melnikov et al. 2005;Inglis et al. 2008;Van Doorsselaere et al. 2011;Su et al. 2012;Chowdhury et al. 2015;Kolotkov et al. 2015;Tian et al. 2016). A model often used to estimate the wave speed is v p =2L/P (Aschwanden 2004), where L is the length of a loop and P is the period of the oscillation. The wave speed provided by this model is often higher than the Alfvén speed and acoustic speed inside the loop. In contrast, the short period oscillation in our simulation is clearly not a global fast standing mode. Turbulence changes the plasma configuration and makes it possible to trap waves in partial loop regions, and this occurred in between the vortical structures set up by the KHI. If we assume that the wave speed equals the acoustic speed, the length of the resonant cavity is about 500 km s −1 ×26 s/2≈6.5 Mm, which is close to the size of the region where we find a high correlation coefficient. As the local acoustic speed is much higher than the local Alfvén speed, this locally trapped oscillation in our simulation is a locally standing fast wave. It is to be noted that this is distinctly different from the usual global standing mode interpretations given to reported quasi-periodic pulsation (QPPs). However, the 2.5D nature of our simulation may influence this finding, as the trapping vortex regions are extended rolls in the invariant perpendicular direction, and this will change in 3D.
We note that 60 km s −1 is not the averaged Alfvén speed over the entire flare loop, but in those local regions where the magnetic field is locally weakened by the KHI turbulence in the far nonlinear stages. In the turbulent apex, some regions have Alfvén speed lower than 100 km s −1 , while others have Alfvén speeds of 600-700 km s −1 . In the loop limbs, where we do not have turbulence, the Alfvén speed is actually as high as 700-800 km s −1 . The temporal evolution of average beta in the red box in Figure 4(a) is shown in Figure 4(c). This panel clearly demonstrates that the local beta increases rapidly when the evaporation flows arrive and then increases gradually due to the turbulent KHI process. The plasma density and the temperature inside the loop in the flare phase are ∼10 10 cm −3 and ∼20 MK, respectively. The magnetic field strength ranges from 20 to 50 G in the nonturbulent regions. In contrast, the density, temperature, and magnetic field strength of a flare loop reported in Nisticò et al. (2017) are ∼10 10 cm −3 , ∼10 MK and 50 G, respectively.

EUV and Spectral Line Broadening
The EUV emission of the flare loop is also studied. Here we synthesize the EUV Fe XXI emission (1354.08 Å) based on the plasma parameters. The local intensity is calculated with where N e is local electron number density and G (T) is the contribution function of this line from the CHIANTI database. The local Fe XXI line spectra are assumed to have a Gaussian profile, where the line width equals the local thermal velocity and the blueshift of the line equals the local upward speed (lines of sight are taken in the negative y-direction). The local spectra are integrated along lines of sight to get spectra as observed. The brightness of the Fe XXI image in the apex clearly shows the turbulent imprint ( Figure 5(a)). The emission is weaker in the regions with pronounced vortical structures (where, on the other hand, the SXR emission is stronger). The signal of the turbulence can also be found in the spectral line profile of EUV lines. Figure 5(c) displays the line profile of Fe XXI for the apex. This profile is broadened by the turbulent motion of the apex plasma. The profile of Fe XXI for the left footpoint ( Figure 5(b)) also demonstrates a broadening, but this is caused by the (mostly upward) speed difference between plasma moving along different field lines rather than by turbulence.

Neupert Effect and HXR Spectra
In addition to turbulence in the flare loop, some interesting aspects of the X-ray emission are as follows. The first phenomenon is the Neupert effect, which predicts a temporal relationship between HXR flux and SXR flux (Neupert 1968;Hudson 1991;Veronig et al. 2005). According to this effect, SXR flux is closely related to the total thermal energy of the flare loop and the HXR flux is closely related to the input rate of the thermal energy. Therefore, the temporal profile of the HXR flux should fit the time derivative of the temporal profile of SXR. There is a trend for the HXR flux to mimic the variation of SXR flux in our simulation (Figure 3), although the SXR flux increases up to t≈3 minutes, while the HXR flux decreases from t≈1.5 minutes.
The second phenomenon is the difference between the spectrum of the HXR looptop source and that of the footpoint sources. Figure 3(b) shows the spectra of the looptop source and that of the footpoint sources. The looptop spectrum has a power-law distribution I(ò)∼ò − γ with γ=δ l +1 for photon energy below 50 keV, which fits the prediction about the relationship between bremsstrahlung photon spectrum and fast electron spectrum provided by Holman et al. (2011). The footpoint HXR spectrum is flatter than the looptop HXR source, as the spectrum of fast electrons becomes flatter and flatter because of collisions when the electrons penetrate downward into the thick chromosphere. The spectral index of the footpoint spectrum below 50 keV is located between δ l and δ l −1, which is a little greater than the value γ=δ l −1 expected by Holman et al. (2011).

Conclusion and Discussion
We simulate the evaporation of chromospheric plasma, the trigger of KHI, and the production of turbulence in a flare loop. The evaporation flow is produced by footpoint heating and the heating rate is calculated based on the penetration and collisions of fast electrons with chromospheric matter. The EUV, SXR, and HXR emission of the flare loop with turbulence inside the apex are studied. Our main conclusions are as follows.
1. The KHI turbulence leads to fluctuations of emission intensity in the images of EUV, SXR, and HXR. SXR light curves show a clear QPP. 2. The KHI turbulence modify the configuration of the apex magnetic field and leads to generation of short period local standing magnetoacoustic waves in cavities, causing the QPP. This is very different from the usual global modes invoked in other interpretations. 3. The spectral profile of the Fe XXI line emitted by the apex plasma is broadened by turbulent motion of the plasma. 4. The temporal profile of the footpoint (and the total) HXR flux closely mimics the profile of the fast electron flux. The SXR flux tends to increase when the footpoint (or total) HXR flux is high and tends to decrease when the footpoint (or total) HXR flux is low. 5. The HXR spectrum of the thick chromospheric source (footpoint source) is flatter than the spectrum of the input fast electrons, while the HXR spectrum of the optically thin coronal source is steeper than the input electron spectrum.
No obvious oscillation is found in the flux of HXR emission. It seems that the HXR flux is not sensitive to the turbulent fluctuations of plasma density. The reported quasi-periodic oscillations in HXR flux (e.g., Ofman & Sui 2006) must therefore be caused by the oscillation of the input fast electron flux. We need to extend our flare loop scenario to include the overlying reconnection site, and can then investigate whether, in combination with turbulence, this can lead to a periodic acceleration of fast electrons. Although turbulence has been simulated in our work, the change of the input fast electron flux/spectrum has not been considered. In addition to periodic particle acceleration due to oscillating electric fields, reconnection could lead to QPPs via other processes. Takasao & Shibata (2016) suggested that the interaction between the reconnection outflow and the magnetic field above the looptop could produce a magnetic tuning fork, two arms of which oscillate periodically. The oscillations in the magnetic tuning fork cause oscillations of the shocks above the looptop, periodically modifying the efficiency of particle acceleration and thus leading to QPPs in nonthermal emission. Our work has as its main limitation the non-self-consistent calculation of the fast electron energy deposition. The input fast electron flux is from an empirical approach (Equations (1) and (2)) rather than a simulation result of particle acceleration. A self-consistent calculation of the fast electron flux can actually be done inside a pure MHD simulation (e.g., Bakke et al. 2018), but then the reconnection point above the flare loop is required. We note that fast electron energy injection is not the only mechanism to transport energy released in magnetic reconnection from corona to chromospheric layers: thermal conduction also can transport sufficient energy to the chromosphere and trigger high-speed evaporation flows. Evaporation flows have been successfully triggered by energy flux transported via thermal conduction in MHD simulations of solar flares reported by Yokoyama & Shibata (1998, 2001, Takasao et al. (2015), Takasao & Shibata (2016), and Cheung et al. (2019). The energy flux injected into the chromosphere owing to thermal conduction in Cheung et al. (2019) reaches a value of 3×10 11 erg cm −2 s −1 , which is higher than the energy flux carried by the fast electrons in our simulation. Such a high injected energy flux also leads to generation of fast evaporation flows of several hundreds km s −1 in their study. Evaporation flows with such a high speed may similarly be subject to KHI and produce looptop turbulence.
Our calculation of SXR and EUV emission is also empirical. MHD models can reproduce the SXR and the EUV well, as these channels are determined by the thermal motion of particles, and examples are found in Xia et al. (2017) and Zhou et al. (2018) to reproduce EUV emission based on 3D simulations of prominences. Finally, our 2.5D simulation neglects the change of plasma parameters in the out-of-plane direction and does not contain information on the depth of the loop in this direction. We artificially assume that the loop has a constant depth in this direction in reproducing the SXR and EUV images. This limitation is unavoidable in 2D simulations and will disappear in our future 3D study.