Ion–neutral Interactions and Nonequilibrium Ionization in the Solar Chromosphere

When assuming LTE, both the ionization fraction (F_ion) and ambipolar diffusion coefficient (η_A) for a given magnetic field strength are uniquely defined by the temperature and density. This is because the ionization degree in LTE is set by Saha–Boltzmann statistics only. This is shown in panel (A) of Figure 3 for F_ion, and in panel (E) for η_A for a magnetic field strength of 30 G. However, in NEI, η_A, calculated from gol_nei, can vary by many orders of magnitude for a given temperature and density (panels (F)–(H), up to ∼10 orders of magnitude), because in this case the ionization fraction also depends on the history of the plasma (panels (B)–(D)).

Zoom In Zoom Out Reset image size

The ionization fraction and ambipolar diffusion have very different values when assuming LTE or NEI. In LTE (calculated from the LTE equation of state):

• 1.  
  F_ion is small (highly neutral) below $\mathrm{log}T=[3.7-4]$ K, and large (highly ionized) above these temperatures (see panel (A) in Figure 3).
• 2.  
  η_A becomes extremely large (up to ∼10¹⁶ m² s⁻¹) in the cold, low-density areas in the chromosphere, and negligible at high temperatures in upper chromosphere (T > 10⁴ K) (see panel (E) in Figure 3 and Martínez-Sykora et al. 2017a).

In NEI (calculated from gol_nei):

• 1.  
  The median of the ionization fraction (panel (B)) is close to fully ionized, and the ambipolar diffusion coefficient (panel (F)) is at least five orders of magnitude smaller in regions with low densities (ρ < 10⁻⁸ kg m⁻³) and low temperatures (T < 10⁴ K) than what is obtained using LTE.
• 2.  
  The ambipolar diffusion coefficient is much larger at high temperatures in the upper chromosphere (T > 10⁴ K) (panels (F) and (H)) as compared to the LTE case, because plasma, in NEI, is not fully ionized in these regions.

Consequently, we conclude that, in order to properly calculate ambipolar diffusion in the solar atmosphere, a proper treatment of the ionization fraction is required. As mentioned above, as a zero-order approximation, the ionization fraction and ambipolar diffusion coefficient are related by Equation (7). One can consider the validity of this relation by comparing the left and right column of Figure 3. However, this relation is imperfect because the ambipolar diffusion coefficient also depends on the number density and the temperature.

One method for estimating the importance of ambipolar diffusion is by dividing it with the Alfvén speed (u_a): L_amb = η_A/u_a. The quantity L_amb is the length scale on which ambipolar diffusion acts. Alfvén waves with wavelength similar to L_amb, or currents perpendicular to the magnetic field (J_⊥) with characteristic spatial scales equal or smaller than L_amb will be efficiently diffused. Figure 4 shows the median and maximum values of L_amb in simulation gol_nei. Ambipolar diffusion can dissipate J_⊥ on length scales as large as some tens of kilometers (up to a megameter) in extended regions in chromospheric and TR plasma. The median of L_amb is above 10 km (10–10³ km) for densities 10⁻⁹ < ρ < 10⁻⁷ kg m⁻³ and low temperatures (T ≤ 3 × 10³ K) as well as at low densities (ρ < 10⁻¹⁰ kg m⁻³) and temperatures T ≤ 10⁴ K. In addition, some plasma reaches large ambipolar diffusivities at densities ρ < 4 × 10⁻¹⁰ kg m⁻³ and up to TR temperatures (see pink–white colors in the bottom panel of Figure 4). Note that the spatial grid spacing is ∼10 km in this simulation, so anything smaller will not be resolved in the numerical simulation. Most of the currents in the simulation are highly confined, i.e., a few tens of kilometers (e.g., Hansteen et al. 2015; Nóbrega-Siverio et al. 2016; Martínez-Sykora et al. 2018).

Zoom In Zoom Out Reset image size

Figure 5 shows Joint Probability Distribution Functions (JPDFs) of the density and temperature in each simulation, computed during a period of 12 minutes of solar time. Below, we discuss the statistical differences in the thermal properties between the two simulations.

