Nonlinear energy transfers in accretion discs MRI turbulenceI-Net vertical field case

# Nonlinear energy transfers in accretion discs MRI turbulence I-Net vertical field case

G. Lesur and P.-Y. Longaretti Department of Applied Mathematics and Theoretical Physics, University of Cambridge, Centre for Mathematical Sciences, Wilberforce Road, Cambridge CB3 0WA, UK
Laboratoire d’Astrophysique, UJF CNRS, BP 53, 38041 Grenoble Cedex 9, France Laboratoire d’Astrophysique, UJF CNRS, BP 53, 38041 Grenoble Cedex 9, France
###### Key Words.:
accretion, accretion disks – instabilities – protoplanetary disks – turbulence

The magnetorotational instability (MRI) is believed to be responsible for most of the angular momentum transport in accretion discs. However, molecular dissipation processes may drastically change the efficiency of MRI turbulence in realistic astrophysical situations. The physical origin of this dependency is still poorly understood as linear and quasi linear theories fail to explain it. In this paper, we look for the link between molecular dissipation processes and MRI transport of angular momentum in non stratified shearing box simulations including a mean vertical field. We show that magnetic helicity is unimportant in the model we consider. We perform a spectral analysis on the simulations tracking energy exchanges in spectral space when turbulence is fully developed. We find that the energy exchanges are essentially direct (from large to small scale) whereas some non linear interactions appear to be non local in spectral space. We speculate that these non local interactions are responsible for the correlation between turbulent transport and molecular dissipation. We argue that this correlation should then disappear when a significant scale separation is achieved and we discuss several methods by which one can test this hypothesis.

## 1 Introduction

The transport of angular momentum in astrophysical discs is a central problem of accretion theory. To explain discs lifetime and accretion rates, it is often assumed that these objects are turbulent. Turbulence is then included in global models using a turbulent viscosity prescription as in the disc model (Shakura & Sunyaev, 1973).

The origin of this turbulence has been the subject of many debates over the past decades. It is now generally assumed that the magnetorotational instability, or shortly MRI (Velikhov, 1959; Chandrasekhar, 1960; Balbus & Hawley, 1991), is responsible for disc turbulence, although hydrodynamic processes might also be at work (Lesur & Ogilvie, 2010; Lesur & Papaloizou, 2010). Although MRI generated turbulence is generally efficient at transporting angular momentum (Hawley et al., 1995), recent results have shown a strong sensitivity of MRI turbulence on small scale dissipation processes (Lesur & Longaretti, 2007; Fromang et al., 2007), and in particular on the magnetic Prandtl number (ratio of microscopic viscosity to resistivity). This effect, called the correlation, casts doubts on the actual efficiency of the MRI in realistic situations since can vary by several order of magnitude in discs (Balbus & Henri, 2008). Several attempts have been made to explain this correlation, either from the linear theory of dissipative MRI modes (Pessah & Chan, 2008) or from the quasi-linear parasitic modes theory (Pessah & Goodman, 2009). However, these approaches were shown to be unsuccessful when compared to high Reynolds number simulations (Longaretti & Lesur, 2010). Instead, Longaretti & Lesur (2010) suggested that the correlation could be due to the nature of the MHD cascade in MRI generated turbulence, in which one might expect inverse cascades and/or non local interaction in spectral space. This kind of process would allow for a direct communication between the injection scales (transport scales) and the largest dissipation scale (either resistive or viscous).

The purpose of the present work is to investigate some of the conjectures presented by Longaretti & Lesur (2010) regarding the nature of the MHD turbulent cascade in accretion discs. To this end, we consider several of the high resolution simulations presented by Longaretti & Lesur (2010) and we analyse the energy exchanges in spectral space. This paper is organized as follows. We describe our model, equations and the spectral analysis we use in section 2. Section 3 is the core of this paper and discusses our numerical results. The implications of these results are presented in the final section.

## 2 Model and spectral analysis

### 2.1 Equations

In the following, we will adopt the shearing box model which accurately represent the local physics of an accretion disc (Hawley et al., 1995; Balbus, 2003; Regev & Umurhan, 2008). The adopted coordinate system is such that and where is the fiducial radius of the shearing box in the disc, and the azimuthal coordinate in the rotating frame. The velocity can be decomposed as a mean velocity plus a fluctuating part where is the local rotation rate and for a Keplerian disc rotation profile. As a simplification, we assume the flow is incompressible, which corresponds to the “small shearing box limit” of Umurhan & Regev (2004). The equations of motion then read

 ∂tv = −v⋅∇v−∇Π+B⋅∇B+qΩx∂yv (1) +2Ωvyex−(2−q)Ωvxey+ν∇2v, ∂tB = −v⋅∇B+B⋅∇v (2) +qΩx∂yB−qΩBxey+η∇2B, ∇⋅v = 0, ∇⋅B = 0,

