Reynolds-number dependence of the dimensionless dissipation rate in homogeneous magnetohydrodynamic turbulence Moritz Linkmann,1,2,∗ Arjun Berera,2 and Erin E. Goldstraw3,2 1Department of Physics & INFN, University of Rome Tor Vergata, Via della Ricerca Scientifica 1, 00133 Rome, Italy 2SUPA, School of Physics and Astronomy, University of Edinburgh, Peter Guthrie Tait Road, EH9 3FD, UK 3School of Mathematics and Statistics, University of St. Andrews, KY16 9SS, UK (Dated: January 5, 2017) ThispaperexaminesthebehaviorofthedimensionlessdissipationrateC forstationaryandnon- ε stationary magnetohydrodynamic (MHD) turbulence in presence of external forces. By combining with previous studies for freely decaying MHD turbulence, we obtain here both the most general 7 model equation for C applicable to homogeneous MHD turbulence and a comprehensive numeri- ε 1 cal study of the Reynolds number dependence of the dimensionless total energy dissipation rate at 0 unitymagneticPrandtlnumber. Wecarryoutaseriesofmediumtohighresolutiondirectnumerical 2 simulations of mechanically forced stationary MHD turbulence in order to verify the predictions of n themodelequationforthestationarycase. Furthermore,questionsofnonuniversalityarediscussed a in terms of the effect of external forces as well as the level of cross- and magnetic helicity. The J measured values of the asymptote C lie between 0.193 (cid:54) C (cid:54) 0.268 for free decay, where ε,∞ ε,∞ 4 the value depends on the initial level of cross- and magnetic helicities. In the stationary case we measure Cε,∞ =0.223. ] n y d - I. INTRODUCTION u l f . The dynamics of conducting fluids is relevant to many areas in geo- and astrophysics as well as in engineering and s c industrial applications. Often the flow is turbulent, and the interaction of the turbulent flow with the magnetic field i leads to considerable complexity. Being a multi-parameter problem, techniques that have been successfully applied s y to turbulence in nonconducting fluids sometimes fail to deliver unambiguous predictions in magnetohydrodynamic h (MHD)turbulence. Thisconcernse.g.thepredictionofinertialrangescalingexponentsbyextensionofKolmogorov’s p arguments[1]toMHD,andconsiderableefforthasbeenputintothefurtherunderstandingofinertialrangecascade(s) [ in MHD turbulence [2–9]. The difficulties are partly due to the many different configurations that can arise in MHD 1 turbulence because of e.g. anisotropy, different levels of vector field correlations, different values of the dissipation v coefficients and different types of external forces, and as such are connected to the question of universality in MHD 6 turbulence[10–21]. Thebehaviorofthe(dimensionless)dissipationrateisrepresentativeofthisproblem,inthesense 8 thattheaforementionedpropertiesofMHDturbulenceinfluencetheenergytransferacrossthescales,i.e.thecascade 0 dynamics [11, 22–26], and thus the amount of energy that is eventually dissipated at the small scales. 1 The behavior of the total dissipation rate in a turbulent non-conducting fluid is a well-studied problem. As such it 0 has been known for a long time that the total dissipation rate in both stationary and freely decaying homogeneous . 1 isotropic turbulence tends to a constant value with increasing Reynolds number following a well-known characteristic 0 curve[27–32]. Forstatisticallysteadyisotropicturbulencethiscurvecanbeapproximatedbythereal-spacestationary 7 energy balance equation, where the asymptote is connected to the maximal inertial flux of kinetic energy [30]. The 1 : correspondingprobleminMHDhasreceivedmuchlessattention,however,recentnumericalresultsforfreelydecaying v MHD turbulence at unity magnetic Prandtl number report similar behavior. Mininni and Pouquet [33] carried out i X directnumericalsimulations(DNSs)offreelydecayinghomogeneousMHDturbulencewithoutameanmagneticfield, showing that the temporal maximum of the total dissipation rate ε(t) became independent of Reynolds number at a r a Taylor-scaleReynoldsnumberR (measuredatthepeakofε(t))ofabout200. DallasandAlexakis[34]measuredthe λ dimensionless dissipation rate C also from DNS data for free decay for random initial fields with strong correlations ε between the velocity field and the current density. Again, it was found that C → const. with increasing Reynolds ε number. Interestingly, a comparison with the data of Ref. [33] showed that the approach to the asymptote was slower than for the data of Ref. [33], suggesting an influence of the level of certain vector field correlations on the approach to the asymptote. A theoretical model for dissipation rate scaling in freely decaying MHD turbulence was put forward recently [35] based on the von K´arm´an-Howarth energy balance equations (vKHE) in terms of ∗ [email protected] 2 Els¨asser fields [36]. For unity magnetic Prandtl number it predicts the dependence of C on a generalized Reynolds ε number R ≡ z−L+/(ν +µ), with z− denoting the root-mean-square value of one Els¨asser field, L+ the integral − scale corresponding to the other Els¨asser field, while ν and µ are the kinematic viscosity and the magnetic resistivity, respectively. The model equation has the following form C D C =C + + +O(R−3) , (1) ε ε,∞ R R2 − − − where C and D are time-dependent coefficients depending on several parameters, which themselves depend on the magnetic,cross-andkinetichelicities. Thepredictionsofthisequationweresubsequentlytestedagainstdataobtained from medium to high resolution DNSs of freely decaying homogeneous MHD turbulence leading to a very good agreement between theory and data. In summary, there is compelling numerical and theoretical evidence for finite dissipation in freely decaying MHD turbulenceatleastforunitymagneticPrandtlnumberPm=ν/µ,whilesofarnosystematicresultsforthestationary casehavebeenreported. InthispaperweextendthederivationcarriedoutinRef.[35]toincludetheeffectsofexternal forces and we present the first systematic study of dissipation rate scaling for stationary MHD turbulence. In order to be able to test the model equation against DNS data for a large range of generalized Reynolds numbers, we concentrate on the case Pm=1. The most general form of Eq. (1) for nonstationary flows with large-scale external forcing is derived, which can be applied to freely decaying and stationary flows by setting the corresponding terms to zero. This generalization of Eq. (1) is the first main result of the paper, it is applicable to both freely decaying and stationary MHD turbulence. It implies that the dissipation rate of total energy is finite in the limit R →∞ in − analogytohydrodynamics,andhighlightsthedependenceofthecoefficientsC andD ontheexternalforces. Assuch, Eq. (1) predicts nonuniversal values of the asymptotic value C of the dimensionless dissipation rate in the infinite ε,∞ Reynolds number limit and of the approach to the asymptote for a variety of MHD flows. The resulting theoretical predictionsforthestationarycasearecomparedtoDNSdataforstationaryMHDturbulenceforthreedifferenttypes of mechanical forcing while the results for the freely decaying case [35] are reviewed for completeness. The DNS data shows good agreement with Eq. (1) and the different forcing schemes have no measurable effect on the values of the coefficients in Eq. (1). The measured values of C lie between 0.193(cid:54)C (cid:54)0.268 for free decay, where the value ε,∞ ε,∞ depends on the initial level of cross- and magnetic helicities. In the stationary case we measure C =0.223. ε,∞ This paper is organized as follows. We begin by reviewing the formulation of the MHD equations in terms of Els¨asser fields in Sec. II where we introduce the basic quantities we aim to study in both formulations of the MHD equations. In Section III we extend the derivation put forward in Ref. [35] to nonstationary MHD turbulence. The model equation is verified against DNS data for statistically steady MHD turbulence and the comparison to data for freely decaying MHD turbulence presented in Ref. [35] is reviewed in Sec. IV, where special emphasis is given to the question of nonuniversality of MHD turbulence in the context of external forces and the level of cross- and magnetic helicities. Our results are summarized and discussed in the context of related work in hydrodynamic and MHD turbulence in Sec. V, where we also outline suggestions for further work. II. THE TOTAL DISSIPATION IN TERMS OF ELSA¨SSER FIELDS InthispaperweconsiderstatisticallyhomogeneousMHDturbulenceintheabsenceofabackgroundmagneticfield. The flow is taken to be incompressible, leading to the following set of coupled partial differential equations 1 1 ∂ u=− ∇P −(u·∇)u+ (∇×b)×b+ν∆u+f , (2) t ρ ρ u ∂ b=(b·∇)u−(u·∇)b+µ∆b+f , (3) t b ∇·u=0 and ∇·b=0 , (4) whereudenotesthevelocityfield,bthemagneticinductionexpressedinAlfv´enunits,ν thekinematicviscosity,µthe magnetic resistivity, P the thermodynamic pressure, f and f are external mechanical and electromagnetic forces, u b whichmaybepresent,andρdenotesthedensitywhichissettounityforconvenience. Equations(2)-(4)areconsidered onathree-dimensionaldomainΩ,whichduetohomogeneitycaneitherbethefullspaceR3 orasubdomain[0,L )3 box withperiodicboundaryconditions. TheMHDequations(2)-(4)canbeformulatedmoresymmetricallyusingEls¨asser variables z± =u±b [37] 1 ∂ z± =− ∇P˜−(z∓·∇)z±+(ν+µ)∆z±+(ν−µ)∆z∓+f± , (5) t ρ ∇·z± =0 , (6) 3 where f± = f ±f and the pressure P˜ consists of the sum of the thermodynamic pressure P and the magnetic u b pressureρ|b|2/2. WhichformulationoftheMHDequationsischosenoftendependsonthephysicalproblem,forsome problems the Els¨asser formalism is technically convenient, while the formulation using the primary fields u and b facilitates physical understanding. The ideal invariants total energy E(t), cross-helicity H (t) and magnetic helicity c H (t) are given in the respective formulations of the MHD equation by m 1(cid:90) 1(cid:90) E(t)= dk (cid:104)|uˆ(k,t)|2+|bˆ(k,t)|2(cid:105)= dk (cid:104)|zˆ+(k,t)|2+|zˆ−(k,t)|2(cid:105) , (7) 2 4 Ω Ω (cid:90) 1(cid:90) H (t)= dk (cid:104)uˆ(k,t)·bˆ(−k,t)(cid:105)= dk (cid:104)|zˆ+(k,t)|2−|zˆ−(k,t)|2(cid:105) , (8) c 4 Ω Ω (cid:90) 1(cid:90) (cid:28)(cid:20)ik (cid:21) (cid:29) H (t)= dk (cid:104)aˆ(k,t)·bˆ(−k,t)(cid:105)= dk ×(zˆ+(k,t)−zˆ−(k,t)) ·(zˆ+(−k,t)−zˆ−(−k,t)) , (9) m 4 k2 Ω Ω with bˆ, uˆ and zˆ± denoting the respective Fourier transforms of the magnetic, velocity and Els¨asser fields, while aˆ is the Fourier transform of the magnetic vector potential a. The angled brackets indicate an ensemble average. Equation (9) is gauge-independent as shown in Appendix A. We now motivate the use of the Els¨asser formulation for the study of the dimensionless dissipation coefficient in MHD. In hydrodynamics, the dimensionless dissipation coefficient C is defined in terms of the Taylor surrogate ε,u expression for the total dissipation rate, U3/L , where U denotes the root-mean-square (rms) value of the velocity u field and L the integral scale defined with respect to the velocity field, as u L C ≡ε u . (10) ε,u kinU3 However, in MHD there are several quantities that may be used to define an MHD analogue to the Taylor surrogate expression, such as the rms value B of the magnetic field, one of the different length scales defined with respect to either b or u, or the total energy. SincethetotaldissipationinMHDturbulenceshouldberelatedtothefluxoftotalenergythroughdifferentscales, one may think of defining a dimensionless dissipation coefficient for MHD in terms of the total energy. However, this wouldleadtoanondimensionalizationofthehydrodynamictransfertermu·(u·∇)uwithamagneticquantity. This can be seen by considering the analog of the von K´arm´an-Howarth energy balance equation in real space [38] stated here for the case of free decay 3 (cid:18)r4 (cid:19) −d E(t)=ε(t)=−∂ (Buu(r,t)+Bbb (r,t))+ ∂ Buuu(r,t)+r4Cbbu (r,t) t t LL LL 2r4 r 6 LLL LLL + 6Cbub(r,t)+ 1 ∂ (cid:0)r4∂ (νBuu(r,t)+µBbb (r,t))(cid:1) , (11) r r4 r r LL LL where Buu, Bbb and Buuu are the longitudinal structure functions, Cbbu the longitudinal correlation function and LL LL LLL LLL Cbub another correlation function. The longitudinal structure and correlation functions are given by Buu(r,t)=(cid:104)(δu (r,t))2(cid:105) , (12) LL L Bbb (r,t)=(cid:104)(δb (r,t))2(cid:105) , (13) LL L Buuu(r,t)=(cid:104)(δu (r,t))3(cid:105) , (14) LLL L Cbbu (r,t)=(cid:104)u (x,t)b (x,t)b (x+r,t)(cid:105) , (15) LLL L L L where r =|r| and v =v·r/r denotes the longitudinal component of a vector field v, that is its component parallel L to the displacement vector r, and r δv (r)=[v(x+r)−v(x)]· , (16) L r its longitudinal increment. The function Cbub is defined through the third-order correlation tensor (cid:16)r r (cid:17) Cbub(r,t)=(cid:104)(u (x)b (x)−b (x)u (x))b (x+r)(cid:105)=Cbub(r,t) jδ − iδ . (17) ij,k i j i j k r ik r jk As can be seen from their respective definitions, the functions Cbbu and Cbub scale with B2U while the function LLL Buuu scales with U3. If Eq. (11) were to be nondimensionalized with respect to the total energy then the purely LLL 4 hydrodynamic term Buuu would be scaled partially by a magnetic quantity. LLL This problem of inconsistent nondimensionalization can be avoided by working with Els¨asser fields, which requires an expression for the total dissipation rate ε(t) in terms of Els¨asser fields. The total rate of energy dissipation in MHD turbulence is given by the sum of Ohmic and viscous dissipation ε(t)=ε (t)+ε (t) , (18) mag kin where (cid:90) ε (t)=µ dk k2(cid:104)|bˆ(k,t)|2(cid:105) , (19) mag Ω (cid:90) ε (t)=ν dk k2(cid:104)|uˆ(k,t)|2(cid:105) . (20) kin Ω Similarly, the total dissipation rate can be decomposed into its respective contributions from the Els¨asser dissipation rates 1(cid:0) (cid:1) ε(t)= ε (t)+ε (t) , (21) 2 + − where the Els¨asser dissipation rates are defined as (cid:90) (cid:90) ε±(t)=ν dk k2(cid:104)|zˆ±(k,t)|2(cid:105)+ν dk k2(cid:104)zˆ±(k,t)·zˆ∓(−k,t)(cid:105) , (22) + − Ω Ω with ν =(ν±µ). The total dissipation rate relates to the sum of the Els¨asser dissipation rates ± ε+(t)+ε−(t)=ε(t)+ε (t)+ε(t)−ε (t)=2ε(t) , (23) Hc Hc where the cross-helicity dissipation rate ε is given by Hc ε (t)= 1(cid:0)ε+(t)−ε−(t)(cid:1) . (24) Hc 2 Sincethispaperisconcernedwithbothstationaryandnonstationaryflows,thetotalenergyinputrateιmustalsobe considered. Similartothedissipationrate,theinputratecanbesplitupintoeitherkineticandmagneticcontributions or the Els¨asser contributions ι±(t) ι(t)=ι (t)+ι (t) (25) mag kin ι(t)= 1(cid:0)ι+(t)+ι−(t)(cid:1) . (26) 2 The latter equation can be rewritten as ι+(t)=ι(t)+ 1(cid:0)ι+(t)−ι−(t)(cid:1)=ι(t)+ι (t) , (27) 2 Hc where ι denotes the input rate of the cross-helicity. Hc III. DERIVATION OF THE EQUATION Since the total dissipation rate can be expressed either in terms of the Els¨asser fields or the primary fields u and b, itshouldbepossibletodescribeitalsobythevKHEforz± [36]. Forthefreelydecayingcasenofurthercomplication arises as the rate of change of total energy, which figures on the left-hand side of the energy balance, equals the total dissipation rate. However, in the more general case the rate of change of the total energy is given by the difference of energy input and dissipation. That is, in the most general case the total energy dissipation rate is given by ε(t)=ι(t)−d E(t) . (28) t 5 Forthestationarycased E(t)=0andoneobtainsε(t)=ι(t). Forthefreelydecayingcaseι(t)=0andthechangein t totalenergyisduetodissipationonly,thatis−d E(t)=ε(t). IntermsofEls¨asservariablesε(t)canalsobeexpressed t as ε(t)=ι(t)−d E(t)=ι(t)−d E±(t)∓d H (t) , (29) t t t c whereE±(t)denotetheEls¨asserenergies. Sincewehaverelatedthetotaldissipationratetotherateofchangeofthe Els¨asser energies, we are now in a position to consider the energy balance equations for z±, which are stated here for the most general case of homogeneous forced nonstationary MHD flows without a mean magnetic field 3 ∂ (cid:18)3r4 (cid:19) −∂ E±(t)+I±(r,t)=− ∂ B±±(r,t)− r C±∓±(r,t) t 4 t LL r4 2 LL,L + 3 ∂ (cid:0)r4∂ (ν+µ)B± (r,t)(cid:1) 4r4 r r LL + 3 ∂ (cid:0)r4∂ (ν−µ)B∓ (r,t)(cid:1) , (30) 4r4 r r LL where I±(r,t) are (scale-dependent) energy input terms and C±∓∓(r,t)=(cid:104)z±(x,t)z∓(x,t)z±(x+r,t)(cid:105) , (31) LL,L L L L B±±(r,t)=(cid:104)(δz±(r,t))2(cid:105) , (32) LL L B±∓(r,t)=(cid:104)δz±(r,t)δz∓(r,t)(cid:105) , (33) LL L L are the third-order longitudinal correlation function and the second-order structure functions of the Els¨asser fields, respectively. As can be seen from the definition, the third-order correlation function scales with (z±)2z∓, where z± denotetherespectivermsvaluesoftheEls¨asserfields. ThispermitsaconsistentnondimensionalizationoftheEls¨asser vKHE using the appropriate quantities defined in terms of Els¨asser variables. As such the complication that arose if the energy balance was written in terms of b and u can be circumvented. This motivates the definition of the dimensionless Els¨asser dissipation rates as ε(t)L (t) C±(t)≡ ± , (34) ε z±(t)2z∓(t) where 3π (cid:90) L (t)= dk k−1(cid:104)|z±(k,t)|2(cid:105) , (35) ± 8E±(t) Ω aretheintegralscalesdefinedwithrespecttoz± [39]. ForbalancedMHDturbulence, i.e.H =0, oneshouldexpect c C+(t)=C−(t), since ε ε E±(t)=2E(t)±2H (t)=2E(t) . (36) c Therefore all quantities defined with respect to the rms fields z+ and z− should be the same in this case. Finally, the dimensionless dissipation rate C (t) is defined as ε ε(t)L (t) ε(t)L (t) C (t)=C+(t)+C−(t)≡ + + − . (37) ε ε ε z+(t)2z−(t) z−(t)2z+(t) Using the definition given in Eq. (34), the Els¨asser energy balance equations (30) can now be consistently nondi- mensionalized. For conciseness the explicit time and spatial dependences are from now on omitted, unless there is a particular point to make. 6 A. Dimensionless von K´arm´an-Howarth equations By introducing the nondimensional variables σ =r/L [12] and non-dimensionalising Eq. (30) as proposed in the ± ± definitions of C± given in Eq. (34) one obtains ε −(cid:0)d E±−I±(cid:1) L± =− 1 ∂ (cid:32)3σ±4CL±L∓,±L(cid:33)− Lz± ∂ 3BL±L± t z±2z∓ σ4 σ± 2z±2z∓ z±2z∓ t 4 ± µ+ν 3 (cid:18) B±±(cid:19) + σ4∂ LL L±z∓4σ±4 ± σ± z±2 ν−µ 3 (cid:18) B±∓(cid:19) + σ4∂ LL . (38) L z±4σ4 ± σ±z±z∓ ± ± Beforeproceedingfurther,thescale-dependentforcingtermontheleft-handsideofthisequationneedstobeanalyzed insomedetailinordertoclarifyitsrelationtotheenergyinputratesιandι±. TheEls¨asserenergyinputI± isgiven by 3 (cid:90) r I±(r)= dr(cid:48)r(cid:48)2(cid:104)z±(x+r(cid:48))·f±(x)(cid:105) . (39) r3 0 Since the energy input rate is given by ι± =(cid:104)z±(x)·f±(x)(cid:105), the correlation function can be expressed as (cid:104)z±(x+r)·f±(x)(cid:105)=ι±φ±(r/L ) , (40) f whereφ± aredimensionlessevenfunctionsofr/L satisfyingφ±(0)=1andL thecharacteristicscaleoftheforcing. f f At scales much smaller than the forcing scale, i.e. for r/L << 1, for suitable types of forces φ±(r/L ) can be f f expanded in a Taylor series [40], leading to the following expression for the energy input 3 (cid:90) r (cid:34) (cid:18) r (cid:19)2 ∂2φ± (cid:12) (cid:32)(cid:18) r (cid:19)4(cid:33)(cid:35) I±(r)= dr(cid:48)r(cid:48)2ι± 1+ (cid:12) +O . (41) r3 0 Lf 2∂(r/Lf)2(cid:12)r/Lf=0 Lf In the limit of infinite Reynolds number the inertial range extends through all wavenumbers, formally implying that L →∞, where Eq. (41) implies I±(r)→ι±. Therefore it should be possible to split the term I±(r) into a constant, f ι±, and a scale-dependent term J±(r), which encodes the additional scale dependence introduced by realistic, finite Reynolds number forcing. For consistency, this scale-dependent term must vanish in the formal limit Re→∞. This can be achieved by writing I±(r) in terms of the correlation of the force and Els¨asser field increments 3 (cid:90) r I±(r)=ι±− dr(cid:48)r(cid:48)2(cid:104)δz±·δf±(cid:105) . (42) 2r3 0 Therefore we define 3 (cid:90) r J±(r)=− dr(cid:48)r(cid:48)2(cid:104)δz±·δf±(cid:105) , (43) 2r3 0 where lim J±(r) = 0. Hence the energy input I±(r) can be expressed as the sum of the scale-independent Re→∞ energy input rate ι± and a scale-dependent term which vanishes in the formal limit Re→∞ I±(r)=ι±+J±(r) , (44) with lim J±(r) = 0. Substitution of Eq. (44) into the nondimensionalized energy balance Eq. (38) leads to Re→∞ the dimensionless version of the Els¨asser vKHE for homogeneous MHD turbulence in the most general case for nonstationary flows at any magnetic Prandtl number ∂ (cid:32)3σ4C±∓±(cid:33) L (cid:18) 3B±± (cid:19) C± =− σ± ± LL,L + ± ±d H −∂ LL −J±∓ι ε σ4 2z±2z∓ z±2z∓ t c t 4 Hc ± 1 3∂ (cid:18) B±±(cid:19) 1 3∂ (cid:18) B±∓(cid:19) + σ± σ4∂ LL + σ± σ4∂ LL , (45) R∓ 2σ±4 ± σ± z±2 R±(cid:48) 2σ±4 ± σ±z±z∓ 7 where R and R(cid:48) denote generalized large-scale Reynolds numbers given by ∓ ± R =z∓L /(ν+µ) and R(cid:48) =z±L /(ν−µ) . (46) ∓ ± ± ± In order to express Eq. (45) more concisely, the following dimensionless functions are defined C±∓± g±∓± = LL,L , (47) z±2z∓ B±± h±± = LL , (48) z±2 B±∓ h±∓ = LL , (49) z±z∓ L H±± = ± ∂ B±± , (50) z±2z∓ t LL L F± = ± J± , (51) z±2z∓ L G± = ± d H , (52) z±2z∓ t c L Q± = ± ι , (53) z±2z∓ Hc such that Eq. (45) can be written as ∂ (cid:18)3σ4 (cid:19) 3 C± =− σ± ±g±∓± ±G±− H±±−F±∓Q± ε σ4 2 4 ± + 3 ∂σ± (cid:0)σ4∂ h±±(cid:1)+ 3 ∂σ± (cid:0)σ4∂ h±∓(cid:1) . (54) R σ4 ± σ± R(cid:48) σ4 ± σ± ∓ ± ∓ ± This equation can be applied to the two simpler cases of freely decaying and stationary MHD turbulence by setting the corresponding terms to zero. For the case of free decay there are no external forces therefore F± = 0, while for the stationary case the terms G± and H± vanish. A further simplification concerns the case Pm=1, that is ν =µ, where the inverse of the generalized Reynolds numbers R(cid:48) vanish. In this case the evolution of C± depends only on ± ε R , and an approximate analysis using asymptotic series is possible. Most numerical results are concerned with this ∓ case due to computational constraints, hence it would be very difficult to test an approximate equation against DNS data if not only Re but also Pm needs to be varied. From now on the magnetic Prandtl number is therefore set to unity, keeping in mind that the analysis could be extended to Pm(cid:54)=1 provided the approximate equation derived in the following section is consistent with DNS data. B. Asymptotic analysis for the case Pm=1 Equation (54) suggests a dependence of C± on 1/R , however, the structure and correlation functions also have ε ∓ a dependence on Reynolds number, which describes their deviation from their respective inertial-range forms. The highest derivative in Eq. (54) is multiplied by the small parameter 1/R , which suggests that this equation may be ∓ viewed as singular perturbation problem amenable to asymptotic analysis [41]. The Els¨asser vKHE was rescaled by the rms values of the Els¨asser fields and the corresponding integral length scales, where the integral scales are by definition the large-scale quantities, the interpretation in hydrodynamics usually being that they represent the size of the largest eddies. As such, the nondimensionalization was carried out with respect to quantities describing the large scales, that is, with respect to ‘outer’ variables. Hence outer asymptotic expansions of the nondimensional structure and correlation functions are considered with respect to the inverse of the (large-scale) generalized Reynolds numbers 1/R . We point out that the case Pm(cid:54)=1 would require expansions in two parameters, where the cases Pm>1 and ∓ Pm<1 must be treated separately due to a sign change in R(cid:48) between the two cases. ± Theformalasymptoticseriesofagenericfunctionf [usedforconcisenessinplaceofthefunctionsontheright-hand side of Eq. (54)] up to second order in 1/R reads ∓ 1 1 f =f + f + f +O(R−3) . (55) 0 R 1 R2 2 ∓ ∓ ∓ 8 AftersubstitutionoftheexpansionsintoEq.(54),collectingtermsofthesameorderin1/R ,onearrivesatequations ∓ describing the behavior of C+ and C− ε ε C± D± C± =C± + + +O(R−3) , (56) ε ε,∞ R R2 ∓ ∓ ∓ up to second order in 1/R , using the coefficients C± , C± and D± defined as ∓ ε,∞ ∂ (cid:18)3σ4 (cid:19) 3 C± =− σ± ±g±∓± ±G±− H±±∓Q± , (57) ε,∞ σ4 2 0 4 0 ± 3∂ (cid:20) (cid:18) g±∓±(cid:19)(cid:21) 3 C± = σ± σ4 ∂ h±±− 1 ∓F±− H±± , (58) σ4 ± σ± 0 2 1 4 1 ± 3∂ (cid:20) (cid:18) g±∓±(cid:19)(cid:21) 3 D± = σ± σ4 ∂ h±±− 2 ∓F±− H±± , (59) σ4 ± σ± 1 2 2 4 2 ± in order to write Eq. (54) in a more concise way. The zero-order term in the expansion of the function F± vanishes, since F± correspondsto thescale-dependent part J± ofthe energyinput whichvanishesin thelimit R →∞, hence ∓ F± =0. According to the definition of C in Eq. (37), the asymptote C is given by 0 ε ε,∞ C =C+ +C− , (60) ε,∞ ε,∞ ε,∞ andusingthedefinitionofthegeneralizedReynoldsnumbers,whichimpliesR =(L /L )(z+/z−)R onecandefine + − + − L z+ C =C++ − C− , (61) L z− + (D is defined analogously), resulting in the following expression for the dimensionless dissipation rate C D C =C + + +O(R−3) . (62) ε ε,∞ R R2 − − − Since the time dependence of the various quantities in this problem has been suppressed for conciseness, it has to be emphasized that Eq. (62) is time dependent, including the Reynolds number R . Equation (62) in conjunction with − eqs. (57)-(59) is the most general asymptotic expression for the Reynolds number dependence of C developed so far. ε It is applicable for freely decaying, stationary and non-stationary MHD turbulence in the presence of external forces, and it may be applied to the corresponding problem in non-conducting fluids by setting b = 0. As such it extends previous results for freely decaying MHD turbulence [35], as well as for the stationary case in homogeneous isotropic turbulence of non-conducting fluids [30]. For nonstationary MHD turbulence at the peak of dissipation the term H±± in Eq. (57) vanishes for constant flux 0 of cross-helicity (that is, d2H = 0), since in the infinite Reynolds number limit the second-order structure function t c will have its inertial range form at all scales. By self-similarity the spatial and temporal dependences of e.g. B++ LL should be separable in the inertial range, that is B++(r,t)∼(ε+(t)r)α , (63) LL for some value α, and ∂ B++ ∼αε+(t)α−1 d ε+rα . (64) t LL t At the peak of dissipation d ε+| =d ε| −d2H =d ε| =0 , (65) t tpeak t tpeak t c t tpeak whichimpliesH++(t )=0. Equation(57)takenfornonstationaryflowsatthepeakofdissipationisthusidentical 0 peak to Eq. (57) for stationary flows, which suggests that at this point in time a nonstationary flow may behave similarly to a stationary flow. We will come back to this point in Sec. IV. Due to selective decay, that is the faster decay of the total energy compared to H and H [25], in most situations one could expect d H to be small compared to ε c m t c in the infinite Reynolds number limit. In this case G± (cid:39)0 and ∂ (cid:18)3σ4 (cid:19) C± (t )=− σ± ±g±∓± , (66) ε,∞ peak σ4 2 0 ± which recovers the inertial-range scaling results of Ref. [36] and reduces to Kolmogorov’s four-fifth law for b=0. 9 C. Relation of C to energy and cross-helicity fluxes ε,∞ In analogy to hydrodynamics, the asymptotes C± should describe the total energy flux, that is the contribution ε,∞ of the cross-helicity flux to the Els¨asser flux should be canceled by the respective terms G± and Q± in Eq. (57). However, since this is not immediately obvious from the derivation, further details are given here. For nonstationary turbulence at the peak of dissipation, Eq. (57) for the asymptotes C± reduces to ε,∞ ∂ (cid:18)3σ4 (cid:19) C± =− σ± ±g±∓± ±G±∓Q± . (67) ε,∞ σ4 2 0 ± The dimensional version of this equation is ∂ (cid:18)3r4 (cid:19) ε=− r C±∓± ±d H ∓ι , (68) r4 2 LL,L t c Hc where it is assumed that the function C±∓± has its inertial range form corresponding to g±∓±. The function C±∓± LL,L 0 LL,L can also be expressed through the Els¨asser increments [36] C±∓± = 1(cid:0)(cid:104)(δz±(r))2δz∓(r)(cid:105)−2(cid:104)z±(x)z±(x)z∓(x+r)(cid:105)(cid:1) , (69) LL,L 4 L L L L L which can be written in terms of the primary fields u and b as 12 C±∓± = (cid:104)(δu (r))3−6b (x)2u (x+r)(cid:105) LL,L 43 L L L 12 ∓ (cid:104)(δb (r))3−6u (x)2b (x+r)(cid:105) , (70) 43 L L L (seee.g.Ref.[36]). ThetwotermsonthefirstlineofEq.(70)arethefluxtermsintheevolutionequationofthetotal energy, while the two terms on last line correspond to the flux terms in the evolution equation of the cross-helicity [36]. Now Eq. (68) can be expressed in terms of the primary fields ∂ (cid:18)3r4 (cid:19) ε=− r C±∓± ±d H ∓ι r4 2 LL,L t c Hc ∂ (cid:18)r4 (cid:19) =− r (cid:104)(δu (r))3−6b (x)2u (x+r)(cid:105) r4 4 L L L ∂ (cid:18)r4 (cid:19) ± r (cid:104)(δb (r))3−6u (x)2b (x+r)(cid:105) ±d H ∓ι r4 4 L L L t c Hc =ε ±ε ±d H ∓ι =ε , (71) T Hc t c Hc T where ε is the flux of total energy and ε the cross-helicity flux, which must equal −d H +ι for nonstationary T Hc t c Hc MHD turbulence. Thus the contribution from the third-order correlator C±∓± resulting in ε is canceled by d H − LL,L Hc t c ι ,or,afternondimensionalization,thecross-helicityfluxε L /[(z±)2z∓]iscanceledbyG±−Q±. Thetwosimpler Hc Hc ± casesoffreelydecayingandstationaryMHDturbulencearerecoveredbysettingeitherQ± =0(freedecay)orG± =0 (stationary case). D. Nonuniversality Since C is a measure of the flux of total energy across different scales in the inertial range, differences for the ε,∞ value of this asymptote should be expected for systems with different initial values for the ideal invariants H and m H . The flux of total energy and thus the asymptote C is an averaged quantity. This implies that cancellations c ε,∞ between forward and inverse fluxes may take place leading on average to a positive value of the flux, that is, forward transfer from the large scales to the small scales. In case of H (cid:54)=0, the value of C should therefore be less than m ε,∞ for H = 0 due to a more pronounced inverse energy transfer in the helical case, the result of which is less average m forward transfer and thus a smaller value of the (average) flux of total energy. For H (cid:54)= 0 the asymptote C is c ε,∞ expected to be smaller than for H = 0, since alignment of u and b weakens the coupling of the two fields in the c 10 inductionequation,leadingtolesstransferofmagneticenergyacrossdifferentscalesandpresumablyalsolesstransfer of kinetic to magnetic energy. Furthermore,fromananalysisofhelicaltriadicinteractionsinidealMHDcarriedoutinRef.[42]itmaybeexpected that high values of cross-helicity have a different effect on the asymptote C , depending on the level of magnetic ε,∞ helicity. The analytical results suggested that the cross-helicity may have an asymmetric effect on the nonlinear transfers in the sense that the self-ordering inverse triadic transfers are less quenched by high levels of H compared c to the forward transfers. The triads contributing to inverse transfers were mainly those where magnetic field modes of like-signed helicity interact, and so for simulations with maximal initial magnetic helicity the dynamics will be dominated by these triads. If the inverse fluxes are less affected by the cross-helicity than the forward fluxes, then the expectationis that for a comparison ofthe valueof C between systems with (i) high H and H , (ii) high H ε,∞ m c m and H =0, (iii) H =0 and high H and finally (iv) H =0 and H =0, the value of C should diminish more c m c m c ε,∞ between cases (i) and (ii) compared to between cases (iii) and (iv). Such a comparison is carried out in Sec. IV using DNS data. As can be seen from Eqs. (57)-(59), the force does not explicitly enter in the asymptote C but does so in the ε,∞ coefficients C and D. Therefore a dependence of C and D, and hence of the approach to the asymptote, on the force may be expected, while C appears to be unaffected by the external force. However, different external forces will ε,∞ lead to different energy transfer scenarios, e.g. mainly dynamo and inertial transfer or mainly conversion of magnetic tokineticenergyduetoastrongLorentzforce,thereforetheasymptotewillbeimplicitlyinfluencedbythat. Inshort, nonuniversal values of C , C and D are expected depending on the level of the ideal invariants and the type of ε,∞ external force. We will address this point in further detail in Secs. IV and V. IV. COMPARISON TO DNS DATA Before comparing Eq. (62) with DNS data the numerical method is briefly outlined. Equations (2)-(4) are solved numerically in a three-dimensional periodic domain of length L =2π using a fully de-aliased pseudospectral MHD box code [43, 44]. Both the initial magnetic and velocity fields are random Gaussian with zero mean with energy spectra given by E (k)=Ak4exp(−k2/(2k )2) , (72) mag,kin 0 whereA(cid:62)0isarealnumberwhichcanbeadjustedaccordingtothedesiredamountofinitialenergy. Thewavenumber k whichlocatesthepeakoftheinitialspectrumistakentobek =5unlessotherwisestated. Nobackgroundmagnetic 0 0 field is imposed. Several series of simulations have been carried out for stationary and freely decaying MHD turbulence. In the case of free decay the dependence of the asymptote on the initial level of the ideal invariants is studied. For the stationary simulations all helicities are initially negligible while the influence of different forcing methods is assessed by applying three different external mechanical forces labeled f , f and f to maintain the simulations in stationary 1 2 3 state, resulting in three different series of stationary DNSs. The forces always act at wavenumbers k (cid:54) k = 2.5, f i.e. at the large scales. The first type of mechanical force f corresponds to the DNS series ND in Tbl. III and is 1 given by fˆ(k,t)=(ι /2E )uˆ(k,t) for 0<|k|<k ; 1 kin f f =0 otherwise, (73) where fˆ(k,t) is the Fourier transform of the forcing and E is the total energy contained in the forcing band. The 1 f second type of mechanical force f , which corresponds to the DNS series HF in Tbl. III is a random δ(t)-correlated 2 process. ItisbasedonadecompositionoftheFouriertransformoftheforceintohelicalmodesandhastheadvantage that the helicity of the force can be adjusted at each wavevector [45], which gives optimal control over the helicity injection. For all simulations using this type of forcing the relative helicity of the force was set to zero. The third type of mechanical force [46] corresponds to the DNS series SF in Tbl. III and is given by   sink z+sink y (cid:88) f f f3 =f0 sinkfx+sinkfz , (74) sink y+sink x kf f f where f is an adjustable constant. This type of force is nonhelical by construction. 0 All three forces have been used in several simulations of stationary homogeneous MHD turbulence. The scheme labeled f was shown by Sahoo et al. [47] to keep the helicities at negligible levels even though zero helicity injection 1