Magnetic Reconnection in Strongly Magnetized Regions of the Low Solar Chromosphere

The HiFi module for partially ionized plasmas includes the electron-impact ionization, radiative recombination, and resonant charge exchange (CX) interactions (Meier 2011; Leake et al. 2012, 2013; Meier & Shumlak 2012). We only consider hydrogen gas in this work. The subscripts "n," "i," and "e" refer to neutrals, ions, and electrons in this work, respectively. The ionization and recombination rates are given by

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{n}^{\mathrm{ion}}\equiv -{n}_{n}{\nu }^{\mathrm{ion}},\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{i}^{\mathrm{rec}}\equiv -{n}_{i}{\nu }^{\mathrm{rec}},\end{eqnarray} \tag{ 2 }$$

with ${{\rm{\Gamma }}}_{i}^{\mathrm{ion}}=-{{\rm{\Gamma }}}_{n}^{\mathrm{ion}}$ and ${{\rm{\Gamma }}}_{i}^{\mathrm{rec}}=-{{\rm{\Gamma }}}_{n}^{\mathrm{rec}}$. The ionization frequency

$$\begin{eqnarray}&&{\nu }^{\mathrm{ion}}=\displaystyle \frac{{n}_{e}A}{X+{\phi }_{\mathrm{ion}}/{T}_{e}^{* }}{\left(\displaystyle \frac{{\phi }_{\mathrm{ion}}}{{T}_{e}^{* }}\right)}^{K}\exp \left(-\displaystyle \frac{{\phi }_{\mathrm{ion}}}{{T}_{e}^{* }}\right)\,{{\rm{m}}}^{3}\,{{\rm{s}}}^{-1}\end{eqnarray} \tag{ 3 }$$

is given by the practical fit from Voronov (1997), using the values $A=2.91\times {10}^{-14}$, K = 0.39, X = 0.232, ${T}_{e}^{* }$ is the electron temperature T_e specified in eV, and the hydrogen ionization potential ${\phi }_{\mathrm{ion}}=13.6\,\mathrm{eV}$. The recombination frequency obtained from Smirnov (2003) is

$$\begin{eqnarray}&&{\nu }^{\mathrm{rec}}=2.6\times {10}^{-19}\displaystyle \frac{{n}_{e}}{\sqrt{{T}_{e}^{* }}}\,{{\rm{m}}}^{3}\,{{\rm{s}}}^{-1}.\end{eqnarray} \tag{ 4 }$$

The CX reaction rate, ${{\rm{\Gamma }}}^{\mathrm{cx}}$, is defined as

$$\begin{eqnarray}&&{{\rm{\Gamma }}}^{\mathrm{cx}}\equiv {\sigma }_{\mathrm{cx}}({V}_{\mathrm{cx}}){n}_{i}{n}_{n}{V}_{\mathrm{cx}},\end{eqnarray} \tag{ 5 }$$

where

$$\begin{eqnarray}&&{V}_{\mathrm{cx}}\equiv \sqrt{\displaystyle \frac{4}{\pi }{V}_{{Ti}}^{2}+\displaystyle \frac{4}{\pi }{V}_{{Tn}}^{2}+{V}_{{in}}^{2}}\end{eqnarray} \tag{ 6 }$$

is the representative speed of the interaction and ${V}_{{in}}^{2}\,\equiv | {{\boldsymbol{V}}}_{i}-{{\boldsymbol{V}}}_{n}{| }^{2}$, with ${{\boldsymbol{V}}}_{\alpha }$ denoting the velocity of species α. The thermal speed of species α is ${V}_{T\alpha }=\sqrt{2{k}_{{\rm{B}}}{T}_{\alpha }/{m}_{\alpha }}$, where ${T}_{\alpha }$ is the temperature, ${m}_{\alpha }$ is the mass of the corresponding particle, and k_B is Boltzmann's constant. The expression for the CX cross section ${\sigma }_{\mathrm{cx}}({V}_{\mathrm{cx}})$ can be found in Meier (2011) and Meier & Shumlak (2012).

The equations solved by the plasma-neutral module of HiFi are the same as those in Murphy & Lukin (2015). The ion and neutral continuity equations are

$$\begin{eqnarray}&&\displaystyle \frac{\partial {n}_{i}}{\partial t}+{\rm{\nabla }}\cdot ({n}_{i}{{\boldsymbol{V}}}_{i})={{\rm{\Gamma }}}_{i}^{\mathrm{rec}}+{{\rm{\Gamma }}}_{i}^{\mathrm{ion}},\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {n}_{n}}{\partial t}+{\rm{\nabla }}\cdot ({n}_{n}{{\boldsymbol{V}}}_{n})={{\rm{\Gamma }}}_{n}^{\mathrm{rec}}+{{\rm{\Gamma }}}_{n}^{\mathrm{ion}}.\end{eqnarray} \tag{ 8 }$$

We assume n_i = n_e in these equations.

The ion and neutral momentum equations are given by

