<figure><img src="image_1.png" /></figure>

# Fluid Energy Cascade Rate and Kinetic Damping: New Insight from 3D Landau-fluid Simulations

R. Ferrand<fnref n="1" />, F. Sahraoui<fnref n="1" />, D. Laveder<fnref n="2" />, T. Passot<fnref n="2" />, P. L. Sulem<fnref n="2" />, and S. Galtier<fnref n="1" />,³

Received 2021 September 7; revised 2021 September 28; accepted 2021 September 29; published 2021 December 15

## Abstract

Using an exact law for incompressible Hall magnetohydrodynamics (HMHD) turbulence, the energy cascade rate is computed from three-dimensional HMHD-CGL (biadiabatic ions and isothermal electrons) and Landau-fluid numerical simulations that feature different intensities of Landau damping over a broad range of wavenumbers, typically $0.05 \lesssim k_\perp d_i \lesssim 100$. Using three sets of cross-scale simulations where turbulence is initiated at large, medium, and small scales, the ability of the fluid energy cascade to “sense” the kinetic Landau damping at different scales is tested. The cascade rate estimated from the exact law and the dissipation calculated directly from the simulation are shown to reflect the role of Landau damping in dissipating energy at all scales, with an emphasis on the kinetic ones. This result provides new prospects on using exact laws for simplified fluid models to analyze dissipation in kinetic simulations and spacecraft observations, and new insights into theoretical description of collisionless magnetized plasmas.

*Unified Astronomy Thesaurus concepts: Space plasmas (1544); Magnetohydrodynamics (1964); Magnetohydrodynamical simulations (1966); Plasma astrophysics (1261)*

## 1. Introduction

The understanding of turbulent astrophysical plasmas remains to date a challenging problem: their chaotic nature and the complexity of the mechanisms at work in such media impose limitations to the methods one can use to study them efficiently. Yet enhancing our understanding of turbulent plasmas would provide the keys to solve a variety of problems related to energy dissipation, particle heating, and acceleration. Examples of systems where these processes are crucial include the solar wind (SW) and planetary magnetospheres (Goldstein et al. 1995; Bruno & Carbone 2005; Matthaeus & Velli 2011; Sahraoui et al. 2020), accretion flows around compact objects (Balbus & Hawley 1998; Quataert & Gruzinov 1999), and fusion devices (Diamond et al. 2005; Garbet 2006; Fujisawa 2021). In the SW, the heating problem is reflected by the slow decline of ion temperature (as a function of the radial distance from the Sun) in comparison with the prediction from the adiabatic expansion model of the wind (Richardson et al. 1995). Turbulence has long been proposed as a way to explain this behavior (Matthaeus et al. 1999), through scale-by-scale transfer of energy (i.e., cascade) toward small (kinetic) scales where dissipation is more effective (Schekochihin et al. 2009). A common tool used to estimate this energy dissipation is the formalism of exact law for fully developed turbulence first introduced by Kolmogorov (1941) to study incompressible neutral fluids. In this formalism, energy is assumed to be injected at large scales at a constant rate per unit volume $\epsilon$, which is assumed to be equal to the rate of cascade to smaller scales and to the rate of dissipation at those scales. Assuming statistical homogeneity and stationarity of the turbulent fields, and the existence of an inertial range in which both forcing and dissipation mechanisms are negligible, the cascade rate $\epsilon$ must remain constant in the inertial range (Kolmogorov 1941; Monin 1959; Antonia et al. 1997). The formalism of exact laws has been extended to (in)compressible magnetized plasmas within various approximations (Politano &

Pouquet 1998; Galtier 2008; Banerjee & Galtier 2013; Andrés & Sahraoui 2017; Andrés et al. 2018; Hellinger et al. 2018; Ferrand et al. 2019).

Exact laws have been used successfully to measure the energy cascade rate in the SW (Smith et al. 2006; Podesta et al. 2007; Sorriso-Valvo et al. 2007; MacBride et al. 2008; Marino et al. 2008; Carbone et al. 2009; Smith et al. 2009; Stawarz et al. 2009; Osman et al. 2011; Coburn et al. 2015; Banerjee et al. 2016; Hadid et al. 2017) and terrestrial magnetosheath (Hadid et al. 2018; Andrés et al. 2019). In those studies, the estimated cascade rate was interpreted as the turbulence energy dissipation rate, and hence used to quantify the amount of plasma heating due to turbulence (Sorriso-Valvo et al. 2007; Carbone et al. 2009; Banerjee et al. 2016). However, as explained above, such an equivalence between injection, cascade, and dissipation rates stems only from the hypothesis underlying exact-law derivation and cannot be demonstrated in spacecraft observations. Indeed, while in numerical simulations the injection, cascade, and dissipation rates can generally be estimated separately and compared to each other as done in this paper, estimating (irreversible) dissipation from spacecraft observation is a challenge and, generally, only the cascade rate, which is directly linked to measurable quantities through the exact law, is accessible (Sorriso-Valvo et al. 2007; Hadid et al. 2017). Thus, in spacecraft data, interpreting the energy cascade rate as the actual dissipation rate is not straightforward. This is particularly true because of the weakly collisional nature of the SW: in such plasmas classical viscous and/or resistive effects are absent, and dissipation is expected to occur via kinetic effects (e.g., Landau and cyclotron resonances; Leamon et al. 1998; Sahraoui et al. 2009, 2010; He et al. 2015; Chen et al. 2019) that are not captured by the usual fluid descriptions of plasmas. A fundamental question arises here: is the *fluid* turbulent cascade rate estimated in simulations and spacecraft observations of space plasmas representative of the actual *kinetic* dissipation in those media? It is the main goal of this paper to address this question,

