Chromospheric Bubbles in Solar Flares

Solar flares are impulsive events with energies reaching up to 10³² erg, caused by the sudden release of free magnetic energy in the corona that is transported into the chromosphere during the impulsive phase of the flare event. The prevailing theory is that this energy transport arises from a propagating beam of accelerated particles driven from the corona (e.g., Brown 1971; Emslie 1978; Holman et al. 2011). Alfvén waves have also been suggested as an alternative energy transport mechanism (Fletcher & Hudson 2008; Kerr et al. 2016; Reep et al. 2018), but this process has not been explored as thoroughly. This Letter focuses on the primary theory of particle beam heating, specifically electron beams.

The majority of nonthermal accelerated particles can penetrate into the dense chromosphere, causing intense heating via Coulomb collisions, resulting in X-ray emission via Bremsstrahlung at the loop foot-points (Hudson et al. 1992; Neidig & Kane 1993; Martínez Oliveros et al. 2012). This rapid heating causes a pressure gradient that drives a flow upward, pushing material to greater geometrical heights in a process known as chromospheric evaporation (Neupert 1968). There are two main types of chromospheric evaporation (Milligan et al. 2006). Gentle evaporation occurs with lower beam fluxes when thermal heating causes chromospheric expansion via an upflow of the order of 10 km s⁻¹, with little evidence of a corresponding downflow. Explosive evaporation occurs when the beam flux heats the chromosphere to coronal temperatures, causing the transition region to initially shift to lower geometrical heights. The intense heating cannot be radiated away sufficiently fast and this results in the expansion of the chromosphere with velocities up to 100 km s⁻¹. The velocity is strong enough to create a shock, and due to the conservation of momentum, a strong downflow can also be present, known as chromospheric condensation (Kosovichev 1986; Hudson 2011). Fisher et al. (1985) suggest that the flux threshold between gentle and explosive evaporation for a 20 keV low-energy cutoff model is ∼10¹⁰ erg cm⁻² s⁻¹. The crossover flux has also been estimated as 2–8 × 10⁹ erg cm⁻² s⁻¹ using both observations and the F-CHROMA grid of numerical models (Sadykov et al. 2019) and has also been estimated with other models (Reep et al. 2015) and observations (Gomory et al. 2016). It has also been proposed that the photosphere can be heated directly by electron beams or even via proton beams that can penetrate deeper into the atmosphere (Švestka 1970; Machado et al. 1978; Aboudarham & Henoux 1986; Procházka et al. 2019).

The parameters of the nonthermal particle beams can be constrained using the X-ray spectrum captured with instruments such as the Reuven Ramaty High Energy Spectroscopic Imager (Lin et al. 2002) or the Fermi Gamma-ray Space Telescope (FERMI; Meegan et al. 2009) during the impulsive phase of a flare. Estimates for the cutoff energy, spectral index, and flux of the accelerated electron beam can be inferred from this X-ray spectrum using the Collisional Thick Target Model, as demonstrated by Petrosian & Chen (2010).

As the majority of the accelerated electrons lose their energy in the chromosphere, this is the part of the solar atmosphere where the majority of the flare radiative output originates (Fletcher et al. 2011). It is therefore imperative to understand how the chromospheric plasma reacts to the dynamics of the magnetic reconnection that mediates the resultant influx of precipitating particles. In this Letter, we describe the new phenomenon of chromospheric bubbles that are related to the rate of chromospheric evaporation and particle deposition in the lower solar atmosphere. Section 2 describes the simulations used, while Section 3 characterizes the bubbles. Section 4 discusses the overall findings and potential future work.

The work presented in this Letter uses a grid of flare models created with the RADYN code (Carlsson & Stein 1992, 1994, 1995) and are part of the F-CHROMA model archive.⁴ RADYN solves the equations of radiative hydrodynamics and equation of state along a one-dimensional atmosphere as described in detail by Allred et al. (2006). It allows for the introduction of energy via various scenarios, such as direct thermal heating (Procházka et al. 2017; Reid et al. 2017), Alfvén wave heating (Kerr et al. 2016), and also via particle beams (Rubio da Costa & Kleint 2017; Procházka et al. 2019) with the modifications created by Abbett & Hawley (1999; see also Allred et al. 2005). The RADYN model atmosphere contains a 6-level hydrogen atom, a 9-level helium atom, and a 6-level calcium atom. It solves the equation of radiative transfer in complete redistribution, which is a good approximation only for nonresonance lines.

The F-CHROMA flare models utilize a Fokker–Planck-type beam (Allred et al. 2015; Daou & Alexander 2016), and use an initial atmosphere with 300 grid points along a VAL3C starting atmosphere. An adaptive grid allows for small-scale dynamic events to be accurately resolved. Allred et al. (2015) compared how the hardness of the applied beam and the low-energy cutoff affects the resultant location of the deposition of the nonthermal electrons for both a Fokker–Planck-type beam and an Emslie beam (Emslie 1978). The F-CHROMA grid also contain the physics of return currents (Holman 2012), which considers how the accelerated electrons produce an electric field that drives a counter-propagating, neutralizing return current. The return current will heat the plasma via Joule heating (van den Oord 1990), but will not largely affect the chromospheric energy deposition.

