Boosting Magnetic Reconnection by Viscosity and Thermal Conduction

Nonlinear evolution of magnetic reconnection is investigated by means of magnetohydrodynamic simulations including uniform resistivity, uniform viscosity, and anisotropic thermal conduction. When viscosity exceeds resistivity (the magnetic Prandtl number Prm>1), the viscous dissipation dominates outflow dynamics and leads to the decrease in the plasma density inside a current sheet. The low-density current sheet supports the excitation of the vortex. The thickness of the vortex is broader than that of the current for Prm>1. The broader vortex flow more efficiently carries the upstream magnetic flux toward the reconnection region, and consequently boosts the reconnection. The reconnection rate increases with viscosity provided that thermal conduction is fast enough to take away the thermal energy increased by the viscous dissipation (the fluid Prandtl number Pr<1). The result suggests the need to control the Prandtl numbers for the reconnection against the conventional resistive model.

efficiently carries the upstream magnetic flux toward the reconnection region, and consequently boosts the reconnection. The reconnection rate increases with viscosity provided that thermal conduction is fast enough to take away the thermal energy increased by the viscous dissipation (the fluid Prandtl number P r < 1). The result suggests the need to control the Prandtl numbers for the reconnection against the conventional resistive model. a) Electronic mail: minoshim@jamstec.go.jp

I. INTRODUCTION
Magnetic reconnection is one of the most fundamental processes in plasma physics, in which stored magnetic energy is rapidly released and converted into kinetic and internal energies through the change of magnetic field topology. It is widely believed to play a major role in explosive phenomena such as magnetospheric substorms and stellar flares.
The reconnection intrinsically contains a hierarchical structure ranging from the fully kinetic scale to the magnetohydrodynamic (MHD) scale. In order to identify the essential physics necessary to model the reconnection, Birn et al. 1 conducted numerical simulations with a variety of codes from kinetic codes to conventional resistive MHD codes. They showed that only the MHD simulation with uniform resistivity fails to trigger fast reconnection, indicating that resistive MHD would be insufficient to model it.
The role of resistive dissipation on the reconnection has been extensively investigated in the framework of MHD. A classical Sweet-Parker model predicts the rate of the reconnection proportional to the square root of resistivity that is too slow to account for observed phenomena. Subsequent studies have demonstrated that the Sweet-Parker type current sheet undergoes secondary instabilities for sufficiently small resistivity (large Lundquist number) [2][3][4][5] .
The resulting reconnection rate seems to be independent of resistivity [6][7][8][9] . Meanwhile, the impact of other dissipation processes should be discussed. We focus on viscosity and heat transfer.
Viscosity might support the nonlinear evolution of the reconnection 10 whereas it suppresses the linear growth 11 . The ratio of kinematic viscosity to resistivity is defined as the magnetic Prandtl number, which relates the dissipation scale of vortex to current. Resistive MHD assumes this number to be zero, meaning that the vortex scale is negligible small compared with the current scale. However, it is not necessarily true in actual plasma environments 12 ; the number can be much larger than unity in a classical Spitzer model for hot tenuous plasmas 13 . Numerical simulations have demonstrated that it affects the nonlinear evolution of MHD phenomena such as small-scale turbulence and dynamo [14][15][16][17] . It may also impact on the reconnection in which small-scale dissipation processes eventually result in large-scale evolution. In the kinetic reconnection composed of collisionless ions and electrons, the reconnection region has a two-scale structure of broad ion diffusion region and narrow electron diffusion region embedded there 18,19 . This structure may be measured as broad vortex and narrow current layers from the viewpoint of MHD, because the momentum and the current are predominantly sustained by ions and electrons, respectively.
Heat transfer is associated with the reconnection. High-energy particles are produced in the vicinity of a reconnection site and stream along a magnetic field line during the collisionless reconnection [20][21][22] . In the solar flare, thermal conduction is effective along a magnetic field line and may affect the evolution of the collisional reconnection [23][24][25] . Including heat transfer increases compressibility that can enhance the reconnection [26][27][28] . The fluid Prandtl number, the ratio of kinematic viscosity to temperature conductivity, is around 10 −3 in the Spitzer model. Therefore, actual plasmas can have the following inequality for the timescale of three diffusion processes, τ heat < τ viscous < τ resistive .
In order to ascertain the effect of viscosity and heat transfer on the nonlinear evolution of the reconnection, we conduct two-dimensional MHD simulations including viscous dissipation and anisotropic thermal conduction as well as resistive dissipation.