• 1.  
  Simulation gol_lte (top panel of Figure 5) shows three horizontal pink-white stripes at log T = 10^3.8, 10⁴, and 10^4.3. These temperatures mark the transition from H i to H ii, He i to He ii, and He ii to He iii with an LTE ionization balance. A large energy gain or loss is needed for gas to pass these temperatures, and gas temperatures consequently tend to cluster there. In the gol_nei simulation, these stripes are not present, owing to the slow ionization/recombination rates. This process is explained in detail in Leenaarts et al. (2011) and Golding et al. (2016). In gol_nei, the low-density upper chromosphere has a large range of temperatures. Note that, for instance, around ρ ≈ 10⁻¹⁰ kg m⁻³, the pink-white patch extends over an order of magnitude in temperature, due to the presence of NEI, and this is where ambipolar diffusion can play a big role. The plasma in this temperature/density range corresponds mainly to the type II spicules and low-lying loops described in Sections 3.3.2–3.3.3.
• 2.  
  Martínez-Sykora et al. (2017a) showed that the gas located in the wake of shocks (cold chromospheric bubbles) is efficiently heated in the gol_lte model. The larger temperatures in the wake of the shocks in the gol_lte simulation (top panel of Figure 5) result from heating through large ambipolar diffusion. The bubbles are cold and have low densities owing to the expansion. This plasma state leads to a large ambipolar diffusion (Equation (6)), and the ambipolar diffusion dissipates and/or advects almost all the currents that are perpendicular to the magnetic field out of the cold bubbles.When hydrogen and helium are treated in NEI (gol_nei, bottom panel of Figure 5), extended regions of the cold chromospheric bubbles reach the minimum temperature set by the ad hoc heating term. Neither recombination nor the formation of H₂ molecules occurs fast enough to counter the cooling caused by adiabatic expansion. Ambipolar diffusion is inefficient, in contrast to the gol_lte model: the ionization degree remains on the order of 10⁻² in NEI, despite the low temperatures, and the ambipolar diffusion is therefore much weaker than in LTE, where the ionization degree can be as low as 10⁻¹⁰. The gol_nei simulation does not reach even lower temperatures because of the ad hoc heating term in the model preventing temperatures from dropping below ∼1800 K.
• 3.  
  The gol_nei simulation has plasma with $-10\lt \mathrm{log}\rho \,\lt -6$, with temperatures higher than reached in the gol_lte model. This is indicated with the yellow contour labeled 2 in the bottom panel of Figure 5. We analyze this population in Section 3.3.2.
• 4.  
  The TR is less sharp in the gol_nei simulation than in the gol_lte. The JPDF shows slightly greater values within the TR in the gol_nei simulation (bottom panel). In addition, the density variation is largest in the gol_nei simulation at low TR temperatures ($4.0\leqslant \mathrm{log}(T)\leqslant 5.0$). These thermal properties are due to the thermal evolution of type II spicules and low-lying loops resulting from NEI and ambipolar diffusion. See the detailed description of type II spicules and low-lying loops in Sections 3.3.2 and 3.3.3.
• 5.  
  Finally, in gol_lte, a few dense (ρ > 10⁻⁶ kg m⁻³) chromospheric plasma elements (JPDF ≈ 10⁻⁷) reach lower temperatures than in the gol_nei simulation. The dense chromosphere with low temperatures in the gol_lte simulation corresponds to a few regions that have accumulated enough magnetic flux in the photosphere to emerge. In such locations, the ambipolar diffusion allows the magnetic flux to expand, carrying dense plasma and producing rather large, cold bubbles (Martínez-Sykora et al. 2008; Tortosa-Andreu & Moreno-Insertis 2009; Ortiz et al. 2014). We do not see enough of such regions of accumulated magnetic flux at the photospheric heights in the gol_nei simulation to produce similar events.

Zoom In Zoom Out Reset image size

We find that most of the fluid is highly ionized in the cold plasma forming during the early phases of the rise of cold bubbles, as well as in thin cold regions within type II spicules and low-lying loops (green, blue, and red arrows in panels (B) and (G) in Figure 2).