<footnote n="1">¹Laboratoire de Physique des Plasmas, CNRS, École polytechnique, Université Paris-Saclay, Sorbonne Université, Observatoire de Paris-Meudon, F-91128 Palaiseau Cedex, France; renaud.ferrand@lpp.polytechnique.fr</footnote>

<footnote n="2">²Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Laboratoire J.L. Lagrange, Boulevard de l'Observatoire, CS 34229, F-06304 Nice Cedex 4, France
³Institut Universitaire de France (IUF), F-75231 Paris Cedex 05, France</footnote>
which impacts the use of fluid models to interpret part of in situ spacecraft observations in the near-Earth space and the theoretical (fluid versus kinetic) modeling of weakly collisional plasmas. In contrast with previous studies based on 2D hybrid particle-in-cell simulations (Hellinger et al. 2018; Bandyopadhyay et al. 2020), the use of 3D LF models give the possibility to isolate the influence of electron and ion Landau damping, neglecting all the other kinetic effects, and is therefore well suited to address the question of interest here.

## 2. Theoretical Model

Although we are dealing with weakly compressible regimes we chose, for the sake of simplicity, to use here the exact law derived by Ferrand et al. (2019) for incompressible Hall magnetohydrodynamics (HMHD; see Section 4 about the use of more general compressible models). Starting from the incompressible HMHD equations, and under the usual assumptions of time stationarity, space homogeneity, and infinite (kinetic and magnetic) Reynolds numbers, one can derive for the energy cascade rate in the inertial range the expression $\epsilon = \epsilon^{\text{MHD}} + \epsilon^{\text{Hall}}$, with

$$ \epsilon^{\text{MHD}} = -\frac{1}{4} \nabla_{\ell} \cdot \langle (|\delta\mathbf{v}|^2 + |\delta\mathbf{b}|^2) \delta\mathbf{v} - 2(\delta\mathbf{v} \cdot \delta\mathbf{b}) \delta\mathbf{b} \rangle, \quad (1) $$

$$ \epsilon^{\text{Hall}} = -\frac{1}{8} d_i \nabla_{\ell} \cdot \langle 2(\delta\mathbf{b} \cdot \delta\mathbf{j}) \delta\mathbf{b} - |\delta\mathbf{b}|^2 \delta\mathbf{j} \rangle, \quad (2) $$

where $\mathbf{v}$, $\mathbf{b} = \mathbf{B}/\sqrt{\mu_0\rho_0}$, and $\mathbf{j} = \nabla \times \mathbf{b}$ are the velocity, magnetic field, and electric current in Alfvén units ($\rho_0$ is the constant mass density) and $d_i$ is the ion inertial length. Fields are taken at points $\mathbf{x}$ and $\mathbf{x}'$ separated by a spatial increment $\ell = \mathbf{x}' - \mathbf{x}$, and the notations $\mathbf{v} \equiv \mathbf{v}(\mathbf{x})$ and $\mathbf{v}' \equiv \mathbf{v}(\mathbf{x}')$ are adopted. We then define the increment operator $\delta$ as $\delta\mathbf{v} = \mathbf{v}' - \mathbf{v}$, and $\nabla_{\ell}$ as the derivative operator with respect to the increment $\ell$.

## 3. Simulation Data

### 3.1. Presentation of the Data

In this study, HMHD-CGL refers to a fluid model with anisotropic ion pressure whose gyrotropic components parallel and perpendicular to the local magnetic field obey nonlinear dynamical equations where the heat fluxes are neglected (biadiabatic approximation introduced by Chew, Goldberger, and Low and known as the CGL approximation; Chew et al. 1956). The electrons are assumed to be isothermal. Differently, the LF model retains the nonlinear dynamics of the parallel and perpendicular pressures and heat fluxes for both the ions and electrons, and involves a closure at the level of the fourth-order moments, consistent with the low-frequency linear kinetic theory (Snyder et al. 1997; Passot & Sulem 2007). The main assumption for modeling Landau damping consists in retaining the imaginary contribution of the plasma response function in the closure relation that expresses the last retained fluid moment of the hierarchy in terms of the lower ones. In Fourier space, this procedure generates factors of the form "$i \operatorname{sgn}(k_z)$" which, in physical space, identifies with the Hilbert transform along the ambient magnetic field (Hammett & Perkins 1990; Hunana et al. 2019). It is then possible to generalize this formulation to take into account magnetic field line distortion, using the convolution form of the Hilbert transform (Snyder et al. 1997). Its approximation in the numerical code is discussed in Passot et al. (2014). In both models, finite ion and

<table><caption>**Table 1**
List of Runs and Their Relevant Parameters, where CGLx and LFx Refer to HMHD-CGL and LF Simulations, Respectively</caption><thead><tr><th>Run</th><th>$k_{\perp,\ell}d_i$</th><th>Resolution</th><th>$\theta$</th><th>$\nu = \eta$</th><th>$\alpha$</th></tr></thead><tbody><tr><td>CGL1</td><td>0.045</td><td>512<sup>3</sup></td><td>83°</td><td>7.35 × 10<sup>-8</sup></td><td>80</td></tr><tr><td>CGL2</td><td>0.045</td><td>512<sup>3</sup></td><td>75°</td><td>7.35 × 10<sup>-8</sup></td><td>10</td></tr><tr><td>CGL3</td><td>0.5</td><td>512<sup>2</sup> × 1024</td><td>75°</td><td>10<sup>-14</sup></td><td>2.5</td></tr><tr><td>CGL4</td><td>0.011</td><td>1024<sup>3</sup></td><td>75°</td><td>3 × 10<sup>-3</sup></td><td>5</td></tr><tr><td>LF1</td><td>0.045</td><td>512<sup>3</sup></td><td>83°</td><td>7.35 × 10<sup>-8</sup></td><td>1</td></tr><tr><td>LF2</td><td>0.045</td><td>512<sup>3</sup></td><td>75°</td><td>7.35 × 10<sup>-8</sup></td><td>1</td></tr><tr><td>LF3</td><td>0.5</td><td>432<sup>3</sup></td><td>75°</td><td>7 × 10<sup>-14</sup></td><td>1.5</td></tr><tr><td>LF4</td><td>0.011</td><td>512<sup>3</sup></td><td>75°</td><td>3 × 10<sup>-3</sup></td><td>2</td></tr></tbody></table>

