text:
# Asymptotic properties of turbulent magnetohydrodynamics

## I Introduction

The nonlinear equations of magnetohydrodynamics (MHD) admit turbulent solutions that can cover a wide range of spatial scales. One expects to see the turbulent regime of MHD as the most probable state of a magnetofluid present in physical systems ranging from the solar corona (Van Ballegooijen 1985; 1986) and solar wind (Burlaga 1991) to galactic magnetic fields (DeYoung 1980; 1992). MHD turbulence in the early universe may have even left a detectable present-day signature at the cosmological scale (Brandenburg et al. 1996).

To study questions such as the asymptotic scaling behavior in magnetic field correlations, the presence or absence of kinetic and magnetic energy cascades, the helicity spectrum and its corresponding cascades, as well as the effect of random perturbations, one needs to consider the evolution of the magnetofluid over a wide range of spatial and temporal scales. In general, turbulence is an essential feature in coming to terms with these problems, whose rigorous treatment calls for an analysis of fully nonlinear MHD.

For the regime that interests us, the large-distance and long time behavior of the fluid velocity and magnetic fields, the equations of MHD lend themselves to an analytical treatment by applying the methods provided by the dynamic renormalization group (RG). The aim of the RG is to describe quantitatively how the dynamics of a system evolves as one changes the length and time scales of the phenomena under study. The RG is an excellent tool for exploring the physics of complicated physical systems that are characterized as having many interacting degrees of freedom and are subject to fluctuations covering many different length scales. RG methods were first used in the study of hydrodynamic (non-magnetic) turbulence by Forster, Nelson and Stephen (1977). They studied the large-distance and long-time properties of velocity correlations generated by the Navier Stokes equations for a randomly stirred incompressible fluid with a variety of Gaussian noise spectra. More recently these methods were used by Yakhot and Orszag (1986) and Dannevik, Yakhot and Orszag (1987) in an attempt to predict several parameters in fully developed hydrodynamic turbulence. A detailed critique of their use of the RG can be found in the paper by Eyink (1994). Some recent reviews of the RG in hydrodynamic turbulence may be found in Frisch (1995) and McComb (1995).

There have been rather fewer applications of the RG to magnetized fluids, probably due to the complexity of the problem and the lack of experimental data for turbulent magnetic fields in laboratory plasmas. Indeed, as indicated above, MHD turbulence is much more likely to be of interest at astrophysical and cosmological scales. The previous works, that apply RG to MHD turbulence for addressing various specific problems, have been carried out by Fournier, Sulem and Pouquet (1982), Longcope and Sudan (1991) and more recently by Camargo and Tasso (1992). A recent survey of MHD turbulence can be found in the monograph by Biskamp (1993).

In this paper, RG will be applied to examine the asymptotic properties of turbulent MHD with a stochastic force term. Our analysis treats the full MHD equations, and thus differs from Longcope and Sudan (1991), who treated a type of ”reduced MHD”. Also, we treat all nonlinear terms as given in the MHD equations on an equal footing, in contrast to Fournier et. al. (1982), who weighted the inertial nonlinearity and Lorentz force differently. In this respect, the starting equations for us are the same as those of Camargo and Tasso (1992), although our treatment differs significantly from their work. In regards to the calculation technique, we will cast the stochastic MHD equations into a functional path integral formalism, in contrast to the direct iterative solution of the equations of motion carried out in (Camargo and Tasso 1992). The functional path integral formalism permits use of useful field theory methods, which are much more compact and efficient than direct iterative expansions. Although in this paper we use this formalism only to implement our one-loop RG analysis, the formalism can be applied much more generally to study stochastic MHD. With respect to results, we differ from all the previous works (Fournier et al. 1982; Longcope and Sudan 1991; Camargo and Tasso 1992). We will examine cascade directions and energy mixing, present solutions of the correlation functions in the scaling regime and implement a RG improved perturbation theory. All of these results are being explored and presented here for the first time. Furthermore, we also will implement a RG analysis when magnetic helicity is present, and to our knowledge no such purely analytic treatment of this problem has been carried out before in the literature.

The paper is organized as follows. In Section II we introduce the stochastic equations of incompressible MHD. These equations then are written in terms of Elsasser variables, which make manifest the symmetry between the velocity and magnetic field for incompressible MHD (Elsasser 1950). This symmetry is useful in checking the correctness of the renormalization group equations that will be derived and for understanding the nature of the correlation functions and energy spectra for the velocity and magnetic fields. In the stochastic equations, the stochastic random forces and currents are intended to model both random initial conditions and subsequent disturbances and fluctuations. The noise therefore serves a two-fold purpose and it summarizes at a coarse-grained level our ignorance about the details of the numerous small-scale phenomena. For example, as applied to the problem of galactic or even cosmic fields, the noise would account for the many diverse particle physics mechanisms (e.g., early phase transitions, vacuum bubble collisions, fluctuating Higgs gradients, superconducting cosmic strings, etc.) that have been proposed for the generation of seed magnetic fields (for a review see eg. Enqvist (1998)). In this particular context then, the noise simulates in part the effects of the particle-physics degrees of freedom at the micro-scale. Although we will not examine models to motivate the noise forces, we emphasize that the noise need not be a stirring force (in the strict sense of turbulence modelling) but also can account, in part, for real physical fluctuations and perturbations intrinsic to the system.

In Section III we turn to the renormalization group analysis of incompressible MHD. First the physical motivations are discussed for implementing the RG techniques to MHD. Then Subsec. III A examines the naive scaling properties of stochastic MHD under independent space and time rescalings. Here, two fundamental scaling exponents are introduced, that encode all the information pertaining to the asymptotic behavior of the solutions of MHD. Subsec. III B describes the two step RG procedure of coarse graining and rescaling. It then explains the relation between scale invariant solutions of the MHD equations and fixed points of the RG analysis. In particular, at a fixed point the MHD equations become scale invariant, and a unique solution exists for the pair of scaling exponents, which are calculable by the RG procedure. To fix ideas, a simple example of the RG procedure is given in the context of linearized MHD. Of course this example does not account for the nonlinear interactions, so crucial for establishing and maintaining turbulence. In Subsec III C the RG analysis is done for the full MHD equations, in particular including the nonlinear terms. As alluded to above, most of the applications of RG to hydrodynamics and MHD require one to iterate the equations of motion to a particular order in the nonlinear coupling. This can result in rather cumbersome expressions after a few iterations. Here, we prefer to exploit the streamlined functional integral and allied diagrammatic methods which allow for the rapid identification of the one-particle irreducible (1PI) diagrams whose renormalization then yields directly the independent RG equations for the parameters appearing in the MHD equations. Since these techniques lie somewhat out of the mainstream of the paper and can be applied more generally to the problem, we relegate them to a series of technical Appendices, where the complete renormalization of the response function, noise spectral function and interaction vertex is carried out in full detail. In Subsec III C mainly the differential RG flow equations are presented, and the calculations leading to these equations are outlined with references to the appropriate places in the Appendices for the various details.

In Sec. IV, the RG flow equations are analyzed with a pair of dimensionless couplings. Then in this section, the fixed points are calculated, their stability properties are elucidated and the associated critical exponents are obtained for each fixed point. All these results are given for any spatial dimension $`d`$ and for a general class of stochastic noise forces, characterized by a single spectral force exponent $`y`$, that ranges from short ($`y<0`$) to long ($`y>0`$) range. In Section V we derive and solve differential equations for the scale dependent viscosity and resistivity, derive the scaling forms obeyed by the velocity and magnetic field correlation functions, and derive both the scaling form and RG improved form of the kinetic and magnetic energy spectra. Here, an expression is said to be renormalization group improved if the bare parameters in that expression are replaced by their scale-dependent forms, typically calculated to some order in perturbation theory. Thus, an improved quantity, whether it be a correlation function, an energy spectrum, etc., is one which combines perturbation theory with the renormalization group (Hochberg et al. 1999a). In Sec. VI, we turn to consider the injection of random magnetic helicity, which is modeled with a power law spectrum, and compute the RG-improved helicity spectrum. In this Section the correlation functions are recomputed, now taking into account the helicity contribution. Then once again, these correlation functions are computed in the RG improved form. An alternative way to understand scaling and the approach to scaling is given in Section VII, which makes use of the Callan-Symanzik equation whose solution yields the improved correlation function automatically. Sec. VIII is concerned with energy and helicity cascades in hydrodynamics and MHD. A methodology for extracting information about cascades from our RG calculations is developed and then implemented. The predictions from our RG analysis are given for cascades and cascade directions at the asymptotically largest length scales. Finally, in Section IX we draw conclusions and summarize our work. Here some discussion is provided, extending from Subsect V C, on the mixing of magnetic and kinetic energy at asymptotic scales.

As mentioned above, the Appendices contain most of the technical calculations of our RG analysis as well as a derivation of the functional path integral formalism that is generally applicable to stochastic MHD. In particular, the Appendices are organized as follows. In Appendix A we set up the dynamic functional for incompressible non-relativistic stochastic MHD and use this to derive the basic elements needed for a Feynman diagram representation of the solutions of the MHD equations and define the response function, noise spectral function and vertex function. In Appendix B we carry out the momentum-shell integration needed to renormalize the response function. From this we obtain the one-loop corrections to the kinematic and magnetic viscosities. The one-loop renormalization of the noise spectral function is carried out in Appendix C. For $`y>2`$, which corresponds to noise spectra from long range down to reasonably short range, we find that there is in fact no one-loop correction to the noise amplitude in the hydrodynamic limit. The one-loop renormalization of the vertex function vanishes identically for all gaussian forcing functions that are temporally white, but for arbitrary spatial correlations, as demonstrated in Appendix D. Some useful identities needed for averaging products of unit vectors over the unit $`d`$-sphere are collected in Appendix E.

## II Stochastic MHD Equations

The equations of motion for viscous and resistive, incompressible, non-relativistic, magnetohydrodynamics are given by (Biskamp 1993)

$$\frac{𝐯}{t}+(𝐯\stackrel{}{})𝐯=\frac{1}{\rho }\stackrel{}{}\left(p+\frac{B^2}{8\pi }\right)+\frac{1}{4\pi \rho }(𝐁\stackrel{}{})𝐁+\nu ^2𝐯+\stackrel{}{\eta }_v,$$
(1)
and

$$\frac{𝐁}{t}=(𝐁\stackrel{}{})𝐯(𝐯\stackrel{}{})𝐁+\frac{c}{4\pi \sigma }^2𝐁+\stackrel{}{\eta }_B,$$
(2)
subject to the constraints

$$\stackrel{}{}𝐯=0,$$
(3)
and

$$\stackrel{}{}𝐁=0,$$
(4)
where $`𝐯`$ denotes the fluid velocity, $`𝐁`$ the magnetic field, $`p`$ the fluid pressure and $`\rho `$ the density. The speed of light is $`c`$ and $`\sigma `$ is the electrical conductivity of the magnetofluid. We allow for the presence of random forces and random currents, denoted by $`\stackrel{}{\eta }_v(\stackrel{}{x},t)`$ and and $`\stackrel{}{\eta }_B(\stackrel{}{x},t)`$, respectively. We assume these forces are Gaussian with zero mean, $`\stackrel{}{\eta }_v=\stackrel{}{\eta }_B=0`$. The random vector forces and currents are specified by their two-point correlation functions, which for convenience we specify directly in Fourier space and in terms of their individual vector components as

$`(\stackrel{}{\eta }_v)_n(\stackrel{}{k},\omega )(\stackrel{}{\eta }_v)_m(\stackrel{}{k}^{},\omega ^{})`$ $`=`$ $`2D_vg(k)(2\pi )^{d+1}𝐏_{nm}(\stackrel{}{k})\delta (\omega +\omega ^{})\delta ^d(\stackrel{}{k}+\stackrel{}{k}^{}),`$ (5)
$`(\stackrel{}{\eta }_B)_n(\stackrel{}{k},\omega )(\stackrel{}{\eta }_B)_m(\stackrel{}{k}^{},\omega ^{})`$ $`=`$ $`2D_Bg(k)(2\pi )^{d+1}𝐏_{nm}(\stackrel{}{k})\delta (\omega +\omega ^{})\delta ^d(\stackrel{}{k}+\stackrel{}{k}^{}),`$ (6)
$`(\stackrel{}{\eta }_v)_n(\stackrel{}{k},\omega )(\stackrel{}{\eta }_B)_m(\stackrel{}{k}^{},\omega ^{})`$ $`=`$ $`0.`$ (7)
Here, $`d`$ is the number of space dimensions, the indices are $`n,m=1,2,\mathrm{},d`$, and the angular brackets denote the average taken over the fluctuations. The presence of both kinetic and magnetic helicities can be taken into account and modeled by means of suitably defined noise terms. We defer the treatment of random helicity to a later section. The transverse projection operator $`𝐏_{nm}(\stackrel{}{k})`$ is needed to ensure the noise is solenoidal, i.e., $`\stackrel{}{k}\stackrel{}{\eta }_v=\stackrel{}{k}\stackrel{}{\eta }_B=0`$, and hence compatable with incompressibility $`(\stackrel{}{}𝐯=0)`$ and the absence of magnetic monopoles $`(\stackrel{}{}𝐁=0)`$. The Fourier transform of the random forces is given by

$$\stackrel{}{\eta }_{v,B}(\stackrel{}{x},t)=\frac{d^dk}{(2\pi )^d}\frac{d\omega }{2\pi }\stackrel{}{\eta }_{v,B}(\stackrel{}{k},\omega )\mathrm{exp}(i\stackrel{}{k}\stackrel{}{x}i\omega t),$$
(8)
and this is the convention used throughout this paper. The magnitudes of the random forces and currents are given by $`D_v0`$ and $`D_B0`$, respectively, and so provide a measure of the size of the fluctuations. The power law function $`g(k)=k^y`$ with exponent $`y`$ is intended to describe the noise spectrum in the inertial range. Eq. (7) states that the random currents and forces are uncorrelated. The noise exponent $`y`$ will be treated as a model parameter. Note depending on the sign of $`y`$, the noise fluctuations are either more correlated at short distances $`(y<0)`$ or at large distances $`(y>0)`$.

We now cast the above coupled vector stochastic partial differential equations into a manifestly symmetric form first by defining $`𝐁𝐁/\sqrt{4\pi \rho }`$, $`pp/\rho `$ and $`p_{}=(p+\frac{1}{2}B^2)`$. Then, following Elsasser (1950), we add and subtract the two equations of motion Eqs. (1) and (2) to obtain

$`{\displaystyle \frac{\stackrel{}{P}}{t}}+(\stackrel{}{Q}\stackrel{}{})\stackrel{}{P}`$ $`=`$ $`\stackrel{}{}p_{}+\gamma _+^2\stackrel{}{P}+\gamma _{}^2\stackrel{}{Q}+\stackrel{}{\eta }_P,`$ (9)
$`{\displaystyle \frac{\stackrel{}{Q}}{t}}+(\stackrel{}{P}\stackrel{}{})\stackrel{}{Q}`$ $`=`$ $`\stackrel{}{}p_{}+\gamma _+^2\stackrel{}{Q}+\gamma _{}^2\stackrel{}{P}+\stackrel{}{\eta }_Q,`$ (10)
where

$`\stackrel{}{P}=\stackrel{}{v}+\stackrel{}{B},`$ $`\stackrel{}{Q}=\stackrel{}{v}\stackrel{}{B},`$ (11)
$`\gamma _+={\displaystyle \frac{1}{2}}(\nu +\nu _B),`$ $`\gamma _{}={\displaystyle \frac{1}{2}}(\nu \nu _B),`$ (12)
$`\stackrel{}{\eta }_P=\stackrel{}{\eta }_v+\stackrel{}{\eta }_B,`$ $`\stackrel{}{\eta }_Q=\stackrel{}{\eta }_v\stackrel{}{\eta }_B.`$ (13)
In terms of these Elsasser variables, the noise statistics is completely determined from Eqs. (5) - (7) together with the definitions Eq. (13) as

$`(\stackrel{}{\eta }_P)_n(\stackrel{}{k},\omega )(\stackrel{}{\eta }_P)_m(\stackrel{}{k}^{},\omega ^{})`$ $`=`$ $`2[D_v+D_B]g(k)(2\pi )^{d+1}𝐏_{nm}(\stackrel{}{k})\delta (\omega +\omega ^{})\delta ^d(\stackrel{}{k}+\stackrel{}{k}^{}),`$ (14)
$`=`$ $`(\stackrel{}{\eta }_Q)_n(\stackrel{}{k},\omega )(\stackrel{}{\eta }_Q)_m(\stackrel{}{k}^{},\omega ^{}),`$ (15)
$`(\stackrel{}{\eta }_P)_n(\stackrel{}{k},\omega )(\stackrel{}{\eta }_Q)_m(\stackrel{}{k}^{},\omega ^{})`$ $`=`$ $`2[D_vD_B]g(k)(2\pi )^{d+1}𝐏_{nm}(\stackrel{}{k})\delta (\omega +\omega ^{})\delta ^d(\stackrel{}{k}+\stackrel{}{k}^{}).`$ (16)
We henceforth write $`\nu _B=\frac{c}{4\pi \sigma }`$, to indicate the magnetic resistivity. This pair of dynamical equations Eqs. (9) and (10) is symmetric under the interchange $`\stackrel{}{P}\stackrel{}{Q}`$ and $`\stackrel{}{\eta }_P\eta _Q`$. Furthermore, in the limit of vanishing magnetic field, $`\stackrel{}{P}=\stackrel{}{Q}`$ and vanishing magnetic noise, $`\stackrel{}{\eta }_P=\stackrel{}{\eta }_Q`$, the pair reduces to two identical copies of the Navier Stokes equation with a random noise source. These symmetries and limits are useful in checking results derived from these equations and can be used as a diagnostic in comparing our calculations with other RG results obtained from MHD and from incompressible hydrodynamics.

For incompressible MHD, we can take the divergence of both equations Eqs. (9) and (10) to eliminate the pressure term $`p_{}`$ since $`\stackrel{}{}𝐯=\stackrel{}{}𝐁=0`$. This yields

$$p_{}=\frac{1}{^2}_n_k(Q_kP_n),$$
(17)
which can be regarded as the equation of state for this system. We note that in solving for $`p_{}`$, there are no boundary terms that contribute to Eq. (17). As in Eqs. (1) and (2), the stochastic random forces and currents can be assumed to be solenoidal, since any longitudinal components can be absorbed into the pressure $`p_{}`$. This simplifies the equations even further. To do this, we introduce the transverse projection operator, which is the same one appearing above in the noise correlation functions

$$𝐏_{jn}=\left(\delta _{jn}_j\frac{1}{^2}_n\right),$$
(18)
for $`j,n=1,2,\mathrm{},d`$. Of course, this projector is a nonlocal operator in coordinate space (as it involves the Green function $`^2`$), but its Fourier transform is local and easier to handle in actual computations. Transforming Eq. (18) using Eq. (8), gives

$$𝐏_{jn}(\stackrel{}{k})=\left(\delta _{jn}\frac{k_jk_n}{k^2}\right).$$
(19)
Now $`p_{}`$ can be eliminated from Eqs. (9) and (10) to give, in coordinate space,

$`{\displaystyle \frac{P_j}{t}}+\lambda _0𝐏_{jn}_l(Q_lP_n)\gamma _+^2P_j\gamma _{}^2Q_j=(\eta _P)_j,`$ (20)
$`{\displaystyle \frac{Q_j}{t}}+\lambda _0𝐏_{jn}_l(P_lQ_n)\gamma _+^2Q_j\gamma _{}^2P_j=(\eta _Q)_j.`$ (21)
This is the final form of the stochastic MHD equations that serves as the starting point for our renormalization group analysis. The presence of the projection operator makes it preferable to express these equations in component form, as we have done here. We introduce a parameter $`\lambda _0=1`$ in front of the nonlinear terms, which is useful for organizing perturbative expansions of the solutions of these equations. This is merely a bookkeeping device, since our true (and dimensionless) perturbation expansion parameter will be defined and calculated below. Finally, note that correlations (and other functions) in the fluid velocity and magnetic fields can be computed in terms of these Elsasser fields, since the two sets of variables are linearly related as

$`\stackrel{}{v}(\stackrel{}{x},t)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[\stackrel{}{P}(\stackrel{}{x},t)+\stackrel{}{Q}(\stackrel{}{x},t)\right],`$ (22)
$`\stackrel{}{B}(\stackrel{}{x},t)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[\stackrel{}{P}(\stackrel{}{x},t)\stackrel{}{Q}(\stackrel{}{x},t)\right].`$ (23)
## III Renormalization Group (RG) Analysis of MHD

The renormalization group (RG) is a systematic method for studying a system at different length scales. One of the key motivations for implementing this method, and the main application in this paper, is to verify a scaling hypothesis. In vague terms a scaling hypothesis is related to the statement that in the asymptotic regime, the correlations in a system all can be expressed in terms of a single suitable length scale. For the renormalization group method, vindication of a scaling hypothesis is related to showing that the underlying dynamical equations are self-similar at different scales, up to an overall rescaling of the space-time coordinates and the dynamical variables of the system. Sometimes self-similar solutions are referred to as fixed points, since the parameters of the system do not change under the RG transformation. For hydrodynamics and magnetohydrodynamics, the scaling hypothesis that has been considered (Forster et al. 1977; Yakhot and Orszag 1986; Dannevik et al. 1987; Fournier et al. 1982; Longcope and Sudan 1991; Camargo and Tasso 1992) and will be considered here, is associated with the asymptotic regime of large length and time separations in the correlation functions.

In this Section, the above statements will be explained in detail and then applied. Subsection III A first will explain the scale transformations of the MHD equations that are important for the renormalization group procedure. Subsection III B then explains the renormalization group in the context of MHD. Finally, Subsection III C presents the results of applying the renormalization group to MHD. This last subsection requires the technically most involved aspect of the calculation in this paper. For ease of presentation, all the technical steps have been given in the Appendices and Subsec. III C simply quotes the results and gives some guidance to the Appendices.

### A Scale Transformation

In preparation for a dynamical renormalization group analysis, and as an important and necessary preliminary step, the stochastic MHD equations Eqs. (20) and (21) are submitted to independent global rescalings of space and time as

$`\stackrel{}{x}`$ $`=`$ $`s\stackrel{}{x}^{},\stackrel{}{}={\displaystyle \frac{1}{s}}\stackrel{}{}^{},`$ (24)
$`t`$ $`=`$ $`s^zt^{}{\displaystyle \frac{}{t}}=s^z{\displaystyle \frac{}{t^{}}},`$ (25)
and the (Elsasser) fields as

$$\stackrel{}{P}(\stackrel{}{x},t)=s^\chi \stackrel{}{P}^{}(\stackrel{}{x}^{},t^{}),\stackrel{}{Q}(\stackrel{}{x},t)=s^\chi \stackrel{}{Q}^{}(\stackrel{}{x}^{},t^{}),$$
(26)
where $`s>1`$. To fix ideas, one can think of this rescaling as examining the original unprimed system in a new primed coordinate system, with the two related by the above rescalings.

The above transformation involves two exponents, the dynamic exponent $`z`$ and the “roughness” exponent $`\chi `$. The choice of the above scale transformation is made because in the renormalization group analysis, the number of independent space and time rescalings one must consider depends on the degree of symmetry and anisotropy of the system under study. As we are considering homogeneous and isotropic magnetohydrodynamic turbulence, this scaling introduces two a-priori independent exponents. Alternatively, for example in reduced (2+1)-dimensional MHD where one has a dominant background magnetic field, say along the $`z`$-direction, one scales the $`(x,y)`$-plane and the perpendicular $`z`$ directions independently (Longcope and Sudan 1991) and this introduces three independent scaling exponents: $`z,\chi ,\xi `$.

An important step in examining self-similar behavior in the renormalization group approach is to examine the original unprimed equations in the primed coordinate system and compare the form of the equations. The precise relevance of this step will be explained in the next subsection. However, at a procedural level, in performing this step, one term in both set of equations always can be made similar by dividing through by an appropriate factor. For us, we will choose this term to be the time derivative terms $`P/t`$, $`Q/t`$. Thus Eqs. (20) and (21) when expressed in the primed coordinates are

$`_t^{}P_j^{}(\stackrel{}{x}^{},t^{})+s^{z+\chi 1}\lambda _0𝐏_{}^{}{}_{jn}{}^{}_l^{}(Q_l^{}P_n^{})`$ $`=`$ $`s^{z2}\gamma _+^{}_{}{}^{}2P_j^{}(\stackrel{}{x}^{},t^{})s^{z2}\gamma _{}^{}_{}{}^{}2Q_j^{}(\stackrel{}{x}^{},t^{})`$ (27)
$`+`$ $`s^{z\chi +\frac{1}{2}(ydz)}\eta _{P_j}(\stackrel{}{x}^{},t^{}),`$ (28)
and similarly for the $`\stackrel{}{Q}`$-equation. There are three points about Eq. (27) to note. First, in arriving at Eq. (27) a factor of $`s^{z+\chi }`$ was divided through the equation so that the time derivative terms in this and Eqs. (20) and (21) are the same. Second, the transformation behavior of the noise functions is fully specified by the space and time rescalings in Eq. (24) and the form of the noise correlations in Eqs. (5) - (7). Third, the MHD equations generically have an implicit short distance cut-off $`\mathrm{\Delta }x>\mathrm{\Delta }x_0`$ or in momentum space $`|𝐤|<\mathrm{\Lambda }1/\mathrm{\Delta }x_0`$. For distances less than this cut-off the hydrodynamic approximation is invalid in properly representing the physics. The continuum description of the fluid breaks down. The scale transformations above also act on this short distance cut-off so that in the primed frame this cut-off becomes $`\mathrm{\Delta }x_0^{}=s\mathrm{\Delta }x_0`$ or, equivalently $`\mathrm{\Lambda }^{}=\mathrm{\Lambda }/s`$. Since $`s>1`$, we see that the cut-off (or lattice spacing) has been scaled up (in real space) towards the infrared.