$$\begin{eqnarray}&&\displaystyle \frac{\partial }{\partial t}({m}_{i}{n}_{i}{{\boldsymbol{V}}}_{i})+{\rm{\nabla }}\cdot ({m}_{i}{n}_{i}{{\boldsymbol{V}}}_{i}{{\boldsymbol{V}}}_{i}+{{\mathbb{P}}}_{i}+{{\mathbb{P}}}_{e})\\ &&\quad =\,{\boldsymbol{J}}\times {\boldsymbol{B}}+{{\boldsymbol{R}}}_{i}^{\mathrm{in}}+{{\rm{\Gamma }}}_{i}^{\mathrm{ion}}{m}_{i}{{\boldsymbol{V}}}_{n}-{{\rm{\Gamma }}}_{n}^{\mathrm{rec}}{m}_{i}{{\boldsymbol{V}}}_{i}\\ &&\qquad +\,{{\rm{\Gamma }}}^{\mathrm{cx}}{m}_{i}({{\boldsymbol{V}}}_{n}-{{\boldsymbol{V}}}_{i})+{{\boldsymbol{R}}}_{{in}}^{\mathrm{cx}}-{{\boldsymbol{R}}}_{{ni}}^{\mathrm{cx}},\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial }{\partial t}({m}_{i}{n}_{n}{{\boldsymbol{V}}}_{n})+{\rm{\nabla }}\cdot ({m}_{i}{n}_{n}{{\boldsymbol{V}}}_{n}{{\boldsymbol{V}}}_{n}+{{\mathbb{P}}}_{n})\\ &&\quad =\,-{{\boldsymbol{R}}}_{i}^{{in}}+{{\rm{\Gamma }}}_{n}^{\mathrm{rec}}{m}_{i}{{\boldsymbol{V}}}_{i}-{{\rm{\Gamma }}}_{i}^{\mathrm{ion}}{m}_{i}{{\boldsymbol{V}}}_{n}\\ &&\qquad +\,{{\rm{\Gamma }}}^{\mathrm{cx}}{m}_{i}({{\boldsymbol{V}}}_{i}-{{\boldsymbol{V}}}_{n})-{{\boldsymbol{R}}}_{{in}}^{\mathrm{cx}}+{{\boldsymbol{R}}}_{{ni}}^{\mathrm{cx}}.\end{eqnarray} \tag{ 10 }$$

Here, the current density is

$$\begin{eqnarray}&&{\boldsymbol{J}}={{en}}_{i}({{\boldsymbol{V}}}_{i}-{{\boldsymbol{V}}}_{e})=\displaystyle \frac{{\rm{\nabla }}\times {\boldsymbol{B}}}{{\mu }_{0}}.\end{eqnarray} \tag{ 11 }$$

Note that the Lorenz force ${\boldsymbol{J}}\times {\boldsymbol{B}}$ only acts on the plasma but not the neutrals. The momentum transfer ${{\boldsymbol{R}}}_{\alpha }^{\alpha \beta }$ is the transfer of momentum to species α due to identity-preserving collisions with species β:

$$\begin{eqnarray}&&{{\boldsymbol{R}}}_{\alpha }^{\alpha \beta }={m}_{\alpha \beta }{n}_{\alpha }{\nu }_{\alpha \beta }({{\boldsymbol{V}}}_{\beta }-{{\boldsymbol{V}}}_{\alpha }),\end{eqnarray} \tag{ 12 }$$

where ${m}_{\alpha \beta }={m}_{\alpha }{m}_{\beta }/({m}_{\alpha }+{m}_{\beta })$. The collision frequency ${\nu }_{\alpha \beta }$ is

$$\begin{eqnarray}&&{\nu }_{\alpha \beta }={n}_{\beta }{{\rm{\Sigma }}}_{\alpha \beta }\sqrt{\displaystyle \frac{8{k}_{{\rm{B}}}{T}_{\alpha \beta }}{\pi {m}_{\alpha \beta }}},\end{eqnarray} \tag{ 13 }$$

with ${T}_{\alpha \beta }=({T}_{\alpha }+{T}_{\beta })/2$. We choose the cross section ${{\rm{\Sigma }}}_{{in}}={{\rm{\Sigma }}}_{{ni}}=5\times {10}^{-19}$ m² as in Leake et al. (2013) and Murphy & Lukin (2015). The momentum transfer from species β to species α due to CX is ${{\boldsymbol{R}}}_{\alpha \beta }^{\mathrm{cx}}$, and the appropriate approximations for these terms between ions and neutrals as presented in Leake et al. (2012) are

$$\begin{eqnarray}&&{{\boldsymbol{R}}}_{{in}}^{\mathrm{cx}}\approx -{m}_{i}{\sigma }_{\mathrm{cx}}({V}_{\mathrm{cx}}){n}_{i}{n}_{n}{{\boldsymbol{V}}}_{{in}}{V}_{{Tn}}^{2}\\ &&\quad \times \,{\left[4\left(\displaystyle \frac{4}{\pi }{V}_{{Ti}}^{2}+{V}_{{in}}^{2}\right)+\displaystyle \frac{9\pi }{4}{V}_{{Tn}}^{2}\right]}^{-1/2}\end{eqnarray} \tag{ 14 }$$

and

$$\begin{eqnarray}&&{{\boldsymbol{R}}}_{{ni}}^{\mathrm{cx}}\approx {m}_{i}{\sigma }_{\mathrm{cx}}({V}_{\mathrm{cx}}){n}_{i}{n}_{n}{{\boldsymbol{V}}}_{{in}}{V}_{{Ti}}^{2}\\ &&\quad \times \,{\left[4\left(\displaystyle \frac{4}{\pi }{V}_{{Tn}}^{2}+{V}_{{in}}^{2}\right)+\displaystyle \frac{9\pi }{4}{V}_{{Ti}}^{2}\right]}^{-1/2}.\end{eqnarray} \tag{ 15 }$$

