Most of the plasma fluid equations have employed the electrical resistivity to generate the magnetic dissipation required for magnetic reconnection to occur in collisionless plasma. However, there has been no clear evidence that such a model is indeed appropriate in the reconnection diffusion region in terms of the kinetic physics. The present study demonstrates that, using a large-scale 3D kinetic simulation and analytical analysis, the spatial distribution of the non-ideal electric field is consistent with the dissipation due to the viscosity rather than the resistivity, when electromagnetic (EM) turbulence is dominant in the electron diffusion region (EDR). The effective viscosity is caused by the EM turbulence that is driven by the flow shear instabilities leading to the electron momentum transport across the EDR. The result suggests a fundamental modification of the fluid equations using the resistivity in the Ohm's law. In contrast, for the 2D current sheet without significant turbulence activity, the non-ideal field profile does not obey the simple form based on the viscosity, so that further investigation is needed for a better description.
I. INTRODUCTION
Magnetic reconnection is an important plasma process that enables fast energy conversion from the magnetic field energy into the plasma kinetic and thermal energies, accompanied often by the topological change in the magnetic field lines. The reconnection process has a significant impact on large-scale dynamics in space, such as solar flares1,2 and magnetospheric substorms.3 On the other hand, the process is driven by the magnetic dissipation that takes place in a kinetic-scale region, so-called the diffusion region, to break the field line connectivity.4 Most of the plasma fluid equations have employed the electrical resistivity to generate the magnetic dissipation required for magnetic reconnection to occur. However, there has been no clear evidence that such a description is appropriate for the diffusion region in collisionless plasma where the Coulomb collision frequencies are negligibly small.
Global responses of magnetic reconnection have mainly been investigated by means of the magnetohydrodynamics (MHD) simulations. Although the MHD simulations are useful in describing large-scale dynamics, the kinetic processes are not explicitly incorporated. The magnetic dissipation in most MHD simulations occurs from ad hoc resistivity that is artificially given in the Ohm's law in the form of , where η is the resistivity, and and are the electric field and current density, respectively. The resistive MHD model enables magnetic reconnection to occur in the simulations, while the large-scale behavior of reconnection and the global responses depend sensitively on the resistivity model.5–8 This fact suggests the potential importance of the dissipation mechanism based on the kinetic physics, not only for magnetic reconnection itself but also for global process in space and astrophysics.
The resistive MHD model relies on the assumption that the resistivity is caused by momentum exchange between ions and electrons. In collisionless plasma, it has been suggested that the momentum exchange can be generated by wave–particle interactions through plasma turbulence.9 However, there have been no convincing physical models of the resistivity that originates from turbulence in the reconnection current layer and is large enough to drive a fast reconnection. It should be also noted that there is no guarantee that is proportional to in the Ohm's law in collisionless plasma. In fact, the momentum exchange can result from drift-type instabilities excited by relative velocity between ions and electrons. Such instabilities tend to be stabilized in the region of high plasma beta10,11 that is typical in the reconnection current layer with anti-parallel configuration.
Recent kinetic simulations in a 3D system have shown that the reconnection current layer can be unstable to flow shear instabilities rather than drift-type instabilities, both in high βe (Ref. 12) and in low βe (Ref. 13) plasmas, where βe is the electron plasma beta. The reconnection layer situation is similar to that of laminar jet flows in the background stationary fluid (the so-called Bickley jet). Such the flows are unstable to flow shear instabilities that work to relax the velocity gradient.14 The shear modes transport the flow momentum across the shear layer, which gives rise to viscosity in the ensemble-averaged shear layer. For the reconnection current layer, the electron flow dominantly contributes to the total current density around the x-line, forming the electron diffusion region (EDR). It has been revealed that the EDR in high βe plasma is unstable to the current sheet shear instability (CSSI)15 and the electron Kelvin–Helmholtz instability (eKHI),16 leading to electromagnetic (EM) turbulence.
The CSSI is an electromagnetic instability excited by the ion and electron flow shears across the current sheet along the inflow direction. The linear dispersion relation is obtained for a current sheet with uniform plasma density and shearing ion and electron flows that is formed typically across the EDR during a fast reconnection. The current sheet structure for the linear analyses differs from the Harris sheet17 that has a large density gradient and is widely used for the initial condition of the simulations. As reconnection proceeds, the Harris sheet plasma is convected away from the EDR and replaced with the upstream low density plasma. Thus, the density profile is almost uniform across the EDR during a quasi-steady reconnection. The electromagnetic perturbation is induced through the electron inertia, while the ions are affected by the perturbed Lorentz force, which determines the typical wavelength intermediate between the ion and electron scales. The CSSI is different from previous flow shear instabilities, such as an electron-ion hybrid instability18 that is purely electrostatic and electron MHD instabilities19–22 that lack the coupling with ions. Meanwhile, the eKHI is excited by the electron flow shear along the outflow direction. Since the ions have no significant flow shear in this direction, the waves hardly couple the ions, so that the mode can be described in the electron MHD framework.19,23
Since these instabilities are mainly excited by the electron flow shears, the resultant turbulence transports the electron flow momentum away from the current sheet center, inducing the magnetic dissipation in the EDR. Our previous study24 showed that the EM turbulence, indeed, gives rise to the viscous dissipation at the reconnection x-line. The corresponding Ohm's law was consistent with the form of rather than at the x-line, where μ is the effective viscosity due to the turbulence. However, the spatial distribution of the dissipation and Ohm's law across the turbulent current layer has not been revealed yet, regardless of its potential impact on the large-scale reconnection process in MHD.
The present study focuses on the spatial distribution across the EDR of the non-ideal electric field generated by the EM turbulence that is responsible for the magnetic dissipation. Using a large-scale particle-in-cell (PIC) simulation in a 3D system, we show that the spatial profile of the non-ideal term is in good agreement with the viscosity term in the Ohm's law, when the EM turbulence is dominant in the EDR. The paper is organized as follows. In Sec. II, the simulation model and initial setup employed in this study are described. Section III shows the results of large-scale 3D PIC simulation and derives the analytical model of the non-ideal electric field based on the flow shear instabilities. The present work is summarized in Sec. IV with discussions on the difference from the EDR without significant turbulence as seen in the 2D simulations.
II. SIMULATION MODEL
We perform a 3D electromagnetic PIC simulation with the adaptive mesh refinement (AMR).25 The AMR technique enables efficient high-resolution simulation of multiscale kinetic processes by subdividing computational cells locally in space and dynamically in time. The initial configuration is provided by a Harris-type current sheet17 with the magnetic field and the number density , where δ is the half width of the current sheet and the subscript b indicates the background parameter. The present study employs and with λi being the ion inertia length based on n0. Note that the background plasma is removed from the current sheet center in order to avoid the stream instabilities due to the drifting current sheet plasma and stationary background plasma. Although there appears a weak pressure imbalance due to the background density profile, the equilibrium is quickly established without any significant modification to the current sheet structure. The mesh refinement is performed in accordance with specific criteria provided by the local electron Debye length λDe and the electron flow velocity Vey defined at the center of each computational cell. Each cell is subdivided when or are satisfied, otherwise the child cells are removed, if any, where ΔL is the size of the cubic computational cells and VA is the Alfvén velocity based on B0 and n0.
The boundary conditions of the simulation domain are periodic in the x and y directions and the conducting wall in the z direction. The ion-to-electron mass ratio and velocity of light are and , respectively, corresponding to , where ωpe and ωce are the electron plasma frequency and cyclotron frequency based on n0 and B0, respectively. The temperature ratios are , and , so that the background ions are colder than the plasma sheet ions. The system size is , which is entirely covered by base-level cells (coarsest cells) with and can be subdivided locally into finer cells up to the dynamic range level (finest level) with . The normalizations of the physical quantities are performed using mi for mass, e for charge, λi for length, and VA for velocity, unless otherwise noted. For comparison, a 2D simulation in the xz plane is carried out with the same initial setup.
III. RESULTS
A. Non-ideal electric field
Magnetic reconnection is initiated with a small perturbation to the magnetic field components Bx and Bz (Ref. 26), which produces an x-line at the center of the xz plane uniformly along the y axis. Once the simulation has been started, the initial density and current profiles are quickly modified to adjust the field perturbation. As a result, a thinner current sheet is formed near the center of the xz plane, so that the tearing instability selectively develops therein. After a long slow development of the tearing mode, a fast reconnection is spontaneously achieved with the rate of , where ER is evaluated by Ey at the instantaneous x-line on the xz plane averaged over the y domain and is normalized to the instantaneous upstream values. As discussed in the previous studies,15,16 the electron current layer formed around the x-line is unstable to the flow shear instabilities, leading to the EM turbulence around the EDR. The turbulence intensity is significantly enhanced due to the formation of magnetic islands that are intermittently ejected from the reconnection current layer. Figure 1(a) shows a 3D view of the turbulent reconnection layer during the fast reconnection. The magnetic islands are locally generated with a finite length in the y direction, which stimulates the turbulence activities in the 3D current layer. The present study considers a macroscopic description of the turbulent current layer in a form that is applicable to the fluid systems such as the MHD where the electron inertia and kinetic effects are not incorporated. For this purpose, we average the current sheet over the entire y domain and evaluate the collective effects of the electron-driven instabilities.
The contributions of the turbulence to the reconnection process are evaluated in the generalized Ohm's law averaged over the y domain (e.g., Refs. 12 and 27),
where with being the average quantity of a variable A, being the deviation from the average, and being the electron pressure tensor. Each term of Eq. (1) is normalized to the upstream values of the reconnection region; that is, is the Alfvén speed, Bu is the magnetic field magnitude, and nu is the plasma density. The difference between LHS and the first term on RHS of Eq. (1), i.e., , represents the non-ideal electric field that is responsible for the magnetic dissipation. The second and third terms on the RHS originate from the average electron inertia,4 while the last three terms are the contributions of the turbulence to the non-ideal electric field. The fourth term, , is caused by electrostatic turbulence, the fifth term, , arises from EM turbulence, and the last term, , results from the eddy viscosity due to turbulence.
The average current density in Fig. 1(b) shows a profile similar to that in a 2D reconnection simulation as discussed in Sec. IV, regardless of the intense fluctuation along the y axis. A thin current layer is generated around the dominant x-line located at , forming the EDR. The super-Alfvénic jets also appear near the downstream edges of the EDR as in the 2D reconnection. Although the vertical width of the EDR is increased due to the turbulence,24 the typical scale of the half width remains the electron scale, no more than , where is the upstream electron inertia length. Thus, the width of the EDR is still much smaller than that of the ion diffusion region. In contrast, the non-ideal electric field in Fig. 1(c) shows persistent negative values over the upstream edges of the EDR in addition to the positive peak near the x-line. These negative ENI regions have not been clearly seen in the 2D reconnection. We find that the ENI distribution is well consistent with that of the EM turbulence term around the EDR [Fig. 1(d)], indicating that the effective magnetic dissipation can be described in terms of the EM turbulence for the current sheet averaged over the y domain.
The more quantitative comparisons between the current density profile and the Ohm's law are performed in Fig. 2, showing the line profiles along the z axis across the x-line as manifested by the magnetic field profiles in Fig. 2(d). One can see that the reconnection electric field almost balances the EM turbulence term at the x-line and the other terms in the generalized Ohm's law (1) are negligibly small [Fig. 2(a)]. Although our AMR-PIC code employs a weak numerical diffusion to suppress grid-scale noise,25 the good consistency between the LHS (thick black curve) of Eq. (1) and the sum of the RHS (light black curve) indicates that the numerical effects are negligible, if any. The non-ideal electric field profile [black curve in Fig. 2(b)] is quantitatively in good agreement with the EM turbulence term (red curve) across the EDR, with clear negative dips around the edges of the EDR. It is interesting to notice that these dips do not appear in the current density profile [Fig. 2(c)]. Thus, the non-ideal field profile in the turbulent current layer cannot be described by the resistivity term in the form of , where the resistivity η is generally positive both in the theoretical9 and numerical6 models.
The formation of the negative dips is closely related to the turbulence activities in the reconnection current layer. Figure 3 shows the time evolutions of the ENI profile across the EDR and wave spectrum of at the location of the instantaneous average x-line together with the reconnection rate. The fast reconnection with starts at [Fig. 3(a)]. One can see that the positive peak of ENI is also increased in association with the fast reconnection [Fig. 3(b)]. However, in the beginning of the fast reconnection, the negative dips in ENI do not appear clearly near the edges of the EDR. Instead, these dips emerge at , which is almost coincident with significant enhancement of the turbulence activities around the x-line [Fig. 3(c)]. One can notice that the enhancement of the power spectrum density of initially occurs at in a low wavenumber range with and spreads to higher wavenumber as time goes on. This suggests that the turbulence is mainly driven by long wavelength modes corresponding to the CSSI and eKHI.
Figure 3(d) compares the z profiles of the EM turbulence term of the generalized Ohm's law before and after the turbulence intensity has been enhanced. When both CSSI and eKHI waves are at a low amplitude without significant cascading at , the EM turbulence term do not show clear negative dips (green curve). As the amplitude increases and the energy cascading proceeds, the dips become more dominant around the upstream edges of the EDR (blue and red curves), indicating close relationship between the dip formation and turbulence activities.
B. Analytical model
Both the eKHI and CSSI are excited through the electron flow shear across the localized electron current layer, i.e., the EDR.15,16 Regardless of the different growth rate, these modes can coexist, since the free energy is different. The CSSI is caused by the flow shear across the EDR along the inflow direction, while the eKHI is excited by that along the outflow direction. The EM turbulence term in the generalized Ohm's law (1) can be rewritten in the following form:
where is assumed for the EM turbulence. The first term on RHS is caused by the flapping motion of the current layer in the z direction, corresponding to the CSSI, while the second term results from the current sheet oscillation in the x direction due to the eKHI. The evolutions of and are given by Faraday's law
where and are assumed because of the long-wavelength properties of the eKHI and CSSI. Although these modes cascade to shorter-wavelength modes due to the nonlinear effects in the simulation [Fig. 3(c)], most wave power remains concentrated around the parent modes expected by the linear wave dispersions. Thus, the long-wavelength approximation is still valid in the turbulent current layer. Based on our previous study,15 may be approximated by the electron inertia term
The convection term on RHS can be neglected in a frame co-moving with the wave in the y direction. Since the other terms are not affected by the frame transformation in the y direction, we may evaluate for the eKHI and CSSI in the separate co-moving frames, so that the convection term can be omitted for each mode.
where τgx and τgz are the growth times of the CSSI and eKHI, respectively, and , and . In Eq. (5), the first derivatives of and in both the x and z directions are neglected near the x-line because of the kink-type modes where both and are expected to be the even functions with respect to the x-line location. Using Eq. (5) in Eq. (2) with an approximation of , we can rewrite the EM turbulence term in the following form:
where and are the phase lags between and and between and , respectively, which are not incorporated in Eq. (5). Taking into account in the localized current layer, we can evaluate the contribution of the EM turbulence to the non-ideal electric field by
where μx and μz are identified as coefficients of the effective viscosities28,29 caused by the turbulence and are given by
Note that each term in Eqs. (6) and (7) is normalized to the upstream values based on VAu, Bu, and nu. From the linear wave analyses,15,16 the growth rates of the eKHI and CSSI may be approximated by and , respectively, where γeKHI and γCSSI are the linear growth rates of the eKHI and CSSI, respectively, normalized to .
We first estimate the phase lags between and to calculate the effective viscosities in Eq. (7). The EM turbulence term in Eq. (1) can be decomposed into a series of the cross-spectrum by using the following relations:
where Ak denotes a Fourier component of A(y), is a complex conjugate of Ak, and represents a Fourier amplitude of Ak. Figure 4 shows the cross-spectra between and Bz [Figs. 4(a) and 4(c)] and between and Bx [Figs. 4(b) and 4(d)] taken at the location of the average x-line. From the spectrum amplitude in Figs. 4(a) and 4(b), one can find that the peak amplitude occurs at for the term and at for the term. The former peak is consistent with the eKHI with (Ref. 16), while the latter peak occurs very close to the typical wavenumber of the CSSI with (Ref. 15), where is the upstream ion inertia length based on at . The corresponding phase lags are taken at the wavenumbers of the peak spectrum amplitudes, that is, for the term and for the term.
Using these representative values for the phase lags in Eq. (7), we obtain a theoretical estimate of the non-ideal electric field due to the EM turbulence in Eq. (6). Figure 5 shows the theoretical profiles (red squares) along the z axis across the EDR, where , and are taken from the simulation. It is found that the theoretical estimate in Fig. 5(a) is in good agreement with the direct simulation result (black curve) in the entire area across the EDR. In particular, the negative dips around the edges of the EDR are clearly reproduced in the theoretical profile, indicating that the modeling based on the flow shear instabilities is appropriate. We also find that the major contribution to the EM turbulence term is from the term [Fig. 5(b)], which is mainly caused by the eKH waves, in the entire area. This is because of the large cross correlation between and . The term [Fig. 5(c)], mainly originated from the CSS waves, is relatively small, but is not negligible with the contribution around 20% of the total non-ideal field in the simulation. The discrepancy between the simulation and theory around in the term may be due to the rough estimate of . From the linear analysis of the CSSI, the peak growth rate is more accurately (Ref. 12), corresponding to , which results in decrease in the peak by a factor of two. Nevertheless, the good agreement in the spatial profile between the simulation and theory, including the negative dips, suggests that the description of the magnetic dissipation based on the viscosity, rather than the resistivity, is appropriate around the EDR, when the EM turbulence is dominant.
IV. SUMMARY AND DISCUSSIONS
The plasma fluid frameworks, such as the MHD, need to incorporate the kinetic effects in the form of a macroscopic model in order to drive magnetic reconnection. However, the appropriate physical model has not been established for collisionless plasma in a form that is applicable to the fluid equations. We have investigated the macroscopic description of the kinetic effects caused by the EM turbulence that is responsible for the magnetic dissipation in the EDR. Following our previous study,24 the present study has considered the spatial distribution of the non-ideal electric field across the EDR and derived a form that is applicable to the fluid systems.
A large-scale 3D PIC simulation shows that the non-ideal electric field almost balances the EM turbulence term in the generalized Ohm's law, Eq. (1), over a broad range across the reconnection layer, when the turbulence activities are significant. It is found that the non-ideal field profile has clear negative dips around the upstream edges of the EDR as well as a positive peak at the x-line. These dips cannot be described by the resistivity term in the form of , since the Jy profile does not show such dips. To achieve a better description of the spatial profile, we derived an analytical form of the EM turbulence term, based on the dispersion relations of the flow shear instabilities. The resultant analytical form suggests that the dissipation should be caused by the viscosity due to the electron momentum transport across the EDR rather than by the resistivity. We found that the analytical profile is in good agreement with the simulation result across the EDR including the negative dips, indicating that the description based on the viscosity terms is appropriate. The present result suggests that a fundamental modification is needed for the plasma fluid frameworks based on the resistivity to generate the magnetic dissipation. The flow shear instabilities in the present simulation are expected to be suppressed in the strong guide-field limit with (Ref. 12), where Bg is the guide-field intensity. Even in such a low βe regime, a flow shear instability may be excited, resulting in the viscous dissipation.13
The origin of the reconnection electric field in collisionless plasma has been extensively investigated for 2D magnetic reconnection, in which there are no plasma waves propagating in the out-of-plane direction. Speiser30 proposed a theoretical model that suggests a finite conductivity can be generated by the momentum transport of the electrons that were accelerated in the EDR and eventually ejected from the EDR. Meanwhile, Vasyliunas4 conjectured from the fluid point of view, that the out-of-plane electric field in the generalized Ohm's law should balance the electron pressure tensor term at the x-line, since the other terms vanish there. The importance of the pressure tensor term at the x-line was later verified in a number of 2D PIC simulations (e.g., Refs. 31–33). The origin of the pressure tensor components in terms of the particle motions was investigated by the test particle simulations34,35 and by analytical analyses.36,37 They have demonstrated that the electron motions about the x-line leading to the Speiser's conductivity indeed give rise to the pressure tensor term to balance the electric field in the generalized Ohm's law.
Figure 6 shows a typical reconnection current layer evolved in the 2D PIC simulation. The 2D current profile in Fig. 6(a) is very similar to that in 3D magnetic reconnection [Fig. 1(b)], that is, a thin current layer is locally generated around the dominant x-line, forming the EDR. Figure 6(b) shows the 1D profile of the generalized Ohm's law (1) across the x-line as manifested by the magnetic field profile [Fig. 6(e)], where the turbulence terms vanish in Eq. (1) with and for a variable A in the 2D simulation. One can see that the electron pressure tensor term balances the reconnection electric field at the x-line. Since the pressure tensor term originates from the electron momentum transport across the flow shear layer, the term is essentially identical to the viscosity term that is proportional to (Ref. 38). In fact, the profile of the pressure tensor term in Fig. 6(b) shows clear negative dips around the edges of the electron current layer, similar to that of the EM turbulence term in 3D reconnection [Fig. 2(a)]. Such the profile cannot be described by the resistivity term in the form of , since the Jy profile has no negative dips [Fig. 6(d)], and is more consistent with the viscosity term as derived for the turbulent current layer in 3D.
However, contrary to the turbulent current layer, the non-ideal electric field in Fig. 6(c) does not show clear negative dips around the upstream edges of the EDR. This is because the negative dips in the pressure tensor term have been almost canceled out due to the electron inertia term [magenta curve in Fig. 6(b)], as pointed out in Ref. 39. The profile in the 2D reconnection is similar to that in the 3D reconnection before the EM turbulence is intensified at [Figs. 3(b) and 3(c)]. Thus, in 2D current sheet without significant turbulence activities, the non-ideal electric field profile cannot be simply expressed in the form proportional to . Further investigation is needed to reveal more appropriate descriptions of the EDR for 2D reconnection in a form that is applicable to the plasma fluid equations.
Finally, it is worthy to note that, even in the 3D turbulent EDR, the local non-ideal field in the turbulence is expressed only by the electron pressure tensor and inertia terms, since there are no other terms in the local Ohm's law. By averaging the local current layer over the y domain, these terms are mostly smoothed out and the dominant collective effect is replaced by the EM turbulence term [Fig. 2(a)]. In other words, the 2D profile represented by Fig. 6(b) is embedded in the 3D turbulent layer whose collective effect is identified as the electron momentum transport across the EDR, leading to the viscous dissipation.
ACKNOWLEDGMENTS
The simulations were carried out by the K computer at the RIKEN Advanced Institute for Computational Science through the HPCI Research project (Nos. hp140129 and hp150123). The data analyses were partly performed on analysis servers at Center for Computational Astrophysics (CfCA), NAOJ. This research was supported by the National Natural Science Foundation of China under Grant Nos. 41874189 and 41821003 and the Natural Sciences and Engineering Research Council of Canada.
AUTHOR DECLARATIONS
Conflict of Interest
The authors have no conflicts to disclose.
Author Contributions
Keizo Fujimoto: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Funding acquisition (lead); Investigation (lead); Methodology (lead); Project administration (lead); Resources (lead); Software (lead); Supervision (lead); Validation (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead). Richard Sydora: Funding acquisition (supporting); Validation (equal); Writing – review & editing (supporting).
DATA AVAILABILITY
The data that support the findings of this study are available from the corresponding author upon reasonable request.