where is the total pressure, the viscosity and the ohmic resistivity. In the following, we impose a mean vertical field which will be conserved during the evolution of the flow thanks to the shearing-sheet boundary conditions. It should be noted that the magnetic field strength is expressed in Alfvén speed for simplicity.

Several dimensionless numbers characterise the equations of motions. In this paper, we will use the following:

• The amplitude of the imposed mean vertical field measured by

 β=(qΩ)2L2zB20 (3)

where is the vertical box size. This definition mimics the usual plasma in vertically stratified discs obeying the vertical hydrostatic equilibrium constraint

• The viscous Elsasser number:

 Λν=B20Ων (4)

which is related to the Reynolds number used in Longaretti & Lesur (2010) by

• The resistive Elsasser number:

 Λη=B20Ωη (5)

connected to the magnetic Reynolds number by a similar relation

• The magnetic Prandtl number

 Pm=νη=ΛηΛν (6)

which compares the amount of viscous and resistive dissipation.

### 2.2 Fourier transform in sheared flows

It is convenient to introduce the shearing frame (, , ) :

 x = x′, y = y′−qΩtx′, z = z′.

Writing the equations of motions in the sheared frame allows us to eliminate the explicit spatial dependency:

 ∂tv = −v⋅∇v−∇Π+B⋅∇B (7) +2Ωvyex−(2−q)Ωvxey+ν∇2v, ∂tB = −v⋅∇B+B⋅∇v (8) −qΩBxey+η∇2B.

where and are now assumed to be functions of so that the nabla operator expression becomes

 ∇=ex(∂x′+qΩt∂y′)+ey∂y′+ez∂z′

In the sheared frame, the shearing sheet boundary conditions make every physical quantity periodic so that can be expanded in Fourier series:

 (9)

This relation defines the time dependent unsheared wave vectors:

 kx = kx′+qΩky′t, (10) ky = ky′, (11) kz = kz′. (12)

As expected, the application of the operator to corresponds to a multiplication of its Fourier components by .

This definition of the unsheared wave vectors with the Fourier decomposition (9) is usually referred to as a “shearing wave” decomposition of a sheared flow. It was first used by Lord Kelvin to study the stability of sheared flows (Thomson, 1887) and later in the astrophysical context by Goldreich & Lynden-Bell (1965) for spiral arms of galaxies.

### 2.3 Shell filter decomposition

Following Frisch (1995) and Alexakis et al. (2007), one defines shell filtered quantities in the unsheared Fourier space. At any given time, a series of linearly spaced shell sizes is defined from to ; by construction is the shell width (for any ). We define the shell-filtered field in shell by:

 XKj=∑Kj−δK/2

It should be noted that since depends on time, the exact number and the distribution of the modes entering the above formula for any given might vary in time, adding an extra complication compared to the homogeneous case of Alexakis et al. (2007). This point is discussed in the Appendices.

### 2.4 Energy transfer equations

The transfers we are interested in relate to the equations of the kinetic and magnetic energies. These shell-restricted, box-averaged equations involve a number of transfer functions that are introduced along with the related equations. We follow here the logic of Alexakis et al. (2007) and extend it to shear flows. The analysis of transfers gives indications about the locality of interactions in Fourier space.

Because of the incompressibility condition, energy transfers (but also stress and magnetic helicity) involve only at most the product of three components of the velocity and magnetic fields and their derivatives. Therefore, couplings in Fourier space depend only on triads of wave-vectors, noted , . The closing condition () furthermore imposes that at least two of the wave vectors are of the same magnitude; the third one can be either much smaller (implying nonlocal couplings through large scales) or of the same magnitude as the other two (implying locality of couplings in Fourier space).

Using the equations of motion in the sheared frame (78), one can derive the equation for the shell filtered energy density:

 ∂tEK = ∑Q[Tvv(Q,K)+Tbv(Q,K)]+Sv,K (13) +qΩIv,K−νDv,K ∂tMK = ∑Q[Tbb(Q,K)+Tvb(Q,K)]+Sb,K (14) +qΩIb,K−ηDb,K

where and . In the above expression, we have defined the transfer functions by

 Tvv(Q,K) = −⟨vK⋅(v⋅∇)vQ⟩, Tbb(Q,K) = −⟨BK⋅(v⋅∇)BQ⟩, Tbv(Q,K) = ⟨vK⋅(B⋅∇)BQ⟩, Tvb(Q,K) = ⟨BK⋅(B⋅∇)vQ⟩,

where denotes an spatial average on the shearing box volume. represents the transfer from energy ‘’ (kinetic or magnetic) from shell to energy ‘’ (kinetic or magnetic) in shell . Note that , so that whatever is taken from one shell of one type is totally transferred to the other shell. Similarly, : there is no effective transfer from a shell to itself, as should be. This justifies the identification of these quantities as shell-to-shell energy transfer functions; in effect, the third member of a triad is only a relay in effective energy exchanges between shells and (for a more detailed discussion, see Verma 2004). In these expressions, , has been used, as well as .