We compare this scaled equation of motion Eq. (27) with its unscaled version Eq. (20). If the scaled equation Eq. (27) is to describe the same dynamics, albeit at the larger scale $`(s>1)`$, then this implies that the coefficients appearing in it must change accordingly. By “same dynamics” we mean of course that the mathematical form of the dynamical equation does not change as one changes the scale of the problem. This implies the following scaling behavior of the parameters appearing in the MHD equations:

$`\lambda _0^{}`$ $`=`$ $`s^{z+\chi 1}\lambda _0,`$ (29)
$`\gamma _+^{}`$ $`=`$ $`s^{z2}\gamma _+,`$ (30)
$`\gamma _{}^{}`$ $`=`$ $`s^{z2}\gamma _{},`$ (31)
$`A^{}`$ $`=`$ $`s^{z2\chi d+y}A,`$ (32)
$`B^{}`$ $`=`$ $`s^{z2\chi d+y}B.`$ (33)
We have defined the noise-amplitude combinations $`A=D_v+D_B`$ and $`B=D_vD_B`$ (Thus, the zero magnetic noise limit $`D_B=0`$ is had by taking $`A=B`$).

It is important to appreciate that Eq. (27) with the definitions Eq. (29) and Eqs. (20) and (21) are simply two expressions of the same equation. In particular no requirement of scale invariance is made. We simply have required that any changes they suffer due to this scaling transformation be absorbable into redefinitions of the (bare) parameters appearing in the equations.

There is one final notational comment in regards to the scale transformation. Due to the central role that this transformation plays in the renormalization group analysis, it is standard to assign all quantities in the dynamical equation a scaling dimension. Typically the momentum coordinates, which scale as $`s^1`$, is referred to as having scaling dimension 1, so that length has scaling dimension $`1`$. Similarly time, the fields, and all the parameters in Eq. (29) have scaling dimension given by minus the exponent to which $`s`$ is raised.

### B Renormalization Group Procedure

The RG transformation is now rather standard fare and is discussed in a number of excellent texts and monographs (Ma 1976; Amit 1978; Zinn-Justin 1996; Cardy 1996). Here we briefly review the RG procedure. The RG is a symmetry transformation that carries a system from one length scale to another through a two step procedure, coarse graining and rescaling. As with any symmetry transformation, such as rotation, translation etc., it is not necessary that the system under examination be invariant to this transformation. However, should this be the case, then is as also the case with any symmetry transformation, important information about the system can be deduced without detailed knowledge of its entire dynamics. In particular, invariance under the RG transformation will imply scale invariance in the system. The existence of such an invariance implies the behavior of correlation functions at different scales can be understood relative to their behavior at one given scale.

Let us review the two step RG procedure. Consider the exact MHD Eqs. (9) and (10) or equivalently Eqs. (20) and (21). As stated earlier, there is an implicit assumption that these equations are only valid down to some distance scale $`x>\mathrm{\Delta }x_0`$ or in momentum space $`k<\mathrm{\Lambda }1/\mathrm{\Delta }x_0`$. This is the smallest length scale at which the MHD dynamics, Eqs. (9) and (10), are meant to be valid. Given this starting point, the RG is implemented in the following two steps:

Coarse-graining: The first step is to solve the dynamics at the very shortest length scales. In Fourier space we integrate the solutions of the stochastic equations of motion over a thin momentum shell $`\mathrm{\Lambda }/s|\stackrel{}{k}|<\mathrm{\Lambda }`$ where $`s=1+\delta `$ and $`0<\delta <<1`$. Physically, this step serves to integrate out the fast or short wavelength modes and thins out the degrees of freedom. This elimination process leads to “running” in the parameters (which depend on $`s=e^{\mathrm{}}`$) appearing in the MHD equations. The short distance cut-off or lattice spacing gives rise to a corresponding momentum cut-off $`\mathrm{\Lambda }/s`$. The resulting MHD equations after this step will be given by

$`_tP_j^<`$ $`+`$ $`\lambda _0^<(\mathrm{})𝐏_{jn}_l(Q_l^<P_n^<)\gamma _+^<(\mathrm{})^2P_j^<`$ (34)
$``$ $`\gamma _{}^<(\mathrm{})^2Q_j^<(\eta _P^<)_j+\mathrm{less}\mathrm{relevant}=0,`$ (35)
and likewise for $`Q_j`$. Clearly in performing this step additional terms will appear in the dynamic equations. These terms have not been explicitly written above, but have been referred to suggestively as ”less relevant”, since we will see later that for the leading asymptotic regime they will not be important. Nevertheless, at this stage if we retain all the terms in the coarse grained equations Eq. (34), including the ”less relevant” terms, the dynamics from these equations is exactly the same as from the original MHD Eqs. (9) and (10), for computing correlations for length scales larger than $`\mathrm{\Delta }x>2\pi s/\mathrm{\Lambda }`$. In particular, correlations computed from Eqs. (9) and (10) or Eqs. (34) will be the same, i.e.

$$P_i(\stackrel{}{x}_1,t_1)P_j(\stackrel{}{x}_2,t_2)=P_i^<(\stackrel{}{x}_1,t_1)P_j^<(\stackrel{}{x}_2,t_2),\mathrm{etc}.$$
(36)
Rescaling: The second step is to represent the above coarse grained Eqs. (34) in the primed coordinates based on the rescalings in Eqs. (24) and (26) with $`s=e^{\mathrm{}}`$. In the primed coordinates, the MHD equations become

$`_t^{}P_j^{}`$ $`+`$ $`\lambda _0^{}(\mathrm{})𝐏_{jn}_l^{}(Q_l^{}P_n^{})\gamma _+^{}(\mathrm{})^2P_j^{}`$ (37)
$``$ $`\gamma _{}^{}(\mathrm{})^2Q_j^{}s^{(yzd)/2}(\eta _P)_j(x^{},t^{})+\mathrm{less}\mathrm{relevant}\mathrm{rescaled}=0,`$ (38)
and likewise for $`Q`$. In these equations

$`\lambda _0^{}(\mathrm{})`$ $`=`$ $`s^{z+\chi 1}\lambda _0^<,`$ (39)
$`\gamma _+^{}(\mathrm{})`$ $`=`$ $`s^{z2}\gamma _+^<,`$ (40)
$`\gamma _{}^{}(\mathrm{})`$ $`=`$ $`s^{z2}\gamma _{}^<,`$ (41)
$`A^{}(\mathrm{})`$ $`=`$ $`s^{z2\chi d+y}A^<,`$ (42)
$`B^{}(\mathrm{})`$ $`=`$ $`s^{z2\chi d+y}B^<.`$ (43)
Note, if we solve Eqs. (37) and examine correlations at a distance $`d=(x_1^{}x_2^{})`$, this is equivalent to solving for the same correlation in the original Eqs. (9) and (10), at distance $`e^{\mathrm{}}d=(x_1x_2)`$.

At this point we can consider Eqs. (37) in place of the original MHD Eqs. (9) and (10), and repeat these two steps iteratively. A very interesting thing can happen upon sufficient iteration, which is a key practical importance of the RG approach. The set of parameters in iterative versions of Eqs. (37) could be the same after a certain number of successive RG transformations. These points in the parameter space are called fixed points. At the fixed point, correlations computed from either of two successive equations clearly are identically the same, since the successive equations are exactly the same. This fact has nontrivial consequences, since due to the scaling relations between the two sets of equations, it also means we can learn about the behavior of correlations at different scales of the original MHD equations. In particular consider the field $`P(\stackrel{}{x},t)`$, $`P^{}(\stackrel{}{x}^{},t^{})`$ (likewise $`Q`$) from two successive iterations of Eqs. (37). If these two successive versions of Eqs. (37) are exactly the same then it is an identity that (dropping the indices)

$$P(\stackrel{}{x}_1,t_1)P(\stackrel{}{x}_2,t_2)=P(\stackrel{}{x}_1^{},t_1^{})P(\stackrel{}{x}_2^{},t_2^{}),$$
(44)
and likewise for $`Q`$. However, based on the rescaling Eqs. (24) and (26) we also have

$$P(x_1,t_1)P(x_2,t_2)=e^{2\chi \mathrm{}}P(e^{\mathrm{}}x_1,e^z\mathrm{}t_1)P(e^{\mathrm{}}x_2,e^z\mathrm{}t_2),$$
(45)
which therefore gives a relation between correlations at two different length/time scales.

To help focus these ideas, let us consider a simple example of a fixed point solution in RG for the case of MHD without the nonlinear terms. In this case, the MHD dynamics from Eqs. (9) and (10) are

$`_tP_j\gamma _+^2P_j\gamma _{}^2Q_j(\eta _P)_j=0,`$ (46)
$`_tQ_j\gamma _+^2Q_j\gamma _{}^2P_j(\eta _Q)_j=0.`$ (47)
Implementing the RG procedure, the first step is coarse graining. Since this is a free (and hence, linear) theory, this step is trivial. We simply remove the highest modes between $`\mathrm{\Lambda }/s<|k|<\mathrm{\Lambda }`$ from the $`P`$ and $`Q`$ fields. The second step then is rescaling which based on Eqs. (24) leads to

$`_t^{}P_j^{}(x^{},t^{})s^{z2}\gamma _+^2P_j^{}s^{z2}\gamma _{}^2Q_j^{}s^{z\chi +\frac{1}{2}(ydz)}\eta _{P_j}(x^{},t^{})=0,`$ (48)
and similarily for $`Q`$. Observe if

$`z`$ $`=`$ $`2`$ (49)
$`\chi `$ $`=`$ $`{\displaystyle \frac{yd+z}{2}},`$ (50)
then both Eqs. (46) and (48) are the same, and so this is a fixed point of the RG.

Using this fact to compute correlations at two different scales in the momentum space correlation functions yields

$`P(k_1,\omega _1)P(k_2,\omega _2)`$ $`=`$ $`e^{2(\chi +d+z)\mathrm{}}P^{}(e^{\mathrm{}}k_1,e^z\mathrm{}\omega _1)P^{}(e^{\mathrm{}}k_2,e^z\mathrm{}\omega _2)`$ (51)
$`=`$ $`e^{2(\chi +d+z)\mathrm{}}P(e^{\mathrm{}}k_1,e^z\mathrm{}\omega _1)P(e^{\mathrm{}}k_2,e^z\mathrm{}\omega _2),`$ (52)
and similarily for $`Q`$. Since this is a free theory, we can check this result of the RG procedure from the exact correlation functions. For this note that the above fixed point Eq. (49) is the same for both linearized MHD Eqs. (46) and (47) or linearized hydrodynamics, which arises from these equations by identifying $`P=Q`$. Also, the exact solutions of both linearized MHD and linearized hydrodynamics have the same basic form, but the latter are much less elaborate. Thus, it is sufficient for our demonstrative purpose here to quote the simpler case of linearized hydrodynamics, for which the exact two-point correlation function is

$$P(k_1,\omega _1)P(k_2,\omega _2)const\times \frac{k_1^y\delta ^d(\stackrel{}{k}_1+\stackrel{}{k}_2)\delta (\omega _1+\omega _2)}{\omega _1^2+(\nu k_1^2)^2},$$
(53)
where $`\nu =\gamma _++\gamma _{}`$ from Eq. (12). The scaling transformations Eqs. (24) and (26)) applied to this imply

$`e^{2(\chi +d+z)\mathrm{}}P(e^{\mathrm{}}k_1,e^z\mathrm{}\omega _1)P(e^{\mathrm{}}k_2,e^z\mathrm{}\omega _2)`$ $``$ $`e^{2(\chi +d+z4ydz)}const\times `$ (55)
$`{\displaystyle \frac{k_1^y\delta ^d(\stackrel{}{k}_1+\stackrel{}{k}_2)\delta (\omega _1+\omega _2)}{\omega _1^2+(\nu k_1^2)^2}}.`$
Comparing the scaling behavior from the RG analysis Eq. (51) and the exact solution, Eq (55), we see they agree.

In summary, returning to the opening statements of this section, finding fixed points in the RG procedure corresponds to identifying scale invariance in the MHD equations with respect to the scale transformations Eqs. (24) and (26). As one practical matter, since the iterative procedure in Eqs. (34) and (37) is focussed on the behavior of the parameters, it is standard and convenient to express the RG procedure through a set of differential equations for the evolution of the parameters, eg. $`d\lambda (\mathrm{})/d\mathrm{},d\gamma _+(\mathrm{})/d\mathrm{}`$, etc. The fixed points then correspond to the solutions where all parameters have zero evolution, that is, $`d\lambda (\mathrm{})/d\mathrm{}=d\gamma _+(\mathrm{})/d\mathrm{}=\mathrm{}=0`$.

### C Differential Renormalization Group Equations for MHD

This subsection presents the RG evolution equations for MHD. There are several steps to this calculation, all of which have been carried out in full detail in the Appendices. Here these steps only will be outlined with reference to the appropriate Appendices for the complete details, and then the results will be given. For solving the equations of stochastic MHD, in Appendix A they have been Fourier transformed and cast into a dynamic functional formalism. This formalism allows a systematic loop-expansion, which in other words is a perturbation expansion in the amplitude of the random fluctuations. The basic building blocks for constructing the perturbative expansion are expressed diagrammatically through the Feynman rules in Fig. 1. Based on this approach, the information contained in the explicit solutions of $`\stackrel{}{P},\stackrel{}{Q}`$ is equivalently contained, but much more compactly, in the response function, noise spectrum and vertex function, whose definitions and diagrammatic representations (expressed to one-loop order) are given in Figs. (2a)-(2c). Thus, we work directly with these objects as they will yield, to a given order in perturbation theory, the corrections to the MHD parameters $`\lambda _0,\gamma _+,\gamma _{},A`$ and $`B`$. As stated in the last subsection, these corrections result from moving to ever larger length scales while taking into account the fluctuations from the smaller scales. We have calculated the one loop corrections to the response function, the noise spectral function and the nonlinear vertex function in the Appendices B, C and D respectively, and from them, have obtained the corresponding one-loop corrections to the parameters $`\lambda _0,\gamma _+,\gamma _{},A`$ and $`B`$ in the hydrodynamic limit (i.e., the limit of large distance and long time). We point out that in arriving at the RG flow equations below, we have taken the limit $`\omega 0`$ at the outset of the calculations as described in the Appendices. This is equivalent to taking the limit $`t\mathrm{}`$, and so serves to “project” out the asymptotic time limit. Finally the second step of the RG procedure is to rescale the coarse grained quantities based on Eqs. (24) and (26)

Applying these steps to the one-loop expressions Eqs. (D17), (C16), (C17), (B22), and (B23), computed in the Appendices, we obtain the following set of one-loop renormalization group equations for nonrelativistic incompressible MHD:

$`{\displaystyle \frac{d\lambda _0}{d\mathrm{}}}`$ $`=`$ $`\left(z+\chi 1\right)\lambda _0,`$ (56)
$`{\displaystyle \frac{dA}{d\mathrm{}}}`$ $`=`$ $`\left(z2\chi d+y\right)A,`$ (57)
$`{\displaystyle \frac{dB}{d\mathrm{}}}`$ $`=`$ $`\left(z2\chi d+y\right)B,`$ (58)
$`{\displaystyle \frac{d(\gamma _++\gamma _{})}{d\mathrm{}}}`$ $`=`$ $`(\gamma _++\gamma _{})\left\{(z2)+{\displaystyle \frac{𝒮}{2}}\lambda _0^2A_d\mathrm{\Lambda }^{dy4}{\displaystyle \frac{(A+B)}{(\gamma _++\gamma _{})^3}}\right\},`$ (59)
$`{\displaystyle \frac{d(\gamma _+\gamma _{})}{d\mathrm{}}}`$ $`=`$ $`(\gamma _+\gamma _{})\left\{(z2)+{\displaystyle \frac{𝒮}{2}}\lambda _0^2A_d\mathrm{\Lambda }^{dy4}{\displaystyle \frac{(AB)}{(\gamma _+\gamma _{})^3}}\right\},`$ (60)
where

$$A_d=\frac{(d^2y4)S_d}{d(d+2)(2\pi )^d},$$
(61)
and $`S_d=\frac{2\pi ^{d/2}}{\mathrm{\Gamma }(d/2)}.`$ Here, $`𝒮=4`$ is the numerical symmetry factor associated to the diagram representing the one-loop correction to the response function (see, Fig. (2a)), which has been calculated in Appendix B. The short-distance or ultraviolet cut-off (or, lattice spacing) is $`\mathrm{\Lambda }`$, below which the hydrodynamic or continuum description of the fluid breaks down. We can set $`\mathrm{\Lambda }=1`$ without loss of generality. In fact, none of the RG results depend explicitly on this cut-off. In Section VII we provide an alternative way to understand scaling which makes implicit use of the short-distance divergences that result from formally taking the limit $`\mathrm{\Lambda }\mathrm{}`$ in the one-loop expressions. It is clear that Eq. (56) describe the scale dependence of the parameters of MHD. Their solution therefore yields the effective parameters $`\lambda _0(\mathrm{}),A(\mathrm{}),B(\mathrm{}),\gamma _+(\mathrm{})`$ and $`\gamma _{}(\mathrm{})`$. Substitution of these scale-dependent parameters into the MHD equations Eqs. (20) and (21) yields the effective dynamical equations corresponding to this scale as parametrized by $`\mathrm{}0`$.

In terms of the dimensionless couplings defined by

$$g_\pm =\frac{\lambda _0^2}{2}A_d\mathrm{\Lambda }^{dy4}\frac{(A\pm B)}{(\gamma _+\pm \gamma _{})^3},$$
(62)
the RG flow can be summarized by means of the two equations

$`{\displaystyle \frac{dg_+}{d\mathrm{}}}`$ $`=`$ $`g_+\left([4d+y]3𝒮g_+\right),`$ (63)
$`{\displaystyle \frac{dg_{}}{d\mathrm{}}}`$ $`=`$ $`g_{}\left([4d+y]3𝒮g_{}\right),`$ (64)
which shows that flow takes place in a two-dimensional coupling or parameter space with axes $`(g_+,g_{})`$. Here, $`𝒮=4`$ is the same numerical symmetry factor as mentioned above.

By shutting off the one-loop corrections (by putting $`\lambda _0=0`$) in these equations Eq. (56), one immediately recovers the naive scaling in Eq. (29). This demonstrates the importance of the nonlinear interaction for scaling behavior. We see that the scale dependence of the viscosity and resistivity is modified from its naive behavior. However, the noise amplitudes and the formal expansion parameter do not get modified at one-loop; see Appendices C and D for further details.

## IV Fixed point solutions and exponents

The solutions of the differential RG equations Eq. (56) (or equivalently Eqs. (63) and (64)) tell us how the parameters appearing in the MHD equations Eqs. (20) and (21) change with length and time scale. At a fixed point, as the name implies, the parameters no longer change under a further RG transformation, and one is instead in the scaling regime, where the dynamics is scale invariant. To find the renormalization group fixed points, we set the differential RG equations to zero, which yields algebraic equations.

The RG equations Eqs. (63) and (64) have the two fixed point solutions:

$`(i)(g_+^{},g_{}^{})`$ $`=`$ $`(0,0)\mathrm{and},`$ (65)
$`(ii)(g_+^{},g_{}^{})`$ $`=`$ $`({\displaystyle \frac{4d+y}{12}},{\displaystyle \frac{4d+y}{12}}).`$ (66)
We wish to emphasize the fact that $`g_\pm `$ provides the bona-fide expansion parameter: it is dimensionless and is less than one whenever $`(4d+y)<12`$. This latter condition is satisfied for all the noise spectra and spatial dimensions considered in this paper. The stability properties of each fixed point is determined by performing a linear analysis about each point. That is, we expand the couplings $`g_\pm =g_\pm ^{}+\delta g_\pm `$ where $`\delta g_\pm `$ is a small fluctuation, and then substitute this expansion into Eqs. (63) and (64) to obtain the RG equations for the linearized perturbations about the fixed points Eqs. (65) and (66). This yields

$`(i){\displaystyle \frac{d\delta g_+}{d\mathrm{}}}`$ $`=`$ $`(4d+y)\delta g_++O(\delta g_+^2),`$ (67)
$`(ii){\displaystyle \frac{d\delta g_+}{d\mathrm{}}}`$ $`=`$ $`(4d+y)\delta g_++O(\delta g_+^2).`$ (68)
Identical equations hold for $`\delta g_{}`$ by simply replacing $`\delta g_+\delta g_{}`$ in these equations. These results tell us that when $`(4d+y)>0`$, the trivial fixed point Eq. (65) is infrared unstable and repulsive (since the perturbation grows for $`\mathrm{}\mathrm{}`$) while the non-trivial fixed point Eq. (66) is infrared stable and attractive. What this means physically is the following. We imagine probing the system at a given initial length and time scale. There is a corresponding set of parameters $`\lambda _0,A,B,\gamma _+,\gamma _{}`$ associated to this scale. Geometrically, this can be represented by a point in (a five-dimensional) parameter space. Then, tracing the effective dynamics of the system at ever larger length and time scales generates a directed path in this space which starts off from this initial point. This flow therefore is always driven away from the repulsive fixed point and towards the attractive fixed point, which therefore acts as a “sink”. Therefore, the asymptotic dynamics is uniquely determined by the attractive fixed point. If $`(4d+y)<0`$, then the stability properties of these two fixed points are interchanged. In this case, all RG flows, irrespective of where they start off in parameter space, inevitably end up in the neighborhood of the trivial fixed point $`\lambda _0=0`$. This implies the asymptotic dynamics is linear (and trivial). We see the stability properties of fixed points depend on the spatial dimension $`d`$ and on the noise spectrum through the noise exponent $`y`$. For the Kolmogorov spectrum $`y=d`$ and for the Iroshnikov-Kraichnan spectrum $`y=d1/4`$, the trivial fixed point Eq. (65) is always repulsive while the non-trivial fixed point Eq. (66) is always attractive in the infrared. Similar conclusions hold for a white noise spectrum, $`y=0`$, provided the dimension $`d<4`$. Finally, the fixed point critical exponents are obtained by substituting the values of the fixed points into the original set of RG equations Eq. (56), which yields the solutions

$`(z,\chi )`$ $`=`$ $`(2,{\displaystyle \frac{1}{2}}(2+yd)),\mathrm{for}\mathrm{fixed}\mathrm{point}(\mathrm{i}),`$ (69)
$`(z,\chi )`$ $`=`$ $`({\displaystyle \frac{1}{3}}(2+dy),{\displaystyle \frac{1}{3}}(1+yd)),\mathrm{for}\mathrm{fixed}\mathrm{point}(\mathrm{ii}).`$ (70)
The exponents for the trivial fixed point are simply the generalization of those of the Edwards-Wilkinson model of surface growth subject to spatially correlated noise (Barabási and Stanley 1995). Note furthermore that the exponents satisfy $`\chi +z=1`$ at the nontrivial fixed point (i.e., when $`\lambda _00`$). This in fact had better be the case as dictated by the renormalization group equation for $`\lambda _0`$ in Eq. (56) which admits a nontrivial fixed point solution iff $`z+\chi 1=0`$. So, this serves as an important consistency check on our calculation of the exponents. Finally, note that for either fixed point, the dynamic exponent is non-zero, $`z0`$.

## V Scaling and the approach to asymptotia

With the fixed point solutions obtained in the previous sections we now turn to the asymptotic scaling behavior (in the infrared) of the velocity and magnetic field correlation functions and energy density spectra.

### A Renormalized scale dependent viscosity and resistivity

The approach to the (nontrivial) fixed point can be quantified and is very useful for carrying out RG-improvement. The effective viscosity $`\nu (k)`$ and magnetic resistivity $`\nu _B(k)`$ at the scale $`k`$ result from all turbulent motions with wave numbers above $`k`$ (i.e., at shorter length scales). To derive the renormalized viscosity and resistivity, we integrate Eqs. (B22) and (B23) over a finite band of momenta from the cut-off $`\mathrm{\Lambda }`$ to $`\mathrm{\Lambda }/s`$, but without rescaling.

From taking first the sum and then the difference of these expressions, we derive the following infinitesimal recursion relations

$`{\displaystyle \frac{d\nu }{d\mathrm{}}}`$ $`=`$ $`{\displaystyle \frac{2A_d(A+B)}{\nu ^2}}\left({\displaystyle \frac{\mathrm{\Lambda }}{e^{\mathrm{}}}}\right)^{(dy4)},`$ (71)
$`{\displaystyle \frac{d\nu _B}{d\mathrm{}}}`$ $`=`$ $`{\displaystyle \frac{2A_d(AB)}{\nu _B^2}}\left({\displaystyle \frac{\mathrm{\Lambda }}{e^{\mathrm{}}}}\right)^{(dy4)},`$ (72)
where we have set $`\lambda _0=1`$ and have used $`s=e^{\mathrm{}}`$. We assume that the noise amplitudes, $`A,B`$ have achieved their asymptotic (i.e., near the fixed point) values.

Integrating Eq. (71) yields the relation

$$\nu ^3(\mathrm{})=\nu ^3(0)+\frac{6A_d(A+B)}{(4+yd)}\left((e^{\mathrm{}}\mathrm{\Lambda })^{dy4}\mathrm{\Lambda }^{dy4}\right),$$
(73)
where $`\nu (0)`$ is the initial value of the viscosity. When a large range of wavenumbers is eliminated this takes a scaling form independent of initial properties and the ultraviolet cut-off $`\mathrm{\Lambda }`$,