In the low to middle chromosphere, the fraction of neutral hydrogen (H i, panel (A) in Figure 6) decreases at a height of z ≈ 1 Mm, which also happens to always be at roughly at the same density of 10⁻⁵ kg m⁻³ (see panel (F) in Figure 2). This is the location where the acoustic waves driven by convective motions at the photosphere become shocks due to the exponential density drop. Shocks appear (in blue in panel (C) of Figure 7) as a result of strong compression and behind the shock fronts, i.e., in the wake, cold bubbles form as the plasma expands (red). Shocks passing through the middle chromosphere ionize hydrogen (panel (A) in Figure 6). When nonequilibrium ionization is accounted for, shocks passing through the chromosphere ionize most of the hydrogen rapidly while the recombination timescales (∼10⁴ s) are much longer than the hydrodynamic timescales (10–300 s). Consequently, the plasma stays ionized despite the drop in temperature in the cold bubbles (red in panel (E) in Figure 6). As a result of this, ${\rho }_{a\hat{I}}{\nu }_{a\hat{I}a^{\prime} 0}$ (see Equation (4)) is high and the coefficient of ambipolar diffusion is extremely low in NEI conditions (panel (I) in Figure 2). Therefore, the effects of ambipolar diffusion in LTE are highly overestimated in the cold bubbles (blue in panel (B) of Figure 7). It is only toward the end of the cold bubble evolution (a few minutes after the shock passes through) that the plasma starts to recombine and the effects of ambipolar diffusion increase in significance.

In NEI, regions with strong expansion in the upper chromosphere and TR appear as red in panel (C) of Figure 7. These regions are cold inside type II spicules and in low-lying loops, but are relatively ionized in hydrogen, while being highly neutral in helium (panels (B) and (C) in Figure 6). As in the middle chromosphere, the long recombination timescales as compared to the timescales of plasma expansion lead to a large percentage of ionized hydrogen compared to what would be found in an atmosphere computed with the assumption of LTE (red in panel (E) in Figure 6). This leads to smaller effects of ambipolar diffusion in these strongly expanding regions than in the LTE case (blue regions in panel (B) of Figure 7).

These expanding regions (cold expanding bubbles and cold thin regions near the TR) are areas termed population 1 in the bottom panel of Figure 5, which will be described in further detail in Section 3.3.1. Statistically speaking, this shows that the medians of the coefficients of ambipolar diffusion η_A (panel (F) of Figure 3) and L_amb (top panel of Figure 4) are small for the plasma contained within the contour labeled 1 (η_A < 10⁸ m² s⁻¹).

The opposite process happens in magneto-acoustic shocks or in the steep temperature rise of the TR (green and yellow contours in Figure 7). The number of neutrals in these regions is larger than expected under the assumption of LTE (blue in panel (E) of Figure 6). This is true for both hydrogen (panel (A) of Figure 6) in chromospheric shocks and for helium (panel (B)) in the upper chromosphere and TR. In these locations, the time between shocks passing is long enough that a large amount of plasma has managed to recombine. In addition, ionization does not happen immediately, and furthermore, the timescales of ambipolar heating are slightly shorter than the ionization timescales. All of these effects added together leads to lower ionization fractions in these regions in NEI than in LTE (blue in panel (E) in Figure 6). Thus, these locations have a larger ambipolar diffusion in NEI than in LTE (red in panel (B) of Figure 7).

The regions with the largest effects of ambipolar diffusion are found in low-lying TR or chromospheric magnetic loops at times long after shocks have passed through (e.g., around x ≈ [54−60] Mm in Figures 6 and 7). In these loops, the density decreases due to the draining of plasma along the loops toward the photosphere. Magneto-acoustic shocks passing through these highly inclined loops are less frequent than in other, more vertical magnetic structures, because the change in effective gravity leads to a lower acoustic cutoff frequency (e.g., Heggland et al. 2007). Consequently, recombination timescales are short compared to the time that elapses between shocks propagating along these loops. In this case, neutrals become important, and due to the low density, the ion–neutral collision frequency is low, resulting in an extremely high ambipolar diffusion. Additionally, these loops have more neutrals, especially in helium (panel (B) in Figure 6) than they would have assuming LTE (blue in panel (E)). The overabundance of neutrals is probably not due to the timescales involved, but rather because under these conditions, the ionization state is better modeled by "coronal equilibrium" conditions, with collisional ionization being balanced by radiative recombination, which predicts more neutrals at a given temperature than what is found in LTE. This state of affairs is detailed for hydrogen in the top panels of Figure 4 of Rutten (2016). Helium behaves in the same manner. A consequence is that the heating timescales due to ambipolar diffusion are much shorter than the ionization timescales. The ionization fraction, mostly in helium, in these regions can be severely overestimated when assuming LTE (see panels (B), (C), and (E) in Figure 6, and Section 3.3.3 for further details), and the ambipolar diffusion is larger in NEI than assuming LTE (red in panel (B) of Figure 7).

