3D turbulent reconnection: Theory, tests, and astrophysical implications

Magnetic fields are ubiquitous in astrophysical systems and critically affect the dynamics and properties of magnetized plasmas over an extended range of scales. These scales are typically much larger than any relevant plasma scales, e.g., the ion inertial length. The magnetic field of the Earth's magnetosphere, the case-study of many in situ measurements, is an exception in the astrophysical context in terms of the involved scales. In what follows, we will discuss magnetic fields at much larger scales.

One textbook concept related to magnetic fields, widely used in astrophysical studies, is that magnetic fields in highly conducting plasmas are nearly perfectly frozen into the fluid and retain their topology for all time (Alfvén, 1942; Parker, 1979). This concept of flux-freezing is at the heart of many theories, e.g., the theory of star formation in the magnetized interstellar medium (ISM).

At the same time, there is ample evidence that magnetic fields in astrophysical systems do change their topology, with solar flares being a classical example of such a process (see, e.g., Parker, 1970; Lovelace, 1976; Priest and Forbes, 2002). The studies on the process enabling such a topology change, i.e., magnetic reconnection, were historically motivated by observations of the solar corona (Innes et al., 1997; Yokoyama and Shibata, 1995; Masuda et al., 1994), which also created a number of misconceptions. First of all, magnetic reconnection is frequently associated only with the processes of conversion of magnetic energy into heat. In addition, this understanding influenced attempts to find very peculiar conditions for flux conservation violation. For instance, in a fundamental study by Priest and Forbes (2002), a number of examples of magnetic configurations that produce fast reconnection are discussed. The attempts to accelerate the reconnection rate for special plasma conditions, e.g., in plasmas with very small collision rates (see, e.g., Shay et al., 1998; Drake, 2001, 2006; Daughton et al., 2006; Uzdensky and Kulsrud, 2006; Bhattacharjee et al., 2003; Zweibel and Yamada, 2009; Yamada et al., 2010), follow the same logic of seeking special circumstances where the magnetic flux-freezing is violated.

We feel that the above treatment of magnetic reconnection is not generic or widely applicable. Solar flares (Sturrock, 1966) are just one vivid example of reconnection activity. In reality, magnetic reconnection is a ubiquitous process taking place in various astrophysical environments, both collisionless and collisional (see Shibata and Magara, 2011). For instance, magnetic reconnection is a part of large-scale dynamo acting in stellar interiors (Parker, 1993; Ossendrijver, 2003). Magnetic reconnection is also required to make the theory of magnetohydrodynamic (MHD) turbulence (Goldreich and Sridhar, 1995) self-consistent (Jafari and Vishniac, 2019b).

The problem of magnetic reconnection is not limited to explaining its typically fast rates. In fact, it is often overlooked that the observations of solar activity indicate that the reconnection rates can significantly vary. This presents serious problems for theories that are based on relating the rate of magnetic reconnection to the peculiar properties of plasmas. These properties do not change, for instance, as the flux gets accumulated prior to a solar flare and gets annihilated during the flare. This may suggest that magnetic reconnection should have a sort of trigger.

It is also worth noting that for decades, a paradoxical situation existed with the treatment of astrophysical magnetic reconnection. On the one hand, the core of the reconnection community claimed numerous restrictions either on magnetic field configurations or the properties of reconnecting plasmas. On the other hand, the practitioners of numerical astrophysics used their MHD codes to simulate various astrophysical settings without worrying about the aforementioned restrictions. At the same time, the question of the actual effect of magnetic reconnection on the validity of numerical studies was somehow avoided.

The absence of an effective communication between the communities is surprising. For instance, for many years, it was considered that in order for reconnection to be fast, it must take place in collisionless environments. However, if collisionless reconnection were the only way to make reconnection rapid, then numerical MHD simulations of many astrophysical processes including those of the interstellar medium (ISM), which is collisional, are in error. Fortunately for numerical practitioners, the observations of the collisional solar photosphere indicate that the reconnection is fast in these environments as well (see Sec. VII A).

In this review, we provide evidence that, in fact, turbulence determines the rate of reconnection in realistic 3D astrophysical systems. The concept of 3D turbulent reconnection was introduced by (Lazarian and Vishniac, 1999, henceforth LV99). This model was followed by several subsequent theoretical studies, for instance, Eyink et al., (2011, henceforth ELV11), Eyink (2015), and Jafari et al. (2018). The LV99 theory has been supported by numerical simulations in nonrelativistic (see Kowal et al., 2009, 2012; Eyink et al., 2013; Beresnyak, 2017; Oishi et al., 2015; Kowal et al., 2017) and relativistic settings (see Takamoto and Lazarian, 2016; Takamoto, 2018). In addition, different pieces of observational evidence, discussed in the review, support the predictions of the 3D turbulent reconnection theory of LV99.

To understand the essence of the LV99 model, one can recall that the problem of magnetic reconnection in astrophysical systems can be viewed as the problem related to the scale disparity. Reconnection occurs on very large scales, while the dissipation processes take place at the smallest plasma scales, which are set by, e.g., resistivity. LV99 theory shows how turbulence, as a multiscale phenomenon, can indeed resolve the problem of scale disparity and other problems that plague the traditional reconnection models.

While the idea of turbulent reconnection is nearly coeval with the ideas of 2D collisionless Hall-MHD reconnection, the latter has been substituted more recently by the ideas of tearing/plasmoid reconnection (Shibata and Tanuma, 2001). In fact, 2D numerical work has already shown that magnetic reconnection can be accelerated significantly by the tearing instability (Loureiro et al., 2007; Lapenta, 2008; Daughton et al., 2009a, 2009b; Bhattacharjee et al., 2009; Cassak et al., 2009).

In a sense, this model has some similarities to the turbulent reconnection. Indeed, it departs from the Sweet-Parker Y-point reconnection but allows instabilities to evolve within the current sheet. We will discuss how in realistic 3D configurations this type of instability transfers to turbulent reconnection.

One should keep in mind that we deal with 3D turbulent reconnection on large scales at which an MHD-like approximation is applicable. In particular, we consider scales much larger than the ion Larmor radius or skin depth, where an ideal Ohm's law is generally agreed to be valid. At the same time, usual tests of reconnection are based on the in situ measurements of magnetospheric reconnection. This type of reconnection occurs on scales comparable to the ion inertial length, and therefore, it is atypical for the large-scale reconnection that occurs in most astrophysical systems. For the reconnection explored at larger scales, e.g., in the solar wind (Gosling, 2012; Lalescu et al., 2015), the measured properties of reconnection are consistent with the general expectations based on the LV99 theory (see Sec. VII B).

Occasionally, the process of magnetic reconnection is understood in a very narrow sense, e.g., magnetic reconnection necessarily in low-β plasmas with the magnetic energy much larger than the thermal energy. Here, we consider turbulent reconnection as a generic process, which is applicable to plasmas of arbitrary β and for arbitrary arrangements of reconnecting magnetic fluxes.

Our review is devoted to the astrophysical reconnection, although the magnetic reconnection is a fundamental physical process, and therefore, our conclusions can be applied also to the laboratory reconnection. We discuss the theory and observations of astrophysical turbulence in Sec. II. The LV99 model of 3D turbulent reconnection and its extensions are described in Secs. III and IV. There we explain how turbulent reconnection and flux-freezing violation in turbulent fluids are inter-related. The numerical tests of turbulent reconnection and the testing of flux freezing violation in turbulent fluids are provided in Sec. V. Testing with both MHD codes and Hall-MHD codes are discussed. The initial tests of selected ideas of turbulent reconnection with 3D kinetic simulations are presented in Sec. VI. The observational evidence supporting turbulent reconnection is briefly discussed in Sec. VII. The magnetic reconnection in a high-Prandtl number medium, reconnection in relativistic fluids, and reconnection in the presence of whistler turbulence are discussed in Sec. VIII. Astrophysical implications of reconnection for turbulent dynamo, star formation, as well as for particle acceleration and gamma ray bursts are covered in Sec. IX. In Sec. X, we discuss the differences between our model of turbulent reconnection and other ideas on the effects of turbulence on reconnection and resistivity. We explain the reasons for our disagreement with (i) the currently popular reconnection-mediated turbulence idea and (ii) the concepts of turbulent resistivity and turbulent ambipolar diffusion. There we also outline the problems of the mean field approach to study reconnection, compare turbulent and tearing reconnection, and explain why reconnection rates are not expected to have a “universal” value of 0.1 Alfvén velocity. We present our summary in Sec. XI.