RADYN simulates a 10 Mm half loop. The spectral index δ ranged between 3 and 8 keV while the low-energy cutoff E_C ranged between 10 and 25 keV in steps of 5 keV. The total beam fluxes chosen were 3.0e + 10 erg cm⁻², 1.00e + 11 erg cm⁻², and 3.0e + 11 erg cm⁻². This results in a total of 72 models. The simulations were run for a total of 50 s, with atmospheric outputs saved every 0.1 s. The beam heating is applied for 20 s with a triangular temporal profile peaking at T = 10 s. The parameters selected for the input beam are listed in Table 1.

Figure 1 shows the deposited electron energy and radiative losses from the electron beams in the simulated grid at a time of T = 10 s. The left panel shows how varying the spectral index has a relatively small effect on the penetration depth of the beam, while varying the low-energy cutoff has a much larger impact. Varying the flux the most significant impact on the temperature increase. If the conditions are just right, the chromosphere will heat up quickly to over 100,000 K. This heating will occur sufficiently deep in the atmosphere to leave a pocket of undisturbed chromosphere sandwiched between the energy deposition and the transition region. As the simulations develop and chromospheric evaporation begins, bubbles of small pockets of chromosphere begin to rise. These bubbles are small regions of chromospheric temperature and increased electron density that also show an increase in mass density relative to the surrounding coronal-type atmosphere. These oddities appear mainly in models with high beam fluxes and softer beams. These small pockets of chromospheric plasma also show an increase in gas pressure at the boundaries at either side of the discontinuity. An example of one of these bubbles can be seen in Figure 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The top panel of Figure 2 clearly shows some "trapped" plasma with chromospheric temperatures within the corona at T = 10–15 s. This bubble is formed when the energy is sufficiently high to cause heating over 100,000 K at chromospheric heights. Importantly, this heating must occur sufficiently deep, and be confined so that it does not heat the initial transition region. We essentially have the creation of a chromospheric pocket of plasma between the transition region and the region of the chromosphere heated by the electron beam. To ensure this is not an artifact of the simulations, the same beam was run with no return current. We have also modified the weighting of the grid points in the velocity domain to check whether models that put less emphasis on small-scale dynamics in the corona could still resolve the bubbles. In all cases, the bubbles created were identical. However, changing the starting atmosphere will change whether a bubble is identified. This is due to the change in the transition region location with respect to the hydrogen column density. The beam will therefore not penetrate at the same location into a different starting atmosphere, which will not trap a portion of the chromosphere. It may very well be that the bubbles can be created in other starting model atmospheres but with a different set of beam parameters.

In the example shown in Figure 2, the initial bubble accelerates upward into the corona initially due to the pressure gradient caused by the beam heating in the chromosphere. However, a shock forms at T = 10 s, which highly accelerates the upper atmosphere, including the bubble. By T = 15 s into the simulation, this velocity is already over 300 km s⁻¹. The value of temperature, density, and electron density are all stationary within the bubble over time, which is essentially the mass motion of the material. Interestingly, due to the gas pressure increasing at the edges of the bubble, the overall width of the bubble decreases over time. For the bubble in Figure 2, 45% of the mass is lost between 10 and 15 s due to redistribution of the bubble mass into the overall wave front. This compression can be most easily seen from T = 5 s to T = 10 s in Figure 2.

Individual well-defined bubbles were identified by investigating the temperature of each atmosphere. The location of the transition region was defined to be the lowest point in the atmosphere where the temperature gradient exceeds 5000 K between grid cells, while having a temperature value of at least 50,000 K. Any point above the transition region that has a temperature below 40,000 K was deemed to be a part of the bubble and logged. This resulted in 21/72 atmospheres with bubbles. Of these 21, all had beam fluxes above 1.0e11 erg cm⁻², and low-energy cutoff values above 15 keV. The most common value of low-energy cutoff being E_c = 15 keV. This value provides the optimal penetration depth to efficiently trap a portion of the chromosphere while accelerating it upward (see Figure 1). The E_c = 10 keV models mainly heated the upper chromosphere and some of the transition region, and did not lead to the creation of bubbles. The E_c = 25 keV models only resulted in bubbles for the highest spectral indices (δ > 5) with the highest levels of beam flux (3.0e11 erg cm⁻²). This is due to the extremely high energy density that needs to be applied to the lower chromosphere in order to heat it sufficiently and create a bubble. Harder beams tend to have a larger spread of energy into the denser portions of the photosphere and lower chromosphere.