In short, the effects of ambipolar diffusion are affected by the choice of assumption for the ionization balance. Assuming one scenario or the other can provide opposite roles for the ambipolar diffusion: in those locations where ambipolar diffusion is large in LTE (cold bubble or inside expanding spicules), it seems unimportant under conditions of NEI; on the other hand, in those locations where the ambipolar diffusion is small in LTE (e.g., shock fronts), with NEI the ambipolar diffusion is large enough to play an important role (as detailed below).

Let us consider the physical processes that play a role in the cold regions that lie along type II spicules or in the cold bubbles formed by rarefaction behind magneto-acoustic shocks. These regions in the gol_nei simulation are visible within the black contours in panels (B)–(D) of Figure 7. In the LTE models described by Martínez-Sykora et al. (2017a), cold bubbles and spicules are heated via ambipolar heating, which is the largest contributor. In those models, the ambipolar diffusion increases drastically, reaching values that rapidly dissipate any current or Aflvénic wave passing through. The ambipolar diffusion becomes extremely large because the ion–neutral collision frequency and ionization degree diminish (panel (D) in Figure 2, and panel (E) in Figure 3) as plasma expansion causes the temperature and density to fall rapidly. Consequently, cold regions behind the shocks are vigorously heated by the ambipolar heating (panel (E) in Figure 2). However, this situation changes once NEI is taken into account (right column in Figure 2). At that point, we find the cold bubbles to be much cooler; see Leenaarts et al. (2011), but cf. our discussion for comparison with the literature. Similarly, elongated thin structures inside type II spicules can be much cooler in the gol_nei simulation than in gol_lte (compare panels (B) and (G) in Figure 2).

Going back to the JPDF maps in Figure 5, we mentioned a very intriguing population in the gol_nei simulation indicated with a solid blue contour labeled 1 in the bottom panel. We added white or black contours in the maps in Figure 7 that correspond to the plasma within this population. There, one can see that this population corresponds to the cold expanding bubbles (see animated Figure 7) as well as regions with strong expansion inside type II spicules. The slope of the contour labeled 1 in the bottom panel of Figure 5 follows ρ ∝ T^5/3. This corresponds to the adiabatic relation between density and temperature. Statistically speaking, in these regions, the plasma behaves adiabatically. This has two causes: (1) ambipolar diffusion is not large enough to dissipate any current in the NEI case as opposed to in LTE, so the entropy sources are minor (see below); and (2) recombination timescales are larger than expansion timescales, at least for the early stages of the lifetime of the bubbles, so recombination does not significantly change the plasma state.

In fact, in these regions, the timescales of the energy losses due to expansion can be as low as $\tfrac{1}{{\rm{\nabla }}\,\cdot \,{\boldsymbol{u}}}\approx 10\,{\rm{s}}$, which is much smaller than recombination timescales (∼10⁴ s). In other words, cooling from expanding regions in simulations assuming LTE has a "reservoir" energy in the recombination energy of the ionized species. In contrast, in the simulation with nonequilibrium ionization/recombination (gol_nei), the cooling due to expansion is not using the "reservoir" energy in recombining the ionized species, due to the rather large recombination timescales. Instead, the temperature drops as dictated by the adiabatic expansion.

As mentioned before, in LTE (simulation gol_lte), the ambipolar diffusion increases abruptly with expansion (panel (D) in Figure 2 and associated animation) because of a decrease in ion–neutral collision frequency and ionization fraction (panel (C)). There, the ambipolar diffusion dissipates, or advects, almost any current or Alfvénic wave traveling through the cold bubbles (panels (D)–(F)).