Observations of the diffuse warm ISM reveal a Kolmogorov spectrum of electron density fluctuations (see Armstrong et al., 1995; Chepurnov and Lazarian, 2010). Similar spectral slopes of supersonic velocity fluctuations are measured using Doppler shifted spectral lines (see Lazarian, 2009, for a review) and through the in situ measurements of the solar wind fluctuations (Leamon et al., 1998). The evidence of turbulence being ubiquitous comes from the nonthermal broadening of spectral lines and measures obtained by other techniques (see Burkhart et al., 2010). This fact is not surprising as magnetized astrophysical plasmas generally have very large Reynolds numbers. Indeed, the length scales involved are large, while the motions of charged particles in the direction perpendicular to magnetic fields are of the order of the Larmor radius. Plasma flows at these high Reynolds numbers are prey to numerous linear and finite-amplitude instabilities. This induces turbulent motions readily.

The precise origin of the plasma turbulence differs from instance to instance. It is sometimes driven by an external energy source, such as supernova explosions in the ISM (Norman and Ferrara, 1996; Ferrière, 2001), merger events and active galactic nuclei (AGN) outflows in the intercluster medium (ICM) (Subramanian et al., 2006; Enßlin and Vogt, 2006; Chandran, 2005), and baroclinic forcing behind shock waves in interstellar clouds. In other cases, turbulence is spontaneous, with available energy released by a rich array of instabilities, such as magnetorotational instability (MRI) in accretion disks (Balbus and Hawley, 1998; Jafari and Vishniac, 2018b) and kink instability of twisted flux tubes in the solar corona (Galsgaard and Nordlund, 1997; Gerrard and Hood, 2003). In all these mentioned cases, turbulence is not driven by reconnection. Nevertheless, we would like to stress that the driving of turbulence through the energy release in the reconnection zone is important for the transfer from laminar to turbulent reconnection. Altogether, whatever its origin is, the signatures of plasma turbulence are seen throughout astrophysical media.

Indeed, observations show that turbulence is ubiquitous in all astrophysical plasmas. The spectrum of electron density fluctuations in the Milky Way is presented in Fig. 1, as one dramatic piece of evidence for cosmic turbulence. Similar spectra are reported in Leamon et al. (1998) and Bale et al. (2005) for solar wind, Padoan et al. (2006) for molecular clouds, and Vogt and Enßlin (2005) for the ICM. New techniques for studying turbulence, e.g., the Velocity Channel Analysis (VCA) and Velocity Coordinate Spectrum (VCS) techniques (Lazarian and Pogosyan, 2000, 2004, 2006), have provided important insights into the velocity spectra of turbulence in molecular clouds (see Padoan et al., 2006, 2010), Galactic, and extragalactic atomic hydrogen [Stanimirović and Lazarian, 2001; Chepurnov et al., 2010, 2015, see also the review by Lazarian (2009), where a compilation of velocity and density spectra obtained with contemporary HI and CO data are presented].

It is generally accepted that an MHD-like description of turbulence is applicable to the motions of magnetized plasmas at sufficiently large scales (e.g., see discussion in Schekochihin et al., 2009). Standard MHD or a multispecies variant is valid at length-scales much larger than the collisional mean-free-path of the plasma. However, even in a nearly collisionless but magnetized plasma, “kinetic MHD” or “kinetic reduced-MHD” equations are satisfied at scales much larger than the ion gyroradius (Kulsrud, 1983; Kunz et al., 2015). Since collisions in those plasmas are too infrequent to bring the particle distributions to Maxwellian, the pressure tensor is anisotropic and given not by a thermodynamic equation of state but instead by 1D kinetic equations along particle guiding line-centers. This subtlety is irrelevant for the incompressible shear-Alfvén wave modes, which nearly uncouple from the kinetic/compressible modes and satisfy the well-known “reduced MHD” (RMHD) equations [for relativistic turbulence (see Takamoto and Lazarian, 2017) the transfer of energy to fast modes is more important]. Shear-Alfvén waves are often energetically dominant, as in the solar wind where they contain about 90% of the total energy of plasma and fields (see Dobrowolny et al., 1980). For this reason, we shall generally assume that MHD is an accurate model at sufficiently large scales or at least a good working approximation within which to explore the effects of turbulence on large-scale reconnection.

We shall invoke below in our discussion of strong MHD turbulence the (Goldreich and Sridhar, 1995, henceforth GS95) theory of Alfvénic turbulence. We believe that the empirical evidence largely supports GS95 over its competitors and gives a good description of MHD turbulence (apart from the phenomenon of small-scale intermittency). However, we wish to avoid here any controversy regarding the “correct” phenomenology. It may be shown following LV99, Appendix D, that the predictions on turbulent reconnection are qualitatively independent of the particular phenomenology adopted and, in particular, do not depend sensitively upon the predicted spectral slopes. Therefore, our use of GS95 theory may be regarded by its persistent skeptics as a mere convenience since it is the simplest model that incorporates the small-scale anisotropy due to a local magnetic field.

MHD turbulence was traditionally viewed as the process dominated by the Alfvén wave interactions. Below, we follow this path to show how to obtain GS95 relations and go beyond them (see also Cho et al., 2003a).

For Alfvénic perturbations, the relative perturbations of velocities and magnetic fields are related in the following way:

where B_l is the perturbation of the magnetic field B at the scale l, B_L is the perturbation of the magnetic field at the injection scale L, while u_l is the velocity fluctuation at the scale l in the turbulent flow with energy injected with velocity u_L.

To understand Alfvén wave damping, it is advantageous to present MHD turbulence as a superposition of colliding Alfvén wave packets. Consider such packets with parallel scales $l||$ and perpendicular scales $l⊥$⁠. The system of reference with respect to the local magnetic field is what is relevant when dealing with MHD turbulence. In the rest of the discussion, we operate with $l‖$ and $l⊥$ which are defined in the frame of wavepackets moving along the local magnetic field. The wave packet collision induces a change

where the first term is the change of the energy of a packet as it interacts with the oppositely moving wave packet. The time of the interaction should be identified with the time of the passage of the given wave packet through the oppositely moving packet of the size $l||$⁠. Thus, the interaction time is given by $Δt∼l||/VA$⁠.

The characteristic rate of the cascade is due to the change in the structure of the oppositely moving wave packet. The latter has the rate $ul/l⊥$⁠. As a result, Eq. (2) gives

The fractional energy change per collision is $ΔE/E$⁠, which measures the strength of the nonlinear interaction

In Eq. (4), f measures the ratio of the rate of shearing of the wave packet $ul/l⊥$ to the rate of the propagation of the wave packet $VA/l||$⁠. Two cases can be identified. For $f≪1$⁠, the shearing rate is much smaller than the propagation rate, and the cascade process is a random walk. This means that

which is the number of steps required for the O(1) cascade of energy, and therefore, the cascade time is

For $ℵ≫1$⁠, the turbulence corresponds to the “weak MHD turbulence.” Naturally, $ℵ$ cannot be smaller than unity. The limiting case of $ℵ≈1$ corresponds to “strong MHD turbulence.”

Consider first the case of weak MHD turbulence, which requires $VA≪uL$ and very high-frequency oscillations. This cascade thus proceeds by resonant 3-wave Alfvénic interactions that satisfy $ω1+ω2=ω3$ in addition to $k1+k2=k3,$ where $k$ and ω denote the wavevector and angular frequency, respectively. Because of the dispersion relation $ω=sVAk‖$ with $s=±1$ for shear-Alfvén waves traveling parallel and antiparallel to magnetic fields and because only counterpropagating waves interact, one finds that $k2,||=0$ and $k1,||=k3,||.$ Thus, resonantly interacting wavepackets have unchanging $k||$⁠, but $k⊥$ increases by distortions through collisions. The decrease in $l⊥$ while $l||$ does not change signifies the increase in the energy change per collision. This eventually makes $ℵ$ of the order of unity and the cascade transitions to strong MHD turbulence. In this case, one gets

which signifies the cascade time being equal to the wave period $∼Δt$⁠. It was hypothesized by Goldreich and Sridhar (1995) that the further decrease in $l⊥$ entails the corresponding decrease in $l||$ to keep Eq. (7) satisfied, which is the “critical balance condition.” The ability to change $l||$ by the nonresonant strong cascade means that the frequencies of interacting waves increase. However, this increase has a natural limit because a too high frequency would restore the resonance condition $ω1+ω2=ω3$⁠, and this would shut off the parallel cascade. Thus, the strong turbulent cascade has a tendency to operate at marginal resonance, where the critical balance relation (7) is satisfied.

The cascade of turbulent energy satisfies the relation (Batchelor, 1953)

which for the hydrodynamic cascade provides

where the relation $tcas≈l/ul$ is used.

For the weak cascade $ℵ≫1$

where Eqs. (8) and (6) are used (see LV99). The isotropic turbulence injection at scale $ℓ‖≃ℓ⊥≃L$ results in the third relation in Eq. (10). Taking into account that $l||$ is constant for weak turbulence, it is easy to see that Eq. (10) provides

which is different from the hydrodynamic $∼l1/3$ scaling [using the relation $kE(k)∼uk2$ it is easy to show that the spectrum of weak turbulence is $Ek,weak∼k⊥−2$ (LV99; Galtier et al., 2000)].