The pressure tensor for species α is ${{\mathbb{P}}}_{\alpha }={P}_{\alpha }{\mathbb{I}}+{\pi }_{\alpha }$, where ${P}_{\alpha }$ is the scalar pressure and ${\pi }_{\alpha }$ is the viscous stress tensor, given by ${\pi }_{\alpha }=-{\xi }_{\alpha }[{\rm{\nabla }}{{\boldsymbol{V}}}_{\alpha }+{({\rm{\nabla }}{{\boldsymbol{V}}}_{\alpha })}^{\top }]$ with ${\xi }_{\alpha }$ as the isotropic dynamic viscosity coefficient.

Combining the electron and ion energy equation together and neglecting terms of order ${({m}_{i}/{m}_{p})}^{1/2}$ and higher, one can get:

$$\begin{eqnarray}&&\displaystyle \frac{\partial }{\partial t}({\varepsilon }_{i}+\displaystyle \frac{{P}_{e}}{\gamma -1})+{\rm{\nabla }}\cdot ({\varepsilon }_{i}{{\boldsymbol{V}}}_{i}+\displaystyle \frac{{P}_{e}{{\boldsymbol{V}}}_{e}}{\gamma -1}\\ &&\quad +\,{{\boldsymbol{V}}}_{i}\cdot {{\mathbb{P}}}_{i}+{{\boldsymbol{V}}}_{e}\cdot {{\mathbb{P}}}_{e}+{{\boldsymbol{h}}}_{i}+{{\boldsymbol{h}}}_{e})\\ &&=\,{\boldsymbol{j}}\cdot {\boldsymbol{E}}+{{\boldsymbol{V}}}_{i}\cdot {{\boldsymbol{R}}}_{i}^{\mathrm{in}}+{Q}_{i}^{\mathrm{in}}-{{\rm{\Gamma }}}_{n}^{\mathrm{rec}}\displaystyle \frac{1}{2}{m}_{i}{{\boldsymbol{V}}}_{i}^{2}-{Q}_{n}^{\mathrm{rec}}\\ &&+\,{{\rm{\Gamma }}}_{i}^{\mathrm{ion}}\left(\displaystyle \frac{1}{2}{m}_{i}{{\boldsymbol{V}}}_{n}^{2}-{\phi }_{\mathrm{eff}}\right)+{Q}_{i}^{\mathrm{ion}}+{{\rm{\Gamma }}}^{\mathrm{cx}}\displaystyle \frac{1}{2}{m}_{i}({{\boldsymbol{V}}}_{n}^{2}-{{\boldsymbol{V}}}_{i}^{2})\\ &&+\,{{\boldsymbol{V}}}_{n}\cdot {{\boldsymbol{R}}}_{{in}}^{\mathrm{cx}}-{{\boldsymbol{V}}}_{i}\cdot {{\boldsymbol{R}}}_{{ni}}^{\mathrm{cx}}+{Q}_{{in}}^{\mathrm{cx}}-{Q}_{{ni}}^{\mathrm{cx}},\end{eqnarray} \tag{ 16 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\varepsilon }_{n}}{\partial t}+{\rm{\nabla }}\cdot ({\varepsilon }_{n}{{\boldsymbol{V}}}_{n}+{{\boldsymbol{V}}}_{n}\cdot {{\mathbb{P}}}_{n}+{{\boldsymbol{h}}}_{n})\\ &&\quad =\,-{{\boldsymbol{V}}}_{n}\cdot {{\boldsymbol{R}}}_{i}^{\mathrm{in}}+{Q}_{n}^{\mathrm{ni}}-{{\rm{\Gamma }}}_{i}^{\mathrm{ion}}\displaystyle \frac{1}{2}{m}_{i}{{\boldsymbol{V}}}_{n}^{2}-{Q}_{i}^{\mathrm{ion}}\\ &&+\,{{\rm{\Gamma }}}_{n}^{\mathrm{rec}}\displaystyle \frac{1}{2}{m}_{i}{{\boldsymbol{V}}}_{i}^{2}+{Q}_{n}^{\mathrm{rec}}+{{\rm{\Gamma }}}^{\mathrm{cx}}\displaystyle \frac{1}{2}{m}_{i}({{\boldsymbol{V}}}_{i}^{2}-{{\boldsymbol{V}}}_{n}^{2})\\ &&+\,{{\boldsymbol{V}}}_{i}\cdot {{\boldsymbol{R}}}_{{ni}}^{\mathrm{cx}}-{{\boldsymbol{V}}}_{n}\cdot {{\boldsymbol{R}}}_{{in}}^{\mathrm{cx}}+{Q}_{{ni}}^{\mathrm{cx}}-{Q}_{{in}}^{\mathrm{cx}}.\end{eqnarray} \tag{ 17 }$$

The internal energy density is ${\varepsilon }_{\alpha }\equiv {m}_{\alpha }{n}_{\alpha }{{\boldsymbol{V}}}_{\alpha }^{2}/2+{P}_{\alpha }/(\gamma -1)$. The term ${{\rm{\Gamma }}}_{i}^{\mathrm{ion}}{\phi }_{\mathrm{eff}}$ represents assumed optically thin radiative losses with ${\phi }_{\mathrm{eff}}=33\,\mathrm{eV}$ (Leake et al. 2012). The ratio of specific heats is denoted by γ. ${Q}_{\alpha }^{\alpha \beta }$ represents the heating of species α due to interaction with species β, which is a combination of frictional heating and a thermal transfer between the two populations: ${Q}_{\alpha }^{\alpha \beta }={{\boldsymbol{R}}}_{\alpha }^{\alpha \beta }\cdot ({{\boldsymbol{V}}}_{\beta }-{{\boldsymbol{V}}}_{\alpha })+3{m}_{\alpha \beta }{n}_{\alpha }{\nu }_{\alpha \beta }({T}_{\beta }-{T}_{\alpha })$. The electron and ion heat conduction terms ${{\boldsymbol{h}}}_{e}$ and ${{\boldsymbol{h}}}_{i}$ are given by