$$\nu (k)=\left(\frac{6A_d(A+B)}{4+yd}\right)^{\frac{1}{3}}k^{(4+yd)/3},$$
(74)
where $`k=e^{\mathrm{}}\mathrm{\Lambda }`$ and $`A+B=D_v>0`$ is the magnitude of the random force fluctuations.

By a similar argument, we find that the $`k0`$ scaling form of the magnetic resistivity is given by

$$\nu _B(k)=\left(\frac{6A_d(AB)}{4+yd}\right)^{\frac{1}{3}}k^{(4+yd)/3},$$
(75)
where $`AB=D_B>0`$ is the magnitude of the random current fluctuations.

Thus we see that the viscosity and magnetic resistivity tend to increase as $`k0`$ (that is, for large length scales) provided $`(4+yd)>0`$, which is the condition that the RG trajectories flow towards the nontrivial fixed point. This trend of increasing effective viscosity has also been reported as a prediction of the RG in pure Navier-Stokes (Dannevik et al. 1987) and in incompressible MHD (Longcope and Sudan 1991; Camargo and Tasso 1992), where in addition there is an increasing effective magnetic resistivity. These scale dependent quantities will be used repeatedly later on when we carry out the renormalization group improvement of the correlation functions and spectral densities.

### B Scaling form of the correlation functions

The scaling properties of the $`\stackrel{}{P}`$-field correlations, the $`\stackrel{}{Q}`$-field correlations as well as the cross correlations provide direct information on the large-distance and long-time behavior of the fluid velocity and magnetic field via the simple linear relations written in Eqs. (22) and (23). These correlations form part of a larger block-matrix structure as given in Eqs. (A34), (A40) and (A36), but we can easily work out the scaling properties of the individual elements of this block array and further note from direct inspection of those expressions that $`P_jP_n=Q_jQ_n`$ and $`P_jQ_n=Q_jP_n`$. So, there are really only two independent types of field correlation functions to consider, whether these are expressed in terms of the Elsasser variables or in terms the velocity and magnetic fields. In this Section, we derive the exact asymptotic scaling form of the correlation functions up to (as of yet) an unknown scaling function. This illustrates the power and limitations of pure scaling arguments. After the inclusion of magnetic helicity, we return to the correlation functions and will determine their scaling functions explicitly by means of RG-improved perturbation theory.

It will prove extremely useful to have the scaling form of these correlations available in both the space and time domain as well as in the Fourier (momentum and frequency) domain. In the scaling regime, the $`\stackrel{}{P}`$-field scales as

$$P_j(\stackrel{}{x},t)=s^\chi P_j^{}(\stackrel{}{x}^{},t^{})=s^\chi P_j^{}(s^1\stackrel{}{x},s^zt)=s^\chi P_j(s^1\stackrel{}{x},s^zt),$$
(76)
and similarly for $`\stackrel{}{Q}`$. The middle expressions are direct consequences of the scale transformations Eqs. (24) and (26) whereas in the scaling regime, the MHD equations are self-similar and then the rightmost equality also holds. Thus in the scaling regime the $`PP`$ correlation function scales as

$`P_j(\stackrel{}{x},t)P_n(0,0)`$ $`=`$ $`s^{2\chi }P_j(s^1\stackrel{}{x},s^zt)P_n(0,0),`$ (77)
$`C_{P_jP_n}(\stackrel{}{x},t)`$ $`=`$ $`s^{2\chi }C_{P_jP_n}(s^1\stackrel{}{x},s^zt).`$ (78)
Here by $`C_{P_jP_n}`$, we mean the correlation function between the $`jth`$ and the $`nth`$ components of the $`P`$-field. Note that because of space and time translational invariance, the correlation for the fields depend only on the spatial separation $`\stackrel{}{x}`$ and temporal separation $`t`$. The above equality relates the correlation function at two different scales, whose separation is measured by $`s>1`$. We emphasize that Eq. (77) is exact, when the system is in the scaling regime, that is to say, in the vicinity of a fixed point.

We now will see what further information can be extracted from the expression Eq. (77). There is an implicit tensor structure in these functions, due to the vector indices, which we make explicit below in the Fourier domain. Consider spatial correlations at equal times (i.e., $`t=0`$) then from Eq. (77) we have that

$$C_{P_jP_n}(\stackrel{}{x},0)=|\stackrel{}{x}|^{2\chi }\widehat{C}_{P_jP_n}(1,0),$$
(79)
whereas temporal correlations at the same spatial point (i.e., $`\stackrel{}{x}=0`$) scale as

$$C_{P_jP_n}(0,t)=t^{\frac{2\chi }{z}}\widehat{C}_{P_jP_n}(0,1),$$
(80)
where we have chosen $`s=|\stackrel{}{x}|`$ in the former and $`s=t^{1/z}`$ in the latter. Here, $`\widehat{C}_{P_jP_n}(1,0),\widehat{C}_{P_jP_n}(0,1)`$ are just constant tensor coefficients. The general scaling form of the correlation function valid for arbitrary spatial and temporal arguments can be written as

$$C_{P_jP_n}(\stackrel{}{x},t)=|\stackrel{}{x}|^{2\chi }\widehat{C}_{P_jP_n}\left(\frac{t}{|\stackrel{}{x}|^z}\right),$$
(81)
where the scaling function $`\widehat{C}`$ obeys the following limits:

$`u0\widehat{C}_{P_jP_n}(u)`$ $``$ $`\widehat{C}_{P_jP_n}(1,0),`$ (82)
$`u\mathrm{}\widehat{C}_{P_jP_n}(u)`$ $``$ $`\widehat{C}_{P_jP_n}(0,1)u^{\frac{2\chi }{z}}.`$ (83)
Note that Eqs. (79) and (80) are just special cases of Eq. (81).

Identical results hold for the $`QQ`$ correlation function by simply replacing the index labels $`PQ`$ from Eq. (77) to Eq. (83). For the cross-correlations $`PQ`$, just replace one of the $`P`$ symbols in Eq. (77) to Eq. (83) by a $`Q`$. This is as far as one can go employing scaling arguments. However, this is enough for determining how the power-law correlations depend on the exponents $`z`$ and $`\chi `$, which are obtained from solving for the RG fixed points.

Fourier transforming Eq. (77) yields equivalent information, but expressed instead in the momentum and frequency domain

$$C_{P_jP_n}(\stackrel{}{k},\omega )=k^{d2\chi z}𝐏_{jn}(\stackrel{}{k})\widehat{C}_{PP}\left(\frac{\omega }{k^z}\right),$$
(84)
with of course, the scaling functions $`\widehat{C}_{PP}=\widehat{C}_{QQ}`$. The scaling of the cross-correlation function involves an a-priori independent scaling function, and thus we have

$$C_{P_jQ_n}(\stackrel{}{k},\omega )=k^{d2\chi z}𝐏_{jn}(\stackrel{}{k})\widehat{C}_{PQ}\left(\frac{\omega }{k^z}\right).$$
(85)
Note that Eqs. (84) and (85) show that the tensor structure of the correlation functions is carried by the projection operator, as was to be expected.

It should be clear from these considerations that the asymptotic behavior of the correlation functions and the quantities derived from them, are all completely determined in terms of the two exponents $`z`$ and $`\chi `$ calculated above. Thus, using Eqs. (22) and (23) to transform back to the physical velocity and magnetic fields, we immediately find in the space and time domain

$`v_j(\stackrel{}{x},t)v_n(0,0)`$ $`=`$ $`|\stackrel{}{x}|^{2\chi }\left\{\widehat{C}_{P_jP_n}\left({\displaystyle \frac{t}{|\stackrel{}{x}|^z}}\right)+\widehat{C}_{P_jQ_n}\left({\displaystyle \frac{t}{|\stackrel{}{x}|^z}}\right)\right\},`$ (86)
$`B_j(\stackrel{}{x},t)B_n(0,0)`$ $`=`$ $`|\stackrel{}{x}|^{2\chi }\left\{\widehat{C}_{P_jP_n}\left({\displaystyle \frac{t}{|\stackrel{}{x}|^z}}\right)\widehat{C}_{P_jQ_n}\left({\displaystyle \frac{t}{|\stackrel{}{x}|^z}}\right)\right\},`$ (87)
$`v_j(\stackrel{}{x},t)B_n(0,0)`$ $`=`$ $`0.`$ (88)
Note there is no net averaged cross-helicity. This is consistent with the fact that no helicity has been injected into the system.

From Eqs. (86) - (88), the spatial dependence of both fluid velocity and magnetic field correlations in the vicinity of the non-trivial IR stable fixed point scale as

$$v_j(\stackrel{}{x},0)v_n(0,0)B_j(\stackrel{}{x},0)B_n(0,0)r^{\frac{2}{3}(1+yd)},$$
(89)
where $`r=|\stackrel{}{x}|`$ and the roughness exponent $`\chi =\frac{1}{3}(1+yd)`$ at the infrared stable fixed point Eq. (66). Note that for the Kolmogorov spectrum, $`y=d`$, and we recover the two-thirds law ($`r^{\frac{2}{3}}`$) for the velocity correlation function in fully developed hydrodynamic turbulence (Frisch 1995) (and we predict that the spatial dependence of the magnetic field correlation function is a power-law obeying a two-thirds law for this noise spectrum). For the Iroshnikov-Kraichnan spectrum, the spatial correlation scales instead as $`r^{\frac{1}{2}}`$. For the temporal correlations, from Eqs. (86) - (88) they scale as

$$v_j(0,t)v_n(0,0)B_j(0,t)B_n(0,0)t^{2\chi /z}=t^{\frac{2(1+yd)}{(2+dy)}},$$
(90)
where the dynamic exponent $`z=\frac{1}{3}(2+dy)`$. For the Kolmogorov spectrum, the temporal correlations scale as $`t`$, whereas for the I-K spectrum, they scale as $`t^{6/7}`$. By marked contrast, a white noise spectrum $`y=0`$ yields spatial correlations that scale as $`r^{\frac{4}{3}}`$ and temporal correlations that scale as $`t^{\frac{4}{5}}`$ in $`d=3`$ dimensions.

To help understand the above results about the correlation functions in the scaling regime, the case of $`d=3`$, $`y=2`$ hydrodynamics will be computed explicitly. The above scaling relations of the correlation functions all apply to hydrodynamics by setting the magnetic field to zero, $`\stackrel{}{B}=0`$. For this case $`4d+y<0`$, thus the trivial fixed point is stable and so it is simple to obtain the exact correlation functions by inverse Fourier transforming Eq. (53) to yield

$`C_{PP}^{(y=2)}(\stackrel{}{x},t)`$ $`=`$ $`const.\times \left[\mathrm{exp}({\displaystyle \frac{|\stackrel{}{x}|^2}{4\nu t}})\right]/(\nu t)^{3/2}`$ (91)
$`=`$ $`const.\times {\displaystyle \frac{1}{|\stackrel{}{x}|^3}}\left[({\displaystyle \frac{|\stackrel{}{x}|^2}{\nu t}})^{3/2}\mathrm{exp}({\displaystyle \frac{|\stackrel{}{x}|^2}{4\nu t}})\right].`$ (92)
Let us verify this solution is consistent with the general scaling forms given above. Since at this trivial fixed point for $`d=3`$, $`y=2`$ from Eq. (69) $`z=2`$, $`\chi =3/2`$, it follows that the above solution is consistent with the general scaling forms Eqs. (77) and (81). Furthermore, from the general scaling form Eq. (80) one expects $`C_{PP}(0,t)t^{3/2}`$ which is consistent with the above solution. For $`t=0`$, the inverse Fourier transform is ill-defined, so comparison is not relevant. This is not inconsistent with the general scaling form Eq. (79) since these general considerations do not guarantee that the function $`C_{PP}(1,0)`$, for example in Eq. (79), necessarily is well defined.

### C Energy spectrum

A quantity of fundamental importance in turbulence research is the energy spectrum which contains, in part, information regarding the flow of energy across different length scales. The energy spectrum also is key in defining and calculating the spectral energy density, which is the energy density per mode or wavenumber. The transfer of energy from smaller to larger (or larger to smaller) length scales as a system evolves in time is known as an inverse (direct) energy cascade. We can investigate whether and under what conditions, randomly forced stochastic MHD in the absence of magnetic helicity, subject to the noise spectrum of Eqs. (5) and (6), will exhibit energy cascades, and more generally, what is the asymptotic behavior of the energy spectrum. To do so, we will begin by focusing attention on the quadratic energy densities for the Elsasser fields and make direct use of the asymptotic fixed point scaling properties for the correlation functions deduced in the previous subsection. After this is done, we then turn to a complete calculation of the spectral energy density functions for both the physical fluid velocity and magnetic field by means of improved perturbation theory using the results obtained in Subsection V A for the renormalized viscosity and magnetic resistivity. This will lead to identical scaling relations as obtained from the “pure” or naive scaling arguments, but with the added advantage of yielding the explicit form of the individual scaling functions, valid to the one-loop order in perturbation theory at which we are working. Moreover, the example will serve to demonstrate how RG-improvement works in practice.

Consider therefore the total energy $`E_{PP}`$ in the Elsasser $`\stackrel{}{P}`$ field averaged over the random fluctuations (this will serve as a useful guide when we come to treat all the fields together below)