II. MODEL
The governing equations are fully-compressible visco-resistive MHD equations coupled with anisotropic thermal conduction, ∂B ∂t where ν, η, and α denote the kinematic viscosity, the resistivity, and the temperature conductivity, I and S the unit and the strain rate tensors, κ the anisotropic thermal conductivity tensor working only along the magnetic field line, ǫ = 10 −6 a small value to avoid the division by zero, and other symbols have their usual meanings. The magnetic permeability is assumed to be unity. The three diffusion coefficients ν, η, α are assumed to be uniform and constant. We use an adiabatic index γ = 5/3.
The equations are advanced based on an operator splitting technique. The ideal MHD part is solved by a nonlinear third-order finite difference scheme coupled with the HLLD approximate Riemann solver 29 . The solenoidal condition for magnetic field is guaranteed by the HLLD Upwind Constrained Transport method 17 . We use a fourth-order central difference in space and a third-order Runge-Kutta time integration for the resistive dissipation term.
The viscous dissipation and the thermal conduction terms are discretized into second-order in space and time, and are updated by the multi-color Gauss-Seidel method.
The initial condition is a Harris current sheet configuration on background stationary plasma, ρ = (ρ 0 − δρ) cosh −2 (y/λ) + δρ, u = 0, B = (B 0 tanh(y/λ), 0, 0), P = (B 2 0 /2)(ρ/ρ 0 ), where δρ = 0.2ρ 0 is the background plasma density and λ is the thickness of the current sheet. The quantities are normalized so that ρ 0 = λ = B 0 = 1. We initiate the reconnection by adding a small perturbation to the flux function δA and a small (1%) uniform random perturbation to u y inside the current sheet. The three diffusion coefficients are expressed as dimensionless numbers of the magnetic and fluid Reynolds numbers, (R m , R e ) = (λV A0 /η, λV A0 /ν), and the magnetic and fluid Prandtl numbers, (P r m , P r) = (ν/η, ν/α), where V A0 = B 0 / ρ 0 is the Alfvén speed. We fix a resistivity value of 10 −3 (R m = 10 3 ) and vary the kinematic viscosity and the temperature conductivity. The rectangular simulation domain [0, 64λ]×[−8λ, 8λ] is resolved with uniformly-spaced grid points of 8, 192 × 2, 048. The symmetric boundary and the conducting wall boundary are set in x and y directions. We confirm that doubling the resolution in space and time does not alter our conclusion.