$$\begin{eqnarray}&&{{\boldsymbol{h}}}_{\alpha }=[{\kappa }_{\parallel ,\alpha }\hat{{\boldsymbol{b}}}\hat{{\boldsymbol{b}}}+{\kappa }_{\perp ,\alpha }({\mathbb{I}}-\hat{{\boldsymbol{b}}}\hat{{\boldsymbol{b}}})]\cdot {\rm{\nabla }}{k}_{{\rm{B}}}{T}_{\alpha },\end{eqnarray} \tag{ 18 }$$

where ${\kappa }_{\parallel ,\alpha }$ and ${\kappa }_{\perp ,\alpha }$ are the conductivity coefficients (Braginskii 1965), which are parallel and perpendicular to the magnetic field direction $\hat{{\boldsymbol{b}}}$, respectively. The isotropic neutral heat conduction is ${{\boldsymbol{h}}}_{n}=-{\kappa }_{n}{\rm{\nabla }}{k}_{{\rm{B}}}{T}_{n}$. The changes in thermal energy of ionized and neutral plasma components due to ionization or recombination are ${Q}_{i}^{\mathrm{ion}}={{\rm{\Gamma }}}_{i}^{\mathrm{ion}}(3/2){k}_{{\rm{B}}}{T}_{n}$ and ${Q}_{n}^{\mathrm{rec}}={{\rm{\Gamma }}}_{n}^{\mathrm{rec}}(3/2){k}_{{\rm{B}}}{T}_{i}$. ${Q}_{\alpha \beta }^{\mathrm{cx}}$ denotes heat flow from species β to species α due to CX (Meier 2011; Meier & Shumlak 2012).

The generalized Ohm's law for this work is given by

$$\begin{eqnarray}&&{\boldsymbol{E}}+{{\boldsymbol{V}}}_{i}\times {\boldsymbol{B}}=\eta {\boldsymbol{J}}+\displaystyle \frac{{\boldsymbol{J}}\times {\boldsymbol{B}}}{{{en}}_{i}}\\ &&\quad -\,\displaystyle \frac{{\rm{\nabla }}\cdot {{\mathbb{P}}}_{e}}{{{en}}_{i}}-\displaystyle \frac{{m}_{e}{\nu }_{{en}}}{e}({{\boldsymbol{V}}}_{i}-{{\boldsymbol{V}}}_{n}).\end{eqnarray} \tag{ 19 }$$

Here, the magnetic diffusion coefficient includes both the ion-electron and electron-neutral collisions and is written as

$$\begin{eqnarray}&&\eta =\displaystyle \frac{{m}_{e}{n}_{e}({\nu }_{{ei}}+{\nu }_{{en}})}{{({{en}}_{e})}^{2}}.\end{eqnarray} \tag{ 20 }$$

We normalize the equations using the characteristic plasma density and magnetic field around the TMR of the low solar chromosphere. The characteristic plasma number density is ${n}_{\star }\,={10}^{21}$ m⁻³, and the characteristic magnetic field is ${B}_{\star }\,=0.05$ T = 500 G. We focus on small-scale magnetic reconnection processes and choose a characteristic length of ${L}_{\star }=100$ m. From these quantities, we derive additional normalizing values to be ${V}_{\star }\equiv {B}_{\star }/\sqrt{{\mu }_{0}{m}_{p}{n}_{\star }}=34.613$ km s⁻¹, ${t}_{\star }\equiv {L}_{\star }/{V}_{\star }={t}_{A}\,=0.0029\,{\rm{s}}$, and ${T}_{\star }\equiv {B}_{\star }^{2}/({k}_{{\rm{B}}}{\mu }_{0}{n}_{\star })=1.441\times {10}^{5}$ K. The initial ionized and neutral fluid densities are set to be uniform with the neutral particle number density of ${n}_{n0}=0.5{n}_{\star }=0.5\,\times {10}^{21}$ m⁻³, and the initial ionization degree is ${f}_{i0}={n}_{i0}/({n}_{i0}\,+{n}_{n0})=0.01 \%$. Thus, the neutral-ion collisional mean free path of the background plasma is ${\lambda }_{{ni}0}=23.74$ m, and the ion inertial length is ${d}_{i0}=0.99$ m. Since the hydrogen gas only has the ground state and the ionized state in this model, the initial temperatures of the ionized and neutral fluids are set to be uniform at ${T}_{i0}={T}_{n0}=8400$ K to keep the ionization degree the same as that around the TMR in solar atmosphere.

The dimensionless magnetic diffusion form is

$$\begin{eqnarray}&&\eta ={\eta }_{{ei}}+{\eta }_{{en}}={\eta }_{{ei}\star }{T}_{e}^{-1.5}+{\eta }_{{en}\star }{({T}_{e}+{T}_{n})}^{1/2}\displaystyle \frac{{n}_{n}}{{n}_{e}},\end{eqnarray} \tag{ 21 }$$