The next terms in the shell energy budget involve energy transfers due to the mean shear and . These terms are singular in time as they correspond to energy fluctuations due to wave entering or leaving the shell as evolves. They can be formally defined by:

 SX,K = ∑k′X∗k′Xk′2δ(t−tk′)ϵk′ = ∑k′qΩkykx(t)|k(t)|X∗k′Xk′2

where is the instant when the wave enters or exits the shell and for an entering/exiting wave (see appendix A-B).

The remaining terms in the shell kinetic and magnetic energy budgets,

 Iv,K = ⟨vx,Kvy,K⟩, Ib,K = −⟨Bx,KBy,K⟩, Dv,K = ⟨(∇×vK)2⟩, Db,K = ⟨(∇×BK)2⟩,

represent the energy injection through the shear and dissipation through viscosity and resistivity in shell .

### 2.5 Energy fluxes in spectral space

Using the transfer function defined above, it is possible to introduce energy flux in Fourier space. In this work, we will use the following fluxes:

 Fv(K0,t) = Kmax∑K=K0∑QTvv(Q,K) (15) Fb(K0,t) = Kmax∑K=K0∑QTbb(Q,K) (16) Fx(K0,t) = Kmax∑K=K0∑QTvb(Q,K)+Tbv(Q,K) (17) Fs,v(K,t) = ∑k′qΩkykx(t)|k(t)|V∗k′⋅Vk′2δ(|k(t)|−K) (18) Fs,b(K,t) = ∑k′qΩkykx(t)|k(t)|B∗k′⋅Bk′2δ(|k(t)|−K) (19)

where , , and are respectively the kinetic, magnetic, exchange and shear fluxes; the shear fluxes are evaluated in . Each of these fluxes computes the amount the energy transferred from shells to shells (i.e. the flux of energy “through” the outer boundary of shell ). The kinetic (respectively magnetic) fluxes compute transfers of kinetic to kinetic (respectively magnetic to magnetic) energy, whereas the exchange flux is a flux of total energy (magnetic plus kinetic) in which magnetic and kinetic energy are constantly transformed into one another. Shear fluxes are singular in time as they are non zero only when a wave enters or leave shells . It should be noted that shear fluxes are statistically non zero only for anisotropic turbulence, as the amplitude of any mode should be statistically equal to the amplitude of the mode in isotropic turbulence.

These fluxes allow one to check the direction of the energy transfers due to the nonlinear terms. Indeed, a direct energy transfer (large to small scale) implies a positive flux with the above definition, whereas an inverse cascade can be characterised by a negative flux.

### 2.6 Numerical method

Equation (1) and (2) are solved using the Snoopy code. Snoopy is a 3D spectral (Fourier) method based on the shearing wave decomposition (10)—(12). The Fourier transforms are evaluated using the FFTW 3 library, with both MPI and OpenMP parallelisation techniques. Nonlinear terms are computed using a pseudo-spectral algorithm (Canuto et al., 1988) and antialiasing is enforced using the “3/2” rule. Time integration is performed by a third order, low storage Runge-Kutta scheme for non-linear terms, whereas an implicit scheme is used for viscous and resistive terms. This spectral scheme uses a periodic remap algorithm in order to continuously follow the smallest wave number of the system in the sheared frame (see Umurhan & Regev 2004 appendix C for a complete description of the periodic remap algorithm). Note that the periodic remap method used in this code is different from the continuous remap method used by Lithwick (2007). Our main motivation to implement a periodic remap is the possibility of using the 3/2 antialiasing rule and power of 2 grid sizes111If the number of point in one direction is a multiple of 2, waves at the Nyquist frequency do not have any properly defined phase. This is not a problem for classical spectral methods or for the periodic remap method since the Nyquist frequency is either in the dissipation range or in the antialiasing dump zone. However, when using a continuous remap method, the Nyquist frequency waves are remapped to large scale waves in physical space, which might lead to non-physical behaviours. for which Fourier transform and parallelisation methods are more efficient. This code or its variant has been used in several context, including the MRI (Longaretti & Lesur, 2010) and the subcritical baroclinic instability (Lesur & Papaloizou, 2010). It is available for download on the author’s website.

## 3 Results

### 3.1 Simulations parameters and averaging procedure

The spectra and transfers presented in this section are all derived from two simulations of an MRI saturated state. These runs corresponds to the and high resolution runs discussed in Longaretti & Lesur (2010). Both runs have a resolution222The quoted resolution now includes the aliasing domain. This convention differs from the one adopted in our previous papers, where only the ”useful” domain was accounted for when quoting resolutions. with a box aspect ratio . We impose a mean vertical field in the box with and (). Each simulation is integrated for 50 orbits starting from random noise and the spectra are averaged from the last 40 orbits to remove any influence of the initial conditions. The two simulations considered in this section only differ by their ohmic resistivity, the run having () whereas the run has ().

The statistical average of any quantity of interest should in principle be computed on different realizations but are evaluated in practice as usual via an ergodic hypothesis:

 ⟨X⟩stat=limT→±∞1T∫T0 X(t)dt≃1T∫T0 X(t)dt≃∑iX(Ti), (20)