III. RESULT
We identify the dominant reconnection site by finding the minimum in x of magnetic flux |B x |dy, and then measure the amount of the magnetic flux reconnected there. The time evolution of the reconnected flux is shown in Figure 1(a). Its slope corresponds to the local reconnection rate. As a global measure of the reconnection, Figure 1(b) shows the time evolution of magnetic energy integrated over a whole domain. We use viscosity values of ν = 0, 10 −3 , 3 × 10 −3 , 10 −2 (P r m = 0, 1, 3, 10), which are indicated as black, green, magenda, and red lines, respectively. The temperature conductivity values are set to α = 0.3 (P r < 1), α = ν (P r = 1), and α = 0, which are shown as solid, dotted, and dashed lines, respectively. In all cases the current sheet gets thin down after the passage of the initial perturbation.
The length L and the thickness δ of the current sheet becomes 20λ and 0.06λ at T = 300 in the resistive case. Subsequently, an elongated Sweet-Parker type current sheet undergoes secondary tearing instability. Secondary plasmoids are observed around x = 0 and x = 30 in Figure 2(a). However, the number of plasmoids is small within the present simulation domain and time, and their impact on the reconnection rate is limited in this case.
The visco-resistive simulations without thermal conduction show that moderate viscosity accelerates the evolution around T = 150 prior to the onset of the secondary instability in comparison with the resistive case (dashed green and magenda lines in Figure 1). The acceleration is weakened in the most viscous case (dashed red lines). A secondary plasmoid is less observed for larger viscosity although the elongated current sheet is formed. Thermal conduction has little effect on the evolution in the resistive case and the visco-resistive cases as long as P r ≥ 1. Their time profiles are almost same as the cases without the conduction (solid black and dotted lines).
Surprisingly, the reconnection is enhanced when the kinematic viscosity is larger than the resistivity (P r m ≥ 1) and the temperature conductivity is further larger than the viscosity (P r < 1). In the case (P r m , P r) = (10, 1/30), the reconnection is initiated slower than the other cases, followed by the explosive onset at T = 180 (solid red lines in Figure 1). Figure   2(c) shows that the spatial structure is considerably different from the previous cases. A single X-type reconnection site is formed around x = 0 and the current is localized there.
The outflow speed reaches the Alfvén speed measured in the inflow region. Similar evolution is observed in all three cases for P r m ≥ 1 and P r < 1. Contrary to the linear theory, the peak reconnection rate increases with the viscosity, 0.005, 0.007, 0.01, 0.013 for ν = 0, 10 −3 , 3×10 −3 , 10 −2 , indicating that viscosity plays a key role in boosting the reconnection. Figure 3(a) shows the one-dimensional distribution across the outflow at x = 4.7 at T = 229 in the visco-resistive case including the conduction. The outflow is predominantly accelerated by the magnetic tension force and is decelerated by the viscous dissipation force. We then approximate the equilibrium in the current sheet as ∂(B x B y )/∂y+ν(∂/∂y)ρ∂u x /∂y = 0 (remaining dissipation terms are relatively small), leading to where quantities with the subscript in(out) are measured in the inflow(outflow) region. By combining it with the induction and continuity equations, the aspect ratio of the reconnection region is dominated by viscosity 30 , where V A,in = B in / ρ in is the Alfvén speed in the inflow region. In the classical Sweet-Parker model, the outflow is decelerated by the inertia force and the aspect ratio is proportional to the inverse square root of resistivity.
The inflow speed in the visco-resistive case follows u in = η/δ similar to the resistive case.
Consequently, we reduce the continuity equation (11) to The plasma can be diluted for P r m > 1 if u out = V A,in holds. Figure 4 supports this relation.
In the nonlinear stage (T > 200), (a) the outflow density decreases with increasing viscosity whereas (b) the outflow velocity u out ∼ V A,in is insensitive to it. The outflow density falls below the inflow density and subsequently the explosive reconnection is triggered. We also observe the density decreasing from the upstream toward the center of the current sheet in Figure 3(b). It rather increases toward the sheet in the resistive case.
The essence of the reconnection is the excitation of a quadrupolar vortex around a reconnection site as well as the dissipation of current. We investigate the excitation mechanism of the enstrophy Q = |ω| 2 /2 where ω = ∇ × u is the vorticity. In the two-dimensional configuration without the out-of-plane magnetic field (that is, ω = ω z e z , j = j z e z ), some vector identities reduce the enstrophy equation to The equation consists of the inertia (advection and compression) term, the Lorentz term, the baroclinic term, and the viscous term. The Lorentz term is reduced to two terms; the first proportional to the gradient of the current, and the second proportional to the gradient of the density. The first Lorentz term is the primary source for the enstrophy under the current sheet configuration since the out-of-plane vector (B · ∇)j is parallel to the quadrupolar vortex. Figure 5  Viscosity has linear and nonlinear stabilizing effects on the reconnection. The viscous dissipation leads to plasma heating via ∂P/∂t = (γ − 1)ρν(S · ∇)u ≃ (γ − 1)ρν|∇ × u| 2 .
Without convection or diffusion, the heated plasma stagnates at the current sheet and inhibits its thinning. The thinning speed of the current sheet is slower than that of the vortex in the visco-resistive case without the conduction (Figure 6(c)). The evolution is saturated when the thickness of the current sheet exceeds the vortex at T = 280. The reconnection rate remains slow in this case.
The simulation indicates that the visco-resistive reconnection can evolve against the inhibition in the presence of thermal conduction. Thermal conduction will decrease the temperature around the reconnection region that leads to the compression of the plasma to satisfy the force balance. It can be accompanied by the upstream magnetic flux provided that the spatial scale of the fluid is broader than that of the magnetic field. It is reasonable to speculate that P r < 1 is a necessary condition to sustain the viscoresistive reconnection since thermal conduction is expected to be fast enough to take away the thermal energy increased by the viscous dissipation. In order to confirm the effect of thermal conduction, Figure 7(a) shows the temperature at the reconnection site as a function of the amount of the reconnected flux. The peak temperature is roughly proportional to the viscosity. When P r m = 1, 3 (green and magenda), thermal conduction is not effective for the cases P r ≥ 1 whereas it successfully decreases the temperature and boosts the reconnection for P r < 1. When P r m = 10 (red), on the other hand, thermal conduction succeeds in decreasing the temperature even with P r = 1. The peak temperature values are 200 for P r = ∞, 40 for P r = 1, and 20 for P r = 1/30. Figure 7(b)-(d) shows the temperature distribution at which the reconnected flux is about the same for these three cases. High temperature plasma is observed in a quite narrow layer in the case (b), indicating a lack of a way to redistribute the thermal energy. The temperature distribution in the case (c) is somewhat diffusive. However, narrow high temperature layer and a large plasmoid remains around x = 10 −30 and x = 0, and the subsequent evolution is similar to the case (b). Thus, thermal conduction with P r = 1 is not fast enough to sustain the visco-resistive reconnection compared with the case (d), in which high temperature region becomes broader toward the downstream.