where ${\eta }_{{ei}\star }=7.457\times {10}^{-6}$ and ${\eta }_{{en}\star }=1.369\times {10}^{-6}$ are normalization constants derived from the characteristic values n_⋆, B_⋆, and ${L}_{\star }$. T_e, T_n, n_e, and n_n are the dimensionless temperatures and number densities for the electron and neutral fluids, respectively. Since the electron and the ion are assumed to be coupled together and only the hydrogen gas is considered in our model, we assume T_i = T_e, n_i = n_e, and the pressure of the ionized component is twice the ion (or electron) pressure, ${P}_{p}={P}_{e}+{P}_{i}=2{P}_{i}$. The initial magnetic diffusion contribution from electron-neutral collisions is about one magnitude higher than that due to electron-ion collisions. However, as is shown below, magnetic diffusion during reconnection becomes dominated by the electron-ion collisions as the ionization degree increases and ${n}_{n}/{n}_{e}$ decreases with time inside the current sheet. We have also calculated the Lundquist number $S={{LV}}_{A}/\eta$ using the magnetic diffusion, the Alfvén speed, and the length of the current sheet in our simulations. Since the length scale is very small, the Lundquist number during the magnetic reconnection process is only around 2000 in this work. In agreement with previous simulations of reconnecting current sheets (e.g., Huang & Bhattacharjee 2010; Leake et al. 2012; Ni et al. 2012, 2013; Comisso et al. 2016; Huang et al. 2017), no plasmoid instability appears in our simulations.

The simulations are initialized with a force-free Harris sheet magnetic equilibrium, where the reconnecting magnetic field is along the y coordinate direction and the guide B-field is in the z direction. Specifically, the initial dimensionless magnetic flux in the z direction is given by

$$\begin{eqnarray}&&{A}_{z0}(y)=-{b}_{p}{\lambda }_{\psi }\mathrm{ln}\left[\cosh \left(\displaystyle \frac{y}{{\lambda }_{\psi }}\right)\right],\end{eqnarray} \tag{ 22 }$$

where b_p is the strength of the the magnetic field and ${\lambda }_{\psi }$ is the initial thickness of the current sheet. The initial magnetic field in the z direction is

$$\begin{eqnarray}&&{B}_{z0}(y)={b}_{p}/\left[\cosh \left(\displaystyle \frac{y}{{\lambda }_{\psi }}\right)\right].\end{eqnarray} \tag{ 23 }$$

We have simulated six cases: Case A, Case A0, Case B, Case C, Case D, and Case E. In order to see how the non-equilibrium ionization–recombination effect impacts the reconnection process, we have eliminated this effect in Case A0. In Case A0, the recombination rate does not depend on the plasma density and temperature as in Equations (2) and (4). Instead, we compel the recombination rate ${{\rm{\Gamma }}}_{i}^{\mathrm{rec}}$ to equal the ionization rate ${{\rm{\Gamma }}}_{n}^{\mathrm{ion}}$ at all times, and the right sides of Equations (7) and (8) are always zero in this case. The non-equilibrium ionization–recombination is included in all the other cases. Except for the non-equilibrium ionization–recombination factor, all the other conditions in Case A0 are the same as in Case A. The electron-neutral collisions and the Hall term are only included in Case B, but the initial conditions in Case A and Case B are identical. The only difference among Case A, Case C, Case D, and Case E is the strength of the initial magnetic field: b_p = 1 in Case A, b_p = 0.2 in Case C, b_p = 2 in Case D, and b_p = 3 in Case E. Therefore, one can calculate the initial plasma β in each case: ${\beta }_{0}=0.058$ in Case A, A0 and B, ${\beta }_{0}=1.46$ in Case C, ${\beta }_{0}=0.0145$ in Case D, and ${\beta }_{0}=0.0064$ in Case E. The differences in initial conditions and evolution equations among these six cases are summarized in Table 1. The reconnection processes are symmetric in both the x and y directions in Cases A, A0, C, D, and E. We therefore only simulate one quarter of the domain ($0\lt x\lt 2$, $0\lt y\lt 1$) in Cases A, A0, C, D, and E. The Hall effect might produce the asymmetries (Lukin & Jardin 2003) and we simulate the full domain ($-2\lt x\lt 2$, $-1\lt y\lt 1$) in Case B. We use the same outer boundary conditions at $| y| =1$ and the same form of the initial electric field perturbations as Murphy & Lukin (2015) to initiate magnetic reconnection in all of the cases in this work. The perturbation electric field is applied for $0\leqslant t\leqslant 1$. The perturbation magnitude is proportional to the value of b_p in each of the cases with an amplitude of $\delta E={10}^{-3}{b}_{p}$. The periodicity of the physical system is imposed in the x direction at $| x| =2$.

Table 1.  The Differences in Initial Conditions and Evolution Equations among the Six Simulation Cases

	Case A	Case A0	Case B	Case C	Case D	Case E
bp	1	1	1	0.2	2	3
with Hall and ${\nu }_{{en}}$	No	No	Yes	No	No	No
with non-equilibrium
ionization and
recombination	Yes	No	Yes	Yes	Yes	Yes

Download table as:  ASCIITypeset image