The bubble lifetimes can be as short as a few seconds or survive to the end of the 50 s simulations. Once bubbles reach the top of the 10 Mm half loop structure, they are reflected due to the top boundary condition. Any physical changes in the atmosphere post-reflection are not considered as realistic as this would imply symmetrical bubbles being created on either side of the loop. The lifetime of a bubble in the simulation ends when the simulation can no longer resolve it due to the minimum resolution between grid points (∼1 km). All bubbles will exponentially shrink with time due to external forces pushing inward. Once a bubble pops, the simulation then redistributes the grid points, inducing some numerical artifacts that appear to look like waves. This artifact will be damped, and the atmosphere will settle to a new equilibrium. This means that the lifetime of a bubble cannot be accurately determined from these simulations. It remains unclear as to what happens to a bubble once it shrinks below 1 km in width.

The upward movement of the bubbles will always have a roughly constant acceleration in line with the acceleration of the upwardly moving wave front from the chromospheric shock. This can range from 4 to 77 km s⁻², depending on the formation height and point along the wave front the bubble first appears, as well as the amplitude of the wave front itself, as the bubble essentially "rides" along this wave front.

Some beams will result in multiple bubbles. Secondary bubbles will generally form only after the initial bubble has progressed into the corona (usually after 10 s), and while beam heating is still being applied. The secondary bubbles are usually formed lower in the atmosphere, due to the atmospheric changes that have occurred from the previous heating event. As such, they will move into the corona with a smaller acceleration as they are formed in denser media causing less upward propagation. Tertiary bubbles also exist in three of the most extreme simulations, with 3.0e11 erg cm⁻² beam fluxes, δ > 6, and E_c = 15 keV.

In order to estimate the observational signatures of the bubbles, we carried out radiative transfer calculations in Hα and Ca ii 8542 Å. Figure 3 shows the line contribution functions along with the corresponding line profiles (shown in green) for Hα and Ca ii 8542 Å for a 3.0e11 erg cm⁻² beam flux with δ = 7 and E_c = 20 keV. The contribution functions shown are described in detail in Carlsson & Stein (1997).

Zoom In Zoom Out Reset image size

A primary bubble exists near the peak of the velocity, at around 6 Mm. This can be seen in the bottom left panels of Figure 3 where the source function (green line) and Planck function (brown line) are shown. The bubble is optically thin with the τ = 1 location at the corresponding Doppler wavelength being in the lower photosphere, with no signs of signal in the resultant intensity. However, the secondary bubble (at 2.6 Mm) does show significant Hα and some Ca ii 8542 Å signal at a Doppler velocity of 140 km s⁻¹. The chromosphere at this location, combined with a strong velocity gradient, is sufficiently dense to cause a large increase in the optical depth (τ_ν). The sudden change into a much denser medium increases the opacity proportionately, as can be seen at a height of 2.6 Mm in the top left panels of both sets of panels of Figure 3.

In this instance, the bubble appears optically thick in the Hα line profile, showing strong absorption features as it appears like a quiet-Sun piece of chromosphere. The Ca ii 8542 Å spectra, however, show the bubble as optically thin, with only a minor change in the line profile. These bubbles have also been recently shown to be dominant in the Lyman α line of hydrogen by Brown et al. (2018).

We have utilized the publicly available F-CHROMA flare model archive hosted by Queen's University Belfast to study the response of the solar atmosphere to electron beam heating. Certain beam parameters lead to the generation of chromospheric bubbles propagating into the corona. These bubbles are confined regions of chromospheric plasma trapped between the transition region and the lower chromosphere that has been heated by the electron beam. The resultant shock wave front accelerates the bubble upward into the corona. The bubbles will compress over time, but can be accelerated to Doppler velocities of up to 500 km s⁻¹. These bubbles will eventually disappear as they compress, releasing an artificial numerical artifact throughout the atmosphere that we do not consider realistic and solely due to the redistribution of the grid points. The bubbles are only formed if sufficient heating occurs in the chromosphere without affecting the transition region, and so require optimal electron beam parameters to manifest. A small number of 1.0e12 erg cm⁻² flux beam models were also created, with all showing evidence of bubbles. Only those models that show explosive evaporation contain bubbles. The bubbles appear strong in the Hα line up to 10 Å from line core. The Ca ii 8542 Å line profiles also show an optically thin remnant of the bubbles. Not all bubbles appear as optically thick in Hα, and indeed some show no noticeable effects on the calculated line profiles. As such, only a small portion of beams will produce these bubbles, and even fewer can be observationally detected. The observational signatures of the bubbles could be detected with an instrument that covers a broad spectral range around chromospheric lines and provides a good signal-to-noise ratio. To our knowledge, no observations currently exist for solar flares that show observational evidence of the bubble. However, these upflowing bubbles may have been observed previously in a flaring dMe star (Gunn et al. 1994) with similar velocities, only assumed to be high-velocity evaporation instead. The large flux of the flare observed aligns with the theory presented in this Letter that stronger fluxes are more likely to produce these bubbles.

The research leading to these results has received funding from the European Community's Seventh Framework Programme (FP7/2007–2013) under grant agreement No. 606862 (F-CHROMA), and from the Research Council of Norway through the Programme for Supercomputing.