$`E_{PP}`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle d^d\stackrel{}{x}\stackrel{}{P}(\stackrel{}{x},t)\stackrel{}{P}(\stackrel{}{x},t)},`$ (93)
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{tr}{\displaystyle d^d\stackrel{}{x}P_j(\stackrel{}{x},t)P_k(\stackrel{}{x},t)},`$ (94)
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{tr}{\displaystyle d^d\stackrel{}{x}\frac{d^d\stackrel{}{k}}{(2\pi )^d}\frac{d^d\stackrel{}{k}^{}}{(2\pi )^d}\mathrm{exp}(i(\stackrel{}{k}+\stackrel{}{k}^{})\stackrel{}{x})P_j(\stackrel{}{k},t)P_k(\stackrel{}{k}^{},t)},`$ (95)
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{tr}{\displaystyle \frac{d^d\stackrel{}{k}}{(2\pi )^d}P_j(\stackrel{}{k},t)P_k(\stackrel{}{k},t)},`$ (96)
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{tr}{\displaystyle \frac{d^d\stackrel{}{k}}{(2\pi )^d}\frac{d\omega }{2\pi }\frac{d\omega ^{}}{2\pi }\mathrm{exp}(i(\omega +\omega ^{})t)P_j(\stackrel{}{k},\omega )P_k(\stackrel{}{k},\omega ^{})},`$ (97)
$`=`$ $`{\displaystyle \frac{V}{2}}\mathrm{tr}{\displaystyle \frac{d^d\stackrel{}{k}}{(2\pi )^d}\frac{d\omega }{2\pi }C_{P_jP_k}(\stackrel{}{k},\omega )},`$ (98)
where we have used $`P_j(\stackrel{}{k},\omega )P_k(\stackrel{}{k}^{},\omega ^{})=(2\pi )^{d+1}\delta ^d(\stackrel{}{k}+\stackrel{}{k}^{})\delta (\omega +\omega ^{})C_{P_jP_k}(\stackrel{}{k},\omega )`$, and the $`d`$-dimensional spatial volume factor is $`V=(2\pi )^d\delta ^d(\stackrel{}{0})`$. From the identity in the last line of Eq. (93), we can immediately identify the energy density $`(E_{PP}/V)`$, which is given by the indicated double integral of the correlation function. The spectral energy density, i.e., the energy density per unit wavenumber, is therefore given by

$$E_{PP}(k)=[\frac{1}{2}S_d/(2\pi )^{d+1}]k^{d1}_{\mathrm{}}^{\mathrm{}}𝑑\omega \mathrm{tr}C_{P_jP_k}(k,\omega ),$$
(99)
since clearly

$$E_{PP}=V_0^{\mathrm{}}𝑑kE_{PP}(k).$$
(100)
By making use of Eq. (84) which holds in the scaling regime, we can discover how the spectral energy density scales as a function of the injected random noise spectrum. Inserting the expression Eq. (84) into Eq. (99) yields

$$E_{PP}(k)=\frac{(d1)S_d}{2(2\pi )^{d+1}}k^\alpha _{\mathrm{}}^{\mathrm{}}𝑑u\widehat{C}_{PP}(u).$$
(101)
This follows after defining the exponent

$$\alpha =12\chi =\frac{5}{3}+\frac{2}{3}(dy),$$
(102)
and making use of the exponent identity $`z+\chi =1`$, which holds at the non-trivial infrared stable fixed point Eq. (66). The steps leading to the derivation of the spectral energy density function $`E_{QQ}(k)`$ for the $`Q`$-field are identical to those outlined above, and in fact, we have $`E_{PP}(k)=E_{QQ}(k)`$, as a consequence of the equality between the scaling functions $`\widehat{C}_{PP}=\widehat{C}_{QQ}`$. We remark at this point that the choice of noise exponent $`y=d`$ formally leads to the Kolmogorov energy density spectrum, since this yields $`\alpha =\frac{5}{3}`$, and the spectral energy density then scales as

$$E_{PP}(k)k^{\frac{5}{3}}.$$
(103)
By way of contrast, at the unstable and trivial fixed point, the spectral density would scale instead as

$$E_{PP}(k)k^3$$
(104)
for $`y=d`$. The difference in the scaling is due to the presence Eq. (103) or absence Eq. (104) of the nonlinear terms in the equations of motion which act to mix distinct modes.

We note of course that the overall magnitude of the spectral density Eq. (101) is determined in part by the scaling function $`\widehat{C}_{PP}`$. It has been common practice in many, if not all, RG analyses of turbulence to deduce only the power law scaling exponent of the physical quantities of interest while leaving the overall proportionality factors undetermined. This could lead to the impression that the RG is incapable of determining these overall factors. In fact, the renormalization group not only yields a determination of the scaling exponents but also provides quantitative information regarding the overall constant factors. This is achieved to a given order in perturbation theory by means of RG-improved perturbation theory (Hochberg et al. 1999a). We will make use of this technique to complete our discussion of the spectral energy density. We now turn to this task, which will be carried out in two steps. First we calculate the energy corresponding to the free-field limit, where the correlation functions are known explicitly, and then we use RG to improve the resulting energy by use of our perturbation theory. This results in the corrected form of the spectral energy, to the same order of perturbation theory at which we are working. We begin again at the level of the Elsasser fields and then at the end of the calculation make use of the linear transformation back to the physical fluid velocity and magnetic fields. The net result will be a one-loop expression for both the fluid velocity and magnetic field spectral energy densities. The calculation of $`E_{PP}(k)`$ carried out in detail above will serve as a useful template for this purpose.

We define an energy function involving both the Elsasser fields $`P`$ and $`Q`$ at once (which is actually a $`2\times 2`$ matrix, as we indicate with the dyad)

$`\stackrel{}{E}`$ $`=`$ $`{\displaystyle \frac{V}{2}}\mathrm{tr}{\displaystyle }{\displaystyle \frac{d^d\stackrel{}{k}}{(2\pi )^d}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\omega }{2\pi }}[\stackrel{}{𝐂}(\stackrel{}{k},\omega )]_{mn},`$ (105)
$`=`$ $`V{\displaystyle _0^{\mathrm{}}}𝑑k\stackrel{}{E}(k),`$ (106)
where the energy spectral density function (also a $`2\times 2`$ matrix) is given by

$$\stackrel{}{E}(k)=[\frac{1}{2}S_d/(2\pi )^{d+1}]k^{d1}_{\mathrm{}}^{\mathrm{}}d\omega \mathrm{tr}[\stackrel{}{𝐂}(\stackrel{}{k},\omega )]_{mn}.$$
(107)
This expression is exact. We can calculate it to one-loop order by means of RG-improved perturbation theory. This proceeds by the two above mentioned steps. That is to say, first we evaluate this spectral density using the free-field correlation function (denoted with the zero subscript)

$$[\underset{0}{\overset{}{𝐂}}(\stackrel{}{k},\omega )]_{mn}=𝐏_{mn}(\stackrel{}{k})𝒢_0(\stackrel{}{k},\omega )D_0(k)𝒢_0(\stackrel{}{k},\omega ),$$
(108)
and then improve the resulting expression by replacing the bare parameters by their one-loop asymptotic scaling forms (in the present case, only the viscosity $`\nu `$ and magnetic resistivity $`\nu _B`$ receive nontrivial corrections at one-loop order, so these are the only MHD parameters that get improved to this order). The free correlation function appearing above is defined and calculated in Eqs. (A40) and (A36), where further details can be found. The frequency integration is straightforwardly performed by means of the residue theorem. From Eqs. (A40) and (A36), we immediately see that there are two simple poles in the lower complex frequency plane and two simple poles in the upper half plane. Closing the contour in the lower half plane and taking the radius of the semicircular contour to infinity yields the desired integral Eq. (107), and we obtain (at zero-loop order, that is, for the free non-interacting theory)

$`\underset{0}{\overset{}{E}}(k)`$ $`=`$ $`\left({\displaystyle \frac{d1}{8}}\right){\displaystyle \frac{S_d}{(2\pi )^d}}k^{d3}D_0(k)\left\{{\displaystyle \frac{1}{\nu }}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{1}{\nu _B}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}`$ (109)
$`=`$ $`\left({\displaystyle \frac{d1}{4}}\right){\displaystyle \frac{S_d}{(2\pi )^d}}k^{3+(dy)}\left\{{\displaystyle \frac{A+B}{\nu }}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{AB}{\nu _B}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}.`$ (110)
The first line holds for arbitrary spatially correlated noise $`D_0(k)`$. We now renormalization group improve this expression by replacing the viscosity and magnetic resistivity by their (one-loop) asymptotic scale dependent forms $`\nu (k)`$ and $`\nu _B(k)`$ calculated above in Eqs. (74) and (75) to obtain the one-loop spectral density

$`\stackrel{}{E}(k)`$ $`=`$ $`({\displaystyle \frac{d1}{4}}){\displaystyle \frac{S_d}{(2\pi )^d}}\left({\displaystyle \frac{6A_d}{4+yd}}\right)^{\frac{1}{3}}k^{\frac{5}{3}+\frac{2}{3}(dy)}`$ (111)
$`\times `$ $`\left\{(A+B)^{\frac{2}{3}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+(AB)^{\frac{2}{3}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}.`$ (112)
This has been evaluated using the specific noise spectrum as written in Eq. (A40). The dimension dependent constant $`A_d`$ is given in Eq. (61). From this expression we can immediately read off the individual spectral energy functions for each Elsasser variable since each entry of this array Eq. (111) involves a quadratic pairing of the Elsasser fields, that is

$$\stackrel{}{E}(k)=\left(\begin{array}{cc}E_{PP}(k)& E_{PQ}(k)\\ E_{QP}(k)& E_{QQ}(k)\end{array}\right).$$
(113)
By means of the transformations given in Eqs. (22) and (23) we readily obtain the spectral energy densities for the physical fluid velocity and magnetic fields,

$`E_V(k)`$ $`=`$ $`\left\{E_{PP}(k)+E_{PQ}(k)\right\},`$ (114)
$`E_B(k)`$ $`=`$ $`\left\{E_{PP}(k)E_{PQ}(k)\right\},`$ (115)
respectively. With these expressions in hand, we can now obtain in detail the relative influence of kinetic and magnetic noise. From Eqs. (111),(114), and (115), we deduce the individual spectral energy densities, namely

$$E_V(k)=(\frac{d1}{2})\frac{S_d}{(2\pi )^d}k^{\frac{5}{3}+\frac{2}{3}(dy)}\left(\frac{6A_d}{4+yd}\right)^{\frac{1}{3}}(A+B)^{\frac{2}{3}},$$
(116)
for the kinetic energy spectrum and

$$E_B(k)=(\frac{d1}{2})\frac{S_d}{(2\pi )^d}k^{\frac{5}{3}+\frac{2}{3}(dy)}\left(\frac{6A_d}{4+yd}\right)^{\frac{1}{3}}(AB)^{\frac{2}{3}},$$
(117)
for the magnetic energy density spectrum. The ratio of the magnetic to kinetic energy spectra is given by

$$\frac{E_B(k)}{E_V(k)}=\left(\frac{D_B}{D_v}\right)^{\frac{2}{3}},$$
(118)
demonstrating that this ratio depends only on the ratio of the magnitudes of magnetic and kinetic noises.

As promised, we have succeeded in calculating the scaling functions for both the fluid velocity and magnetic field spectral functions. To one-loop, these results may be read off from Eqs. (116) and (117). Moreover, the difference in the two spectral functions is due solely from the dependence on the noise amplitudes. Thus for example, in the limit of vanishing magnetic noise, we have $`A=B`$. In this limit, the magnetic resistivity Eq. (75) does not renormalize and moreover, the magnetic spectral energy density Eq. (117) vanishes identically. Conversely, in the limit of zero kinetic noise, the fluid viscosity Eq. (74) would not renormalize and the kinetic energy density Eq. (116) would vanish identically. Another observation is that when both kinetic and magnetic noises are present, both spectral functions scale in the same way.

The scaling behavior of the energy spectrum has a special interest due to the significance of the Kolmogorov spectrum as well as others, such as the Iroshnikov-Kraichnan spectrum. In hydrodynamics, for a freely decaying turbulent fluid, which in particular means without an external driving force, Kolmogorov gave general arguments that the energy spectrum would behave as $`k^{5/3}`$. For forced hydrodynamics and/or MHD, it is clear from Eq. (111) that the energy spectrum can have a range of behavior, with the $`5/3`$ form not being unique. In particular, for a specific choice, namely $`y=d`$, the Kolmogorov spectrum emerges. However, since the validity of the present RG analysis only is in the asymptotic region, this spectral behavior only occurs in the $`k0`$ limit. This differs from the region relevant for freely decaying turbulence, which extends from some intermediate value of $`k`$ up to the higher $`k`$ regime. Nevertheless, an interesting hypothesis can be motivated based on Kolmogorov’s universality arguments for this spectrum. If the energy spectrum in the forced case has the same behavior as found for freely decaying turbulence, then the asymptotic limit of other hydrodynamic quantities, such as the Prandtl number, Reynolds number, and the viscosity also should be the same in both cases. Following this hypothesis, it would allow the stochastic hydrodynamic equations for $`y=d`$ to be used as a tool for studying various properties about freely decaying turbulence. Similarly, other types of noise (e.g., Iroshnikov-Kraichnan) could be studied, by choosing the noise exponent $`y`$ accordingly. In this paper, this direction of reasoning will not be developed or applied, but it is worth noting it here.

## VI Magnetic Helicity

Up to now we have been considering the effects of both random forces and random currents in turbulent MHD and the observable consequences these lead to at large spatial separations and for long times. One can also assess the influence of magnetic helicity on the asymptotic behavior of the magnetized turbulent fluid. This is especially interesting since helicity can arise naturally in rotating turbulent systems such as are encountered in astrophysics (e.g., in accretion disks. see Papaloizou and Lin 1995) and galactic dynamics. In MHD, the magnetic helicity is an independent quadratic invariant conserved in the limit of nondissipative turbulence. In three dimensions ($`d=3`$), it is given by

$$H=\frac{1}{2}d^3\stackrel{}{x}𝐀𝐁,$$
(119)
where $`𝐀`$ is the vector potential. In 2d MHD, the corresponding quadratic invariant is given instead by the square of the magnetic potential (Biskamp 1993). For hydrodynamics, the analogous quantity of interest is the kinetic helicity, but in nondissipative MHD, this is not a conserved quantity.

In much the same way as was done for the energy spectral functions derived above, one can derive a spectral density for magnetic helicity and use this to study the helicity spectrum at very large spatial and temporal scales. We will do so for $`d=3`$ dimensions. First, we note that Eq. (119) is a gauge-invariant quantity. Thus, by means of the relation between magnetic field and vector potential $`𝐁=\times 𝐀`$, in the Coulomb gauge we can write

$$𝐀(\stackrel{}{x})=\frac{1}{4\pi }d^3\stackrel{}{y}\frac{\times 𝐁(\stackrel{}{y})}{|\stackrel{}{x}\stackrel{}{y}|}.$$
(120)
With this, the helicity can be expressed entirely in terms of the magnetic field as

$$H=\frac{1}{8\pi }d^3\stackrel{}{x}d^3\stackrel{}{y}\frac{(\times 𝐁(\stackrel{}{y}))}{|\stackrel{}{x}\stackrel{}{y}|}𝐁(\stackrel{}{x}).$$
(121)
By Fourier transforming the magnetic field and the vector potential Eq. (120) to wavenumber space, it allows the helicity to be expressed as

$$H=\frac{1}{8\pi }\frac{d^3\stackrel{}{k}}{(2\pi )^3}k^2[i\stackrel{}{k}\times \stackrel{}{B}(\stackrel{}{k},t)]\stackrel{}{B}(\stackrel{}{k},t).$$
(122)
We now define a helicity density (helicity per unit spatial volume) and an associated helicity spectral density. For this, first Fourier transform the integrand in Eq. (122) to frequency variables and take the stochastic average of $`H`$ over the random fluctuations. This allows the averaged magnetic helicity to be expressed as

$$H=V_0^{\mathrm{}}𝑑k(k),$$
(123)
where $`V=(2\pi )^3\delta ^3(\stackrel{}{0})`$ is once again a spatial volume factor and the helicity spectral density is given by

$$(k)=\frac{i}{8\pi }\frac{d\mathrm{\Omega }_3}{(2\pi )^3}_{\mathrm{}}^{\mathrm{}}\frac{d\omega }{2\pi }ϵ_{jnm}k_n\{C_{P_mP_j}(\stackrel{}{k},\omega )C_{P_mQ_j}(\stackrel{}{k},\omega )\}.$$
(124)
Here, $`ϵ_{ijm}`$ is the fully antisymmetric tensor density in three dimensions. The indicated integrations are taken over the unit sphere and frequency, respectively. In arriving at this final expression, we have averaged over all the sources of random fluctuations. We also point out that $``$ is a real quantity; this can be checked directly in Eq. (122) by using $`\stackrel{}{B}(\stackrel{}{k},t)^{}=\stackrel{}{B}(\stackrel{}{k},t)`$ together with the symmetry properties of $`ϵ_{ijm}`$. In arriving at this final expression, we have written the magnetic field correlation function in terms of Elsasser field correlations $`(C_{PP},C_{PQ})`$.

It is clear that the magnetic helicity spectral density (as well as the total helicity) vanishes identically in a system subject strictly to non-helical forces and currents. This is because the resulting magnetic field correlation function is proportional to a symmetric tensor, $`<B_m(\stackrel{}{k},t)B_i(\stackrel{}{k},t)>𝐏_{mi}(\stackrel{}{k})`$. The same symmetry characteristics are of course shared by the correlation functions for the Elsasser variables. It is furthermore clear from Eq. (124) that to have a net averaged helicity, the correlation function must contain an antisymmetric contribution, and that this helical contribution must in fact be proportional to the epsilon tensor. In other words, to include helicity one must “build” an antisymmetric contribution to the correlation function using only the the metric, or Kronecker delta $`\delta _{ij}`$, products of the wavevector $`\stackrel{}{k}_m`$, and epsilon itself $`ϵ_{lnk}`$. This we will do dynamically by injecting an initial spectrum of random magnetic helicity into the MHD equations by means of a straightforward extension of the noise term. We will then track this random helicity dynamically by means of the renormalization group applied to helical MHD and compute the asymptotic limit of the helicity spectrum.

The presence of random helicity in three dimensions can be incorporated in the system by extending the noise spectrum Eq. (A33) as

$$[\stackrel{}{𝚪}(\stackrel{}{k},\omega )]_{mn}=D_0(k)𝐏_{mn}(\stackrel{}{k})iF_0(k)ϵ_{mnl}k_l,$$
(125)
where the nonhelical contribution arising from the random forces and currents is that given in Eq. (A33) and the second term involving the epsilon tensor is due entirely to helicity. In higher ($`d>3`$) spatial dimensions, the “extra” indices on the corresponding epsilon tensor would be contracted by additional factors of the momentum vector. The noise function (matrix) $`D_0(k)`$ is specified in Eq. (A40). The noise function (matrix) associated with the helical contribution can be modeled by taking the following initial power law spectrum for random helicity,

$$F_0(k)=2h_Bk^w\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right],$$
(126)
where $`h_B>0`$ is the amplitude of the helical noise and the number $`w`$ is an exponent (analogous to the exponent $`y`$ characterizing nonhelical noise). Both are free parameters. We have written the helical noise directly in the Elsasser field basis Eq. (11). In terms of the physical fields, this source Eq. (126) represents adding a helical noise term to the right hand side of Eq. (6), leaving both Eqs. (5) and (7) unaltered. One can add a random kinetic helicity to the right hand side of Eq. (5) if one wishes, by a simple modification of Eq. (126). To our knowledge, in contrast to the kinds of scaling arguments leading to the derivation of the Kolmogorov spectrum, there are no arguments that can shed light on what kind of helicity spectrum one might expect to see for example in a freely decaying turbulent helical fluid. As such, there are no specifically interesting values of $`w`$ of which to take note.

In a series of Appendices, we have set up and carried out the renormalization of the response, noise and vertex functions associated with the MHD equations in the absence of helicity, $`h_B0`$. Now that we have modified the noise spectrum in Eq. (125), we must repeat the renormalization group calculations taking into account the new helical contribution. If we refer to the elementary Feynman diagrams in Fig. (1a)-(1c), we immediately see that the only elementary diagram to be modified is the correlation function, Fig (1b). On the other hand, the bare response and vertex functions are unaffected by the presence of Eq. (126). When computing the loop diagrams, rather than Eq. (A40), one must use

$$[\underset{0}{\overset{}{𝐂}}(\stackrel{}{k},\omega )]_{mn}=𝒞_0(\stackrel{}{k},\omega )𝐏_{mn}(\stackrel{}{k})i𝒞_0^F(\stackrel{}{k},\omega )ϵ_{mnl}k_l,$$
(127)
where $`𝒞_0^F`$ is obtained from $`𝒞_0`$ in Eq. (A40) by replacing $`D_0(k)`$ by $`F_0(k)`$. In terms of the diagrammatic language, the response function diagram, Fig. (2a), and vertex function diagram, Fig. (2c), are unaffected, while the “blob” in the diagram Fig. (2b) for the correlation function is modified by the addition (insertion) of the helical noise. All the calculations in the Appendices are to be carried out with this modified diagram and setting $`d=3`$, as we are considering helical noise in three dimensions. The details and mechanics of the calculations are very similar to those presented in the Appendices. The net result of this procedure simply is that at one loop order the random helicity in Eq. (125) does not renormalize the response function, the noise spectrum nor the vertex function in the hydrodynamic limit, i.e., for $`k0`$ and for $`t\mathrm{}`$. In the case of the response and vertex functions, the angular integrations are responsible for yielding null results in the hydrodynamic limit. The non-helical noise spectrum does not renormalize for $`w>1`$, which is a mild condition we impose on this exponent. This means that the original set of renormalization group equations Eq. (56) (as well as their solutions, fixed points and stability properties) remain unaffected (at least at one-loop) by the presence of random magnetic helicity.

Independently from this, observe that the helicity density never contributes to the system’s total energy Eq. (105) (since helicity is traceless). Conversely, the energy density never contributes to the helicity (the energy is linked to a purely symmetric tensor). Nevertheless, the helicity has an impact on the dynamics of the system, generating its own contribution to the velocity and magnetic field correlations, as is clear from Eq. (127). We come back to this important point after calculating the helicity spectral density by means of RG-improved perturbation theory at one-loop.

### A Helicity spectral density

In analogy with the energy spectrum, we can calculate the helicity spectrum for both, non-interacting and interacting dynamics. To this end, we define the following helicity (matrix) function (in three dimensions)

$`\stackrel{}{H}`$ $`=`$ $`{\displaystyle \frac{V}{8\pi }}i{\displaystyle }{\displaystyle \frac{d^3\stackrel{}{k}}{(2\pi )^3}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\omega }{2\pi }}k^2ϵ_{jmn}k_j[\stackrel{}{𝐂}(\stackrel{}{k},\omega )]_{mn},`$ (128)
$`=`$ $`V{\displaystyle _0^{\mathrm{}}}𝑑k\stackrel{}{}(k)`$ (129)
where the helicity spectral density function (also a $`2\times 2`$ matrix) is given by

$$\stackrel{}{}(k)=\frac{i}{8\pi }\frac{d\mathrm{\Omega }_3}{(2\pi )^3}_{\mathrm{}}^{\mathrm{}}\frac{d\omega }{2\pi }ϵ_{jnm}k_n[\stackrel{}{𝐂}(\stackrel{}{k},\omega )]_{mj}.$$
(130)
We need to know the correlation function (for the Elsasser fields) $`[\stackrel{}{𝐂}(\stackrel{}{k},\omega )]_{mj}`$ in order to determine the helicity spectrum. We will obtain this spectrum by means of (1-loop) RG-improved perturbation theory. So, we first evaluate the spectrum Eq. (130) in the non-interacting field limit and then improve the resultant expression in order to obtain the corrected spectrum, which is valid at one-loop. For the first step, we simply replace $`[\stackrel{}{𝐂}(\stackrel{}{k},\omega )]_{mj}`$ by its free-field limit as given in Eq. (127). Note in comparing Eq. (127) with Eqs. (128) and (130), we see that only the antisymmetric part of the full correlation function contributes to the helicity and to the helicity spectrum, which should come as no surprise. Carrying out this replacement yields for the zeroth-order helicity spectrum

$$\underset{0}{\overset{}{}}(k)=\frac{S_3}{4\pi (2\pi )^3}k^2_{\mathrm{}}^{\mathrm{}}\frac{d\omega }{2\pi }𝒞_0^F(\stackrel{}{k},\omega ),$$
(131)
where the factor $`S_3`$ results from the integration over angles (in $`d=3`$), and we must integrate

$$𝒞_0^F(\stackrel{}{k},\omega )=𝒢_0(\stackrel{}{k},\omega )F_0(k)𝒢_0(\stackrel{}{k},\omega )$$
(132)
over the frequency. This involves the identical frequency integration performed for the calculation of the energy spectrum, in arriving at Eq. (109). Thus, direct use can be made of that step, without any further calculation, to obtain

$`\underset{0}{\overset{}{}}(k)`$ $`=`$ $`{\displaystyle \frac{S_3}{(2\pi )^3}}F_0(k)\left\{{\displaystyle \frac{1}{\nu }}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{1}{\nu _B}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}`$ (133)
$`=`$ $`{\displaystyle \frac{S_3h_B}{4\pi (2\pi )^3}}{\displaystyle \frac{k^w}{\nu _B}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right],`$ (134)
for the zeroth-order value of the spectrum (which is also the spectrum corresponding to the trivial fixed point). The matrix products are worked out making use of Eq. (126). The RG-improved helicity spectrum is therefore given by

$`\stackrel{}{}(k)`$ $`=`$ $`\left(\begin{array}{cc}_{PP}(k)& _{PQ}(k)\\ _{QP}(k)& _{QQ}(k)\end{array}\right)`$ (135)
$`=`$ $`{\displaystyle \frac{h_B}{4\pi }}{\displaystyle \frac{S_3}{(2\pi )^3}}\left({\displaystyle \frac{6A_3(AB)}{1+y}}\right)^{\frac{1}{3}}k^{w+\frac{1}{3}(1+y)}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right],`$ (136)
which indicates this helicity spectral matrix refers to the Elsasser $`P,Q`$ field basis. In terms of the physical fields, the magnetic helicity spectral density is given by

$`_M(k)`$ $`=`$ $`_{PP}(k)_{PQ}(k)`$ (137)
$`=`$ $`{\displaystyle \frac{h_B}{2\pi }}{\displaystyle \frac{S_3}{(2\pi )^3}}\left({\displaystyle \frac{6A_3D_B}{1+y}}\right)^{\frac{1}{3}}k^{w+\frac{1}{3}(1+y)}.`$ (138)
We note that the other linear combination yields the kinetic helicity spectrum, but this is identically zero, $`_V=_{PP}+_{PQ}=0.`$ This is consistent with the fact that no kinetic helicity has been injected into the dynamics. We see the difference in scaling of the helical spectrum as determined from the fixed points: apart from numerical prefactors, the magnetic helicity scales as

$$_M(k)k^w,$$
(139)
for the free, non-interacting theory (trivial fixed point) while it scales as

$$_M(k)k^{w+\frac{1}{3}(1+y)},$$
(140)
at the IR stable nontrivial fixed point (interacting theory).

### B Renormalization Group Improved Correlation Functions

Just as the general scaling analysis used for obtaining the general form of the energy density spectra is rigorously complemented by improved perturbation theory, the scaling behavior derived above for the field correlation functions can be complemented by renormalization group improvement. Thus, the scaling form of the correlations displayed in Eqs. (84) and (85) is correct, as it stands. However, this only makes explicit the scaling exponent, while lumping together with the scaling function the as yet unknown proportionality constants. To obtain the latter, a bit more work is required. Once again, renormalization group improved perturbation theory leads directly to the answer.

Having introduced the magnetic helicity, we opt for improving the full correlation function Eq. (127), which contains both kinetic and helicity effects, since the same basic steps are needed for improving each contribution separately. Improvement of a particular quantity consists, as we have already seen, in replacing the bare viscosity and magnetic resistivity by their (one-loop) running expressions Eqs. (74) and (75).

From the structure in Eq. (127), we see that this is tantamount to improving the product of response functions Eq. (A36) that appear in both the $`𝒞_0`$ and $`𝒞_0^F`$ contributions, as can be checked in Eq. (A40). Since both noise matrices commute with the product of response function matrices, the problem reduces to improving the following matrix function

$`𝒢_0(\stackrel{}{k},\omega )𝒢_0(\stackrel{}{k},\omega )`$ $`=`$ $`f(\nu ,\nu _B;k,\omega )\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right],`$ (141)
$`f(\nu ,\nu _B;k,\omega )`$ $`=`$ $`{\displaystyle \frac{\gamma _+k^4}{(\gamma _++\gamma _{})^2(\gamma _+\gamma _{})^2}}C\left({\displaystyle \frac{\omega }{k^2}}\right),`$ (142)
where the dimensionless (in the free-field limit) scaling function $`C`$ is explicitly given by

$$C\left(\frac{\omega }{k^2}\right)=\frac{\frac{\omega ^2}{(\gamma _+k^2)^2}+\frac{\gamma _{}^2}{\gamma _+^2}}{\left(1+\frac{\omega ^2}{[k^2(\gamma _++\gamma _{})]^2}\right)\left(1+\frac{\omega ^2}{[k^2(\gamma _+\gamma _{})]^2}\right)}.$$
(143)
We make no special distinction between functions of $`\gamma _+,\gamma _{}`$ and functions of $`\nu ,\nu _B`$, since these viscosity and resistivity parameters are linearly related via Eq. (12), and are freely interchangeable.

From Eq. (141) and the commutativity of the block matrices, we can express Eq. (127) as

$`[\underset{0}{\overset{}{𝐂}}(\stackrel{}{k},\omega )]_{mn}`$ $`=`$ $`\left[\begin{array}{cc}PP& PQ\\ QP& QQ\end{array}\right]_{mn},`$ (144)
$`=`$ $`f(\nu ,\nu _B;k,\omega )\left\{2k^y\left[\begin{array}{cc}A& B\\ B& A\end{array}\right]𝐏_{mn}(\stackrel{}{k})+2h_Bk^w\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]ϵ_{mnl}k_l\right\}.`$ (145)
Of course, the term proportional to $`h_B`$ holds only in $`d=3`$, while the remainder of this expression can be evaluated in any spatial dimension. This yields the individual correlations in terms of the Elsasser fields, but we can skip immediately to the quadratic correlations involving the fluid velocity and magnetic fields by means of the linear relations given in Eq. (11). From this, it implies

$`v_iv_j`$ $`=`$ $`f(\nu ,\nu _B;k,\omega )(A+B)k^y\left\{𝐏_{ij}(\stackrel{}{k})+{\displaystyle \frac{h_B}{A+B}}k^{yw}ϵ_{ijm}k_m\right\},`$ (146)
$`B_iB_j`$ $`=`$ $`f(\nu ,\nu _B;k,\omega )(AB)k^y\left\{𝐏_{ij}(\stackrel{}{k})+{\displaystyle \frac{h_B}{AB}}k^{yw}ϵ_{ijm}k_m\right\},`$ (147)
$`v_iB_j`$ $`=`$ $`f(\nu ,\nu _B;k,\omega )h_Bk^wϵ_{ijm}k_m.`$ (148)
The improved correlations are obtained immediately upon making the replacements mentioned above. We see that this operation affects only the overall factor $`f(\nu ,\nu _B;k,\omega )`$, which is common to each individual correlation function. As such, this is the only function that needs to be improved, with the correlations in Eqs. (146), (147), and (148) improved by replacing

$$f(\nu ,\nu _B;k,\omega )f(\nu (k),\nu _B(k);k,\omega ).$$
(149)
From Eqs. (141), (74) and (75), we find that

$`f(\nu (k),\nu _B(k);k,\omega )`$ $`=`$ $`{\displaystyle \frac{\left(\nu (k)+\nu _B(k)\right)^2k^4}{4\nu (k)^2\nu _B(k)^2}}C\left({\displaystyle \frac{\omega }{k^z}}\right),`$ (150)
$`=`$ $`{\displaystyle \frac{(a+b)^2}{4(a^2b^2)}}k^{4+\frac{2}{3}(4+yd)}C\left({\displaystyle \frac{\omega }{k^z}}\right),`$ (151)
where the constants

$$a=\left(\frac{6A_d(A+B)}{4+yd}\right)^{\frac{1}{3}},b=\left(\frac{6A_d(AB)}{4+yd}\right)^{\frac{1}{3}}.$$
(152)
We point out that the scaling function $`C`$ is also modified by the running of the viscosity and resistivity as should be clear from Eq. (143), and this is accounted for by the exponent $`z`$.

Substituting Eq. (150) into Eqs. (146), (147), and (148) gives the complete structure of the correlation functions correct at one-loop order. Note that the net exponent $`4+\frac{2}{3}(4+yd)yd2\chi z`$, so that Eqs. (146) and (147) reproduce the same scaling (for $`h_B=0`$, i.e. zero helicity) as we found in Eqs. (84) and (85). However, we also explicitly obtain the associated scaling functions. We see that the presence of random magnetic helicity modifies these correlations, since the helical contribution scales as $`k^{yw+1}`$ relative to the non-helical part. Furthermore, the relative strength of this term goes as the ratio of the amplitudes of the helical to the non-helical fluctuations: $`h_B/(A+B),h_B/(AB)`$.

## VII The Callan-Symanzik equation

In this section, we derive the scaling behavior of the correlation functions in the Fourier domain from a distinct but entirely complementary point-of-view provided by the methods of field theory. In this approach, the renormalization of the ultraviolet (short-distance) divergences of the stochastic field theory for MHD formulated in Appendix A leads to certain partial differential equations, frequently referred to as Callan-Symanzik equations<sup>*</sup><sup>*</sup>*Equations of this type are typically referred to as Callan-Symanzik equations, though strictly speaking, this name refers to a specific type of renormalization group equation resulting from so-called mass-subtraction (in quantum field theory). A mass-independent subtraction leads to a different equation, similar to the one actually derived here, and might be more properly called a t’Hooft-Weinberg equation. A clear distinction between Gell-Mann-Low, Callan-Symanzik and t’Hooft-Weinberg renormalization group equations is drawn by (Gross 1981). for the one-particle irreducible (1PI) Green functions. The exact solutions of these equations provide considerable information that goes well beyond perturbation theory. These solutions tell us how the 1PI functions scale near the fixed points of the theory and permit one to quantify the approach to the fixed points. Moreover, the concept and mechanics of RG-improvement is given a rigorous basis by means of the Callan-Symanzik equation. We therefore devote some space in the present section to an application of these ideas to stochastic MHD field theory. In what follows, we consider non-helical MHD. The extension to the helical case is straightforward but will not be carried out here.

To begin, we first tabulate the canonical dimensions (in units of momentum) of all the parameters and fields appearing in the generating functional for MHD. These are easily obtained by noting that the action $`𝒜`$ in the dynamic functional,

$$Z=[𝒟\stackrel{}{\sigma }][𝒟𝚽]\mathrm{exp}(𝒜),$$
(153)
must be dimensionless: $`[𝒜]=0.`$ We recall that we have employed a short-hand notation for the pair of Elsasser fields: $`𝚽=(\stackrel{}{P},\stackrel{}{Q})`$, and the conjugate fields $`\stackrel{}{\sigma }`$ are defined in Appendix A. We can read off the explicit structure of the action $`𝒜`$ directly from Eqs. (A17), (A22), (A23) and (A33). From this zero-dimensionality condition we therefore find that (where $`[X]=d_X`$ denotes the dimension of $`X`$ in units of momentum)

$$[\gamma _\pm ]=d_\gamma ,[\omega ]=d_\omega =2+d_\gamma ,$$
(154)
$$[A]=[B]=d_A=d_B,[\stackrel{}{\mathrm{\Gamma }}]=d_\mathrm{\Gamma }=d_Ay,$$
(155)
$$[\lambda ]=d_\lambda =\frac{(4d+y)}{2}+\frac{3}{2}d_\gamma \frac{1}{2}d_A,[g_\pm ]=0.$$
(156)
The Elsasser and conjugate fields have dimensions

$$[𝚽]=d_\mathrm{\Phi }=3\frac{d}{2}\frac{3}{2}d_\gamma +\frac{1}{2}d_\mathrm{\Gamma },[\stackrel{}{\sigma }]=d_\sigma =1\frac{d}{2}\frac{1}{2}d_\gamma \frac{1}{2}d_\mathrm{\Gamma },$$
(157)
respectively. Note that the canonical dimension of the fields depends on the noise exponent $`y`$ since $`d_\mathrm{\Gamma }=d_Ay`$. Finally, the two 1PI functions of particular interest to us are the response function and the noise spectral function. The response function has dimension

$$[\stackrel{}{\mathrm{\Gamma }_{11}}]=d_{11}=[(\underset{0}{\overset{1}{\stackrel{}{𝐆}}})_{ij}]=(2+d_\gamma ),$$
(158)
while the non-helical part of the noise spectral function has dimension

$$[\stackrel{}{\mathrm{\Gamma }_{02}}]=d_{02}=[(\stackrel{}{𝚪})_{ij}]=d_Ay.$$
(159)
We remark that it is customary to denote the 1PI Green functions by a capital $`\mathrm{\Gamma }`$, which is not to be confused with the block array of noise correlation functions. These 1PI functions are built up out of certain numbers of conjugate fields and physical fields: the indices on the first two terms on left hand side of Eq. (158) indicate that one (1) conjugate and one (1) physical field are involved, while in Eq. (159), no conjugate field (0) and two (2) physical fields are involved. We suppress the latin indices here to avoid clutter. Note all canonical dimensions can be expressed in terms of the dimension of the viscosity (and resistivity) $`d_\gamma `$, the noise amplitude $`d_A=d_B`$, and exponent $`y`$.

The one-loop corrections to the MHD parameters calculated in the Appendices have been cut-off with a momentum or wavenumber regulator $`\mathrm{\Lambda }`$. Recall, this actually is a physical cut-off in momentum. It is associated with the smallest length scale below which the continuum description of the fluid breaks down and must be replaced by an atomic/molecular description. If this short-distance cutoff is taken to zero ($`\mathrm{\Lambda }^10`$), i.e. take $`\mathrm{\Lambda }\mathrm{}`$, it would result in short distance or ultraviolet divergences showing up in the one-loop response and correlation functions that are calculated in the Appendices. To systematically remove these ultraviolet divergences, renormalization constants would have to be introduced in the parameters and fields that appear in the above action $`𝒜`$. This may be implemented by introducing renormalized parameters (denoted by $`Z_a,Z_\mathrm{\Phi },Z_\sigma `$) and renormalized fields (denoted by an $`R`$) as (Gross 1981; Amit 1978)

$`\xi _a`$ $`=`$ $`Z_a^1\xi _a^R,`$ (160)
$`𝚽`$ $`=`$ $`Z_\mathrm{\Phi }^{1/2}𝚽^R,`$ (161)
$`\stackrel{}{\sigma }`$ $`=`$ $`Z_\sigma ^{1/2}\stackrel{}{\sigma }^R.`$ (162)
Here the symbol $`\xi _a`$ is a shorthand notation denoting the collection of all parameters that appear in the dynamical equations, $`\xi _a=(A,B,\gamma _+,\gamma _{},\lambda )`$, with the label running from $`a=1,..,5`$. We will be able to draw some general conclusions about the asymptotic properties of the correlation functions without having to explicitly calculate these renormalization parameters.

In the following, we first consider an arbitrary 1PI function, and then specialize to our response and noise functions to obtain information regarding the correlation function. Keep in mind that for the MHD theory considered here, these functions all have the structure of two-by-two block arrays multiplied by a $`d`$-by-$`d`$ tensor factor. We suppress the tensor indices in what follows. The renormalization of the 1PI functions containing $`N`$ factors of the physical field and $`\stackrel{~}{N}`$ factors of the conjugate field is given by

$$\underset{N,\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}(\stackrel{}{k},\omega ;\{\xi _a^R\},g;\mu )=Z_\mathrm{\Phi }^{N/2}Z_\sigma ^{\stackrel{~}{N}/2}\underset{N,\stackrel{~}{N}}{\overset{}{\mathrm{\Gamma }}}(\stackrel{}{k},\omega ;\{\xi _a\},g_0;\mathrm{\Lambda }).$$
(163)
In arriving at this expression, we have expressed the bare fields in terms of the renormalized fields by means of Eqs. (161) and (162). The label $`R`$ denotes a renormalized quantity. Most importantly, the renormalization requires introducing a finite but arbitrary momentum scale $`\mu `$, which remains in the finite, renormalized expressions. As it stands, this equation merely expresses the fact that the function in question is renormalizable. In other words, the UV divergences can be consistently subtracted out. This subtraction is done at some arbitrarily chosen scale $`\mu `$, and this fact must be recorded in the resulting finite renormalized function. Using the fact that the bare, unrenormalized 1PI function does not depend on this arbitrary momentum scale $`\mu `$ (it depends only on the ultraviolet cut-off $`\mathrm{\Lambda }`$), we can easily derive a differential equation for the renormalized 1PI function as follows:

$$\left(\mu \frac{d}{d\mu }\right)_0\underset{N,\stackrel{~}{N}}{\overset{}{\mathrm{\Gamma }}}(\stackrel{}{k},\omega ;\{\xi _a\},g_0;\mathrm{\Lambda })=0.$$
(164)
This implies (see footnote 1)

$$\left(\mu \frac{}{\mu }+\beta (g)\frac{}{g}+\underset{a}{}\mathrm{\Delta }_a\xi _a^R\frac{}{\xi _a^R}\frac{N}{2}\gamma \frac{\stackrel{~}{N}}{2}\stackrel{~}{\gamma }\right)\underset{N,\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}(\stackrel{}{k},\omega ;\{\xi _a^R\},g;\mu )=0,$$
(165)
where the coefficient functions appearing in this differential equation are defined by

$`\mathrm{\Delta }_a`$ $`=`$ $`\left(\mu {\displaystyle \frac{\mathrm{ln}Z_a}{\mu }}\right)_0,`$ (166)
$`\beta (g)`$ $`=`$ $`\left(\mu {\displaystyle \frac{g}{\mu }}\right)_0,`$ (167)
$`\gamma `$ $`=`$ $`\left(\mu {\displaystyle \frac{\mathrm{ln}Z_\mathrm{\Phi }}{\mu }}\right)_0,\stackrel{~}{\gamma }=\left(\mu {\displaystyle \frac{\mathrm{ln}Z_\sigma }{\mu }}\right)_0.`$ (168)
The zero-subscript (0) indicates the derivatives are taken while holding fixed the bare parameters This equation Eq. (165) indicates that a change in the arbitrary but finite subtraction scale $`\mu `$ is compensated for by corresponding changes in the physical parameters. As Eq. (165) is a linear first-order partial differential equation, it can be readily solved by means of characteristics (Amit 1978). The exact solution is given by

$$\underset{N,\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}(\stackrel{}{k},\omega ;\{\xi _a^R\},g;\mu )=\mathrm{exp}\left(_1^b\frac{du}{u}[\frac{N}{2}\gamma (u)+\frac{\stackrel{~}{N}}{2}\stackrel{~}{\gamma }(u)]\right)\underset{N,\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}(\stackrel{}{k},\omega ;\{\xi _a^R(b)\},g(b);b\mu ),$$
(169)
where the equations for the characteristics are

$`\beta (g(b))`$ $`=`$ $`b{\displaystyle \frac{g(b)}{b}},`$ (170)
$`\mathrm{\Delta }_a(g(b))`$ $`=`$ $`b{\displaystyle \frac{\mathrm{ln}\xi _a^R(b)}{b}}.`$ (171)
These equations are subject to the boundary conditions at $`b=1`$, $`g(1)=g`$ and $`\xi _a^R(1)=\xi _a^R`$. The solution can be readily verified by differentiating Eq. (169) on both sides with respect to $`b`$ and using the chain rule. The set of characteristics define a ruled (hyper) surface on which $`\underset{N\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}`$ is defined and $`b`$ parametrizes the characteristic curves on this surface. It is related to the RG flow parameter via $`b=e^{\mathrm{}}`$. This dimensionless parameter $`b`$ can be related to ratios of momentum (or length) scales, as we will see below. The physical content of the solution Eq. (169) is as follows. A change in momentum scale (as measured by $`b`$) is accompanied by a change in the coupling $`g(b)`$ and in the other parameters $`\xi _a^R(b)`$. Furthermore, the vertex function picks up an additional overall exponential factor depending on the two anomalous dimensions $`\gamma `$ and $`\stackrel{~}{\gamma }`$ that are defined and calculated in Eq. (168). In particular, we see that the anomalous dimensions depend on the renormalization of the physical and conjugate fields (wavefunction renormalization). These dimensions control in turn, the approach to the fixed point, as we will see below. In a nutshell, this equation relates the 1PI function at two different scales. We now use this fact to deduce how the correlation functions depend on scale.

Define dimensionless scaling functions ($`d_{N\stackrel{~}{N}}`$ is the canonical dimension of the 1PI function) $`\underset{N\stackrel{~}{N}}{\overset{}{F}}`$ as

$$\underset{N\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}(\stackrel{}{k},\omega ;\{\xi _a^R\},g;\mu )=\mu ^{d_{N\stackrel{~}{N}}}\underset{N\stackrel{~}{N}}{\overset{}{F}}(\frac{k}{\mu },\frac{\omega }{\mu ^{2+d_\gamma }},\{\frac{\xi _a^R}{\mu ^{d_{\xi _a}}}\};g),$$
(172)
so that the dimension of $`\underset{N\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}`$ is carried entirely by an appropriate power of $`\mu `$. The individual arguments of $`\underset{N\stackrel{~}{N}}{\overset{}{F}}`$ are rendered dimensionless since each is expressed as a dimensionless quotient, by making use of the dimensions listed above in Eq. (154). Here, $`d_{\xi _a}`$ denotes one of $`d_A,d_B,d_\lambda ,d_{\gamma _+},d_\gamma _{}`$. Then from the solution Eq. (169)

$`\underset{N,\stackrel{~}{N}}{\overset{R}{\stackrel{}{\mathrm{\Gamma }}}}(\stackrel{}{k},\omega ;\{\xi _a^R\},g;\mu )`$ $`=`$ $`\left(b\mu \right)^{d_{N\stackrel{~}{N}}}\mathrm{exp}\left({\displaystyle _1^b}{\displaystyle \frac{dx}{x}}[{\displaystyle \frac{N}{2}}\gamma (x)+{\displaystyle \frac{\stackrel{~}{N}}{2}}\stackrel{~}{\gamma }(x)]\right)`$ (173)
$`\times `$ $`\underset{N\stackrel{~}{N}}{\overset{}{F}}({\displaystyle \frac{k}{b\mu }},{\displaystyle \frac{\omega }{(b\mu )^{2+d_\gamma }}},\{{\displaystyle \frac{\xi _a^R(b)}{(b\mu )^{d_{\xi _a}}}}\};g(b)).`$ (174)
At this point, we have succeeded in writing each renormalized 1PI Green function as a product involving a dimensionful power of the momentum scale $`\mu `$ times an exponential factor times a dimensionless scaling function. Now, the correlation function for the Elsasser fields $`P,Q`$, which recall is not 1PI, is given by (see Eq. (A40)) the block-matrix product $`\stackrel{}{C}=(\underset{11}{\overset{}{\mathrm{\Gamma }}})^1\underset{02}{\overset{}{\mathrm{\Gamma }}}(\underset{11}{\overset{}{\mathrm{\Gamma }}})^1`$. Thus, from Eq. (173)

$$\stackrel{}{C}(\stackrel{}{k},\omega ;\{\xi _a^R\},g;\mu )=\left(b\mu \right)^{d_Ay2d_\gamma 4}\mathrm{exp}\left(+_1^b\frac{dx}{x}\gamma (x)\right)\widehat{C}(\frac{k}{b\mu },\frac{\omega }{(b\mu )^{2+d_\gamma }},\{\frac{\xi _a^R(b)}{(b\mu )^{d_{\xi _a}}}\};g(b)),$$
(175)
where $`\widehat{C}=(\underset{11}{\overset{}{F}})^1\underset{02}{\overset{}{F}}(\underset{11}{\overset{}{F}})^1`$ and from Eqs. (158) and (159) $`d_{02}2d_{11}=d_Ay42d_\gamma `$. Recall that $`\widehat{C}`$ is a $`2\times 2`$ block array with a pair of $`d`$-dimensional indices. Note that in Eq. (175) the dependence on the conjugate field anomalous dimension $`\stackrel{~}{\gamma }`$ has dropped out.

It is worth commenting that provided the characteristic equations Eqs. (170) and (171) are solved, Eq. (175) yields the exact solution for the correlation functions in terms of the parameter $`b`$ which is a measure of the momentum scale (or inverse length scale) at which one observes the system. The exact solution is known if explicit knowledge is available of the anomalous dimension $`\gamma (x)`$ and the solutions of the above characteristic equations. However, typically these are known only up to some order in perturbation theory, which usually is to a given number of loops. Thus, the solution of the Callan-Symanzik equation is only known to the same number of loops. The important utility of the solution lies in the fact that it tells us how quantities behave away from the fixed point and it places the results obtained from naive scaling arguments on a solid footing. In general, the beta function Eq. (167) will not vanish for arbitrary $`b`$ (i.e., far away from fixed points) nor will the coupling $`g(b)`$ be constant. As such, the expression Eq. (175) tells how the fields are correlated in general, that is, everywhere in parameter space, and not just in the vicinity of the fixed points. Of course, if the system happens to be at a fixed point, $`gg^{}=\mathrm{const}.`$, $`\gamma \gamma ^{}=\mathrm{const}.`$, etc. Choosing $`b=k/\mu `$ and restoring the block labels and tensor indices gives for the $`PP`$ entry of the block array

$$C_{P_jP_n}(\stackrel{}{k},\omega )=k^{d_Ay2d_\gamma 4}\left(\frac{k}{\mu }\right)^\gamma ^{}𝐏_{jn}(\stackrel{}{k})\widehat{C}_{PP}(\frac{\omega }{k^{2+d_\gamma }},\{\frac{\xi _a^R(k)}{k^{d_{\xi _a}}}\};g^{}),$$
(176)
and similarly for the $`PQ`$ and $`QQ`$ entries. At a fixed point, i.e. $`\mathrm{\Delta }_a^{}=\mathrm{\Delta }_a(g^{})`$,

$$\xi _a^R(b)=\xi _a^R(1)b^{\mathrm{\Delta }_a^{}},\xi _a^R(k)=\xi _a^R(1)\left(\frac{k}{\mu }\right)^{\mathrm{\Delta }_a^{}},$$
(177)
which follows from evaluating Eq. (171) at a fixed point. In particular we see that the correlation function, calculated up to $`n`$-loops, is given by the free correlation function but expressed in terms of the scale dependent (running) couplings and parameters (calculated up to $`n`$-loop order). This is nothing but the prescription leading to RG-improved quantities, but here it is a natural and automatic property of the solutions of the Callan-Symanzik equation.

Comparing (176) to Eqs. (84) and (85) yields the exponent relations

$$z=2+d_\gamma ,\chi =1\frac{d}{2}+\frac{1}{2}(d_\gamma d_A+y\gamma ^{}).$$
(178)
These relations allow the various “engineering” dimensions to be related with the fixed point exponents calculated earlier. Moreover, we see how the concept of RG-improvement used throughout this paper appears as a rigorous feature of the solutions of the Callan-Symanzik equation. The equality in Eq. (176) shows that the correlation function in the vicinity of the fixed point is a power law multiplied by the scaling function. Furthermore, the scaling function is evaluated with the running or scale-dependent values of the parameters, that are computed in the neighborhood of that fixed point, as exactly indicated by Eq. (177). Note that if there is any nontrivial wavefunction renormalization (indicated by a $`\gamma 0`$), it is automatically subsumed in the exponent $`\chi `$ as shown in Eq. (178). Since we have calculated both $`z`$ and $`\chi `$ by independent means, there is no need for a separate calculation of $`\gamma ^{}`$ if we are interested only in what is happening at the fixed point.

### A Approach to the fixed point

Observe that the exponent (also called the anomalous dimension) $`\gamma (g)`$ controls the approach to the scaling regime, which is reflected in the pure power law behavior of the correlation function at the fixed point. In general, it is unlikely the dynamical system will be sitting exactly at the fixed point. Therefore one should study the form of the correlation function for values of the coupling that are in the basin of attraction of the fixed point and see how this function depends on arbitrary initial values of the coupling.

This can be done systematically by returning to the exponential factor in Eq. (175) and carefully expanding this about the (non-trivial) fixed point. First make a convenient change of variables as

$$_1^b\frac{dx}{x}\gamma (g(x))=_0^l𝑑l\gamma (g(l))=_{g(0)}^{g(l)}\frac{\gamma (y)}{\beta (y)}𝑑y,$$
(179)
where $`b=e^l`$. Note, this $`l`$ is not the same as the RG flow parameter $`\mathrm{}`$ in Eq. (56), but they are related as $`\mathrm{ln}(b)=\mathrm{}`$ (Frey and Tauber 1994). As already pointed out, $`b`$ parametrizes the sense of “flow” along the characteristic curves that arise in the solution of the CS equation and does not automatically single out either the infrared (decreasing $`b`$) or the ultraviolet (increasing $`b`$) limit. The connection to these physical limits is made once we assign infrared and ultraviolet flow “directions” to the characteristics. The one-loop beta function can be read off directly from Eq. (63) and is given by $`\beta (g)=12g(g^{}g)`$. Solving for the coupling by integrating the differential equation Eq. (63) gives

$$\frac{g(\mathrm{})}{g^{}}=\frac{\frac{g(0)}{g^{}g(0)}\mathrm{exp}(ϵ\mathrm{})}{1+\frac{g(0)}{g^{}g(0)}\mathrm{exp}(ϵ\mathrm{})},$$
(180)
where $`ϵ=4d+y>0`$. The infrared stable non-trivial fixed point is reached by taking $`b=\frac{k}{\mu }0`$ or $`l\mathrm{}`$. In terms of the Wilsonian flow parameter $`\mathrm{}`$, this corresponds to the limit $`\mathrm{}\mathrm{}`$. Note this fixed point is reached for all initial values $`g(0)`$ of the coupling.

Next, Taylor expand the anomalous dimension about the fixed point as

$$\gamma (g)=\gamma ^{}+\gamma ^{}(gg^{})+O(gg^{})^2.$$
(181)
With this, evaluating the integral

$$_{g(0)}^{g(l)}\frac{\gamma (g)}{\beta (g)}𝑑g=\gamma ^{}l\frac{\gamma ^{}}{12}\mathrm{ln}\left(\frac{g(l)}{g(0)}\right),$$
(182)
it implies the exponential factor contributing to the correlation function Eq. (175) works out to be

$`\mathrm{exp}\left({\displaystyle _1^b}{\displaystyle \frac{dx}{x}}\gamma (g(x))\right)`$ $`=`$ $`b^\gamma ^{}\left({\displaystyle \frac{g(b)}{g(1)}}\right)^{\frac{\gamma ^{}}{12}},`$ (183)
$`=`$ $`\left({\displaystyle \frac{k}{\mu }}\right)^\gamma ^{}\left({\displaystyle \frac{\frac{g^{}}{g^{}g}\left(\frac{k}{\mu }\right)^ϵ}{1+\frac{g}{g^{}g}\left(\frac{k}{\mu }\right)^ϵ}}\right)^{\frac{\gamma ^{}}{12}}.`$ (184)
Recall that the $`k`$-dependent factor with exponent $`\gamma ^{}`$ is already subsumed into the definition of the fixed point “roughness” exponent $`\chi `$, as indicated in Eq. (178). However the second term, which depends on the exponents $`\gamma ^{}`$ and $`ϵ`$, clearly is a correction to the fixed point (power law) form of the correlation function written in Eq. (176). Depending on the algebraic sign of $`\gamma ^{}`$, this extra term, which controls the approach to the fixed point, can either augment or diminish the energy spectrum with respect to its asymptotic fixed point value Eq. (103). Similar consideration holds at the nontrivial fixed point for the asymptotic form of the helicity spectrum in Eq. (140).

Some comments are in order about the region away from the fixed point. First, any amplification or diminution in the spectra densities with respect to their steady state fixed point values, which is brought about by this extra term, is a statement about how the densities evolve as they approach the nontrivial fixed point. It is a conceivable possibility that numerical simulations, due to their spatial and temporal limitations, actually may be seeing some sufficiently long-lived non-steady state of this type. This point will be elaborated upon in Sec. IX. However, we hasten to point out that if the wavefunction renormalization constant turns out to be equal to unity, $`Z_\mathrm{\Phi }=1`$, then the anomalous dimension Eq. (168) is identically zero, $`\gamma =0`$. When this is the case, then the expontial factor in (175) reduces to unity, and any corrections to pure power law scaling will be due entirely to the scaling function itself. For incompressible hydrodynamics subject to random forces, it is known that $`Z=1`$ (DeDominicis and Martin 1979). This in turn can be shown to be a consequence of the underlying Galilean invariance of the Navier-Stokes equation, which is maintained for all spatially correlated random forces. Since the equations of nonrelativistic MHD are also Galilean invariant, we suspect that the anomalous dimension Eq. (168) will also be vanishing. This issue will not be discussed further in this paper.

## VIII Cascade Directions

One important question of general interest in the study of hydrodynamic systems is how the nonlinear terms redistribute energy and, if present, helicity amongst the various scales. In particular, if energy and/or helicity is introduced into the system at a given scale, it then will be distributed by the subsequent evolution of the hydrodynamic or MHD equations, either to smaller (direct or normal cascade) or larger (inverse cascade) length scales. In general, this discussion of cascades is restricted to the inertial range, which corresponds to length scales sufficiently large that they are not significantly dissipated by the viscosity term.

It is interesting to reflect on the history of the subject of cascades. The earliest expectation, based on simulations and general arguments of the Navier-Stokes equation in $`d=3`$, was that a direct cascade (also suggestively called a normal cascade) would occur from hydrodynamic evolution. Subsequently in (Kraichnan 1967), for the Navier-Stokes equation in $`d=2`$, based on general observations about the two conserved constants of motion, enstrophy and kinetic energy per unit mass, it was argued that an inverse cascade of energy should occur. For MHD, these arguments of inverse cascade are not applicable in either $`d=2,3`$ and so a direct cascade of energy is expected (Biskamp 1993; Kraichnan 1973). On the other hand for MHD, considerations based on conservation of magnetic helicity in $`d=3`$ and its counterpart, the square of the magnetic potential in $`d=2`$ showed that an inverse cascade is expected of the respective quantities (Frisch et al. 1975). All the above expectations are based on compelling, but not rigorous, arguments. Verification of these expectations mainly have come from computer simulations directly of the hydrodynamic/MHD equations or indirectly through some approximations to these equations (Brandenburg et al. 1996; Lilly 1969; Frisch and Sulem 1984; Sommeria 1986; Pouquet et al. 1976; Pouquet and Patterson 1978; Pouquet 1978). Clearly this approach has its limitations in the interesting regime of high Reynolds number, due to finite simulation time and system size restrictions.

Returning from this digression back to our RG analysis, cascades clearly are associated with time dependent processes. On the other hand, the RG analysis in this paper examined steady state solutions of MHD with time translationally invariant force terms and only treated the leading asymptotic behavior at large distance. From this analysis, information about the energy spectrum can be obtained for the largest length scale. In particular, this means length scales beyond the largest characteristic scale at which energy is input from the external force. As such, the RG analysis of this paper can determine presence or absence of inverse cascade, depending whether the infra-red limit of the energy spectrum is enhanced or not. This is a useful application for RG, since generally some direct cascade always occurs and the interesting question is whether there also is some inverse cascade.

In order to extract information about cascades from this framework, a few points must be appreciated. The force term is continuosly inputing energy and/or helicity into the system with some prescribed spectrum. A cascade for such a case therefore is associated with the difference in the steady state magnetic and velocity spectra directly created with the force (in the absence of mode mixing) versus the steady state spectra that emerges from the subsequent MHD evolution. For example, if the spectra associated with the force is less enhanced in the ultraviolet versus the resulting spectra after MHD evolution, then the MHD dynamics is directing energy flow to smaller length scales thus exhibiting a direct cascade. Similarily, if the input spectra associated with the force is less enhanced in the infrared versus the resulting spectra from MHD, then MHD is exhibiting an inverse cascade.

To proceed with this approach, the implied magnetic and velocity spectra associated with the force terms must be determined. The point here is Eqs. (5) - (7) give the spectra of the force, but what we seek is the associated spectra this implies for the magnetic and velocity fields. To determine these, MHD with the force terms, but without the nonlinear terms, must be examined. These linearized MHD equations express the corresponding magnetic and velocity spectra for given input force spectrum, before there is any mode mixing. Equivalently stated, these are the associated input magnetic and velocity spectra for the given force, prior to any nonlinear MHD evolution.

The procedure then for determining information about cascades is to compare for a given force term the corresponding magnetic and velocity spectra from the linearized MHD equations versus those from the full nonlinear MHD equations. Depending on how the ultraviolet and infrared behavior of the latter changes compared to the former, cascade directions can be deduced.

Following this criteria, let us determine the cascade directions predicted at the largest length scales from our RG analysis. First note that the spectra for linearized MHD equivalently stated is the one also for the trivial fixed point. Therefore one fact immediately follows from this criteria. In the region of stability of the trivial fixed point, the RG analysis predicts no cascade whatsoever, at the largest length scales.

Turning next to energy spectra in MHD, from Eqs. (101) and (102) $`E(k)k^\alpha `$ with $`\alpha _{\mathrm{linear}}=3+dy`$ and $`\alpha _{\mathrm{nontrivial}}=5/3+(2/3)(dy)`$, so that $`\alpha _{\mathrm{linear}}()\alpha _{\mathrm{nontrivial}}`$ in the region $`4d+y()0`$. Thus we see that in the stability region of the nontrivial fixed point, the input energy spectra $`k^{\alpha _{\mathrm{linear}}}`$ is more enhanced in the infrared and less in the ultraviolet versus the resulting spectra from MHD evolution $`k^{\alpha _{\mathrm{nontrivial}}}`$. Thus our RG analysis of MHD predicts direct energy cascades at the largest length scales. Since the above case involves the entire stability region of the nontrivial fixed point, this analysis finds at the largest length scales, no prediction of an inverse cascade of energy in MHD for any $`d`$.

The RG calculations in this paper easily can be reduced from MHD to hydrodynamics by identifying the fields $`P=Q`$. The outcome of this is the predictions of the fixed points for the hydrodynamic case are found to be exactly the same as given in Sec. IV. As such, our analysis also expects direct but no inverse energy cascades at the largest length scales in hydrodynamics for any $`d`$.

It is interesting to note that the change from direct cascade to inverse cascade in the RG analysis occurs right at the boundary of stability for the nontrivial fixed point. In particular the unstable region of the nontrivial fixed point would predict inverse cascade. In the Conclusion we will discuss interpretations of this observation.

The same analysis as above also can be done for helicity. The helicity spectra given in Subsect. VI A are $`_M(k)k^\beta `$ with $`\beta _{\mathrm{linear}}=w`$ from Eq. (139) and $`\beta _{\mathrm{nontrivial}}=w+(1/3)(1+y)`$ from Eq. (140). The region of direct helicity cascade is $`\beta _{\mathrm{linear}}<\beta _{\mathrm{nontrivial}}`$, which implies $`1+y>0`$. Since this is for $`d=3`$, note this last condition also corresponds to the stability region of the nontrivial fixed point, i.e., $`4d+y=1+y`$. So once again, RG predicts a direct magnetic helicity cascade in the entire stability region of the nontrivial fixed point and no cascade in the stability region of the trivial fixed point. Furthermore, once again the change from absence to presence of inverse cascade for the nontrivial fixed point is at its boundary of stability.

## IX Conclusion

In this paper a Renormalization Group (RG) analysis has been performed for MHD with a stochastic force term. Such an analysis enables the study of highly turbulent MHD in the asymptotic regime of large length and time scales. In the format of the RG approach, the central problem is to determine the effective MHD equations at large scales that emerge once all short distance physics is integrated out of the initial fundamental MHD equations. The analysis was carried out for a general class of force functions Eq. (5) - (7) and for arbitrary spatial dimension $`d`$. The outcome of this analysis, as detailed in Sec. IV, showed that MHD for this class of force functions, has two types of asymptotic behavior, depending on $`d`$ and the nature of the force. In one regime, the system becomes a free theory in the asymptotic or hydrodynamic limit which is represented by a trivial fixed point in Sec. IV , and in the other regime the system is a nonlinear interacting theory at the largest scales, which is represented by the nontrivial fixed point in Sec. IV. For both regimes, correlations were computed in Subsec. V B of the magnetic and velocity fields at large length and time scales and in Subsec. V C of the energy spectrum . In addition, the viscosity coefficients in the effective long distance theory were computed in Subsec. V A. Helical force also was treated in Sec. VI. Up to the one-loop analysis in this paper, the fixed point structure was unchanged with respect to that in Sec. IV. The helicity spectrum at both fixed points was given in Subsec. VI A. It would be interesting to test the RG predictions against simulations for all these quantities.

Another interesting direction examined was mixing between velocity and magnetic energy at large scales. The results of the one-loop analysis in Sec. V C indicate that in fact at the largest scales magnetic and velocity energy will not mix. Thus, for example, if only the velocity force was nonzero, then no magnetic energy would be created at the largest scales and visa-versa. These results of course are only concerned with the asymptotic large scales. Clearly the nonlinear interactions in the MHD equations generally should mix velocity and magnetic energy and this is seen in simulations (Pouquet et al. 1976; Pouquet and Patterson 1978), which suggest an approximate equipartition of energy. However, the RG analysis indicates that at increasingly larger length scales this mixing diminishes. It would be interesting to see if such trends also are seen in simulations.

The RG analysis of this paper addressed the steady state solutions of MHD with time translationally invariant stochastic force terms. As shown in Sec. VIII, information about cascade directions also can be deduced from this framework. To our knowledge, this is the first time the RG has been used to predict cascade directions in turbulence. In that section, we examined directions for energy cascade in both hydrodynamics and MHD and helicity cascade in MHD. In the stability region of the trivial fixed point, in all cases and in any $`d`$ the one loop RG analysis showed no cascade at asymptotic scales. On the other hand, in the stability region of the nontrivial fixed point, the RG analysis predicts for the energy in all cases and for all $`d`$, direct or normal cascades and for helicity in MHD also a direct cascade. For energy cascades in $`d=2,3`$ for MHD and $`d=3`$ for hydrodynamics, the results of the RG analysis found here are consistent with general expectations. However, for energy cascade in $`d=2`$ hydrodynamics and helicity cascade in $`d=2,3`$ MHD, some general arguments (Kraichnan 1967; Frisch et al. 1975) and simulations (Lilly 1969; Frisch and Sulem 1984; Sommeria 1986; Pouquet at al. 1976; Pouquet and Patterson 1978; Pouquet 1978) suggest an inverse cascade. There are no rigorous proofs in any of these cases. Nevertheless, the cases where the RG analysis differs from conventional expectations deserves further discussion. Below we identify some possible sources for the discrepancy and suggest ways they might be resolved. However first, let us consider the random force terms. The class of random force functions we consider are reasonably general in terms of covering a broad range of scales, and we do not believe the special form of the force is the main cause for this discrepancy. These have been taken to be power laws in wavenumber, characterized by a free exponent $`y`$ that can either be positive, zero or negative. This is enough to cover most cases of spatially correlated noise and white uncorrelated noise (if $`y=0`$).

Returning to the discrepancies, there are several possible reasons that come to mind for these disagreements between our RG analysis versus the various simulations and general arguments. We will first list them and and then address each in turn. The main possibilities are: (a) intrinsic limitations of perturbative RG applied to MHD; (b) limitations of the one-loop RG analysis which may be corrected by higher loop calculations; (c) the RG analysis is correct and the inverse cascades suggested by simulations are a metastable effect in the unstable region of the nontrivial fixed point.

Concerning point (a), we remark that it has long been known and criticized that since the formal (bare and unrenormalized) perturbation expansion parameter $`\lambda _0=1`$ in the Navier-Stokes equation is equal to unity, the (primitive) perturbation series based upon it are of dubious utility (see, e.g. discussion in McComb 1995, and the references therein). Analogous primitive perturbation series constructed for MHD would not escape this criticism either. In response, we point out that in the loop expansion carried out here, we have defined and expanded our perturbation series in powers of the amplitude of the noise fluctuations (which could be large or small). This is because the noise amplitude serves to count the number of loops in the series (represented by a graphical expansion). An expansion in $`\lambda _0`$ is not the same as an expansion in loops. The parameter $`\lambda _0`$ is only a strict bookkeeping device, or “flag”, allowing us to keep track of the nonlinear adjective mode-mixing terms. Note moreover, that $`\lambda _0`$ is a dimensionful parameter, and it makes no physical sense to expand in terms of it. Meaningful perturbation series must ultimately refer to dimensionless expansion parameters or couplings. The actual perturbative expansion parameter for MHD is not $`\lambda _0`$, but rather, $`g_\pm `$, which as we have calculated, is less than unity at both the fixed points. Thus, it is even possible to have reasonably large fluctuations and maintain a small expansion parameter, (see Eq. (62)). This coupling $`g_\pm `$ is dimensionless, and provides a bona-fide expansion parameter. So, we do not believe there to be any inherent internal inconsistency in making use of a perturbative loop calculation of hydrodynamics or MHD, provided we correctly recognize that the true (and small) expansion parameter is not $`\lambda _0`$ but rather Eq. (62).

Turning to point (b), this has to do with possible limitations of the one-loop analysis itself. The number of loops considered is intimately related to the magnitude of the stochastic fluctuations. If this number is small, then one reasonably expects a low order loop expansion should be able to capture the salient physics in the problem. However, if the fluctuations are large, one may have to push the calculation to higher loops. Once again however, the noise amplitude is not a dimensionless parameter, and what counts loops, properly speaking, is the dimensionless ratio written in Eq. (62). Note that large effective viscosity and/or magnetic resistivity can counteract the presence of large stochastic fluctuations. The number of terms or graphs in the loop expansion grows geometrically with the number of loops, so ample physical justification would be needed before taking the calculational leap. A tempting and quite interesting, though exceedingly nontrivial, proposal would be to try and go beyond the loop expansion by using nonperturbative methods to make predictions of cascade directions, following Wilson’s ideas for the exact renormalization group (Wilson and Kogut 1974). In any case, neither of these options have yet been pursued.

Finally, coming to the point raised in (c), it may very well be the case that the RG analysis is correct in all its physical predictions despite the one-loop order, and there is no actual discrepancy with the numerical results. The key to the resolution is in recognizing that the numerical simulations are always limited by the finite duration times of the computer runs and the finite size of the lattice, whereas the RG analysis here refers always to the asymptotically largest spatial and temporal scales, which are evidently not accessible to present day computer simulations. So, we are comparing phenomena (cascade directions) associated to very different spatial and temporal regimes. In other words, the inverse cascade might not be an asymptotic phenomena, but is linked instead to intermediate scale physics, and this is evidently what the simulations are seeing. One possible explanation is that simulations are seeing a low frequency mode that temporarily shifts energy to larger length scales but eventually shifts the energy to smaller length scales. The RG analysis in this paper would ignore such modes. Let us elaborate further on this suggestion. Both inverse and direct cascades are possible, but different cascade directions would come to predominate at different time scales. The complete dynamical evolution may occur in two phases: an earlier inverse cascade phase followed by a later direct cascade phase. If true, this two-phase behavior would be an example of a “cross-over” phenomenon, familiar from critical dynamics. Of course, we do not see any evidence for cross-over from our one-loop RG calculation, and only a higher-order (or preferably, non-perturbative) calculation could begin to shed some light on this possibility (see, point (b) above). To get a further hold on the problem, one could either (1) run the simulations for ever increasing times and larger lattices, or (2) study the RG away from the vicinity of the fixed points. In the latter case, we are “backing away” from the leading asymptotic limit and the terms in the equations of motion denoted as “less relevant” in Eqs. (34) and (37) would now have to be included in the RG analysis.

## Acknowledgments

We thank Axel Brandenburg and Kari Enqvist for helpful discussions. We also thank W.D. McComb for a useful discussion after the paper was completed. A.B. is funded by the United Kingdom Particle Physics and Astronomy Research Council (PPARC) and D.H. is supported in part by Spain’s National Aerospace Technology Institute (INTA) and the Centro de Astrobiología (CAB). A.B. acknowledges the kind hospitality provided him by the CAB in Madrid and D.H. acknowledges the Department of Physics of the University of Edinburgh for the hospitality provided him during a reciprocal visit. Travel to and from our respective home institutions between Madrid and Edinburgh has been made possible with funds provided in part by CAB (Spain) and PPARC (U.K.).

## A Dynamic generating functional for stochastic MHD

The purpose of this Appendix is to set up and apply a powerful and extremely useful dynamic functional integral technique that leads to an efficient, systematic and organized perturbative loop-expansion for the solutions of the fully nonlinear magnetohydrodynamic equations Eqs. (9) and (10), when the latter are subject to random forces and fluctuations. The ability to cast stochastic differential equations in terms of functional integrals opens the door to the calculation of many physically important quantities including the so-called effective action and effective potential (Hochberg et al. 1999b). These quantities can be used to investigate the onset of pattern formation and pattern selection and can be used to identify symmetry breaking states of the system. In this paper, we will use functional methods to deduce the Feynman rules needed to construct a compact graphical representation of the perturbative solutions of the coupled nonlinear MHD equations. These rules permit one to identify a distinguished class of diagrams, the so-called one-particle irreducible (1PI) graphs. The renormalization of these diagrams to a given number of loops will be used to derive the differential renormalization group equations Eq. (56), that govern the scale dependence of the parameters appearing in the MHD equations. Additional details regarding the application of functional and path integral methods to stochastic differential equations can be found in (Martin et al. 1973; DeDominicis and Peliti 1978; Rivers 1987; Zinn-Justin 1996; Hochberg et al. 1999b) and references therein. In the following, vector quantities will generally be denoted in boldface (or sometimes with an arrow) and matrix operators indicated by means of the dyad (double-arrow) notation.

Let $`𝐐_{soln},𝐏_{soln}`$ denote a general solution of the coupled stochastic MHD equations Eqs. (9) and (10). We pass immediately to a functional or path-integral representation for an arbitrary function $`F`$ of this solution-pair by means of the trivial identity

$$F(𝐐_{soln},𝐏_{soln})=[𝒟𝐏][𝒟𝐐]F(𝐐,𝐏)\delta [𝐏𝐏_{soln}]\delta [𝐐𝐐_{soln}],$$
(A1)
where the path integrals are taken over all field configurations. The product of delta-functionals selects out the particular field configurations $`𝐏=𝐏_{soln}`$ and $`𝐐=𝐐_{soln}`$ corresponding to the solutions of the MHD equations. They constrain the functional integrals which otherwise freely sum over the totality of all $`𝐐`$ and $`𝐏`$ field configurations.

The solutions of the MHD equations depend of course on the random noise sources, so that $`𝐏_{soln}=𝐏_{soln}[\stackrel{}{\eta }_P,\stackrel{}{\eta }_Q]`$, and $`𝐐_{soln}=𝐐_{soln}[\stackrel{}{\eta }_P,\stackrel{}{\eta }_Q]`$. The solutions themselves can be regarded as functionals of the noise source. Ultimately, we are only interested in the solutions of MHD averaged over the random forces and random currents. This statistical average is denoted by angular brackets and is given by

$$F(𝐐,𝐏)=[𝒟\stackrel{}{\eta }]F(𝐐_{soln},𝐏_{soln})𝒫[\stackrel{}{\eta }],$$
(A2)
where the noise distribution functional $`𝒫`$ is assumed to be given (i.e. it is a Gaussian functional for Gaussian noise) and is normalized to unity

$$[𝒟\stackrel{}{\eta }_P][𝒟\stackrel{}{\eta }_Q]𝒫[\stackrel{}{\eta }_P,\stackrel{}{\eta }_Q]=[𝒟\stackrel{}{\eta }]𝒫[\stackrel{}{\eta }]=1.$$
(A3)
We abbreviate the $`2d`$-component noise vector field by $`\stackrel{}{\eta }=(\stackrel{}{\eta }_P,\stackrel{}{\eta }_Q)`$, where $`\stackrel{}{\eta }_P`$ and $`\stackrel{}{\eta }_Q`$ each have $`d`$-components in a $`d`$-dimensional space. The goal and objective of the following sequence of manipulations consists in choosing a particular form for $`F`$ and explicitly performing all the functional integrals encountered along the way. The resulting functional then is used to generate all the statistically averaged products of solutions of the nonlinear MHD equations to any desired order in the expansion parameter $`\lambda _0`$. The Feynman rules (see Fig. 1) can be read off directly from this final functional.

We need to carry out a change of variables in Eq. (A1). To see how this works, consider first a simple example. Let $`\stackrel{}{f}(\stackrel{}{x})=(f_1(x,y),f_2(x,y))`$ be a two-component vector-valued function of a vector variable $`\stackrel{}{x}=(x,y)`$ and suppose this has a unique root at $`\stackrel{}{x}_0`$, $`\stackrel{}{f}(\stackrel{}{x}_0)=0.`$ Now make use of the functional generalization of the identity

$`\delta (\stackrel{}{f}(\stackrel{}{x}))=\delta (f_1(x,y))\delta (f_2(x,y))`$ $`=`$ $`{\displaystyle \frac{\delta (\stackrel{}{x}\stackrel{}{x}_0)}{det|\frac{f_{1,2}(x,y)}{(x,y)}|}},`$ (A4)
where

$`{\displaystyle \frac{f_{1,2}(x,y)}{(x,y)}}`$ $`=`$ $`\left(\begin{array}{cc}\frac{f_1(x,y)}{x}& \frac{f_1(x,y)}{y}\\ \frac{f_2(x,y)}{x}& \frac{f_2(x,y)}{y}\end{array}\right),`$ (A5)
and

$`\delta (\stackrel{}{x}\stackrel{}{x}_0)`$ $`=`$ $`\delta (xx_0)\delta (yy_0),`$ (A6)
to express the product of delta functionals in Eq. (A1) as

$`\delta [𝐏𝐏_{soln}]\delta [𝐐𝐐_{soln}]`$ $`=`$ $`𝒥\delta [_tP_j+\lambda _0𝐏_{jn}_l(Q_lP_n)\gamma _+^2P_j\gamma _{}^2Q_j(\eta _P)_j]`$ (A7)
$`\times `$ $`\delta [_tQ_j+\lambda _0𝐏_{jn}_l(P_lQ_n)\gamma _+^2Q_j\gamma _{}^2P_j(\eta _Q)_j].`$ (A8)
Here, the Jacobian functional determinant $`𝒥`$ is the determinant of the following $`2\times 2`$-block array:

$$\left(\begin{array}{cc}(_t\gamma _+^2)\delta _{jk}+\lambda _0𝐏_{jk}_l(Q_l)& \gamma _{}^2\delta _{jk}+\lambda _0𝐏_{jn}_k(P_n)\\ \gamma _{}^2\delta _{jk}+\lambda _0𝐏_{jn}_k(Q_n)& (_t\gamma _+^2)\delta _{jk}+\lambda _0𝐏_{jk}_l(P_l)\end{array}\right)\delta (\stackrel{}{x}\stackrel{}{x}^{})\delta (tt^{}).$$
(A9)
When the identity Eq. (A7) is substituted into Eqs. (A1) and (A2), the stochastic MHD equations Eqs. (20) and (21) appear explicitly under the functional integrals as constraints, and we see that the stochastic dynamics is built in from the outset. The next step is to replace each of delta-functionals appearing on the left-hand side of Eq. (A7) with its (functional) Fourier integral representation. That is, for any (vector-valued) space-time field $`𝚽(\stackrel{}{x},t)`$,

$$\delta [𝚽]=[𝒟\stackrel{}{\sigma }]\mathrm{exp}\left(id^d\stackrel{}{x}𝑑t\stackrel{}{\sigma }𝚽\right),$$
(A10)
where $`\stackrel{}{\sigma }`$ is the Fourier conjugate field, and has the same number of components as $`𝚽`$. Also, in the above, the dot product indicates a sum over the repeated indices, so that $`\stackrel{}{\sigma }𝚽=\sigma _j𝚽_j`$, for $`j=1,2,\mathrm{},d`$. Replacing each delta-functional by its corresponding Fourier integral introduces two conjugate field variables into the stochastic average Eq. (A2). Since there are two Elsasser fields, namely $`𝐏`$ and $`𝐐`$, each requiring a separate conjugate field, we distinguish each individual conjugate field by a superscript: $`\stackrel{}{\sigma }=(\stackrel{}{\sigma }^P,\stackrel{}{\sigma }^Q)`$.

At this stage what results is a multiple functional integral representation for the stochastic average of an arbitrary function $`F`$, which depends on the exact solutions of the original pair of MHD equations. That is

$`F(𝐐,𝐏)`$ $`=`$ $`{\displaystyle [𝒟\stackrel{}{\eta }][𝒟𝐏][𝒟𝐐][𝒟\stackrel{}{\sigma }^P][𝒟\stackrel{}{\sigma }^Q]F(𝐐,𝐏)𝒥𝒫[\stackrel{}{\eta }]}`$ (A11)
$`\times `$ $`\mathrm{exp}\left(i{\displaystyle d^d\stackrel{}{x}𝑑t\sigma _j^P[_tP_j+\lambda _0𝐏_{jn}_l(Q_lP_n)\gamma _+^2P_j\gamma _{}^2Q_j(\eta _P)_j]}\right)`$ (A12)
$`\times `$ $`\mathrm{exp}\left(i{\displaystyle d^d\stackrel{}{x}𝑑t\sigma _j^Q[_tQ_j+\lambda _0𝐏_{jn}_l(P_lQ_n)\gamma _+^2Q_j\gamma _{}^2P_j(\eta _Q)_j]}\right).`$ (A13)
Following the standard practice in both quantum and statistical field theory, now choose a convenient form for $`F`$. In particular, let $`𝐉^P,𝐉^Q,𝐉^{\sigma ^P},𝐉^{\sigma ^Q}`$ denote a set of four arbitrary vector source functions, one source for each of the Elsasser and conjugate fields, respectively. Then let

$$F=\mathrm{exp}(𝐏𝐉{}_{}{}^{𝐏}+𝐐𝐉{}_{}{}^{𝐐}+\stackrel{}{\sigma }^P𝐉{}_{}{}^{\sigma ^P}+\stackrel{}{\sigma }^Q𝐉{}_{}{}^{\sigma ^Q}),$$
(A14)
and define the generating functional $`Z`$ by

$$Z[𝐉^P,𝐉^Q,𝐉^{\sigma ^P},𝐉^{\sigma ^Q}]\mathrm{exp}(𝐏𝐉{}_{}{}^{𝐏}+𝐐𝐉{}_{}{}^{𝐐}+\stackrel{}{\sigma }^P𝐉{}_{}{}^{\sigma ^P}+\stackrel{}{\sigma }^Q𝐉{}_{}{}^{\sigma ^Q}).$$
(A15)
We have suppressed the omnipresent integrations over space and time to avoid clutter (they can be restored at any time). The source functions allow one to factor the non-linear (in this case, cubic) interaction terms out from under the functional integral, where they now appear as functional derivative operators acting on a simpler and exactly calculable generating functional. Using the functional generalization of the identity

$$𝑑ug(u)\mathrm{exp}(uJ)=g\left(\frac{}{J}\right)𝑑u\mathrm{exp}(uJ),$$
(A16)
it then follows from Eqs. (A11), (A14) and (A15), that

$`Z[𝐉^P,𝐉^Q,𝐉^{\sigma ^P},𝐉^{\sigma ^Q}]`$ $`=`$ $`\mathrm{exp}\left(i{\displaystyle d^d\stackrel{}{x}𝑑t\lambda _0\frac{\delta }{\delta J_j^{\sigma ^P}}𝐏_{jn}_l\left(\frac{\delta }{\delta J_l^Q}\frac{\delta }{\delta J_n^P}\right)}\right)`$ (A17)
$`\times `$ $`\mathrm{exp}\left(i{\displaystyle d^d\stackrel{}{x}𝑑t\lambda _0\frac{\delta }{\delta J_j^{\sigma ^Q}}𝐏_{jn}_l\left(\frac{\delta }{\delta J_l^P}\frac{\delta }{\delta J_n^Q}\right)}\right)Z_0[𝐉^P,𝐉^Q,𝐉^{\sigma ^P},𝐉^{\sigma ^Q}],`$ (A18)
where $`Z_0`$ is the functional $`Z`$ evaluated at zero coupling: $`\lambda _0=0`$. The functional derivatives have the property that, e.g.,

$$\frac{\delta J_k(\stackrel{}{x},t)}{\delta J_l(\stackrel{}{x}^{},t^{})}=\delta _{kl}\delta (\stackrel{}{x}\stackrel{}{x}^{})\delta (tt^{}).$$
(A19)
Before continuing, special mention should be made regarding the treatment of the Jacobian functional determinant $`𝒥`$ that appears in $`Z`$. The determinant can be exponentiated by means of the identity $`detM=\mathrm{exp}(tr\mathrm{ln}M)`$, and it can be shown that the argument of the exponential is proportional to the formally divergent object $`\delta ^d(0)`$. This happens to be an ultraviolet or short-distance divergence, and when renormalized by dimensional continuation, yields zero, $`\delta ^d(0)=0`$ (Zinn-Justin 1996). Thus, in this regularization scheme, the determinant renormalizes to yield a finite and field independent constant (unity), and so contributes nothing to the normalized correlation functions calculated from $`Z`$. We henceforth drop these determinants for the remainder of this discussion (by setting $`𝒥1`$). If one wants to investigate the short-distance or ultraviolet limit of turbulent MHD however, then it is preferable to keep these Jacobian determinants intact. For example, this is important when calculating the effective action and effective potential associated with Eq. (A17). The retention of these determinants can be done with ghost fields, as described at length in Hochberg et al. (1999b), where a number of details pertaining to the calculation of Jacobian determinants can be found. To our knowledge, the calculation of $`𝒥`$ for MHD is an open problem.

Up to this stage, we have assumed a completely general noise distribution functional $`𝒫`$. For the remainder of this section, and throughout this paper, only Gaussian noise is considered. For Gaussian noise, $`Z_0`$ can be computed exactly and in closed form. In fact, we consider translationally invariant Gaussian noise, that is

$$𝒫[\stackrel{}{\eta }]=𝒩\mathrm{exp}\left(\frac{1}{2}d^d\stackrel{}{x}𝑑td^d\stackrel{}{x}^{}𝑑t^{}\eta _i(\stackrel{}{x},t)\mathrm{\Gamma }_{ij}^1(\stackrel{}{x}\stackrel{}{x}^{},tt^{})\eta _j(\stackrel{}{x}^{},t^{})\right),$$
(A20)
where $`𝒩`$ is a normalization constant and $`\stackrel{}{𝚪}=\stackrel{T}{\stackrel{}{𝚪}}`$ is the symmetric covariance matrix, which will be specified below. By substituting this form of $`𝒫`$ into $`Z_0`$, and completing the square in the conjugate field variable, which is a $`2d`$-component vector field $`\stackrel{}{\sigma }=(\stackrel{}{\sigma }^P,\stackrel{}{\sigma }^Q)`$, exact integration can be done immediately over the noise fieldGaussian functional integration is treated in any number of excellent texts devoted to (quantum) field theory. We recommend the excellent coverage given in (Ramond 1981; Rivers 1987; Zinn-Justin 1996)..

Following this procedure, completing the square yields

$$\frac{1}{2}\stackrel{}{\eta }\stackrel{1}{\stackrel{}{𝚪}}\stackrel{}{\eta }i\stackrel{}{\sigma }\stackrel{}{\eta }=\frac{1}{2}[\stackrel{}{\eta }+i\stackrel{}{𝚪}\stackrel{}{\sigma }]^T\stackrel{1}{\stackrel{}{𝚪}}[\stackrel{}{\eta }+i\stackrel{}{𝚪}\stackrel{}{\sigma }]\frac{1}{2}\stackrel{}{\sigma }\stackrel{}{𝚪}\stackrel{}{\sigma },$$
(A21)
where $`T`$ denotes transpose. Next, the functional integration over the noise leads to

$$Z_0[𝐉,𝐉^\sigma ]=[𝒟\stackrel{}{\sigma }][𝒟𝚽]\mathrm{exp}(\frac{1}{2}\stackrel{}{\sigma }\stackrel{}{𝚪}\stackrel{}{\sigma }+i\stackrel{}{\sigma }\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽+𝚽𝐉+\stackrel{}{\sigma }𝐉^\sigma ),$$
(A22)
where we have used the fact that the noise distribution function $`𝒫`$ is normalized to unity Eq. (A3). We have streamlined the notation even further by defining two-component fields for the Elsasser fields, their sources and their conjugate fields respectively as $`𝚽=(𝐏,𝐐)`$, $`𝐉=(𝐉^P,𝐉^Q)`$, and $`𝐉^\sigma =(𝐉^{\sigma ^P},𝐉^{\sigma ^Q})`$. The (free) inverse response function in the space and time domain is defined to be

$$[\underset{0}{\overset{1}{\stackrel{}{𝐆}}}]_{ij}=\left(\begin{array}{cc}_t\gamma _+^2& \gamma _{}^2\\ \gamma _{}^2& _t\gamma _+^2\end{array}\right)\delta _{ij}.$$
(A23)
In Eq. (A22) we have a quadratic form in the conjugate field, so we can integrate over this field exactly by completing the square as

$`{\displaystyle \frac{1}{2}}\stackrel{}{\sigma }\stackrel{}{𝚪}\stackrel{}{\sigma }+i\stackrel{}{\sigma }\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽`$ $`+`$ $`\stackrel{}{\sigma }𝐉^\sigma =`$ (A24)
$``$ $`{\displaystyle \frac{1}{2}}(\stackrel{}{\sigma }i\stackrel{1}{\stackrel{}{𝚪}}[\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma ])^T\stackrel{}{𝚪}(\stackrel{}{\sigma }i\stackrel{1}{\stackrel{}{𝚪}}[\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma ])`$ (A25)
$``$ $`{\displaystyle \frac{1}{2}}(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma )^T\stackrel{1}{\stackrel{}{𝚪}}(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma ).`$ (A26)
Integrating over the conjugate fields yields (up to an overall and unimportant multiplicative constant)

$$Z_0[𝐉,𝐉^\sigma ][𝒟𝚽]\mathrm{exp}(\frac{1}{2}(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma )^T\stackrel{1}{\stackrel{}{𝚪}}(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma )+𝚽𝐉).$$
(A27)
Finally, by means of the identity

$`(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}`$ $``$ $`𝚽i𝐉^\sigma )^T\stackrel{}{𝚪}^1(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}𝚽i𝐉^\sigma )=`$ (A29)
$`(𝚽i\underset{0}{\overset{}{𝐆}}𝐉^\sigma )^T(\underset{0}{\overset{1}{\stackrel{}{𝐆}}}\stackrel{1}{\stackrel{}{𝚪}}\underset{0}{\overset{1}{\stackrel{}{𝐆}}})(𝚽i\underset{0}{\overset{}{𝐆}}𝐉^\sigma ),`$
we are led to make the constant shift change-of-variable $`𝚽^{}=𝚽i\underset{0}{\overset{}{𝐆}}𝐉^\sigma `$, integrate over the remaining field $`𝚽^{}`$, and note that $`[𝒟𝚽^{}]=[𝒟𝚽]`$. Thus, up to an overall irrelevant constant prefactor, the free $`(\lambda _0=0)`$ generating functional is given by

$$Z_0[𝐉^P,𝐉^Q,𝐉^{\sigma ^P},𝐉^{\sigma ^Q}]=\mathrm{exp}(i𝐉^\sigma \underset{0}{\overset{}{𝐆}}𝐉)\mathrm{exp}(\frac{1}{2}𝐉[\underset{0}{\overset{}{𝐆}}\stackrel{}{𝚪}\underset{0}{\overset{}{𝐆}}]𝐉).$$
(A30)
This can be expressed in Fourier space as

$`Z_0[𝐉^P,𝐉^Q`$ , $`𝐉^{\sigma ^P},𝐉^{\sigma ^Q}]=\mathrm{exp}(i{\displaystyle }{\displaystyle \frac{d^d\stackrel{}{k}d\omega }{(2\pi )^{d+1}}}𝐉^\sigma (\stackrel{}{k},\omega )\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )𝐉(\stackrel{}{k},\omega ))`$ (A31)
$`\times `$ $`\mathrm{exp}({\displaystyle \frac{1}{2}}{\displaystyle }{\displaystyle \frac{d^d\stackrel{}{k}d\omega }{(2\pi )^{d+1}}}𝐉(\stackrel{}{k},\omega )\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )\stackrel{}{𝚪}(\stackrel{}{k},\omega )\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )𝐉(\stackrel{}{k},\omega )).`$ (A32)
From the free generating functional, it is straightforward to read off the bare response and bare correlation functions (or by taking the appropriate functional derivatives of $`Z_0`$ with respect to the sources). In fact, the whole point in defining $`Z`$ the way we did back in Eq. (A15) was so we could compute the statistical average of products of fields and conjugate fields. As we see, functionally differentiating $`Z`$ with respect to any combination of sources brings down the corresponding factors of fields under the functional integral. Then, setting all the sources to zero after differentiating yields the correlation and response functions associated with this stochastic field theory. Let us now specify the explicit form for the noise or covariance matrix as

$$[\stackrel{}{𝚪}(\stackrel{}{k},\omega )]_{mn}=D_0(k)𝐏_{mn}(\stackrel{}{k}),$$
(A33)
which is the product of a $`2\times 2`$-block noise matrix $`D_0(k)`$ times the transverse projection operator. The presence of this projector enforces the solenoidal character of the noise, and appears as well in the original noise correlation functions above Eqs. (5) and (C17). The noise matrix can be used to model an arbitrary spatially correlated noise. To handle both spatial and temporal correlations in the noise, one would write $`D_0(k,\omega )`$. In this paper, we consider spatially correlated noise as specified by $`2\times 2`$-block $`D_0(k)`$ given below in Eq. (A40). The effects of random helicity (both kinetic and magnetic) can be studied by an appropriate and simple extension of Eq. (A33); this is introduced and analyzed in the main text.

Thus, the bare response function in Fourier space is computed from

$`\sigma _m(\stackrel{}{p}_1,\omega _1)\mathrm{\Phi }_n(\stackrel{}{p}_2,\omega _2)_0`$ $`=`$ $`{\displaystyle \frac{(2\pi )^{2(d+1)}}{Z_0[0]}}{\displaystyle \frac{\delta ^2}{\delta J_m^\sigma (\stackrel{}{p}_1,\omega _1)\delta J_n(\stackrel{}{p}_2,\omega _2)}}|_{J^\sigma =J=0}Z_0[𝐉^\sigma ,𝐉]`$ (A34)
$`=`$ $`i(2\pi )^{d+1}\delta ^d(\stackrel{}{p}_1+\stackrel{}{p}_2)\delta (\omega _1+\omega _2)[\underset{0}{\overset{}{𝐆}}(\stackrel{}{p}_2,\omega _2)]_{mn}.`$ (A35)
The bare response function obtained by calculating the inverse of Eq. (A23) is

$`[\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )]_{mn}`$ $`=`$ $`𝒢_0(\stackrel{}{k},\omega )\delta _{mn},`$ (A36)
$`𝒢_0(\stackrel{}{k},\omega )`$ $`=`$ $`{\displaystyle \frac{1}{[i\omega +k^2(\gamma _++\gamma _{})][i\omega +k^2(\gamma _+\gamma _{})]}}\left(\begin{array}{cc}i\omega +\gamma _+k^2& \gamma _{}k^2\\ \gamma _{}k^2& i\omega +\gamma _+k^2\end{array}\right),`$ (A37)
where the convention is that momentum and frequency are taken to flow from the conjugate field to the physical field (otherwise, one must flip the sign of $`\omega `$). Next, the bare correlation function (in Fourier space) is computed from

$`\mathrm{\Phi }_m(\stackrel{}{p}_1,\omega _1)\mathrm{\Phi }_n(\stackrel{}{p}_2,\omega _2)_0`$ $`=`$ $`{\displaystyle \frac{(2\pi )^{2(d+1)}}{Z_0[0]}}{\displaystyle \frac{\delta ^2}{\delta J_m(\stackrel{}{p}_1,\omega _1)\delta J_n(\stackrel{}{p}_2,\omega _2)}}|_{J^\sigma =J=0}Z_0[𝐉^\sigma ,𝐉],`$ (A38)
$`=`$ $`(2\pi )^{d+1}\delta ^d(\stackrel{}{p}_1+\stackrel{}{p}_2)\delta (\omega _1+\omega _2)[\underset{0}{\overset{}{𝐂}}(\stackrel{}{p}_2,\omega _2)]_{mn},`$ (A39)
where,

$`[\underset{0}{\overset{}{𝐂}}(\stackrel{}{k},\omega )]_{mn}`$ $`=`$ $`𝒞_0(\stackrel{}{k},\omega )𝐏_{mn}(\stackrel{}{k}),`$ (A40)
$`𝒞_0(\stackrel{}{k},\omega )`$ $`=`$ $`𝒢_0(\stackrel{}{k},\omega )D_0(k)𝒢_0(\stackrel{}{k},\omega ),`$ (A41)
$`D_0(k)`$ $`=`$ $`\left(\begin{array}{cc}2k^yA& 2k^yB\\ 2k^yB& 2k^yA\end{array}\right).`$ (A42)
To complete the specification of the full generating functional $`Z`$, we must also express the non-linear interaction vertex in Fourier space. This is obtained by Fourier-transforming the cubic terms appearing in Eq. (A11) or in Eq. (A17). Note there are two kinds of cubic interaction, one corresponding to the $`\sigma ^P𝐏𝐐`$ vertex and the other corresponding to the $`\sigma ^Q𝐏𝐐`$ vertex (the distinguishing feature is the conjugate field). Of course, there is only one (bare) cubic coupling parameter $`\lambda _0`$ (which we set to unity at the end of the calculation). This means corrections only to one of the two vertices needs to be considered, even though both bare vertices are needed in building up loop-expansions. Taking the convention that all momentum and frequency flow inward to the vertex, the bare vertex corresponding to the cubic interaction $`\sigma _j^P(\stackrel{}{k}_1,\omega _1)P_n(\stackrel{}{k}_2,\omega _2)Q_l(\stackrel{}{k}_3,\omega _3)`$ as well as the other cubic interaction $`\sigma _j^Q(\stackrel{}{k}_1,\omega _1)Q_n(\stackrel{}{k}_2,\omega _2)P_l(\stackrel{}{k}_3,\omega _3)`$ is

$$i\lambda _0𝐏_{l,jn}(\stackrel{}{k}_1)=i\lambda _0(\stackrel{}{k}_1)_l𝐏_{jn}(\stackrel{}{k}_1)=i\lambda _0(\stackrel{}{k}_1)_l\left(\delta _{jn}\frac{(\stackrel{}{k}_1)_j(\stackrel{}{k}_1)_n}{(\stackrel{}{k}_1)^2}\right).$$
(A43)
The vertex is independent of the external frequencies. Also, here it has been written so that it depends only on the momentum carried by the conjugate field(s), since momentum conservation implies $`\stackrel{}{k}_1+\stackrel{}{k}_2+\stackrel{}{k}_3=0`$.

In summary, in this Appendix the stochastic PDEs of MHD Eqs. (20) and (21) have been expressed in terms of a functional integral representation. From this representation, the bare response function Eq. (A36), correlation functions Eq. (A40), and the bare non-linear interaction vertex Eq. (A43) have been identified and calculated. As explained earlier, this set of functions constitute the basic elements (see Fig. 1) upon which a diagrammatic expansion (see Fig. 2) of the generating functional can be systematically developed.

## B Response function renormalization

This Appendix carries out the one-loop renormalization of the response function. Transcribing the diagrammatic equation in Fig. (2a) for the one-loop corrected response function into the corresponding analytic expression by means of the Feynman rules in Fig. 1 yields

$$[\stackrel{}{𝐆}(\stackrel{}{k},\omega )]_{mn}=[\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )]_{mn}𝒮\lambda _0^2[\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )]_{m\overline{m}}[𝐏_{l,\overline{n}s}(\stackrel{}{k})+𝐏_{s,\overline{n}l}(\stackrel{}{k})]I_{\overline{m}ls}(\stackrel{}{k},\omega )[\underset{0}{\overset{}{𝐆}}(\stackrel{}{k},\omega )]_{\overline{n}n},$$
(B1)
where the symmetry factor for the one loop diagram $`𝒮=4`$ and the loop-integral is given by

$$I_{\overline{m}ls}(\stackrel{}{k},\omega )=^>\frac{d^d\stackrel{}{q}}{(2\pi )^d}_{\mathrm{}}^{\mathrm{}}\frac{d\mathrm{\Omega }}{2\pi }[𝐏_{\overline{m},lr}(\stackrel{}{k}\stackrel{}{q})+𝐏_{r,l\overline{m}}(\stackrel{}{k}\stackrel{}{q})]𝒢_0(\stackrel{}{k}\stackrel{}{q},\omega \mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝐏_{rs}(\stackrel{}{q}).$$
(B2)
The precise meaning of the momentum-shell integration will be given below. Using the factorization of the response function into the product of a block matrix times a Kronecker delta as given in Eqs. (A36) and (B1), the following equation is obtained for the one-loop correction to the inverse response function,

$$𝒢^1(\stackrel{}{k},\omega )\delta _{mn}=𝒢_0^1(\stackrel{}{k},\omega )\delta _{mn}+𝒮\lambda _0^2[𝐏_{l,ns}(\stackrel{}{k})+𝐏_{s,nl}(\stackrel{}{k})]I_{mls}(\stackrel{}{k},\omega )+O(\lambda _0^4),$$
(B3)
where the renormalized $`2\times 2`$ inverse block is parametrized as

$$𝒢^1(\stackrel{}{k},\omega )=\left(\begin{array}{cc}i\omega +\gamma _+^<k^2& \gamma _{}^<k^2\\ \gamma _{}^<k^2& i\omega +\gamma _+^<k^2\end{array}\right).$$
(B4)
This structure follows from Fourier transforming Eq. (A23). The infrared (i.e., large-distance and long-time) renormalizability of the response function implies that the renormalized or corrected response function Eq. (B4) must have the identical mathematical structure as its bare counterpart $`𝒢_0^1`$. It therefore has the same frequency and momentum dependence as the bare response function and contains the same number of parameters. Here $`\gamma _\pm ^<=\gamma _\pm ^<(A,B,\gamma _+,\gamma _{},\lambda _0)`$ are computable functions of the indicated bare parameters which will be determined below.

Thus the renormalization of the response function yields the individual renormalization of the two independent viscosities $`\gamma _+,\gamma _{}`$. The renormalization of the two viscosities requires some amount of careful calculation. We will go through the major points step by step. To proceed we need to explicitly carry out the internal frequency integration that appears in the loop integral Eq. (B2). From Eq. (B4), we see that the viscosities $`\gamma _\pm `$ are the coefficients of the $`k^2`$ terms, so we can set the external frequency $`\omega `$ to zero from the outset. This frequency integral is straightforwardly calculated by means of the calculus of residues. We close the contour in the lower half plane, where there are two simple poles, and obtain

$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}𝒢_0(`$ $`\stackrel{}{k}`$ $`\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })=`$ (B6)
$`{\displaystyle \frac{D_0(q)}{8q^4}}\left(1+{\displaystyle \frac{\stackrel{}{q}\stackrel{}{k}}{q^2}}+O(k^2)\right)\left\{{\displaystyle \frac{1}{(\gamma _+\gamma _{})^2}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{1}{(\gamma _++\gamma _{})^2}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\},`$
$`=`$ $`{\displaystyle \frac{1}{4}}q^{y4}\left(1+{\displaystyle \frac{\stackrel{}{q}\stackrel{}{k}}{q^2}}+O(k^2)\right)\left\{{\displaystyle \frac{AB}{(\gamma _+\gamma _{})^2}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{A+B}{(\gamma _++\gamma _{})^2}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}.`$ (B7)
The result displayed in the first line holds for all spatially correlated Gaussian noise $`D_0(q)`$. We should point out this was evaluated assuming a positive viscosity, $`\nu >0`$, and a positive resistivity, $`\nu _B>0`$. A change in sign will change the location of the poles and correspondingly modify the frequency integration. The second line results from substituting the particular form for the noise spectrum used in this paper, Eq. (A40). If one wants to treat noise with temporal correlations, i.e., noise spectra of the form $`D_0(q,\mathrm{\Omega })`$, then the frequency integral must be evaluated in a case-by-case fashion, and will generally involve both poles and branch cuts.