All six cases have been tested using a lower resolution and a higher resolution. The highest number of grid elements in Cases A, A0, C, D, and E is m_x = 64 elements in the x direction and m_y = 64 elements in the y direction. The highest number of grid elements in Case B is m_x = 128 and m_y = 128. We use sixth-order basis functions for all simulations, resulting in an effective total grid size $({M}_{x},{M}_{y})\,=6({m}_{x},{m}_{y})$. Grid-packing is used to concentrate mesh in the reconnection region. Therefore, the mesh-packing along the y direction is concentrated in a thin region near y = 0. Note that the quantities shown in the figures in this work are in dimensionless units, except for temperature and velocity, which can be used to compare with the temperature and velocity derived from a model of the line-forming process and observations of spectral line profiles.

We have simulated magnetic reconnection with different initial magnetic field strengths, resulting in different plasma β. How the plasma β affects the reconnection process in initially weakly ionized plasmas around TMR will be presented in this subsection.

First, we compare the reconnection process in Case A with initial ${\beta }_{0}=0.058$ and Case C with initial ${\beta }_{0}=1.46$. Figure 1 shows the current density J_z and ionization degree f_i in one quarter of the domain at three different times in Case A and Case C, respectively. The current sheet lengths in Case A at t = 4.024 and in Case C at t = 20.052 are the same, as are those shown in Figures 1(b) and (e), and those shown in Figures 1(c) and (f). As presented in Figure 1, about 40% of the neutral particles are eventually ionized inside the current sheet region in Case A. In contrast, while the maximum ionization degree in Case C increased during the reconnection process by a factor of 300 from $0.01 \%$ to 3%, the plasmas are still weakly ionized even during the later stage of the reconnection process. Figure 3 shows that the neutral fluid can be fully ionized inside the current sheet when the magnetic field is strong enough, as in Case D and Case E.

Zoom In Zoom Out Reset image size

Figure 2 shows the profiles of ion and neutral temperatures across the current sheet in Cases A and C at the same three pairs of times as in Figure 1. It is apparent that (1) in both cases, the ion and neutral temperatures are nearly equal throughout the evolution due to rapid thermal exchange between the plasma components; (2) the initial temperature increase observed in Case A due to the Joule heating within the broad Harris current sheet prior to the onset of the reconnection process is later moderated by radiative cooling; and (3) the peak temperature within the narrow reconnection current sheet is maintained at about $1.6\times {10}^{4}$ K in both cases. As shown in Figure 3, the plasma inside the current sheet can only be heated to above 2 × 10⁴ K when the plasma becomes fully ionized, and the highest temperature can reach around $4.6\times {10}^{4}$ K in Case E.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We have calculated the time evolution of the total radiated energy ${Q}_{\mathrm{rad}}={\int }_{0}^{1}{\int }_{0}^{2}{L}_{\mathrm{rad}}{dxdy}$, the Joule heating ${Q}_{\mathrm{Joule}}\,={\int }_{0}^{1}{\int }_{0}^{2}\eta {J}^{2}{dxdy}$, the frictional heating between ions and neutral particles Q_in and the viscous heating of ions and neutral particles Q_vis in the whole simulation domain. Together, these four terms represent all sources and sinks of plasma thermal energy within the domain. In all the simulations described here, both Q_in and Q_vis are found to be several magnitudes smaller than Joule heating Q_Joule, (not shown). The time evolution of these quantities for Case A is plotted in Figure 4(a), showing that the radiated thermal energy strongly increases once the reconnection current sheet is formed. This is consistent with our interpretation of the temperature evolution shown in Figure 2(a). Furthermore, it is shown that most of the total generated thermal energy ${Q}_{\mathrm{Joule}}+{Q}_{{in}}+{Q}_{\mathrm{vis}}$ is radiated during the rapid ionization–recombination cycle in Case A. In Case C, which has a high plasma β and weak magnetic field, the radiated thermal energy increases much slower with time during the reconnection process, as shown in Figure 4(b). However, a large amount of the thermal energy is similarly radiated during the later stage of the reconnection process.

Zoom In Zoom Out Reset image size

From our simulations results, we find that the ionized and neutral fluids are coupled well in the reconnection outflow regions. This is consistent with the previous results by Leake et al. (2012, 2013), and those by Murphy & Lukin (2015). These prior studies in the high β regime also showed that the neutral and ionized fluid components decouple upstream of the reconnection site on scales smaller than ${\lambda }_{{ni}}$. As shown in Figure 5, similar behavior is observed in our Case C simulation, but the neutral and ion inflows are observed to be well-coupled in the steady-state reconnection phase in Case A, as well as in the other low β Cases: B, D, and E (not shown). We note that for the low β plasma in Case A, there is a positive peak value of ${V}_{{iy}}\mbox{--}{V}_{{ny}}$ in the upstream of the reconnection site; the peak value of ${V}_{{iy}}\mbox{--}{V}_{{ny}}$ is at around y = 0.125 at t = 20.126, which means the speed of the inflowing neutral fluid is even greater than that of the ion fluid at this position.

Zoom In Zoom Out Reset image size

One possible explanation for the observed differences between the high β and low β cases lies in the degree of plasma ionization self-consistently formed within the reconnection current sheets. We note that for two plasma elements with the same total atom particle density, a plasma element that is 50% ionized would have the neutral-ion mean free path be 50 times shorter than that which is 1% ionized. Thus, the decoupling between ion and neutral fluid inflows observed in the high β simulations that remain weakly ionized can be suppressed in the low β simulations where the plasma becomes strongly ionized.