IV. DISCUSSION
Based on the visco-resistive MHD simulation coupled with anisotropic thermal conduction, we propose the viscosity-dominated reconnection model in Figure 8. Viscosity and thermal conduction can be a key to boost the reconnection in the MHD regime. The reconnection rate is found to increase with viscosity within the explored range provided that thermal conduction is fast enough. However, the reconnection rate in the present model is still slower than that in kinetic models 31,32 . The Hall-MHD model is thought to be a minimal model for fast reconnection 33 . Terasawa 34 carried out a linear analysis of the resistive tearing instability including the Hall effect, and argued that the thickness of the vortex is broader than the current and the broad vortex flow is expected to enhance the reconnection. One of differences between the two models is the presence or absence of dissipation. The Hall effect is purely dispersive mode without dissipation. It should more efficiently convert magnetic energy (current) into kinetic energy (vortex) than the viscosity-dominated reconnection.
The inflow speed is characterized by the whistler rather than the Alfvén speed in the kinetic regime, u in ∼ V A,in d i /L, where d i is the ion inertia length [35][36][37] . If we assume u out ∼ V A,in and relate the ratio V A,in /u in to the aspect ratio of the reconnection region, the viscositydominated model (eq. (12)) gives the effective viscosity for the kinetic reconnection as It also gives the effective resistivity as η eff ∼ u in δ ∼ V A,in d i (δ/L). The effective magnetic Prandtl number P r m,eff ∼ d i /δ may be larger than unity in the kinetic reconnection in which the thickness of the current sheet gets thinner than the ion inertia length (down to the electron inertia length). Transition from the slow resistive MHD to the fast Hall-MHD reconnection can be observed when the current sheet thickness falls below the ion inertia length 38 .
The decreasing density distribution from the upstream toward the current sheet is not a situation observed only in the viscosity-dominated reconnection. It has been extensively studied that ad hoc localized resistivity triggers fast reconnection in the MHD regime, the so-called Petscheck-type reconnection. The decreasing density distribution is observed in the Petscheck-type reconnection due to fast magnetosonic rarefaction waves emanated from the reconnection site 39 . The rarefaction wave drives the upstream plasma toward the reconnection site since ∇ · u ∼ −(u y /ρ)∂ρ/∂y > 0. The dilatation in the upstream is also seen in the viscosity-dominated reconnection. Furthermore, thermal conduction facilitates the Petscheck-type reconnection 24 . Observation of the solar flare supports these theoretical models 40 .
Currently, it remains unclear whether the viscosity-dominated reconnection is controlled by diffusion coefficients (η, ν, α) or their ratio (P r m , P r) (or both), because we have fixed the resistivity value. Subsequent works will investigate the dependence on resistivity (Lundquist number), which is a key parameter to classify the reconnection dynamics 41 . We may anticipate that the Prandtl numbers control the dynamics to some extent because P r m > 1 leads to the two-scale structure of the reconnection region and P r < 1 is required to sustain the boost. This implies that one cannot ignore viscosity and thermal conduction whenever they exceed resistivity, even if their absolute values are small. The condition P r m ≫ 1 and