In contrast to this, when modeling NEI (simulation gol_nei), ambipolar diffusion is nonexistent in these regions (Section 3.2.1), and thus the ambipolar heating there is negligible (panel (J) in Figure 2). Figure 8 shows the histogram of the ambipolar diffusion (top) and of the ambipolar heating (bottom) for the gol_nei simulation within the three selected regions: 1 (associated with the adiabatic expanding plasma; solid line), 2 (associated with the shocks; dashed line) and 3 (associated with low-density plasma; dotted line) shown in the bottom panel of Figure 5. Figure 8 confirms that the ambipolar diffusion is small and the ambipolar heating is negligible in the expanding cold bubbles and cold regions within the spicules (solid lines).

Zoom In Zoom Out Reset image size

In simulation gol_nei, the Joule heating in cold bubbles is basically negligible, not only due to the rather low ambipolar diffusion, but also because most of the current and shocks are ahead of the cold bubble (panels (C)–(D) of Figure 7). In spicules, the region that corresponds to the population 1 (white contours in panel (A) of Figure 7) is highly confined along the spicule. Similarly, in these thin regions, the Joule heating is lower than in other locations within the spicule, since fewer and smaller currents are present.

The second population of interest in this work has the opposite behavior. This population is indicated with the dashed yellow contour labeled 2 in the bottom panel of Figure 5 for simulation gol_nei. The plasma comprising this population is localized mostly in the shock fronts or boundaries and/or footpoints of type II spicules (yellow contours in Figure 7 and corresponding animation).

As mentioned above, in magneto-acoustic shocks and at sharp transitions with enough ambipolar diffusion and large currents, the ambipolar diffusion allows the dissipation of magnetic energy into heat. Therefore, ambipolar heating in the gol_nei simulation contributes to heating these locations (panel (J) in Figure 2 and corresponding animation). In fact, the ambipolar heating (dashed line in bottom panel of Figure 8) is up to an order of magnitude larger than inside the expanding cold bubbles (solid). The heat in shocks and those locations in spicules with large ambipolar diffusion comes not only from ambipolar heating, but also artificial diffusion, viscous heating, and compression. The timescales of all these heating processes can reach a similar order of magnitude (∼1 minute).

It is interesting to consider what this rapid heating in spicules means for a long-standing open question about the rapid apparent motions associated with spicules visible in TR lines. There are conflicting observations between type II spicules (typically observed spectrally with lines shifts showing line-of-sight velocities of up to ∼100 km s⁻¹), and apparent fast motions (the so-called network jets that show apparent speeds in the plane of the sky of a few hundred of km s⁻¹) (Beck et al. 2013; Tian et al. 2014; Rouppe van der Voort et al. 2015). De Pontieu et al. (2017) have shown that network "jets" do not necessarily involve very high plasma flows but can be caused by rapidly propagating heating fronts associated with ambipolar dissipation of currents traveling along type II spicules. The quoted work assumed LTE. In the current work, we relax the LTE assumption and investigate NEI conditions. The NEI calculations do not fundamentally change this picture. We again find that the extremely large currents in spicules (see bottom panel of Figure 7) are dissipated via ambipolar diffusion. The extremely large currents explain how this population (dashed line versus solid line in Figure 8) can be subject to much larger ambipolar heating than what is found in cold expanding bubbles, despite the relatively smaller ambipolar diffusion in spicules. For instance, the example shown in panel (J) in Figure 2 reveals large heating occurring along the boundaries of (and inside of) the spicules. In other words, the NEI effects mostly affect the spatial distribution and (local) amount of heating, but do not not impact how fast a propagating current or Alfvén wave is propagating or dissipated. As in shock fronts, artificial diffusion, viscous heating, and compression also play an important role. In addition, spicules have temperature gradients high enough that thermal conduction contributes substantially to their thermal properties. The importance of each of these source terms varies within the spicules and evolution.