For “sub-Alfvénic turbulence” injected isotropically at length scale L with u_L < V_A, the transition to the strong regime, i.e., $ℵ≈1,$ happens at the scale (LV99)

Therefore, weak turbulence has a limited, i.e., $[L,LMA2]$⁠, inertial range and transits at l_trans into strong turbulence. The velocity at the transition follows from the $ℵ≈1$ condition given by Eqs. (5) and (4)

The scaling for the strong turbulence range at $l<ltrans$ can be easily obtained. Indeed, the turbulence cascades over one wave period, which according to Eq. (23), are equal to $l⊥/ul$⁠. Substituting the latter in Eq. (8), one gets

The latter is analogous to the hydrodynamic Kolmogorov cascade in the direction perpendicular to the local direction of the magnetic field. This strong MHD turbulence cascade starts at l_trans, and it has the injection velocity given by Eq. (13). This provides the another way to get the scaling relations for velocity in strong Alfvénic turbulence (LV99)

which can be rewritten in terms of the injection velocity u_L as

Substituting this in Eq. (7), we get the relationship between the parallel and perpendicular scales of the eddies (LV99)

The relations Eqs. (17) and (15) reduce to the GS95 scaling for trans-Alfvénic turbulence if $MA≡1$⁠.

In the opposite case of “super-Alfvénic turbulence,” i.e., for $uL>VA,$ the scales close to the injection scale are essentially hydrodynamic as the influence of magnetic forces is minimal. Thus, the velocity is Kolmogorov

The cascade changes its nature at the scale

for which the turbulent velocity is equal to the Alfvén velocity (Lazarian and Pogosyan, 2016). The cascade rate for $l<lA$ is

This case likewise reduces to trans-Alfvénic turbulence but with l_A acting as the injection scale. At scales $ℓ<lA$

and

The wave description of strong MHD turbulence is very productive in the case when the turbulence is imbalanced, i.e., the flow of wave energy from one direction along field lines significantly exceeds the flow from the other direction. An example of imbalanced MHD theory that agrees with numerical simulation is the theory by Beresnyak and Lazarian (2008) (see also Beresnyak and Lazarian, 2010).

In our above discussion, we revealed the analogy between the Alfvénic turbulence and Kolmogorov turbulence. This analogy can be extended in the case of strong Alfvénic turbulence. In fact, the GS95 relations can be naturally obtained if one assumes that such turbulence presents a collection of eddies that are mixing the conducting fluid perpendicular to the direction of the magnetic field present at the location of the eddies. This description can be justified in the presence of fast reconnection that we are describing in the review, and this supports our claim that the MHD turbulence and turbulent magnetic reconnection are intrinsically connected.

As it was described in LV99, the reconnection of magnetic fields associated with an eddy takes place within the eddy turnover time (see also Sec. III). In this situation, the random driving preferentially excites the motions that induce minimal bending of magnetic field lines, i.e., the hydrodynamic type motions of eddies that mix magnetized fluid perpendicular to the direction of the magnetic field. This is the picture of MHD turbulence that follows from the LV99 description. The key concept in the eddy description is that magnetic perturbations are measured in terms of the local direction of the magnetic field, not the mean field. This concept is missing in GS95, and this causes a lot of confusion with different authors attempting to test the GS95 scale-dependent anisotropy in the frame of the mean magnetic field. In fact, the latter anisotropy only exists in the local frame of the eddies, while the anisotropy in the reference frame of the mean magnetic field is scale-independent and is determined by the anisotropy of the largest eddies (see Cho et al., 2003a).

Incidentally, within the description suggested in LV99, it is clear that the gradients of both velocities and magnetic field are expected to be perpendicular to the local magnetic field direction. Therefore, by tracing the directions of the gradients, one can trace the magnetic field direction in turbulent media. This idea is at the foundations of the new technique of magnetic field tracing that employs gradients of intensities within spectroscopic channel maps (Lazarian and Yuen, 2018b), gradients of synchrotron intensity (Lazarian et al., 2018), and gradients of synchrotron polarization (Lazarian and Yuen, 2018a). The magnetic field tracing by the new technique provides results in good agreement with polarization data (see Hu et al., 2019), which is an indirect way of confirming both the ubiquitous presence of the MHD turbulence in the interstellar medium and a successful confirmation of the theoretical predictions.

Adopting the picture of eddies perpendicular to the local magnetic field of the eddies, one can claim that in terms of perpendicular motions, we expect to have the Komogorov picture, i.e., $uℓ∼l⊥1/3$⁠. The critical balance relations in the picture of the eddies also follow naturally. Indeed, as the eddy provides the mixing of magnetic field lines in the perpendicular direction over the time scale $l⊥/uℓ$⁠, it sends an Alfvén wave with the same period. This is exactly

where $l||$ and $l⊥$ are eddy axes parallel and perpendicular to the “local” direction of the magnetic field and the magnetic field of the eddy. Note that this requires a more sophisticated way to analyze numerical simulations. For instance, a way of determining the local magnetic field in numerical simulations using structure functions is discussed in the study by Cho and Vishniac (2000). In fact, the notion of the local magnetic field is the essential part of the modern understanding of Alfvénic turbulence, and it was added to the MHD turbulence theory by the later studies (LV99, Cho and Vishniac, 2000; Maron and Goldreich, 2001). This notion is natural within the picture of LV99 turbulent reconnection. Indeed, if one considers the aforementioned eddy picture, the use of the local magnetic field means that the eddies are only affected by the magnetic field around them and do not know anything about the mean magnetic field obtained by averaging over the entire turbulent volume.

It is assumed in GS95 that the energy is injected at the scale L with the injection velocity equal to the Alfvén velocity. If the injection velocity u_L is smaller than $VA,$ then the eddies at a small scale get more elongated. For $MA<1$⁠, the eddy description is valid only for the strong turbulence, i.e., starting with the scales less than $ltrans=LMA2$⁠. Taking l_trans as the injection scale of the turbulence and the velocity at this scale $VAMA2$ [see Eq. (13)] as the injection velocity, one can repeat the arguments about the correspondence of the time scales of parallel and perpendicular motions to obtain the relations given by Eqs. (15) and (17).

The compressible MHD turbulence can be approximated by the superposition of 3 cascades of basic MHD modes, namely, Alfvénic, slow, and fast. The original idea of the decomposition of low amplitude MHD perturbations into the modes can be traced to (Dobrowolny et al., 1980), while the practical realization for the actual trans-Alfvénic turbulence was first performed in (Cho and Lazarian, 2002). We, however, will concentrate on the Alfvénic part of the MHD turbulent cascade, as this component is most important for the theory of reconnection that we consider in the review. In nonrelativistic MHD turbulence, the backreaction of slow and fast magnetosonic modes (Cho and Lazarian, 2002, 2003; Kowal and Lazarian, 2010) on the Alfvénic cascade is insignificant (Cho and Lazarian, 2002; Goldreich and Sridhar, 1995; Lithwick and Goldreich, 2001). Therefore, in nonrelativistic MHD turbulence, the above relations obtained for Alfvénic turbulence are not much affected by the presence of compressible modes.

In the rest of this section, we discuss how the partial ionization and relativistic effects change the picture of Alfvénic turbulence.

In partially ionized astrophysical plasmas in, e.g., the early universe, cold interstellar phases, protoplanetary disks, and solar atmosphere, the MHD turbulence is subjected to the damping effects due to the presence of neutrals. Since we are interested in the turbulent reconnection associated with Alfvénic turbulence, here we focus on the damping of Alfvénic turbulence.

The damping effects in a partially ionized medium mainly arise from the collisional friction between ions and neutrals, i.e., ion-neutral collisional damping (IN) and the viscosity in neutral fluid, i.e., neutral viscous damping (NV) (Lithwick and Goldreich, 2001; Lazarian et al., 2004). Their relative importance depends on the turbulence and plasma parameters (Xu et al., 2015, 2016; Xu and Lazarian, 2017b). In super-Alfvénic turbulence with the turbulent energy larger than the magnetic energy, IN dominates over NV if the condition

is satisfied, where ξ_n is the neutral fraction, ν_ni is the neutral-ion collisional frequency, ν_n is the neutral viscosity, L is the driving scale of turbulence, and u_L is the turbulent velocity at L. In sub-Alfvénic turbulence with the turbulent energy smaller than the magnetic energy, under the condition

IN is more important than NV. Here, $MA=uL/VA$ is the Alfvén Mach number and V_A is the Alfvén speed. An illustration of the above condition in sub-Alfvénic turbulence can be found in Fig. 2, where the typical parameters of the interstellar turbulence are used (Xu and Lazarian, 2017b). Xu et al. (2016) found that in the warm neutral medium of the ISM and the solar chromosphere, NV plays a significant role in damping MHD turbulence.