Panels (c) and (f) of Figure 5 show the terms contributing to the pressure balance across the reconnection current sheet during the quasi-steady-state phase in Cases A and C, respectively. We note that in both cases neither neutral pressure nor guide-field B_z act to balance the magnetic field pressure from the reconnecting B-field component within the current sheet. As a result, as in the previous work by Leake et al. (2013), the decrease in the magnetic field pressure of the reconnecting field has to be balanced by a corresponding increase in the ionized fluid density and pressure. This observation points to a conclusion that the ionization fraction within a reconnecting current sheet in a partially ionized plasma is primarily controlled by force-balance requirements (compression) rather than by a plasma temperature increase due to Ohmic heating.

The maximum values of the inflow plasma velocity are almost the same in Cases A, C, D, and E with different plasma β, but the outflow velocity is higher in the lower plasma β case. The maximum V_ix in Case E is about 16 km s⁻¹ and it is about 7 km s⁻¹ in Case C.

We have used the same method as Leake et al. (2013) to calculate the reconnection rate, ${M}_{\mathrm{sim}}={\eta }^{* }{j}_{\max }/({V}_{A}^{* }{B}_{\mathrm{up}})$, where j_max is the maximum value of the out-of-plane current density J_z, located at $(x,y)=(0,0)$ in all the simulations in this work. B_up is B_x evaluated at $(0,{\delta }_{\mathrm{sim}})$, where ${\delta }_{\mathrm{sim}}$ is defined as the half-width at half-maximum in J_z. ${V}_{A}^{* }$ is the relevant Alfvén velocity defined using B_up and the total number density n^* at the location of j_max. ${\eta }^{* }$ is the magnetic diffusion coefficient defined by in Equation (21) at the location of j_max, where the electron-neutral collisions are only included in Case B. Figure 6(a) shows the time evolution of the reconnection rates in Cases A, B, C, and D, and Figure 6(b) shows the time evolution of the normalized magnetic diffusion coefficient ${\eta }^{* }$ at the X-point. We note that each of the reconnection rate curves includes both the formation phase of the reconnection region and the later quasi-steady-state phase; and the provided reconnection rate measure is most meaningful in the latter well-developed phase of each of the simulations. The measured reconnection rate in Case C with high β is several times higher than that in the other lower β cases. The maximum of the reconnection rate reaches above 0.1 in Case C, while the maximum rate in Case A, Case B, and Case D is only around 0.025. We do not presently have a quantitative theory to explain the measured reconnection rates; this result is consistent with the conjecture that decoupling of ionized and neutral fluid inflows, as observed in Case C but not in Cases A, B, or D, can accelerate the reconnection rate within a reconnection current sheet. Figure 6(b) further shows that the magnetic diffusion coefficient ${\eta }^{* }$ at the reconnection X-point has approximately the same value of 2 × 10⁻⁴ in Cases A, B, and C during the quasi-steady-state reconnection phase. Thus, it has no contribution to the difference in the reconnection rates among these three simulations and further provides evidence that omitting the contribution of electron-neutral collisions to the electrical resistivity in Case A, as compared to the Case B, did not significantly modify the outcome. ${\eta }^{* }$ in Case D drops to a lower value because of the heating at the reconnection X-point during the later stage of the reconnection process.

Zoom In Zoom Out Reset image size

Figures 7(a) and (b) show the four terms in Equation (7) contributing to $\partial {n}_{i}/\partial t$ in Cases A and C. The values of the four contributions are the average values inside the current sheet domain at each time. Figure 7(c) shows the variations of the half length and width of the currents versus times in Cases A and C. The half length of the current sheet abruptly drops to 0.5–0.6 during a very short period in both Cases A and C, and the half-width of the sheet gradually drops to around 0.015. In each case, the time when the measured length of the current sheet abruptly drops corresponds to the formation of a fully nonlinear reconnection current sheet and onset of the self-consistent reconnection process. In Case A, the ionization is much larger than the other three terms prior to the formation of a reconnection current sheet, but the contribution of all four terms substantially increases with the onset of reconnection. The value of $-{{\rm{\Gamma }}}_{i}^{\mathrm{rec}}+\partial ({n}_{i}{V}_{{ix}})/\partial x$ gradually gets close to the value of $-\partial ({n}_{i}{V}_{{iy}})/\partial y+{{\rm{\Gamma }}}_{i}^{\mathrm{ion}}$ during the later stage of the reconnection process, and the current sheet ion density reaches an approximate steady-state in Case A. In Case C, the ionization rate is always the maximum among the four components, as shown in Figure 7(b). However, we should also note that the total value of contributing terms to $\partial {n}_{i}/\partial t$ in Case A is about two orders of magnitude larger than that in Case C. The imbalance between $-{{\rm{\Gamma }}}_{i}^{\mathrm{rec}}+\partial ({n}_{i}{V}_{{ix}})/\partial x$ and $-\partial ({n}_{i}{V}_{{iy}})/\partial y+{{\rm{\Gamma }}}_{i}^{\mathrm{ion}}$ in Case C is small relative to the neutral fluid density. In both Leake et al. (2012, 2013) studies of null-point reconnection, the ionization rate ${{\rm{\Gamma }}}_{i}^{\mathrm{ion}}$ did not play an important dynamical role and it was always the minimum among the four components, which is significantly different from the simulation results presented here. The ionization rate ${{\rm{\Gamma }}}_{i}^{\mathrm{ion}}$ always plays an important role in the whole reconnection process in all the cases studied in this work. From Figures 7 (a) and (b), we can see that the ionization rate ${{\rm{\Gamma }}}_{i}^{\mathrm{ion}}$ is always higher than the recombination rate $-{{\rm{\Gamma }}}_{i}^{\mathrm{rec}}$ in both Cases A and C. The plasma β inside the current sheet in Case C is still smaller than that in previous work by Leake et al. (2013) without a guide field. Therefore, the plasma β and the magnetic field structures appear to determine whether ionization will be dynamically important inside the current sheet region.