Substituting this result back into Eq. (B2), we turn next to the remaining and final integration over the internal loop momentum. From Eq. (B4) the renormalization of the parameters $`\gamma _\pm `$ requires that we expand the one-loop contribution up to second order in the external momentum. The projection operator in Eq. (B3) is linear in $`\stackrel{}{k}`$, so the loop integral only need be expanded to first order in external momentum. It is important to note that this integral depends on the external momentum not only via the integrand but also through its domain of integration. This is because, when integrating over the momentum shell, we must ensure that all combinations of internal and external momenta appearing in Eq. (B2) remain within the band or shell of momenta to be integrated over. In the case of the present integral, this means we must integrate over the intersection of the domains $`\mathrm{\Lambda }/s|\stackrel{}{q}|\mathrm{\Lambda }`$ and $`\mathrm{\Lambda }/s|\stackrel{}{k}\stackrel{}{q}|\mathrm{\Lambda }`$, where $`s>1`$ is a measure of the fraction of modes to be eliminated in the coarse-graining. To first order in $`\stackrel{}{k}`$, the second inequality can be written as $`\mathrm{\Lambda }/s+k\mathrm{cos}\theta <q<\mathrm{\Lambda }+k\mathrm{cos}\theta `$, where $`\theta `$ is the angle between $`\stackrel{}{k}`$ and $`\stackrel{}{q}`$. There are two cases to consider: when $`\mathrm{cos}\theta >0`$, the intersection of the two intervals can be expressed as the difference $`[\mathrm{\Lambda }/s,\mathrm{\Lambda }][\mathrm{\Lambda }/s,\mathrm{\Lambda }/s+k\mathrm{cos}\theta ]`$, and when $`\mathrm{cos}\theta <0`$, the intersection can be written as $`[\mathrm{\Lambda }/s,\mathrm{\Lambda }][\mathrm{\Lambda }+k\mathrm{cos}\theta ,\mathrm{\Lambda }]`$. This means the complete momentum-shell integration valid up to $`O(k^2)`$ can be written as