where are a sufficiently large number of instants of flow snapshots.

The spectra are averaged in the spherical shells introduced in (2.3). The shells are defined so that and . This means that some power is present in the shell , as it contains large scale horizontal waves with no vertical structure. The shell-integrated spectra and transfers obtained by this procedure are then averaged in time over 40 instantaneous snapshots (1 snapshot per orbit). For simplicity, we have renormalized the wavevectors so that on all the plots in this section. Note that shells are incomplete in the direction since the resolution per scale height is smaller in that direction. This is not a problem since these shells are in the dissipative range and high modes are weaker than the equivalent high modes due to the anisotropy of MRI turbulence (see section 3.2). Note that, with the procedure we have used, one can reconstruct the box averaged quantities by summing the spectra over the integers K.

The shear transfer terms and are computed in a special way. Indeed, one cannot compute for a given shell and snapshot time numerically because of the functions. Instead, we introduce a shell-averaged flux:

 Fcs(K0,t)=1δK∫K0+δK/2K0−δK/2dKFs(K,t). (21)

As depend only on (the turbulence is statistically stationary) and varies little with on scales of the order of , one has . One can therefore use in the averaging procedure described above to estimate .

The numerical flux we obtain is then averaged over time following the procedure described above.

### 3.2 Spectra and energy injection

We first present the energy spectra on Fig. 1 for and . The standard deviation, measured from 40 instantaneous snapshots, is shown as a shaded region on these spectra. This dispersion is due to temporal fluctuations of the turbulence intensity. The most obvious feature observed in these spectra is the presence of a spectrum for the kinetic energy; the traditional Kolmogorov scaling appears to be excluded in the run, but it cannot be strictly excluded in the run. A power-law was also found in zero net flux MRI turbulence (Fromang, 2010), although the spectrum shape differs, our spectra being exempt of any “bump” at intermediate scale. As our runs do not yet resolve the inertial range of the turbulent cascade, these apparent spectral shapes require some comment. The presence of a spectrum is usually related to the theoretical argument of Iroshnikov (1963) and Kraichnan (1965) (or shortly IK). However MRI turbulence is not strongly magnetized, and therefore falls outside the domain of validity of the IK phenomenology. Moreover, the magnetic field spectrum does not follow any well-defined power law, as expected from the wide and overlapping injection (see below) and dissipation spectra, indicating that the spectrum we get is not an inertial spectrum. We are therefore forced to conclude that although the kinetic spectrum looks like an IK or Kolmogorov spectrum, it is not described by the IK or Kolmogorov phenomenologies, nor by recent extensions (Boldyrev, 2005).

Changing the magnetic Prandtl number does not change the power-law index for the kinetic energy. We note however 2 major effects: the overall spectra amplitudes are reduced and the dissipation scales move to larger scale as one reduces . These two effects are expected since it is known that smaller turbulence is associated with a smaller transport efficiency and therefore a weaker injection of energy in the cascade. This effect is confirmed by the injection spectra (Fig. 2) which are significantly reduced at smaller .

We note that the energy injection peaks at the largest scale of the box, although injection still exist at . Therefore, although a power-law spectrum is found for , this spectrum cannot be described as an “inertial range” since energy is still injected at these intermediate scales.

We present on Fig. 3 bidimensional spectra of magnetic energy for (kinetic spectra are not shown as they share essentially the same properties). These spectra were obtained averaging 3D energy spectra over 40 orbits and taking the average in the , , or directions. We first note a strong anisotropy in the plane which indicates that trailing shearing waves () have more energy than leading shearing waves (). As we will see below, this results in non zero shear transfer terms.

Looking at the aspect ratio of the energy contours, we see that turbulence is slightly less anisotropic at large than at small (the contours are less “elongated” at large ), although complete isotropy is not yet reached in this simulation. Let us however remark that the spectral truncation (due to the finite resolution) tends to deform the contours at large , which might accelerate the return to isotropy. One should therefore perform higher resolution runs (or at least double ) in order to confirm this return to isotropy. In principle, one would expect a return to isotropy at small scales if the nonlinear transfer terms dominate all the other terms (injection, body forces) at large enough . However this is not always the case (e.g. in the presence of a strong mean magnetic field).

The spectrum shows that turbulence is essentially isotropic at large in that plane. For , we find a slight anisotropy where modes with are favoured. This is probably a result of large scale MRI unstable modes which all have in the presence of a mean vertical field. Note that this anisotropy disappears very quickly as one moves to larger .

### 3.3 Magnetic helicity and cross helicity

The presence of kinetic and/or magnetic helicity in MHD turbulence is often invoked to explain large scale dynamo action. Indeed, it is known that an inverse cascade of magnetic helicity can appear in fully developed helical MHD turbulence (Frisch et al., 1975), potentially leading to a buildup of magnetic energy at large scale. Some kinematic effects often used in mean field dynamos, like the effect (Moffatt, 1978), also lead to the generation to large scale helical fields (Brandenburg & Subramanian, 2005). Magnetic helicity has therefore been suggested as a possible driving mechanism (or at least a tracer) of disc dynamos (Blackman, 2010). Moreover, several authors have tried to link magnetic helicity conservation and magnetic helicity flux to the saturation properties of the MRI (Vishniac, 2009; Käpylä & Korpi, 2010).

