# Local cascade and dissipation in incompressible Hall magnetohydrodynamic
turbulence: the Coarse-Graining approach

D. Manzini  Laboratoire de Physique des Plasmas (LPP), CNRS, École
Polytechnique, Sorbonne Université, Université Paris-Saclay, Observatoire de
Paris, 91120 Palaiseau, France  Dipartimento di Fisica E.Fermi, University of
Pisa, Italy  F. Sahraoui  Laboratoire de Physique des Plasmas (LPP), CNRS,
École Polytechnique, Sorbonne Université, Université Paris-Saclay,
Observatoire de Paris, 91120 Palaiseau, France  F. Califano  Dipartimento di
Fisica E.Fermi, University of Pisa, Italy  R. Ferrand  Laboratoire de Physique
des Plasmas (LPP), CNRS, École Polytechnique, Sorbonne Université, Université
Paris-Saclay, Observatoire de Paris, 91120 Palaiseau, France

###### Abstract

We derive the coarse-graining (CG) equations of incompressible Hall
Magnetohydrodynamics (HMHD) turbulence to investigate the local (in space)
energy cascade rate as a function of the filtering scale $\ell$. First, the CG
equations are space averaged to obtain the analytical expression of the mean
cascade rate. Its application to 3 dimensional (3D) simulations of (weakly
compressible) HMHD shows a cascade rate consistent with the value of the mean
dissipation rate in the simulations and with the classical estimates based on
the “third-order” law. Furthermore, we developed an anisotropic version of CG
that allows us to study the magnitude of the cascade rate along different
directions with respect to the mean magnetic field. Its implementation on the
numerical data with moderate background magnetic field shows a weaker cascade
along the magnetic field than in the perpendicular plane, while an isotropic
cascade is recovered in the absence of a background field. The strength of the
CG approach is further revealed when considering the local-in-space energy
transfer, which is shown theoretically and numerically to match at a given
position $\bm{x}$, when locally averaged over a neighboring region, the
(quasi-)local dissipation. Prospects of exploiting this new model to
investigate local dissipation in spacecraft data are discussed.

##  I Introduction

Turbulence plays a key role in space and astrophysical plasmas as, for
instance, it mediates energy conversion stored on large-scale fields into
particle heating and/or acceleration at smaller scales. The standard theory of
turbulence for hydrodynamics predicts an energy cascade from large scales,
where it is injected, to the small scales where it is dissipated by viscosity
stemming from particle collisions at the microscopic level (Batchelor, 1953;
Frisch, 1995; Monin and Yaglom, 2013). On the other hand, the low density and
high temperature conditions of most heliospheric plasmas make them nearly
collisionless. In those plasmas a turbulent cascade and particle heating are
frequently observed, however, the precise mechanisms by which the turbulent
fluctuations of the electromagnetic fields and plasma flow are damped still
elude our full understanding (Bruno and Carbone, 2013; Matthaeus and Velli,
2011; Goldstein  _et al._ , 2015; Sahraoui  _et al._ , 2020).

A key step in answering these fundamental questions is to identify and
characterize the regions of plasma involved in intense cross-scale energy
transfers. A popular tool that has been widely used in turbulence studies, in
particular those based on spacecraft observations, is the so-called “third-
order” law: a statistical relation that links the mean energy cascade rate
(equal in the formalism to the rate of energy injection and dissipation) to
the turbulent fluctuations at a given scale. The theoretical models used range
from Incompressible MHD (Politano and Pouquet, 1998) to more complex systems
that involve density fluctuations and/or small (sub-ion) scale effects
(Banerjee and Galtier, 2013; Andrés and Sahraoui, 2017; Andrés  _et al._ ,
2018; Hellinger  _et al._ , 2018; Ferrand  _et al._ , 2019, 2021a). All these
studies have greatly helped to gain deeper insight into the turbulence
dynamics in a variety of heliospheric plasmas, including the solar wind (SW)
(Smith  _et al._ , 2006; Sorriso-Valvo  _et al._ , 2007; Marino  _et al._ ,
2008; Stawarz  _et al._ , 2010; Osman  _et al._ , 2011; MacBride  _et al._ ,
2008; Coburn  _et al._ , 2015; Banerjee and Galtier, 2016; Hadid  _et al._ ,
2017) and planetay magnetospheres (Hadid  _et al._ , 2018; Andrés  _et al._ ,
2020; Sorriso-Valvo  _et al._ , 2019). The “third-order” law, albeit rigorous
and derived under fairly non-restrictive hypotheses, does however require
ensemble averages computed as time and/or space average when applied to
simulations or spacecraft data under the assumption of ergodicity (Frisch,
1995). As such, those laws fail to describe cross-scale energy transfer in
localized regions of space. To overcome this shortcoming some heuristic tools
have been proposed such as the Local Energy Transfer (LET) (Sorriso-Valvo  _et
al._ , 2019), which relaxes the statistical average used in the “third-order”
law, or the Partial Variance Index (PVI) (Greco  _et al._ , 2008; Chasapis
_et al._ , 2015) used to localize regions of space with large magnetic shear,
a proxy to identify regions of strong electric current. However, those tools
lack a solid theoretical foundation, which is mandatory to justify their use
as a means to measure energy rates. The present work fills this gap by
providing a robust theoretical model based on filtering (or coarse-graining)
the Hall-MHD (HMHD) equations that retains spatial locality while allowing one
to recover results consistent with the “third-order” law once spatially
averaged. We furthermore show analytically and numerically that the local (in
space) energy cascade rate across a scale $\bm{\ell}$ is a good proxy to
measure local dissipation within limited regions of space.

Note that other local theories of turbulence based on the concept of inertial
dissipation have been proposed in recent years to include the role of
discontinuities in dissipating energy, which are not rigorously accounted for
in the “third-order” law formalism (Galtier, 2018; Dubrulle, 2019). Those
models have been recently used to compute dissipation within discontinuities
observed in spacecraft data (Kuzzay  _et al._ , 2019; David and Galtier, 2021;
David  _et al._ , 2022)

##  II Coarse-Grained Incompressible HMHD equations

We introduce the notion of a coarse grained (CG) measurement of a field $f$
with a scale resolution $\ell$. We choose a kernel function $G$ normalized to
one $\int d\bm{r}G(\bm{r})=1$, centered $\int d\bm{r}\bm{r}G(\bm{r})=0$ and
with variance of order unity $\int d\bm{r}|\bm{r}|^{2}G(\bm{r})\simeq 1$. For
any given scale $\ell$ we define $G_{\ell}(\bm{x})=\ell^{-3}G(\bm{x}/\ell)$ so
that the normalization is preserved but the variance is now of order
$\ell^{2}$. The coarse-graining operation is defined as

$\bar{f}_{\ell}(\bm{x})=\int d\bm{r}G_{\ell}(\bm{r})f(\bm{x}+\bm{r}).$ (1)

It represents a local average of $f$ on a spatial region of radius $\sim\ell$
centered around the point $\bm{x}$. This convolution smooths out the
fluctuations with scale smaller than $\ell$, and gives a coarser
representation of the field, hence the name. Coarse-Graining a field $f$ at a
given scale $\ell$ individuates two different quantities: the large-scale
field $\bar{f}_{\ell}$, which in virtue of the smoothing operation retains
only the scales $>\ell$ (in $\mathbf{k}$-space the wave-vectors
$|\bm{k}|\lesssim 1/\ell$) and the _un-resolved_ or _sub-scale_ field
$f^{\prime}=f-\bar{f}_{\ell}$ which accounts for scales $<\ell$ (wave-vectors
$\bm{k}\gtrsim 1/\ell$). The term “unresolved” comes from the fact that when
coarse-graining the fields we choose to resolve fluctuations down to the scale
$\ell$ only. To keep the notation simple we will omit the subscript (filtering
scale) $\ell$ when not strictly necessary and denote $\bar{f}_{\ell}$ simply
as $\bar{f}$. The quality of the filtering in $\bm{k}$ space depends on the
choice of the filtering function $G_{\ell}$. For instance, a sharp spectral
filter such as $G_{\ell}^{sp}=\Pi_{i=1}^{3}sin(\pi x_{i}/\ell)/(\pi x_{i})$
allows one to clearly separate between scales but it looses the spatial
locality, while a box filter $G^{box}_{\ell}=\ell^{-3},\;|x_{i}|<\ell/2$,
$G^{box}_{\ell}=0,\;|x_{i}|\geq\ell/2$ is local in space, but does not allow
for unambiguous separation between scales. Other filters with intermediate
properties can be defined such as the Gaussian filter
$G^{Gss}_{\ell}=(2/\pi\ell^{2})^{3/2}exp\\{-2|\bm{x}|^{2}/\ell^{2}\\}$ used in
this work (Meneveau and Katz, 2000; Eyink and Aluie, 2009).
We start from the incompressible Hall-MHD (IHMHD) equations normalized to
Alfvén units:

$\begin{split}&\partial_{t}\bm{u}=-(\bm{u}\cdot\nabla)\bm{u}+(\bm{b}\cdot\nabla)\bm{b}-\nabla
P+\bm{d}_{\nu}+\bm{f}\\\
&\partial_{t}\bm{b}=\nabla\times(\bm{u}\times\bm{b})-d_{i}\nabla\times(\bm{j}\times\bm{b})+\bm{d}_{\eta}\\\
&\nabla\cdot\bm{u}=0\end{split}$ (2)

where $\bm{b}={\bf B}/\sqrt{\mu_{0}\rho_{0}}$ is the (scale dependent) Alfvén
speed, $\rho_{0}$ the mean plasma density, ${\bm{j}=\nabla\times\bm{b}}$,
$d_{i}$ is the ion inertial length, $P=p/\rho_{0}+|\bm{b}|^{2}/2$ is the total
pressure, $\bm{d}_{\nu},\bm{d}_{\eta}$ are the velocity and magnetic field
dissipation terms, respectively, and $\bm{f}$ is an external force injecting
energy at large scales. The CG operation is a convolution and therefore
commutes with space and time derivatives. The equations filtered at scale
$\ell$ are readily obtained by convolving equations (2) with the filtering
kernel $G_{\ell}$

$\displaystyle\partial_{t}\bar{\bm{u}}=-(\bar{\bm{u}}\cdot\nabla)\bar{\bm{u}}+(\bar{\bm{b}}\cdot\nabla)\bar{\bm{b}}-\nabla\bar{P}+\bar{\bm{d}_{\nu}}-\nabla\cdot\bm{\tau}+\bar{\bm{f}}$
(3a)
$\displaystyle\partial_{t}\bar{\bm{b}}=\nabla\times\left[(\bar{\bm{u}}-d_{i}\bar{\bm{j}})\times\bar{\bm{b}}\right]+\bar{\bm{d}}_{\eta}+\nabla\times(\mathbfcal{E}_{\text{MHD}}+d_{i}\mathbfcal{E}_{\text{Hall}})$
(3b)

where we introduced the notations $\tau(f,g)=\overline{fg}-\bar{f}\bar{g}$
and, for the sake of readability, we define the second order tensor
$\tau_{ij}=\tau(u_{i},u_{j})-\tau(b_{i},b_{j})=\overline{u_{i}u_{j}}-\bar{u}_{i}\bar{u}_{j}-\left(\overline{b_{i}b_{j}}-\bar{b}_{i}\bar{b}_{j}\right)$.