Zoom In Zoom Out Reset image size

As presented in the last paragraph in the above section for Case A and Case C, the ionization rate is always higher than the recombination rate in Cases A, B, C, D, and E. In this subsection, we will present the numerical results in Case A0, in which have set the recombination rate to equal the ionization rate. By comparing Cases A and A0, we demonstrate how the non-equilibrium ionization–recombination effect impacts magnetic reconnection in initially weakly ionized plasmas around the solar TMR.

Figure 1(c) presents the distributions of the current density J_z and ionization degree f_i at t = 5.441, t = 13.13, and t = 23.381 in Case A0. One can see that the maximum current density during the later stage is higher than that in Case A. However, the maximum ionization degree in Case A0 can only reach around $0.08 \%$, which is much lower than that in Case A. The ionization rate increases extremely slowly with time in Case A0. Figure 2(c) shows the profiles of ion and neutral temperatures across the current sheet in Cases A0 at the same three time instances as in Figure 1(c). Both the plasma and neutral temperatures reach very high values inside the current sheet. The maximum plasma temperature is above 5 × 10⁴ K during the later stage of the reconnection process, which is much higher than that in Case A. These results are consistent with the previous one-fluid MHD simulations (Ni et al. 2016), where the plasmas were heated from 4200 K to very high temperatures (about $8\times {10}^{4}$ K) with reconnection magnetic fields of 500 G at around the TMR region in that work. However, it is not realistic to have such high temperatures and weakly ionized hydrogen plasmas, which is a reflection of the unphysical nature of the ionization–recombination balance imposed in this calculation. Figure 4(c) presents the total generated thermal energy and radiated energy inside the simulation domain in Case A0. Joule heating is also the dominant term to generate the thermal energy in this case, and most of the generated thermal energy is also radiated away.

The black dashed line in Figure 6(a) represents the time-dependent reconnection rate in Case A0. The method for calculating the reconnection rate is the same as that presented in the above section. Though the maximum current density in Case A0 is higher than that in Case A, the magnetic reconnection rate in Case A0 is more than two times lower than that in Case A. The much higher plasma temperature at the reconnection X-point in Case A0 results in a much lower magnetic diffusion ${\eta }^{* }$, which is the main reason for a lower reconnection rate.

By comparing the numerical results in Case A and A0, one can conclude that the non-equilibrium ionization–recombination factor leads to much lower temperature increases inside the reconnection region. This factor also leads to a faster reconnection rate in our low Lundquist number simulations. Therefore, it is very important and necessary to include the non-equilibrium ionization–recombination for studying the reconnection events around the solar TMR region, especially for answering the questions about how high the plasma temperature can be increased and if the plasma temperatures inside the reconnection region are high enough to produce Si iv emissions from the IRIS bombs around the solar TMR region.

We next describe the roles of the Hall effect and the magnetic diffusion contributed by the electron-neutral collisions. We include the Hall effect and electron-neutral collisions in Case B, and we performed calculations in the whole domain. We can then compare the evolution of each variable in the reconnection process in Cases A and B. Figure 8 shows the current density J_z, the ionization degree f_i and the magnetic field in z direction B_z in the whole domain at three different times in Case B. Comparing Figures 1(a)–(c) and 8(a)–(f), we find that current density distributions in Cases A and B are very similar and the maximum current density in Case B is slightly higher than in Case A. The maximum ionization degree can also reach about 40% in Case B. The electron-neutral collision can quickly heat the plasma temperature to above 4 × 10⁴ K at the beginning, but then the plasma temperature drops fast to low values below 2 × 10⁴ K. As shown in Figure 9, the maximum temperature in Case B eventually drops to around $1.5\times {10}^{4}$ K. Figures 9(b) and (c) show the ion inflow velocity V_iy and the difference in speeds of between ion and neutral inflow ${V}_{{iy}}\mbox{--}{V}_{{ny}}$ at x = 0 along the y direction at three different times. Comparing Figures 5(a) and (b) with Figures 9(b) and (c), we note that the plasma inflow velocity distributions are also very similar in Cases A and B; the maximum of V_iy in Case B is slightly higher than that in Case A. The plasma and the neutral inflows are also coupled well in the reconnection upstream region in Case B.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As shown in Figure 8, the Hall effect does not result in the obvious asymmetry in either the x or y directions for any of the variables. Tilting of the current sheet can only occur on the scales of the ion inertial length d_i, and d_i decreases to a very small value, even smaller than the grid size during the reconnection process in our simulations. Therefore, a tilted current sheet as shown by Lukin & Jardin (2003) does not appear in Case B in this work. From Figure 6(a), we see that the Hall effect just slightly enhances the reconnection rate. Figure 6(b) shows that the electron-neutral collision only strongly affects the magnetic diffusion at the very beginning and the magnetic diffusion coefficient contributed by electron-neutral collisions quickly drops to a low value before the reconnection process starts. Values of ${\eta }^{* }$ are eventually very close to one another in Cases A and B. Though the Hall effect increases the reconnection rate slightly and it enhances the inflow velocity and the maximum current density in the same fashion more apparently in Case B than in Case A, including the Hall effect and electron-neutral collisions does not significantly change the temporal evolution and the spatial structure of the reconnection process.