$$^>\frac{d^d\stackrel{}{q}}{(2\pi )^d}=𝑑\mathrm{\Omega }_d\left(_{\mathrm{\Lambda }/s}^\mathrm{\Lambda }_{\mathrm{\Lambda }/s}^{\mathrm{\Lambda }/s+k\mathrm{cos}\theta }_{\mathrm{\Lambda }+k\mathrm{cos}\theta }^\mathrm{\Lambda }\right)\frac{dqq^{d1}}{(2\pi )^d}+O(k^2).$$
(B8)
Furthermore, note it is convenient to carry out the integrations over the momentum modulus before working out the angular integrations over the element of solid angle $`d\mathrm{\Omega }_d`$.

From Eqs. (B6) and (B2), the remaining integration to be carried out after the frequency integration is performed is

$`I`$ $`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}[𝐏_{m,lr}(\stackrel{}{k}\stackrel{}{q})+𝐏_{r,lm}(\stackrel{}{k}\stackrel{}{q})]𝐏_{rs}(\stackrel{}{q})q^{y4}\left(1+{\displaystyle \frac{\stackrel{}{k}\stackrel{}{q}}{q^2}}+O(k^2)\right),`$ (B9)
$`=`$ $`I_1+I_2+I_3,`$ (B10)
where we break the full integral $`I`$ down into somewhat smaller and more manageable pieces as