**Note.** The ratio of the longitudinal to transverse box sizes is given by $\tan(\theta)$.

electron Larmor radius corrections are neglected, thus reducing the kinetic effects to Landau damping. Ohm's law includes the Hall term and the electron pressure contribution. Turbulence is forced with counter-propagating kinetic Alfvén waves (KAWs) making an angle $\theta$ with the ambient magnetic field, at the largest scales of the simulation domain. This corresponds to transverse wavenumbers $k_{\perp,\ell}$, whose values are summarized in Table 1. The amplitudes obey a Langevin equation, with an oscillation frequency given by the KAW linear dispersion relation (TenBarge et al. 2014). We also introduce two thresholds in order to constrain the sum of perpendicular kinetic and magnetic energies to stay within a certain range. Small-scale dissipation is ensured by the hyperviscosity and hyperdiffusivity terms in the velocity and induction equations, of the form $\mathbf{d}_{\nu} = \nu(\Delta_{\perp} + \alpha\partial_{z}^{2})^{4}\mathbf{v}$ and $\mathbf{d}_{\eta} = \eta(\Delta_{\perp} + \alpha\partial_{z}^{2})^{4}\mathbf{b}$, with $\alpha$ being an anisotropy coefficient.

In all the simulations, $\beta_i = 1$ and the ion and electron pressures are taken to be isotropic and equal initially. The other parameters are reported in Table 1. The simulations are performed using a desaliased spectral code (at 2/3 of the maximum wavenumber) with a third-order Runge-Kutta scheme for time stepping.

Two different propagation angles for the KAWs driven at the largest scale of each LF simulation were chosen, hence tuning Landau damping to two different levels (Kobayashi et al. 2017). This can be seen in Figure 1, which compares the linear dispersion relation and damping rate of the KAWs: the higher the propagation angle, the lower the damping rate at a given scale. Note that, while changing the angle, we do not change the amplitude of the fluctuations at the driving scale (i.e., $\omega_{NL}$ remains constant), and so the nonlinear parameter $\chi = \omega_{NL}/\omega_{L}$ also varies: when the angle decreases, $k_{\parallel}$ increases and so does $\omega_{L}$, thus $\chi$ is reduced. As the ratio $\gamma/\omega_{L}$ is approximately constant for high oblique angles $\theta$ (e.g., Figure 7 in Sahraoui et al. 2012), the strength of the Landau damping relative to the cascade rate $\gamma/\omega_{NL} = (1/\chi)(\gamma/\omega_{L})$ thus increases as the angle decreases.

### 3.2. Energy Balance and Time Stationarity

Let $\mathcal{E}_{\text{tot}}(t)$ be the total energy of the system at time $t$, $\mathcal{I}_t(t)$ the injection rate due to the external forcing on the perpendicular velocity components, and $\mathcal{D}_h(t)$ the total dissipation rate due to the hyperviscous and hyperdiffusive terms. Since Landau damping does not affect the total energy balance, total energy conservation implies

$$ \frac{d}{dt} \mathcal{E}_{\text{tot}}(t) = \mathcal{I}_t(t) - \mathcal{D}_h(t). \quad (3) $$
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Frequency $\omega$ and damping rate $\gamma$ (normalized by the ion gyrofrequency $\omega_{ci}$) of KAWs vs. the normalized transverse wavenumber $k_{\perp}d_i$ (where $d_i$ refers to the ion inertial length) for the LF model at the two propagation angles $\theta = 83^{\circ}$ and $75^{\circ}$ used in driving the simulations ($\beta_i = 1$, $T_i = T_e$).</figcaption></figure>

Denoting by $\mathcal{E}_{\text{int}}$ and $\mathcal{E}_{\parallel}$ the parts of the total energy associated with the pressure components (internal energy) and the parallel velocity and magnetic field components entering the kinetic and magnetic parts, we can write

$$ \frac{d}{dt} \mathcal{E}_{\text{tot}}(t) - \frac{d}{dt} \mathcal{E}_{\text{int}}(t) - \frac{d}{dt} \mathcal{E}_{\parallel}(t) \equiv \frac{d}{dt} \mathcal{E}_{\perp}(t) \approx 0, \quad (4) $$

where $\mathcal{E}_{\perp}$ is the sum of the perpendicular kinetic and magnetic energies, a quantity bound to remain nearly constant by the forcing procedure.

Because of computational constraints, the time evolution of the different energy components is computed for low-resolution (LR) simulations analogs of runs CGL3 and LF3 and shown in Figure 2 (injection and dissipation rates needed to perform this extra study were not output at a high-enough frequency in the high-resolution simulations). From this figure it is conspicuous that the time evolution of total energy, injection, and hyperdissipation is consistent with the energy conservation (3). Moreover, one can see the driving procedure at play in keeping the perpendicular energy $\mathcal{E}_{\perp}$ roughly constant. Its time stationarity is in practice established when the hyperdissipation rate has reached a constant value. When comparing CGL3-LR with LF3-LR, one notices that run LF3-LR requires a larger injection rate to maintain the same level of turbulence on the magnetic and perpendicular velocity than in run CGL3-LR, since Landau damping efficiently converts a part of the injected energy into internal energy. This is evidenced by the dashed green curve in Figure 2 (bottom), which shows that the increase of the internal energy is consistent with the heating by heat fluxes. Moreover, the hyperdissipation rate is lower on run LF3-LR, suggesting that part of the cascading energy is taken by Landau damping, as will be evidenced in the next section.