In this work, we define the magnetic helicity , where is the vector potential of the fluctuation (this expression is gauge-invariant in the shearing box); denotes a volume average. We show on Fig. 4 (left) the spectrum of relative helicity for the run. As it can be seen, the relative helicity is less than 1% for all scales of these simulations. Moreover, this quantity is strongly fluctuating and its sign is not well defined333Note that the absolute value of the relative magnetic helicity is plotted on Fig. 4.

These results tend to indicate that magnetic helicity is dynamically unimportant in the unstratified simulations presented here and that MRI saturation is not related to a magnetic quenching effect due to magnetic helicity accumulated at large scale. This was to be expected in the first place as unstratified shearing boxes (both with and without mean field) are mirror symmetric. Note however that this picture might change when stratification is included.

An other quantity of interest which might play a role in the MHD turbulence cascade is the cross helicity (see e.g. Perez & Boldyrev 2010). When cross helicity is non zero, the energy of Alfvén waves traveling along and against the mean field are not equal. For this reason, turbulence with cross helicity is often called imbalanced turbulence. Locally imbalanced turbulence is often observed in strong MHD turbulence, where the guide field is weaker than the turbulent fluctuations. To check whether non stratified MRI turbulence was imbalanced, we have computed cross helicity spectra of our simulations [Fig. 4 (right)]. As for the magnetic helicity, we find that the relative cross helicity is small () and highly fluctuating at all scales. This tends to indicate that MRI turbulence is not imbalanced in our setup. As for magnetic helicity, this result was to be expected because of the mirror symmetry properties of the non stratified shearing box. The absence of any significant cross helicity also shows that energy spectra in Elsässer variables are equal and proportional to the kinetic plus magnetic energy spectrum.

### 3.4 Energy fluxes

The energy fluxes (15)-(19) allow one to check the average direction of the energy flux in spectral space. To explain the dependence of the turbulent transport on , several authors (e.g. Lesur & Longaretti 2007; Fromang et al. 2007) have suggested that an inverse cascade driven by resistive and viscous scales might be at work. Since magnetic helicity is irrelevant to this problem, only the kinetic, magnetic, exchange and shear fluxes are important for the non stratified shearing box and should be checked for an inverse cascade.

We present the energy fluxes at and on Fig. 5. Standard deviations are shown for kinetic and magnetic fluxes as shaded regions. These deviations are computed following the procedure described in section 3.2. We always find positive fluxes, meaning that the non linear transfers are forward or direct (from large to small scales) on average. However, at larger scale, the standard deviation may allow for an inverse cascade of kinetic energy in the run. This indicates that, although the kinetic cascade is direct on (time) average over most of the spectrum, inverse cascades can sometimes be observed on the largest scales. This inverse cascade of kinetic energy could be an explanation to the large scale hydrodynamic structures which are observed in several MRI turbulence simulations, such as vortices (Fromang & Nelson, 2005) and zonal flows (Johansen et al., 2009).

We also observe that the energy cascade is dominated by the magnetic and exchange fluxes down to the resistive scale, the kinetic flux and shear fluxes being almost negligible. Note also that the shear fluxes are always positive. This is due to the anisotropy of MRI turbulence in which shearing waves have statistically a larger amplitude than leading waves (see section 3.2).

In the case, the kinetic flux becomes dominant at subresistive scale (), indicating that the cascade becomes essentially hydrodynamic below the resistive scale, as expected for small MHD turbulence (see also Fig. 6). Moreover, the kinetic flux dominates the kinetic shear transfer term at least for the larger , which indicates that, as far as the non linear transfers are concerned, the cascade is close to isotropic at small scales, as noted in section 3.2. Finally, note that none of the fluxes reach a plateau at intermediate scales, which would be expected in the presence of an inertial range. This indicates that the kinetic spectrum found in Fig. 1 is not properly speaking an inertial spectrum.

### 3.5 Energy transfer locality

To test the locality of the energy transfer in spectral space, we have plotted the transfer functions in the case for several values of : at the injection scale (), in the intermediate range () and in the resistive range (). We first plot the kinetic to kinetic and magnetic to magnetic transfers on Fig. 7. The transfers and obtained at all scales perfectly illustrate local energy exchanges. Energy is taken from wavenumbers slightly smaller than and is transferred to wavenumbers slightly larger than , except (not surprisingly) for . As expected from the energy flux, we also find that the cascade is direct, with energy going from small to large wavenumbers, Finally, note that the transfers are always much smaller than the transfers above the dissipation range, illustrating the magnetically dominated energy transfer described above.