When the damping effects are significant, MHD turbulence is suppressed. The damping scale where the MHD turbulence is dissipated can be determined by comparing the turbulent energy cascade rate with the damping rate. Different from the damping scale of linear MHD waves, the damping scale of MHD turbulence has parallel and perpendicular components, $kdam,||$ and $kdam,⊥$⁠, due to the anisotropy of MHD turbulence. In the frame of local magnetic field, i.e., the magnetic field averaged over the length scale of interest, magnetic fields are mixed by turbulent motions in the perpendicular direction, and the perturbed magnetic fields have wavelike motions in the parallel direction. For strong MHD turbulence, the lifetime of waves is equal to the eddy turnover time. As the turbulence anisotropy increases toward small length scales, when the damping scale of MHD turbulence is much smaller than L, it can be approximated by its perpendicular component in terms of magnitude.

Under the damping effects of both IN and NV, a general expression of the damping scale is (Xu and Lazarian, 2017b)

for super-Alfvénic turbulence and

for sub-Alfvénic turbulence, where V_Ai is the Alfvén speed of the ionized gas, ν_in is the ion-neutral collisional frequency, and $lA=LMA−3$ in super-Alfvénic turbulence.

When IN is the dominant damping effect, the damping is related to the coupling state between ions and neutrals. The most severe damping arises in the weakly coupled regime in a weakly ionized medium, where neutrals decouple from ions, but ions are still coupled with neutrals due to their frequent collisions with surrounding neutrals. In this situation, the (perpendicular) damping scale for super- and sub-Alfvénic turbulence can be simplified as (Xu and Lazarian, 2017b)

and

Its relation to the neutral-ion decoupling scale is

for both super- and sub-Alfvénic turbulence. The difference between the decoupling and damping scales has been confirmed by two-fluid MHD simulations (Burkhart et al., 2015). The above relation shows that the decoupling of neutrals from ions plays a key role in the IN damping of MHD turbulence.

After neutrals decouple from ions, the hydrodynamic turbulence carried by neutrals is only subject to viscous damping, and its NV damping scale is smaller than the IN damping scale of MHD turbulence in ions. This differential damping of neutrals and ions can account for the observed linewidth difference of the coexistent neutrals and ions in molecular clouds (Li and Houde, 2008; Xu et al., 2015).

When NV dominates over IN, the damping scale is approximately (Xu and Lazarian, 2017b)

for super-Alfvénic turbulence and

for sub-Alfvénic turbulence. MHD turbulence is damped at the NV damping scale. On the length scales smaller than the NV damping scale and larger than the IN damping scale, magnetic fluctuations persist, and the magnetic energy has a spectral form $k−1$⁠. This “new regime of MHD turbulence” in the subviscous range was first analytically predicted by Lazarian et al. (2004) and further confirmed by MHD simulations (Cho et al., 2002a, 2003b). The $k−1$ spectrum of magnetic energy in the subviscous range has also been observed in dynamo simulations in the nonlinear stage of dynamo (Haugen et al., 2004; Schekochihin et al., 2004).

Due to its numerical and theoretical simplicity, MHD turbulence in the relativistic force-free regime has been studied first. Relativistic force-free formalism can be used for a system, such as the magnetosphere of a pulsar or a black hole, in which the magnetic energy density is much larger than that of matter.

Scaling relations for relativistic Alfvénic MHD turbulence were first derived by Thompson and Blaes (1998). The predictions of these theories in terms of the spectral slope and anisotropy coincide with those in GS95 theory. These relations were first numerically tested by Cho (2005) who performed a numerical simulation of a decaying relativistic force-free MHD turbulence with a numerical resolution of 512³ and calculated the energy spectrum and anisotropy of eddy structures. After t = 3, the energy spectrum with $E∼k5/3$ was shown to decrease in amplitude without changing its slope. The anisotropy was also measured and was shown to be consistent with the GS95 expectations, i.e., $k||∼k⊥2/3$⁠. Thus, in terms of Alfvénic turbulence, the correspondence between the nonrelativistic theory and its relativistic counterpart is excellent. Interestingly enough, the relativistic simulations of imbalanced Alfvénic turbulence, i.e., the turbulence with the excess of the energy flux moving in one direction, performed by Cho and Lazarian (2014) provided results in agreement with the theory of nonrelativistic imbalanced turbulence in Beresnyak and Lazarian (2008).

The situation is somewhat different in the case of relativistic compressible theory. The initial studies of relativistic MHD turbulence (see Zhang et al., 2009; Inoue et al., 2011; Beckwith and Stone, 2011; Zrake and MacFadyen, 2012, 2013; Garrison and Nguyen, 2015) have not revealed much differences between the nonrelativistic turbulence and relativistic turbulence. However, more recent studies by Takamoto and Lazarian (2016) and (2017) provided the decomposition of the turbulent motions into Alfvén, slow, and fast modes. The comparison of the obtained results with the analogous decomposition into modes in the nonrelativistic case in the study by Cho and Lazarian (2002) and (2003) and Kowal and Lazarian (2010) revealed significant differences in terms of the coupling of the compressible and incompressible motions in the two cases. These differences are illustrated in Fig. 3 where the transfer of energy between Alfvén and fast modes is shown in the relativistic case. Note that the corresponding transfer of energy is suppressed for nonrelativistic compressible turbulence. In terms of anisotropy, this energy transfer makes fast modes more similar to Alfvén ones. We find that this effect of energy transfer is important for describing turbulent relativistic reconnection.

For our discussion of turbulent reconnection, it is important that turbulence is ubiquitous in astrophysical environments and that Alfvénic modes usually represent the dominant component of the MHD turbulent cascade. These modes are present in collisional and collisionless environments. The turbulence in most astrophysical settings has the sources different from the magnetic reconnection. However, it is difficult to understand the properties of turbulent fluids if magnetic reconnection does not happen on the dynamical time of the fluid motions.

The partial ionization of plasmas changes the properties of MHD turbulence increasing the scale at which MHD turbulence is damped. The relativistic MHD turbulence and nonrelativistic MHD turbulences have significant similarities, although they show a difference in the coupling of fast and Alfvén modes. A more detailed description of the theory and implications of MHD turbulence can be found in a book by Beresnyak and Lazarian (2019).

The 3D model of turbulent reconnection proposed in LV99 provides a natural generalization of the classical Sweet-Parker model (Parker, 1957; Sweet, 1958). Thus, we remind the reader first the basic properties of this model.

In the Sweet-Parker model, two regions with uniform “laminar” magnetic fields are separated by a thin current sheet. The magnetic fields are frozen in through the two regions, and the violation of the flux freezing takes place over a thin slot δ, with the latter being determined by Ohmic diffusion. Thus, the speed of reconnection is determined roughly by the resistivity divided by the sheet thickness, i.e.,

Obviously, the smaller the Δ value, the larger the reconnection speed. Small Δ, however, prevents plasmas from leaving the reconnection region.

The Sweet-Parker model deals with the “steady state reconnection,” and therefore, the plasma in the diffusion region must be ejected from the edge of the current sheet at the Alfvén speed, V_A. As a result, the reconnection speed is limited by the conservation of the mass condition

where L_x is the lateral extend of the current sheet.

The conditions imposed by Eqs. (33) and (34) act in the opposing directions as far as the thickness of the reconnection region Δ is concerned. As a result, the Sweet-Parker reconnection cannot be fast as it faces contradictory requirements. Therefore, the compromise speed of reconnection is the Sweet-Parker reconnection speed that is reduced from the Alfvén speed by the square root of the Lundquist number, $S≡LxVA/η$⁠, i.e.,

The Lundquist number in astrophysical conditions is enormous, and it can be 10¹⁶ or larger. In this situation, the Sweet-Parker reconnection rate is absolutely negligible.

Due to the disparity between L_x and Δ, the Sweet-Parker model cannot explain most of the observed astrophysical reconnection events, e.g., solar flares. A way to deal with this problem was suggested by Petschek (1964), who speculated that in special conditions, the reconnection configurations may have magnetic field lines converging to the reconnection zone at a sharp angle making L_x and Δ comparable. The corresponding configuration of the magnetic field is known as X-point reconnection as opposed to the Sweet-Parker extended current sheet Y-type reconnection. The research exploring X-point reconnection was the dominant trend in the attempts of obtaining fast reconnection for decades. However, the major problem of such configurations is to preserve the X-point configuration over astrophysically large scales.

The direction taken by LV99 was to modify the Y-type reconnection in order to make Δ macroscopic and even astrophysically large. Naturally, this is impossible with plasma effects. However, LV99 proposed that turbulence can induce magnetic field wandering that can do the job. The corresponding model of magnetic reconnection is illustrated in Fig. 4.

Similar to the Sweet-Parker model, the turbulent model developed in LV99 deals with a generic configuration, which arises naturally as magnetic flux tubes crossing each other make their way one through another. The outflow thickness Δ is determined by the large-scale magnetic field wandering. This wandering depends on the level of turbulence. As a result, the rate of reconnection in the LV99 model is determined by the level of MHD turbulence. The LV99 reconnection rate is fast, i.e., its rate does not depend on the fluid/plasma resistivity. However, the rates of turbulent reconnection can change significantly depending on the level of turbulence. This is an important prediction of LV99 theory which makes it different from the theories that, for instance, appeal to plasma effects for describing fast reconnection.