### 4. Calculation of the Energy Cascade Rate

Using Equations (1)–(2), we compute the transverse energy cascade rate for each simulation as a function of the perpendicular increment by averaging over all spatial positions in the simulation box and different increment vectors, at a time for which the simulations reached a stationary state. Increment vectors $\ell = (\ell_{\perp}, \ell_{\parallel})$ are selected following the angle averaging method of Taylor et al. (2003), and only the increments forming an angle of at least $45^{\circ}$ with the parallel direction are retained. As already mentioned, the exact law used is the one for incompressible HMHD, whereas the simulations are weakly compressible. A comparison (not shown) with a full compressible HMHD exact law (Andrés et al.

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Low-resolution runs of CGL3 (CGL3-LR, top) and LF3 (LF3-LR, bottom): time evolution of the total energy injected in the system $\int_0^t \mathcal{I}_t dt$ (solid red line), total energy $\mathcal{E}_{\text{tot}}$ (solid blue), internal energy $\mathcal{E}_{\text{int}}$ (solid green), perpendicular energy $\mathcal{E}_{\perp}$ (solid magenta, roughly constant), and the time-integrated hyperdissipation $\int_0^t \mathcal{D}_h dt$ (dashed red). The piece of the dashed green curve (right) starting at $t = 250$, whose vertical position is arbitrary, displays the heating due to heat fluxes, which is consistent with the increase of internal energy.</figcaption></figure>

2018) only showed slight change ($\lesssim 10\%$ in the inertial range) of the cascade rate with respect to the current estimate from the incompressible model. The transverse cascade rate is then averaged over all increments of an equal value of $\ell_{\perp}$. The transverse hyperdissipation is computed in Fourier space as