Due to the relatively large ionization timescales (a few tens of seconds) in these regions of the gol_nei simulation, heating goes directly to increasing the plasma temperature instead of ionizing the plasma. Consequently, magneto-acoustic shocks are much hotter in the gol_nei simulation. In spicules, we find both: (1) cold regions with low ambipolar diffusion, low currents, large expansion, and large recombination timescales; and (2) hot regions with high ambipolar diffusion and large currents. These regions experience strong a temperature increase due to long ionization timescales. As a result of these two properties, spicules in gol_nei have large temperature variations within them. Note that, roughly when the spicule has fully formed, at its boundaries, the ambipolar diffusion can get as high as 10¹⁰ m² s⁻¹ (panel (A) in animated Figure 7). The large density and temperature variations within spicules may explain the apparently conflicting observations of temperatures and speeds in spicules (compare Beck et al. 2013; Tian et al. 2014; Rouppe van der Voort et al. 2015).

Another important aspect is that type II spicules seem to impact the associated coronal loops (Henriques et al. 2016; De Pontieu et al. 2017). Due to NEI effects, the amount of current (which, in our models, causes heating of coronal plasma) that pierces into the corona is larger than in the LTE case, and its role in the corona seems to be larger than what was previously reported by Martínez-Sykora et al. (2018). We find that the corona has on average 10% more current density in the gol_nei simulation than in the gol_lte simulation. A proper calculation of synthetic observables from this simulation and comparison with observations is required for further studies on the apparent motions and impact of spicules on the corona.

As mentioned above, in NEI, low-lying loops with either chromospheric or TR temperatures show the largest ambipolar diffusion (η_A > 10⁸ m² s⁻¹, Section 3.2.2). These are regions of relatively low density (10⁻¹² < ρ < 10⁻⁹ kg m⁻³), with temperatures T < 10⁵ K (labeled as region 3 and shown with dashed–dotted green contours in Figures 5 and 7, and with dotted lines in Figure 8).

Since the ambipolar diffusion is larger in NEI than assuming LTE in low-lying loops (Section 3.2.2), the heating is also larger in NEI than in LTE. In NEI, heating from the ambipolar diffusion is more than an order of magnitude larger in these loops than in the other features described in the previous sections (dotted line in bottom panel of Figure 8). This is because the ambipolar diffusion coefficient is several orders of magnitude larger in these loops (dotted line in top panel of Figure 8) than in the other features (dashed or solid lines). Consequently, almost any Alfvénic wave or current is dissipated there.

These low-lying loops share similarities with two different types of loops that are observed in the solar atmosphere: (1) loops that appear to outline the canopy and connect plage regions of opposite polarity in active regions, typically clearly visible in chromospheric (e.g., Hα) or TR (e.g., Si iv 1400 Å) diagnostics; (2) so-called Unresolved Fine Structure (UFS, Feldman 1983; Hansteen et al. 2014).

Our results suggest that "canopy loops" may be composed of plasma with a large range of temperatures due to NEI effects, "adiabatic" expansions, and ambipolar heating, with a mix of hot and cool strands right next to each other. In our simulations, these low-lying loops initially are filled with cool plasma due to shocks passing through from one footpoint all the way to the other footpoint (see animated Figures 2, 6, and 7). After the cold loop has formed, the fluid recombines, ambipolar diffusion increases, and the loop is heated up to TR temperatures via current dissipation from the ambipolar diffusion. This could explain the presence of chromospheric and TR loop-like structures that both appear to outline the canopy and occur in close proximity of one another.

Similarly, observations suggest that UFS loops seem to be formed by sporadic heating instead of cooling from hot loops (Hansteen et al. 2014; Brooks et al. 2016). Three-dimensional radiative MHD models (excluding ambipolar diffusion) show that the cycling of TR plasma starts with being first rapidly heated from chromospheric temperatures, to later drain and cool down in small magnetic loops (Guerreiro et al. 2013). Our model suggests that the ambipolar diffusion may play an important role in the thermal evolution of these low-lying loops: the ambipolar diffusion dissipates Alfvénic waves or currents into thermal energy. Since the ambipolar heating timescales are shorter or comparable to ionization timescales, most of the heating goes into heating the plasma instead of ionizing the plasma. In addition, the ambipolar diffusion is important all the way up to TR temperatures, even above helium ionization temperatures in LTE, owing to the large role of NEI as shown in Figures 2, 3, and 6.