We next plot the exchange transfers and on Fig. 8. We note in this case that the scales involved in each transfer are much broader. In particular, measured at the resistive scale () has contribution coming from all scales, including the largest injection scales. This effect can also be seen in , which exhibit a very long tail toward large , down to the resistive scale. Comparing directly these transfers to and show that these terms are highly non local. In turns, this indicates that the exchange flux computed in the previous section is non local. The results in the case are not shown here as they are very similar to the case.

We note that despite the non locality of the energy exchanges, the overall cascade direction is still forward, confirming our previous interpretation regarding the exchange flux. We also remark that the shear transfer terms are local by definition as they transfer energy to neighbouring shells.

## 4 Discussion

In this paper, we have described some properties of the turbulent cascade found in incompressible MRI turbulence. We have shown that compared to isotropic MHD turbulence, the presence of a mean shear led to several new transfer terms and introduced a source of anisotropy. We have computed the effect of each non linear term and found that all the terms contribute to a direct cascade of energy (from large to small scales) but some terms involved non local transfers in Fourier space. This non locality is due to the Lorentz force and to the magnetic stretching term of the induction equation (combined here in the exchange transfer term). We have also shown that magnetic helicity, although non zero, was totally negligible and should not play any role in the behaviour of MRI turbulence.

The presence of non local transfer terms in the MRI turbulent cascade is the most important finding of this work. It indicates that in principle, the large scales – responsible for the transport – can directly interact with the small dissipative ones through non linear terms. This direct interaction could of course explain the correlation observed between and the turbulent transport of angular momentum (Longaretti & Lesur, 2010). However, it should be pointed out that nonlocal transfers were already found in isotropic MHD turbulence (e.g. Alexakis et al., 2007). Therefore, MRI turbulence has nothing special regarding the nature of these non linear transfers.

Although some nonlinear terms are found to directly connect injection and dissipation scales in current simulations, one might wonder if this would be true in a more realistic setup where the injection and dissipation scales are separated by a wide range of scales (typically ). In other words, what is the maximum scale separation these terms can connect? A partial answer to this question is given by Aluie & Eyink (2010). In order to describe their result, let us define the structure functions:

 δvl,p=⟨|v(x+l)−v(x)|p⟩

In the inertial range the structure function depends only on and , where is the structure function index of order . It is then possible to derive an upper bound to the nonlinear transfer terms thanks to the Hölder inequality. Applying this procedure to the non local transfer , Aluie & Eyink (2010) found

 |Tub(Q,K)|<(const)Q1−ζu3/3K−2ζb3/3 (22)

where and are dyadic (octave) wavenumbers and . Similar terms can be obtained for and . If one assumes Iroshnikov-Kraichnan theory, one has . On the contrary, considering Goldreich-Sridhar (GS) phenomenology, which should be valid for MRI turbulence, one gets and . In all these cases, (22) indicates that the non-locality of these transfer terms cannot extend over several decades, with a typical scaling for GS turbulence ( being the usual turbulence energy injection rate).

We therefore conclude that the non locality in Fourier space is somewhat relative. Although and are non local compared to and , these terms should be local when one considers transfers over several decades. Unfortunately, separating the injection scale from the dissipative scales by several decades is numerically difficult. It is even harder for MRI turbulence since the injection term is rather broad in spectral space compared to forced turbulence. Assuming the injection and dissipation scales both spread over one decade in Fourier space, one typically needs simulations to get a 2 decades inertial range in which non local transfers are significantly reduced. This kind of resolution is for the moment out of reach of the best computational facilities.

Nevertheless, we can conjecture that if the correlation is effectively due to the non local transfers, then it should vanish when the injection and dissipation scales are well separated as it is the case in some accretion discs. Although this conclusion looks rather reassuring for the relevance of today simulations regarding small scale dissipation, it does not tell us what the asymptotic value of is in this limit, nor how MRI turbulence behaves when the scale separation is not achieved, a situation which probably occurs in the inner regions of protoplanetary discs where is not very large.

## Appendix A Shearing waves approach to the shear transfer term

The shell-filter decomposition can be properly defined using a projector operator on the field in the sheared frame:

 [ΠKj(F)](x′,t)=∑k′∈Σj(t)Fk′(t)exp(ik′⋅x′)

where is the shell containing all the shearing waves with a norm between and :

 Σj(t)={Kj−δK/2<|k′+qΩk′ytex|≤Kj+δK/2} (23)

Our notation indicates that the projected is function of space and time. Note also that is real for real fields .

As it can be seen, the waves included in change with time. This is to be expected as shearing waves should in principle enter and exit shells as the one defined above (see Fig. 9). As a result, the projector operator has an explicit time dependence which leads to non trivial transfer effects. The above projector operator can be written using Heaviside function :

 [ΠKj(F)](x′,t) = ∑k′Θ(|k(t)|−Kj+δK/2) (24) ×Θ(−|k(t)|+Kj+δK/2) ×Fk′(t)exp(ik′⋅x′).

where we have defined as a function of as in (10)—(12).