$`I_1`$ $`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}q^{y4}\left(q_m𝐏_{ls}(\stackrel{}{q})\right)+O(k^2),`$ (B11)
$`I_2`$ $`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}q^{y4}\left(q_m𝐏_{ls}(\stackrel{}{q})\right){\displaystyle \frac{\stackrel{}{k}\stackrel{}{q}}{q^2}}+O(k^2),`$ (B12)
$`I_3`$ $`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}q^{y4}\left(k_m𝐏_{ls}(\stackrel{}{q})(k_s\stackrel{}{k}\stackrel{}{q}{\displaystyle \frac{q_s}{q^2}}){\displaystyle \frac{q_mq_l}{q^2}}+k_r𝐏_{lm}(\stackrel{}{q})𝐏_{rs}(\stackrel{}{q})\right)+O(k^2).`$ (B13)
The following projection operator product expansions proved to be useful in arriving at these expressions:

$`𝐏_{a,bc}(\stackrel{}{k}\stackrel{}{q})𝐏_{dc}(\stackrel{}{q})`$ $`=`$ $`(k_aq_a)\left(\delta _{bd}{\displaystyle \frac{q_bq_d}{q^2}}\right){\displaystyle \frac{q_aq_b}{q^2}}\left(k_d{\displaystyle \frac{\stackrel{}{q}\stackrel{}{k}q_d}{q^2}}\right)+O(k^2),`$ (B14)
$`𝐏_{a,bc}(\stackrel{}{k}\stackrel{}{q})𝐏_{ad}(\stackrel{}{q})`$ $`=`$ $`k_a𝐏_{bc}(\stackrel{}{q})𝐏_{ad}(\stackrel{}{q})+O(k^2).`$ (B15)
The integrand of $`I_1`$ is independent of $`k`$ and the contribution from the interval $`[\mathrm{\Lambda }/s,\mathrm{\Lambda }]`$ vanishes upon angular averaging. Only the latter two contributions in Eq. (B8) will contribute to the angular averaging and we find that

$$I_1=[\left(\frac{\mathrm{\Lambda }}{s}\right)^{dy4}\mathrm{\Lambda }^{dy4}]\left\{\frac{1}{2}\frac{S_d}{d(2\pi )^d}k_m\delta _{ls}\frac{1}{2}\frac{S_d}{d(d+2)(2\pi )^d}(k_m\delta _{ls}+k_l\delta _{ms}+k_s\delta _{ml})\right\},$$
(B16)
where $`S_d`$ is defined in Eq. (E2). The integrand in $`I_2`$ is already linear in $`k`$, so it suffices to work out the contribution from the interval $`[\mathrm{\Lambda }/s,\mathrm{\Lambda }]`$. This yields

$$I_2=\frac{2}{(dy4)}I_1.$$
(B17)
Lastly, the integrand in $`I_3`$ is also linear in $`k`$. Integrating over $`[\mathrm{\Lambda }/s,\mathrm{\Lambda }]`$ and averaging over the angles yields

$`I_3={\displaystyle \frac{[\mathrm{\Lambda }^{dy4}(\frac{\mathrm{\Lambda }}{s})^{dy4}]}{(dy4)}}({\displaystyle \frac{(d1)S_d}{d(2\pi )^d}}k_m\delta _{ls}`$ $`+`$ $`{\displaystyle \frac{(d3)S_d}{d(2\pi )^d}}k_s\delta _{lm}`$ (B18)
$`+`$ $`{\displaystyle \frac{2S_d}{d(d+2)(2\pi )^d}}(k_m\delta _{ls}+k_l\delta _{ms}+k_s\delta _{ml})).`$ (B19)
Various identities needed for the angular integrations are collected in Eq. (E2). This completes the calculation of the one-loop integral. Combining all the results from Eqs. (B9), (B11), (B12), (B13), (B16) ,(B17), and (B18), we obtain the one-loop matrix equation

$`𝒢^1(\stackrel{}{k},0)\delta _{mn}=𝒢_0^1(\stackrel{}{k},0)\delta _{mn}`$ $`+`$ $`𝒮{\displaystyle \frac{\lambda _0^2}{4}}{\displaystyle \frac{(d^2y4)S_d}{d(d+2)(dy4)(2\pi )^d}}[\mathrm{\Lambda }^{dy4}\left({\displaystyle \frac{\mathrm{\Lambda }}{s}}\right)^{dy4}]`$ (B21)
$`\left\{{\displaystyle \frac{AB}{(\gamma _+\gamma _{})^2}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{A+B}{(\gamma _++\gamma _{})^2}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}\times k^2𝐏_{mn}(\stackrel{}{k}).`$
We can project out individual equations for each of the two viscosities by multiplying Eq. (B21) through by $`𝐏_{ml}(\stackrel{}{k})`$ and using the idempotency of this projection operator, $`𝐏_{mn}(\stackrel{}{k})𝐏_{ml}(\stackrel{}{k})=𝐏_{nl}(\stackrel{}{k})`$. Thus we find that each viscosity renormalizes as

$`\gamma _+^<`$ $`=`$ $`\gamma _++𝒮{\displaystyle \frac{\lambda _0^2}{4}}[\mathrm{\Lambda }^{dy4}({\displaystyle \frac{\mathrm{\Lambda }}{s}})^{dy4}]{\displaystyle \frac{(d^2y4)S_d}{d(d+2)(dy4)(2\pi )^d}}\left\{{\displaystyle \frac{(AB)}{(\gamma _+\gamma _{})^2}}+{\displaystyle \frac{(A+B)}{(\gamma _++\gamma _{})^2}}\right\},`$ (B22)
$`\gamma _{}^<`$ $`=`$ $`\gamma _{}+𝒮{\displaystyle \frac{\lambda _0^2}{4}}[\mathrm{\Lambda }^{dy4}({\displaystyle \frac{\mathrm{\Lambda }}{s}})^{dy4}]{\displaystyle \frac{(d^2y4)S_d}{d(d+2)(dy4)(2\pi )^d}}\left\{{\displaystyle \frac{(AB)}{(\gamma _+\gamma _{})^2}}+{\displaystyle \frac{(A+B)}{(\gamma _++\gamma _{})^2}}\right\}.`$ (B23)
Individual equations for the sum $`(\gamma _++\gamma _{})`$ and difference $`(\gamma _+\gamma _{})`$ can be easily obtained by taking the sum and difference, respectively, of Eqs. (B22) and (B23). This yields the individual renormalizations of the fluid viscosity $`\nu `$ and magnetic resistivity $`\nu _B`$.