These equations describe the dynamics of the CG (large scale) fields and
closely resemble the HMHD equations (2). The difference lies in the presence
of additional contributions stemming from the filtering of the nonlinear
terms. These quantities represent the action of the “unresolved” scales
($<\ell$) on the filtered fields. In particular, in equation (3a) we find the
divergence of the subscale Reynolds and Maxwell stress tensors
$\tau_{ij}=\tau(u_{i},u_{j})-\tau(b_{i},b_{j})$, while in equation (3b) we
find the curl of the subscale electric field in the MHD limit
$\mathbfcal{E}_{\text{MHD}}=\overline{\bm{u}\times\bm{b}}-\bar{\bm{u}}\times\bar{\bm{b}}$
and the correction due to the Hall term
$\mathbfcal{E}_{\text{Hall}}=-\left(\overline{\bm{j}\times\bm{b}}-\bar{\bm{j}}\times\bar{\bm{b}}\right)$.

Multiplying equation (3a) by $\bar{\bm{u}}$ and equation (3b) by
$\bar{\bm{b}}$ we obtain the time evolution of the large scale kinetic and
magnetic energy densities:

$\displaystyle\begin{split}\partial_{t}\frac{|\bar{\bm{u}}|^{2}}{2}=&-\nabla\cdot\left[\frac{|\bar{\bm{u}}|^{2}}{2}\bar{\bm{u}}+\bar{P}\bar{\bm{u}}-(\bar{\bm{u}}\cdot\bar{\bm{b}})\bar{\bm{b}}+\bm{\tau}\cdot\bar{\bm{u}}\right]\\\
-&\bar{\bm{b}}\cdot(\bar{\bm{b}}\cdot\nabla)\bar{\bm{u}}+\bar{\bm{u}}\cdot\bar{\bm{d}_{\nu}}-\pi^{u}+\bar{\bm{u}}\cdot\bar{\bm{f}}\end{split}$
(4a)
$\displaystyle\begin{split}\partial_{t}\frac{|\bar{\bm{b}}|^{2}}{2}=&-\nabla\left[\frac{|\bar{\bm{b}}|^{2}}{2}\bar{\bm{u}}+d_{i}(\bar{\bm{j}}\times\bar{\bm{b}})\times\bar{\bm{b}}-(\mathbfcal{E}_{\text{MHD}}+d_{i}\mathbfcal{E}_{\text{Hall}})\times\bar{\bm{b}}\right]\\\
&+\bar{\bm{b}}\cdot(\bar{\bm{b}}\cdot\nabla)\bar{\bm{u}}+\bar{\bm{b}}\cdot\bar{\bm{d}_{\eta}}-\pi^{b,\text{MHD}}-d_{i}\pi^{b,\text{Hall}}\end{split}$
(4b)

where we introduced the quantities:

$\pi^{u}=-\partial_{j}\bar{u}_{i}\tau_{ij}\quad\pi^{b,\text{MHD}}=-\bar{\bm{j}}\cdot\mathbfcal{E}_{\text{MHD}}\quad\pi^{b,\text{Hall}}=-\bar{\bm{j}}\cdot\mathbfcal{E}_{\text{Hall}}.$
(5)

These terms, the study of which is the main focus of this work, are the local
(in space) energy transfers across the scale $\ell$. They appear as a sink in
the large scale energy equations (4) and a source in the small scales ones
(see Appendix B).

While these equations allow us to analyze separately the magnetic and kinetic
energy cascades, here we rather focus on the study of the cascade of the total
energy (note that in the energy balance equation the internal energy is not
included as it is a conserved quantity in incompressible pressure-isotropic
flows (Simon and Sahraoui, 2021)). Summing equations (4a)-(4b) and separating
the cascade rate into its MHD component,
$\pi^{\text{MHD}}=\pi^{u}+\pi^{b,\text{MHD}}$, and Hall component
$\pi^{b,Hall}$, hereafter denoted simply as $\pi^{\text{Hall}}$, we find:

$\begin{split}&\partial_{t}\left(\frac{|\bar{\bm{u}}|^{2}+|\bar{\bm{b}}|^{2}}{2}\right)\\\
+&\nabla\cdot\left[\frac{|\bar{\bm{u}}|^{2}+|\bar{\bm{b}}|^{2}}{2}\bar{\bm{u}}+\bar{P}\bar{\bm{u}}-(\bar{\bm{u}}\cdot\bar{\bm{b}})\bar{\bm{b}}+d_{i}(\bar{\bm{j}}\times\bar{\bm{b}})\times\bar{\bm{b}}\right]\\\
+&\nabla\cdot\left[\bm{\tau}\cdot\bar{\bm{u}}-(\mathbfcal{E}_{I}+d_{i}\mathbfcal{E}_{H})\times\bar{\bm{b}}\right]\\\
-&\bar{\bm{u}}\cdot\bar{\bm{d}_{\nu}}-\bar{\bm{b}}\cdot\bar{\bm{d}_{\eta}}-\bar{\bm{u}}\cdot\bar{\bm{f}}=\\\
-&\pi^{\text{MHD}}-d_{i}\pi^{\text{Hall}}\end{split}$ (6)

This equation is the starting point of our study and the first rigorous result
of this paper, which extends previous results Aluie (2017) to IHMHD. Equation
(6) shows that, at a given position $\bm{x}$, the time variation of the large
scale energy density is a balance between the spatial advection due to large
and small scale fields (the second and third lines, respectively), the effects
of dissipation and forcing at large scales (fourth line) and, lastly, the
energy transfer across the scale $\ell$, namely $\pi^{\text{MHD}}$ and
$d_{i}\pi^{\text{Hall}}$. We recall that in equation (6) the filtering scale
$\ell$ can be varied to gauge the magnitude of each term as a function of
scale and its choice individuates two separate range of the spectrum: the
resolved large scales (corresponding to wavevectors $|\bm{k}|\lesssim 1/\ell$)
and the “unresolved” small scales ($|\bm{k}|\gtrsim 1/\ell$).

##  III Space Integration and the Cascade Rate

Equation (6) describes the temporal evolution of the large scale energy
density. By performing a spatial average over the whole domain, denoted
$\langle\rangle$, and assuming no net flux at the boundaries, we recover the
following expression for the temporal evolution of the mean large scale energy
$\bar{E}=\left\langle|\bar{\bm{u}}|^{2}+|\bar{\bm{b}}|^{2}\right\rangle$/2,

$\partial_{t}\bar{E}=-\Pi+\langle\bar{\bm{u}}\cdot\bar{\bm{d}}_{\nu}+\bar{\bm{b}}\cdot\bar{\bm{d}}_{\eta}\rangle+\langle\bar{\bm{u}}\cdot\bar{\bm{f}}\rangle$
(7)

where we introduced
$\Pi=\Pi^{\text{IMHD}}+d_{i}\Pi^{\text{Hall}}=\langle\pi^{\text{IMHD}}\rangle+d_{i}\langle\pi^{\text{Hall}}\rangle$.

All the quantities involved in Eq.(7) are functions of the filtering scale
$\ell$ only, as the spatial dependence is lost when averaging over the spatial
domain. We see that the large scale energy is affected by three processes: the
forcing mechanism that injects energy via the term
$\langle\bar{\bm{u}}\cdot\bar{\bm{f}}\rangle$, the effect of dissipation at
scales larger than the filtering scale $\ell$ given by the term
$\langle\bar{\bm{u}}\cdot\bar{\bm{d}}_{\nu}+\bar{\bm{b}}\cdot\bar{\bm{d}}_{\eta}\rangle$
and, lastly, the large scale energy transfer due to nonlinearities given by
$\Pi(\ell)=-\partial_{t}\bar{E}|_{\text{NL}}$, where we denote with
$\partial_{t}(\cdot)|_{\text{NL}}$ the rate of change due to non linear
processes. $\Pi(\ell)$ stands for the cross-scale energy transfer rate, that
is: the amount of energy flowing from the large-scale, resolved quantities (
$>\ell$) to the small-scales “unresolved” ones ($<\ell$). In Fourier space it
writes

$\Pi(\ell)=-\frac{\partial}{\partial t}\left[(2\pi)^{3}\int
d\bm{k}\left(\frac{|\hat{\bm{u}}(\bm{k})|^{2}+|\hat{\bm{b}}(\bm{k})|^{2}}{2}\right)|\hat{G}_{\ell}(\bm{k})|^{2}\right]_{\text{NL}}$
(8)

where the squared modulus of the filtering function
$|\hat{G}_{\ell}(\bm{k})|^{2}$ plays the role of a low-pass filter. Eq.(8)
shows that $\Pi$ closely resembles the formal definition of the cascade rate
across the wave-vector $K=1/\ell$ (Frisch, 1995):

$-\frac{\partial}{\partial t}\left[(2\pi)^{3}\int_{|\bm{k}|\leq
K}d\bm{k}\left(\frac{|\hat{\bm{u}}(\bm{k})|^{2}+|\hat{\bm{b}}(\bm{k})|^{2}}{2}\right)\right]_{\text{NL}}$
(9)

if one replaces the sharp cut-off in ${\bf k}$-space with a gentler slope
given by the shape of $|\hat{G}_{\ell}(\bm{k})|^{2}$. In line with this view,
for a given (short) interval $\Delta t$, the variation of the large scale
energy due to nonlinear processes is $-\Pi(\ell)\Delta t$. This energy is
directly transferred to small scales $<\ell$ so that a positive $\Pi(\ell)$ is
the signature of a direct cascade with energy flowing from large to small
scales, while a negative value yields an inverse cascade (i.e., from small to
large scales).

It is our understanding that the CG quantity $\Pi$ represents an extension
over the current state of the art “third-order” laws that, under a set of
assumptions, provide an estimate of the cascade rate, often denoted
$\varepsilon$. Indeed “third order” laws rely heavily on the assumption of the
existence of an inertial range and on the ergodicity hypothesis to compute
ensemble averages as time and/or space averages, while equation (6) holds
without the need of those assumptions. Thus, it allows us to estimate the
transfer rate at any scale, i.e. not necessarily in the inertial range, and
regardless of how important is dissipation (Eyink, 2018). This remark might
prove to be useful in particular in collisionless plasmas where dissipation
(via, e.g., Landau damping) may occur at all scales, which would question the
very existence of the inertial range (Ferrand  _et al._ , 2021b).

Nevertheless, if Kolmogorov hypotheses are satisfied we expect the two
quantities $\Pi$ and $\varepsilon$ to converge to the same value. This can be
readily explained by recalling the formal definition of $\varepsilon$ as the
time evolution of the auto-correlation function due to nonlinear terms
(Frisch, 1995):

$\varepsilon(\bm{l})=-\frac{1}{2}\partial_{t}\langle\bm{u}(\bm{x})\cdot\bm{u}(\bm{x}+\bm{l})+\bm{b}(\bm{x})\cdot\bm{b}(\bm{x}+\bm{l})\rangle|_{\text{NL}}$
(10)

where $\langle\rangle$ denotes an ensemble average, which is computed as a
space average under the assumption of ergodicity. Using Parseval theorem we
can recover from the definition of $\varepsilon$, the expression:

$\varepsilon(\bm{l})=-\frac{\partial}{\partial t}\left[(2\pi)^{3}\int
d\bm{q}\left(\frac{|\hat{\bm{u}}(\bm{k})^{2}|+|\hat{\bm{b}}(\bm{k})^{2}|}{2}\right)e^{-i\bm{q}\cdot\bm{l}}\right]_{\text{NL}}$
(11)

and performing a Fourier transform we obtain:

$\int
d\bm{k}\varepsilon(\bm{l})e^{i\bm{k}\cdot\bm{l}}=-(2\pi)^{3}\frac{\partial}{\partial
t}\left[\frac{|\hat{\bm{u}}(\bm{k})^{2}|+|\hat{\bm{b}}(\bm{k})^{2}|}{2}\right]_{\text{NL}}.$
(12)

Lastly, substitution of relation (12) in the definition of $\Pi$ given by
Eq.(8) yields:

$\Pi(\ell)=\int
d\bm{l}^{\prime}\varepsilon(\bm{l}^{\prime})\phi_{\ell}(\bm{l}^{\prime})$ (13)

where the function $\phi_{\ell}(\bm{l}^{\prime})=(2\pi)^{-3}\int
d\bm{k}|\hat{G}_{\ell}(\bm{k})|^{2}e^{i\bm{k}\cdot\bm{l}^{\prime}}$ is
directly related to the shape of the filtering function $G_{\ell}$. For the
Gaussian filter used in this work, we have
$\phi_{\ell}(\bm{l}^{\prime})=(\pi\ell^{2})^{-3/2}\exp\\{-|\bm{l}^{\prime}|^{2}/\ell^{2}\\}$
so that the two quantities $\Pi$ and $\varepsilon$ differ by a Gaussian
smoothing operation. It appears therefore that when Kolmogorov hypotheses are
verified (and as such $\varepsilon$ is constant in the inertial range) the two
quantities have the same value and $\varepsilon$ is a good estimate of the
cascade rate $\Pi$.

###  III.1 Numerical Validation

In this section we will apply the theoretical results obtained above to
simulation data that feature 3D freely-decaying weakly compressible HMHD
turbulence (see Ferrand  _et al._ (2022)). The two simulations were performed
using the Fourier pseudo-spectral code GHOST (Gómez  _et al._ , 2005; Mininni
_et al._ , 2011) on a $N=1024^{3}$ grid, spanning a real space cubic domain of
side $L_{0}=50\times 2\pi d_{i}$ and a grid spacing of $\Delta\sim 0.05\times
2\pi d_{i}$. The two simulations were performed with two different values of
the background field, namely $B_{0}=0$ (hereafter, run I) and $B_{0}=2$ (run
II), which is very convenient to test the anisotropic CG approach that we
introduce further below. In both runs dissipation is implemented via viscous
and resistive terms and the dimensionless viscosity and magnetic diffusivity
are taken to be equal $\nu=\eta=3\times 10^{-4}$. The data cubes were taken
when turbulence was deemed to have reached a fully developed state.
First, we estimate all terms of equation (7) using the simulation data of run
I ($B_{0}=0$). To improve the graphical representation we add and subtract to
the left side of equation (7) $-\partial_{t}E$, which by definition in freely
decaying turbulence is $\partial_{t}E=-\varepsilon_{diss}$, so that we can
write:

$\varepsilon_{diss}=-\partial_{t}(E-\bar{E})+\Pi-\langle\bar{u}\cdot\bar{\bm{d}_{\nu}}+\bar{b}\cdot\bar{\bm{d}_{\eta}}\rangle$
(14)

where we did not include the term
$\langle\bar{\bm{u}}\cdot\bar{\bm{f}}\rangle$ as there is no forcing in our
simulations. The different terms of Eq. (14) are shown in Fig.1.

Figure 1: Evaluation of the terms in Eq. (14) for run I ($B=0$) normalized to
the mean dissipation rate. Dashed lines indicate negative values while the
dotted horizontal line shows the reference value
$\varepsilon/\varepsilon_{diss}=1$. The black line represents the sum of all
the terms.

We observe that the cascade rate (blue curve) in the central region
$(0.3\lesssim\ell\lesssim 3)$ is representative of the dissipation rate. As
$\ell$ tends to zero the effects of dissipation (red curve) on the large
scales, _i.e._ larger than the filtering scale $\ell$, become comparable to
the total dissipation rate. Furthermore, the contribution of the Hall term to
the energy cascade starts to increase around $\ell\sim d_{i}$ however, it
remains negligible compared to the Ideal MHD contribution. This could be due
to the small scale separation between the dissipative scales and $d_{i}$
(Ferrand  _et al._ , 2022). The black line, which is the sum of the four terms
displayed, remains constant within $10\%$ of the total dissipation rate for
all values of $\ell$.

###  III.2 Link with the “third-order” law of IHMHD

We now turn to the assertion that $\Pi$ matches the classical ”third-order”
law theory of $\varepsilon$ when the latter is applicable. The standard
approach to obtain an expression for $\varepsilon$ is to derive a generalized
von Kármán–Howarth (vKH) equation for IHMHD turbulence (Banerjee and Galtier,
2017; Hellinger  _et al._ , 2018; Ferrand  _et al._ , 2019), which describes
the time evolution of the spatial auto-correlation function
$R_{E}=\langle\bm{u}(\bm{x})\cdot\bm{u}(\bm{x}+\bm{l})+\bm{b}(\bm{x})\cdot\bm{b}(\bm{x}+\bm{l})\rangle/2$
where now the brackets $\langle\rangle$ denote formally an _ensemble average_.
This dynamical equation has contribution from both linear processes, such as
the forcing or dissipation mechanisms and nonlinear ones. The latter are of
particular interest as they provide the cross-scale interactions needed to
sustain the energy cascade. The key quantity under study is the contribution
to the rate of change $\partial_{t}R_{E}$ stemming from nonlinear processes.
This quantity is denoted as
$\varepsilon(\bm{l})=-\partial_{t}R_{E}(\bm{l})|_{\text{NL}}$. In the present
work we use the vKH equation of IHMHD turbulence derived in Banerjee and
Galtier (2017) (hereafter BG17). Denoting quantities evaluated at
$\bm{x}+\bm{l}$ with a prime (_e.g._ $\bm{u}(\bm{x}+\bm{l})=\bm{u}^{\prime}$)
and defining the field increments as
$\delta\bm{u}=\bm{u}(\bm{x}+\bm{l})-\bm{u}(\bm{x})=\bm{u}^{\prime}-\bm{u}$,
the BG17 law reads:

$\begin{split}&\partial_{t}R_{E}(\bm{l})=\\\
&-\frac{1}{2}\langle\delta[\bm{u}\times\bm{w}+\bm{j}\times\bm{b}]\cdot\delta\bm{u}\rangle-\frac{1}{2}\langle\delta[(\bm{u}-d_{i}\bm{j})\times\bm{b}]\cdot\delta\bm{j}\rangle\\\
&+\frac{1}{2}\langle(\bm{u}^{\prime}\cdot\bm{d}_{\nu}+\bm{u}\cdot\bm{d}^{\prime}_{\nu})+(\bm{u}^{\prime}\cdot\bm{d}_{\eta}+\bm{u}\cdot\bm{d}^{\prime}_{\eta})\rangle\\\
&+\frac{1}{2}\langle\bm{u}\cdot\bm{f}^{\prime}+\bm{u}^{\prime}\cdot\bm{f}\rangle\end{split}$
(15)

where we can identify the effects of non linearities (second line) which can
be split into the MHD component from the Hall one (proportional to $d_{i}$),
$\partial_{t}R_{E}|_{NL}=-[\varepsilon_{\text{IMHD}}+d_{i}\varepsilon_{\text{Hall}}]$,
in addition to the dissipation and forcing terms (third and fourth line,
respectively), which are simply denoted $\mathcal{D},\mathcal{F}$ in the
following.

We evaluated all the terms in equation (15) using the data from run I, but the
term $\mathcal{F}$ that does not apply here because of the free-decay nature
of our simulations. Furthermore to improve the graphical representation we add
and subtract to the left side of (15) $-\partial_{t}E$ which by definition, in
freely decaying turbulence simulations, corresponds to the total dissipation
rate $\partial_{t}E=-\varepsilon_{diss}$. We can re-write equation (15) in a
more compact form:

$\varepsilon_{diss}=\partial_{t}(R_{E}-E)+\varepsilon_{\text{IMHD}}+d_{i}\varepsilon_{\text{Hall}}-\mathcal{D}$
(16)

The right hand side terms of this equation are plotted in Fig.2 where the
spatial lag $\bm{l}$ is directly related to the scale
$\bm{\ell}=\bm{l}/(2\pi)$.

Figure 2: Numerical verification of the BG17 Law, the energy rates normalized
to the mean dissipation rate. Negative values are represented with dashed
lines, while the dotted horizontal line signals
$\varepsilon/\varepsilon_{diss}=1$. The black line represents the sum of all
the terms.

The results in Fig.2 show that the sum of the different terms remains constant
at all scales. We observe a good resemblance between the cascade rates given
by our CG model (see Fig.1) and those given by the “third-order” law of IHMHD
(Fig.2). This is better emphasized when looking at Fig.3 where the various
quantities are compared.

Figure 3: Comparison between the cascade rates $\Pi$ and $\varepsilon$ (BG17
law) for run I. The negative values are shown as dashed lines for
$\varepsilon$ and as empty markers for $\Pi$.

##  IV Anisotropic Cascade Rate

In the presence of a background magnetic field $B_{0}\neq 0$ plasma turbulence
becomes anisotropic (Shebalin  _et al._ , 1983; Matthaeus  _et al._ , 1984;
Goldreich and Sridhar, 1995; Schekochihin  _et al._ , 2009), that is energy
cascades preferentially in the direction perpendicular to the mean magnetic
field. Turbulence anisotropy has been essentially investigated by looking at
the energy spectra in $(k_{\parallel},k_{\perp})$ or, equivalently, at the
second-order structure functions (Cho and Vishniac, 2000; Sahraoui  _et al._ ,
2006; Meyrand and Galtier, 2013; Sahraoui  _et al._ , 2010; Chen  _et al._ ,
2010). Here we are interested in analyzing directly the anisotropy of the
cascade rate. It is therefore mandatory to extend the CG approach beyond an
isotropic treatment. Instead of using a spherically symmetric filtering kernel
$G_{\ell}(\bm{x})$ we can define a more general filter $G_{\bm{\ell}}(\bm{x})$
that has different characteristic widths $\ell_{x},\ell_{y},\ell_{z}$ in the
three real space directions. The CG quantity $\bar{f}=f\ast
G_{\bm{\ell}}(\bm{x})$ is low-pass filtered in an anisotropic way aimed at
retaining mainly wave-vectors $k_{x}\lesssim 1/\ell_{x},k_{y}\lesssim
1/\ell_{y},k_{z}\lesssim 1/\ell_{z}$. In this way we can highlight possible
presence of plasma anisotropy. In the limit case when the filtering scales
along two directions go to zero, _e.g._ $\ell_{x},\ell_{z}\to 0$, computing
the quantity $\Pi$ allows one to recover the rate of energy flowing from
$k_{y}\lesssim 1/\ell_{y}$ to $k_{y}\gtrsim 1/\ell_{y}$, thus recovering the
1-D cascade rate in the $y$ direction. Continuing our example, we choose a
filtering function
$G_{\bm{\ell}}(\bm{x})=\delta(x)\delta(z)\psi_{\ell_{y}}(y)$, where $\delta$
and $\psi_{\ell_{y}}$ are respectively the Dirac and a 1-D filtering kernel
with a characteristic width $\ell_{y}$. Then, the quantity
$\Pi=-\partial_{t}\bar{E}|_{\text{NL}}$ can be written as:

$\Pi_{\bm{\ell}}=-\frac{\partial}{\partial t}\left[(2\pi)^{3}\int
dk_{x}dk_{z}\int
dk_{y}\left(\frac{|\hat{\bm{u}}(\bm{k})|^{2}+|\hat{\bm{b}}(\bm{k})|^{2}}{2}\right)|\hat{\psi}_{\ell_{y}}(k_{y})|^{2}\right]_{\text{NL}}$

corresponding to the one-dimensional cascade rate along direction $y$.

The designed scheme is implemented on the two sets of simulations data
obtained from run I and run II and the results are shown in Fig.4.

Figure 4: Anisotropic cascade rate for run I ($B_{0}=0$, top) and run I
($B_{0}=2$, bottom) along the three space directions $x,y,z$. A weaker cascade
is evidenced in the $z$ direction when $B_{0}=2$.

We found that in the absence of a strong mean field the cascade rate is almost
isotropic at most scales, with only weak predominance of the cascade in the
direction $x$ at the smallest scales. In the presence of a background field
(run II, $B_{0}=2$) the cascade rate in the parallel direction is weaker in
comparison with the two perpendicular directions, which look overall similar
(but at the largest scales). In both runs, the violation of gyrotropy at the
largest scales is likely to be a residual effect of the initial modes used to
inject energy at largest scales of the simulation box.

A global picture of how energy flows in the $(\ell_{\parallel},\ell_{\perp})$
plane can be obtained by computing the cascade rate
$\Pi(\ell_{\perp},\ell_{\parallel})$, which measures the amount of energy that
goes across wavevectors $|\bm{k}_{\perp}|\sim
1/\ell_{\perp},k_{\parallel}=1/\ell_{\parallel}$. This is achieved by using
the filtering function
$G_{\ell_{\perp},\ell_{\parallel}}=\varphi_{\ell_{\perp}}(x,y)\psi_{\ell_{\parallel}}(z)$.
Implementation of this procedure on the same simulations data as above yields
the results plotted in Fig.5. Run I ($B_{0}=0$, Fig.5(a) ) shows an almost
isotropic behaviour with only a slightly weaker cascade rate in the $z$
direction, while in run II ($B_{0}=2$, Fig.5(b)) the cascade develops
preferentially in the perpendicular plane with a very weak dependence on
$\ell_{\parallel}$.

Figure 5: 2D Cascade rate $\Pi(\ell_{\perp},\ell_{\parallel})$ normalized to
the dissipation rate for run I ($B_{0}=0$, left) and run II ($B_{0}=2$,
right).

The results of Fig.5 showing the anisotropic energy transfers when $B_{0}\neq
0$, are in agreement with theoretical expectation for MHD and HMHD turbulence.
However, to the best of our knowledge, this is the first time that a full 3D
anisotropic cascade rate is directly revealed, while previous studies dealt
with energy spectral anisotropy Meyrand and Galtier (2013); Cho and Vishniac
(2000), which is rather a consequence of the anisotropic cascade shown in
Fig.5. It is worth noticing that a similar plot of anisotropic cascade is
obtained experimentally in rotating fluid turbulence (Lamriben  _et al._ ,
2011)

##  V The cascade equation for the small scale energy

It is possible to derive an equation complementary to Eq. (6) that describes
the time evolution of the energy density contained in scales smaller than the
filtering scale $\bm{\ell}$. This sub-scale energy density is defined as
$h_{\ell}=\left[\tau_{\ell}(u_{i},u_{i})+\tau_{\ell}(b_{i},b_{i})\right]/2$
(see Aluie (2017)) and its time evolution is given by (see Appendix B):

$\begin{split}\partial_{t}h_{\ell}=&-\nabla\cdot\mathbfcal{J}^{\text{LS}}-\nabla\cdot\mathbfcal{J}^{\text{SS}}-\nabla\cdot\mathbfcal{J}^{\text{p}}\\\
&+\pi_{\ell}+\tau_{\ell}(u_{i},d_{\nu,i})+\tau_{\ell}(u_{i},d_{\eta,i})+\tau_{\ell}(u_{i},f_{i})\end{split}$
(17)

where all quantities are function of both the position $\bm{x}$ and the scale
$\ell$.
In the first line we find
$\nabla\cdot\mathbfcal{J}^{\text{LS}},\nabla\cdot\mathbfcal{J}^{\text{SS}},\nabla\cdot\mathbfcal{J}^{\text{p}}$
(whose expressions are given in Appendix B), which describe the spatial
transport of $h_{\ell}$ due to large and small scale fields and pressure
interactions, respectively. The local in space cascade rate
$\pi_{\ell}(\bm{x})$ appears as a sink in the large scale equation (6) and as
a source here, its full expression is given by:

$\begin{split}\pi_{\ell}(\bm{x})=&-\partial_{j}\bar{u}_{i}(\tau_{\ell}(u_{i},u_{j})-\tau_{\ell}(b_{i},b_{j}))-\epsilon_{\alpha\beta\gamma}\bar{j}_{\alpha}\tau(u_{\beta},b_{\gamma})\\\
&+d_{i}\varepsilon_{\alpha\beta\gamma}\bar{j}_{\alpha}\tau(j_{\beta},b_{\gamma})\end{split}$
(18)

where we introduced the Levi-Civita tensor $\varepsilon_{\alpha\beta\gamma}$,
with the usual summation rule over repeated indices. We stress here that
$\pi_{\ell}(\bm{x})$ is the only term able to exchange energy between the
large and small scales at a given position $\bm{x}$. The other terms instead
are associated either to spatial transport or to forcing/dissipation. These
last two processes enter the small scale energy equation as a difference of
filtered terms, _e.g._
$\tau_{\ell}(u_{i},d_{\nu,i})=\overline{(\bm{u}\cdot\bm{d}_{\nu})}_{\ell}-\bar{\bm{u}}_{\ell}\cdot\bar{\bm{d}_{\nu}}_{\ell}$.
To aid in the physical interpretation we show (see Appendix A) that by
averaging over a spatial region of characteristic size $L\gg\ell$ we recover
the relation

$\langle\tau(f,g)\rangle_{L}\approx\langle f^{\prime}g^{\prime}\rangle_{L}$

where $f^{\prime}=f-\bar{f},g^{\prime}=g-\bar{g}$ are the “unresolved”
(subscale) fluctuations. In this view, we can interpret the
forcing/dissipation terms in equation (17) as the contributions to these
processes coming from scales $<\ell$. For this reason we can write:

$\begin{array}[]{c}\langle\tau_{\ell}(u_{i},f_{i})\rangle_{L}\approx\langle\bm{u}^{\prime}\cdot\bm{f}^{\prime}\rangle_{L}\approx
0\\\
\langle\tau_{\ell}(u_{i},d_{\nu,i})+\tau_{\ell}(u_{i},d_{\eta,i})\rangle_{L}\approx\langle\bm{u}^{\prime}\cdot\bm{d}_{\nu}^{\prime}+\bm{b}^{\prime}\cdot\bm{d}_{\eta}^{\prime}\rangle_{L}\end{array}$
(19)

where in the first line we stated that the forcing injects energy at large
scales only and in the second line we recover the contrtibution of dissipation
due to scales smaller than $\ell$.

##  VI Cascade and Local Dissipation

One of the main results of the Kolmogorov theory of turbulence is that the
mean cascade rate (_i.e._ averaged over the whole simulation box) is
representative of the dissipation rate. We want to show that this holds even
when averaging over much smaller spatial regions. In particular, we will prove
that at small filtering scales $\ell$, and under some assumptions, the local
cascade rate and the local dissipation rate match _quasi-locally_. In other
words, the amount of energy cascading across the (small) scale $\ell$ at
position $\bm{x}$ will be eventually dissipated at close locations. Choosing a
region of characteristic size $L$, we will derive the smallest $L$ for which
the local cascade rate and the local dissipation match when spatially averaged
over such region.
More precisely, we will prove that if the region size $L$ satisfies the two
inequalities $L\gg\ell$,$L\gg\ell B_{0}/\delta b$ (the latter becoming
redundant in strong turbulence with $\delta b/B_{0}>1$) then the following
relation holds:

$\langle\pi_{\ell}(\bm{x})\rangle_{L}\approx-\langle\bm{u}^{\prime}\cdot\bm{d}^{\prime}_{\nu}+\bm{b}^{\prime}\cdot\bm{d}^{\prime}_{\eta}\rangle_{L}=\langle\varepsilon_{diss}^{<\ell}(\bm{x})\rangle_{L}$
(20)

where on the right side of Eq. (20) we obtain the _local_ dissipation due to
scales smaller than $\ell$ averaged on a region of size $L$. We stress that
the size of the region is not fixed and can be varied at will. In general we
expect the agreement to get better as $L\to L_{box}$, eventually recovering
the “global” result of Section III, however, we will show that this relation
holds for regions of smaller size, effectively allowing us to study
(quasi-)locally the processes of cascade and dissipation.

The starting point in deriving equation (20) is the small scales
($\lesssim\ell$) equation (17). For each position $\bm{x}$, we average over a
region $\chi(\bm{x},L)$, centered on $\bm{x}$, of characteristic size $L$ and
of volume $V(L)\sim L^{3}$. For each quantity $q$ we can write

$\langle q\rangle_{L}=V(L)^{-1}\int_{\chi(\bm{x},L)}q(\bm{r})d\bm{r}.$

This operation (_de-facto_ a new CG operation) is linear and as such it
commutes with space and time derivatives. We can therefore write:

$\partial_{t}\langle
h_{\ell}\rangle_{L}=-\nabla\cdot\langle\mathbfcal{J}\rangle_{L}+\langle\pi_{\ell}\rangle_{L}+\langle\tau_{\ell}(u_{i},d_{\nu,i})+\tau_{\ell}(b_{i},d_{\eta,i})\rangle_{L}$
(21)

where $\mathbfcal{J}=\mathbfcal{J}^{LS}+\mathbfcal{J}^{SS}+\mathbfcal{J}^{p}$.
We can use this equation to study the characteristic time scales associated
with the different terms.
We recall that there are two distinct scales: $\ell$ is the filtering scale
across which energy cascades; $L$ is the size of the spatial region upon which
we compute the spatial average. The filtering scale $\ell$ helps us to
introduce two characteristic velocities: the large (at scales $>\ell$) and
small ($<\ell$) scale velocity $\bar{u},u^{\prime}$, respectively. The same
filtering can be applied for the magnetic field with the only difference that
the mean field cannot be removed by a Galilean transformation. However, we can
decompose the characteristic large scale field as
$\bar{b}=B_{0}+\delta\bar{b}$. The small scale field $\bm{b}^{\prime}$is not
affected by $B_{0}$. We furthermore assume the following Alfvénic ordering:

$\bar{u}\sim\delta\bar{b},\quad u^{\prime}\sim b^{\prime}$

and denote these quantities as $\bar{U},U^{\prime}$, which apply indistinctly
to $u$ and $b$. Furthermore, we use
$\bar{\bm{j}}\sim\delta\bar{b}/\ell,\bm{j}^{\prime}\sim b^{\prime}/\ell$ and,
more importantly, $\tau(f,g)\sim f^{\prime}g^{\prime},\tau(f,g,h)\sim
f^{\prime}g^{\prime}h^{\prime}$ (see Appendix A for the justification of the
latter).

###  VI.1 Nonlinear cascade

Let us consider $\pi_{\ell}$, the local in space cascade rate. By spatially
averaging $\pi_{\ell}$, see Eq. (18), over a region of size $L$ we get:

$\begin{split}\langle\pi_{\ell}\rangle_{L}=&-\langle\partial_{j}\bar{u}_{i}(\tau_{\ell}(u_{i},u_{j})-\tau_{\ell}(b_{i},b_{j}))+\varepsilon_{\alpha\beta\gamma}\bar{j}_{\alpha}\tau(u_{\beta},b_{\gamma})\rangle_{L}\\\
&+d_{i}\langle\varepsilon_{\alpha\beta\gamma}\bar{j}_{\alpha}\tau(j_{\beta},b_{\gamma})\rangle_{L}\end{split}$
(22)

It is straightforward to show that the order of
$\langle\pi_{\ell}(\bm{x})\rangle_{L}$, denoted
$[\langle\pi_{\ell}\rangle_{L}]$ is given by:

$[\langle\pi_{\ell}\rangle_{L}]=\frac{\bar{U}}{\ell}U^{\prime
2}+\frac{d_{i}}{\ell}\frac{\bar{U}}{\ell}U^{\prime 2}$

so that the two characteristic times of the cascade can be found by computing
$[\langle h_{\ell}\rangle_{L}]/[\langle\pi_{\ell}\rangle_{L}]$:

$t_{\text{NL}}^{\text{MHD}}=\frac{\ell}{\bar{U}},\qquad
t_{\text{NL}}^{Hall}=t_{\text{NL}}^{\text{MHD}}\frac{\ell}{d_{i}}$

These quantities describe how fast energy cascades across scale $\ell$. In
particular we recover the _eddy turn-over time_ at scale $\ell$ as the
characteristic time of nonlinear cascade.

###  VI.2 Large scale spatial transport

The terms governing the spatial transport of $\langle h_{\ell}\rangle_{L}$ due
to the large scale fields read (see Appendix B):

$\begin{split}&\nabla\cdot\langle\mathbfcal{J}^{\text{LS}}\rangle_{L}=\nabla\cdot\langle
h_{\ell}\bar{u}_{j}-\tau(u_{i},b_{i})\bar{b}_{j}\rangle_{L}\\\
&d_{i}\nabla\cdot\langle\bar{b}_{j}\tau(b_{i},j_{i})+\bar{j}_{i}\tau(b_{i},b_{j})-\bar{j}_{j}\tau(b_{i},b_{i})-\bar{b}_{i}\tau(j_{j},b_{i})\rangle_{L}\end{split}$
(23)

The divergence acts on quantities averaged over a region of size $L$.
Therefore the characteristic scale of the spatial derivatives is $1/L$ as all
fluctuations with smaller scales have been removed. With this in mind we can
proceed to derive the characteristic time by which each of these processes
extract or bring energy inside this region of size $L$:

$t_{\text{LS}}=\frac{L}{\bar{U}},\quad t_{A}=\frac{L}{B_{0}},\quad
t_{\text{LS}}^{Hall}=t_{\text{LS}}\frac{\ell}{d_{i}},\quad
t_{\text{Hall}}=t_{A}\frac{\ell}{d_{i}}$

Alongside the characteristic time of energy advection by the large scale flow
$t_{\text{LS}}$ we recognize the linear Alfvén time $t_{A}$ that represents
the energy transport due to propagating Alfvén waves at scale $L$.
Additionally, faster modes that have a time scale $t_{Hall}=t_{A}\ell/d_{i}$
exist as the scales approach $d_{i}$, which can be identified as whistler
waves with a dispersion relation $\omega\sim k_{\parallel}kB_{0}d_{i}$
(Sahraoui  _et al._ , 2007).

###  VI.3 Small scale transport

The effect of the small scale fluctuations in the spatial transport of
$\langle h_{\ell}\rangle_{L}$ is given by (see Appendix B):

$\begin{split}\nabla\cdot\langle\mathbfcal{J}^{\text{SS}}\rangle_{L}&=\nabla\cdot\left\langle\frac{\tau(b_{i},b_{i},u_{j})+\tau(u_{i},u_{i},u_{j})}{2}-\tau(u_{i},b_{i},b_{j})\right\rangle_{L}\\\
&+d_{i}\nabla\cdot\left\langle\tau(b_{i},b_{j},j_{i})-\tau(j_{j},b_{i},b_{i})\right\rangle_{L}\end{split}$
(24)

The corresponding characteristic time scales read:

$t_{\text{SS}}=\frac{L}{U^{\prime}},\quad
t_{\text{SS}}^{\text{Hall}}=t_{\text{SS}}\frac{\ell}{d_{i}}$

###  VI.4 Pressure term

The last term that we need to analyze is the spatial transport one,
$\nabla\cdot\langle\tau(u,P)\rangle_{L}$, involving the total plasma pressure.
We recall that in the incompressible HMHD model pressure is completely
determined by $\bm{u},\bm{b}$: taking the divergence of (2), pressure is found
by solving the Poisson equation:

$\nabla^{2}P=\partial_{i}\partial_{j}(u_{i}u_{j}-b_{i}b_{j})$

the equation for $P^{\prime}=P-\bar{P}$ is easily derived, whose ordering is
given by

$[P^{\prime}]=U^{\prime 2}+\bar{U}U^{\prime}$ (25)

which yields the following ordering of the transport term

$[\nabla\cdot\langle\tau(u,P)\rangle_{L}]=\frac{1}{L}P^{\prime}U^{\prime}\approx\frac{U^{\prime
2}U^{\prime}}{L}+\frac{\bar{U}U^{\prime 2}}{L}$ (26)

and the corresponding characteristic time scales

$t_{\text{SS}}=\frac{L}{U^{\prime}},\quad t_{\text{LS}}=\frac{L}{\bar{U}}$

both of which were already found in the analysis of the large and small scale
transport terms.

###  VI.5 The fastest dynamical process

We want to show that we can choose the region size $L$ so that the cross-scale
energy transfer in the region is faster than spatial transport across the
region borders. To do so we analyze the time-scales associated with different
processes. For the sake of simplicity the comparison is limited to the MHD
range where the Hall term can be neglected. Including the Hall term will not
modify any of the conclusion of this study based on the observation that the
Hall term modifies the time scales of the nonlinear cascade and transport by
the same (small scale) factor $\ell/d_{i}$. We consider the nonlinear time
scale $t_{\text{NL}}^{\text{MHD}}$ as a reference to obtain:

$\frac{t_{\text{NL}}^{\text{MHD}}}{t_{\text{LS}}}=\frac{\ell}{L};\quad\frac{t_{\text{NL}}^{\text{MHD}}}{t_{\text{SS}}}=\frac{\ell}{L}\frac{U^{\prime}}{\bar{U}};\quad\frac{t_{\text{NL}}^{\text{MHD}}}{t_{\text{A}}}=\frac{\ell}{L}\frac{B_{0}}{\bar{U}}\sim\frac{\ell}{L}\frac{B_{0}}{\delta\bar{b}}$
(27)

we notice that when $\ell/L\ll 1$ and $U^{\prime}<\bar{U}$ (the latter is
justified by a power-law decline of the fluctuations in MHD turbulence) the
first two conditions in relations (27) yield the ordering
$t_{\text{NL}}^{\text{MHD}}\ll t_{\text{LS}},t_{\text{SS}}$.

If we further assume high amplitude fluctuations with respect to the
background field, _i.e._ , $\delta b/B_{0}\gtrsim 1$, the condition $\ell/L\ll
1$ automatically implies the ordering $t_{\text{NL}}^{\text{MHD}}\ll
t_{\text{A}}$. This should not be confused with the critical balance ordering
(Goldreich and Sridhar, 1995) as the quantity $t_{A}$ is the time that an
Alfvén wave takes to propagate over a distance $L$: as $L$ can be varied at
will, $t_{A}$ can be made as large or as small as needed. In the general case
the stronger the mean magnetic field compared to the fluctuations, the larger
we should choose $L$ to maintain $t_{\text{NL}}^{\text{MHD}}\ll
t_{\text{A}}\rightarrow L\gg\ell B_{0}/\delta b$ for a given $\ell$.

To summarize, under the conditions $\ell/L\ll 1$ and $\ell/L\ll\delta
b/B_{0}$, the latter becoming redundant in strong turbulence,
$t_{\text{NL}}^{\text{MHD}}$ becomes the fastest time of the averaged
dynamics. The interpretation is the following: considering a spatial region of
characteristic size $L$, a certain amount of energy cascades across scale
$\ell$ inside this region on time scales much faster than what it takes for
the same amount of energy to be spatially transported across the region
surface by large and small scale processes including linear (Alfvén) waves.

For the aforementioned reason, in the averaged small scale equation the
spatial transport terms are slower than the nonlinear cascade. At a small
filtering scale $\ell$ we anticipate the averaged cascade to be balanced by
the averaged dissipation: as energy cascades in a given region it does not
have time to be spatially transported outside the region before it is
dissipated (note that this would imply some kind of _short-time_ stationarity
of the small scale, averaged, energy $\langle h_{\ell}\rangle_{L}$).
Therefore, in equation (21) we can expect the two fast processes to match:

$\langle\pi_{\ell}(\bm{x})\rangle_{L}\approx-\langle\tau_{\ell}(u_{i},d_{\nu,i})+\tau_{\ell}(b_{i},d_{\eta,i})\rangle_{L}\approx-\langle\bm{u}^{\prime}\cdot\bm{d}^{\prime}_{\nu}+\bm{b}^{\prime}\cdot\bm{d}^{\prime}_{\eta}\rangle_{L}$
(28)

This equation shows, as stated at the beginning of this section, that the
(quasi-)local cascade matches the contribution to the (quasi-)local
dissipation coming from scales $\lesssim\ell$.
The result is general and does not require the existence of an inertial range
as it does not involve the total dissipation, but only a part of it. However,
if we assume the dissipation rate due to the large scales (_i.e._
$\gtrsim\ell$) to be negligible, it can be further refined to obtain

$\langle\pi_{\ell}(\bm{x})\rangle_{L}\approx-\langle\bm{u}^{\prime}\cdot\bm{d}^{\prime}_{\nu}+\bm{b}^{\prime}\cdot\bm{d}^{\prime}_{\eta}\rangle_{L}\approx-\langle\bm{u}\cdot\bm{d}_{\nu}+\bm{b}\cdot\bm{d}_{\eta}\rangle_{L}$
(29)

which shows the approximate (quasi-)local balance between the cascade rate and
total dissipation. This result will now be tested numerically.

###  VI.6 Numerical Validation

Figure 6: 2D plots of the cascade rate $\pi_{\ell}$ across scale
$\ell=0.3d_{i}$ (a) and dissipation (b) based on the data from run I. Figures
(c,d) are obtained by computing for each point $\bm{x}$ the average on a cubic
region of size $L=1.4d_{i}$ centered on $\bm{x}$. Figure 7: Ratio
$R(\bm{x},L)=-\langle\pi_{\ell}\rangle/\langle\bm{u}\cdot\bm{d}_{\nu}+\bm{b}\cdot\bm{d}_{\eta}\rangle$
between the local cascade rate and the local dissipation rate averaged on a
cubic region of size $L$ centered around a given point $\bm{x}$. The results
for three different points on the plane $z=25\times 2\pi d_{i}$ in run I (a-c)
and run II (d-f) are shown. The shaded region individuates the interval
$0.5<R<2$.

A first glimpse on the local (in space) behavior of the cascade and
dissipation at the filtering scale $\ell=0.3d_{i}$ is given in Fig. (6(a) and
(c)) based on the data from run I. Overall the two quantities exhibit similar
patterns, indicating an approximate balance between them, although the extrema
of the cascade rate can be larger locally by a factor $\sim 2$. We also
observe that cascade rate $\pi_{\ell}(x)$ presents both positive and negative
values. This is a clear sign that the nonlinear interactions work in both
directions, bringing energy from large to small scales and _vice-versa_ ,
albeit at this small filtering scale $\ell=0.3d_{i}$ most of the energy is
going towards the small scales and $\pi$ is mostly positive. Nevertheless,
when integrated over a cubic region of size $L=1.4d_{i}$ centered on each
point $\bm{x}$ (panel 6(b)), we recover a nearly positive flux, i.e.
$\langle\pi_{\ell}(\bm{x})\rangle_{L}>0$, a sign that on average the quasi-
local turbulent cascade is direct since energy is carried towards the small
scales. The same observations can be made about dissipation (panels 6(c) and
(d)), which is indeed positive definite only when no energy leakage through
the boundaries is assumed.

A thorough analysis can be done by pinpointing individual locations in space
$\bm{x}_{0}$ that correspond to intense local transfers on which we can test
the balance between the energy cascade and dissipation given by relation (29).
The chosen locations are indicated in Fig.7 by labels (a)-(c) for run I and
(d)-(f) for run II. Varying $L$, the size of the region over which the
integration is performed, we are able to gauge $R(\bm{x}_{0},L)$, the ratio
between $\langle\pi_{\ell}\rangle_{L}$ and
$-\langle\bm{u}\cdot\bm{d}_{\nu}+\bm{b}\cdot\bm{d}_{\eta}\rangle_{L}$. The
results in Fig.7 show that the ratio between the two terms is close to 1 not
only when $L=50d_{i}=L_{\text{box}}$ but it remains constant for a very large
range of scales all the way down to $\ell\sim 8d_{i}$. This is consistent with
relation (29) which holds when $\ell/L\ll 1$. It is however remarkable to
notice that even for smaller values $\ell\sim L$ the ratio between the two
terms remains comparable with 1. In particular, Fig.7 shows that $0.5\leq
R\leq 2$ at almost all scales for all the three points under study in run I.
The same study for run II shows an overall good agreement between the local
cascade and dissipation magnitude. However, the matching is less good than in
run I, as the ratio $R(\bm{x}_{0},L)$ departs from the reference value $1$
(i.e., perfect balance) at larger size box $L$ for run II than run I. This is
likely to be caused by the presence of the mean magnetic field $B_{0}=2$ in
run II, which introduces (large scale) Alfvén waves that spatially transport
energy on comparable time scales than those of the nonlinear cascade and
dissipation as discussed in Section VI.5. Said differently, for run I the time
scale ratio $t_{\text{NL}}^{\text{MHD}}/t_{\text{A}}\sim\ell
B_{0}/(L\langle\delta\bar{b}\rangle_{L})\sim 0$ since $B_{0}=0$ (i.e.,
$t_{A}\to\infty$) regardless of how large is the ratio $\ell/L$. This could
explain the very good matching between the cascade and dissipation rates in
run I even for very small integrated regions, i.e., $\ell/L\sim 1$ as observed
in Fig.7 (a)-(c). In run II, with $B_{0}=2$, the linear transport time scale
$t_{A}$ becomes finite and can be comparable to $t_{\text{NL}}^{\text{MHD}}$
even for small ratios of $\ell/L$. In the limit case of a very strong mean
field, $\delta b/B_{0}\ll 1$ the linear Alfvén time will become faster than
the nonlinear time at all scales, and relation (29) would no longer hold.

left

Figure 8: The normalized distributions of the ratio
$-\langle\pi_{\ell}(\mathbf{x})\rangle_{L}/\langle\mathbf{u}\cdot\mathbf{d}_{\nu}+\mathbf{b}\cdot\mathbf{d}_{\eta}\rangle_{L}$
(c-e) computed at scale $\ell=0.3d_{i}$ and plane $z=25\times 2\pi d_{i}$ in
the regions of intense cross-scale transfer marked by crosses for runs I (a)
and II (b). Three different sizes $L$ of the integrated regions are
considered. Filled, colored values correspond to run I ($B_{0}=0$), while
black, hatched values to run II $B_{0}=2$.

To obtain a more complete picture of the balance between the quasi-local
cascade and dissipation we performed the same study over a larger sample of
locations shown in Fig 8. For both runs we consider some of the local maxima
of $\pi_{\ell}(\bm{x})$, $\ell=0.3d_{i}$ in the plane $z=25\times 2\pi d_{i}$,
for each of these points the ratio $R(\bm{x}_{0},L)$ is computed and its
histogram is shown for each size $L$ of the integrated box. As expected, for
both runs we observe that as $L$ decreases the width of the histogram (i.e.,
dispersion) increases, and its mean tends to shift to larger values. This
indicates that, on average, the local cascade tends to slightly dominate the
local dissipation. We also observe that dispersion around the reference value
$R(\bm{x}_{0},L)\sim 1$ (perfect balance) is more prominent in run II than run
I, which confirms the role of the background magnetic field $B_{0}$ in
shortening the linear transport time, thus competing with the pair cascade-
dissipation discussed above.

##  VII Conclusions

In this paper we have derived the _large_ and _small_ scales CG equations for
incompressible HMHD and showed that the quantity $\pi_{\ell}(\bm{x})$ is the
local (at position $\bm{x}$) energy transfer across the scale $\ell$. This
quantity, when averaged over the simulation domain allows us to recover
results consistent with the so-called “third-order” laws. It represents an
improvement over the current state-of-the-art in this field as it does not
rely on the assumptions of ergodicity or the existence of the inertial range.
The model is further generalized to account for anisotropic cascade in the
presence of a background magnetic field. The major strength of this new theory
is to provide a means of estimating local energy dissipation in turbulent
plasmas. This is shown to be achieved at a cost of fairly loose assumptions,
the main ones being the scale separation between the filtering scale $\ell$
and the integration region size $L$, moderate-to-large amplitude of the
turbulent fluctuations with respect to the background field $B_{0}$. This last
condition is expected in strong turbulence. These two assumptions are required
to minimize the role of the spatial transport across the region of
integration, enforcing thus the balance between local nonlinear cascade and
dissipation.

The theory was tested successfully on two simulations featuring different
intensities of the background field $B_{0}$. In agreement with the time scale
estimate of the various processes, a very good local balance between cascade
and dissipation is found in run with $B_{0}=0$ even when $L\sim\ell$, while a
moderate imbalance between the two is observed at larger values $L$ for run
with $B_{0}=2$. We conjecture that this behaviour is due to the role of Alfvén
waves that spatially transport energy across the integration region on time
scales comparable to those of the cascade and dissipation.

An immediate application of this novel theory is to estimate energy
dissipation in localized structures (e.g., reconnecting current sheets)
frequently reported in numerical simulations (Arró  _et al._ , 2020) and
spacecraft observations in the near-Earth space (Chasapis  _et al._ , 2015;
Huang  _et al._ , 2018). Indeed, even if the theory is based on a fluid (HMHD)
model, the fact that $\pi_{\ell}(\bm{x})$ is shown to reflect local
dissipation (within the aforementioned assumptions) regardless of the explicit
form of the dissipation operators makes it particularly relevant to
collisionless plasmas where the damping of the fluid and electromagnetic
fluctuations is believed to originate from kinetic processes that would show
up in the moment equation as complex damping operators
$\bm{d}_{\nu},\bm{d}_{\eta}$.

The present theory can be extended in various directions. For instance to the
two-fluid model where electron inertia can be accounted for (particularly
relevant in reconnection studies and small scale turbulence, see for instance
(Faganello  _et al._ , 2009) and references therein), or by incorporating
compressible effects (Eyink, 2018; Eyink and Drivas, 2018). Regarding its
implementation on simulations and spacecraft data, an integration over an
anisotropic domain $L_{\perp}\ll L_{\parallel}$ would be more appropriate to
capture the difference in time scale of the various processes along the two
directions. Lastly, in numerical simulation an integration over the actual 3D
volume of the structures (e.g., current sheets, plasmoids) would be a
significant improvement with respect to the cubic box adopted here.

## Acknowledgement

The authors thank Prof. P. Mininni and the developers of the GHOST code for
providing the code used to run the simulations presented in this work.
This work was granted access to the HPC resources of CINES under allocation
2021 A0090407714 made by GENCI. It is supported by the CNRS/CONICET
Laboratoire International Associé (LIA) MAGNETO.

##  Appendix A On the generalized central moments

When studying turbulence using the statistical approach the original field $u$
is decomposed into its (ensemble) average part $\langle\bm{u}\rangle$ and a
fluctuations $\bm{u}^{\prime}$ whose sum gives the original quantity
$\bm{u}=\langle\bm{u}\rangle+\bm{u}^{\prime}$. In this framework we can derive
relations concerning the central moments of the fluctuations $\langle
u^{\prime}_{i}u^{\prime}_{j}\rangle,\langle
u^{\prime}_{i}u^{\prime}_{j}u^{\prime}_{k}\rangle$,…. It can be shown (Frisch,
1995; Monin and Yaglom, 2013) that:

$\begin{split}&\langle u^{\prime}_{i}u^{\prime}_{j}\rangle=\langle
u_{i}u_{j}\rangle-\langle u_{i}\rangle\langle u_{j}\rangle\\\ &\langle
u^{\prime}_{i}u^{\prime}_{j}u^{\prime}_{k}\rangle=\langle
u_{i}u_{j}u_{k}\rangle-\langle u_{i}\rangle\langle
u_{j}^{\prime}u_{k}^{\prime}\rangle-\langle u_{j}\rangle\langle
u_{i}^{\prime}u_{k}^{\prime}\rangle\\\ &\qquad\qquad\qquad-\langle
u_{k}\rangle\langle u_{i}^{\prime}u_{j}^{\prime}\rangle-\langle
u_{i}\rangle\langle u_{j}\rangle\langle u_{k}\rangle\\\ &\langle
u^{\prime}_{i}u^{\prime}_{j}u^{\prime}_{k}u^{\prime}_{l}\rangle=\dots\end{split}$
(30)

When using the filtering approach (CG), we decompose the field into its large
scale $\bar{\bm{u}}$ and “unresolved” $\bm{u}-\bar{\bm{u}}$ contributions,
replacing the ensemble average $\langle\cdot\rangle$ with the CG operation. To
recover a set of relations equivalent to (30) the seminal work of M.Germano
Germano (1992) introduces the _generalized central moments_
$\tau(u_{i},u_{j}),\tau(u_{i},u_{j},u_{k}),\dots$, formally defined as:

$\begin{split}&\tau(u_{i},u_{j})=\overline{u_{i}u_{j}}-\bar{u}_{i}\bar{u}_{j}\\\
&\tau(u_{i},u_{j},u_{k})=\overline{u_{i}u_{j}u_{k}}-\bar{u}_{i}\tau(u_{j},u_{k})-\bar{u}_{j}\tau(u_{i},u_{k})\\\
&\qquad\qquad\qquad-\bar{u}_{k}\tau(u_{i},u_{j})-\bar{u}_{i}\bar{u}_{j}\bar{u}_{k}\\\
&\tau(u_{i},u_{j},u_{k},u_{l})=\dots\end{split}$ (31)

which are employed in this work and allow us to recover simple filtered
equations. To aid in the physical interpretation some useful relations can be
derived to link the generalized central moments to the field fluctuations
(Vreman  _et al._ , 1994). We start with:

$\begin{split}&\tau(f,g)=\overline{fg}-\bar{f}\bar{g}=\left(\overline{fg}-\bar{f}\bar{g}\right)+\left(-\bar{f}\bar{g}+\bar{f}\bar{g}\right)\\\
=&\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})f(\bm{x}+\bm{r})g(\bm{x}+\bm{r})-\bar{f}(\bm{x})\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})g(\bm{x}+\bm{r})+\\\ -&\bar{g}(\bm{x})\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})f(\bm{x}+\bm{r})+\bar{f}(\bm{x})\bar{g}(\bm{x})\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})\\\ =&\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})\left[f(\bm{x}+\bm{r})-\bar{f}(\bm{x})\right]g(\bm{x}+\bm{r})-\bar{g}(\bm{x})\times\\\
\times&\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})\left[f(\bm{x}+\bm{r})-\bar{f}(\bm{x})\right]=\\\
=&\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})\left[f(\bm{x}+\bm{r})-\bar{f}(\bm{x})\right]\left[g(\bm{x}+\bm{r})-\bar{g}(\bm{x})\right]\end{split}$
(32)

Performing a second order Taylor expansion:

$f(\bm{x}+\bm{r})=f(\bm{x})+r_{j}\frac{\partial f}{\partial
x_{j}}+\frac{r_{i}r_{j}}{2}\frac{\partial^{2}f}{\partial x_{i}\partial
x_{j}}+o(|\bm{r}|^{2})$

and denoting the turbulent, subscale fluctuations with $f^{\prime}=f-\bar{f}$,
we have:

$\begin{split}\tau(f,g)\sim&f^{\prime}(\bm{x})g^{\prime}(\bm{x})+\left[g^{\prime}\nabla
f+f^{\prime}\nabla g\right]_{\bm{x}}\cdot\int d\bm{r}\bm{r}G_{\ell}(\bm{r})\\\
+&\left[\frac{\partial f}{\partial x_{i}}\frac{\partial g}{\partial
x_{j}}+\frac{\partial f}{\partial x_{j}}\frac{\partial g}{\partial
x_{i}}\right]\int d\bm{r}G_{\ell}(\bm{r})r_{i}r_{j}\\\
+&\left[g^{\prime}\frac{\partial^{2}f}{2\partial x_{i}\partial
x_{j}}+f^{\prime}\frac{\partial^{2}g}{2\partial x_{i}\partial
x_{j}}\right]\int d\bm{r}G_{\ell}(\bm{r})r_{i}r_{j}\end{split}$

Since $\int d\bm{r}G_{\ell}(\bm{r})\bm{r}=0$, and the second order centered
moments matrix $\int d\bm{r}G_{\ell}(\bm{r})r_{i}r_{j}\sim\ell^{2}\delta_{ij}$
(with the Gaussian filter used), the relation $\tau(f,g)\sim
f^{\prime}g^{\prime}$ is valid up to the second order in $\ell/\Lambda$. Here
$\Lambda^{-1}\sim\min\\{\nabla f/f,\nabla g/g\\}$ is the characteristic scale
of the field gradients. The parameter $\ell/\Lambda$ is however not guaranteed
to be small as the fields have a Fourier spectrum that extends all the way
down to the dissipation scales. Therefore, an even more precise relation can
be obtained by performing another CG operation, or, in general, by averaging
over a spatial region of size $L\gg\ell$:

$\langle\tau(f,g)\rangle_{L}=\langle
f^{\prime}g^{\prime}\rangle_{L}+\mathcal{O}\left((l/\Lambda^{\prime})^{2}\right)$
(33)

where $\Lambda^{\prime}$, the gradient characteristic scale of the averaged
fields, is of order $L$ since the spatial average removes fluctuations at
smaller scales.

A similar result can be derived for the third order generalized moment defined
as:

$\tau(f,g,h)=\overline{fgh}-\bar{f}\bar{g}\bar{h}-\bar{f}\tau(g,h)-\bar{g}\tau(f,h)-\bar{h}\tau(f,g)$
(34)

a direct computation, denoting $\bm{x}^{\prime}=\bm{x}+\bm{r}$ shows:

$\begin{split}&\int
d\bm{r}G_{\ell}(\bm{x}^{\prime})\left[f(\bm{x}^{\prime})-\bar{f}(\bm{x})\right]\left[g(\bm{x}^{\prime})-\bar{g}(\bm{x})\right]\left[h(\bm{x}^{\prime})-\bar{h}(\bm{x})\right]=\\\
=&\overline{fgh}-\bar{f}\bar{g}\bar{h}+\,[-\bar{f}\overline{gh}+\bar{f}\bar{g}\bar{h}]+[-\bar{g}\overline{fh}+\bar{g}\bar{f}\bar{h}]+[-\overline{fg}\bar{h}-\bar{f}\bar{g}\bar{h}]=\\\
=&\overline{fgh}-\bar{f}\bar{g}\bar{h}\qquad-\bar{f}\tau(g,h)\qquad-\bar{g}\tau(f,g)\qquad-\bar{h}\tau(f,g)=\\\
&=\tau(f,g,h)\end{split}$ (35)

and, as for the second order generalized central moment, a Taylor expansion
followed by a spatial average leads to

$\langle\tau(f,g,h)\rangle_{L}=\langle
f^{\prime}g^{\prime}h^{\prime}\rangle_{L}+\mathcal{O}((\ell/L)^{2})$ (36)

which is precisely the relation used in the main text.

##  Appendix B DERIVATION OF THE SMALL SCALE EQUATIONS

The energy contained at sub-filtering scales $<\ell$ is defined as
$h_{\ell}=\left[\tau_{\ell}(u_{i},u_{i})+\tau_{\ell}(b_{i},b_{i})\right]/2$
and is a positive quantity at every point in space in virtue of relation (32),
_e.g._ $\tau(u_{i},u_{i})=\int
d\bm{r}G_{\ell}(\bm{x}+\bm{r})|u_{i}(\bm{x}+\bm{r})-\bar{u}_{i}(\bm{x})|^{2}>0$.
Furthermore when integrating over the spatial domain we recover the total sub-
scale energy $\int d^{3}\bm{x}h_{\ell}(\bm{x})=\int
d^{3}x\left[(|\bm{u}|^{2}-|\bar{\bm{u}}_{\ell}|^{2})+(|\bm{b}|^{2}-|\bar{\bm{b}}_{\ell}|^{2})\right]/2$,
(see Aluie (2017)).
To compute the evolution of $h_{\ell}$ we filter on a scale $\ell$ the energy
equation and we subtract the equation for the large scale energy
$(|\bar{\bm{u}}|^{2}+|\bar{\bm{b}}|^{2})/2$. We readily obtain for kinetic and
magnetic energy densities:

$\begin{split}\partial_{t}\frac{\tau(u_{i},u_{i})}{2}=&-\partial_{j}\left[\frac{\tau(u_{i},u_{i})}{2}\bar{u}_{j}+\frac{\tau(u_{i},u_{i},u_{j})}{2}+\tau(u_{j},P)-\tau(u_{i},b_{i},b_{j})-\tau(u_{i},b_{j})\bar{b}_{i}-\tau(u_{i},b_{i})\bar{b}_{j}\right]\\\
&+\pi^{u}+(\bar{b_{i}}\bar{b}_{j}\partial\bar{u}_{j}-\overline{b_{i}b_{j}\partial_{j}u_{i}})+\tau(u_{i},d_{\nu,i})+\tau(u_{i},f)\end{split}$
(37)

$\begin{split}\partial_{t}\frac{\tau(b_{i},b_{i})}{2}=&-\partial_{j}\left[\frac{\tau(b_{i},b_{i})}{2}\bar{u}_{j}+\frac{\tau(b_{i},b_{i},u_{j})}{2}+\bar{b}_{i}\tau(u_{i},b_{j})\right]\\\
&-d_{i}\partial_{j}\left[\tau(b_{i},b_{j},j_{i})-\tau(j_{j},b_{i},b_{i})+\bar{b}_{j}\tau(b_{i},j_{i})+\bar{j}_{i}\tau(b_{i},b_{j})-\bar{j}_{j}\tau(b_{i},b_{i})-\bar{b}_{i}\tau(j_{j},b_{i})\right]\\\
&+\pi^{b}+(\overline{b_{i}b_{j}\partial_{j}u_{i}}-\bar{b_{i}}\bar{b}_{j}\partial\bar{u}_{j})+\tau(b_{i},d_{\eta,i})\end{split}$
(38)

And taking the sum of the two equations we obtain

$\begin{split}\partial_{t}\frac{\tau(u_{i},u_{i})+\tau(b_{i},b_{i})}{2}=&-\partial_{j}\left[\frac{\tau(u_{i},u_{i})+\tau(b_{i},b_{i})}{2}\bar{u}_{j}+\tau(u_{j},P)-\tau(u_{i},b_{i})\bar{b}_{j}\right]-\partial_{j}\left[\frac{\tau(b_{i},b_{i},u_{j})+\tau(u_{i},u_{i},u_{j})}{2}-\tau(u_{i},b_{i},b_{j})\right]\\\
&-d_{i}\partial_{j}\left[\tau(b_{i},b_{j},j_{i})-\tau(j_{j},b_{i},b_{i})+\bar{b}_{j}\tau(b_{i},j_{i})+\bar{j}_{i}\tau(b_{i},b_{j})-\bar{j}_{j}\tau(b_{i},b_{i})-\bar{b}_{i}\tau(j_{j},b_{i})\right]\\\
&+\pi^{u}+\pi^{b}+\tau(b_{i},d_{\eta,i})+\tau(u_{i},d_{\nu,i})+\tau(u_{i},f)\end{split}$
(39)

We can group the contribution of the spatial transport of the small scale
energy $h_{\ell}=\left[\tau(u_{i},u_{i})+\tau(b_{i},b_{i})\right]/2$ due to
the large scale fields:

$\begin{split}&\nabla\cdot\mathbfcal{J}^{\text{LS}}=\nabla\cdot\left[h_{\ell}\bar{u}_{j}-\tau(u_{i},b_{i})\bar{b}_{j}\right]+\\\
+&d_{i}\nabla\cdot\left[\bar{b}_{j}\tau(b_{i},j_{i})+\bar{j}_{i}\tau(b_{i},b_{j})-\bar{j}_{j}\tau(b_{i},b_{i})-\bar{b}_{i}\tau(j_{j},b_{i})\right]\end{split}$

Similarly, the spatial transport of $h_{\ell}$ due to the subscale fields is
given by:

$\begin{split}&\nabla\cdot\mathbfcal{J}^{\text{SS}}=\nabla\cdot\left[\frac{\tau(b_{i},b_{i},u_{j})+\tau(u_{i},u_{i},u_{j})}{2}-\tau(u_{i},b_{i},b_{j})\right]+\\\
+&d_{i}\nabla\cdot\left[\tau(b_{i},b_{j},j_{i})-\tau(j_{j},b_{i},b_{i})\right]\end{split}$

and lastly the spatial transport due to pressure interactions:

$\nabla\cdot\mathbfcal{J}^{\text{p}}=\nabla\cdot\tau(u_{j},P)$

so that equation (39) can be written in compact form as:

$\begin{split}\partial_{t}h_{\ell}=&-\nabla\cdot\left[\mathbfcal{J}^{LS}+\mathbfcal{J}^{SS}+\mathbfcal{J}^{p}\right]\\\
+&\pi_{\ell}+\tau(b_{i},d_{\eta,i})+\tau(u_{i},d_{\nu,i})+\tau(u_{i},f)\end{split}$
(40)