To quantify the reconnection rate in the LV99 model, one should use the scaling relations for Alfvénic turbulence from Sec. II. In fact, it is reconnection at low Alfvén Mach numbers which corresponds to the classical representation of magnetic reconnection flux tubes. This in fact required LV99 to generalize the GS95 MHD theory for arbitrary Alfvén Mach numbers.

As the first step consider the problem of magnetic field wandering, as this magnetic wandering determines the thickness of outflow Δ. We may note that magnetic wandering determined in LV99 is not only important for magnetic reconnection. In fact, it is a fundamental property of the turbulent magnetic field, and it is important for different problems, e.g., problems of heat transfer (Narayan and Medvedev, 2001, Lazarian, 2006) and cosmic ray propagation (Yan and Lazarian, 2008, Lazarian and Yan, 2014).

Using the MHD turbulence relations, one can describe the diffusion of magnetic field lines. Indeed, a bundle of field lines confined within a region of width y at some particular point spreads out perpendicular to the mean magnetic field direction as magnetic field lines are traced in either direction. The rate of field line spread is given by (see more in Sec. III B 3)

where $λ||−1≈ℓ||−1, ℓ||$ is the parallel scale. The corresponding transversal scale, $ℓ⊥$⁠, is $∼⟨y2⟩1/2$⁠, and x is the distance measured along the magnetic field direction. As a result, using Eq. (21), one gets

where the substitution $⟨y2⟩1/2$ for $ℓ⊥$ is used. This expression for the magnetic field dispersion is applicable only when y is small enough for us to use the strong turbulence scaling relations, i.e., when $⟨y2⟩<Li2(uL/VA)4$⁠.

For scales less than L_i, it follows from Eq. (37) that the mean squared magnetic field wandering is equal to

which provides a superdiffusive dependence of $⟨y2⟩∼s3$⁠, which was one of the major findings in LV99 as far as the properties of the magnetic field in turbulence are concerned. Due to the above scaling, the probability of the wandering magnetic field line to return back to the reconnection layer is very small.

The dependence of $⟨y2⟩∼MA4$ is significant as it shows that if the turbulent kinetic energy is a small portion of magnetic energy of the fluid, the magnetic field lines are nearly straight. This provides a smooth transition to the slow reconnection in the absence of turbulence and, as we discuss later, to the flux freezing being a fair approximation to the low M_A fluids.

For describing magnetic field wandering over scales larger than the turbulence injection scale, i.e., $s≫Li$⁠, magnetic field wandering obeys the usual random walk scaling. The corresponding random walk step is $s/Li$ with the mean squared displacement per step equal to $Li2(uL/VA)4$⁠. As a result

which provides the other limit for the magnetic field wandering. While the diffusive regime of magnetic field lines given by Eq. (39) was well known in the cosmic ray literature (see Schlickeiser, 2002), the superdiffusive one was introduced in LV99 [using the relation $kE(k)∼uk2$ it is easy to show that the spectrum of weak turbulence is $Ek,weak∼k⊥−2$ (LV99; Galtier et al., 2000)].

Using Eqs. (38) and (39), one can describe the thickness of the outflow for the reconnection of the antiparallel magnetic fields in turbulent media, i.e., in the setting generalizing Sweet-Parker reconnection without the guide field for turbulent media. Measuring s along the horizontal x-axis along the magnetic field reversal and combining Eqs. (38) and (39), it is possible to derive the thickness of the outflow Δ in the aforementioned two regimes. Then, using the mass conservation Eq. (52), it is straightforward to obtain

where $VAMA2$ is proportional to the turbulent eddy speed. As we discussed earlier, the obtained reconnection rate varies depending on Alfvén Mach number M_A and for $MA∼1$ can represent a large fraction of the Alfvén speed in the case when L_i and L_x are not too comparable.

Similar to the Sweet-Parker picture, the contribution of the guide field does not change the dynamics of magnetic reconnection. Figure 4 presents a XY plane of the cross section of the reconnection region perpendicular to the guide field B_z shared by two fluxes. The guide field should be ejected together with plasma from the reconnection region. Indeed, both the random velocity and magnetic field can be decomposed into components in the direction of the guide field and perpendicular to it. The field wondering in the direction perpendicular to the guide field is affected by the velocity perturbations in the same direction. The guide field does not oppose to the corresponding bending of the magnetic field. The outflow velocity is, however, determined only by the Alfvén velocity corresponding to the perpendicular component of the magnetic field. In fact, for the sake of clarity, it could be right to change M_A for $MA,XY$⁠, where the latter denotes the Alfvén Mach number corresponding to the Alfvén velocity in the plane perpendicular to the guide field. This notation of the Alfvén Mach number is not common, and we will not use it therefore.

It was communicated to us by Andrey Beresnyak that the independence of the reconnection rate on the guide field can be also justified from the point of the Reduced MHD (RMHD) approach (see Beresnyak, 2012). Formally, the RMHD is applicable for describing incompressible fluids when the perturbations are much smaller than V_A. However, in practical cases, RMHD is well applicable for many collisionless systems, e.g., solar wind (see Schekochihin et al., 2009), and its practical validity is good even for trans-Alfvénic turbulence (Beresnyak et al., 2009). Within the RMHD, the dynamics of turbulent motions is not changed if the ratio $B0/l||$ stays the same. If no constraints are imposed on the parallel scale of perturbations $l||$⁠, they increase with the increase in the guide field B₀ as was demonstrated in numerical simulations by Beresnyak (2017). The dynamics in terms of perpendicular perturbations $l⊥$ and field wondering that these perturbations entailed increasing the reconnection outflow stayed the same.

For distances along the magnetic field less than the injection scale L_i using the definition of the cascading rate for strong MHD turbulence given by Eq. (14), we can rewrite the expression for magnetic field wandering [see Eq. (38)] as

from which it is clear that the increase in V_A can be naturally compensated by the increase in the parallel scale along the magnetic field.

For the scales larger than the injection scale, Eq. (39) can be rewritten as

where the expression for the cascading of the weak MHD cascade given by Eq. (10) is used. It is also obvious that the simultaneous increase in the parallel scale of both injection and turbulent motions and the Alfvén velocity does not change the magnetic field wandering.

For the formal applicability of RMHD, one should assume that the guide field B_z is much stronger than the reconnecting magnetic field B_x. The relevant setting corresponds to the reconnection of two flux tubes intersecting at a small angle α. The corresponding setting was discussed in LV99 for the reconnection induced by individual eddies within MHD turbulence, and it was shown there that the reconnection happens within one eddy turnover time. Here, we consider a more generic situation of reconnecting flux tubes with an arbitrary turbulence driving.

In the above setting, $Bx≈Bz sin θ/2$⁠, and the perturbations at the injection scale L_i result in the perturbations in the XY plane with the scale $Li,XY∼Liθ/2$⁠. Similarly, the path x measured along the magnetic flux tubes projects into the path $xXY≈xθ/2$ in the XY plane. If one considers the dynamics of magnetic field lines in the XY-plane, the speed of propagation of Alfvénic perturbations along the x-axis is $VAθ/2$⁠, and therefore, in this approximation, it is obvious that the growth of the outflow region Δ as seen in the XY-plane does not depend on the angle θ between reconnecting fluxes, i.e., magnetic reconnection proceeds independently of the value of guide field but depends only on the B_x field and the properties of turbulence in the XY-plane.

In numerical experiments described in Sec. V, the properties of driving and the value of the reconnecting component of the magnetic field were fixed in the XY-plane. Therefore, the magnetic field wandering for such settings can be described by Eqs. (41) and (42) with V_A arising only from B_x and L_i being the injection scale of turbulence in the XY-plane. Correspondingly, Eq. (40) should use $MA≈δBx/Bx$⁠, where $δBx$ is a turbulent perturbation of B_x at the injection scale L_i.

The LV99 theory of turbulent reconnection has provided a number of predictions very different from the contemporary theories of fast magnetic reconnection. First of all, the turbulent reconnection was identified as a generic process that takes place everywhere in the magnetized conducting fluid. A violation of flux freezing implicitly follows from that. Second, LV99 does not prescribe the given reconnection rate but provides an expression that shows that the speed of reconnection varies depending on the level of turbulence. This can explain the variability of the energy release in reconnection events observed in, e.g., solar flares. Indeed, many of the observation processes require that the rate of reconnection changes. In LV99 theory, this is controlled by the level of turbulence that can vary. For instance, the low rate of turbulent reconnection can explain how magnetic flux can be accumulated prior to a solar flare. Third, the independence on the plasma parameters makes the turbulent reconnection a generic reconnection process in the astrophysical settings.

While the earlier theories of turbulent reconnection were focused on X-point reconnection, LV99 showed that Y-type reconnection can be fast. The tearing reconnection that was developed later (see Sec. X C) presented another version of the fast Y-type reconnection. Similar to LV99, it departed from the regular structures prescribed by Sweet-Parker and Petschek reconnection models.

Being Y-type fast reconnection, LV99 necessarily involved a significant volume of magnetized fluid in the reconnection process. This volume reconnection changes the nature of the energetic particle acceleration as it is discussed in Sec. IX C.

Since the seminal work by Alfvén (1942), the “flux freezing principle” has become a powerful tool to estimate the solutions of MHD equations in a diverse range of problems (Parker, 1979; Kulsrud, 2005). Besides its utility in plasma physics and in astrophysics, it has found many applications, for example, in explaining the low angular momentum of stars and the spiral structure of magnetic field lines in the solar wind. Despite its success in addressing diverse problems, however, there have been also fatal failures (see Burlaga et al., 1982; Khabarova and Obridko, 2012; Richardson et al., 2013; Eyink, 2015). For instance, assuming flux freezing, the magnetic topology in a high-conductivity magnetized fluid should not change rapidly, whereas the observations of solar flares and coronal mass ejections (CMEs) indicate otherwise. Also, in star formation, the magnetic pressure of the in-falling magnetized material should be strong enough to prevent gravitational collapse altogether, if magnetic flux-freezing held with good accuracy. Finally, if the magnetic field were frozen into the fluid in small-scale astrophysical dynamos, the tangled structure of the generated magnetic field would quench its further growth.

In order to account for these difficulties with magnetic flux-freezing, one might find it tempting to appeal to additional effects besides collisional resistivity. For instance, one might consider additional terms, such as Hall effect and electron pressure anisotropy, in the generalized Ohm's law (electron momentum equation) (Burch et al., 2016; Cassak, 2016). This does not answer the question, unfortunately, how magnetic reconnection events can be fast in fluids with high rates of collisions. Indeed, both observations and theory imply fast magnetic reconnection in the solar photosphere and other environments such as collisional parts of accretion disks, which cannot be explained by appealing to collisionless mechanisms (Chae et al., 2010; Singh et al., 2011). The problem persists even for weakly collisional systems where an ideal MHD-type description should be valid at length-scales much larger than the ion gyroradius, above which flux-freezing is often assumed to hold. However, observations indicate that magnetic structures with length scales much larger than ion gyroradius reconnect very rapidly, e.g., at the heliospheric current sheet (HCS) in the solar wind (Gosling, 2012). It is unclear how collisionless mechanisms at scales below the ion cyclotron radius can lead to fast magnetic reconnection on much larger scales. This is the essence of the “scale problem” in magnetic reconnection (Ji et al., 2019).

While the necessity of fast magnetic reconnection was widely accepted, most of the studies assumed that the fast reconnection happens for very special configurations of magnetic fields. This way the community could reconcile at least some of the observational evidence of fast reconnection and the concept of magnetic flux freezing. LV99, as discussed previously, identified magnetic reconnection as an intrinsic property of the MHD turbulent cascade, implying the violation of the flux freezing principle (see Vishniac and Lazarian, 1999, for an explicit statement about the flux freezing in turbulent fluids).

Below we provide a description of the Richardson dispersion, the phenomenon that manifests a gross violation of the flux freezing in magnetized turbulent fluids. This numerically confirmed phenomenon (see Sec. V B) provides another way of rederiving the LV99 rate of turbulent reconnection. A more rigorous approach to the flux freezing violation is discussed in Sec. IV.

Another way to understand the violation of flux freezing in turbulent fluids is related to the concept of Richardson dispersion. This process initially introduced in hydrodynamic turbulence carries over to MHD turbulence. Richardson (1926) empirically analyzed the dispersion of the volcanic ashes during the eruption and formulated the law that can be derived using hydrodynamic turbulence theory. Note, however, that the corresponding theory was formulated much later by Kolmogorov (1941). Using this theory, one can predict that the separation between two particles obeys the equation $d/dt[l(t)]∼v(l)∼αl1/3$⁠, where α is proportional to a cubical root of the energy cascading rate. The solution of this equation $l2∼t3$ describes the Richardson dispersion, i.e., the mean square separation between particles increases in proportion to the cube of time. The accelerated growth of the separation is due to the fact that the larger the separation, the larger the eddies that induce the separation. In Kolmogorov turbulence, the larger eddies have larger velocity dispersion. An interesting feature of this solution is that it provides the type of fast separation even if the initial separation of particles is zero. Formally, this means the violation of Laplacian determinism. Mathematically, the above paradox is resolved by accounting to the fact that the turbulent field is not differentiable (the Kolmogorov velocity field is Holder continuous, i.e., $|v(r1)−v(r2)|≲C|r1−r2|1/3.$⁠), and therefore, it is not unexpected that the initial value problem does not have a unique solution. In any realistic physical settings, the turbulence is damped at a nonzero scale l_min, which makes the paradox less vivid. At the scales less than l_min, the separation of points is not stochastic, and it follows the Lyapunov growth (see Lazarian, 2006). However, at the scales larger l_min but smaller than the injection scale of turbulence L_max, i.e., over the inertial range of turbulence, the explosive growth of separation is still present.

In terms of motions perpendicular to the local direction of the magnetic field, the scaling of velocities in the GS95 picture corresponds to the Kolmogorov one. Therefore, we do expect the plasma particles to separate in accordance with the Richardson law over the inertial range. The situation is illustrated in Fig. 5. The magnetic field is shown to be perpendicular to the plane of the picture. A cloud of test particles undergoes an ordinary diffusion for $a≫Lmax$⁠, undergoing Lyapunov separation of field lines for $a<lmin$ and exhibiting Richardson dispersion over the inertial range $[Lmax,lmin]$⁠. Obviously, the Richardson-type dynamics is not possible without magnetic field lines disentangling themselves via fast reconnection. This corresponds to the notion in LV99 that turbulent reconnection is required for turbulent magnetized fluid to preserve the fluid-type behavior. We note that if reconnection were slow, the intersecting magnetic field flux tubes would form unresolved knots. This inevitably would arrest the fluid motions, making the fluid with a dynamically important magnetic field Jello or feltlike substance. The spectrum of magnetic fields would be very modified in this case with a significant increase in energy at small scales. Only acoustic type waves can propagate in this substance and definitely no present-day numerical simulations of MHD type can describe it. This picture of slow reconnection in turbulent media contradicts both numerical simulations and observations of interstellar turbulence (see Secs. II A and II B). The turbulent reconnection, on the contrary, predicts the free evolution of eddies that produce the Kolmogorov-type cascade for eddies perpendicular to the local direction of the magnetic field as discussed in Sec. II B.

The Richardson dispersion is a numerically proven phenomenon (see Sec. V B), and it was used in ELV11 to rederive the LV99 rates of reconnection.

Consider first the Sweet-Parker reconnection. The Ohmic diffusion induces the mean square displacement across the reconnection layer in a time t given by

where $λ=c2/4πσ$ is the magnetic diffusivity. As in our previous treatment, the magnetized plasma is ejected out of the sides of the reconnection layer of length L_x at a velocity of order V_A. Thus, the time that the lines can spend in the resistive layer is the Alfvén crossing time $tA=Lx/VA$⁠. Therefore, the reconnected field lines are separated by a distance

where S is the Lundquist number. Combining Eqs. (34) and (44), one obtains the Sweet-Parker result, $vrec=VA/S$⁠.

The difference with the turbulent case is that instead of Ohmic diffusion, the Richardson dispersion dominates (ELV11). In this case, the mean squared separation of particles follows the Richardson law, i.e., $⟨|x1(t)−x2(t)|2⟩≈ϵt3$⁠, where t is the time and ϵ is the energy cascading rate. For sub-Alfvénic turbulence, $ϵ≈uL4/(VALi)$ (see LV99), and therefore, one can write

where it is assumed that $L<Li$⁠. Combining Eqs. (34) and (45), it is easy to obtain

in the limit of $L<Li$⁠. Similar considerations allowed ELV11 to rederive the LV99 expression for $L>Li$⁠, which differs from Eq. (46) by the change in the power 1/2 to $−1/2$⁠.

Magnetic fields produce a new effect compared to the Richardson dispersion in hydrodynamics. In magnetized fluids, one can trace magnetic line separation. This effect has been discussed for decades in the cosmic ray literature as the cause of the perpendicular diffusion of cosmic rays in the galactic magnetic field (Jokipii, 1973). It is worth mentioning that the correct quantitative description of the effect was obtained only in LV99. Earlier papers assumed the diffusive separation of magnetic field lines, which is the effect taking place only if magnetic field lines are separated by more than the turbulence injection scale L. It was shown in LV99 that at scales less than L, the magnetic field lines separate superdiffusively. This changes a lot both the propagation and the acceleration of cosmic rays (see Lazarian and Yan, 2014).

Dealing with magnetic reconnection in Sec. III, we have discussed the separation of magnetic field lines. Figure 6 illustrates the spread of magnetic field lines in the perpendicular direction in the situation when magnetic field lines are traced by particles moving along them. This figure confirms the prediction for the magnetic field line separation in Eq. (38), i.e., the squared separation between the magnetic field lines increases in proportion to the cube of the distance measured along the magnetic field lines (see Fig. 6).

The other important prediction in Eq. (38) is the change in the square of the separation between the magnetic field lines, which increases with the forth power of the Alfvén Mach number, i.e., $∼MA4$⁠. The successful testing of this prediction is demonstrated in Fig. 7, which is taken from the study by Xu and Yan (2013). The asterisks correspond to the calculations obtained using the parallel particle velocity. This way of calculating the separation provides a better description of magnetic field wandering. The fitted dependence is $MA−3.8$⁠, which is very different from the $MA2$ dependence traditionally accepted in the cosmic ray literature.

We note that the magnetic field tracing with particles is not the only way to numerically explore the magnetic field separation. In fact, the analytical predictions of magnetic field line superdiffusion were also explored numerically by the direct tracing of magnetic field lines first with low resolution by Lazarian et al. (2004) and with higher resolution in Beresnyak (2013).

Taken together, these numerical testing provides a solid confirmation of the predictions given by Eq. (38). It is significant that the turbulent magnetic field wandering induced by strong MHD turbulence provides the separation of magnetic field lines $δlines2∼s3$⁠, where the distance s along the instantaneous realization of a magnetic field line. Thus, for sufficiently large s, all parts of the volume of magnetized fluid get connected. In other words, the entire volume becomes accessible to particles moving along magnetic field lines. This is important for a wide variety of astrophysical processes, e.g., the cosmic ray propagation and acceleration (Yan and Lazarian, 2008; Lazarian and Yan, 2014) and heat propagation (Narayan and Medvedev, 2001; Lazarian, 2006).

We would like to stress again that magnetic wandering cannot be understood using the notion of magnetic field lines preserving their identify in turbulent flows. Magnetic field lines wandering in the turbulent magnetized flow constantly reconnect, inducing the exchange of plasmas and magnetic field. In view of this, it is interesting to recall the considerations on the nontrivial nature of astrophysical magnetic field lines that can be found in a prophetic book by Parker (1979). There it was mentioned that a magnetic field line may be traced through a significant part of the galaxy. We add that this is only possible through constant reforming of the magnetic field lines via reconnection.

We note that there is a formal analogy between the Richardson diffusion that predicts that the mean squared separation between particles $δparticles2$ that grow with time as t³ and the mean squared separation between the magnetic field lines $δlines2$ that grow with the distance s as s³. Therefore, rather than using a new term to describe the latter dependence theoretically predicted in LV99, it was decided in ELV11 to term this new type of dispersion as “the Richardson diffusion in space,” while keeping the term “Richardson dispersion” for the original time-dependent process empirically discovered in Richardson (1926).

The derivation of the LV99 reconnection rate assumes that all magnetic energy is converted to the kinetic energy of the outflow. This may not be true when M_A is large. In this case, the energy dissipation in the reconnection layer becomes non-negligible, and the decrease in the outflow velocity is expected (see ELV11).

The corresponding effect of turbulent dissipation can be estimated from steady-state energy balance in the reconnection layer

where the kinetic energy of the outflow is balanced against the magnetic energy transported into the layer minus the energy dissipated by turbulence. The turbulent dissipation is given by $ε=uL4/VALi$ for sub-Alfvénic turbulence (Kraichnan, 1965). Dividing (47) by $Δ=Lxvrec/vout$⁠, one obtains

The corresponding solutions are easy to get by introducing the ratios $f=vout/VA$ and $r=2MA4(Lx/Li)$ which are the measures of the outflow speed and the energy dissipated by turbulence in units of the available magnetic energy. This gives

When r = 0, the only solution of (49) with f > 0 is $f=1,$ which is the LV99 estimate $vout=VA$ for $MA≪1$⁠. For the larger values of r, $f≃1−(r/2)$⁠, in agreement with the formula $f=(1−r)1/2$ that follows from Eq. (65) in ELV11. This implies a decrease in v_out as compared to V_A. Note that formula (49) can be used when r is not too large. The largest value of r for which a positive, real f exists is $rmax=2/27≈0.385$⁠. For r_max, one gets $fmin=1/3≈0.577$⁠. This implies that the LV99 approach is limited to sufficiently small M_A. For $Lx≃Li$⁠, one may consider values of M_A up to 0.662. Therefore, the original LV99 approach is accurate for $MA≲1$⁠.

While the accuracy of LV99 expressions can be decreased due to the additional dissipation of the turbulent energy, the predicted phenomenon of fast reconnection arising from the turbulence persists for any M_A. In fact, one can see that the effect of the reduced outflow velocity “increases” the reconnection rate. The reason is that field-lines now spend a time $Lx/vout$ exiting from the reconnection layer, greater than that assumed in LV99 by a factor of $1/f$⁠. This implies a thicker reconnection layer Δ due to the longer time-interval of Richardson dispersion in the layer, greater than the LV99 estimate by a factor of $(1/f)3/2$⁠. The net reconnection speed $vrec=voutΔ/Lx$ is thus larger by a factor of $(1/f)1/2$⁠. The increased width Δ more than offsets the reduced outflow velocity v_out. At even larger M_A, magnetic fields at a large scale are weak that they are easily bent and twisted by the fluid turbulent motions. However, at the scales smaller than l_A, we expect the original LV99 treatment to be accurate again.

The model of turbulent reconnection suggested in LV99 departed from the earlier attempts to explain astrophysical magnetic reconnection as something that required special conditions. On the contrary, it suggested that reconnection takes place everywhere in realistically turbulent fluids. Moreover, the transition to turbulent reconnection is inevitable as the 3D outflow from the reconnection region gets turbulent for sufficiently large Lundquist numbers.

The model identified turbulent magnetic field wandering as the process that opens up the reconnection region outflow making reconnection fast, i.e., independent of fluid resistivity or microscopic plasma effects. The derived dependences for the field wandering were successfully tested with the numerical simulations.

The model provided analytical expressions relating the reconnection rate with the level of turbulence. It also provided the relationships between the thickness of the reconnection layer and the level of turbulence. This explains that the rate of magnetic reconnection may vary significantly and also provides ways for quantitative numerical and observational testing of the model.

In Sec. III B, we discussed the issue of flux freezing on the basis of the simplified treatment of Richardson dispersion, treating the latter more like an empirical fact and using the corresponding scaling relations to rederive the LV99 reconnection rates. However, the phenomenon of Richardson dispersion presents a much more fundamental property of magnetic fields in turbulent fluids. Below we discuss the violation of flux freezing on the basis of a more rigorous theoretical treatment.

The conventional flux-freezing principle has not only problems with explaining observations. It also encounters serious theoretical difficulties. This can be exemplified by naively extending the concept to turbulent MHD flows with both small resistivity and small viscosity. To understand the subtlety of such singular limits, consider the problem of incompressible fluid turbulence. As is well known, the kinetic energy dissipation rate, in a fluid with viscosity ν and velocity field $u$⁠, does not vanish in the limit of vanishing viscosity [e.g., see Sreenivasan (1984, 1998)]. On the contrary, $limν→0ν|∇u|2>0$⁠, which implies a “dissipative anomaly.” This in turn requires that velocity gradients must blow up, with $∇u→∞$ in the limit $ν→0$⁠. In fact, as Onsager pointed out in a prescient work (see Sec. IV B below), the velocity field in a turbulent fluid must develop Hölder-type singularities in the limit as viscosity tends to zero (Onsager, 1949; Eyink and Sreenivasan, 2006; Eyink, 2018). These singularities render all hydrodynamic equations containing velocity gradients naively undefined. Extending Onsager's analysis to MHD plasmas, turbulent solutions are found to suffer several dissipative anomalies that require diverging gradients of all MHD fields, e.g., velocity field, magnetic field, and density (see, e.g., Mininni and Pouquet, 2009; Caflisch et al., 1997; Eyink and Drivas, 2018). For example, magnetic-field gradients blow up in the limit as magnetic diffusivity η tends to zero, with $∇B→∞$ as $η→0$⁠. Because velocity fields in particular become nonsmooth and singular as $ν→0$⁠, the magnetic field lines cannot freeze into the fluid in a standard deterministic sense. Instead, a nondeterministic frozen-in property known as “stochastic flux freezing” (Sec. IV C) has been established by Eyink (2009) and (2011), which is intrinsically random in nature due to the roughness of the advecting velocity field (Eyink, 2010) and is associated with the “spontaneous stochasticity” of Lagrangian particle trajectories (Bernard et al., 1998; Eyink et al., 2011, 2013).