One next defines the energy within a shell where denotes a volume average. The energy equation then reads

 dEFKjdt = ⟨ΠKj(F)∂ΠKj(F)∂t⟩ (25) =

The first term on the right handside leads to the terms obtained in isotropic turbulence (Alexakis et al., 2007) and introduced in Eqs. (13) and (14) along with the injection terms and . The second one however is due to shearing waves entering and exiting the shells. Using (24), it is possible to obtain an exact (though singular) expression of the operator time derivative:

 ∂t[ΠKj](F) = ∑k′qΩkykx(t)|k(t)| (26) ×[δ(|k(t)|−Kj+δK/2)−δ(−|k(t)|+Kj+δK/2)] ×Fk′(t)exp(ik′⋅x′).

This expression can be interpreted easily. As an example, let us consider waves with . Then, the first function represents waves entering the shell, the second delta represents waves leaving the same shell, and the factor in front of the delta functions quantifies the “flux” of waves going through a shell boundary. This interpretation is similar to the phenomenological picture one can have of waves traveling through a fixed shell in the unsheared Fourier space (Fig. 9).

We then deduce from (25)

 dEXKjdt = ⟨ΠKj(F)ΠKj(∂tF)⟩ +∑k′qΩkykx(t)|k(t)|F∗k′Fk′2 ×[δ(|k(t)|−Kj+δK/2)−δ(−|k(t)|+Kj+δK/2)].

where we have used the property444This can be shown differentiating the relation . As expected, this expression shows two contributions to the energy evolution inside a shell: a volume contribution and a surface contribution equal to the energy of the waves entering and leaving the shell. Introducing the equation of motions in sheared space (7)—(8) in the above relation leads to the energy equations (13)—(14).

Note finally that using the relation

 ∂t|k(t)|=qΩkykx(t)|k(t)|, (28)

one can write (A) as

 dEXKjdt=⟨ΠKj(F)ΠKj(∂tF)⟩+∑k′F∗k′Fk′2δ(t−tk′)ϵk′ (29)

where is the instant when the wave enters or exits the shell and for an entering/exiting wave. This somewhat simpler expression has the same interpretation as (A).

## Appendix B Unsheared Fourier transform approach to the shear transfer term

It is possible to understand the origin of the shear transfer term (26) starting from the equations of motion in unsheared coordinates (1)—(2) with an appropriate use of continuous Fourier transforms. Shear-periodic functions are not absolutely integrable, but this difficulty can be circumvented because their continuous Fourier transform is well-defined as a distribution. To demonstrate this point, let us consider a 2D infinite medium in which a field obeys the model equation

 ∂tF(x,t)−qΩx∂yF(x,t)=0. (30)

Let us introduce the unsheared Fourier transform:

 ~F(k,t) = 1(2π)3∭dxF(x,t)exp(−ik⋅x) F(x,t) = ∭dk~F(k,t)exp(ik⋅x). (31)

If obeys the shearing sheet boundary conditions, the solution of the equation is a Fourier series of the form (9) with Fourier coefficients . Its Fourier transform in the unsheared spectral space is then:

 ~Fp(k,t)=∑k′δ(k−k0(t))Fk′

where and where we have defined the “Fourier velocity” , evaluated in . By construction is solution of:

 ∂t~F(k)+qΩky∂kx~F(k)=0, (32)

which is the Fourier transform of our model equation555The Fourier transform of (30) gives (32) except for an extra term that vanishes once the Fourier series expression of is used.. The left handside of this equation can be interpreted as a comoving derivative of in Fourier space with the Fourier velocity . This equation tells us that the amplitude of the waves is constant when one moves at velocity in Fourier space consistently with the form of the solution . It can also be interpreted as constant amplitude shearing waves, as expected.

It is then possible to introduce the projector operator, now time-independent as it is defined in unsheared coordinates:

 [ΠKj(F)](x,t) = ∭dkΘ(|k|−Kj+δK/2) Ê ×Θ(−|k|+Kj+δK/2) (33) ×~F(k,t)exp(ik⋅x).

As previously, the shell energy () time variation follows from:

 dEFKjdt = (34) = −⟨ΠKj(F)ΠKj(∇k⋅(Vk~F))⟩

where (32) has been used in the last equality as well as . Because the volume average selects the zero frequency contributions, this leads us back to (A) with the help of the relation

 ∫dkxΘ(|k|−C)∂kx~F(k) = −∫dkxkx|k|δ(|k|−C)~F(k).

and with the use of and (B).

In the spirit of the integral expressions used in this appendix, Eq. (A) can be recast in integral form by introducing the energy density in Fourier space defined by

 EF(k,t)=∑k′δ(k−k0(t))Fk′F∗k′2. (35)

With this definition,

 EFKj=∭dkΘ(|k|−K−)Θ(−|k|+K+) EF(k,t) (36)

while

 dEFKjdt = (37) = −\oiint∂KjdkVk⋅n EF(k,t) (38)