## C Noise spectrum renormalization

From the diagrammatic representation of the one-loop correction to the noise spectral function in Fig. (2b), (the noise spectral function is obtained by amputating the external legs from the correlation function) which is a $`2\times 2`$ block matrix, we obtain the following expression for the one-loop noise spectrum in terms of the bare response, correlation and vertex functions,

$$D(k)𝐏_{rs}(\stackrel{}{k})=D_0(k)𝐏_{rs}(\stackrel{}{k})+𝒮\lambda _0^2[𝐏_{l,rn}(\stackrel{}{k})+𝐏_{n,rl}(\stackrel{}{k})]I_{lm,nj}(\stackrel{}{k},\omega )[𝐏_{m,sj}(\stackrel{}{k})+𝐏_{j,sm}(\stackrel{}{k})],$$
(C1)
where the symmetry factor in this case is $`𝒮=2`$ and the associated one-loop integral is given by

$`I_{lm,nj}(\stackrel{}{k},\omega )`$ $`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}[\underset{0}{\overset{}{𝐂}}(\stackrel{}{k}\stackrel{}{q},\omega \mathrm{\Omega })]_{lm}[\underset{0}{\overset{}{𝐂}}(\stackrel{}{q},\mathrm{\Omega })]_{nj}`$ (C2)
$`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}𝐏_{lm}(\stackrel{}{k}\stackrel{}{q})𝐏_{nj}(\stackrel{}{q})𝒞_0(\stackrel{}{k}\stackrel{}{q},\omega \mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega }),`$ (C3)
$`𝒞_0(\stackrel{}{k}\stackrel{}{q},\omega \mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })`$ $`=`$ $`𝒢_0(\stackrel{}{k}\stackrel{}{q},\omega \mathrm{\Omega })D_0(\stackrel{}{k}\stackrel{}{q})𝒢_0(\stackrel{}{q}\stackrel{}{k},\mathrm{\Omega }\omega )`$ (C4)
$`\times `$ $`𝒢_0(\stackrel{}{q},\mathrm{\Omega })D_0(q)𝒢_0(\stackrel{}{q},\mathrm{\Omega }).`$ (C5)
The infrared renormalizability of the noise spectral function implies that the renormalized spectral function must have the identical mathematical structure as its bare counterpart $`D_0(k)`$ in Eq. (A40). It must therefore have the same frequency and momentum dependence as the bare noise spectral function and contain the same number of parameters,

$$D(k)=\left(\begin{array}{cc}2k^yA^<& 2k^yB^<\\ 2k^yB^<& 2k^yA^<\end{array}\right).$$
(C6)
Here $`A^<=A^<(A,B,\gamma _+,\gamma _{},\lambda _0)`$ and $`B^<=B^<(A,B,\gamma _+,\gamma _{},\lambda _0)`$ are calculable in terms of the bare parameters and will be determined below.

We are interested in evaluating this loop integral in the hydrodynamic limit (i.e., for small $`k`$ and small $`\omega `$). However, there is already an overall quadratic factor of $`k^2`$ multiplying the loop integral arising from the product of individual projection operator factors appearing in Eq. (C1). As can be seen from Taylor expanding about $`k=0`$, the loop integral itself can only contribute a constant plus positive integer powers of external momentum $`\stackrel{}{k}`$ in the hydrodynamic limit. Thus, the complete one-loop correction in the long-wavelength regime is a power series in $`k`$, and it starts off at quadratic order. For the type of noise of interest to us, we now show that this one loop correction is irrelevant in the hydrodynamic regime. In fact, by setting $`\stackrel{}{k}`$ and $`\omega `$ to zero inside the loop-integral Eq. (C2), it gives

$`I_{lm,nj}(0,0)`$ $`=`$ $`{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}𝐏_{lm}(\stackrel{}{q})𝐏_{nj}(\stackrel{}{q})𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })`$ (C7)
$`=`$ $`{\displaystyle 𝑑\mathrm{\Omega }_d𝐏_{lm}(\stackrel{}{q})𝐏_{nj}(\stackrel{}{q})_{\mathrm{\Lambda }/s}^\mathrm{\Lambda }\frac{q^{d1}dq}{(2\pi )^d}_{\mathrm{}}^{\mathrm{}}\frac{d\mathrm{\Omega }}{2\pi }𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })}`$ (C8)
$`=`$ $`{\displaystyle 𝑑\mathrm{\Omega }_d𝐏_{lm}(\stackrel{}{q})𝐏_{nj}(\stackrel{}{q})_{\mathrm{\Lambda }/s}^\mathrm{\Lambda }\frac{q^{d1}dq}{(2\pi )^d}(q)}.`$ (C9)
Here

$$(q)=\left(\begin{array}{cc}M_1(q)& M_2(q)\\ M_2(q)& M_1(q)\end{array}\right),$$
(C10)
is the $`2\times 2`$ block-matrix function of the modulus $`q=|\stackrel{}{q}|`$, which results from integrating the product of the bare correlation functions in Eq. (C2) over the internal loop frequency. This can be straightforwardly calculated by means of residues (there are four simple poles in both the upper and lower half-planes, respectively), but we need not do so for the purposes of the present discussion. Note that the integration over the sphere factors out from the modulus and frequency integrations. This is because the product of projection operators depends only on angles, while the product of the correlation block-matrices depends only on frequency and momentum modulus. Integrating the product of the projection operators that appear in Eq. (C7) over the unit sphere using the identities in Eq. (E2) yields the result

$$𝑑\mathrm{\Omega }_d𝐏_{lm}(\stackrel{}{q})𝐏_{nj}(\stackrel{}{q})=\frac{S_d}{d(d+2)}\left((d^23)\delta _{lm}\delta _{nj}+[\delta _{ln}\delta _{mj}+\delta _{lj}\delta _{nm}]\right).$$
(C11)
After using the properties of the projection operators to further simplify the final expression, it leads to

$`[𝐏_{l,rn}(\stackrel{}{k})+𝐏_{n,rl}(\stackrel{}{k})]I_{lm,nj}(0,0)[𝐏_{m,sj}(\stackrel{}{k})`$ $`+`$ $`𝐏_{j,sm}(\stackrel{}{k})]=`$ (C13)
$`2k^2\times 𝐏_{rs}(\stackrel{}{k})\left({\displaystyle \frac{S_d(d^22)}{d(d+2)}}\right){\displaystyle _{\mathrm{\Lambda }/s}^\mathrm{\Lambda }}{\displaystyle \frac{q^{d1}dq}{(2\pi )^d}}(q).`$
The one-loop correction is of order $`O(k^2)`$ and is proportional to the same projector as the bare noise spectrum. Inserting the above line into Eq. (C1) gives

$$D(k)=D_0(k)+𝒮\lambda _0^2\mathrm{\hspace{0.17em}2}k^2\times \left(\frac{S_d(d^22)}{d(d+2)(2\pi )^d}\right)_{\mathrm{\Lambda }/s}^\mathrm{\Lambda }q^{d1}𝑑q(q)+O(k^3).$$
(C14)
Referring to the form of the noise spectrum matrix in Eq. (A40), we see that the one-loop correction to the noise spectrum is irrelevant in the long-wavelength limit provided $`y>2`$ (The case $`y=2`$ is akin to Model A of Forster, Nelson and Stephen 1977). The Taylor expansion of $`I_{lm,nj}(\stackrel{}{k},0)`$ about $`\stackrel{}{k}=0`$ (taking into account contributions from the domain of integration as well, see discussion surrounding Eq. (B8) ) will just introduce positive powers of the external momentum. Thus we conclude that at one-loop as $`k0`$ and for $`y>2`$

$`D(k)`$ $`=`$ $`D_0(k)`$ (C15)
$`A^<`$ $`=`$ $`A,`$ (C16)
$`\mathrm{and}B^<`$ $`=`$ $`B.`$ (C17)
Note that for noise spectra having $`y<2`$, the one-loop correction will induce relevant operators in the hydrodynamic limit. For example, if $`y=3`$, the one-loop correction goes as $`O(k^2)`$ and in fact dominates over the the initial bare spectrum in the infrared. To be able to renormalize consistently in these other situations requires adding the appropriate dominant relevant operators to the noise spectrum. For the case $`y=2`$, one needs only complete the calculation sketched in Eq. (C7) and insert this into Eq. (C14). In this case, one would find that both $`A`$ and $`B`$ receive nonzero corrections in the hydrodynamic limit.

## D Vertex renormalization

The diagrammatic expansion for the one-loop correction to the vertex function, or coupling constant, is depicted in Fig. (2c). As can be seen from this figure, there are three distinct one-loop diagrams contributing to the correction. Momentum and frequency conservation implies that the vertex function can depend at most on two independent external momenta and two independent external frequencies. With the momentum and frequency flow assignment taken as shown in Fig. (2c), where momentum and frequency enter via the conjugate field and exit via the two physical fields, the equation for the one-loop vertex correction is

$$\lambda ^<𝐏_{l,nm}(\stackrel{}{k})=\lambda _0𝐏_{l,nm}(\stackrel{}{k})𝒮\lambda _0^3[𝐏_{j,nt}(\stackrel{}{k})+𝐏_{t,nj}(\stackrel{}{k})]I_{lmjt}(\stackrel{}{k},\stackrel{}{k}_2,\stackrel{}{k}\stackrel{}{k}_2;\omega ,\omega _2,\omega \omega _2),$$
(D1)
where the diagrammatic symmetry factor is $`𝒮=4`$. Here $`\lambda ^<=\lambda ^<(A,B,\gamma _+,\gamma _{},\lambda _0)`$ is the renormalized coupling parameter, which in principle can depend on all the bare parameters appearing in the equations of motion. The loop-integral representing the sum of the three types of (amputated) triangle graphs can be written as

$`I_{lmjt}`$ $`(`$ $`\stackrel{}{k},\stackrel{}{k}_2,\stackrel{}{k}\stackrel{}{k}_2;\omega ,\omega _2,\omega \omega _2)=`$ (D3)
$`\mathrm{tr}{\displaystyle ^>}{\displaystyle \frac{d^d\stackrel{}{q}}{(2\pi )^d}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}\{𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝐏_{jk}(\stackrel{}{q})[𝐏_{k,rl}(\stackrel{}{k}_2\stackrel{}{q})+𝐏_{l,rk}(\stackrel{}{k}_2\stackrel{}{q})]𝒢_0(\stackrel{}{q}\stackrel{}{k}_2,\mathrm{\Omega }\omega _2)\delta _{rs}`$
$`\times `$ $`[𝐏_{s,pm}(\stackrel{}{k}\stackrel{}{q})+𝐏_{m,ps}(\stackrel{}{k}\stackrel{}{q})]𝒢_0(\stackrel{}{q}\stackrel{}{k},\mathrm{\Omega }\omega )\delta _{tp}`$ (D4)
$`+`$ $`𝒢_0(\stackrel{}{q},\mathrm{\Omega })\delta _{jk}[𝐏_{l,kr}(\stackrel{}{q})+𝐏_{r,kl}(\stackrel{}{q})]𝒞_0(\stackrel{}{q}\stackrel{}{k}_2,\mathrm{\Omega }\omega _2)𝐏_{rs}(\stackrel{}{q}\stackrel{}{k}_2)`$ (D5)
$`\times `$ $`[𝐏_{s,pm}(\stackrel{}{k}\stackrel{}{q})+𝐏_{m,ps}(\stackrel{}{k}\stackrel{}{q})]𝒢_0(\stackrel{}{q}\stackrel{}{k},\mathrm{\Omega }\omega )\delta _{tp}`$ (D6)
$`+`$ $`𝒢_0(\stackrel{}{q},\mathrm{\Omega })\delta _{jk}[𝐏_{l,kr}(\stackrel{}{q})+𝐏_{r,kl}(\stackrel{}{q})]𝒢_0(\stackrel{}{k}_2\stackrel{}{q},\omega _2\mathrm{\Omega })\delta _{rs}[𝐏_{p,sm}(\stackrel{}{q}\stackrel{}{k}_2)+𝐏_{m,sp}(\stackrel{}{q}\stackrel{}{k}_2)]`$ (D7)
$`\times `$ $`𝒞_0(\stackrel{}{k}\stackrel{}{q},\omega \mathrm{\Omega })𝐏_{tp}(\stackrel{}{k}\stackrel{}{q})\},`$ (D8)
which follows from translating the diagrams into their associated mathematical expressions using the basic elements derived above for the bare response Eq. (A36), bare correlation Eq. (A40) and bare vertex functions Eq. (A43). The trace (tr) is taken over the product of the $`2\times 2`$ block matrices. We are interested in evaluating this integral in the hydrodynamic limit (i.e., for small $`k`$ and small $`\omega `$). As there is already a factor linear in $`\stackrel{}{k}`$ multiplying the loop integral, and the bare vertex is itself linear in momentum, we can immediately take the limit $`\stackrel{}{k}=\stackrel{}{k}_20`$ and $`\omega =\omega _20`$ inside the loop integral right from the outset and focus attention on $`I_{lmjt}(0,0,0;0,0,0)`$. This results in a tremendous simplification. The properties of the projection operators allow one to further reduce products such as

$$[𝐏_{l,jr}(\stackrel{}{q})+𝐏_{r,jl}(\stackrel{}{q})][𝐏_{p,rm}(\stackrel{}{q})+𝐏_{m,rp}(\stackrel{}{q})]𝐏_{tp}(\stackrel{}{q})=q_lq_m𝐏_{jt}(\stackrel{}{q})=q^2\widehat{n}_l\widehat{n}_m𝐏_{jt}(\stackrel{}{q}),$$
(D9)
and likewise for the other products that one encounters in Eq. (D3).

The angular integrations immediately can be computed by using

$$𝑑\mathrm{\Omega }_d\widehat{n}_l\widehat{n}_m𝐏_{jt}(\stackrel{}{q})=\frac{S_d}{d(d+2)}\left((d+1)\delta _{jt}\delta _{lm}[\delta _{lj}\delta _{mt}+\delta _{lt}\delta _{jm}]\right),$$
(D10)
with the result that

$`\lambda ^<𝐏_{l,nm}(\stackrel{}{k})`$ $`=`$ $`\lambda _0𝐏_{l,nm}(\stackrel{}{k})+\lambda _0^3\left({\displaystyle \frac{2S_d}{d(d+2)}}\right)[𝐏_{l,nm}(\stackrel{}{k})+𝐏_{m,nl}(\stackrel{}{k})]\times `$ (D11)
$`\mathrm{tr}`$ $`{\displaystyle ^>}{\displaystyle \frac{q^{d+1}dq}{(2\pi )^d}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}\{𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })`$ (D12)
$`+`$ $`𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })\}.`$ (D13)
The frequency integral can be performed by the method of residues. The integrand has both simple and double poles in the complex frequency plane; the integration contours can be closed in either the upper or lower half-plane. The first and third factors appearing under the frequency integral have simple poles in the upper and lower half-plane, respectively. These separate contributions can be neatly combined to yield

$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}(𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })`$ $`+`$ $`𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega }))`$ (D14)
$`=`$ $`{\displaystyle \frac{D_0(q)}{8q^6}}\left\{{\displaystyle \frac{1}{(\gamma _+\gamma _{})^3}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{1}{(\gamma _++\gamma _{})^3}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}.`$ (D15)
On the other hand, the middle factor has double poles in both the upper and lower half-planes and its contribution to the frequency integral turns out to yield

$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}𝒢_0(\stackrel{}{q},\mathrm{\Omega })𝒞_0(\stackrel{}{q},\mathrm{\Omega })𝒢_0(\stackrel{}{q},\mathrm{\Omega })={\displaystyle \frac{D_0(q)}{8q^6}}\left\{{\displaystyle \frac{1}{(\gamma _+\gamma _{})^3}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]+{\displaystyle \frac{1}{(\gamma _++\gamma _{})^3}}\left[\begin{array}{cc}1& 1\\ 1& 1\end{array}\right]\right\}.`$ (D16)
Subtracting Eq. (D16) from Eq. (D14) as per Eq. (D11) yields a null result for the frequency integral that holds for all spatially correlated Gaussian noise $`D_0(q)`$. Thus the vertex coupling parameter does not renormalize in the long-wavelength limit and we have that

$$\lambda ^<=\lambda _0.$$
(D17)
This non-renormalization of the nonlinear coupling parameter holds also for the Navier-Stokes equation driven by arbitrary spatially correlated Gaussian noise and is ascribed to the Galilean invariance of the stochastic NS equation (Forster et al. 1977). The importance of Galilean invariance for the Burgers equation (and the KPZ) equation is further addressed in Medina et al. (1989). We suspect that the non-renormalization of the nonlinearity holds to all orders in the loop expansion for stochastic MHD. To establish this one could investigate the renormalization of the vertex coupling via suitable Ward identities, as was done for the case of the KPZ equation (Frey and Tauber 1994). We leave that problem for a separate discussion elsewhere. Temporally correlated noise breaks the Galilean invariance and can ( and usually does) lead to non-trivial renormalization of $`\lambda _0`$.

## E Integrating over the d-sphere

The various angular integrations and angular averages that are needed in the renormalization group calculations are collected here. The area element of the unit sphere in $`d`$-dimensions is

$$d\mathrm{\Omega }_d=d\varphi \mathrm{sin}\theta _1d\theta _1\mathrm{sin}^2\theta _2d\theta _2\mathrm{}\mathrm{sin}^{(d2)}\theta _{d2}d\theta _{d2},$$
(E1)
where $`0\varphi <2\pi `$ and $`0\theta _j<\pi `$, for $`j=1,2,\mathrm{}(d2)`$. Let $`S_d`$ represent the area of the unit $`d`$-sphere, $`\widehat{n}_j`$ denote a unit vector in the $`j`$-th direction, and $`\mathrm{\Gamma }`$ denote the Gamma function. Then

$`{\displaystyle 𝑑\mathrm{\Omega }_d}`$ $`=`$ $`{\displaystyle \frac{2\pi ^{d/2}}{\mathrm{\Gamma }(\frac{d}{2})}}=S_d,`$ (E2)
$`{\displaystyle 𝑑\mathrm{\Omega }_d\widehat{n}_i\widehat{n}_j}`$ $`=`$ $`{\displaystyle \frac{S_d}{d}}\delta _{ij},`$ (E3)
$`{\displaystyle 𝑑\mathrm{\Omega }_d\widehat{n}_i\widehat{n}_j\widehat{n}_n\widehat{n}_m}`$ $`=`$ $`{\displaystyle \frac{S_d}{d(d+2)}}\left(\delta _{ij}\delta _{mn}+\delta _{im}\delta _{jn}+\delta _{in}\delta _{jm}\right).`$ (E4)
The angular average of the product of an odd number of unit vectors over the unit sphere vanishes identically.

Amit, D.J. 1978 Field Theory, the Renormalization Group, and Critical Phenomena. McGraw-Hill.

Barabási A.-L., and Stanley, H.E. 1995 Fractal Concepts in Surface Growth. Cambridge.

Biskamp, D. 1993 Nonlinear Magnetohydrodynamics. Cambridge.

Brandenburg, A., Enqvist, K. and Olesen, P. 1996 Large-scale magnetic fields from hydromagnetic turbulence in the very early universe. Phys. Rev. D54, 1291-1300.

Burlaga, L.F. 1991 Intermittent turbulence in the solar-wind. J. Geophys. Res. 96 (A4), 5847-5851.

Camargo, S.J. and Tasso, H. 1992 Renormalization group in magnetohydrodynamic turbulence. Phys. Fluids B4 (5), 1199-1212.

Cardy, J. 1996 Scaling and Renormalization in Statistical Physics. Cambridge.

Dannevik, W.P., Yakhot, V. and Orszag, S.A. 1987 Analytical theories of turbulence and the $`ϵ`$ expansion. Phys. Fluids 30 (7), 2021-2029.

De Dominicis, C. and Peliti, L. 1978 Field-theory renormalization and critical dynamics above $`T_c`$: helium, antiferromagnets and liquid-gas systems. Phys. Rev. B18, 353-376.

DeDominicis, C. and Martin, P.C. 1979 Energy spectra of certain randomly-stirred fluids. Phys. Rev. A 19, 419-422.

DeYoung, D.S. 1980 Turbulent Generation of Magnetic Fields in Extended Extragalactic Radio Sources. Ap. J. 241, 81-97.

DeYoung, D.S. 1992 Galaxy-Driven Turbulence and the Growth of Intracluster Magnetic Fields. Ap. J. 386, 464-472.

Elsasser, W.M. 1950 The hydromagnetic equations. Phys. Rev. 79, 183.

Enqvist, K. 1998 Primordial Magnetic Fields. Int. J. Mod. Phys. D7, 331-350.

Eyink, G.L. 1994 The renormalization group method in statistical hydrodynamics. Phys. Fluids 6 (9), 3063-3078.

Forster, D., Nelson, D.R. and Stephen, M.J. 1977 Large-distance and long-time properties of a randomly stirred fluid. Phys. Rev. A16, 732-749.

Fournier, J.-D., Sulem, P.-L. and Pouquet, A. 1982 Infrared properties of forced magnetohydrodynamic turbulence. J. Phys. A: Math. Gen. 15, 1393-1420.

Frey, E. and Tauber, U.C. 1994 Two-loop renormalization-group analysis of the Burgers- Kardar-Parisi-Zhang equation. Phys. Rev. E50, 1024-1044.

Frisch, U., Pouquet, A., Leorat, J. and Mazure, A. 1975 Possibility of an inverse cascade of magnetic helicity in magnetohydrodynamic turbulence. J. Fluid Mech. 68, 769-778.

Frisch, U. and Sulem, P.L. 1984 Numerical simulation of the inverse cascade in two-dimensional turbulence. Phys. Fluids 27, 1921-1923.

Frisch, U. 1995 Turbulence. Cambridge.

Gross, D.J. 1981 in Methods in Field Theory, Les Houches Summer School 1975 Session XXVIII, edited by R. Balian and J. Zinn-Justin. North-Holland/World Scientific

Hochberg, D., Pérez-Mercader, J., Molina-París, C. and Visser, M. 1999a Renormalization group improving the effective action: a review. Int. J. Mod. Phys. A14, 1485-1521.

Hochberg, D., Molina-París, C., Pérez-Mercader, J. and Visser, M. 1999b Effective action for stochastic partial differential equations. Phys. Rev. E60, 6343-6360.

Kraichnan, R.H. 1967 Inertial ranges in two-dimensional turbulence. Phys. Fluids 10, 1417-1423.

Kraichnan, R.H. 1973 Helical turbulence and absolute equilibrium. J. Fluid Mech. 59, 745-752.

Lilly, D. 1969 Numerical simulation of two-dimensional turbulence. Phys. Fluids Suppl. 12 II, 240-249.

Longcope, D.W. and Sudan, R.N. 1991 Renormalization group analysis of reduced magnetohydrodynamics with application to subgrid modeling. Phys. Fluids B3 (8), 1945-1962.

Ma, S.-K. 1976 Modern Theory of Critical Phenomena. Addison-Wesley

Martin, P.C., Siggia, E.D. and Rose, H.A. 1973 Statistical dynamics of classical systems. Phys. Rev. A8, 423-437.

McComb, W.D. 1995 The Physics of Fluid Turbulence. Oxford.

Medina, E., Hwa, T., Kardar, M. and Zhang, Y-C. 1989 Burgers equation with correlated noise: Renormalization-group analysis and applications to directed polymers and interface growth. Phys. Rev. A39, 3053-3075.

Papaloizou, J.C.B and Lin, D.N.C. 1995 Theory of accretion disks I: angular momentum transport processes. Annu. Rev. Astron. Astrophys. 33, 505-540.

Pouquet, A., Frisch, U. and Leorat J. 1976 Strong MHD helical turbulence and the nonlinear dynamo effect. J. Fluid Mech. 77, 321-354.

Pouquet, A. 1978 On two-dimensional magnetohydrodynamic turbulence. J. Fluid Mech. 88, 1-16.

Pouquet, A. and Patterson, G.S. 1978 Numerical simulation of helical magnetohydrodynamic turbulence. J. Fluid Mech. 85, 305-323.

Ramond, P. 1981 Field Theory: A Modern Primer Benjamin-Cummings, Reading, Massachusetts.

Rivers, R.J. 1987 Path integral methods in quantum field theory. Cambridge.

Sommeria, J. 1986 Experimental study of the two-dimensional inverse energy cascade in a square box. J. Fluid Mech. 170, 139-168.

Van Ballegooijen, A.A. 1985 Electric currents in the solar corona and the existence of magnetostatic equilibrium. Ap. J. 298, 421-430.

Van Ballegooijen, A.A. 1986 Cascade of magnetic energy as a mechanism of coronal heating. Ap. J. 311, 1001-1014.

Wilson K.G. and Kogut, J. 1974 The renormalization group and the $`ϵ`$ expansion. Physics Reports 12C, 75-200.

Yakhot, V. and Orszag, S.A. 1986 Renormalization-group analysis of turbulence. Phys. Rev. Lett. 57, 1722-1724.

Zinn–Justin, J. 1996 Quantum field theory and critical phenomena. Oxford.