$$ \varepsilon^{\text{diss}}(\ell_{\perp}) = \int_{0}^{k_{\perp}} dk'_{\perp} \int_{0}^{k'_{\perp}} k'^{8}_{\perp} (\eta |\mathbf{b}(\mathbf{k}')|^2 + \nu |\mathbf{v}(\mathbf{k}')|^2) k'_{\perp} d\theta' dk'_{z}, \quad (5) $$

where we use $\ell_{\perp} = \pi/k_{\perp}$.

For simulations forced at intermediate scales, whose results are reported in Figure 3, the behavior of the MHD and Hall contributions to the energy cascade rate are similar, with the latter rising up at sub-ion scales, then dominating the former at about the ion inertial length. The total energy cascade rate is roughly constant on more than one decade of scales in the simulations with $\theta = 83^{\circ}$, in particular in CGL1, which demonstrates the existence of an inertial range. To highlight the effect of Landau damping on the cascade rate, we compare in Figure 4 the cascade rates from the LF simulations with $\theta = 75^{\circ}$ and $\theta = 83^{\circ}$, normalized to the corresponding ones from the CGL simulations. We observe a stronger decrease (by up to a factor 5) in the normalized cascade rate at small scales
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Energy cascade rate $\epsilon$ and its ideal MHD and Hall components together with the transverse hyperdissipation computed for runs CGL1 (top) and LF1 (bottom). Plain lines represent positive values and dashed lines negative values.</figcaption></figure>

for LF2 ($\theta = 75^\circ$), i.e., for the simulation with the strongest Landau damping, than for LF1 ($\theta = 83^\circ$) for which the normalized cascade rate remains nearly constant at all scales. This result clearly relates the enhancement of Landau damping at kinetic scales to the decline of the energy cascade rate at these scales. We note also the consistency between the (transverse) hyperdissipation and the cascade rate at the smallest scale of the simulation box (Figure 3).

We complement our study with the cascade rates estimated from simulations forced at even smaller scales (LF3 and CGL3) with $\theta = 75^\circ$ and reported in Figure 5. Simulation LF3 exhibits a strong decrease in $\epsilon^{MHD}$ partially compensated by a quick rising of the Hall component, giving no clear inertial range, in contrast to CGL3, which still behaves similarly to the simulations forced at intermediate scales. As shown below, this effect may be attributed to the fact that Landau dissipation reaches high levels at the sub-ion scales of LF3, whereas CGL3 contains no dissipation mechanisms other than hyperviscosity and hyperdiffusivity, which are bound to act only at the smallest scales. Note that the sudden changes of sign observed at large scales in some components of the cascade rates in Figures 3 and 5 are likely to be due to the proximity of the forcing. Those observed at small scales for the MHD component of run CGL3 would result from numerical errors in the calculation of $\epsilon^{MHD}$ given its very small magnitude at those scales.

To obtain a full picture as to how Landau damping affects the energy cascade rate, we performed simulations forced at large scales (LF4 and CGL4). Combining the runs CGL2-3-4 and LF2-3-4 we construct a multiscale energy cascade rate over nearly three decades of scales that highlights the effect of Landau damping on it. As the simulations were run at different scales, the amplitude of the forcing was changed to ensure that each simulation reaches a fully turbulent state. Therefore, we renormalized the cascade rate $\epsilon$ obtained from the different

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Ratios of the energy cascade rate computed for LF simulations over the one for CGL simulations for a driving wave angle $\theta = 83^\circ$ (black) and $\theta = 75^\circ$ (red).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Same as in Figure 3, but for runs CGL3 (top) and LF3 (bottom).</figcaption></figure>

simulations to match the one of intermediate runs CGL2 and LF2, while taking care to discard the smallest scales of intermediate- and large-scale forcing cascade rates to ensure that hyperviscosity is not acting at intermediate scales of the reconstructed energy cascade. Figure 6 shows the full energy cascade rate for CGL and LF runs for the driving wave angle $\theta = 75^\circ$. CGL runs exhibit an almost constant energy cascade rate over two and a half decades of scales, whereas $\epsilon_{LF}$ decreases steadily over scales and reaches its minimum value at the smallest ones, confirming that the behavior already observed in Figure 4 remains valid over a broader range of scales.

# 5. Influence of Landau Dissipation

## 5.1. Heating due to Heat Fluxes

In the wake of the previous results an important question arises: can the drop in the energy cascade rate for LF runs be directly connected to Landau damping? For this purpose, we calculate the heating due to heat fluxes in the presence of Landau damping. For each species, the pressure equations with
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Energy cascade rates reconstructed with CGL2-3-4 runs and LF2-3-4 runs. The ranges spanned by each simulation are delimited by the black dotted lines. A slight irregularity is observed on the green curve at the transition between CGL2 and CGL3, which is caused by an insufficient overlap of the cascade rates at these scales.</figcaption></figure>

the Hall term and the gyrotropic heat fluxes read

$$ \frac{d}{dt} \ln \left( \frac{p_{\parallel} |B|^2}{\rho^3} \right) = - \frac{2c}{|B|} \hat{\mathbf{b}} \cdot \nabla \times E_H - \frac{1}{p_{\parallel}} \left( -2q_{\perp} \nabla \cdot \hat{\mathbf{b}} + \nabla \cdot (q_{\perp} \hat{\mathbf{b}}) \right), \quad (6) $$

$$ \frac{d}{dt} \ln \left( \frac{p_{\perp}}{\rho |B|} \right) = \frac{c}{|B|} \hat{\mathbf{b}} \cdot \nabla \times E_H - \frac{1}{p_{\perp}} \left( q_{\perp} \nabla \cdot \hat{\mathbf{b}} + \nabla \cdot (q_{\perp} \hat{\mathbf{b}}) \right). \quad (7) $$

We define the parallel, perpendicular, and total entropies per unit mass

$$ s_{\parallel} = \frac{c_V}{3} \ln \left( \frac{p_{\parallel} |B|^2}{\rho^3} \right), \quad s_{\perp} = \frac{2c_V}{3} \ln \left( \frac{p_{\perp}}{\rho |B|} \right), \quad (8) $$

$$ s = s_{\parallel} + s_{\perp} = \frac{c_V}{3} \ln \left( \frac{p_{\parallel} p_{\perp}^2}{\rho^5} \right), \quad (9) $$

where $c_V$ is the specific heat at constant volume. Denoting by $e$ the internal energy per unit mass, the internal energy per unit volume reads $E \equiv \rho e = p_{\perp} + \frac{1}{2} p_{\parallel} = \frac{3}{2} n T$ where $T = \frac{1}{3}(2T_{\perp} + T_{\parallel})$. From $e = c_V T$, one gets $c_V = \frac{3}{2m}$ (the Boltzmann constant is included in the definition of temperature). The total entropy then obeys

$$ \begin{aligned} \partial_t(\rho s) + \nabla \cdot \left( \rho s \mathbf{u} + \left( \frac{q_{\perp}}{T_{\perp}} + \frac{q_{\parallel}}{2T_{\parallel}} \right) \hat{\mathbf{b}} \right) & \\ = \left( \frac{1}{T_{\parallel}} - \frac{1}{T_{\perp}} \right) q_{\perp} \nabla \cdot \hat{\mathbf{b}} & \\ - \left( \frac{q_{\perp}}{T_{\perp}} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\perp} + \frac{q_{\parallel}}{2T_{\parallel}} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\parallel} \right) & \end{aligned} \quad (10) $$

From the form of the right-hand side of Equations (6)–(7), we can conclude that the rates of change of the parallel and perpendicular entropies per unit mass ($s_{\parallel}^p$ and $s_{\perp}^p$, respectively) associated with a production (or destruction) and excluding transport or exchanges between the parallel and perpendicular

directions (see, e.g., Hazeltine et al. 2013) are given by

$$ \frac{d}{dt} s_{\parallel}^{p} = \frac{1}{\rho T_{\parallel}} q_{\perp} \nabla \cdot \hat{\mathbf{b}} - \frac{q_{\parallel}}{2 \rho T_{\parallel}} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\parallel}, \quad (11) $$

$$ \frac{d}{dt} s_{\perp}^{p} = - \frac{1}{\rho T_{\perp}} q_{\perp} \nabla \cdot \hat{\mathbf{b}} - \frac{q_{\perp}}{\rho T_{\perp}} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\perp}. \quad (12) $$

The associated rates of heat production per unit mass are related by $dQ_{\parallel}/dt = T_{\parallel} ds_{\parallel}^p/dt$ and $dQ_{\perp}/dt = T_{\perp} ds_{\perp}^p/dt$. We thus get for the total heat production $Q = Q_{\parallel} + Q_{\perp}$

$$ \partial_t(\rho Q) + \nabla \cdot (\rho Q \mathbf{u}) = -\frac{q_{\parallel}}{2} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\parallel} - q_{\perp} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\perp}. \quad (13) $$

The global heating is thus given by