In these relations, has been defined; the second integration is performed on the surface of the shell and is the normal to this surface. The last expression has an explicit flux form.

This approach can be applied to the original shearing sheet MHD equations; the time dependence due to the shear term will produce the shear flux contribution just computed.

###### Acknowledgements.
GL acknowledges support by STFC. PYL acknowledges the hospitality of the Isaac Newton Institute and of the DAMTP in Cambridge where parts of this research has been conducted. GL thanks S. Fromang for his comments on the initial version of the manuscript. This work was granted access to the HPC resources of IDRIS under allocation x2009042231 made by GENCI (Grand Equipement National de Calcul Intensif).

## References

• Alexakis et al. (2007) Alexakis, A., Mininni, P. D., & Pouquet, A. 2007, New Journal of Physics, 9, 298
• Aluie & Eyink (2010) Aluie, H. & Eyink, G. L. 2010, Physical Review Letters, 104, 081101
• Balbus (2003) Balbus, S. A. 2003, ARA&A, 41, 555
• Balbus & Hawley (1991) Balbus, S. A. & Hawley, J. F. 1991, ApJ, 376, 214
• Balbus & Henri (2008) Balbus, S. A. & Henri, P. 2008, ApJ, 674, 408
• Blackman (2010) Blackman, E. G. 2010, Astronomische Nachrichten, 331, 101
• Boldyrev (2005) Boldyrev, S. 2005, ApJ, 626, L37
• Brandenburg & Subramanian (2005) Brandenburg, A. & Subramanian, K. 2005, Phys. Rep., 417, 1
• Canuto et al. (1988) Canuto, C., Hussaini, M. Y., Quarteroni, A., & Zang, T. A. 1988, Spectral methods in fluid dynamics (Springer)
• Chandrasekhar (1960) Chandrasekhar, S. 1960, Proc. Nat. Acad. Sci., 46, 253
• Frisch (1995) Frisch, U. 1995, Turbulence (Cambridge University Press)
• Frisch et al. (1975) Frisch, U., Pouquet, A., Leorat, J., & Mazure, A. 1975, Journal of Fluid Mechanics, 68, 769
• Fromang (2010) Fromang, S. 2010, A&A, 514, L5+
• Fromang & Nelson (2005) Fromang, S. & Nelson, R. P. 2005, MNRAS, 364, L81
• Fromang et al. (2007) Fromang, S., Papaloizou, J., Lesur, G., & Heinemann, T. 2007, A&A, 476, 1123
• Goldreich & Lynden-Bell (1965) Goldreich, P. & Lynden-Bell, D. 1965, MNRAS, 130, 125
• Hawley et al. (1995) Hawley, J. F., Gammie, C. F., & Balbus, S. A. 1995, ApJ, 440, 742
• Iroshnikov (1963) Iroshnikov, P. S. 1963, AZh, 40, 742
• Johansen et al. (2009) Johansen, A., Youdin, A., & Klahr, H. 2009, ApJ, 697, 1269
• Käpylä & Korpi (2010) Käpylä, P. J. & Korpi, M. J. 2010, ArXiv e-prints
• Kraichnan (1965) Kraichnan, R. H. 1965, Phys. Fluids, 8, 1385
• Lesur & Longaretti (2007) Lesur, G. & Longaretti, P.-Y. 2007, MNRAS, 378, 1471
• Lesur & Ogilvie (2010) Lesur, G. & Ogilvie, G. I. 2010, MNRAS, 404, L64
• Lesur & Papaloizou (2010) Lesur, G. & Papaloizou, J. C. B. 2010, A&A, 513, A60+
• Lithwick (2007) Lithwick, Y. 2007, ApJ, 670, 789
• Longaretti & Lesur (2010) Longaretti, P. & Lesur, G. 2010, A&A, 516, A51+
• Moffatt (1978) Moffatt, H. K. 1978, Magnetic field generation in electrically conducting fluids (Cambridge University Press)
• Perez & Boldyrev (2010) Perez, J. C. & Boldyrev, S. 2010, Physics of Plasmas, 17, 055903
• Pessah & Chan (2008) Pessah, M. E. & Chan, C.-k. 2008, ApJ, 684, 498
• Pessah & Goodman (2009) Pessah, M. E. & Goodman, J. 2009, ApJ, 698, L72
• Regev & Umurhan (2008) Regev, O. & Umurhan, O. M. 2008, A&A, 481, 21
• Shakura & Sunyaev (1973) Shakura, N. I. & Sunyaev, R. A. 1973, A&A, 24, 337
• Thomson (1887) Thomson, W. 1887, Phil. Mag., 24, 188
• Umurhan & Regev (2004) Umurhan, O. M. & Regev, O. 2004, A&A, 427, 855
• Velikhov (1959) Velikhov, E. P. 1959, Sov. Phys.-JETP, 36, 1398
• Verma (2004) Verma, M. K. 2004, Phys. Rep, 401, 229
• Vishniac (2009) Vishniac, E. T. 2009, ApJ, 696, 1021
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters