Multi-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer
UUnder consideration for publication in J. Plasma Phys. Multi-scale dynamics of magnetic flux tubesand inverse magnetic energy transfer
Muni Zhou , Nuno F. Loureiro , Dmitri A. Uzdensky Plasma Science and Fusion Center , Massachusetts Institute of Technology , Cambridge , MA02139 , USA Center for Integrated Plasma Studies , Physics Department , UCB-390 , University of Colorado , Boulder , CO 80309 , USA(Received xx; revised xx; accepted xx)
We report on an analytical and numerical study of the dynamics of a three-dimensionalarray of identical magnetic flux tubes in the reduced-magnetohydrodynamic descriptionof the plasma. We propose that the long-time evolution of this system is dictated by flux-tube mergers and that such mergers are dynamically constrained by the conservation ofthe pertinent (ideal) invariants, viz. the magnetic potential and axial fluxes of each tube.We also propose that in the direction perpendicular to the merging plane, flux tubesevolve in critically-balanced fashion. These notions allow us to construct an analyticalmodel for how quantities such as the magnetic energy and the energy-containing scaleevolve as functions of time. Of particular importance is the conclusion that, like itstwo-dimensional counterpart, this system exhibits an inverse transfer of magnetic energythat terminates only at the system scale. We perform direct numerical simulations thatconfirm these predictions and reveal other interesting aspects of the evolution of thesystem. We find, for example, that the early time evolution is characterized by a sharpdecay of the initial magnetic energy, which we attribute to the ubiquitous formation ofcurrent sheets. We also show that a quantitatively similar inverse transfer of magneticenergy is observed when the initial condition is a random, small-scale magnetic seed field.
1. Introduction
Inverse transfer of magnetic energy from small to large scales in highly conducting plas-mas is an important topic in modern plasma physics and astrophysics. It is particularlyimportant for understanding magnetogenesis (Kulsrud & Zweibel 2008), i.e., the originof large-scale magnetic fields that are ubiquitous in space and astrophysical systems.Examples of such systems include planets, such as Earth and Jupiter, main-sequence starsincluding the Sun, neutron stars and white dwarfs, accretion disks around black holes, theinterstellar medium (ISM) in galaxies including the Milky Way, and, on even larger scales,the hot gas in the intracluster medium (ICM). While the leading paradigm for the originof large-scale fields invokes turbulent magnetohydrodynamic (MHD) α Ω dynamo (Parker1955), driven by turbulence coupled with large-scale differential rotation (Goldreich &Lynden-Bell 1965), successful dynamo action needs a seed magnetic field of sufficientstrength and sufficiently large coherence length to overcome resistive (or other nonideal)losses. In some natural systems, however, the seed field may originate at very small,kinetic scales and then has to traverse many orders of magnitude in coherence length toreach the desired, extremely large, astronomical scales.An important example of this scenario arises when the initial field is generatedvia the Weibel instability (Weibel 1959; Fried 1959) — for example, in a collisionless a r X i v : . [ a s t r o - ph . H E ] J un Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky shock (Medvedev & Loeb 1999). This is relevant to systems as diverse as gamma-ray-burst (GRB) jets, supernova explosions, and large-scale accretion shocks at the outskirtsof galaxies during early stages of galaxy formation (Gruzinov 2001). The initial field inthese situations can be conceptualized as small-scale (on the order of a few electron or ionskin depths) elongated magnetic flux ropes associated with current filaments (Silva et al. et al. β ∼ ), they might dissipate rather quickly dueto their small scales. On the other hand, current filaments are prone to coalesce and thusgradually transfer energy to larger scales. Therefore, it is not a priori clear whether suchelectron-scale (or ion-scale) fields survive on long time-scales, how quick and effective theinverse magnetic-energy cascade arising from their mergers might be, and what fraction ofthe initial magnetic energy the inverse cascade can deliver to large-enough scales, wherethe field can be picked up and amplified by the ambient turbulence (e.g., in the contextof Galactic magnetogenesis). An estimation of such scales can be found in Appendix A.In the context of GRB prompt emission and afterglow, an important additional questionis whether the Weibel-produced field in a relativistic shock can survive long enough toexplain the observed powerful synchrotron emission (Medvedev & Loeb 1999; Silva et al. et al. et al. et al. et al. et al. et al. et al. et al. et al. et al. et al. et al. et al. a , b ). These studies consider only afew interacting flux ropes and thus cannot be viewed as statistical investigations of ahierarchical, multi-stage coalescence process. However, in recent years, several theoreticalstudies have used statistical approaches to understand systems of merging structuresin various astrophysical contexts (Gruzinov 2001; Medvedev et al. et al. et al. et al. a , b ; Zrake & Arons 2017)and dynamics of coronal magnetic loops (Uzdensky & Goodman 2008). Despite all theseefforts, important questions remain to be answered. For example, what are the underlyingphysical mechanisms of the merging process, and what are the merging statistics in 3Dsystems? On the experimental side, there have been detailed investigations of the 2Ddynamics during the merger of pairs of annular flux ropes (Yamada et al. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer et al. et al. et al. et al. et al. et al. et al. et al. et al. ˜ t − , where ˜ t is time normalized tothe (appropriately defined) reconnection timescale; and the correlation length of thefield grows as ˜ t / . Our 2D hierarchical model was directly verified by 2D reduced-magnetohydrodynamics (RMHD) simulations that we carried out.In the present paper, we generalize our theory to 3D systems (Section 2 and Ap-pendix B), which we successfully test with 3D RMHD simulations (Section 3). InAppendix A, we apply the scaling laws of magnetic field evolution obtained in thiswork to the problem of galactic magnetogenesis, and address the relevance of small-scalemagnetic seed fields to the galactic dynamo. Additionally, we demonstrate (Appendix C)that reconnection-enabled inverse transfer remains a robust phenomenon when the initialcondition is a random magnetic field, rather than a structured array of magnetic fluxtubes. Conclusions from this study and a discussion of some implications of our work canbe found in Section 4.
2. Statistical theory of three-dimensional flux-tube mergers
Three-dimensional hierarchical model
In this section, we introduce a 3D hierarchical model to describe the merger ofmagnetic flux tubes. Our goal is to obtain explicit expressions for the evolution ofquantities such as magnetic energy and the coherence length of magnetic field in boththe merging plane and the direction perpendicular to it, as functions of time. Similarto its 2D counterpart (Z19), the 3D model is constructed in the framework of resistiveincompressible RMHD, although the key ideas should remain valid in some collisionlessplasma descriptions.In our RMHD system, a constant strong background magnetic guide field B g = B g ˆ z is assumed, and the xy plane is referred to as the perpendicular plane. We consideran ensemble of magnetic flux tubes of the same size and with the same magnetic fieldstrength as the initial condition. The magnetic field in each flux tube can be decomposedinto the constant guide field B g and the perpendicular magnetic field B ⊥ , correspondingto an electric current parallel or anti-parallel to the guide field. For simplicity, each flux Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky tube is assumed to be of cylindrical shape, with radius R (in the xy plane) and length l inthe z direction. The lengths R and l represent the coherence lengths of the perpendicularmagnetic field, B ⊥ , in the directions perpendicular and parallel, respectively, to theguide field. In the RMHD approximation, B ⊥ and R are asymptotically smaller than B g and l respectively. The magnetic field lines in each flux tube thus have eitherright- or left-handed helical shapes depending on the direction of the current, withasymptotically small pitch angles. The flux tubes are assumed to be randomly distributedin the system and volume-filling, with half of them carrying currents parallel to theguide field, and the other half anti-parallel. Each flux tube is a system of concentriccylindrical magnetic flux surfaces nested around the tube’s single magnetic axis alignedwith the guide field. It contains axial magnetic flux Ψ a = B g πR and poloidal (i.e.,azimuthal around the magnetic axis) magnetic flux Ψ p = B tube ⊥ Rl , where B tube ⊥ is thecharacteristic perpendicular magnetic field of each flux tube. In the perpendicular xy plane, the magnetic potential of a flux tube is defined as ψ tube = B tube ⊥ R , similar toits 2D counterpart. Thus, Ψ p = ψ tube l . Fig. 1 schematically represents the simplifiedgeometry of a single tube as we have just described.The magnetic energy associated with the perpendicular magnetic field contained in aflux tube is (cid:15) M (cid:39) πR l ( B tube ⊥ ) / (8 π ) = ( ψ tube ) l/ . Denoting by N the number densityof flux tubes per unit volume [ N ∼ / ( πR l ) ], the total magnetic energy density (of theperpendicular magnetic field) of the system is, therefore, E M = (cid:15) M N = ( B tube ⊥ ) / (8 π ) .In addition, we denote the number density of flux tubes per unit area in the xy plane as N ⊥ , where N ⊥ ∼ / ( πR ) and N ⊥ = N l .A non-zero local magnetic helicity exists associated with each flux tube, defined as h = (cid:82) dV tube A · B , where B = ∇ × A . In the RMHD limit, however, this expressionsimplifies as follows: h = (cid:90) dV tube A · B = (cid:90) dV tube ( A z B g + A ⊥ B tube ⊥ ) ≈ − B g (cid:90) dV tube ψ tube ≈ − ( B g πR )( B tube ⊥ Rl ) = − Ψ a Ψ p , (2.1)represented by the product of axial flux and poloidal flux. In the derivation, the approx-imation is made based on the RMHD ordering A ⊥ (cid:28) A z and B tube ⊥ (cid:28) B g .When two flux tubes with parallel axial currents are close to each other, they areunstable to the coalescence instability and approach each other (Finn & Kaw 1977),merging into a wider flux tube via reconnection of their poloidal fluxes. For simplicity,we assume the merging process to happen in discrete stages, denoted by index n . Weassume the mergers to proceed in hierarchical fashion and, therefore, take all flux tubesto be identical at each stage. The index n is thus added as a subscript to the abovequantities to indicate the generation of mergers to which they belong. As the flux tubesmerge, the above quantities evolve; i.e., they are functions of n and so, implicitly, of time.We also assume mergers to be binary ( N ⊥ ,n +1 = N ⊥ ,n / ) and governed by basicconservation laws [similar arguments can be found in Medvedev et al. (2005); Polomarov et al. (2008); Lyutikov et al. (2017 a ); Zrake & Arons (2017)].First, the merger of two flux tubes should conserve the axial magnetic flux: Ψ a,n +1 =2Ψ a,n . Since the axial magnetic field B g is a constant, it follows that the area of theperpendicular cross section should be conserved and, therefore, R n +1 = √ R n . (2.2)Second, the perpendicular magnetic potential associated with each flux tube is not ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer R , length l , perpendicular magnetic field B tube ⊥ , axialflux Ψ a and poloidal flux Ψ p .affected by reconnection and, as in the two-dimensional case (Z19), should remain thesame: ψ tube n +1 = ψ tube n . Combined with Eq. (2.2), we obtain the evolution of perpendicularmagnetic field and its energy density: B tube ⊥ ,n +1 = B tube ⊥ ,n / √ , E M,n +1 = E M,n / . (2.3)We assume that the parallel coherence length of each flux tube is determined by thecritical balance condition, i.e., the statement that the pertinent parallel and perpendiculartimescales should roughly be the same (Goldreich & Sridhar 1995). Parallel dynamics isdictated by the propagation of shear Alfvén waves along the magnetic field. Therefore,the parallel timescale is τ (cid:107) ∼ l/V A , with V A the Alfvén speed computed with the guidefield. In the perpendicular plane, the relevant timescale is τ ⊥ ∼ R/v A , with v A the Alfvénspeed computed with B tube ⊥ . That this is the relevant perpendicular timescale to considerremains true even during mergers, since merging tubes are not static, but rather movein the perpendicular plane at roughly the ambient v A . Note also that another relevantperpendicular timescale is the reconnection time — that is, the time it takes to merge twotubes through reconnection. This is always slower than Alfvénic, implying that criticalbalance can always establish itself dynamically.We thus find τ (cid:107) ∼ τ ⊥ ⇒ l/B g ∼ R/B tube ⊥ . (2.4)Plugging Eqs. (2.2-2.3) into this relation, we obtain the evolution of coherence length ofmagnetic field in the parallel direction: l n +1 = 2 l n . (2.5)That is, l ∼ R , i.e., the magnetic flux tubes grow faster in the parallel direction thanon perpendicular planes. Knowing the magnetic potential remains the same during themerging process, the poloidal flux ( Ψ p = ψ tube l ) increases with the length of the fluxtubes: Ψ p,n +1 = 2Ψ p,n . As a result of the increasing flux contained in each generation offlux tubes as they merge and grow, the magnetic helicity associated with each flux tubealso increases: h n +1 = 4 h n .Our next step in order to obtain a continuous-time description of the evolution is to Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky compute the lifetime of each generation; this entails a discussion of the reconnectionprocess mediating the mergers.The Lundquist number corresponding to the merging tubes at any given n -th stage isdefined as S n ≡ R n v A,n /η , where η is the (constant) magnetic diffusivity. Importantly, S ∝ RB tube ⊥ /η ∝ ψ tube /η . Therefore, and as we pointed out in our earlier work (Z19),the constancy of the magnetic potential ψ tube during the mergers implies the same forthe Lundquist number: S n +1 = S n . In particular, in MHD, if S (cid:46) , reconnectionproceeds in the Sweet-Parker (SP) regime (Sweet 1958; Parker 1957). In this regime, thedimensionless reconnection rate, defined as β rec ≡ v rec /v A , scales with the Lundquistnumber as β rec (cid:39) S − / . In contrast, if S (cid:38) , then magnetic reconnection pro-ceeds in the plasmoid-dominated regime (Biskamp 1986; Loureiro et al. et al. et al. et al. et al. β rec (cid:39) . .Therefore, we reach the non-trivial conclusion that the normalized reconnection rate, β rec , remains unchanged throughout the evolution, and the merging process remains inthe same reconnection regime (SP or plasmoid-dominated) in which it started initially.We can now proceed to calculate the lifetime of each generation. We note first that,by definition, the time it takes to reconnect (merge) two flux tubes is τ rec = β − R/v A = β − τ ⊥ . Since β rec (cid:28) , the time scale for the n -th stage is approximately the reconnectiontime at this stage: τ n ≈ τ rec , n ∝ β − R n /B tube ⊥ ,n . (2.6)Consequently, τ n +1 = 2 τ n . The time taken to reach the n -th stage is thus: t n = n − (cid:88) k =0 τ k = τ n − (cid:88) k =0 k ≈ τ n , n (cid:29) , (2.7)where τ = β − R /B tube ⊥ , , R and B tube ⊥ , are the radius and perpendicular magnetic fieldof the initial flux tubes respectively. Therefore we can express the index n as a functionof time: n = log t/τ . The explicit formulae for the different quantities as a function of n can be obtained from the recursive relations derived earlier: B tube ⊥ ,n = B tube ⊥ , − n/ , E M,n = E M, − n , (2.8) k ⊥ ,n = k ⊥ , − n/ , R n = R n/ , l n = l n , (2.9) N ⊥ ,n = N ⊥ , − n . (2.10)By using n = log t/τ , we obtain the continuous time evolution of the relevant quantities: B tube ⊥ = B tube ⊥ , ˜ t − / , E M = E M, ˜ t − , (2.11) k ⊥ = k ⊥ , ˜ t − / , R = R ˜ t / , l = l ˜ t, (2.12) N ⊥ = N ⊥ , ˜ t − , (2.13)where ˜ t ≡ t/τ is time normalized to the initial reconnection timescale and k ⊥ ≡ π/R is the perpendicular wavenumber. We note that the aspect ratio of a flux tube l/R ∼ ˜ t / . Therefore, as flux tubes merge and grow, they become progressively more slender.However, they always satisfy the critical balance and remain (marginally) stable to kink-type instabilities (Appendix B). We also note that l = l ˜ t = l ( v A, β rec /R ) t . Combinedwith the critical balance scaling of initial flux tubes l /V A ∼ R /v A, , where V A is theAlfvén speed based on the guide field, we obtain l = V A β rec t . We can thus consider V A β rec as the characteristic speed for the growth of tube’s length. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer Self-similar properties and magnetic power spectra
The power-law time dependencies [Eqs. (2.11)–(2.13)] suggest that the system evolvesin a self-similar manner. We note that the perpendicular dynamics in the 3D RMHDsystem is the same as that in its 2D counterpart (Z19). Therefore, the self-similarproperties and the behavior of perpendicular magnetic power spectrum should also beidentical between the 2D and 3D cases. We discuss these properties in more detail inthis section, for the completeness of this paper. In the perpendicular plane, the growinglength scale and decreasing field strength [Eq.(2.11) and Eq.(2.12)] can be representedby a dynamical renormalization, which we show as follows. With the power-law timeevolution, quantities at any given two moments of time ˜ t and ˜ t (cid:48) can be related througha scaling factor λ = (˜ t (cid:48) / ˜ t ) / k ⊥ (˜ t ) = λk ⊥ (˜ t (cid:48) ) , B tube ⊥ (˜ t ) = λB tube ⊥ (˜ t (cid:48) ) , ˜ t = λ − ˜ t (cid:48) . (2.14)That is, the time evolution of the system is equivalent to performing the above scaletransformation. In a system consisting of volume-filling flux tubes, the perpendicularmagnetic field B ⊥ follows the same scaling as B tube ⊥ . That is, B ⊥ (˜ t ) = λ B ⊥ (˜ t (cid:48) ) .Let us see how this self-similarity affects the evolution of perpendicular magnetic powerspectrum, which we define as M ( k ⊥ , ˜ t ) ≡ π πk ⊥ (2 π ) (cid:90) d r e i k ⊥ · r (cid:104) B ⊥ ( x , ˜ t ) · B ⊥ ( x + r , ˜ t ) (cid:105) , (2.15)where the angle brackets represent a 3D average over all x , and r is a vector on theperpendicular ( x - y ) plane. The spectrum is also averaged over the z dimension, which isnot explicitly shown in Eq. (2.15). Using Eq. (2.14), the spectra at different times arerelated as M [ k ⊥ ( t (cid:48) ) , ˜ t (cid:48) ] = M ( k ⊥ /λ, λ ˜ t ) = 18 π π ( k ⊥ /λ )(2 π ) (cid:90) d ( λ r ) e i ( k ⊥ /λ ) · ( λ r ) (cid:104) B ⊥ ( λ x , λ ˜ t ) · B ⊥ [ λ ( x + r ) , λ ˜ t ] (cid:105) = λ − π πk ⊥ (2 π ) (cid:90) d r e i k ⊥ · r (cid:104) B ⊥ ( x , ˜ t ) · B ⊥ ( x + r , ˜ t ) (cid:105) = λ − M ( k ⊥ , ˜ t ) . (2.16)It can be shown that the perpendicular magnetic power spectrum satisfying the relationEq. (2.16) has the self-similar solution (Olesen 1997) M ( k ⊥ , ˜ t ) = ˜ t − / ¯ M ( k ⊥ ˜ t / ) = k ⊥ F ( k ⊥ ˜ t / ) , (2.17)where ¯ M and F are scaling functions of the variable k ⊥ ˜ t / .The solution Eq. (2.17) represents the pattern of energy decay and scale growth,determined by the physics of the merger events. However, it does not contain anyinformation about energy as a function of scale. Assuming the case of a power-law Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky spectrum, the scaling function has to satisfy ¯ M ( k ⊥ ˜ t / ) ∝ ( k ⊥ ˜ t / ) − γ where − γ is theindex of the spectrum. Eq. (2.17) can thus be further expressed as M ( k ⊥ , ˜ t ) ∝ t − / ( k ⊥ ˜ t / ) − γ ∝ ˜ t − α k − γ ⊥ , (2.18)where α = γ + 1 .
3. Numerical study
We now proceed to test the predictions of the previous section by means of directnumerical simulations with the code
Viriato (Loureiro et al. et al. ∂ψ∂t + [ φ, ψ ] = V A ∂φ∂z + η ∇ ⊥ ψ, (3.1) ∂∂t ∇ ⊥ φ + [ φ, ∇ ⊥ φ ] = V A ∂∂z ∇ ⊥ ψ + [ ψ, ∇ ⊥ ψ ] + ν ∇ ⊥ φ, (3.2)where V A is the Alfvén speed based on the (constant and uniform) guide-field, V A = B g / √ πρ , and [ A, B ] ≡ ∂ x A∂ y B − ∂ y A∂ x B . The flux function ψ is defined as ψ ≡− A z / √ πρ , where A z is the axial component of the magnetic vector potential. Thestream function φ is φ ≡ cϕ/B g where ϕ is the electrostatic potential. The magnetic field(measured in velocity units) and velocity are thus related to ψ and φ as: B = ˆ z × ∇ ψ + ˆ z B g u ⊥ = ˆ z × ∇ φ. (3.3)The viscosity ν is taken equal to the magnetic diffusivity η in all simulations.3.1. Simulation setup In Viriato , the length scales on the perpendicular plane and in the parallel directionare normalized to (arbitrary) reference lengths L ⊥ and L (cid:107) respectively ( L ⊥ (cid:28) L (cid:107) inthe RMHD limit). Times are normalized to the parallel Alfvén time τ A = L (cid:107) /V A . Inthe following description, all quantities are given in dimensionless form. The simulationdomain is a box of dimensions L x × L y × L z , where L x = L y in all the simulationspresented. We set L x = 2 π , and show the results from simulations with L z = 4 L x . Amore elongated simulation box allows more generations of mergers before the length l offlux tubes in the z direction reaches the length of the box. In order to adequately resolvethe dynamics along z , and considering the computational cost of the simulation, we findthe box with L z = 4 L x to be optimal. The box has periodic boundary conditions in alldirections.The initial condition is a straightforward extension to 3D of our previous 2Dstudy (Z19). There is no initial flow, φ ( x, y, z, t = 0) = 0 , and the magnetic flux tubesare invariant in the z direction: ψ ( x, y, z, t = 0) = ψ cos( k x ) cos( k y ) . It is easy tocheck that this initial condition satisfies the dissipationless version of Eqs. (3.1-3.2) andis therefore an (ideal) equilibrium. We thus have a k × k static array of guide-field-aligned magnetic flux tubes of alternating polarities. In all runs we set k = 8 and so theinitial radius of each flux tube is R = L x / k = π/ . We further fix ψ k = 1 to set thestrength of the initial in-plane magnetic field B ⊥ , ≡ ψ /R = 2 /π . A spatially-randomperturbation with small (linear) amplitude is added to trigger the dynamic evolution ofthe system. This perturbation breaks the axial symmetry of the system. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer H M = (cid:82) dV A · B ≈ B g (cid:82) dV ψ . We note that while the local magnetic helicity associated with each flux tubeis finite [Eq. (2.1)], the net magnetic helicity of the entire setup is zero. Since H M isan ideal invariant, it should remain approximately zero throughout the evolution of thesystem. Therefore, our simulations are relevant to the study of inverse energy transferin non-helical turbulent systems, a question briefly alluded to in the Introduction andfurther discussed in Appendix C.We set ν = η = 10 − in the simulation, corresponding to the Lundquist number ofthe initial flux tubes S ≡ R v A, /η = 1250 . The perpendicular plane is (spectrally)resolved with grid points. As per our discussion of Section 2.1, this value of S implies that we expect the reconnection of two flux tubes to fall in the Sweet-Parker(SP) regime; the thickness of the SP current layers during the first generation of mergers, δ SP, ≈ S − / R is only marginally resolved with about one cell. But we note that asthe system evolves and the flux tubes grow to larger scales, the thickness of the currentsheets increases (the Lundquist number is expected to remain constant, as we have arguedin Section 2.1, but the radius of the tubes progressively increases) and the simulationbecomes progressively better resolved. Naturally, it would be desirable to investigateinitial conditions characterized by even higher values of S , such as to access the plasmoid-mediated reconnection regime. Unfortunately, the computational requirements presentedby such simulations are too great: we estimate that achieving S (cid:38) with the samenumber of initial flux tubes ( k = 8 ) requires a numerical resolution of at least 3000grid points in each perpendicular dimension. This is too demanding even if the parallelresolution requirements were to remain unchanged.The parallel length of the box is resolved with 512 grid points. A seventh-order upwindscheme is used to advect the fields in this direction. To check the convergence with respectto parallel resolution, we performed a simulation with the same physical parameters butwith grid points instead (thus doubling the resolution in the z -direction), andevolved the system for the first few generations of mergers (during which the simulationis least resolved). The results are indistinguishable from the fiducial run with × grid points.Finally, we also perform a run in which the resistivity and viscosity are replacedwith hyper-dissipation. Specifically, the dissipative terms in Eq. (3.1) and Eq. (3.2)become η H ∇ ⊥ ψ and ν H ∇ ⊥ φ , respectively. Since the timestep with which the equationsare integrated is constantly recomputed by Viriato to maximize the efficiency of thecalculation, it follows that the values of η H and ν H are also continuously adjusted in thecode to ensure that the hyper-dissipation is reduced to the minimum possible for a givenspatial resolution (Loureiro et al. η H and ν H as constants. The corresponding effective hyper-Lundquist numberis S H, ≡ R v A, /η H ≈ × . 3.2. Visualizations
We begin the analysis of the results of the numerical simulations by showing visual-izations at different times of the run with regular dissipation, S = 1250 . The hyper-dissipative run looks qualitatively similar, except that structures are, of course, sharper.In Fig. 2 we show the vorticity ω z ≡ ∇ ⊥ φ and current density J z ≡ ∇ ⊥ ψ at varioustimes. The top panel shows ω z on the xy plane at z = 0 [i.e., ω z ( x, y, z = 0) ], and themiddle panel shows J z on that same plane [i.e., J z ( x, y, z = 0) ]. The bottom panel shows0 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky J z on the xz plane at y = 0 [i.e., J z ( x, y = 0 , z ) ], plotted at the same times as the middlepanels.These plots suggest a two-stage evolution of the system. The length of the initial fluxtubes ( π ) is much longer than the length that critical balance would dictate (about . ) and, therefore, the flux tubes are not causally connected from one end of the boxto the other. When reconnection events start at different positions along z -axis (i.e., ondifferent xy planes), the flux tubes break in the z -direction (kink-type instabilities arenot expected in our RMHD system, as we discuss in Appendix B). Due to this unstableset-up of elongated flux tubes, the system first enters a transient stage of developingturbulence. At this stage, finer structures appear and the system becomes fully turbulent(see the snapshot at t = 11 τ A ). After the transient stage, the system relaxes to a morenatural state, where the initial setup is “washed out” by the developing turbulence. Thesystem then enters a prolonged stage of decaying turbulence, during which progressivelylarger flux tubes emerge (both in the xy and in xz planes). The length of the flux tubesgrows faster than the radius, consistent with our prediction in Eq. (2.12) (This can beseen in Fig. 2, and a quantitative demonstration is shown in Section 3.6).In Fig. 3 we show current sheets, defined as the regions where | J z ( x, y, z ) | > J rms and J rms is the root mean square of J z . The increasing complexity at the transient stage ismanifested by this plot: the breaking of flux tubes leads to increasing number of currentsheets, indicating more structures and reconnection sites. This phenomenon is also shownquantitatively by the probability distribution function of current density square J z inFig. 7, which we will describe in Section 3.4.3.3. Magnetic topology
In this subsection we introduce a method to study the topology of the magnetic field;specifically, we aim to quantify the number of flux tubes on the perpendicular planes soas to make contact with our predictions of Section 2.Any flux tube crossing a given xy plane will leave a local maximum or minimum pointof the flux function on the plane (called the O-point). That is, the cross section of theflux tube on the plane has a structure of a magnetic island. Similarly, the current layersbetween two merging flux tubes will leave a saddle point on the flux function of the plane(called the X-point). Both the O-points and the X-points are the critical (or stationary)points of a scalar field — in this case, the flux function ψ ( x, y ) on the plane. In principle,on each xy plane, the number of O-points should be the same as that of X-points. Theyhelp characterize in a quantitative way the topology of the magnetic field and can beidentified using the technique which we describe next.We identify the X/O points using the method of Hessian matrix, following previouswork (Servidio et al. xy plane, we first identify all the criticalpoints of the magnetic potential, defined by ∇ ⊥ ψ = 0 , where the perpendicular magneticfield B ⊥ = 0 . Then we numerically calculate the second-order partial derivatives of theflux function ψ at the critical points, constructing the Hessian matrix H = (cid:18) ∂ xx ψ ∂ xy ψ∂ yx ψ ∂ yy ψ (cid:19) (3.4)and computing its determinant and eigenvalues. If H is positive definite at point ( x, y ) ,then the point is an O-point with local minimum; and vice-versa : a point with negativedefinite H is an O-point with local maximum. If H has one positive and one negativeeigenvalues, such point is a saddle point, called X-point. Rarely, the determinant of H iszero; the locations where this happens are called degenerate critical points.We apply this diagnostic to time slices throughout the simulation, thus obtaining the ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer z / rms - - /2 0 /2-- /20/2 t = 0 A - - /2 0 /2 t = 11 A - - /2 0 /2 t = 65 A - - /2 0 /2 t = 214 A J z / J rms - - /2 0 /2-- /20/2 t = 0 A - - /2 0 /2 t = 11 A - - /2 0 /2 t = 65 A - - /2 0 /2 t = 214 A - - /2 0 /2-4-2024 - - /2 0 /2 - - /2 0 /2 - - /2 0 /2 Figure 2: Visualizations for the run with S = 1250 : vorticity density ω z (normalized byits root mean square value) at various times on the xy planes at z = 0 (top); currentdensity J z (normalized by its root mean square value) at various times on the xy planesat z = 0 (middle) and on the xz planes at y = 0 (bottom).2 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky
Figure 3: Current density contours selected according to | J z | > J rms at various times forthe run with S = 1250 . ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer J z (colors) and magnetic flux ψ (contours) on the xy plane at z = L z , at t = 100 τ A for the run with S = 1250 . Green circles (‘ ◦ ’ symbols) identifythe positions of O-points; Green crosses (‘ × ’ symbols) mark the positions of X-points.
10 100 t / A N X-pointsO-pointsindex -1
Figure 5: Time evolution of the number of X and O points, N ⊥ , for S = 1250 . The datapoints show the values of the number of X/O-points averaged over eight xy planes foreach timeslice. The shadowed area shows the variation between different xy planes. A t − power law, predicted in Section 2, is shown for reference.time evolution of the number of X/O-points, N ⊥ ( t ) . For each time slice, we examine eight xy planes spread evenly along the z -direction. An example is shown in Fig. 4. Lastly,the number of X/O-points is averaged over the eight planes (all xy planes should bestatistically equivalent).Fig. 5 shows the time evolution of N ⊥ for the run with S = 1250 . The data points arethe number of X/O-points averaged over eight xy planes and the shadowed area showsthe variation between different xy planes. Initially there are ×
16 = 256
X/O points,as specified by our initial condition. During the transient stage (between ∼ τ A − τ A ),the number of X/O points first increases by roughly a factor of 2, and then rapidlydecreases. During this stage, the increase in the number of X/O-points is correlated tothe deformation of flux tubes and the formation of current sheets (Fig. 2 and Fig. 3).4 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky
This indicates the formation of reconnection sites at multiple places and leads to rapidmagnetic energy dissipation, as we will discuss in Section 3.4.After the transient stage, the decay of N ⊥ slows down and scales with time as N ⊥ ∼ t − , as it predicted by our hierarchical model [Eq. (2.13)]. This suggests that the late-timeevolution of the system is indeed governed by the merger of magnetic flux tubes.3.4. Energy decay
From the visualizations and the number of magnetic structures, we can see that thesystem first goes through a transient stage of developing turbulence. This is followedby a prolonged stage of decaying turbulence, with the formation of progressively largermagnetic structures. The behavior of magnetic energy evolution, as we will now describe,is consistent with this physical picture.Fig. 6 shows the time evolution of magnetic and kinetic energies for the run with S = 1250 . Rapid magnetic energy decay occurs during the transient stage (around 6 τ A to τ A ), corresponding to the increase, followed by rapid decrease, of the number of X/Opoints (Fig. 5). Also, it corresponds to the stage, as we observe from the visualizations,of increasing complexity with the breaking of flux tubes and the emergence of finerstructures. Therefore, we attribute the rapid magnetic energy dissipation during thisstage to the formation of a large amount of reconnection sites with localized currentsheets among the twisted and deformed structures in this complex system. After thisstage ( t (cid:38) τ A ), the magnetic energy decay clearly follows a t − scaling, as predictedby our hierarchical model, Eq. (2.11). It coincides with the stage when the number ofX/O points exhibits the same scaling, N ⊥ ∼ t − (Fig. 5), also in accordance with ourpredictions, Eq. (2.13). These numerical confirmations of our predictions suggest thatindeed the key process underlying the evolution of our system at this stage is the mergerof flux tubes, which appears to proceed as we envisionsed in Section 2.Further evidence to support the claims on enhanced energy dissipation during thetransient stage is provided in Fig. 7, which shows the compensated probability distri-bution function (PDF) of current density square, J z , normalized to its mean square, J , at different moments of time [in resistive RMHD with constant resistivity, themagnetic energy dissipation (ohmic heating) rate per unit volume is proportional to J z ].We arbitrarily choose one snapshot at the transient stage ( t = 11 τ A ) and two from theself-similar stage ( t = 65 τ A and t = 214 τ A ). The PDF of the two later snapshots (duringthe self-similar stage) are similar, while the PDF at t = 11 τ A has a clear tail at thehigh current density end. This suggests that, during the transient stage, there are morelocalized current sheets facilitating the rapid dissipation of magnetic energy.The kinetic energy, starting from zero, rapidly increases as the flux tubes start tointeract. It then traces the magnetic energy and undergoes the transient rapid-decaystage and t − self-similar stage. Comparing the vorticity density ω z and the currentdensity J z on the same xy planes at the same moments of time (Fig. 2), we find that theflow in the system is dominated by the reconnection outflows with the Alfvén speed basedon the perpendicular magnetic field. Therefore, when the magnetic energy decays throughthe merger of flux tubes, the kinetic energy follows magnetic energy and decays as t − .During the self-similar stage, a small oscillating exchange between magnetic energy andkinetic energy is observed, which we will discuss in Section 3.5.The rapid magnetic energy decay that we observe in the transient stage bears sim-ilarities to that observed in the numerical study of magnetized jets from active galax-ies (O’Neill et al. et al. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer t / A E M E k index -1 Figure 6: Time evolution of total magnetic energy ( E M ) and kinetic energy ( E K ) for S = 1250 . A t − power law is shown for reference, in accordance with our prediction,Eq. (2.11). J z / J rms J z F ( J z ) A A A Figure 7: Compensated probability distribution function (PDF) of J z /J for S = 1250 at three moments of time. The tail at large J z /J for the PDF at t = 11 τ A indicatesa greater number of localized current sheets in the transient stage.the dissipation of magnetic energy. In our system with a strong guide field, we do notexpect kink-type instabilities (see discussion in Appendix B). However, the configurationof a large ensemble of flux tubes has more degrees of freedom to form reconnection sitesvia their twisting and deformation.3.5. Spectra and related quantities
We now proceed to investigate the time evolution of the magnetic and kinetic energyspectra. In order to obtain the clearest spectra possible, we analyze the run with hyper-dissipation instead of normal (Laplacian) diffusivity. For reference, the time evolutionof energies for this run is shown in Fig. 8. Comparison with Fig. 6 shows qualitativelysimilar behavior between the two runs, although two differences stand out: (i) in thehyper-dissipative run the kinetic energy remains about a factor of 2 below the magneticenergy, whereas in the S = 1250 resistive run they have comparable magnitudes inthe self-similar part of the evolution (i.e., for t/τ A (cid:38) ); and (ii) the amplitude of the6 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky t / A E M E k index -1 Figure 8: Time evolution of total magnetic ( E M ) and kinetic ( E K ) energies for the runwith hyper-dissipation.oscillations in the late evolution of the energies is significantly reduced in the hyper-dissipative case. We will first discuss the energy spectra and then turn to the reasons forthese differences at the end of this subsection.Generally speaking, the magnetic energy in the system is contained not only inmagnetic flux tubes but also in Alfvén wave packets launched by the convective motion offlux tubes in the perpendicular plane at the scale of flux tubes. Magnetic energy containedin flux tubes is transferred to larger scales through the flux-tube mergers, while that con-tained in Alfvén wave packets cascades to smaller scales through nonlinear interactions.In the Alfvén wave packets, the kinetic energy is in equipartition with the magneticenergy, cascading to smaller scales. A k − / spectrum is expected for both the magneticand the kinetic energy as a result of the turbulent direct cascade (Boldyrev 2006). Thereis another ingredient in the system that can potentially contribute significantly to themagnetic and kinetic energy spectrum — the thin current sheets between merging fluxtubes. The sharp magnetic field reversals at the current sheets are expected to yield a k − magnetic energy spectrum (Burgers 1948). The Alfvénic outflows from reconnectionlayers, whose spatial profile is derived in Loureiro et al. (2013), are expected to yield a flatkinetic energy spectrum. That is, the thin current sheets in the system can potentiallysteepen the turbulent magnetic energy spectrum, and flatten the turbulent kinetic energyspectrum.With these considerations borne in mind, let us analyze the numerical results. Theperpendicular magnetic energy spectrum M ( k ⊥ , t ) at different moments of time is shownin Fig. 9 (left panel). The initial spectrum peaks at the scale corresponding to the crosssection of the initial flux tubes. As the flux tubes merge, the peak of the spectrum shiftsto larger scales, until it becomes limited by the box size. An inertial range appears tothe right of the peak, with a power-law exponent γ ≈ at earlier times, and γ ≈ / at later times. The kinetic spectrum, shown in Fig. 10, does not have a steady powerlaw. At later times, while the magnetic spectrum peaks at k ⊥ = 4 , corresponding to thesize of islands, the kinetic spectrum peaks at k ⊥ = 2 , corresponding to the large-scalemotion of the islands. Such motion injects kinetic energy into the system, potentiallyyielding a power-law inertial range at smaller scales, which we tentatively fit with k − / ⊥ .However, at earlier times, the kinetic spectrum follows a shallower power-law ∼ k − ⊥ . ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer k M ( k , t ) k k k / k , max ( t )10 M ( k , t ) / M m a x ( t ) t = 0 A t = 10 A t = 20 A t = 30 A t = 40 A t = 60 A t = 80 A Figure 9: Raw (left) and normalized (right) perpendicular magnetic power spectra, forsimulation with hyper-dissipation. Dotted and dashed lines indicate reference slopes of k − ⊥ and k − / ⊥ , respectively. k K ( k , t ) k k t = 10 A t = 20 A t = 30 A t = 40 A t = 60 A t = 80 A Figure 10: Perpendicular kinetic power spectra for the simulation with hyper-dissipation. k − ⊥ and k − / ⊥ slopes are shown for reference.While we do not have an explanation for the − index, we think it is consistent with theobservation that at earlier times the system has many localized thin current sheets, wherethe kinetic energy is dominated by the (Alfvénic) outflows. The outflows have a spatialprofile yielding a flat spectrum as we described above, and act as a source of kineticenergy at all scales. As the flux tubes merge and grow to larger scales, the thickness ofthe current sheets also grows and the magnetic reversal becomes less sharp, resulting ina weaker flattening effect on the spectrum. The kinetic spectrum thus becomes steeperat later times. At early times (before the turbulence is fully developed), the magneticspectrum follows a k − ⊥ spectrum resulting from the sharp field reversals (Burgers 1948).At later times, as the turbulent cascade becomes stronger, the k − ⊥ spectrum is replacedby a shallower k − / ⊥ turbulent spectrum. Interestingly, in our 2D study, a k − ⊥ magneticenergy spectrum was maintained throughout the evolution [Fig. 4 in Z19]. The reasonfor this difference is that, in the 2D system, the magnetic energy significantly dominatesthe kinetic energy, making the system less turbulent than its 3D counterpart. The directturbulent cascade is too weak to form a k − / spectrum, and therefore, the total magneticenergy spectrum in our 2D system is mostly determined by the current sheets, yieldinga k − ⊥ spectrum.Our analytical model of Section 2 leads to the prediction that the long-time evolution8 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky of the magnetic energy spectrum should be self-similar, Eq. (2.18). To test this prediction,we show in the right panel of Fig. 9 the spectra normalized to their instantaneous maxi-mum value, M max ( t ) , as a function of wavenumber normalized to the values k ⊥ , max ( t ) atwhich M max ( t ) are achieved. We observe that, except for the transient-stage spectrum at t = 10 τ A , the spectra thus normalized approximately collapse onto the same distribution,and that this collapse becomes progressively better as time increases. Mathematically,this behavior is expressed as M ( k ⊥ , t ) = M max ( t ) ¯ M ( k ⊥ /k ⊥ , max ) , where ¯ M ( k ⊥ /k ⊥ , max ) is the static universal distribution to which all the curves collapse after normalization.We also observe that after the transient stage, the time evolution of k ⊥ , max and M max approximately follow the power-law-in-time behavior: k ⊥ , max ∝ t − / , and M max ∝ t − / — see the second and third panels in Fig. 11. Thus, the magnetic energy spectrumcan be expressed as M ( k ⊥ , t ) ∝ t − / ¯ M ( k ⊥ t / ) , in agreement with the prediction ofEq. (2.18). Furthermore, since the scaling function is itself a power law in the inertialrange, ¯ M ( k ⊥ t / ) ∝ k − / ⊥ ( γ = 3 / ), according to our prediction, the time evolution ofspectral magnetic energy density at any fixed wavenumber should behave as M k ⊥ ( t ) ∝ t − α , where α = ( γ + 1) / / . This is indeed observed in the bottom panel of Fig. 11,where after the transient stage, M k ⊥ ( t ) (for moderate wavenumbers) approaches the t − / power law. At the largest scale the energy keeps building up, explicitly showingthe inverse transfer of magnetic energy. From the measured scaling laws described above,we conclude that the self-similar solution of the magnetic power spectrum, based on ourhierarchical analytical model, is confirmed by our numerical results.Complementary to the time evolution of k ⊥ , max , the growth of the perpendicularcoherence length of our system is illustrated by computing the magnetic energy containingscale, defined as ξ M ( t ) ≡ (cid:82) πk − ⊥ M ( k ⊥ , t ) dk ⊥ (cid:82) M ( k ⊥ , t ) dk ⊥ . (3.5)As shown in the top panel of Fig. 11, after the transient stage, ξ M increases with timeroughly following the t / scaling, consistent with our model, until the later stage ( t (cid:38) τ A ) during which the parallel length of flux tubes becomes limited by the length ofthe box.To characterize the system’s dynamics along the guide field, we show in Fig. 12 themagnetic power spectrum as a function of k z . A broad energy distribution develops oncethe flux tubes break in the z direction. As flux tubes merge, magnetic energy at high k z rapidly decreases, indicating the increasing coherence length in z direction. This isconsistent with the measurement of magnetic structure functions (Fig. 13), showing therapid growth of the characteristic flux tube, as we explain in Section 3.6.Finally, let us address the issue of the oscillatory energy exchange between the magneticfield and the flow. From the time evolution of magnetic energy density at differentwavenumbers k ⊥ (Fig. 11, bottom panel), we see that energy oscillation mostly happens atthe smallest k ⊥ , corresponding to system-size magnetic flux tubes. These tubes cannotmerge and grow into larger scales due to the system-size limit (in a sense, they havereached a final ‘equilibrium’ consisting of two flux tubes of opposite polarities). Instead,they exchange energy with the flow — possibly both in the form of Alfvén waves, andbouncing motions: background flows can push these flux tubes against each other, leadingto their (incompressible) deformation. This increases magnetic field line tension andresults in a restitution force leading to an oscillation whose amplitude, therefore, dependson the ratio of the kinetic-to-magnetic energy. Since the magnetic energy is larger inthe hyper-dissipative case than in the Laplacian resistivity case (simply because hyper- ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer M t k , m a x t M m a x ( t ) t t / A M k ( t ) t k=1k=3 k=5k=10 k=15k=20 Figure 11: Time evolution of ξ M (top), k ⊥ , max (second), M max (third) and M k ⊥ (forselected values of k ; bottom), for the simulation with hyper-dissipation.resistivity leads to better conservation of magnetic energy), the hyper-dissipation casehas smaller amplitude oscillations.3.6. Magnetic structure functions
In order to obtain statistical information on the growth of flux tubes, and to test theassumption of critical balance, we calculate the second-order structure function of themagnetic field B ( x , t ) , defined as S (2) B ( δ x , t ) ≡ (cid:104)| B ( x , t ) − B ( x , t ) | (cid:105) x , (3.6)where δ x ≡ x − x and the angle brackets indicate the average over the position x in the simulation domain. We define the local mean field as B m ≡ [ B ( x ) + B ( x )] / and the local unit magnetic field vector ˆ b m ≡ B m / | B m | represents its direction. Thesedefinitions ensure that the structure function is computed with respect to the local meanmagnetic field [as it should, if it is to be used to gauge whether the turbulence we observeis critically balanced (Cho & Vishniac 2000; Mallet et al. δx ⊥ , δx (cid:107) ) can thus be specified by δx ⊥ ≡ | ˆ b m × δ x | , δx (cid:107) ≡ ˆ b m · δ x . (3.7)The correlation lengths of the magnetic field parallel and perpendicular to the localmean field, which we interpret as the statistical length and radius of the flux tubes0 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky k z + 1 M ( k z ) t = 0 A t = 20 A t = 40 A t = 60 A t = 80 A Figure 12: Parallel (to the guide field B g ) magnetic power spectra at different momentsof time for the simulation with hyper-dissipation. x x . . . . t=20.7 A x . . . . t=40.7 A x . . . . . . t=60.1 A x . . . . . . . t=80.5 A Figure 13: Contours of normalized (to its maximum value) magnetic structure function S (2) B ( δx ⊥ , δx (cid:107) ) at various times for the run with hyper-dissipation.respectively, are thus expressed by the contours of the structure function S (2) B ( δx ⊥ , δx (cid:107) ) on the ( δx ⊥ , δx (cid:107) ) plane. This is shown in Fig. 13 (normalized to the instantaneousmaximum value of the structure function). We then define the characteristic length l andthe radius R of the flux tubes using the intersection of the contour line S (2) B /S (2)B , max = 0 . (roughly the largest value at which the contour lines are not affected by noise) withthe δx (cid:107) axis and the δx ⊥ axis, respectively. Note that because of periodic boundaryconditions, the largest values of δx ⊥ and δx (cid:107) , and hence maximum radius and lengthof a flux tube that can be attained in our simulations, are roughly L x / √ √ π and L z / π , respectively. The growth of the typical flux tube can be clearly observed asthe system evolves with time. The growth of parallel length is faster than that in theperpendicular direction, which agrees with our discussion in Section 2 that the aspectratio of a flux tube increases with time ( l/R ∼ ˜ t / ).In order to characterize the time evolution of magnetic flux tubes, we take eightsnapshots between t = 10 τ A and t = 80 τ A (roughly every τ A ). For each snapshot, ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer t / A l l = 0.14 t / A R / B l l = R / B Figure 14: Time evolution of l [left panel, dashed line represents a t power law, Eq. (2.5)]and scaling of l versus R/ ¯ B ⊥ [right panel, dashed line represents the critical balancescaing, Eq. (2.4)] for the run with hyper-dissipation, where l and R are the length andradius of the characteristic flux tube respectively.we compute the values of l and R . The time evolution of R (not shown) roughly followsa t / scaling (similar to the time evolution of ξ M shown in Figure 11), whereas thetime evolution of l scales linearly with time (left panel of Figure 14), in accordancewith the prediction of Eq. (2.12). We also compute the root-mean-square value of theperpendicular magnetic field, ¯ B ⊥ , and test if it satisfies the critical balance scaling, l/ ˆ B g ∼ R/ ¯ B ⊥ , where ˆ B g = B g L ⊥ /L (cid:107) = 1 is the guide field normalized in the sameway as used to compute the structure function. We plot l versus R/ ¯ B ⊥ for these eightsnapshots in the right panel of Figure 14. As seen, the critical balance condition isindeed verified (except for the last snapshot, whereupon the parallel length of flux tubeshas reached the length of the simulation box, and the critical balance scaling is thusexpected to break). This confirms our assumption made in Section 2 that in our systemof merging flux tubes, the aspect ratio of the flux tubes, l/R , is determined by criticalbalance. 3.7. Energy transfer function
In this subsection we use shell-to-shell energy transfer functions to quantify the transferbetween magnetic and kinetic energy at different scales. This enables us to identify directand inverse energy transfer in our simulations. This technique has been used in studiesof MHD and kinetic turbulence [e.g., Alexakis et al. (2005); Told et al. (2015)], and wemade modifications for its use in the RMHD system.We define the shell-filtered variables as B K ⊥ ( x ⊥ , z ) = (cid:88) K Figure 15: Averaged energy transfer from shell Q to shell K at t = 10 τ A (first stage) and t = 60 τ A (second stage) for the run with hyper-dissipation. Top: magnetic-to-magnetictransfer T bb ( Q, K ) . Bottom: kinetic-to-magnetic transfer T ub ( Q, K ) .interested in the time evolution of magnetic energy in each shell K : ∂ t (cid:28)(cid:90) d x ⊥ 12 ( B K ⊥ ) (cid:29) z = (cid:42)(cid:90) d x ⊥ (cid:88) Q (cid:104) B K ⊥ · ( B ⊥ · ∇ ⊥ ) u Q ⊥ − B K ⊥ · ( u ⊥ · ∇ ⊥ ) B Q ⊥ (cid:105) − η |∇ ⊥ B K | (cid:43) z . (3.10)The transfer of magnetic energy from shell Q to shell K is T bb ( Q, K ) = −(cid:104) (cid:90) d x ⊥ B K ⊥ · ( u ⊥ · ∇ ⊥ ) B Q ⊥ (cid:105) z , (3.11)and the transfer of kinetic energy in shell Q to magnetic energy in shell K is T ub ( Q, K ) = (cid:104) (cid:90) d x ⊥ B K ⊥ · ( B ⊥ · ∇ ⊥ ) u Q ⊥ (cid:105) z . (3.12)We perform the energy transfer function diagnostic on the results of the run withhyper-dissipation. We arbitrarily choose one snapshot during the first stage of developingturbulence ( t = 10 τ A ) and another during the second stage of decaying turbulence ( t =60 τ A ) to illustrate the behaviour of shell-to-shell energy transfer, as shown in Figure 15.To guide the reader, we note that red color below (above) the diagonal in these plotsdenotes an inverse (direct) energy transfer; and vice-versa for the blue color. We firstfocus on the contour plot of T bb ( Q, K ) at t = 10 τ A , shown in the top-left panel. At thismoment, the magnetic energy spectrum shown in Figure 9 peaks at k ⊥ = 8 . There is ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer 23a strong direct energy transfer from shells Q = 7 , to shells K > Q , indicating thenon-local direct transfer characterizing the formation of thin current sheets (consistentwith the k − ⊥ magnetic energy spectrum as we explained earlier, shown in the left panelof Figure 9). There is also a weak but statistically significant inverse transfer of magneticenergy from shells Q = 7 , to shells K = 5 − (corresponding to the peak of magneticspectrum). This suggests that the inverse transfer happens at the scale of flux-tubes andis caused by flux-tube merger. Due to the development of turbulence, T bb is dominatedby the direct transfer from small- Q shells to large- K shells overall. At this stage, thereis also a transfer of magnetic energy to kinetic energy close to the diagonal, as shown inthe contour plot of T ub in the bottom-left panel, consistent with the formation of plasmaflow.At t = 60 τ A , flux tubes have reached larger scales (with a peak at around k ⊥ = 4 ).Shells K, Q > are dominated by the strong direct transfer between a wide rangeof scales. However, although relatively weak in strength, there is also a finite inversemagnetic energy transfer from shells Q > to shells K = 1 − shown by the redhorizontal band along the lower boundary of the contour plot of T bb in the top-rightpanel. From the T ub plot in the bottom-right panel, we observe a mild transfer frommagnetic energy in shells K = 1 , to kinetic energy in shells Q > . For K > , kineticenergy in shells Q is transferred to magnetic energy in shells K if Q < K , while for Q > K , the transfer between kinetic and magnetic energy happens in both directions.In summary, the energy transfer diagnostics paint a complex picture that is consistentwith the analytical model and numerical results discussed earlier. Both inverse and directenergy transfer occurs: the former driven by flux-tube coalescence; the latter by theturbulent dynamics that arises. While the direct energy transfer is dominant, it is theinverse transfer that progressively delivers energy to larger and larger scales until reachingthe simulation domain size. 4. Conclusion In this work, we investigate the inverse transfer of magnetic energy in 3D MHDturbulence, and propose the merger of magnetic flux tubes through magnetic reconnectionas its underlying physical mechanism. An analytical model is developed in the 3D RMHDframework, based on the hierarchical merger of flux tubes. Perpendicular to the strongguide field, the dynamics are found to be identical to those in the analogous 2D system.The perpendicular magnetic field energy decays as ˜ t − and the perpendicular correlationlength grows as ˜ t / , where ˜ t is the time normalized to the reconnection time scale.Parallel to the guide field, the dynamics are determined by critical balance, setting theaspect ratio of the flux tubes. The power-law time evolution of physical quantities leadsto the self-similar evolution of the perpendicular magnetic spectrum M ( k ⊥ , ˜ t ) ∝ ˜ t − α k − γ ⊥ ,where α = γ + 1 .Our analytical model is confirmed by direct numerical solution of 3D incompressibleRMHD equations. The initial equilibrium (without plasma flow) is set up with anensemble of volume-filling flux tubes of alternating polarities aligned with the guide field.Given a small perturbation, the system undergoes a two-stage evolution. The system firstenters a transient stage of developing turbulence, during which the flux tubes break inthe parallel direction, leading to the development of turbulence and enhanced magneticenergy dissipation. A prolonged self-similar stage of decaying turbulence ensues, duringwhich the evolution of various quantities is captured by our analytical model and theassumption of critical balance is confirmed by the measurement of 2nd-order magneticstructure function.4 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky Our 3D system shares similarities and differences with its 2D counterpart. The presenceof a strong guide field (underlying the RMHD formalism that we employ in our simula-tions) implies that the merger of flux tubes in the perpendicular planes is identical in 3Dand in 2D. That is, in both cases the mergers are well described by rules stemming fromthe conservation of mass and poloidal flux (or magnetic potential). A novel aspect of ourwork is the observation that, in the 3D case, critical balance determines the relationshipbetween perpendicular and parallel flux-tube dynamics. The 3D system tends to have alarger kinetic-to-magnetic energy ratio than the 2D version (for the same initial Lundquistnumber S ), mainly due to the flow arising in the initial stage of developing turbulenceand the Alfvén waves propagating along the guide field. This renders the 3D system moreturbulent than 2D. The k − magnetic energy spectrum appears in 2D systems, causedby the sharp field reversals at thin current sheets. Although such current sheets are alsopresent in the 3D simulations, this effect is masked by a k − / magnetic spectrum causedby the direct turbulent cascade.This work lies exclusively in the (reduced) MHD framework. In many realistic as-trophysical situations the magnetic field generation and early-stage evolution are mostlikely happening in collisionless plasmas without background magnetic fields. One maythus question the applicability of our results to such situations. While we do think thatdirect applicability may not be warranted, we also note that our analytical model forstatistical flux tube mergers stems from the conservation of quantities that we wouldexpect to remain conserved in the collisionless regime in the absence of backgroundmagnetic fields. We thus think that the ideas we have proposed here should form thebasis of the corresponding collisionless model, which we defer to future study.One may also wonder if the scaling laws of inverse energy transfer — and magneticreconnection as the process that enables such inverse transfer — in systems with astrong guide field remain valid in isotropic systems (i.e., systems without a guide field),as investigated by Brandenburg et al. (2015) and Zrake (2014). This is the subject ofongoing investigations (Bhat et al. Acknowledgements. This work was supported by NSF CAREER award No. 1654168(MZ and NFL), NSF grants AST-1411879 and AST-1806084, and NASA ATP grantsNNX16AB28G and NNX17AK57G (DAU). The authors thank Alexander Schekochi-hin, Annick Pouquet, Pallavi Bhat, Abhay Ram, Ryan White and Greg Werner forinsightful discussions. This research used resources of the MIT-PSFC partition of theEngaging cluster at the MGHPCC facility, funded by DOE award No. DE-FG02-91-ER54109 and the National Energy Research Scientific Computing Center (NERSC)–a USDepartment of Energy Office of Science User Facility operated under Contract No. DE-AC02-05CH11231. This work also used the Extreme Science and Engineering DiscoveryEnvironment (XSEDE), which is supported by National Science Foundation Grant No.ACI-1548562. This work used the XSEDE supercomputer Stampede2 at the TexasAdvanced Computer Center (TACC) through allocation No. TG-PHY140041 (Towns et al. Appendix A. Interaction between inverse magnetic transfer andambient turbulence in the context of galacticmagnetogenesis In this appendix we estimate the scale where the magnetic seed fields can be pickedup and amplified by the turbulent dynamo. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer L and cascade to smaller scales. Magnetic seed fieldsare generated at the small scale l seed and inverse transfer to larger scales.We consider an unmagnetized hydro-turbulent astrophysical system where there is anasymptotically large scale separation between the injection of kinetic energy at largescales and the generation of magnetic seed fields at small scales. A Kolmogorov cascadeof kinetic energy is assumed, resulting in a k − / inertial range of kinetic energy spectrum E K ( k ) , starting from the large driving scale L down to the intermediate scales we areinterested in. The evolution of magnetic seed fields is assumed to follow the scalinglaws we derived in Section 2 [Eq.(2.11) and Eq.(2.12)]. We consider the seed fields tobe generated at a small scale l seed , and described by a time-depending characteristicwavenumber k (˜ t ) ∼ ˜ t − / as they merge and grow into larger scales. The magneticenergy density of seed fields associated with its instantaneous characteristic wavenumberevolves as E M (˜ t ) ∼ ˜ t − , therefore, the magnetic energy density scales as E M ( k ) ∼ k .Fig. 16 shows a cartoon of the situation we have in mind: the kinetic and magneticenergy spectra resulting from the direct cascade and inverse transfer, respectively. Kineticenergy dominates at the large scales while magnetic energy dominates at small scales.In this appendix, we provide an estimation of the critical scale l c where the kinetic andmagnetic energy are comparable. At such scale, the ambient turbulent flow is dynamicallyimportant to the evolution of magnetic seed fields, and can act as a dynamo to amplifythe magnetic field.The kinetic energy at the critical scale l c , E K ( l c ) , and at the driving scale L , E K ( L ) ,are related by the E K ( k ) ∼ k − / spectrum, while the magnetic energy at l c , E M ( l c ) ,and at the initial scale l seed , E M ( l seed ) , are related by the E M ( k ) ∼ k spectrum: E K ( l c ) = E k ( L ) (cid:18) l c L (cid:19) / , E M ( l c ) = E M ( l seed ) (cid:18) l c l seed (cid:19) − . (A 1)The kinetic energy at the driving scale L and magnetic energy at the initial scale l seed can both be related to the plasma thermal pressure. At the driving scale L , we have E K ( L ) ∼ M P , where M ∼ is the Mach number of the turbulence. At the seedfield generation scale l seed , we have E M ( l seed ) ∼ (cid:15) B P , where (cid:15) B is a parameter of thesaturation amplitude of the instability, and (cid:15) B ∼ . if we consider the magnetic seedfields generated by the Weibel instability (Silva et al. Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky The critical scale l c calculated from the condition E K ( l c ) ∼ E M ( l c ) is: l c ∼ l / L / (cid:16) (cid:15) B M (cid:17) / , (A 2)and the magnetic energy at this critical scale l c is: E M ( l c ) ∼ (cid:18) l seed L (cid:19) / P (cid:15) / B M / . (A 3)From Eq. (A 2) we learn that the critical scale l c is close to the geometric mean betweenthe driving large scale L and the microscopic seed-field scale l seed . We also note thatthe dependence of l c on the equipartition parameter (cid:15) B , which is determined by themechanisms of seed field generation, is very weak. The dependence of l c on the Machnumber M , which is determined by various astrophysical events injecting kinetic energyinto the system, is also weak. Therefore, we think this estimation of the critical scale l c [Eq. (A 2)] should be a robust result, since it does not strongly depend on the astrophysicalprocesses that generate the turbulence and magnetic seed fields, but on the genericfeatures of the turbulent energy cascade, and the scaling laws for inverse energy transfer[Eq.(2.11) and Eq.(2.12)], which are derived based on the fundamental conservation laws.The scale of magnetic seed fields l seed generated by the Weibel instability should be ofthe order of electron (or ion) skin depth d e (or d i ). For a typical galaxy, L ∼ parsec ∼ cm, d e ∼ cm and d i ∼ cm (Schekochihin et al. l seed ∼ d e , we obtain l c ∼ cm ∼ d e (if l seed ∼ d i , we obtain l c ∼ cm ∼ . AU). Thestrength of magnetic seed fields B ( l seed ) ∼ µG [about 3 µG magnetic field correspondsto equipartition with the plasma pressure, where electron density n e ∼ cm − andtemperature T e ∼ eV are assumed in the interstellar medium (ISM)]. We thus obtainthe strength of magnetic field at the critical scale B ( l c ) ∼ B ( l seed )( l seed /l c ) ∼ − Gfor l seed ∼ d e [ B ( l c ) ∼ − G for l seed ∼ d i ]. We note that the strength of this magneticfield is significantly larger than that generated by the Biermann battery effect, while itsscale is deep in the MHD range (since l c is approximately the geometric mean of thedriving scale and plasma skin depth).Our estimation is relevant for the problem of Galactic dynamo. Imagine we start witha completely unmagnetized and hydro-turbulent ISM, and a shock is launched throughthe ISM by an astrophysical event. The shock will produce microscopic filaments withtypical radius of cm [on the order of ion skin depth (Spitkovsky 2008; Kato & Takabe2008)], and with the amplitude of magnetic seed fields of the order of microgauss. Thesemagnetic filaments would at first undergo a self evolution, without being affected by theturbulent flow at larger scales. During this stage, filaments will quickly merge with eachother through magnetic reconnection, increasing their length scale and decreasing thestrength of magnetic field. This merging process would allow the filaments to deliver themagnetic field of the strength of − G to the scale of cm, where the turbulentflow is dynamically important. This can be viewed as the starting point for the turbulentMHD dynamo. That is, the above gives an estimate for the scale and strength of theinitial seed field for the Galactic dynamo problem. Appendix B. Current sheet stability to ideal kink-type modes Our treatment and discussion of the dynamical evolution of an ensemble of flux tubesmade no reference to the possible role of ideal kink-type instabilities. In part this isjustified by the fact that, in RMHD, a flux tube (more precisely, a screw pinch) witha circular cross-section is stable to the ideal internal kink mode (Hazeltine & Meiss ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer ψ ( x ) and no flow, φ ( x ) = 0 .Consider, therefore, the ideal RMHD system of equations, and linearize them aroundthe equilibrium ψ ( x ) , whose exact functional shape need not be specified for this proof.We obtain: γψ − ik y φ ∂ x ψ = ik z V A φ , (B 1) γ ( ∂ x − k y ) φ = [ ψ , ∇ ⊥ ψ ] + [ ψ , ∇ ⊥ ψ ] + ik z V A ( ∂ x − k y ) ψ . (B 2)The combination of these two relations yields the following equation for the perturbedstream function φ : (cid:2) γ + ( k z V A + k y ∂ x ψ ) (cid:3) ( ∂ x − k y ) φ + 2 k y ( k z V A + k y ∂ x ψ ) ∂ x ψ ∂ x φ = 0 . (B 3)It can be shown that, as expected, Eq. (B 3) is self-adjoint, since it is derived from theideal MHD equations possessing the self-adjoint properties. That is, it can be written inthe form L φ ≡ ∂ x [ f ( x ) ∂ x φ ] + g ( x ) φ = 0 , (B 4)where L is a self-adjoint operator and f ( x ) = γ + ( k z V A + k y ∂ x ψ ) , (B 5) g ( x ) = − k y f ( x ) . (B 6)Due to the properties of self-adjoint operators, the eigenvalues γ can only be eitherreal or imaginary. If γ < , the modes are purely oscillatory, and if γ > , the modesare purely growing or damped. The self-adjoint operator L satisfies (cid:82) ∞−∞ φ ∗ L φ dx = 0 .Therefore Eq. (B 3) yields: (cid:90) ∞−∞ φ ∗ (cid:2) ∂ x ( f ( x ) ∂ x φ ) − k y f ( x ) φ (cid:3) dx = 0 . (B 7)Using the boundary condition φ → as x → ±∞ , we obtain: (cid:90) ∞−∞ (cid:0) f ( x ) ∂ x φ ∗ ∂ x φ + k y f ( x ) φ ∗ φ (cid:1) dx = 0 . (B 8)Equivalently, Eq. (B 8) can be written as an expression for γ : γ = − (cid:82) ∞−∞ ( k z V A + k y ∂ x ψ ) (cid:0) ∂ x φ ∗ ∂ x φ + k y φ ∗ φ (cid:1) dx (cid:82) ∞−∞ (cid:0) ∂ x φ ∗ ∂ x φ + k y φ ∗ φ (cid:1) dx . (B 9)This expression makes the fact that γ < explicit for any nontrivial solution φ ( x ) . I.e.,any admissible modes are purely oscillatory: essentially, shear Alfvén waves propagatingin a nontrivial magnetic field. The specific values of γ will depend on the functionalform of ψ ( x ) , but not their sign.Based on this analysis, we conclude that two representative structures of the systemthat we study in this paper — cylindrical flux tubes and the current sheets betweenthem — are stable to ideal kink-type modes. This is not, of course, a proof that thekink instability plays no role in the evolution of our system. However, in addition, wecall attention to the fact that, at least at the twiddle-algebra level, the critical balance8 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky condition, which our system does seem to obey rather well, is the condition for kinkstability. Altogether, we think these are reasonably convincing arguments for ignoringideal kink-type modes in our analysis. Appendix C. Inverse energy transfer of random magnetic fields The inverse transfer of magnetic energy and the formation of large-scale structures in asystem consisting of volume-filling flux tubes can be explained by our hierarchical modelbased on the merger of flux tubes, as we demonstrate in the main part of this paper.However, it is possible that this fundamental physical process — merger of magneticstructures through magnetic reconnection — can be a more general phenomenon. It iswell known that in 3D MHD turbulence an inverse helicity cascade and the correspond-ing inverse energy transfer occurs due to the conservation of (non-zero) net magnetichelicity (Pouquet 1978; Christensson et al. et al. (2015) and Zrake (2014) numerically demonstrated that an inverse energy transfer existseven when the net magnetic helicity is zero, in 3D decaying non-relativistic and relativisticturbulence respectively. A careful parameter study of this problem can be found inReppin & Banerjee (2017). Bhat et al. (2020) argue that this non-helical inverse transferphenomenon might be explained by the merger of magnetic structures in a qualitativelysimilar way to the process that we have described in this paper and in Z19.To further explore and substantiate this idea – and, more generally, to demonstratethat the inverse transfer that we observe is a generic phenomenon, and not something thatdepends on the specific form of the initial magnetic field — we present in this appendixa simulation which differs from the ones analyzed in the main body of the paper inthat the initial condition is instead a random magnetic field. Specifically, we initializea Gaussian-random magnetic field, with the perpendicular magnetic spectrum narrowlypeaked around k ⊥ , = 20 . In the parallel direction we impose sinusoidal perturbationsusing the modes with the 16 longest wavelengths and random phases. The simulation isperformed with × grid points, and we have used hyper-dissipation rather thanregular (Laplacian) resistivity and viscosity.From the visualizations of current density shown in Fig. 17, we see the emergence oflarger structures from the initial small-scale random magnetic field. The behaviors of themagnetic energy decay and of the perpendicular magnetic spectrum are similar to theruns with the flux-tube setup, and are consistent with the predictions of our hierarchicalmodel. The time evolution of magnetic energy (Fig. 18) shows the t − scaling. In Fig. 19we show the time evolution of the perpendicular magnetic power spectrum M ( k ⊥ , t ) . Asthe system evolves, a k − / ⊥ inertial range is developed, and the spectra are observed tomigrate to larger scales (Fig. 19, top panel), demonstrating the inverse energy transfer.When normalized as described in Section 3.5, the spectra roughly overlap (Fig. 19 bottompanel), indicating the self-similar evolution of the system, as predicted by our hierarchicalmodel in Section 2.2.The similarities between this simulation and those reported in the main section of thispaper strongly support the notion that non-helical inverse transfer enabled by magneticreconnection is a generic phenomenon. REFERENCESAlexakis, A., Mininni, P. D. & Pouquet, A. Physical Review E (4), 046301. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer J z / J z , rms - - /2 0 /2-- /20/2 t = 0 A - - /2 0 /2 t = 200 A - - /2 0 /2 t = 402 A - - /2 0 /2 t = 599 A Figure 17: Current density J z (normalized by its root mean square value) at various timeson the xy plane at z = 0 for a run starting with a random magnetic field. t / A E M E k index -1 Figure 18: Time evolution of total magnetic ( E M ) and kinetic ( E K ) energies for the runstarting with random magnetic field. k M ( k , t ) k k / k , max ( t )10 M ( k , t ) / M m a x ( t ) t = 0 A t = 100 A t = 200 A t = 300 A t = 402 A t = 500 A t = 598 A Figure 19: Raw (top) and normalized (bottom) perpendicular magnetic power spectra, fora simulation starting with random magnetic fields. A k − / ⊥ slope is shown for reference.0 Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky Alves, E. P., Zrake, J. & Fiuza, F. Physical Review Letters (24), 245101. Bhat, P., Zhou, M. & Loureiro, N. F. In preparation . Bhattacharjee, A., Huang, Y., Yang, H. & Rogers, B. Physics of Plasmas (11),112102. Biskamp, D. The Physics of Fluids , 1520–1531. Boldyrev, S. (11),115002. Brandenburg, A., Kahniashvili, T. & Tevzadze, A. G. Physical Review Letters (7), 075001. Burgers, J. M. Advancesin applied mechanics , , vol. 1, pp. 171–199. Elsevier. Cho, J. & Vishniac, E. T. The Astrophysical Journal (1), 273–282. Christensson, M., Hindmarsh, M. & Brandenburg, A. Physical Review E (5), 056405. Dmitruk, P. & Gómez, D. O. TheAstrophysical Journal Letters , L63–L66. Drake, J. F., Opher, M., Swisdak, M. & Chamoun, J. N. The Astrophysical Journal , 963–974. East, W. E., Zrake, J., Yuan, Y. & Blandford, R. D. Physical Review Letters (9), 095002. Einaudi, G. & Velli, M. Physics of Plasmas , 4146–4153. Fermo, R. L., Drake, J. F. & Swisdak, M. Physics of Plasmas (1), 010702. Finn, J. M. & Kaw, P. K. The Physics ofFluids (1), 72–78. Fried, B. D. The Physics of Fluids (3), 337–337. Furno, I., Intrator, T. P., Hemsing, E. W., Hsu, S. C., Abbate, S., Ricci, P.& Lapenta, G. Physics of Plasmas (5), 055702. Gekelman, W., De Haas, T., Daughton, W., Van Compernolle, B., Intrator, T. &Vincena, S. Physical Review Letters (23), 235101. Gekelman, W, Lawrence, E & Van Compernolle, B The Astrophysical Journal (2), 131. Goldreich, P. & Lynden-Bell, D. Monthly Notices of the Royal Astronomical Society , 125. Goldreich, P. & Sridhar, S. The Astrophysical Journal , 763–775. Gruzinov, A. The Astrophysical Journal Letters (1), L15. Hazeltine, R. D. & Meiss, J. D. Plasma confinement . Courier Corporation. Hu, Q., Chen, Y. & le Roux, J. Journal of Physics: Conference Series , 012005. Huang, Y. & Bhattacharjee, A. Physics ofPlasmas (6), 062104. Huntington, C. M., Fiuza, F., Ross, J. S., Zylstra, A. B., Drake, R. P., Froula,D. H., Gregori, G., Kugland, N. L., Kuranz, C. C., Levy, M. C. & others ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer Observation of magnetic field generation via the Weibel instability in interpenetratingplasma flows. Nature Physics (2), 173. Intrator, T. P., Sun, X., Dorf, L., Sears, J. A., Feng, Y., Weber, T. E. & Swan,H. O. PlasmaPhysics and Controlled Fusion (12), 124005. Jara-Almonte, J., Ji, H., Yamada, M., Yoo, J. & Fox, W. Physical ReviewLetters (9), 095001. Kadomtsev, B. B. & Pogutse, O. P. Soviet physics, Journal of Experimental and Theoretical Physics , 575–590. Kato, T. N. Physics of Plasmas (8), 080705. Kato, T. N. & Takabe, H. (2), L93. Katz, B., Keshet, U. & Waxman, E. The AstrophysicalJournal (1), 375. Khabarova, O., Zank, G. P., Li, G., le Roux, J. A., Webb, G. M., Dosch, A. &Malandraki, O. E. The Astrophysical Journal (2), 181. Klimchuk, J. A., Patsourakos, S. & Cargill, P. J. The Astrophysical Journal (2), 1351–1362. Kulsrud, R. M. & Zweibel, E. G. Reports onProgress in Physics (4), 046901. Lapenta, G. PhysicalReview Letters (23), 235001. Lazarian, A. & Opher, M. The Astrophysical Journal , 8–21. Le Roux, J. A., Zank, G. P., Webb, G. M. & Khabarova, O. The Astrophysical Journal (2), 112. Loureiro, N. F., Dorland, W., Fazendeiro, L., Kanekar, A., Mallet, A., Vilelas,M. S. & Zocco, A. Computer Physics Communications , 45–63. Loureiro, N. F., Samtaney, R., Schekochihin, A. A. & Uzdensky, D. A. Physicsof Plasmas (4), 042303. Loureiro, N. F., Schekochihin, A. A. & Cowley, S. C. Physics of Plasmas (10), 100703–100703. Loureiro, N. F., Schekochihin, A. A. & Uzdensky, D. A. Physical Review E (1), 013102. Loureiro, N. F. & Uzdensky, D. A. Plasma Physics and Controlled Fusion (1),014021. Lyutikov, M., Sironi, L., Komissarov, S. S. & Porth, O. a Particle acceleration inrelativistic magnetic flux-merging events. Journal of Plasma Physics (6). Lyutikov, M., Sironi, L., Komissarov, S. S. & Porth, O. b Particle acceleration inrelativistic magnetic flux-merging events. Journal of Plasma Physics (6), 635830602. Mallet, A., Schekochihin, A. A., Chandran, B. D. G., Chen, C. H. K., Horbury,T. S., Wicks, R. T. & Greenan, C. C. Monthly Notices of the RoyalAstronomical Society (2), 2130–2139. Matthaeus, W. H. & Lamkin, S. L. Physics of Fluids , 2513–2534. McComas, D. J., Allegrini, F., Bochsler, P., Bzowski, M., Christian, E. R., Crew,G. B., DeMajistre, R., Fahr, H., Fichtner, H., Frisch, P. C., Funsten, H. O., Muni Zhou, Nuno F. Loureiro and Dmitri A. Uzdensky Fuselier, S. A., Gloeckler, G., Gruntman, M., Heerikhuisen, J., Izmodenov,V., Janzen, P., Knappenberger, P., Krimigis, S., Kucharek, H., Lee, M.,Livadiotis, G., Livi, S., MacDowall, R. J., Mitchell, D., Möbius, E., Moore,T., Pogorelov, N. V., Reisenfeld, D., Roelof, E., Saul, L., Schwadron, N. A.,Valek, P. W., Vanderspek, R., Wurz, P. & Zank, G. P. Science , 959. Medvedev, M. V., Fiore, M., Fonseca, R. A., Silva, L. O. & Mori, W. B. The AstrophysicalJournal Letters (2), L75. Medvedev, M. V. & Loeb, A. The Astrophysical Journal (2), 697. Olesen, P. Physics Letters B (3-4),321–325. O’Neill, S. M., Beckwith, K. & Begelman, M. C. MonthlyNotices of the Royal Astronomical Society (2), 1436–1452. Opher, M., Drake, J. F., Swisdak, M., Schoeffler, K. M., Richardson, J. D.,Decker, R. B. & Toth, G. The Astrophysical Journal , 71. Parker, E. N. The Astrophysical Journal , 293. Parker, E. N. Journal of Geophysical Research , 509–520. Parker, E. N. TheAstrophysical Journal , 635–647. Polomarov, O., Kaganovich, I. & Shvets, G. PhysicsReview Letter , 175001. Pouquet, A. Journal of FluidMechanics (1), 1–16. Reppin, J. & Banerjee, R. Physical Review E (5), 053105. Ruyer, C. & Fiuza, F. Physics Review Letter , 245002. Samtaney, R., Loureiro, N. F., Uzdensky, D. A., Schekochihin, A. A. & Cowley,S. C. Physical ReviewLetters (10), 105004. Schekochihin, A. A., Cowley, S. C., Dorland, W., Hammett, G. W., Howes, G. G.,Quataert, E. & Tatsuno, T. The Astrophysical JournalSupplement Series (1), 310. Sears, J., Feng, Y., Intrator, T. P., Weber, T. E. & Swan, H. O. Plasma Physics and Controlled Fusion (9), 095022. Servidio, S., Matthaeus, W. H., Shay, M. A., Cassak, P. A. & Dmitruk, P. PhysicalReview Letters (11), 115003. Servidio, S., Matthaeus, W. H., Shay, M. A., Dmitruk, P., Cassak, P. A. & Wan,M. Physics of Plasmas (3), 032315. Silva, L. O., Fonseca, R. A., Tonge, J. W., Dawson, J. M., Mori, W. B. & Medvedev,M. V. The Astrophysical Journal Letters (1), L121. Spitkovsky, A. The Astrophysical Journal Letters (1), L39. Stone, E. C., Cummings, A. C., McDonald, F. B., Heikkila, B. C., Lal, N. & Webber,W. R. Science , 2017–2020. ulti-scale dynamics of magnetic flux tubes and inverse magnetic energy transfer Stone, E. C., Cummings, A. C., McDonald, F. B., Heikkila, B. C., Lal, N. & Webber,W. R. Nature , 71–74. Strauss, H. R. The Physics of Fluids (1), 134–140. Sweet, P. A. Electromagnetic Phenomenain Cosmical Physics (ed. B. Lehnert), IAU Symposium , vol. 6, p. 123. Told, D., Jenko, F., TenBarge, J. M., Howes, G. G. & Hammett, G.W. Physicalreview letters (2), 025003. Towns, J., Cockerill, T., Dahan, M., Foster, I., Gaither, K., Grimshaw, A.,Hazlewood, V., Lathrop, S., Lifka, D., Peterson, G. D. & others Computing in Science & Engineering (5), 62–74. Uzdensky, D. A. & Goodman, J. The Astrophysical Journal , 608–629. Uzdensky, D. A., Loureiro, N. F. & Schekochihin, A. A. Physical Review Letters (23), 235002. Weibel, E. S. Physical Review Letters (3), 83. Yamada, M., Ji, H., Hsu, S., Carter, T., Kulsrud, R., Bretz, N., Jobes, F., Ono,Y. & Perkins, F. Physics of Plasmas (5), 1936–1944. Yuan, Y., Nalewajko, K., Zrake, J., East, W. E. & Blandford, R. D. The Astrophysical Journal (2), 92. Zank, G. P., Hunana, P., Mostafavi, P., Le Roux, J. A., Li, Gang, Webb,G. M., Khabarova, O., Cummings, A., Stone, E. & Decker, R. The Astrophysical Journal (2), 137. Zank, G. P. & Matthaeus, W. H. Journal of Plasma Physics , 85–100. Zhou, M., Bhat, P., Loureiro, N. F. & Uzdensky, D. A. PhysicalReview Research , 012004. Zrake, J. TheAstrophysical Journal Letters (2), L26. Zrake, J. & Arons, J. TheAstrophysical Journal847