$$ H = - \int \left( \frac{q_{\parallel}}{2} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\parallel} + q_{\perp} (\hat{\mathbf{b}} \cdot \nabla) \ln T_{\perp} \right) d^3x, \quad (14) $$

where $q_{\parallel}$ and $q_{\perp}$ are the heat fluxes obtained from the integration of the model closed at the level of the fourth-rank moments.

We can define a spectral density for the heating rate $H$ (also referred to as a cospectrum) in the form

$$ H(\mathbf{k}) = -\frac{1}{2} \left( \frac{1}{2} \mathcal{F}\{q_{\parallel}\}(-\mathbf{k}) \mathcal{F}\{(\hat{\mathbf{b}} \cdot \nabla) \ln T_{\parallel}\}(\mathbf{k}) + \mathcal{F}\{q_{\perp}\}(-\mathbf{k}) \mathcal{F}\{(\hat{\mathbf{b}} \cdot \nabla) \ln T_{\perp}\}(\mathbf{k}) + \text{c.c.} \right), \quad (15) $$

where $\mathcal{F}$ denotes the Fourier transform.

A few remarks can be made here:

1. In all the simulations we have performed, the volume-integrated heat production is observed to be positive, but its pointwise value can be negative in relatively small regions of space. This contrasts with the (semi)collisional regime where the heat fluxes obey Fourier laws of the form $q = -\kappa (\hat{\mathbf{b}} \cdot \nabla) T$, making the heat production positive everywhere in space.

2. Inserting in Equation (14) the quantities $q_{\perp}$ and $q_{\parallel}$ obtained by the integration of the dynamical equation for the heat fluxes results in taking into account in the heating rate contributions originating from the heat flux present when a quasi-normal closure is implemented (i.e., where the fourth-rank cumulants are taken to be equal to zero, thus making the Landau damping disappear). In the present simulations, this contribution does not exceed 15% of the total heating rate. In order to only deal with the heat flux originating from the Landau damping, it would be necessary to define a conserved entropy for the quasi-normal closure and evaluate its rate of change due to the introduction of Landau damping. This is left for future work as it is not straightforward.

3. More importantly, this heating rate takes into account the Landau damping on all the waves present in the simulations, including the magnetosonic waves. At this level, it appears difficult to separate the contributions of the KAW and to evaluate their dissipation by Landau damping. Nevertheless, these magnetosonic waves get dissipated at large scales; thus, at small enough scales the estimated heating rate mostly results from Landau damping of KAWs and it becomes possible to compare it to the cascading energy. This particularity is also the reason why Landau damping appears to be acting at all scales in all the results presented above, even in simulations forced at large scales.
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Spectral densities of the heating rate $D_L(k_⊥)$ (red) and of the magnetic (blue) and kinetic (green) hyperdissipation as functions of the transverse wavenumber $k_⊥$ for run LF3. A straight line of slope 5.2 is supplemented, for comparison with the scale variation of the magnetic hyperdissipation.</figcaption></figure>

The fact that Landau damping is present at all scales in the simulation can be seen by estimating the spectral density of the total heating rate at a given wavenumber $k_⊥$, $D_L(k_⊥) = \int k_⊥H(k)dk_⊥dθ$, where $H(k)$ is the sum of the spectral densities given by Equation (15) for both the ions and the electrons. This spectral density is represented in Figure 7 along with the densities of hyperdissipation and hyperdiffusivity. One clearly sees that the heating rate due to the presence of heat fluxes dominates hyperdissipation over a broad range of scales due to the dissipation of KAWs and magnetosonic modes, the two becoming comparable only at the smallest scales (note that the magnetic hyperdissipation dominates at small scales over the kinetic one).

## 5.2. Dissipation due to Landau Damping

The energy $E_⊥(t)$ of the (quasi-incompressible) KAWs that cascade toward small scales, and which is the subject of our study, obeys

$$ \frac{d}{dt} E_⊥(t) = \mathcal{I}^C(t) - \mathcal{D}_L^C(t) - \mathcal{D}_h^C(t), \quad (16) $$

where $\mathcal{I}^C$ is the part of the injection rate that contributes to the KAW cascade (the other part is transferred to magnetosonic modes that are dominantly dissipated at large scales), while $\mathcal{D}_L^C$ and $\mathcal{D}_h^C$ are the parts of the Landau and hyperviscous (and hyperdiffusive) dissipation that affect the cascading modes. Using cylindrical coordinates and assuming time stationarity, one can write the integrated energy balance at each Fourier mode as (adopting roman scripts for spectral densities):

$$ \epsilon(k_{\perp}) = \int_{0}^{k_{\perp}} \{I^{C}(k'_{\perp}) - D_{L}^{C}(k'_{\perp}) - D_{h}^{C}(k'_{\perp})\} dk'_{\perp}. \quad (17) $$

Considering two wavenumbers $k_{⊥1}$ and $k_{⊥2}$ large enough so that the forcing (which is concentrated at large scales) leads to $\int_{0}^{k_{⊥1}} I^C(k'_{⊥})dk'_{⊥} = \int_{0}^{k_{⊥2}} I^C(k'_{⊥})dk'_{⊥} = \mathcal{I}^C$, yet small enough for

<figure><img src="image_9.png" /><figcaption>**Figure 8.** Energy cascade rate $\epsilon$ (red) and transverse hyperdissipation (violet) for run LF3. The orange and brown curves show the same $\epsilon$ corrected by Landau damping integrated between $\ell_⊥$ and a reference scale $\ell_⊥ = 0.1\pi d_i$ and $\ell_⊥ = 0.2\pi d_i$, respectively.</figcaption></figure>

hyperviscous dissipation to be negligible, one obtains

$$ \epsilon(k_{\perp 1}) - \epsilon(k_{\perp 2}) = \int_{k_{\perp 1}}^{k_{\perp 2}} D_L^C(k'_\perp) dk'_\perp \lesssim \int_{k_{\perp 1}}^{k_{\perp 2}} D_L(k'_\perp) dk'_\perp. \quad (18) $$

The inequality draws closer to an equality for values of $k_{⊥}$ large enough so that all magnetosonic modes have been dissipated.

Equation (18) can be used to estimate a correction to the energy cascade rate that would take into account the energy lost due to Landau damping. We do so for run LF3: using this equation we add to the transfer rate the cumulative Landau dissipation between an arbitrary scale (chosen however to be not too large nor too small) and the running (smaller) scale $l_⊥$. Two of these resulting corrected rates $\epsilon^{corr}$ are shown in Figure 8. They appear to be almost constant, and as such they behave very similarly to the transfer rate of run CGL3 (Figure 5). The slight increase of $\epsilon^{corr}$ toward small scales probably reflects the (weak) contribution of some remaining magnetosonic waves to the calculated Landau damping. This clearly demonstrates that the energy lost along the cascade due to Landau damping is well captured by the decline of the (fluid) cascade rate at the corresponding scales.

A complementary estimate of energy dissipation can be done in Fourier space by also taking into account hyperdissipation. Indeed, assuming stationarity, one can also derive that

$$ \epsilon(k_{\perp}) = \mathcal{I}^C - \mathcal{D}_L^C + \int_{k_{\perp}}^{\infty} D_L^C(k'_{\perp}) dk'_{\perp} \\ - \int_{0}^{k_{\perp}} D_h^C(k'_{\perp}) dk'_{\perp} = \int_{k_{\perp}}^{\infty} \{D_h^C(k'_{\perp}) + D_L^C(k'_{\perp})\} dk'_{\perp}. \quad (19) $$

Equation (19) indicates that, as expected, the rate of energy transfer at the wavenumber $k_⊥$ identifies with the sum of the rates of Landau and hyperdissipation beyond this wavenumber. One can compare the second right-hand-side term of this equation to the energy cascade rate $\epsilon(k_⊥)$ obtained from the IHMHD exact law, as displayed in Figure 9. The difference between the two curves, which is especially significant at large scales, is due to the fact that the estimation of the dissipation includes the Landau damping of magnetosonic modes, whereas the cascade rate considers only incompressible modes. At smaller scales, however, where magnetosonic modes have already been dissipated, the dissipation and cascade rates decrease parallel to each other: this indicates that, at scales not
<figure><img src="image_10.png" /><figcaption>**Figure 9.** Energy cascade rate $\epsilon(k_\perp)$ (black line) together with Landau and hyperdissipation (red line) computed with Equation (19) for run LF3.</figcaption></figure>

yet affected by hyperdissipation, the decay of $\epsilon(k_\perp)$ in a spectral interval identifies with Landau dissipation within this interval.

Figures 8 and 9 clearly demonstrate that, through the cascade, the energy lost due to Landau damping is well captured by the decline of the (fluid) cascade rate at the corresponding scales. Note that a similar decline of the fluid cascade rate at kinetic scales was reported in 2D hybrid PIC simulations and spacecraft observations in the SW and magnetosheath (Hellinger et al. 2018; Bandyopadhyay et al. 2020). Also, Sorriso-Valvo et al. (2019) found a correlation between enhancement of a proxy of the local cascade rate and the signatures of wave-particle interactions in data from the Magnetospheric MultiScale (MMS) mission.

## 6. Conclusion

In this study, we tackle a fundamental question about the ability of fluid exact laws to reflect the presence of kinetic (Landau) damping. By constructing multiscale energy cascade and dissipation rates using the HMHD model on a variety of turbulence simulations bearing different intensities of Landau damping, we showed that the presence of Landau damping at small (kinetic) scales is reflected by the steady decline of the energy cascade rate at the same scales, which was found to be comparable to the effective Landau dissipation at those scales. By demonstrating the ability of a fluid exact law to provide a correct estimate of kinetic dissipation in the sub-ion range of numerical simulations, this work provides a means to evaluate the amount of energy that is dissipated into particle heating in spacecraft data: the decline of the cascade rate allows one to evaluate the kinetic dissipation as a function of scale. This should help in investigating (at least partially) a longstanding problem in astrophysical plasmas about energy partition between ions and electrons (Kawazura et al. 2019), which are generally heated at different scales.

The study presented in this paper only makes use of Landau damping. It would be interesting in future works to extend these conclusions to a broad variety of kinetic effects and to test them on more general simulations of the SW, featuring a plasma turbulence driven by other types of waves than slightly perturbed KAWs. It is also important to stress that, even if the oversimplified (yet fully nonlinear) fluid models of turbulence can provide good estimates of the amount of energy that is dissipated into particle heating, they do not specify how this dissipation occurs. The answer to this question and those

related to the fate of energy when handed to the plasma particles requires a kinetic treatment.

This work was granted access to the HPC resources of CINES/IDRIS under the allocation A0060407042. Part of the computations have also been done on the “Mesocentre SIGAMM” machine, hosted by Observatoire de la Côte d’Azur.

## ORCID iDs

R. Ferrand https://orcid.org/0000-0002-0044-3405

P. L. Sulem https://orcid.org/0000-0002-6231-5126

S. Galtier https://orcid.org/0000-0001-8685-9497

## References

Andrés, N., Galtier, S., & Sahraoui, F. 2018, PhRvE, 97, 013204

Andrés, N., & Sahraoui, F. 2017, PhRvE, 96, 053205

Andrés, N., Sahraoui, F., Galtier, S., et al. 2019, PhRvL, 123, 245101

Antonia, R., Ould-Rouis, M., Anselmet, F., & Zhu, Y. 1997, JFM, 332, 395

Balbus, S. A., & Hawley, J. F. 1998, RvMP, 70, 1

Bandyopadhyay, R., Sorriso-Valvo, L., Chasapis, A., et al. 2020, PhRvL, 124, 225101

Banerjee, S., & Galtier, S. 2013, PhRvE, 87, 013019

Banerjee, S., Hadid, L., Sahraoui, F., & Galtier, S. 2016, ApJL, 829, L27

Bruno, R., & Carbone, V. 2005, LRSP, 2, 4

Carbone, V., Marino, R., Sorriso-Valvo, L., Noullez, A., & Bruno, R. 2009, PhRvL, 103, 061102

Chen, C. H. K., Klein, K. G., & Howes, G. G. 2019, NatCo, 10, 740

Chew, G. F., Goldberger, M. L., & Low, F. E. 1956, RSPSA, 236, 112

Coburn, J., Forman, M., Smith, C., Vasquez, B., & Stawarz, J. 2015, RSPTA, 373, 20140150

Diamond, P. H., Itoh, S. I., & Hahm, T. S. 2005, PPCF, 47, R35

Ferrand, R., Galtier, S., Sahraoui, F., et al. 2019, ApJ, 881, 50

Fujisawa, A. 2021, PJAB, 97, 103

Galtier, S. 2008, PhRvE, 77, 015302

Garbet, X. 2006, CRPhy, 7, 573

Goldstein, M. L., Roberts, D. A., & Matthaeus, W. 1995, ARA&A, 33, 283

Hadid, L., Sahraoui, F., & Galtier, S. 2017, ApJ, 838, 9

Hadid, L. Z., Sahraoui, F., Galtier, S., & Huang, S. Y. 2018, PhRvL, 120, 055102

Hammett, G. W., & Perkins, F. W. 1990, PhRvL, 64, 3019

Hazeltine, R. D., Mahajan, S. M., & Morrison, P. J. 2013, PhPl, 20, 022506

He, J., Wang, L., Tu, C., Marsch, E., & Zong, Q. 2015, ApJL, 800, L31

Hellinger, P., Verdini, A., Landi, S., Franci, L., & Matteini, L. 2018, ApJ, 857, L19

Hunana, P., Tenerani, A., Zank, G. P., et al. 2019, JPIPh, 85, 205850603

Kawazura, Y., Barnes, M., & Schekochihin, A. 2019, PNAS, 116, 771

Kobayashi, S., Sahraoui, F., Passot, T., et al. 2017, ApJ, 839, 122

Kolmogorov, A. 1941, DoSSR, 32, 16

Leamon, R. J., Smith, C. W., Ness, N. F., Matthaeus, W. H., & Wong, H. K. 1998, JGR, 103, 4775

MacBride, B., Smith, C., & Forman, M. 2008, ApJ, 679, 1644

Marino, R., Sorriso-Valvo, L., Carbone, V., et al. 2008, ApJ, 677, L71

Matthaeus, W. H., & Velli, M. 2011, SSRv, 160, 145

Matthaeus, W. H., Zank, G. P., Smith, C. W., & Oughton, S. 1999, PhRvL, 82, 3444

Monin, A. S. 1959, DoSSR, 125, 515

Osman, K. T., Wan, M., Matthaeus, W. H., Weygand, J. M., & Dasso, S. 2011, PhRvL, 107, 165001

Passot, T., Henri, P., Laveder, D., & Sulem, P.-L. 2014, EPJD, 68, 207

Passot, T., & Sulem, P. L. 2007, PhPl, 14, 082502

Podesta, J. J., Forman, M. A., & Smith, C. W. 2007, PhPl, 14, 092305

Politano, H., & Pouquet, A. 1998, PhRvE, 57, 21

Quataert, E., & Gruzinov, A. 1999, ApJ, 520, 248

Richardson, J. D., Paularena, K. I., Lazarus, A. J., & Belcher, J. W. 1995, GeoRL, 22, 325

Sahraoui, F., Belmont, G., & Goldstein, M. L. 2012, ApJ, 748, 100

Sahraoui, F., Goldstein, M. L., Belmont, G., Canu, P., & Rezeau, L. 2010, PhRvL, 105, 131101

Sahraoui, F., Goldstein, M. L., Robert, P., & Khotyaintsev, Y. V. 2009, PhRvL, 102, 231102

Sahraoui, F., Hadid, L., & Huang, S. 2020, RvMPP, 4, 4

Schekochihin, A., Cowley, S. C., Dorland, W., et al. 2009, ApJS, 182, 310
Smith, C. W., Hamilton, K., Vasquez, B. J., & Leamon, R. J. 2006, ApJL, 645, L85

Smith, C. W., Stawarz, J. E., Vasquez, B. J., Forman, M. A., & MacBride, B. T. 2009, PhRvL, 103, 201101

Snyder, P. B., Hammett, G. W., & Dorland, W. 1997, PhPI, 4, 3974

Sorriso-Valvo, L., Catapano, F., Retinò, A., et al. 2019, PhRvL, 122, 035102

Sorriso-Valvo, L., Marino, R., Carbone, V., et al. 2007, PhRvL, 99, 115001

Stawarz, J. E., Smith, C. W., Vasquez, B. J., Forman, M. A., & MacBride, B. T. 2009, ApJ, 697, 1119

Taylor, M. A., Kurien, S., & Eyink, G. L. 2003, PhRvE, 68, 026310

TenBarge, J., Howes, G., Dorland, W., & Hammett, G. 2014, CoPhC, 185, 578