The latter concept is crucial and must be briefly explained. In essence, spontaneous stochasticity means that stochastic Lagrangian particle trajectories remain nonunique and indeterministic (random) for turbulent flows even in the limit of vanishing viscosity and molecular diffusivity (Bernard et al., 1998). This phenomenon is directly related to the Richardson superballistic behavior (Richardson, 1926) (see Sec. III B). In a turbulent flow, therefore, infinitely many magnetic field lines are advected to each point in space even in the ideal limit and must be averaged to obtain the resultant magnetic field (Eyink, 2011; Eyink et al., 2013). Due to Richardson dispersion, in fact, even a tiny breakdown of magnetic flux freezing at small scales will be explosively amplified into larger scales in the turbulent inertial range (Eyink, 2015, see also Sec. IV D below). The rates of such rapid amplification are nonvanishing in the ideal limit (Eyink, 2011; Eyink et al., 2011). This breakdown of conventional magnetic flux freezing has already been verified in numerical simulations of resistive MHD turbulence at high conductivity (e.g., Eyink et al., 2013).

To appreciate the importance of (asymptotic) singularities in MHD, we may note that without the differentiability of fields, even fundamental concepts such as the notion of field lines needs to be reconsidered. Taking magnetic field $B(x)$ as an example, its “field line” through point $x$ is conventionally defined as the integral curve $ξ(s;x)$ whose tangent vector at any point is parallel to the vector field at that point, so that, parameterized with arc length s

If the magnetic field satisfies the following inequality for some constant C > 0

with h = 1 (Lipschitz condition), then there exists a unique field-line $ξ(s;x)$ passing through point $x$⁠; however, if this condition holds merely with $0<h<1$ (Hölder continuity), then there are generally infinitely many such integral curves (see Jafari and Vishniac, 2019a, 2019b; Jafari et al., 2019b). If the Lipschitz condition is not satisfied and hence the uniqueness theorem cannot be applied for a given magnetic field, what would a magnetic field line mean after all? When resistivity is tiny but nonzero, then $B$ is differentiable, and indeed, Eq. (50) has unique solutions. However, these solutions exhibit such extreme sensitivity to initial conditions that the distance $|ξ(s;x′)−ξ(s;x)|$ between field lines through neighboring points becomes “independent of” $x′−x$ at distances $s≫|x′−x|$ along the lines (LV99) (we provide more detail on this effect in Sec. III B). Such explosive separation of magnetic field lines is not observed if the field lines exhibit standard “deterministic chaos,” where the memory of the initial separation is preserved. Exactly analogous difficulties plague the concepts of “Lagrangian fluid particle” or “fluid trajectory,” which are ill-conditioned and asymptotically nonunique if the fluid velocity $v(x,t)$ becomes spatially non-Lipschitz as $ν→0.$

The singularities associated with dissipative anomalies lead to subtle effects whose consistent description requires great care; otherwise, contradictory conclusions are obtained. Similar difficulties were encountered already many decades ago in quantum field theory and condensed matter physics, which could be alleviated by mathematical methods collectively known as regularization and renormalization. As we review in Sec. IV B, such methods can be applied also to plasma turbulence by exploiting a spatial coarse-graining operation, which removes divergences and leads to meaningful equations. Physically speaking, “coarse-graining” simply means averaging quantities in a volume in space. Indeed, what is measured in the laboratory as magnetic or velocity fields at a point in space and time is actually always an average over a small volume around that point during a time interval. Paradoxes are avoided by developing a theory based upon such quantities that can be measured in principle rather than idealized pointwise fields that are strictly unobservable. Such renormalization introduces an “observation scale” or “resolution scale” which is entirely arbitrary. However, obviously, no objective physical fact can depend upon this arbitrary scale. In the language of quantum field theory and statistical physics, this principle is called “renormalization-group invariance” (Gross, 1976).

We illustrate Onsager's renormalization-type arguments in one of the simplest models with realistic astrophysical applications and the standard MHD equations for a compressible, magnetized, single-component fluid with mass density $ρ,$ velocity $v,$ magnetic field $B,$ and internal energy density u

Here, $I$ is the identity tensor and $p=p(ρ,u)$ is the thermodynamic equation of state for the pressure and

with shear viscosity $μ,$ bulk viscosity $ζ,$ thermal conductivity $κ,$ and electrical conductivity σ (or magnetic diffusivity $η=c2/4πσ$⁠). See Landau et al. (2013), Chap. VIII, Secs. 65 and 66.

The local balance of mechanical energy for the above model is given by the following equation:

Just as for the case of incompressible, neutral fluids studied by Onsager, turbulent solutions are characterized by “dissipative anomalies” (see, e.g., Mininni and Pouquet, 2009) in which

Similar anomalous dissipation must then appear also in the balances of internal energy and entropy. Such anomalies obviously require diverging space-gradients $|∇v|, |∇B|→∞$ for $μ, ζ, κ, η→0$ or a “violet catastrophe” in the words of Onsager (1945).

To obtain a well-defined dynamical description in this singular ideal limit, one must regularize these UV divergences. For example, one may employ a coarse-graining or “block-spin” regularization, defined as

and likewise for $v,$$ρ,$u, etc., with $Gℓ(r)=ℓ−3G(r/ℓ)$ for some positive, smooth, rapidly decaying kernel $G(ρ)$ normalized as $∫d3ρ G(ρ)=1.$ Coarse-grained gradients such as $∇B¯ℓ(x,t)$ necessarily remain finite as $μ,ζ,κ,η→0,$ and their UV divergences are removed. It should be emphasized again that such a coarse-graining is purely passive and corresponds to observing the flow “with spectacles off” so that eddies of size $<ℓ$ cannot be resolved. The basic principle of renormalization-group invariance is that this regularization scale $ℓ$ is completely arbitrary, and no objective physics can depend upon it. Because the coarse-graining operation commutes with space-time derivatives, the effective equations of motion for eddies of size $>ℓ$ are exactly the same as the original MHD equations with just an additional overbar $()¯ℓ$⁠.

However, even when there are dissipative anomalies in the microscopic or “bare” balance equations for conserved quantities, the dissipative transport terms such as $∇×(η∇×B)¯ℓ$ in the coarse-grained MHD equations can be shown to vanish for $μ,ζ,κ,η≪1$ at fixed $ℓ$ or for $ℓ≫1$ at fixed $μ,ζ,κ,η$⁠. These dissipative terms are thus “irrelevant” in the renormalization-group sense. It follows, when $μ,ζ,κ,η≪1,$ that there exists an “inertial-range” of length-scales $ℓ$ at which the ideal MHD equations hold in the “coarse-grained sense,” i.e.,

Taking the limit $μ,ζ,κ,η→0$⁠, these coarse-grained equations hold for all $ℓ>0,$ which is equivalent to the standard mathematical notion of a “weak”/“distributional” solution or what Onsager called a “more general description” than standard differential equations.

The validity of the ideal MHD equations in this generalized sense in not sufficient, however, to guarantee that naive balances of mechanical energy, internal energy, and entropy will hold. To show this clearly, the above ideal equation can be cosmetically rewritten by introducing the “density-weighted Favre-average”

and “pth-order cumulants” $τ¯ℓ(f1,f2,…,fp)$ and $τ̃ℓ(f1,f2,…,fp)$ for the two types of coarse-graining averages. With these definitions, the effective equations of motion for scales $>ℓ$ become

and

Here, we have introduced the “turbulent stress tensor”

“turbulent electromotive force (emf)”

“turbulent internal energy”

and “turbulent enthalpy transport”

with $h=u+p$ the enthalpy per volume.

Although Eqs. (55)–(58) are exactly equivalent to the coarse-grained ideal MHD equations (51)–(54), they appear “nonideal.” Integrating out the small-scale eddies has led to “renormalized” equations with new terms $Tℓ,$ $εℓ,$ $δu¯ℓ,$ and $Jℓh$⁠, which would not appear in naive ideal MHD equations for the resolved fields $ρ¯ℓ,$ $ṽℓ,$ $B¯ℓ,$ and $u¯ℓ$ at length-scales $>ℓ.$ The precise condition for the validity of the ideal equations (51)–(54) at length-scale $ℓ$ is that the standard dissipative transport terms are negligible compared with these new terms. For example, the ideal magnetic induction equation (57) is valid precisely when

One of the Onsager key insights was that contributions from coarse-graining such as $εℓ$ can be expressed entirely in terms of “space-increments” (see Eyink and Aluie, 2006; Eyink and Drivas, 2018), such as

It follows that

In that case, the condition (60) becomes

which precisely defines the inertial-range of scales $ℓ$ for which the ideal induction equation is valid. This condition can be restated as $ℓ≫ℓη$⁠, where $ℓη$ is the resistive dissipation length specified by $ℓηδv(ℓη)≃η.$ Similar conditions guarantee the validity of all of the other ideal equations as well (see Eyink and Drivas, 2018).

From the coarse-grained ideal MHD equations rewritten in this manner, one obtains a balance equation for “resolved mechanical-energy”

where the renormalization has introduced new apparently nonideal terms. In particular, the right side of this balance contains the “turbulent energy flux”