## References

  * Batchelor (1953) G. Batchelor,  _The Theory of Homogeneous Turbulence_ (Cambridge University Press, 1953).
  * Frisch (1995) U. Frisch,  _Turbulence: The Legacy of A.N. Kolmogorov_ (CAMBRIDGE UNIVERSITY PRESS, 1995).
  * Monin and Yaglom (2013) A. S. Monin and A. M. Yaglom,  _Statistical Fluid Mechanics, Volume II: Mechanics of Turbulence_ (Courier Corporation, 2013) google-Books-ID: 6xPEAgAAQBAJ.
  * Bruno and Carbone (2013) R. Bruno and V. Carbone, Living Rev. Solar Phys. 10 (2013), 10.12942/lrsp-2013-2.
  * Matthaeus and Velli (2011) W. H. Matthaeus and M. Velli, Space Sci Rev 160, 145 (2011).
  * Goldstein  _et al._ (2015) M. L. Goldstein, R. T. Wicks, S. Perri, and F. Sahraoui, Philosophical Transactions. Series A, Mathematical, Physical, and Engineering Sciences 373, 20140147 (2015).
  * Sahraoui  _et al._ (2020) F. Sahraoui, L. Hadid, and S. Huang, Rev. Mod. Plasma Phys. 4, 4 (2020).
  * Politano and Pouquet (1998) H. Politano and A. Pouquet, Geophys. Res. Lett. 25, 273 (1998).
  * Banerjee and Galtier (2013) S. Banerjee and S. Galtier, Phys. Rev. E 87, 013019 (2013).
  * Andrés and Sahraoui (2017) N. Andrés and F. Sahraoui, Phys. Rev. E 96, 053205 (2017).
  * Andrés  _et al._ (2018) N. Andrés, F. Sahraoui, S. Galtier, L. Z. Hadid, P. Dmitruk, and P. D. Mininni, J. Plasma Phys. 84, 21 (2018).
  * Hellinger  _et al._ (2018) P. Hellinger, A. Verdini, S. Landi, L. Franci, and L. Matteini, The Astrophysical Journal 857, L19 (2018), publisher: American Astronomical Society.
  * Ferrand  _et al._ (2019) R. Ferrand, S. Galtier, F. Sahraoui, R. Meyrand, N. Andrés, and S. Banerjee, ApJ 881, 50 (2019).
  * Ferrand  _et al._ (2021a) R. Ferrand, S. Galtier, and F. Sahraoui, Journal of Plasma Physics 87, 905870220 (2021a), publisher: Cambridge University Press.
  * Smith  _et al._ (2006) C. W. Smith, K. Hamilton, B. J. Vasquez, and R. J. Leamon, ApJ 645, L85 (2006), publisher: IOP Publishing.
  * Sorriso-Valvo  _et al._ (2007) L. Sorriso-Valvo, R. Marino, V. Carbone, A. Noullez, F. Lepreti, P. Veltri, R. Bruno, B. Bavassano, and E. Pietropaolo, Phys. Rev. Lett. 99, 115001 (2007), publisher: American Physical Society.
  * Marino  _et al._ (2008) R. Marino, L. Sorriso-Valvo, V. Carbone, A. Noullez, R. Bruno, and B. Bavassano, ApJ 677, L71 (2008), publisher: IOP Publishing.
  * Stawarz  _et al._ (2010) J. E. Stawarz, C. W. Smith, B. J. Vasquez, M. A. Forman, and B. T. MacBride, ApJ 713, 920 (2010).
  * Osman  _et al._ (2011) K. T. Osman, M. Wan, W. H. Matthaeus, J. M. Weygand, and S. Dasso, Phys. Rev. Lett. 107, 165001 (2011).
  * MacBride  _et al._ (2008) B. T. MacBride, C. W. Smith, and M. A. Forman, ApJ 679, 1644 (2008), publisher: IOP Publishing.
  * Coburn  _et al._ (2015) J. T. Coburn, M. A. Forman, C. W. Smith, B. J. Vasquez, and J. E. Stawarz, Phil. Trans. R. Soc. A. 373, 20140150 (2015).
  * Banerjee and Galtier (2016) S. Banerjee and S. Galtier,  50, 015501 (2016), publisher: IOP Publishing.
  * Hadid  _et al._ (2017) L. Z. Hadid, F. Sahraoui, and S. Galtier, ApJ 838, 9 (2017), arXiv: 1612.02150.
  * Hadid  _et al._ (2018) L. Hadid, F. Sahraoui, S. Galtier, and S. Huang, Phys. Rev. Lett. 120, 055102 (2018).
  * Andrés  _et al._ (2020) N. Andrés, N. Romanelli, L. Z. Hadid, F. Sahraoui, G. DiBraccio, and J. Halekas, The Astrophysical Journal 902, 134 (2020), publisher: American Astronomical Society.
  * Sorriso-Valvo  _et al._ (2019) L. Sorriso-Valvo, F. Catapano, A. Retinò, O. Le Contel, D. Perrone, O. W. Roberts, J. T. Coburn, V. Panebianco, F. Valentini, S. Perri, A. Greco, F. Malara, V. Carbone, P. Veltri, O. Pezzi, F. Fraternale, F. Di Mare, R. Marino, B. Giles, T. E. Moore, C. T. Russell, R. B. Torbert, J. L. Burch, and Y. V. Khotyaintsev, Phys. Rev. Lett. 122, 035102 (2019), publisher: American Physical Society.
  * Greco  _et al._ (2008) A. Greco, P. Chuychai, W. H. Matthaeus, S. Servidio, and P. Dmitruk, Geophysical Research Letters 35 (2008), 10.1029/2008GL035454.
  * Chasapis  _et al._ (2015) A. Chasapis, A. Retinò, F. Sahraoui, A. Vaivads, Y. V. Khotyaintsev, D. Sundkvist, A. Greco, L. Sorriso-Valvo, and P. Canu, The Astrophysical Journal 804, L1 (2015), publisher: American Astronomical Society.
  * Galtier (2018) S. Galtier, J. Phys. A: Math. Theor. 51, 205501 (2018), publisher: IOP Publishing.
  * Dubrulle (2019) B. Dubrulle, J. Fluid Mech. 867, P1 (2019).
  * Kuzzay  _et al._ (2019) D. Kuzzay, O. Alexandrova, and L. Matteini, Phys. Rev. E 99, 053202 (2019).
  * David and Galtier (2021) V. David and S. Galtier, Physical Review E 103, 063217 (2021), publisher: American Physical Society.
  * David  _et al._ (2022) V. David, S. Galtier, F. Sahraoui, and L. Z. Hadid, arXiv:2201.02377 [astro-ph, physics:physics] (2022), arXiv: 2201.02377.
  * Meneveau and Katz (2000) C. Meneveau and J. Katz, Annual Review of Fluid Mechanics 32, 1 (2000), aDS Bibcode: 2000AnRFM..32….1M.
  * Eyink and Aluie (2009) G. L. Eyink and H. Aluie, Physics of Fluids 21, 115107 (2009), publisher: American Institute of Physics.
  * Simon and Sahraoui (2021) P. Simon and F. Sahraoui, ApJ , 9 (2021).
  * Aluie (2017) H. Aluie, New Journal of Physics 19, 025008 (2017), publisher: IOP Publishing.
  * Eyink (2018) G. L. Eyink, Phys. Rev. X 8, 041020 (2018).
  * Ferrand  _et al._ (2021b) R. Ferrand, F. Sahraoui, D. Laveder, T. Passot, P. L. Sulem, and S. Galtier, arXiv:2109.03123 [astro-ph, physics:physics] (2021b), arXiv: 2109.03123.
  * Ferrand  _et al._ (2022) R. Ferrand, F. Sahraoui, S. Galtier, N. Andrés, P. Mininni, and P. Dmitruk, arXiv:2201.10819 [physics] (2022), arXiv: 2201.10819.
  * Gómez  _et al._ (2005) D. O. Gómez, P. D. Mininni, and P. Dmitruk, Advances in Space Research Fundamentals of Space Environment Science, 35, 899 (2005).
  * Mininni  _et al._ (2011) P. D. Mininni, D. Rosenberg, R. Reddy, and A. Pouquet, Parallel Computing 37, 316 (2011).
  * Banerjee and Galtier (2017) S. Banerjee and S. Galtier, J. Phys. A: Math. Theor. 50, 015501 (2017).
  * Shebalin  _et al._ (1983) J. V. Shebalin, W. H. Matthaeus, and D. Montgomery, Journal of Plasma Physics 29, 525 (1983), publisher: Cambridge University Press.
  * Matthaeus  _et al._ (1984) W. H. Matthaeus, J. J. Ambrosiano, and M. L. Goldstein, Phys. Rev. Lett. 53, 1449 (1984).
  * Goldreich and Sridhar (1995) P. Goldreich and S. Sridhar, ApJ 438, 763 (1995).
  * Schekochihin  _et al._ (2009) A. A. Schekochihin, S. C. Cowley, W. Dorland, G. W. Hammett, G. G. Howes, E. Quataert, and T. Tatsuno, ApJS 182, 310 (2009).
  * Cho and Vishniac (2000) J. Cho and E. T. Vishniac, The Astrophysical Journal 539, 273 (2000), publisher: IOP Publishing.
  * Sahraoui  _et al._ (2006) F. Sahraoui, G. Belmont, L. Rezeau, N. Cornilleau-Wehrlin, J. L. Pinçon, and A. Balogh, Phys. Rev. Lett. 96, 075002 (2006).
  * Meyrand and Galtier (2013) R. Meyrand and S. Galtier, Phys. Rev. Lett. 111, 264501 (2013).
  * Sahraoui  _et al._ (2010) F. Sahraoui, M. L. Goldstein, G. Belmont, P. Canu, and L. Rezeau, Phys. Rev. Lett. 105, 131101 (2010).
  * Chen  _et al._ (2010) C. H. K. Chen, T. S. Horbury, A. A. Schekochihin, R. T. Wicks, O. Alexandrova, and J. Mitchell, Physical Review Letters 104, 255002 (2010), publisher: American Physical Society.
  * Lamriben  _et al._ (2011) C. Lamriben, P.-P. Cortet, and F. Moisy, Physical Review Letters 107, 024503 (2011), publisher: American Physical Society.
  * Sahraoui  _et al._ (2007) F. Sahraoui, S. Galtier, and G. Belmont, J. Plasma Phys. 73, 723 (2007).
  * Arró  _et al._ (2020) G. Arró, F. Califano, and G. Lapenta, Astronomy & Astrophysics 642, A45 (2020), publisher: EDP Sciences.
  * Huang  _et al._ (2018) Z. Huang, C. Mou, H. Fu, L. Deng, B. Li, and L. Xia, The Astrophysical Journal 853, L26 (2018), publisher: American Astronomical Society.
  * Faganello  _et al._ (2009) M. Faganello, F. Califano, and F. Pegoraro, New Journal of Physics 11, 063008 (2009), publisher: IOP Publishing.
  * Eyink and Drivas (2018) G. L. Eyink and T. D. Drivas, Phys. Rev. X 8, 011022 (2018).
  * Germano (1992) M. Germano, Journal of Fluid Mechanics 238, 325–336 (1992).
  * Vreman  _et al._ (1994) B. Vreman, B. Geurts, and H. Kuerten, Journal of Fluid Mechanics 278, 351 (1994), publisher: Cambridge University Press.