eurm10 msam10 msam10

# The Fundamental Parameters of Astrophysical Plasma Turbulence and its
Dissipation: Nonrelativistic Limit

Gregory G. Howes  Department of Physics and Astronomy, University of Iowa,
Iowa City IA 54224, USA

(?; revised ?; accepted ?. - To be entered by editorial office)

###### Abstract

A specific set of dimensionless plasma and turbulence parameters is introduced
to characterize the nature of turbulence and its dissipation in weakly
collisional space and astrophysical plasmas. Key considerations are discussed
for the development of predictive models of the turbulent plasma heating that
characterize the partitioning of dissipated turbulent energy between the ion
and electron species and between the perpendicular and parallel degrees of
freedom for each species. Identifying the kinetic physical mechanisms that
govern the damping of the turbulent fluctuations is a critical first step in
constructing such turbulent heating models. A set of ten general plasma and
turbulence parameters are defined, and reasonable approximations along with
the exploitation of existing scaling theories for magnetohydrodynamic
turbulence are used to reduce this general set of ten parameters to just three
parameters in the isotropic temperature case. A critical step forward in this
study is to identify the dependence of all of the proposed kinetic mechanisms
for turbulent damping in terms of the same set of fundamental plasma and
turbulence parameters. Analytical estimations of the scaling of each damping
mechanism on these fundamental parameters are presented, and this information
is synthesized to produce the first phase diagram for the turbulent damping
mechanisms as a function of driving scale and ion plasma beta.

##  1 Introduction

Turbulence is a fundamental, yet incompletely understood, process in space and
astrophysical plasmas that mediates the transfer of the energy of chaotic
plasma flows and electromagnetic fields into the energy of the plasma
particles, either as heat of the plasma species or as acceleration of a small
population of particles (Howes, 2017). The 2013 NRC Heliophysics Decadal
Survey (National Research Council, 2013) identifies plasma turbulence as a
ubiquitous phenomenon occurring both within the heliosphere and throughout the
Universe. Predicting the heating or acceleration of the different plasma
species by turbulence, based on the observable turbulence and plasma
parameters at large scales, constitutes one of the grand challenge problems in
heliophysics and astrophysics.

Turbulent plasma heating plays a key role in governing the flow of energy in
the heliosphere, impacting the mesoscale and macroscale evolution of key
heliospheric environments comprising the coupled solar-terrestrial system and
connecting the solar corona, solar wind, and Earth’s magnetosphere. The
accuracy of ongoing efforts (Mikić et al., 1999; Lionello et al., 2009; van
der Holst et al., 2014; Adhikari et al., 2020) to model globally the flow of
energy from the Sun through the interplanetary medium to the Earth, to the
other planets, and on to the boundary of the heliosphere with the surrounding
interstellar medium (Opher et al., 2020) relies on the availability of
prescriptive models of the turbulent plasma heating (Howes, 2010; Chandran et
al., 2011; Rowan et al., 2017).

Such turbulent heating prescriptions are also critical for the interpretation
of remote astronomical observations of emissions from black hole accretion
disks, such as the ground-breaking observations by the Event Horizon Telescope
of the supermassive black holes at the center of M87 (Event Horizon Telescope
Collaboration et al., 2019) and at Sagittarius A∗ in the Milky Way (Event
Horizon Telescope Collaboration et al., 2022), In these cases, alternative
turbulent heating prescriptions (Howes, 2010; Rowan et al., 2017) have been
shown to yield drastically different predictions for the emitted radiation
(Chael et al., 2019).

###  1.1 How Does Turbulence Mediate Energy Transport and Plasma Heating?

Figure 1:  (a) Schematic diagram of the turbulent magnetic-energy spectrum in
the solar wind, depicting the local transfer of energy from large scales
(small wavenumbers $k$) to small scales (large $k$) through a turbulent
cascade (black arrows). Instabilities may alter this energy flow by nonlocally
transporting energy directly to small, kinetic scales (red arrow). (b) In
wavevector space, the local transfer of the Alfvénic turbulent cascade follows
the anisotropic, critically balanced path, here plotted quantitatively for the
case of the _GS95_ scaling theory.

The transport of energy mediated by plasma turbulence is depicted in Fig. 1,
where (a) a plot of the turbulent magnetic energy wavenumber spectrum shows
that energy injected into the turbulent cascade at the outer scale
($k_{0}\rho_{i}=10^{-4}$) is transferred nonlinearly without loss through the
inertial range to the transition into the dissipation range at ion scales
($k\rho_{i}\sim 1$). At the ion scales and smaller ($k\rho_{i}\gtrsim 1$),
poorly understood mechanisms can remove energy from the turbulence,
transferring energy to the ions (blue arrows) or electrons (magenta arrows).
In Fig. 1(b), modern plasma turbulence theory [(Goldreich & Sridhar, 1995,
hereafter _GS95_) or (Boldyrev, 2006, hereafter _B06_)] predicts that the
turbulent cascade generates a scale-dependent anisotropy in wavevector space
(blue shaded region); isotropy at the driving scale ($k_{\perp
0}\rho_{i}=k_{\parallel 0}\rho_{i}=10^{-4}$) develops into a significant
anisotropy $k_{\perp}\gg k_{\parallel}$ as the cascade transfers energy to
smaller scales111 _Parallel_ and _perpendicular_ here are defined relative to
the local mean magnetic field direction.. In addition, kinetic instabilities
can generate unstable turbulent fluctuations with wavevectors $k\rho_{i}\sim
1$ (red ovals in both panels), which can mediate the nonlocal transfer of
energy from the large-scale motions directly to the ion kinetic scales. I
explicitly define “turbulence” here as the collection of physical mechanisms
involved in the conversion of the energy of large-scale plasma flows and
magnetic fields into plasma heat or the energy of accelerated particles,
including the nonlinear scale-to-scale energy transfer of the turbulent
cascade, the mechanisms damping the turbulence at small scales, and any
possible instabilities that may generate turbulent fluctuations or govern
nonlocal energy transfer.

###  1.2 What Do We Want to Predict?

A long-term goal of the heliophysics community is to characterize the
partitioning of turbulent energy in terms of observable plasma and turbulence
parameters, producing a predictive heating model that can be used in numerical
modeling of the global evolution of the heliosphere, thereby advancing our
capability for predictive modeling of the coupled solar-terrestrial system.
Understanding how the multiscale physics of plasma turbulence couples the
large-scale plasma conditions and evolution to the microphysical heating and
particle energization is essential to creating predictive heating models.
Unfortunately, current numerical codes cannot simulate the vast dynamic range
of scales inherent to the multiscale nature of heliospheric plasma turbulence
while simultaneously capturing the kinetic physics governing the conversion of
turbulent energy to some form of particle energy. Developing turbulent heating
models that capture this coupling is one promising avenue for formulating
predictive models at a computationally feasible cost. To develop such models,
it is essential to establish a thorough understanding of the microphysical
kinetic processes that govern the energization of particles.

In a statistically steady-state turbulent cascade in a weakly collisional
astrophysical plasma, the energy injected into the turbulence at the outer
scale is transferred locally in scale down to increasingly smaller scales
until reaching the kinetic length scales at which damping mechanisms can
collisonlessly remove energy from the turbulence and energize the different
particle species, as depicted in Fig. 1(a). This local cascade of energy is
mediated by nonlinear interactions among the turbulent fluctuations, such as
between counterpropagating Alfvén wave packets, often denoted “Alfvén wave
collisions” (Kraichnan, 1965; Sridhar & Goldreich, 1994; Goldreich & Sridhar,
1995; Ng & Bhattacharjee, 1996; Galtier et al., 2000; Howes et al., 2012b;
Howes & Nielson, 2013; Nielson et al., 2013; Howes et al., 2013; Drake et al.,
2013; Howes, 2016; Verniero et al., 2018; Verniero & Howes, 2018; TenBarge et
al., 2021; Ripperda et al., 2021). Nonlocal energy transfer mechanisms—such as
kinetic instabilities driven by the large-scale dynamics that can cause the
plasma conditions to deviate from local thermodynamic equilibrium—can also
directly generate turbulent fluctuations at the small kinetic length scales,
as illustrated in Fig. 1(b). Together, these mechanisms lead to a turbulent
energy density cascade rate $\varepsilon$ from the large scales down to the
small scales of the turbulence. This nonlinear transfer of the turbulent
energy density is ultimately damped at small scales by kinetic plasma physics
mechanisms, converting the energy density of the small-scale turbulent plasma
flows and magnetic fields into particle energy with a statistically steady-
state total plasma energy density heating rate $Q\sim\varepsilon$.

Quantity | Symbol | Bounded
---|---|---
Turbulent energy density cascade rate | $\varepsilon$ |
Plasma energy density heating rate | $Q$ |
Ion-to-Electron heating ratio | $Q_{i}/Q_{e}$ | $Q_{i}/Q$
Anisotropic ion energy density heating ratio | $Q_{\perp,i}/Q_{\parallel,i}$ | $Q_{\perp,i}/Q_{i}$
Anisotropic electron energy density heating ratio | $Q_{\perp,e}/Q_{\parallel,e}$ | $Q_{\perp,e}/Q_{e}$
Table 1: The key input quantities (top section) and deliverables (bottom
section) for predictive turbulent heating models. Bounded quantities are
normalized such that their values vary only over the interval $[0,1]$.

For a single-ion-species plasma, such as a fully ionized hydrogenic plasma of
protons and electrons, the primary deliverable of a turbulent heating model is
the determination of the partitioning of the energy by the kinetic turbulent
damping mechanisms between the ions and electrons, $Q_{i}/Q_{e}$. Unlike in
the case of strongly collisional plasmas, in which the dissipated turbulent
energy would be rapidly equilibrated between ion and electron species by
collisions, under the weakly collisional plasma conditions typical of space
and astrophysical plasmas, the differential heating of the plasma species can
impact the thermodynamic evolution of the plasma on mesoscales and
macroscales.

Furthermore, energy transferred to a given plasma species may lead to
energization of the parallel degree of freedom or the perpendicular degrees of
freedom relative to the direction of the local magnetic field. Thus, secondary
deliverables of a turbulent heating model are measures of the anisotropic ion
energy density heating rate $Q_{\perp,i}/Q_{\parallel,i}$ and the anisotropic
electron energy density heating rate $Q_{\perp,e}/Q_{\parallel,e}$. In a
strongly collisional plasma, the anisotropic heating for each species would be
rapidly isotropized by collisions, but weakly collisional plasmas allow the
parallel and perpendicular222Note that, under typical astrophysical plasma
conditions in which the turbulent fluctuation frequencies are not much faster
than the species cyclotron frequencies, $\omega\lesssim\Omega_{s}$, the two
perpendicular degrees of freedom are rapidly equilibrated by the gyromotion of
the particles about the magnetic field, leading to gyrotropic plasma
conditions with $T_{\perp,1}\sim T_{\perp,2}$. degrees of freedom to be
energized independently (Kawazura et al., 2019).

Different kinetic particle energization mechanisms may energize different
degrees of freedom— _e.g._ , Landau damping energizes particles parallel to
the magnetic field, but cyclotron damping energizes particles perpendicular to
the magnetic field—so determining $Q_{\perp,i}/Q_{\parallel,i}$ and
$Q_{\perp,e}/Q_{\parallel,e}$ first requires one to identify which kinetic
physical mechanisms govern the damping of the turbulent fluctuations as a
function of the plasma and turbulence parameters. Therefore, identifying the
contributing turbulent damping mechanisms as a function of the plasma and
turbulence parameters is an additional deliverable of a predictive turbulent
heating model. Such information about the dominant turbulent damping
mechanisms, as a function of the key dimensionless parameters describing the
turbulent plasma, can be concisely presented in a _phase diagram_ for the
dissipation of kinetic plasma turbulence, analogous to that developed for the
behavior of magnetic reconnection in astrophysical plasmas (Ji & Daughton,
2011). A first attempt at a phase diagram for the turbulent damping mechanisms
in weakly collisional space and astrophysical plasmas is presented here in
Fig. 11.

The key input quantities for a turbulent heating model ($\varepsilon$, $Q$)
and the deliverables ($Q_{i}/Q_{e}$, $Q_{\perp,i}/Q_{\parallel,i}$,
$Q_{\perp,e}/Q_{\parallel,e}$) for a turbulent heating model are summarized in
Table 1. Note that for ion cyclotron damping, which energizes only the
perpendicular degrees of freedom in the ion velocity distribution function,
the anisotropic ion heating ratio $Q_{\perp,i}/Q_{\parallel,i}$ would be
formally infinite. To avoid these inelegant infinities, we may choose to use
alternative _bounded_ versions of the heating ratios in the lower part of
Table 1 that vary only over the range $[0,1]$. These bounded quantities are
related to the original quantities by
$Q_{i}/Q=(Q_{i}/Q_{e})/[1+(Q_{i}/Q_{e})]$ for the species partition and
$Q_{\perp,s}/Q_{s}=(Q_{\perp,s}/Q_{\parallel,s})/[1+(Q_{\perp,s}/Q_{\parallel,s})]$
for the anisotropic heating ratio for each species.

It is worthwhile noting that the determination of the partitioning of the
turbulent energy between ions and electrons $Q_{i}/Q_{e}$ can be simplified by
the properties of the turbulent cascade that are evident in Fig. 1(a). The
inertial range is defined as the range of scales over which the turbulent
energy density is transferred to smaller scales without loss, leading to a
constant energy density cascade rate $\varepsilon$ throughout this
range333This standard definition of the inertial range from fluid turbulence
theory neglects the possibility of a nonlocal transfer of energy that may
arise in weakly collisional plasmas.. Therefore, the turbulent cascade rate
$\varepsilon$ at the onset of damping at ion kinetic scales at $k\rho_{i}\sim
1$ is the same as that determined at the outer scale of the turbulence at
$k_{0}\rho_{i}\ll 1$. The collisionless wave-particle interactions that serve
to remove energy from the turbulent fluctuations and transfer that energy to
the _ions_ generally play a significant role only within the range of length
scales similar to the ion kinetic length scales, $k\rho_{i}\sim 1$: (i) at
larger scales $k\rho_{i}\ll 1$, the collisionless damping rate with the ions
$\gamma_{i}$ is typically negligible in comparison to the turbulent cascade
rate $\gamma_{i}\ll\omega_{nl}$; and (ii) at smaller scales $k\rho_{i}\gg 1$,
the gyromotion of the ions averages out the effect of the small-scale
electromagnetic fluctuations, so the ions decouple from the turbulent fields
and do not exchange significant energy. Whatever turbulent energy manages to
cascade beyond the ion scales to smaller scales at $k\rho_{i}\gg 1$ will
ultimately be deposited with the electrons. Thus, the task of determining
$Q_{i}/Q_{e}$ boils down to determining what fraction of the turbulent energy
is damped on the ions as the cascade progresses through the ion kinetic scales
at $k\rho_{i}\sim 1$.

Ultimately, the kinetic mechanisms that serve to damp the turbulent
fluctuations at kinetic length scales depend on the nature of the turbulent
fluctuations at those scales. Therefore, a viable strategy for developing a
useful phase diagram for plasma turbulence is to use turbulence scaling
theories (_e.g._ , _GS95_ or _B06_) to predict the properties of the turbulent
fluctuations upon reaching the ion and electron kinetic scales from the plasma
and turbulence parameters at large scales. The primary aim of this
investigation is to establish a theoretical framework upon which to base such
predictions of the dominant turbulent damping mechanisms as a function of the
plasma and turbulence parameters.

###  1.3 What is the Approach for Developing a Predictive Turbulent Heating
Model?

The Buckingham Pi Theorem (Buckingham, 1914) uses dimensional analysis
(Barenblatt, 1996) to determine the minimum number of dimensionless parameters
upon which the physical behavior of a system depends. This approach is the key
principle underlying the development of predictive models of turbulent plasma
heating in terms of the fundamental dimensionless parameters of kinetic plasma
turbulence. For sufficiently simple systems, in which there are not multiple
physical quantities of importance that have the same dimensions, one can
define a unique set of fundamental dimensionless parameters. But for more
complicated systems, in which there may be more than one relevant physical
quantity with the same dimensions, such as multiple characteristic length
scales, the determination of the fundamental dimensionless parameters is
generally non-unique. The case of kinetic turbulence in a weakly collisional
plasma certainly falls into this more complicated category. In that situation,
it is essential to use physical insight into the system’s properties and
behavior to define a particular set of dimensionless parameters that is most
useful in characterizing the behavior of the system.

The goal of this paper is to define a particular set of fundamental
dimensionless parameters that can be used to characterize the properties of
kinetic plasma turbulence and its damping mechanisms. Modern scaling theories
for anisotropic plasma turbulence can then be expressed in terms of these
fundamental parameters to determine the properties of the turbulent
fluctuations at the kinetic length scales at which kinetic physical mechanisms
can act to remove energy from the turbulent fluctuations and energize the
plasma particles. How those kinetic damping mechanisms depend on the plasma
parameters and properties of the turbulent fluctuations can then be used to
predict the dominant mechanisms for the damping of turbulence as a function of
the dimensionless parameters. All of this information can then be synthesized
to generate improved turbulent heating prescriptions for comparison to direct
numerical simulations of turbulent dissipation and for use in global modeling
efforts for space and astrophysical plasma systems, such as the heliosphere or
supermassive black hole accretion disks.

##  2 The Fundamental Dimensionless Parameters of Kinetic Plasma Turbulence

The first step in developing a predictive model for plasma turbulence in
weakly collisional space and astrophysical plasmas is to identify the key
dimensionless parameters upon which the turbulent energy cascade and its
dissipation depend. In this paper, the focus is to characterize kinetic plasma
turbulence in the non-relativistic limit. There are two categories of
governing parameters: _Plasma Parameters_ and _Turbulence Parameters_. Modern
scaling theories for anisotropic plasma turbulence can be employed to reduce
the number of parameters needed to characterize the turbulence, as described
in §3.1. The reduced set of dimensionless parameters can then be used to
develop useful phase diagrams that identify the dominant physical mechanisms
of turbulent dissipation and which of these key parameters has the strongest
impact on how particles are energized by the turbulence.

###  2.1 Plasma Parameters

The _Plasma Parameters_ play a critical role in determining the nature of the
turbulent fluctuations and their interaction with the underlying particle
velocity distributions. For simplicity of notation, we assume here a fully
ionized, hydrogenic plasma of protons and electrons with bi-Maxwellian
equilibrium velocity distributions characterized by separate parallel and
perpendicular temperatures444Note that we absorb the Boltzmann constant into
the temperature throughout this paper, giving temperature in units of energy.
for each species, $T_{\parallel,s}$ and $T_{\perp,s}$. In this idealized case,
the equilibrium ion and electron number densities are equal, $n_{i}=n_{e}$.

The first and most important plasma parameter is the parallel ion plasma beta,
$\beta_{\parallel,i}=8\pi n_{i}T_{\parallel i}/B^{2}$, which characterizes the
ratio of parallel thermal pressure to the magnetic pressure in the plasma and
dominantly controls the phase speeds of different linear wave modes in a
magnetized plasma. Note that with definitions for the parallel ion thermal
velocity $v_{t\parallel,i}^{2}=2T_{\parallel,i}/m_{i}$ and the Alfvén velocity
$v_{A}^{2}=B^{2}/(4\pi n_{i}m_{i})$, the parallel ion plasma beta can be
alternatively expressed as
$\beta_{\parallel,i}=v_{t\parallel,i}^{2}/v_{A}^{2}$. The second parameter is
the ion-to-electron parallel temperature ratio
$\tau_{\parallel}=T_{\parallel,i}/T_{\parallel,e}$. The third and forth
parameters describe the ion and electron temperature anisotropies,
$A_{i}=T_{\perp,i}/T_{\parallel,i}$ and $A_{e}=T_{\perp,e}/T_{\parallel,e}$.
The fifth parameter characterizes the collisionality of the plasma by
comparing the mean free path scale for electron-electron collisions relative
to the parallel driving scale $k_{\parallel 0}\lambda_{mfp,e}$.

Although for a hydrogenic plasma of protons and electrons the ion-to-electron
mass ratio $\mu=m_{i}/m_{e}$ has a fixed physical value of 1836, this ratio is
included as an additional independent parameter in the calculations below
since a reduced mass ratio is often used in numerical simulations of plasma
turbulence for reasons of computational efficiency. Note also that, once $\mu$
and $\tau$ are specified, the single dimensionless parameter $k_{\parallel
0}\lambda_{mfp,e}$ is sufficient to specify all possible charged particle
collision rates, $\nu_{ee}$, $\nu_{ei}$, $\nu_{ii}$, and $\nu_{ie}$.

Therefore, $(\beta_{\parallel,i},\tau_{\parallel},A_{i},A_{e},k_{\parallel
0}\lambda_{mfp,e})$ are five key dimensionless plasma parameters that govern
the properties of plasma turbulence, along with the fixed mass ratio $\mu$.
These Plasma Parameters are summarized in Table 2.

Plasma Parameters
---
Parameter | Symbol and Definition
Parallel Ion Plasma Beta | $\beta_{\parallel,i}=8\pi n_{i}T_{\parallel i}/B^{2}$
Ion-to-Electron Parallel Temperature Ratio | $\tau_{\parallel}=T_{\parallel,i}/T_{\parallel,e}$
Ion Temperature Anisotropy | $A_{i}=T_{\perp,i}/T_{\parallel,i}$
Electron Temperature Anisotropy | $A_{e}=T_{\perp,e}/T_{\parallel,e}$
Plasma Collisionality | $k_{\parallel 0}\lambda_{mfp,e}$
Ion-to-Electron Mass Ratio | $\mu=m_{i}/m_{e}$
Turbulence Parameters
Parameter | Symbol and Definition
Perpendicular Driving Scale | $k_{\perp 0}\rho_{i}$
Driving Scale Anisotropy | $k_{\parallel 0}/k_{\perp 0}$
Nonlinear Parameter | $\chi_{0}=k_{\perp 0}\delta B_{\perp 0}\theta_{0}/(k_{\parallel 0}B_{0})$
Imbalance Parameter | $Z_{0}^{+}/Z_{0}^{-}$
Driving Compressibility | $E_{comp}/E_{inc}$
Table 2: For the general case of homogeneous turbulence in a fully ionized,
hydrogenic plasma with bi-Maxwellian equilibrium velocity distributions, these
are the fundamental dimensionless parameters of kinetic plasma turbulence,
including Plasma Parameters and Turbulence Parameters. Note that the Boltzmann
constant is absorbed to give temperature in units of energy.

Note that if the plasma contains multiple ion species, possibly each with
multiple charge states, the number of parameters increases dramatically,
requiring a minimum of five new dimensionless parameters for each new species
$s$: (i) mass ratio $m_{s}/m_{i}$, (ii) number density ratio $n_{s}/n_{e}$,
(iii) charge ratio $q_{s}/q_{e}$, (iv) parallel temperature ratio
$T_{\parallel,s}/T_{\parallel,i}$, and (v) temperature anisotropy
$T_{\perp,s}/T_{\parallel,s}$. Partially ionized plasmas will similarly
require additional dimensionless parameters to be characterized. Since the
intention here is _not_ to treat kinetic plasma turbulence in full generality,
but rather to describe how one can use the Pi Theorem and turbulence scaling
theories to minimize the number of relevant dimensionless parameters, the
remainder of this paper will focus on the case of a fully ionized, hydrogenic
plasma of protons and electrons with bi-Maxwellian equilibrium velocity
distributions. Refinements of turbulent heating models to include minor ions,
such as $\alpha$-particles that comprise a small fraction of ions in space and
astrophysical plasmas, can be constructed following the same principles
presented here.

###  2.2 Turbulence Parameters

The _Turbulence Parameters_ characterize the scale, amplitude, and properties
of the turbulent driving though five additional parameters. For simplicity,
the focus here is restricted to only the case of turbulence in a spatially
homogeneous system in which the scale length of gradients in the equilibrium
quantities are much larger than the outer scale of the system. Dimensionless
ratios of length scales are notationally simplified by characterizing the
length-scales of the turbulent fluctuations in terms of their perpendicular
and parallel wavenumbers. Below quantities at the outer scale of the
turbulence, commonly denoted the driving scale or energy injection scale, are
denoted by the subscript “0”.

The driving of the turbulence can be characterized by three dimensionless
parameters. The first parameter is the perpendicular wavenumber of the
driving-scale fluctuations normalized by the ion Larmor radius, $k_{\perp
0}\rho_{i}$. The second parameter is the wavevector anisotropy of the
fluctuations at the driving scale, $k_{\parallel 0}/k_{\perp 0}$. The third
parameter is the the nonlinearity parameter that characterizes the amplitude
of the driving, $\chi_{0}=(k_{\perp 0}\delta B_{\perp
0}\theta_{0})/(k_{\parallel 0}B_{0})$, where $\delta B_{\perp 0}$ is the
amplitude of the perpendicular magnetic field fluctuations, $B_{0}$ is the
equilibrium magnetic field magnitude, and $\theta_{0}$ is the angle of dynamic
alignment at the driving scale in radians (Boldyrev, 2006). In the next
section, we demonstrate how modern scaling theories for anisotropic MHD
turbulence can be used to combine these three dimensionless parameters into a
single parameter, the _isotropic driving wavenumber_ , $k_{0}\rho_{i}$. Note
also that, in principle, the dynamic alignment angle $\theta_{0}$ at driving
scale is an additional dimensionless parameter, but dynamic alignment is
typically understood to be a phenomenon that develops self-consistently due to
the turbulence as it cascades to smaller scales, so we neglect this parameter
at the outerscale, generally taking $\theta_{0}\sim 1$.

The nature of the turbulent fluctuations at the driving scale also plays an
important role in determining the ultimate fate of the turbulent energy. The
fourth turbulence parameter is the imbalance parameter, measuring the ratio of
the amplitudes in upward-to-downward Elsasser fields, $Z_{0}^{+}/Z_{0}^{-}$,
where the Elsasser fields at the outer scale are defined by
$Z_{0}^{\pm}=\delta\mathbf{U}_{\perp,0}\pm\delta B_{\perp,0}/\sqrt{4\pi
n_{i}m_{i}}$ (Elsasser, 1950), and where $\delta\mathbf{U}_{\perp,0}$ and
$\delta B_{\perp,0}$ are the perpendicular perturbed plasma flow velocity and
perpendicular perturbed magnetic field of the turblent fluctuations at the
driving scale. Physically, Alfvén wave packets propagating up the mean
magnetic field have $Z^{+}=0$ and $Z^{-}\neq 0$, and those propagating down
the field have $Z^{+}\neq 0$ and $Z^{-}=0$. The Elsasser fields effectively
characterize the wave energy flux in each direction along the magnetic
field555Note that the imbalance of counterpropagating wave energy fluxes in
plasma turbulence has alternatively been characterized in the literature using
the cross helicity (Dobrowolny et al., 1980; Marsch & Tu, 1990; Matthaeus et
al., 1982; Perez & Boldyrev, 2009, 2010) or the normalized Poynting flux (Hnat
et al., 2003; Di Mare & Howes, 2024). of the typically dominant Alfvénic
fluctuations in heliospheric plasma turbulence (Tu & Marsch, 1995;
Schekochihin et al., 2009; Bruno & Carbone, 2013). Since the nonlinear
interactions that mediate the turbulent cascade in Alfvénic turbulence at MHD
scales arise only between counterpropagating and perpendicularly polarized
Alfvén waves (Kraichnan, 1965; Sridhar & Goldreich, 1994; Goldreich & Sridhar,
1995; Ng & Bhattacharjee, 1996; Galtier et al., 2000; Howes & Nielson, 2013;
Howes:2014c), the imbalance plays an important role in characterizing the
turbulence. The importance of $Z_{0}^{+}/Z_{0}^{-}$ under sufficiently
imbalanced conditions, such as occurs in the inner heliosphere where wave
energy fluxes are dominantly anti-sunward, has been demonstrated clearly with
the recent discovery of new phenomenon known as the _helicity barrier_ in
plasma turbulence (Squire et al., 2022, 2023).

A fifth turbulence parameter is the ratio of the energy in compressible
fluctuations to that in incompressible fluctuations, $E_{comp}/E_{inc}$. Since
a magnetized plasma supports two distinct compressible wave modes, the fast
and slow magnetosonic waves (differing by whether the thermal and magnetic
pressure perturbations associated with the wave act in conjunction or in
opposition), one could further characterize the compressible component of the
turbulent fluctuations by the additional ratio of the energy in the fast
magnetosonic fluctuations to the total energy of compressible fluctuations
$E_{F}/E_{comp}$, where the total compressible energy is the sum of the
contrbutions from fast-mode and slow-mode fluctuations,
$E_{comp}=E_{F}+E_{s}$. The ratio $E_{F}/E_{comp}$ is not included here as one
of the fundamental parameters based on a study of the inertial range in solar
wind turbulence at 1 AU that showed negligible energy in fast-mode mode
fluctuations, or $E_{F}/E_{comp}\ll 1$ (Howes et al., 2012a). In significant
energy exists in fast-mode fluctuations in other turbulent astrophysical
plasmas, $E_{F}/E_{comp}$ may need to be included.

Therefore, $(k_{\perp 0}\rho_{i},k_{\parallel 0}/k_{\perp
0},\chi_{0},Z_{0}^{+}/Z_{0}^{-},E_{comp}/E_{inc})$ are five key dimensionless
turbulence parameters that govern the properties of plasma turbulence. These
Turbulence Parameters are summarized in Table 2.

If the restriction to spatially homogeneous turbulent plasmas is relaxed, the
large-scale context of the turbulent plasma can also impact the nonlinear
dynamics of the turbulence and the evolution of the turbulent fluctuations. In
the solar wind, the spherical expansion of the solar wind governs the
evolution of the plasma equilibrium temperatures with heliocentric radius,
possibly triggering kinetic temperature anisotropy instabilities that can
mediate nonlocal transfer of energy to small scales, as depicted in Fig. 1(b).
In rotating plasmas, such as found in the solar convection zone, the
turbulence may be driven with helical motions that can initiate the growth of
magnetic field energy through a dynamo effect (Parker, 1955; Glatzmaier, 1984,
1985; Parker, 1993; Ossendrijver, 2003). Differential rotation in
astrophysical accretion disks, which can trigger plasma turbulence through the
magnetorotational instability (Balbus & Hawley, 1991; Hawley & Balbus, 1991;
Balbus & Hawley, 1998) has been recently shown to generate turbulence with a
dominant fraction of its energy in compressible fluctuations (Kawazura et al.,
2022). Finally, kinetic instabilities arising in the foot and ramp regions of
collisionless shocks (Brown et al., 2023) can generate turbulent fluctuations
that are swept into the downstream region. In all of these astrophysically
relevant cases, additional dimensionless turbulence parameters will be
required to characterize properly the turbulent fluctuations and their
evolution.

##  3 Scaling Theories for Anisotropic Plasma Turbulence

The key to predicting the evolution of the turbulent fluctuations as energy is
cascaded to ever smaller scales, and therefore also to developing a useful
phase diagram for plasma turbulence, is to use turbulence scaling theories to
predict the properties of the turbulence upon reaching the ion kinetic scales
from the characteristics of the turbulence at large scales. In particular, we
will include the two prominent theories for the anisotropic cascade of
Alfvénic turbulence in MHD plasmas, the Goldreich-Sridhar 1995 (_GS95_)
formulation for weak and strong MHD turbulence (Sridhar & Goldreich, 1994;
Goldreich & Sridhar, 1995), and the modification of this theory by Boldyrev
(_B06_) that accounts for a weakening of the nonlinearity by dynamic alignment
(Boldyrev, 2006). Note that since these are MHD theories of plasma
turbulence—where strong collisionality is assumed in the standard MHD
approximation (Kulsrud, 1983; Gurnett & Bhattacharjee, 2017)—the equilibrium
species temperatures are assumed to be anisotropic and equal, corresponding to
$\tau_{\parallel}=A_{i}=A_{e}=1$. For the case of weakly collisional plasma
turbulence that is relevant to most space and astrophysical plasma
environments, these theories may need to be modified appropriately, but such
development of refined turbulence theories is beyond the scope of the
treatment here.

The particular scalings of the parallel wavenumber, intermediate wavenumber,
alignment angle, and magnetic field perturbation as a function of the
perpendicular wavenumber are given by

$k_{\parallel}\propto k_{\perp}^{2/(3+\alpha)}$ (1) $k_{i}\propto
k_{\perp}^{3/(3+\alpha)}$ (2) $\theta\propto k_{\perp}^{-\alpha/(3+\alpha)}$
(3) $\delta B_{\perp}\propto k_{\perp}^{-1/(3+\alpha)}$ (4)

where $\alpha=0$ corresponds to the _GS95_ theory and $\alpha=1$ corresponds
to the _B06_ theory. Here the turbulent fluctuations as a function of scale
are characterized by the three spatial components of the wave vector relative
to the local magnetic field direction: (i) the parallel wavenumber
$k_{\parallel}$; and (ii) the intermediate wavenumber $k_{i}$ and (iii) the
perpendicular wavenumber $k_{\perp}$, which together describe the anisotropy
of the fluctuations in the plane perpendicular to the local magnetic field,
$k_{i}\ll k_{\perp}$, arising from the dynamic alignment of the fluctuations
(Boldyrev, 2006; Mason et al., 2006; Boldyrev et al., 2009; Perez et al.,
2012) and predicting the development of current sheets at small scales. These
scaling theories predict the scaling of $k_{\parallel}$, $k_{i}$, the angle of
dynamic alignment $\theta$, and the amplitude of the perpendicular magnetic
field fluctuations $\delta B_{\perp}$ as a function of the perpendicular
wavenumber $k_{\perp}$.

These theories yield specific predictions for the perpendicular magnetic
energy spectrum of turbulence in the MHD regime, $k_{\perp}\rho_{i}\gg 1$, and
for the scale-dependent anisotropy of the turbulence, a qualitative feature
that is well-established in MHD turbulence (Cho & Vishniac, 2000; Maron &
Goldreich, 2001; Cho et al., 2002; Cho & Lazarian, 2003; Mason et al., 2006;
Boldyrev et al., 2009; Narita et al., 2010; Sahraoui et al., 2010; Perez et
al., 2012; Roberts et al., 2015). The perpendicular magnetic energy spectrum
for the two theories is given by

$E_{B_{\perp}}(k_{\perp})\propto\frac{(\delta
B_{\perp})^{2}}{k_{\perp}}\propto
k_{\perp}^{-(5+\alpha)/(3+\alpha)}=\left\\{\begin{array}[]{cc}k_{\perp}^{-5/3}&\alpha=0\\\
k_{\perp}^{-3/2}&\alpha=1\\\ \end{array}\right..$ (5)

The parallel-to-perpendicular wavevector anisotropy is given by

$\frac{k_{\parallel}}{k_{\perp}}\propto
k_{\perp}^{-(1+\alpha)/(3+\alpha)}=\left\\{\begin{array}[]{cc}k_{\perp}^{-1/3}&\alpha=0\\\
k_{\perp}^{-1/2}&\alpha=1\\\ \end{array}\right.,$ (6)

demonstrating that both theories predict that the turbulent fluctuations
develop a scale-dependent anisotropy, whereby the fluctuations become
progressively more elongated along the magnetic field as they cascade to
smaller scale, eventually yielding $k_{\parallel}/k_{\perp}\ll 1$ at
sufficiently small scales. This $k_{\parallel}/k_{\perp}$ anisotropy follows
from the imposing of a critical balance between the linear and nonlinear
timescales of the turbulence (Goldreich & Sridhar, 1995; Mallet et al., 2015).
The intermediate-to-perpendicular wavevector anisotropy is given by

$\frac{k_{i}}{k_{\perp}}\propto
k_{\perp}^{-\alpha/(3+\alpha)}=\left\\{\begin{array}[]{cc}k_{\perp}^{0}&\alpha=0\\\
k_{\perp}^{-1/4}&\alpha=1\\\ \end{array}\right..$ (7)

Here, the _GS95_ theory predicts no development of anisotropy in the
perpendicular plane, at odds with the findings of current sheets at small
perpendicular scales in numerical simulations (Matthaeus & Montgomery, 1980;
Meneguzzi et al., 1981; Uritsky et al., 2010; Zhdankin et al., 2012, 2013) and
spacecraft observations (Borovsky, 2008; Osman et al., 2011). The _B06_ theory
predicts the development of anisotropy in the perpendicular plane,
corresponding to the development of current sheets with $k_{i}/k_{\perp}\ll 1$
at sufficiently small scales, potentially triggering tearing instabilities
(Furth et al., 1963; Zocco & Schekochihin, 2011) that enable magnetic
reconnection to disrupt the turbulent cascade (Boldyrev & Loureiro, 2017;
Loureiro & Boldyrev, 2017b, a; Mallet et al., 2017b, a; Walker et al., 2018).

It is worthwhile emphasizing that the idea of critical balance in the
turbulent cascade does _not_ imply that all of the turbulent energy is
concentrated at the parallel wavenumber $k_{\parallel}^{cb}$ determined by
critical balance. Rather, $k_{\parallel}^{cb}$ is the upper limit of parallel
wavenumber, and thus generally also implies an upper limit on the linear wave
frequency (since $\omega\propto k_{\parallel}$ for $\omega\ll\Omega_{i}$).
Rather, the turbulent fluctuation power extends over the range $0\leq
k_{\parallel}\leq k_{\parallel}^{cb}$ (Goldreich & Sridhar, 1995; Maron &
Goldreich, 2001). Accounting for how the turbulent energy is distributed over
this $k_{\parallel}$ range would lead to order unity or less quantitative
changes in the heating rates. Since the goal of predictive turbulent heating
models is generally to calculate the lowest-order, order-of-magnitude
quantities for species heating rates, any effects due to the distribution of
turbulent energy over $0\leq k_{\parallel}\leq k_{\parallel}^{cb}$ are likely
to yield only higher-order corrections to these predictions; thus, the effect
of this range of turbulent energy in $k_{\parallel}$ will be neglected in the
first generation of turbulent heating models addressed here.

In closing, note that these theories for MHD turbulence are largely focused on
the transport of energy in Alfvénic fluctuations, so if significant energy is
transferred nonlinearly into other (non-Alfvénic ) wave modes, this may
significantly alter the energy transport from the predictions shown here, with
potentially qualitatively different predictions for the damping mechanisms of
these non-Alfvénic turbulent fluctuations.

###  3.1 Exploiting Turbulent Scaling Theories to Reduce the Number of
Fundamental Dimensionless Parameters

Given the large number of fundamental dimensionless parameters that
characterize weakly collisional plasma turbulence, shown in Table 2, it is
imperative to reduce the number of important parameters to make progress in
developing predictive models of turbulent heating. Adopting idealized
approximations, such as maintaining isotropic equilbrium velocity
distributions for the plasma, enables one to eliminate some of the parameters
in Table 2. Another approach is to exploit the MHD turbulent scaling theories
described in §3 to combine several parameters into a new single parameter, an
approach presented here.

Figure 2: Diagram illustrating how to determine the isotropic driving
wavenumber $k_{0}\rho_{i}$ in terms of the perpendicular driving wavenumber
$k_{\perp 0}\rho_{i}$, the driving anisotropy $k_{\parallel 0}/k_{\perp 0}$,
and the nonlinearity parameter of the driving, $\chi_{0}=k_{\perp 0}\delta
B_{\perp}\theta_{0}/(k_{\parallel 0}B_{0})$.

In general, the driving of the plasma turbulence depends on at least three
independent parameters: (i) the perpendicular wavenumber at which the
turbulence is driven $k_{\perp 0}\rho_{i}$, (ii) the anisotropy of the driven
fluctuations with respect to the local magnetic field direction $k_{\parallel
0}/k_{\perp 0}$, and (iii) the nonlinearity parameter which characterizes the
strength of the driving, $\chi_{0}\equiv k_{\perp 0}\delta B_{\perp
0}\theta_{0}/(k_{\parallel 0}B_{0})$. Recall that quantities evaluated at the
driving scale are given the subscript “0”. Fortunately, for turbulence driven
at large scales (relative to the ion kinetic scales), as is typical in most
space or astrophysical plasma environments, MHD scaling theory can be used to
predict the properties of the turbulence upon reaching the kinetic scales,
where kinetic plasma processes act to dissipate the turbulence. In fact, as
detailed here, it is possible to combine these three dimensionless parameters
into a new, single dimensionless parameter: the _isotropic driving wavenumber_
, $k_{0}\rho_{i}$.

MHD turbulent scaling theories can be used, for a given set of driving
parameters ($k_{\perp 0}\rho_{i}$, $k_{\parallel 0}/k_{\perp 0}$, $\chi_{0}$),
to predict the properties of the turbulent fluctuations (_e.g._ , anisotropy,
amplitude) upon reaching the ion kinetic scales, here taken to be the
perpendicular wavenumber $k_{\perp}\rho_{i}=1$. The isotropic driving
wavenumber $k_{0}\rho_{i}$ is the unique wavenumber that produces the same
conditions at $k_{\perp}\rho_{i}=1$ as if the turbulence were driven
isotropically at that wavenumber with $k_{\perp 0}=k_{\parallel 0}\equiv
k_{0}$ and in a condition of critical balance $\chi_{0}=1$.

The key procedure for determining the isotropic driving wavenumber
$k_{0}\rho_{i}$ is illustrated in Fig. 2. Panel (a) shows the perpendicular
wavenumber magnetic energy spectrum $E_{B_{\perp}}(k_{\perp})$ driven weakly
with $\chi_{0}=0.1$ at perpendicular wavenumber $k_{\perp 0}\rho_{i}=10^{-5}$,
where a weak turbulent cascade leads to a steep weak MHD turbulence spectrum
with $E_{B_{\perp}}\propto k_{\perp}^{-2}$ (Sridhar & Goldreich, 1994; Galtier
et al., 2000; Boldyrev & Perez, 2009), followed by a flatter strong MHD
turbulence spectrum with $E_{B_{\perp}}\propto
k_{\perp}^{-(5+\alpha)/(3+\alpha)}$ (Goldreich & Sridhar, 1995; Boldyrev,
2006). The cascade reaches the perpendicular ion kinetic scale
$k_{\perp}\rho_{i}=1$, and at higher wavenumbers the cascade continues as
strong kinetic Alfvén wave (KAW) turbulence with a significantly steeper
magnetic energy spectrum with $E_{B_{\perp}}\propto k_{\perp}^{-2.8}$
(Sahraoui et al., 2009; Kiyani et al., 2009; Alexandrova et al., 2009; Chen et
al., 2010; Howes et al., 2011; Sahraoui et al., 2013). When the collisionless
damping of the KAW turbulent cascade begins to diminish the rate of nonlinear
energy cascade, the spectrum begins to fall off exponentially in the weakly
dissipating KAW turbulence (WDKT) regime (Howes et al., 2011), ultimately
terminating the turbulent cascade at scales $k_{\perp}\rho_{e}\gtrsim 1$.

In Fig. 2(b) is presented the characteristic path of the anisotropic cascade
of energy through $(k_{\perp},k_{\parallel})$ wavevector space for the
different regimes of the turbulent cascade. There is no parallel cascade in
the weak MHD turbulence regime (Sridhar & Goldreich, 1994; Galtier et al.,
2000), followed by a scale-dependent anisotropy given by $k_{\parallel}\propto
k_{\perp}^{2/(3+\alpha)}$ in the strong MHD turbulence regime (Goldreich &
Sridhar, 1995; Boldyrev, 2006). The turbulence transitions to the strong KAW
regime at $k_{\perp}\rho_{i}\sim 1$, with the anisotropy becoming stronger at
$k_{\perp}\rho_{i}>1$, scaling as $k_{\parallel}\propto k_{\perp}^{1/3}$ (Cho
& Lazarian, 2004; Howes et al., 2008a, 2011). As the turbulence weakens due to
collisionless damping of turbulent fluctuations in the WDKT regime, the
parallel cascade ceases once again (Howes et al., 2011).

The diagrams in Fig. 2 provide guidance to understand how MHD turbulence
scaling theories can be used to combine ($k_{\perp 0}\rho_{i}$, $k_{\parallel
0}/k_{\perp 0}$, $\chi_{0}$) into a single new parameter, the isotropic
driving wavenumber, $k_{0}\rho_{i}$. If the turbulence is driven at critical
balance ($\chi_{0}=1$) and isotropically ($k_{\parallel 0}/k_{\perp 0}=1$),
then the isotropic driving wavenumber is simply the same as the perpendicular
driving wavenumber, $k_{0}\rho_{i}=k_{\perp 0}\rho_{i}$. Therefore, it is
necessary to use MHD turbulence scaling theory to determine $k_{0}\rho_{i}$
only in the case that the turbulence is driven weakly $\chi_{0}<1$ or
anisotropically $k_{\parallel 0}/k_{\perp 0}\neq 1$.

For weakly driven turbulence with $\chi_{0}<1$, MHD turbulence theory predicts
that the parallel cascade is inhibited and that the perpendicular wavenumber
magnetic energy spectrum of the turbulence scales as
$E_{B_{\perp}}\propto(\delta B_{\perp})^{2}/k_{\perp}\propto k_{\perp}^{-2}$.
Therefore, in terms of the driving amplitude $\delta B_{\perp 0}$ and
perpendicular driving wavenumber $k_{\perp 0}$, the amplitude of the
turbulence as a function of $k_{\perp}$ is given by $\delta
B_{\perp}(k_{\perp})=\delta B_{\perp 0}(k_{\perp 0}/k_{\perp})^{1/2}$. Given
the lack of parallel cascade in weak MHD turbulence, we fix
$k_{\parallel}=k_{\parallel 0}$, and obtain the following scaling666Here we
have assumed that the dynamic alignment angle $\theta_{0}$ also does not
evolve as a function of $k_{\perp}$ in weak MHD turbulence, as there is no
prediction in the literature for the scaling of dynamic alignment in weak
turbulence. for the nonlinearity parameter as a function of $k_{\perp}$,

$\chi(k_{\perp})=\chi_{0}\left(\frac{k_{\perp}\rho_{i}}{k_{\perp
0}\rho_{i}}\right)^{1/2}.$ (8)

MHD turbulence theory predicts that the nonlinearity parameter will increase
through the weak turbulent cascade until it reaches a value of unity,
$\chi\sim 1$, at which point the turbulence has reached a state of critical
balance in the strong turbulence regime. From that point, theory predicts that
an anisotropic cascade will ensue (transferring energy to smaller scales more
rapidly in the perpendicular than the parallel direction), with the parallel
and perpendicular cascades balanced to maintain a nonlinearity parameter of
unity, $\chi\sim 1$. Therefore, the key point in the cascade is perpendicular
wavenumber at which the turbulence first reaches $\chi=1$, given by

$\left.k_{\perp}\rho_{i}\right|_{\chi=1}=\frac{1}{\chi_{0}^{2}}k_{\perp
0}\rho_{i}$ (9)

An example of this weakly turbulent, strictly perpendicular cascade in
$(k_{\perp},k_{\parallel})$ space is depicted in Fig. 2(b), where turbulence
driven weakly with $\chi_{0}=0.1$, cascades to $k_{\perp}\rho_{i}=100k_{\perp
0}\rho_{i}$, at which point the turbulence becomes strong, with $\chi=1$.

To obtain the isotropic driving wavenumber $k_{0}\rho_{i}$ from this point in
wavevector space $(\left.k_{\perp}\rho_{i}\right|_{\chi=1},k_{\parallel 0})$,
we need to extend the cascade backward along the critical balance line (solid
black and red dashed lines) until it intersects with isotropy
$k_{\parallel}=k_{\perp}$ (dotted line). From the MHD scaling theories in §3,
for a critically balanced, strong turbulent cascade that is driven strongly
and isotropically at $k_{0}\rho_{i}$, the anisotropy as a function of
$k_{\perp}$ scales as

$\frac{k_{\parallel}}{k_{\perp}}=\left(\frac{k_{0}\rho_{i}}{k_{\perp}\rho_{i}}\right)^{\frac{1+\alpha}{3+\alpha}}.$
(10)

The anisotropy at the point
$(\left.k_{\perp}\rho_{i}\right|_{\chi=1},k_{\parallel 0})$ is given by
$k_{\parallel}/k_{\perp}=(k_{\parallel 0}/k_{\perp 0})\chi_{0}^{2}$.
Substituting this expression for $k_{\parallel}/k_{\perp}$ on the left-hand
side of (10) and substituting (9) into the denominator on the right-hand side
of (10), we can solve for the isotropic driving wavenumber $k_{0}\rho_{i}$ ,
obtaining the result

$k_{0}\rho_{i}=k_{\perp 0}\rho_{i}\left(\frac{k_{\parallel 0}}{k_{\perp
0}}\right)^{\frac{3+\alpha}{1+\alpha}}\chi_{0}^{\frac{4}{1+\alpha}}$ (11)

In Fig. 2(b), the isotropic driving wavenumber $k_{0}\rho_{i}$ (red circle) is
the shown where the critical balance line (solid black and red dashed lines)
intersects the isotropy line, $k_{\parallel}=k_{\perp}$ (dotted line). Note
that this reduction of the three driving parameters to $k_{0}\rho_{i}$ is
valid only in the limit that the turbulent cascade becomes strong before it
reaches the ion kinetic scales at $k_{\perp}\rho_{i}\sim 1$.

The relation (11) provides the means to combine the three turbulence driving
parameters ($k_{\perp 0}\rho_{i}$, $k_{\parallel 0}/k_{\perp 0}$, $\chi_{0}$)
into the single new effective parameter, the isotropic driving wavenumber
$k_{0}\rho_{i}$. The properties of the turbulent fluctuations at perpendicular
scales smaller than the point where the turbulent cascade becomes strong, or
$k_{\perp}\rho_{i}>\left.k_{\perp}\rho_{i}\right|_{\chi=1}$, will be the same
for the original case with driving parameterized by ($k_{\perp 0}\rho_{i}$,
$k_{\parallel 0}/k_{\perp 0}$, $\chi_{0}$) as for turbulence driven strongly
and isotropically at $k_{0}\rho_{i}$. In the development of a turbulent
heating model, it is the properties of the turbulent fluctuations at the
small, kinetic length scales that will determine the kinetic damping mechanism
and the resulting particle energization, so reducing three turbulence driving
parameters to a single parameter substantially reduces the dimensionality of
the parameter space to be modeled.

It is worthwhile noting that the evolution of the dynamic alignment angle with
perpendicular wavenumber $\theta(k_{\perp})$ is not well modeled by
$k_{0}\rho_{i}$. The development of current sheets at very small scales due to
dynamic alignment, and the possibility for the interruption of the turbulent
cascade by magnetic reconnection (Boldyrev & Loureiro, 2017; Loureiro &
Boldyrev, 2017b, a; Mallet et al., 2017b, a; Walker et al., 2018), will depend
on this scaling of $\theta(k_{\perp})$. But, from a practical point of view,
this issue in modeling $\theta(k_{\perp})$ while reducing ($k_{\perp
0}\rho_{i}$, $k_{\parallel 0}/k_{\perp 0}$, $\chi_{0}$) to $k_{0}\rho_{i}$ is
likely to be a problem only in rare cases. First, the scaling of dynamic
alignment angle $\theta\propto k_{\perp}^{-1/4}$ in _B06_ is rather weak,
scaling as the $-1/4$ power, so only when the MHD inertial range spans many
orders of magnitude in the perpendicular length scale will the development of
current sheets likely have a significant impact the turbulent dynamics.
Second, dynamic alignment appears to arise only in strong MHD turbulence, so
it seems reasonable to take $\theta_{0}\sim 1$ radian (roughly equivalent to
isotropy in the perpendicular plane) at the perpendicular scale
$\left.k_{\perp}\rho_{i}\right|_{\chi=1}$ where the turbulence becomes strong.
One may then estimate the evolution of that dynamic alignment angle
$\theta(k_{\perp})$ from that point for the cases in which it impacts the
evolution. Finally, although the development of current sheets in plasma
turbulence is well established (Matthaeus & Montgomery, 1980; Meneguzzi et
al., 1981; Borovsky, 2008; Uritsky et al., 2010; Osman et al., 2011; Zhdankin
et al., 2012, 2013), it is not clear that the scaling theory in _B06_ indeed
properly describes this development; an alternative explanation of the
development of current sheets due to from Alfvén wave collisions exists (Howes
& Nielson, 2013; Howes, 2015, 2016; Verniero et al., 2018), although a
specific scaling for this alternative mechanism has not yet been proposed.

###  3.2 Efficiently Modeling the Small-Scale End of the Turbulent Inertial
Range

The multiscale physics of astrophysical plasma turbulence presents a
particular challenge for numerical modeling that covers the entire dynamic
range of the turbulent cascade while simultaneously capturing the kinetic
mechanisms the govern the damping of the turbulence under weakly collisional
plasma conditions, in particular since physics of plasma turbulence is a
inherently three-dimensional (Howes, 2015). For the development of turbulent
heating models, it is critical to investigate the physical mechanisms that
remove energy from the turbulence and energize the plasma particles.
Fortunately, since the inertial range of turbulence is defined as the range of
scales over which the turbulent energy cascades to smaller scales but the
dissipation is negligible, as depicted in Fig. 1(a), it is not strictly
necessary to model the the entire inertial range in order to capture the
physics of the turbulent dissipation.

Figure 3: Modeling the small-scale end of the turbulent inertial range using
MHD scaling relations, choosing _B06_ scaling with $\alpha=1$,
$m_{i}/m_{e}=100$, $k_{0}\rho_{i}=10^{-4}$, $k_{\perp D}\rho_{i}=10^{-2}$.

One less computationally costly approach to capture numerically the mechanisms
of turbulent damping is to choose to model directly only the small-scale end
of the inertial range, preserving as much of the resolved dynamic range of the
simulation as possible to capture the dissipation mechanisms. To do so, one
takes advantage of the scaling relations for turbulence in the MHD regime
$k_{\perp}\rho_{i}\ll 1$, summarized in §3, to “drive” the turbulence at the
domain scale of the simulation with turbulent fluctuations that have the
properties dictated by the scaling theories.

As an example of how to set up simulations in this way, consider turbulence
physically driven strongly and isotropically at an isotropic driving
wavenumber of $k_{0}\rho_{i}=10^{-4}$, as depicted in Fig. 3. For a hydrogenic
plasma of protons and electrons with a realistic mass ratio
$m_{i}/m_{e}=1836$, the turbulent cascade extends at least down to the
perpendicular scales of the electron Larmor radius $k_{\perp}\rho_{e}\sim 1$
(Sahraoui et al., 2010; Howes et al., 2011; TenBarge et al., 2013; Kiyani et
al., 2015). Thus, a dynamic range of
$L_{o}/\rho_{e}=k_{\perp\mbox{max}}/k_{0}=(k_{\perp\mbox{max}}\rho_{e})/(k_{0}\rho_{i})(\rho_{i}/\rho_{e})\propto(k_{0}\rho_{i})^{-1}(T_{i}/T_{e})^{1/2}(m_{i}/m_{e})^{1/2}\sim
4.3\times 10^{5}$ in each of the three spatial dimensions is necessary when
$T_{i}/T_{e}\sim 1$, well beyond the capabilities of any kinetic simulation
code. If the simulation code can resolve a dynamic range of $10^{3}$ in each
spatial dimension, we can choose a perpendicular domain scale in the middle of
the physical inertial range at $k_{\perp D}\rho_{i}=10^{-2}$ and use a reduced
mass ratio $m_{i}/m_{e}=100$ so that $k_{\perp}\rho_{e}=1$ corresponds to
$k_{\perp}\rho_{i}=10$ for $T_{i}/T_{e}=1$. Thus, the resolved perpendicular
length scales of the simulation span $10^{-2}\leq k_{\perp}\rho_{i}\leq 10$,
corresponding to $10^{-3}\leq k_{\perp}\rho_{e}\leq 1$ in terms of the
electron Larmor radius $\rho_{e}$, as shown in Fig. 3.

To drive the turbulence in the middle of the turbulent inertial range at
$k_{\perp D}\rho_{i}=10^{-2}$, we use the scaling relations (1)–(4) to
determine the appropriate turbulent scale-dependent anisotropy and amplitude
of the fluctuations at the numerical domain scale for a _physical_ inertial
range driven strongly ($\chi_{0}=1$) and isotropically ($k_{\perp
0}=k_{\parallel 0}=k_{i0}\equiv k_{0}$ and $\theta_{0}=1$ rad) at
$k_{0}\rho_{i}=10^{-4}$. For this example, we will adopt the _B06_ scaling,
choosing $\alpha=1$. The wavevector anisotropy and amplitude of the
fluctuations at the domain scale can be expressed using scaling relations
(1)–(4) to yield

$k_{\parallel D}\rho_{i}=k_{0}\rho_{i}\left(\frac{k_{\perp
D}\rho_{i}}{k_{0}\rho_{i}}\right)^{1/2}$ (12)
$k_{iD}\rho_{i}=k_{0}\rho_{i}\left(\frac{k_{\perp
D}\rho_{i}}{k_{0}\rho_{i}}\right)^{3/4}$ (13)
$\theta_{D}=\theta_{0}\left(\frac{k_{\perp
D}\rho_{i}}{k_{0}\rho_{i}}\right)^{-1/4}$ (14) $\frac{\delta B_{\perp
D}}{B_{0}}=\left(\frac{k_{\perp D}\rho_{i}}{k_{0}\rho_{i}}\right)^{-1/4}$ (15)

The properties of the domain scale fluctuations at $k_{\perp
D}\rho_{i}=10^{-2}$ determined by these scaling relations are illustrated in
Fig. 3 for a turbulent cascade parameterized by dimensionless parameters
$k_{0}\rho_{i}=10^{-4}$ and $\theta_{0}=1$ rad. The normalized value of energy
spectrum at $k_{\perp D}\rho_{i}=10^{-2}$ is
$(k_{0}/B_{0})^{2}E_{B_{\perp}}(k_{\perp})=10^{-3}$, consistent with the
solution $\delta B_{\perp D}/B_{0}=0.32$ of (15) at the domain scale, given by
the red circle in Fig. 3(a). The three-dimensional anisotropy of the
wavevector of the turbulent fluctuations at the perpendicular domain scale
$k_{\perp D}\rho_{i}=10^{-2}$ is characterized by: (i) $k_{\parallel
D}\rho_{i}=10^{-3}$ from (12), given by the blue circle in Fig. 3(b); and (ii)
$k_{iD}\rho_{i}=3.2\times 10^{-3}$ from (13), given by the green circle in
Fig. 3(b). The dynamic alignment of the velocity and magnetic field
fluctuations at the domain scale is given by $\theta_{D}=0.32$ rad from (14).
Note that for these parameters, the scalings lead to a constant nonlinearity
parameter as a function of perpendicular wavenumber
$\chi(k_{\perp})=\chi_{0}=1$ at all scales $k_{\perp}\rho_{i}<1$, as dictated
by the conjecture of critical balance in strong turbulence (Goldreich &
Sridhar, 1995; Mallet et al., 2015). These scaling calculations determine the
nature of the turbulent fluctuations at the numerical domain scale, in the
middle of the physical inertial range being modeled.

Since the large-scale of the resolved dynamic range of the simulation at
$k_{\perp D}\rho_{i}$ is constantly being fed turbulent cascade energy
mediated by nonlinear interactions among turbulent fluctuations at scales
slightly larger than the domain scale with $k_{\perp}<k_{\perp D}$, it is
imperative that these turbulence simulations are driven, rather than simply
simulating the decay of the energy turbulent fluctuations initialized at the
domain scale. The nonlinear interactions among perpendicularly polarized,
counterpropagating Alfvén waves drive the cascade of energy in Alfvénic plasma
turbulence (Kraichnan, 1965; Sridhar & Goldreich, 1994; Goldreich & Sridhar,
1995; Ng & Bhattacharjee, 1996; Galtier et al., 2000; Howes et al., 2012b;
Howes & Nielson, 2013; Nielson et al., 2013; Howes et al., 2013; Drake et al.,
2013; Howes, 2016; Verniero et al., 2018; Verniero & Howes, 2018). Exploiting
this physical insight, driving counterpropagating Alfvén waves, polarized in
both dimensions perpendicular to the equilibrium magnetic field, has been
shown to generate successfully a steady-state plasma turbulent cascade (Howes
et al., 2008b, 2011; TenBarge et al., 2013; Told et al., 2015; Verniero et
al., 2018; Verniero & Howes, 2018). A finite-time correlated driving, such as
an antenna driven by a Langevin equation (TenBarge et al., 2014), mimics the
unsteady nature of the turbulent cascade from larger scales expected of
turbulent driving at a chosen scale within the turbulent inertial range.

To set up the simulation to describe optimally the domain scale fluctuations,
it has become common practice to utilize a numerical domain that is elongated
along the direction of the equilibrium magnetic field, consistent with the
anisotropic scaling of turbulent fluctuations at small scales,
$k_{\parallel}/k_{\perp}\ll 1$. For the example shown in Fig. 3, individual
Alfvénic fluctuations are predicted to have a perpendicular plane anisotropy
given by $k_{iD}/k_{\perp D}=0.32$. But since these Alfvénic fluctuations can
be polarized in either of the directions in the perpendicular plane, it is
necessary for both perpendicular directions to be large enough to resolve the
intermediate scale (the larger of the two perpendicular scales since
$k_{iD}/k_{\perp D}<1$). Thus, the simulation domain should be defined by
$L_{\parallel}\times L_{i}^{2}$, where
$L_{\parallel}=2\pi\rho_{i}(k_{\parallel D}\rho_{i})^{-1}$ and
$L_{i}=2\pi\rho_{i}(k_{iD}\rho_{i})^{-1}$. Within this domain,
counterpropagating Alfvénic fluctuations should be driven, polarized in both
perpendicular directions, with wavevectors characterized by $(k_{\perp
D},k_{iD},k_{\parallel D})$, with dynamic alignment characterized by
$\theta_{D}$, and with amplitudes $\delta B_{\perp D}$. This approach enables
the efficient simulation of the small-scale end of the inertial range
consistent with the scaling properties of MHD turbulence, enabling kinetic
simulations to focus on the physical mechanisms damping the turbulent
fluctuations.

Plasma Parameters
---
Parameter | Symbol and Definition
Ion Plasma Beta | $\beta_{i}=8\pi n_{i}T_{i}/B^{2}$
Ion-to-Electron Temperature Ratio | $\tau=T_{i}/T_{e}$
Ion-to-Electron Mass Ratio | $\mu=m_{i}/m_{e}$
Turbulence Parameters
Parameter | Symbol and Definition
Isotropic Driving Wavenumber | $k_{0}\rho_{i}$
Table 3: _Isotropic Temperature Case_ : Reduced set of fundamental
dimensionless parameters for isotropic Maxwellian equilibrium velocity
distributions with $A_{i}=A_{e}=1$, collisionless conditions $k_{\parallel
0}\lambda_{mfp,e}\gg 1$, strong turbulence assumed at the ion kinetic scales,
balanced turbulence $Z_{0}^{+}/Z_{0}^{-}\sim 1$, and negligible energy in
compressible fluctuations $E_{comp}/E_{inc}\ll 1$. This results in a reduced
set of plasma and turbulence parameters $(\beta_{i},\tau;k_{0}\rho_{i})$.

###  3.3 A Reduced Parameter Space for Kinetic Plasma Turbulence

Although the physics of plasma turbulence, its dissipation, and the resulting
heating of the plasma species formally depends on all ten of the general
plasma and turbulence parameters
$(\beta_{\parallel,i},\tau_{\parallel},A_{i},A_{e},k_{\parallel
0}\lambda_{mfp,e};k_{\perp 0}\rho_{i},k_{\parallel 0}/k_{\perp
0},\chi_{0},Z_{0}^{+}/Z_{0}^{-},E_{comp}/E_{inc})$ listed in Table 2, focusing
initially on the development of completely general turbulent heating models on
a ten-dimensional parameter space is unlikely to be successful. Instead, first
developing predictive turbulent heating models over a reduced parameter space
is likely a better strategy.

Note that kinetic temperature anisotropy instabilities generally bound the
temperatures in the vicinty of isotropy, although if these instabilities are
triggered, at $\beta_{i}\gg 1$, they can impact both the linear wave physics
at play in the turbulent dynamics of critically balanced turbulence (Squire et
al., 2016, 2017b, 2017a) and mediate the nonlocal transfer of energy directly
from the large driving scales to the kinetic length scales (Arzamasskiy et
al., 2023). As an initial exploration of Alfvénic turbulence and its kinetic
damping mechanisms, however, we choose to focus on the simplified limit of
isotropic temperatures, modifying our approach as needed for cases that
violate these conditions.

For many of the calculations in this study, a reduced isotropic temperature
case is adopted, assuming that the equilibrium velocity distribution for each
species remains well described as an isotropic Maxwellian distribution. Note
that kinetic temperature anisotropy instabilities generally bound the species
temperature anisotropies in the vicinity of isotropy $A_{s}\sim 1$ (Kasper et
al., 2002; Hellinger et al., 2006; Štverák et al., 2008; Bale et al., 2009;
Štverák et al., 2009), although if these instabilities are triggered, they can
impact both the turbulent dynamics and energy transfer to particles (Squire et
al., 2016, 2017b, 2017a; Arzamasskiy et al., 2023). In this isotropic
temperature case, the temperature anisotropies for ions and electrons are both
set to unity, $A_{i}=1$ and $A_{e}=1$. In addition, the ion plasma beta and
ion-to-electron temperature ratio are redefined in terms of these isotropic
temperatures, yielding the isotropic parameters $\beta_{i}=8\pi
n_{i}T_{i}/B^{2}$ and $\tau=T_{i}/T_{e}$. In addition, for this reduced
parameter space, we will also assume that the turbulent inertial range is
sufficiently large, or the driving is sufficiently strong, that the dynamics
at the perpendicular ion kinetic length scales $k_{\perp}\rho_{i}\sim 1$ are
in a state of strong plasma turbulence with a nonlinearity parameter $\chi\sim
1$. Under these conditions, we can replace the three turbulent driving
parameters $(k_{\perp 0}\rho_{i},k_{\parallel 0}/k_{\perp 0},\chi_{0})$ by the
single isotropic driving wavenumber $k_{0}\rho_{i}$, as detailed in §3.1.

We will also make two other simplfying assumptions for the isotropic
temperature case: We will assume that the turbulence is balanced, with equal
wave energy fluxes up and down the local magnetic field
$Z_{0}^{+}/Z_{0}^{-}\sim 1$, which is a necessary condition to apply the MHD
turbulence scaling relations given by (1)–(4). Additionally, based on direct
spacecraft measurements of turbulence in the heliosphere which find that
incompressible Alfvénic fluctuations appear to dominate turbulent dynamics in
space plasmas (Belcher & Davis, 1971; Tu & Marsch, 1995; Alexandrova et al.,
2008; Bruno & Carbone, 2013), as well as theoretical considerations that
suggest that the dynamics of Alfvén waves dominate the turbulent cascade
(Lithwick & Goldreich, 2001; Maron & Goldreich, 2001; Cho & Lazarian, 2003;
Schekochihin et al., 2009; Howes, 2015), we will also assume a negligible
amount of energy in compressible turbulent fluctuations, $E_{comp}/E_{inc}\ll
1$.

Therefore, we propose the _isotropic temperature case_ , a reduced parameter
space for weakly collisional plasma turbulence in the limit of strong,
balanced, incompressible turbulence, with a governing set of three
dimensionless parameters, $(\beta_{i},\tau;k_{0}\rho_{i})$, as shown in Table
3. Although the isotropic temperature case eliminates a number of important
physical effects such as temperature anisotropy instabilities, imbalanced
turbulent wave energy fluxes, and compressible fluctuations, it nonetheless
represents an excellent starting point for the development of refined
turbulent heating models. Such models can produce theoretical predictions that
can be compared to kinetic numerical simulations and spacecraft observations
of weakly collisional plasma turbulence, exemplified by the phase diagram for
turbulent damping mechanisms shown in Fig. 11.

In terms of the fundamental parameters of the isotropic temperature case, the
turbulent scaling relations (1)–(4) simplify to the following easy-to-use
forms

$k_{\parallel}\rho_{i}=(k_{0}\rho_{i})^{(1+\alpha)/(3+\alpha)}(k_{\perp}\rho_{i})^{2/(3+\alpha)}$
(16)
$k_{i}\rho_{i}=(k_{0}\rho_{i})^{\alpha/(3+\alpha)}(k_{\perp}\rho_{i})^{3/(3+\alpha)}$
(17)
$\theta=\theta_{0}\left(\frac{k_{0}\rho_{i}}{k_{\perp}\rho_{i}}\right)^{\alpha/(3+\alpha)}$
(18) $\frac{\delta
B_{\perp}}{B_{0}}=\left(\frac{k_{0}\rho_{i}}{k_{\perp}\rho_{i}}\right)^{1/(3+\alpha)}$
(19)

##  4 Review of Previous Turbulent Heating Models

Several previous studies have constructed turbulent heating models that
attempt to quantify how the energy removed from the turbulent cascade is
partitioned between protons and electrons and between parallel and
perpendicular degrees of freedom. Here we briefly summarize these existing
turbulent heating models, including the dissipation mechanisms considered, the
assumptions of how that energy is removed from the turbulence, and that stated
limitations of each model.

One of the earliest efforts to estimate the proton and electron energization
rates by turbulence due to collisionless damping mechanisms was carried out by
Quataert & Gruzinov (1999) (hereafter _QG99_), following earlier results by
each of these authors to determine the collisionless damping rates under
astrophysically relevant conditions (Quataert, 1998; Gruzinov, 1998). This
study adopted the _GS95_ scaling for dominantly Alfvénic turbulence at
perpendicular scales restricted to the MHD limit, $k_{\perp}\rho_{i}\lesssim
1$. The turbulent fluctuations upon reaching the end of the MHD limit at the
perpendicular scale of the ion Larmor radius $k_{\perp}\rho_{i}\sim 1$ are
predicted to be significantly anisotropic with $k_{\perp}\gg k_{\parallel}$;
since the MHD Alfvén wave frequency $\omega=k_{\parallel}v_{A}$ is
proportional to $k_{\parallel}$, this yields turbulent fluctuations at this
perpendicular ion scale with low frequencies relative to the ion cyclotron
frequency, $\omega\ll\Omega_{i}$. Thus, all collisionless damping is provided
by the $n=0$ Landau resonance, specifically Landau damping and transit-time
damping. In addition, this model assumed that all turbulent energy at a given
perpendicular wavenumber $k_{\perp}\rho_{i}$ was concentrated at the parallel
wavenumber given by the conjecture of critical balance (Goldreich & Sridhar,
1995; Quataert & Gruzinov, 1999). They constructed a simple cascade model that
balanced the collisionless damping by ions and electrons given by the Vlasov-
Maxwell linear dispersion relation (Stix, 1992) over the MHD range of scales,
$k_{\perp}\rho_{i}\lesssim 1$. Beyond these scales at
$k_{\perp}\rho_{i}\gtrsim 1$ in the kinetic regime, uncertainty about the
nature of the turbulence at these sub-ion length scales led them to assume
simply that any energy reaching $k_{\perp}\rho_{i}\gtrsim 1$ ultimately heated
the electrons. They calculated the fraction of the turbulent energy that heats
the electrons, $Q_{e}/Q$, as function of total plasma beta
$\beta=\beta_{i}+\beta_{e}=\beta_{i}(1+T_{i}/T_{e})$, yielding a numerical
result for $Q_{e}/Q(\beta)$. Their resulting energy partition depended
strongly on a constant assumed in the model that characterized the balance of
the nonlinear turbulent cascade rate to the linear collisionless damping rate;
uncertainty in the value of this constant (by a factor from $1/4$ to $4$) led
the resulting prediction of $Q_{e}/Q(\beta)$ to vary by more than one order-
of-magnitude over the range of $3\lesssim\beta\lesssim 100$.

The Howes 2010 (hereafter _H10_) turbulent heating model (Howes, 2010) was
based on a turbulent cascade model (Howes et al., 2008a) that extended the
_GS95_ turbulence scaling from MHD scales at $k_{\perp}\rho_{i}\lesssim 1$
into the kinetic regime at $k_{\perp}\rho_{i}\gtrsim 1$, generalizing the
conjecture of critical balance to account for the dispersive nature of kinetic
Alfvén waves at these sub-ion length scales (Biskmap et al., 1999; Cho &
Lazarian, 2004; Krishan & Mahajan, 2004; Shaikh & Zank, 2005; Schekochihin et
al., 2009). The _H10_ model was based on three primary assumptions: (i) the
Kolmogorov hypothesis that the energy cascade is determined by local
interactions (Kolmogorov, 1941); (ii) the conjecture that turbulence maintains
a state of critical balance at all scales (Goldreich & Sridhar, 1995); and
(iii) the speculation that linear collisionless damping rates are applicable
even in the presence of the strong nonlinear interactions that mediate the
turbulent cascade to small scales. The casade model used by _H10_ balanced at
each perpendicular wavenumber $k_{\perp}\rho_{i}$ the nonlinear energy cascade
rate given by the extended turbulent scaling theories with the linear
collisionless damping rates at that value of $k_{\perp}\rho_{i}$. Similar to
_QG99_ , it assumed all turbulent energy at a given value of
$k_{\perp}\rho_{i}$ resided at the parallel wavenumber governed by the
extended critical balance. The _H10_ model assumed turbulence under
nonrelativistic conditions $v_{te}/c\ll 1$ with isotropic ion and electron
temperatures, strictly Alfvénic turbulence with $E_{comp}/E_{inc}=0$, balance
turbulence with $Z_{0}^{+}/Z_{0}^{-}=1$, and strongly and isotropically driven
turbulence with $\chi_{0}=1$ and $k_{\parallel 0}/k_{\perp 0}=1$ (the driving
was therefore characterized by the isotropic driving wavenumber,
$k_{0}\rho_{i}$). Furthermore, the same as _QG99_ , the model assumed that the
inertial range of the turbulence is sufficiently large that the anisotropic
fluctuations predicted by the turbulent scaling relations remain at a low
enough frequency relative to the ion cyclotron frequency $\omega\ll\Omega_{i}$
that the cyclotron resonance plays a negligible role; the _H10_ model
calculated an estimate of the maximum value of $k_{0}\rho_{i}$ as a function
of $(\beta_{i},T_{i}/T_{e})$ for this assumption to be satisfied. Thus, the
_H10_ turbulent heating model predicted the damping of the turbulent cascade
due to Landau damping and transit-time damping over the full range of the
turbulent cascade, including MHD scales at $k_{\perp}\rho_{i}\lesssim 1$ and
kinetic scales at $k_{\perp}\rho_{i}\gtrsim 1$. The _H10_ model predicted the
partitioning of energy between ions and electrons as a function of the ion
plasma beta $\beta_{i}$ and the ion-to-electron temperature ratio
$T_{i}/T_{e}$, providing a simple analytical fit to the numerically calculated
predictions, $Q_{i}/Q_{e}(\beta_{i},T_{i}/T_{e})$. In short, the model found
that $Q_{i}/Q_{e}$ is a montonically increasing function of $\beta_{i}$ with
weak dependence on $T_{i}/T_{e}$, where the transition from $Q_{i}/Q_{e}<1$ to
$Q_{i}/Q_{e}>1$ happens at approximately $\beta_{i}\sim 1$.

The Chandran 2011 (hereafter _C11_) turbulent heating model (Chandran et al.,
2011) looked to incorporate the recently quantified damping rate due to
nonlinear ion stochastic heating (Chandran et al., 2010; Chandran, 2010) into
a model that also accounted for ion and electron Landau damping and transit-
time damping. Similar to _QG99_ and _H10_ , it assumed a sufficiently large
inertial range that no significant cyclotron damping would occur since the
anisotropic cascade would lead to turbulent fluctuations with low frequencies,
$\omega\ll\Omega_{i}$. Focusing also on Alfvénic turbulence that transitions
from an MHD Alfvén wave cascade at $k_{\perp}\rho_{i}\lesssim 1$ and to a
kinetic Alfvén wave cascade at $k_{\perp}\rho_{i}\gtrsim 1$, the _C11_ model
assumed that dissipation occured dominantly in two distinct wavenumber ranges:
(i) at $k_{\perp}\rho_{i}\sim 1$, where collisionless damping via the Landau
($n=0$) resonance can lead to parallel energization of the ions and electrons
and also where nonlinear stochastic heating can lead to perpendicular
energization of the ions; and (ii) at $k_{\perp}\rho_{i}\gg 1$, where all
remaining energy in the turbulent cascade is ultimately deposited with the
electrons. The model used fitting forms for the Landau-resonant collisionless
damping rates for the ions and electrons at $k_{\perp}\rho_{i}=1$ over
$10^{-3}<\beta_{i}<10$ and $1<T_{i}/T_{e}<5$. It also included an empirically
derived formula for the ion stochastic heating rate at $k_{\perp}\rho_{i}=1$
(Chandran et al., 2010). This empirical prescription for the ion stochastic
heating depended on two dimensionless constants: one governed the scaling of
the stochastic heating rate relative to the turbulent frequencies; and the
other, which appeared as a factor in the argument of the exponential function,
effectively established a threshold amplitude below which ion stochastic
heating was negligible. Thus, the _C11_ model calculated how the turbulent
cascade rate at $k_{\perp}\rho_{i}=1$ is partitioned among ion Landau and
transit-time damping, electron Landau and transit-time damping, ion stochastic
heating, and nonlinear transfer to smaller scales with $k_{\perp}\rho_{i}>1$.
The model provided a partitioning of the total turbulent heating rate into
parallel ion heating by ion Landau and transit-time damping
$Q_{\parallel,i}/Q$, perpendicular ion stochastic heating $Q_{\perp,i}/Q$, and
electron heating $Q_{e}/Q$, each as a function of ion plasma beta $\beta_{i}$,
the ion-to-electron temperature ratio $T_{i}/T_{e}$, and the amplitude of the
turbulent fluctuations at the ion Larmor radius scale determined using the
_B06_ scaling for the amplitude in terms of the driving scale, as given here
by (19), indicating that the amplitude dependence can be characterized by
$k_{0}\rho_{i}$.

Motivated by the desire to understand the partitioning of dissipated turbulent
energy between ions and electrons in hot collisionless accretion flows, the
Rowan 2017 (hereafter _R17_) heating model (Rowan et al., 2017) took a
fundamentally different approach from that adopted in the _QG99_ , _H10_ , and
_C11_ turbulent heating models by focusing instead in the heating arising from
transrelativistic magnetic reconnection. In this case, the protons are
nonrelativistic, but the electrons can be ultrarelativistic. Utilizing a large
suite of 2D particle-in-cell (PIC) simulations of anti-parallel reconnection,
the _R17_ model studied the irreversible particle energization over a range of
ion plasma beta $10^{-4}\lesssim\beta_{i}\leq 2$ and magnetization
$0.1\leq\sigma_{w}\leq 10$, where $\sigma_{w}$ is defined as the ratio of the
magnetic energy density to the enthalpy density. Starting with a Harris
current sheet (Harris, 1962) of the antiparallel magnetic field (given by a
$\tanh$ profile) with width $a\sim 40d_{e}$, they carefully separated the
irreversible heating associated with an increase in the entropy from the
reversible heating associated with adiabatic compression. From their synthesis
of the results of a large suite of simulations varying different parameters,
they fit a formula to their numerical results that provided to electron-to-
total heating ratio $Q_{e}/Q$ as a function of the ion plasma beta $\beta_{i}$
and magnetization $\sigma_{w}$.

The numerical study by Kawazura et al. (2019) (hereafter _K19_) used a large
suite of hybrid gyrokinetic simulations to determine directly the ion-to-
electron heating ratio $Q_{i}/Q_{e}$ over a broad range of ion plasma beta
$0.1\leq\beta_{i}\leq 100$ and ion-to-electron temperature ratio $0.05\leq
T_{i}/T_{e}\leq 100$. Using a hybrid code (Kawazura & Barnes, 2018) that
evolves the ions gyrokinetically and the electrons as an isothermal fluid,
this study adopted many of the same assumptions as the _H10_ model: (i) the
anisotropic ($k_{\perp}\gg k_{\parallel}$) Alfvénic turbulence was driven in a
statistically balanced manner with $Z_{0}^{+}/Z_{0}^{-}=1$; (ii) the
equilibrium ion and electron temperatures were assumed to be isotropic; (iii)
the gyrokinetic approximation eliminated the physics of the ion cyclotron
resonances, leaving the Landau resonances and collisionless magnetic
reconnection in the large-guide-field limit as potential damping mechanisms;
and (iv) the strategy explained in §3.2 was exploited to model numerically
only the small-scale end of the inertial range, where the isotropic driving
wavenumber $k_{0}\rho_{i}$, characterizing the physical driving scale, was
assumed to be sufficiently small that the amplitudes of the turbulent
fluctuations at $k_{\perp}\rho_{i}\sim 1$ were small enough that ion
stochastic heating was inhibited (Chandran et al., 2010, 2011). The _K19_
model provides a numerically validated prescription for $Q_{i}/Q_{e}$ as a
function of $\beta_{i}$ and $T_{i}/T_{e}$, presenting a simple analytical
formula that fit the numerical results well. The numerical results largely
validated the qualitative predictions of the _H10_ turbulent heating model
that $Q_{i}/Q_{e}$ is a monotonic function of $\beta_{i}$ with little
dependence on $T_{i}/T_{e}$, with a couple of minor quantitative differences:
the _K19_ model found a ceiling of $Q_{i}/Q_{e}\simeq 30$ at $\beta_{i}\gg 1$,
cutting off at lower values than the _H10_ model; and the numerical results
did not find the drop in the heating rate $Q_{i}/Q_{e}\ll 0.1$ predicted by
the _H10_ model at $\beta_{i}\ll 1$.

Subsequently, using the same hybrid gyrokinetic code, Kawazura et al. (2020)
explored the effect of compressible driving $E_{comp}/E_{inc}>0$ on the ion-
to-electron heating ratio $Q_{i}/Q_{e}$ (hereafter _K20_). The turbulence was
driven at MHD scales $k_{\perp}\rho_{i}\ll 1$ with an adjustable mixture of
incompressible Alfvénic fluctuations and compressible slow magnetosonic
fluctuations (fast magnetosonic wave modes are ordered out of gyrokinetic
simulations by the perpendicular pressure balance that is maintained under the
gyrokinetic approximation (Howes et al., 2006)). At $\beta_{i}\ll 1$, they
found that $Q_{i}/Q_{e}=E_{comp}/E_{inc}$, confirming a previous theoretical
prediction that compressive energy strictly heats ions at low ion plasma beta
$\beta_{i}\ll 1$ (Schekochihin et al., 2019). Surprisingly, the _K20_ model
also found that $Q_{i}/Q_{e}=E_{comp}/E_{inc}$ when the compressible driving
dominates, $E_{comp}/E_{inc}\gg 1$. A simple fitting formula was constructed
that simply added the term $E_{comp}/E_{inc}$ to the previous prescription for
$Q_{i}/Q_{e}$ due to Alfvénic turbulence in the _K19_ model. Thus, the _K20_
model provided a simple analytical formula for
$Q_{i}/Q_{e}(\beta_{i},T_{i}/T_{e},E_{comp}/E_{inc})$.

The application of the _K20_ model to understand the ion-to-electron heating
ratio in accretion disks where the magnetorotational instability (MRI) (Balbus
& Hawley, 1991; Hawley & Balbus, 1991; Balbus & Hawley, 1998) drives the
turbulence requires a determination of the ratio of the compressible-to-
incompressible driving $E_{comp}/E_{inc}$. The application of a rotating
reduced MHD (RRMHD) model to simulate the collisional MRI turbulence threaded
by a nearly aximuthal magnetic field in a 3D pseudo-spectral, shearing reduced
MHD code (Kawazura, 2022) found a ratio of compressible-to-incompressible
turbulent energy of $2\lesssim E_{comp}/E_{inc}\lesssim 2.5$ (Kawazura et al.,
2022), providing a solid foundation for the application of the _K20_ model to
astrophysical accretion disks.

##  5 Dependence of Dissipation Mechanisms on Fundamental Parameters

The development of turbulent heating models for weakly collisional space and
astrophysical plasmas requires two fundamental steps: (i) identifying the
mechanisms that play a role in the damping of the turbulence; and (ii)
determining upon which of the fundamental dimensionless parameters of
turbulence each mechanism depends and approximating a functional form for
those dependencies.

It is worthwhile noting an important subtlety regarding irreversible plasma
heating in the thermodynamic picture of turbulent dissipation in weakly
collisional plasmas (Howes, 2017). In the strongly collisional limit of fluid
plasma turbulence, the dissipation of turbulence is fundamentally collisional
on a microscopic level. For example, in an MHD plasma, the physical mechanisms
of viscosity and resistivity arise in the limit of a non-negligible
collisional mean free path of the plasma particles. These two dissipative
mechanisms may be derived rigorously beginning with kinetic theory through a
hierarchy of moment equations derived in the limit of small mean free path by
the Chapman-Enskog procedure (Chapman & Cowling, 1970) for neutral gases, or
an analogous procedure for plasma systems (Spitzer, 1962; Grad, 1963;
Braginskii, 1965). Because viscosity and resistivity are ultimately
collisional in nature, the resulting energy transfer from the turbulent
fluctuations to the particles is irreversible, leading to an increase in the
entropy of the system.

Under the weakly collisional conditions typical of many space and
astrophysical plasmas, the collisionless removal of energy from the turbulent
fluctuations and the ultimate conversion of that removed energy into plasma
heat is a two-step process (Howes et al., 2006; Howes, 2008; Schekochihin et
al., 2009; Navarro et al., 2016; Howes, 2017; Howes et al., 2018). First,
collisionless damping mechanisms remove energy from the turbulent
fluctuations, transferring that energy into nonthermal energy contained in
non-Maxwellian fluctuations in the particle velocity distributions, an energy
transfer that is reversible in principle. That free energy in the fluctuations
of the velocity distribution about the equilibrium may be transferred to
smaller scales in velocity space— _e.g._ , via linear phase mixing by the
ballistic term of the Boltzmann equation (Landau, 1946; Hammett et al., 1992;
Snyder et al., 1997; Schekochihin et al., 2009) or via nonlinear phase mixing,
a process recently denoted the entropy cascade (Dorland & Hammett, 1993;
Schekochihin et al., 2009; Tatsuno et al., 2009; Plunk et al., 2010; Plunk &
Tatsuno, 2011; Kawamori, 2013)—until it reaches sufficiently small velocity-
space scales that an arbitrarily weak collisionality is sufficient to smooth
out those fluctuations, providing an irreversible conversion of that energy to
plasma heat and thereby increasing the entropy of the system (Howes, 2008;
Schekochihin et al., 2009).

The subtleties of energy conversion and plasma heating in weakly collisional
plasmas that are not in a state of local thermodynamic equilibrium remains a
topic of vigorous investigation, with many recent developments that will not
be reviewed in detail here (Schekochihin et al., 2016; Parker et al., 2016;
Yang et al., 2016, 2017a, 2017b; Howes et al., 2018; Kawazura et al., 2019;
Meyrand et al., 2019, 2021; Barbhuiya & Cassak, 2022; Cassak & Barbhuiya,
2022; Cassak et al., 2022; Squire et al., 2022; Cassak et al., 2023). In this
paper, we use the term _turbulent dissipation_ to mean the entire process of
removing energy from the turbulent fluctuations and subsequently thermalizing
that energy via collisions to realize irreversible heating of the plasma and
an increase in the entropy of the system. When we refer to only the first step
of the process, the collisionless removal of energy from the turbulent
fluctuations, we denote those specific processes as _turbulent damping
mechanisms_.

A primary aim of this study is to provide a unified framework for
characterizing the dependence of the turbulence and its dissipation on a
common set of dimensionless parameters. The first step pursued here is to
identify the fundamental dimensionless parameters upon which each proposed
turbulent damping mechanism primarily depends, where the set of dimensionless
parameters for the general case is presented in Table 2 and for the isotropic
temperature case is in Table 3. The key parameter dependencies found here are
summarized for each proposed damping mechanism in Table 4. Although, for any
single particular mechanism, different choices of dimensionless parameters may
be more natural (_e.g._ , taking $\beta_{e}$ instead of $\beta_{i}$ in the
case of magnetic reconnection), the particular sets of dimensionless
parameters in Tables 2 and 3 are proposed here as a useful framework for
describing _all_ of the proposed turbulent damping mechanisms, thereby
simplifying the task of detemining how different mechanisms compete with each
other throughout the full parameter space.

Category | Abbrev | Mechanism | Key Parameters
---|---|---|---
Resonant | iLD | Ion Landau Damping | $\beta_{i}$
| eLD | Electron Landau Damping |  $\beta_{i},T_{i}/T_{e}$, ($\mu$)
| iTTD | Ion Transit-Time Damping | $\beta_{i}$
| eTTD | Electron Transit-Time Damping |  $\beta_{i},T_{i}/T_{e}$, ($\mu$)
| iCD | Ion Cyclotron Damping | $k_{0}\rho_{i},\beta_{i}$
Non- | iSH | Ion Stochastic Heating | $k_{0}\rho_{i},\beta_{i}$
Resonant | iMP | Ion Magnetic Pumping | $E_{comp}/E_{inc},\beta_{i},T_{i}/T_{e},k_{\parallel 0}\lambda_{mfp,e},$
|  |  |  $k_{\parallel 0}/k_{\perp 0},\chi_{0}$
| eMP | Electron Magnetic Pumping | $\beta_{i},T_{i}/T_{e},k_{\parallel 0}\lambda_{mfp,e},k_{\parallel 0}/k_{\perp 0},\chi_{0},$
|  |  |  $E_{comp}/E_{inc}$, ($\mu$)
| iVH | Ion Kinetic Viscous Heating | $\beta_{i},k_{\perp 0}/k_{\parallel 0},\chi_{0}$
Coherent | iRXN | Ion Magnetic Reconnection |  $k_{0}\rho_{i},\beta_{i},T_{i}/T_{e}$, ($\mu$)
Structures | eRXN | Electron Magnetic Reconnection |  $k_{0}\rho_{i},\beta_{i},T_{i}/T_{e}$, ($\mu$)
Table 4: Proposed turbulent damping mechanisms and their key dimensionless
parameter dependencies.

###  5.1 Proposed Damping Mechanisms for Kinetic Turbulence

In weakly collisional space and astrophysical plasmas, the proposed kinetic
damping mechanisms that remove energy from the turbulent fluctuations and
transfer it to the particle species fall into three broad categories: (i)
resonant wave-particle interactions, (ii) non-resonant wave-particle
interactions, and (iii) damping occurring in coherent structures.

Resonant mechanisms for the damping of the turbulent fluctuations in a
magnetized plasma transfer energy to particles in regions of the velocity
distribution that satisfy the resonant condition

$\omega-k_{\parallel}v_{\parallel}-n\Omega_{s}=0,$ (20)

where $\Omega_{s}$ is the cyclotron frequency for species $s$ and $|n|>1$
indicates the harmonics of the cyclotron frequency. The $n=0$ resonance is
often denoted the _Landau resonance_ and includes two damping mechanisms: (i)
Landau damping is mediated by the parallel component of the electric field
$E_{\parallel}$ doing work on the charged particles (Landau, 1946; Leamon et
al., 1998b, a; Quataert, 1998; Leamon et al., 1999; Quataert & Gruzinov, 1999;
Leamon et al., 2000; Howes et al., 2008a; Schekochihin et al., 2009; TenBarge
& Howes, 2013; Howes, 2015; Li et al., 2016; Chen et al., 2019; Afshari et
al., 2021); and (ii) transit-time damping, also known as Barnes damping, is
mediated by the magnetic mirror force associated with magnetic field magnitude
fluctuations doing work on the magnetic moment of the charged particle
gyromotion (Barnes, 1966; Quataert, 1998; Quataert & Gruzinov, 1999; Howes et
al., 2008a; Huang et al., 2024). Both of these Landau-resonant collisionless
damping mechanisms energize particles in the degree of freedom parallel to the
local magnetic field. The $n\neq 0$ resonances yield cyclotron damping of the
turbulent flucutations by $E_{\perp}$, energizing particles in the two
perpendicular degrees of freedom (Coleman, 1968; Denskat et al., 1983;
Isenberg & Hollweg, 1983; Goldstein et al., 1994; Leamon et al., 1998a; Gary,
1999; Isenberg et al., 2001; Hollweg & Isenberg, 2002; Isenberg & Vasquez,
2019; Afshari et al., 2023).

Non-resonant damping mechanisms can remove energy from turbulent fluctuations
through wave-particle interactions that do not require satisfying a particular
resonance condition. Stochastic heating can arise when sufficiently large
amplitude turbulent electromagnetic fluctuations can effectively scatter the
otherwise organized gyromotion of particles, leading to a heating of the
particles in the perpendicular degrees of freedom (Johnson & Cheng, 2001; Chen
et al., 2001; White et al., 2002; Voitenko & Goossens, 2004; Bourouaine et
al., 2008; Chandran et al., 2010; Chandran, 2010; Chandran et al., 2011;
Bourouaine & Chandran, 2013; Chandran et al., 2013; Klein & Chandran, 2016;
Arzamasskiy et al., 2019; Martinović et al., 2020; Cerri et al., 2021).
Magnetic pumping energizes particles through a combination of two effects: (i)
particles undergo the reversible, double adiabatic evolution of their
perpendicular and parallel velocities in a magnetic field with a time varying
magnitude; and (ii) those particles are scattered by collisions during this
evolution, introducing irreversibility and leading to a small net transfer of
energy to the particles over a full oscillation cycle (Spitzer & Witten, 1951;
Berger et al., 1958; Lichko et al., 2017; Lichko & Egedal, 2020; Montag &
Howes, 2022). One final proposed non-resonant damping mechanism is the
recently identified kinetic “viscous” heating mediated by the effective
collisionality associated with temperature anisotropy instabilities driven by
large-scale compressible fluctuations in high beta plasmas (Arzamasskiy et
al., 2023).

The final category involves the removal of turbulent fluctuation energy in
coherent structures, such as particle energization arising from collisionless
magnetic reconnection in current sheets that are found to arise naturally in
plasma turbulence (Ambrosiano et al., 1988; Dmitruk et al., 2004; Markovskii &
Vasquez, 2011; Matthaeus & Velli, 2011; Osman et al., 2011; Servidio et al.,
2011; Osman et al., 2012b, a; Wan et al., 2012; Karimabadi et al., 2013;
Zhdankin et al., 2013; Dalena et al., 2014; Osman et al., 2014a, b; Zhdankin
et al., 2015b, a; Mallet et al., 2017b, a; Loureiro & Boldyrev, 2017b).

Each of the kinetic damping mechanisms in these three categories (with helpful
abbreviations) are summarized in Table 4. Also included in Table 4 are the
primary dimensionless parameters upon which each damping mechanism depends, as
explored in more detail in the following subsections.

###  5.2 A Competition between Turbulent Cascade Rates and Damping Rates

One of the most challenging aspects of modeling the damping of plasma
turbulence is that the net energy transfer to each species arises from the
competition between the turbulent cascade rate and the kinetic damping rates
as a function of scale, as depicted by the arrows in Fig. 1(a). Nonlinear
interactions among turbulent fluctuations—denoted Alfvén wave collisions
(Howes & Nielson, 2013) in the case of Alfvénic turbulence—mediate the
transfer of energy from large to small scales that is dominantly local in
scale space (Howes et al., 2011; Told et al., 2015), where energy is
transferred locally from a wavenumber $k$ to a wavenumber $2k$, then on to
$4k$, and so on. This is indicated by the black curved arrows in the figure.
Because the physics of the cascade is dominantly local in nature (Howes et
al., 2011; Told et al., 2015), scaling theory implies that the energy transfer
rate at a given wavenumber $k$ is determined by the conditions at that scale.
Note that, when analyzed in Fourier space, the energy transfer to small scales
in the turbulent cascade is the sum of the nonlinear interactions among many
triads of plane-wave modes, with a wide spread of energy transfer rates about
zero; the net energy transfer given by this sum at a given scale is much
smaller in amplitude than the spread of individual three-wave interaction
rates, and the sum typically indicates a forward cascade of energy to small
scales (Coburn et al., 2014, 2015). The turbulent cascade rate $\varepsilon$
used here specifically refers to this net energy rate obtained by integrating
over all contributing nonlinear interactions involving that scale.

Figure 4: Diagram of the balance of the energy transfer rates in a steady-
state turbulent cascade at a local wavenumber $k_{*}$, where the following
terms sum to zero: (i) the turbulent nonlinear energy transfer rate from lower
wavenumbers $\varepsilon(k^{-})^{(k^{-}<k_{*})}$; (ii) the turbulent nonlinear
energy transfer rate to higher wavenumbers
$\varepsilon(k_{*})^{(k_{*}<k^{+})}$; and (iii) the local (in scale) kinetic
damping rates of the turbulent energy, due to the sum of ion energization
rates $Q_{i}(k_{*})$ (blue) and electron energization rates $Q_{e}(k_{*})$
(magenta).

Competing with that turbulent cascade rate at each scale is the damping rate
due to collisionless interactions between the electromagnetic fields and the
plasma particles, as shown by the blue arrows for the transfer of energy to
ions and the magenta arrows for the transfer of energy to electrons in Fig.
1(a). The net damping rate of the entire turbulent cascade, resulting from the
sum of the energization rates for each particle species $s$, can be computed
by integration over all scales of the turbulent cascade. In a steady-state, at
each wavenumber $k_{*}$, the sum of three terms must balance to zero, as
illustrated in Fig. 4: (i) the energy transfer from lower wavenumbers
$k^{-}<k_{*}$ (larger scales), denoted by
$\varepsilon(k^{-})^{(k^{-}<k_{*})}$; (ii) the energy transfer to higher
wavenumbers $k^{+}>k_{*}$ (smaller scales), denoted by
$\varepsilon(k_{*})^{(k_{*}<k^{+})}$; and (iii) the kinetic damping of the
fluctations at wavenumber $k_{*}$, given by $Q(k_{*})$, which is the sum of
the local (in scale) ion damping $Q_{i}(k_{*})$ (blue) and electron damping
$Q_{e}(k_{*})$ (magenta). Complicating the determination of the partitioning
of damped turbulent energy between ion and electrons $Q_{i}/Q_{e}$ is the fact
that the kinetic damping rates due to ions and electrons often overlap as a
function of scale (Howes et al., 2011, 2011; Told et al., 2015).

Figure 5: Weakened cascade model results for a turbulent cascade driven weakly
at $k_{\perp 0}\rho_{i}=k_{\parallel}\rho_{i}=10^{-5}$ with nonlinearity
parameter $\chi(k_{\perp 0})=0.1$ in a plasma with isotropic equilibrium
Maxwellian velocity distributions with $T_{i}/T_{e}=9$ and $\beta_{i}=1$: (a)
One-dimensional perpendicular magnetic energy spectrum
$E_{B_{\perp}}(k_{\perp})$, (b) anisotropic cascade through wave vector space
given by , $k_{\parallel}(k_{\perp})$, (c) energy transfer rates to the ions
$Q_{i}(k_{\perp})$ (blue) and electrons $Q_{e}(k_{\perp})$ (magenta), and (d)
integrated energy damping rates by ions
$\int_{0}^{k_{\perp}}dk^{\prime}_{\perp}Q_{i}(k^{\prime}_{\perp})$ (blue) and
by electrons
$\int_{0}^{k_{\perp}}dk^{\prime}_{\perp}Q_{e}(k^{\prime}_{\perp})$ (magenta),
all plotted as a function of the perpendicular wavenumber $k_{\perp}\rho_{i}$.

Determining this balance of cascade and damping rates as a function of scale,
and ultimately using that to determine the partitioning of turbulent energy
among species and degrees of freedom—parameterized by $Q_{i}/Q_{e}$,
$Q_{\perp,i}/Q_{i}$, and $Q_{\perp,e}/Q_{e}$—has been attempted in empirical
turbulent cascade models (Pao, 1965; Howes et al., 2008a; Podesta et al.,
2010; Howes et al., 2011; Zhao et al., 2013; Schreiner & Saur, 2017). An
example of the modeling of the turbulent energy cascade and its kinetic
dissipation is presented in Fig. 5 from the results of the _weakened cascade
model_ by Howes et al. (2011). In addition to balancing the physics of the
local cascade rate $\varepsilon(k_{\perp})$ with the kinetic damping rates due
to ions $Q_{i}(k_{\perp})$ and electrons $Q_{e}(k_{\perp})$ via the Landau
resonance, this cascade model also includes the contributions to the local
cascade rate from nonlocal interactions by large-scale shearing and small-
scale diffusion. The model results shown here in Fig. 5 are from the same
calculation as presented in Fig. 8 of Howes et al. (2011), where we reproduce
(a) the one-dimensional perpendicular magnetic energy spectrum
$E_{B_{\perp}}(k_{\perp})$ and (b) the anisotropic cascade through wave vector
space characterized by the characteristic parallel wavenumber as a function of
the perpendicular wavenumber, $k_{\parallel}(k_{\perp})$. This model for a
plasma with $\beta_{i}=1$ and isotropic equilibrium Maxwellian velocity
distributions with $T_{i}/T_{e}=9$ is driven weakly with $\chi(k_{\perp
0})=0.1$ at $k_{\perp 0}\rho_{i}=k_{\parallel}\rho_{i}=10^{-5}$, taking the
model Kolmogorov constants to be $C_{1}=1.4$ and $C_{2}=1.0$.

In Fig. 5(c) is plotted the rate of energy transfer to the ions
$Q_{i}(k_{\perp})$ (blue) and to the electrons $Q_{e}(k_{\perp})$ (magenta) as
a function of $k_{\perp}\rho_{i}$, both normalized by the turbulent cascade
rate at the driving scale $\varepsilon_{0}$. This shows that the ion kinetic
damping occurs at the ion kinetic scales $k_{\perp}\rho_{i}\sim 1$, whereas
the electron damping rises monotonically from the ion scales, peaking around
$k_{\perp}\rho_{i}\sim 30$ for this case. To obtain the net partitioning of
energy between ions and electrons $Q_{i}/Q_{e}$, it is necessary to integrate
these energy damping rates over the entire cascade, as shown in Fig. 5(d): the
turbulent energy cascade rate $\varepsilon(k_{\perp})$ (black) is ultimately
terminated at $k_{\perp}\rho_{i}\gg 1$ by the integrated energy damping rates
by ions $\int_{0}^{k_{\perp}}dk^{\prime}_{\perp}Q_{i}(k^{\prime}_{\perp})$
(blue) and by electrons
$\int_{0}^{k_{\perp}}dk^{\prime}_{\perp}Q_{e}(k^{\prime}_{\perp})$ (magenta).
Note that this turbulent cascade model determines the solutions presented in
Fig. 5 by computing the balance of the turbulent cascade rates with the local
kinetic damping rates, as depicted in Fig. 4. Integration of those damping
rates over the entire cascade yields the predicted partitioning of dissipated
turbulent energy between ions and electrons $Q_{i}/Q_{e}$.

The energy transfer rates $Q_{i}(k_{\perp})$ and $Q_{e}(k_{\perp})$ that
remove energy from the turbulent fluctuations and transfer it to the plasma
particles in turbulent cascade models require as input the kinetic damping
rates due to ions $\gamma_{i}(k_{\perp})$ and electrons
$\gamma_{e}(k_{\perp})$. In the case of the weakened cascade model shown in
Fig. 5, the damping rates used are collisionless damping rates via the Landau
($n=0$) resonance for the Alfvén wave mode (and its extension as a kinetic
Alfvén wave at scales $k_{\perp}\rho_{i}\gtrsim 1$) from the Vlasov-Maxwell
linear dispersion relation (Howes et al., 2011). Previous studies have called
into question whether resonant collisionless wave-particle interactions, such
as Landau damping, can effectively damp turbulent fluctuations in the presence
of the strong nonlinear interactions that mediate the turbulent cascade
(Plunk, 2013; Schekochihin et al., 2016). Using the field-particle correlation
technique (Klein & Howes, 2016; Howes et al., 2017; Klein et al., 2017),
analysis of kinetic simulations (Klein et al., 2017; Howes et al., 2018; Klein
et al., 2020) and spacecraft observations (Chen et al., 2019; Afshari et al.,
2021, 2023) of plasma turbulence find clear velocity-space signatures of
Landau damping. These findings indicate indeed that resonant wave-particle
interactions do play a role, and possibly a dominant role (Afshari et al.,
2021), in the damping of space and astrophysical plasma turbulence, settling
the question of whether Landau damping can effectively damp strong plasma
turbulence.

In quantifying how a particular damping mechanism competes with the nonlinear
turbulent cascade, the key dimensionless quantity to evaluate is the ratio of
linear damping rate due to nonlinear cascade rate $\gamma/\omega_{nl}$ (Howes
et al., 2008a, 2011). Adopting the conjecture of critical balance between
linear and nonlinear timescales in strong plasma turbulence
$\omega\sim\omega_{nl}$ (Goldreich & Sridhar, 1995; Howes et al., 2008a, 2011;
Mallet et al., 2015), the key dimensionless measure of the importance of that
particular mechanism is the ratio $\gamma/\omega$, which can be estimated
using the linear theory for each of the different damping mechanisms.

###  5.3 The Relevance of Linear Physics to Strong Plasma Turbulence

The concept of critical balance in plasma turbulence (Goldreich & Sridhar,
1995; Mallet et al., 2015)—which suggests that the nonlinear timescales of the
turbulent cascade remain in approximate balance with the timescales of the
linear wave physics, $\omega_{nl}\sim\omega$—implies that the physics of the
linear plasma response remains valid in strong plasma turbulence. But, perhaps
in analogy with the case of incompressible hydrodynamic turbulence in which no
linear response exists (Howes, 2015), the validity of linear physics has been
called into question in the presence of strong plasma turbulence (Matthaeus et
al., 2014). Although there exist numerous counterexamples in which linear
physics properties have been shown to be relevant to strong plasma turbulence
(Maron & Goldreich, 2001; Cho & Lazarian, 2003; Alexandrova et al., 2008;
Howes et al., 2008b, a; Svidzinski et al., 2009; Sahraoui et al., 2010; Hunana
et al., 2011; Howes et al., 2011; TenBarge et al., 2012; Chen et al., 2013),
the question has persisted. A more recent careful study of kinetic numerical
simulations and spacecraft observations has shown clearly that both small-
amplitude and large-amplitude fluctuations, including variations interpreted
to be coherent structures, preserve the linear properties of Alfvénic
fluctuations (Grošelj et al., 2019), hopefully providing sufficient evidence
to establish that relevance of linear physics to strong plasma turbulence.

The key to predicting qualitatively which kinetic damping mechanisms are
likely to dominate for a given instance of turbulence is to determine how the
various collisionless damping mechanisms in §5.1 depend on the fundamental
dimensionless plasma and turbulence parameters. For a fully ionized proton and
electron plasma with anisotropic (bi-Maxwellian) equilibrium velocity
distribution functions (Stix, 1992; Quataert, 1998; Quataert & Gruzinov, 1999;
Swanson, 2003; Klein & Howes, 2015), the complex eigenfrequencies determined
by the linear Vlasov-Maxwell dispersion relation depend on the dimensionless
parameters that characterize the nature of the turbulent fluctuations through
the typical wave vector components in a plane-wave decomposition of the
turbulent fluctuations $k_{\perp}\rho_{i}$ and $k_{\parallel}\rho_{i}$ and the
plasma parameters through the parallel ion plasma beta $\beta_{\parallel,i}$,
the parallel ion-to-electron temperature ratio
$\tau_{\parallel}=T_{\parallel,i}/T_{\parallel,e}$, the temperature
anisotropies of each species $A_{i}=T_{\perp,i}/T_{\parallel,i}$ and
$A_{e}=T_{\perp,e}/T_{\parallel,e}$, and the ratio $v_{t\parallel,i}/c$. Thus,
these dependencies can be expressed through the general relation

$\omega=\omega_{VM}^{aniso}(k_{\perp}\rho_{i},k_{\parallel}\rho_{i},\beta_{\parallel,i},\tau_{\parallel},A_{i},A_{e},v_{ti}/c).$
(21)

In the limit of isotropic equilibrium velocity distributions for the ions and
electrons—$A_{i}=1$, $A_{e}=1$, and $\tau=T_{i}/T_{e}$—this set of
dependencies is reduced to

$\omega=\omega_{VM}(k_{\perp}\rho_{i},k_{\parallel}\rho_{i},\beta_{i},\tau,v_{ti}/c).$
(22)

For both the isotropic and anisotropic equilibirum temperatures, in the non-
relativistic limit $v_{ti}/c\ll 1$ typical of many space and astrophysical
plasmas of interest, the linear Vlasov-Maxwell dispersion relation is
practically independent of $v_{ti}/c$, so this parameter may effectively be
dropped for the case of non-relativistic plasma turbulence.

In the limit of a sufficiently large inertial range of turbulence, the MHD
turbulence scaling parameters given by (1)–(4) predict that turbulent
flucuations at the small-scale end of the inertial range, where the kinetic
damping mechanisms are believed to become significant, will be anisotropic
with $k_{\parallel}\ll k_{\perp}$ and small amplitude relative to the
equilibrium magnetic field $|\delta\mathbf{B}|\ll B_{0}$ (Howes et al., 2006,
2008a; Schekochihin et al., 2009). In this limit, the turbulent frequencies of
the incompressible Alfvénic fluctuations directly observed to dominate
turbulence in space plasmas (Belcher & Davis, 1971; Tu & Marsch, 1995;
Alexandrova et al., 2008; Howes et al., 2012a; Bruno & Carbone, 2013) have
frequencies much smaller than the ion cyclotron frequency,
$\omega/\Omega_{i}=k_{\parallel}d_{i}\ll 1$, where
$d_{i}=v_{A}/\Omega_{i}=\rho_{i}/\sqrt{\beta_{i}}$ is the ion inertial length.
In this limit, the turbulent fluctuations at the small-scale end of the
inertial range satisfy the conditions for the applicabilty of gyrokinetic
theory (Rutherford & Frieman, 1968; Taylor & Hastie, 1968; Catto, 1978;
Antonsen & Lane, 1980; Catto et al., 1981; Frieman & Chen, 1982; Dubin et al.,
1983; Hahm et al., 1988; Brizard, 1992; Howes et al., 2006; Brizard & Hahm,
2007; Schekochihin et al., 2009). Gyrokinetic theory can be rigorously derived
from plasma kinetic theory under the limits of low frequency
($\omega/\Omega_{i}\ll 1$) and anisotropic $k_{\parallel}\ll k_{\perp}$
fluctuations in the non-relativisitc limit $v_{ti}/c\ll 1$. In the gyrokinetic
approximation, the physics of the the cyclotron resonances and of the fast
magnetosonic fluctuations—and their kinetic extension to whistler fluctuations
(Howes et al., 2014)—are ordered out of the system. Gyrokinetics retains the
physics of the $\mathbf{E}\times\mathbf{B}$ nonlinearity that mediates the
turbulent cascade to small scales, full finite-Larmor-radius effects, and the
collisionless damping associated with the $n=0$ Landau resonance. In this
reduced limit, the complex eigenfrequencies of the linear gyrokinetic
dispersion relation for isotropic equilibrium velocity distributions are
linearly dependent on the parallel wavenumber, and have the significantly
reduced dimensionless parameter dependency given by (Howes et al., 2006)

$\overline{\omega}\equiv\frac{\omega}{k_{\parallel}v_{A}}=\overline{\omega}_{GK}(k_{\perp}\rho_{i},\beta_{i},\tau).$
(23)

Although the ion-to-electron mass ratio $\mu=m_{i}/m_{e}$ is not a parameter
that varies in the physical world, many numerical studies of plasma turbulence
choose to employ a reduced mass ratio for computational efficiency. Therefore,
we choose to include the $\mu$ dependence of different proposed turbulent
damping mechanisms in our calculations below to determine how the use of a
reduced mass ratio may impact the partitioning of turbulent energy removed by
different damping mechanisms. Note also that it is important to employ a
sufficiently large mass ratio to achieve a separation between the ion and
electron scales, with a previous study suggesting $m_{i}/m_{e}\geq 32$ is
necessary to model the distinct ion and electron responses to the
electromagnetic fluctuations due to this scale separation (Howes et al.,
2018).

###  5.4 Landau Damping and Transit-Time Damping

To explore the two Landau-resonant collisionless damping mechansims—Landau
damping and transit-time damping—we take $n=0$ in the resonance condition in
(20), obtaining the condition that particles with parallel velocities equal to
the parallel phase velocity of the waves,
$v_{\parallel}=\omega/k_{\parallel}$, may resonantly exchange energy with the
electromagnetic waves.

For simplicity in this analysis, we will assume the isotropic temperature case
for a proton-electron plasma with linear wave properties given by the
isotropic linear Vlasov-Maxwell dispersion relation (22); extending these
calculations to allow for anisotropic temperatures can be done simply
following the same logical procedure and using the linear Vlasov-Maxwell
dispersion relation for anisotropic temperatures given by (21). Furthermore,
justified by observations in solar wind turbulence that incompessible Alfvénic
fluctuations energetically dominate the turbulence (Belcher & Davis, 1971; Tu
& Marsch, 1995; Alexandrova et al., 2008; Bruno & Carbone, 2013), we will
focus on the linear collisionless damping rates of the Alfvénic fluctautions,
including both Alfvén waves in the MHD limit at $k_{\perp}\rho_{i}\ll 1$, and
kinetic Alfvén waves at $k_{\perp}\rho_{i}\gtrsim 1$. Consideration of the
damping of compressible turbulent fluctuations, specifically fast and slow
magnetosonic fluctuations at $k_{\perp}\rho_{i}\ll 1$ (Howes et al., 2012a),
is easily performed using the same procedure with the linear dispersion
relations for these alternative wave modes777Although, it is important to keep
in mind that the recent identified nonlinear effect of anti-phase-mixing in
turbulent plasmas can suppress the collisionless damping predicted by linear
kinetic theory for the compressible wave modes in the inertial range of scales
at $k_{\perp}\rho_{i}\ll 1$ (Schekochihin et al., 2016; Parker et al., 2016)..

Normalizing the resonant parallel velocity by the thermal velocity for both
species in turn and converting the result into our dimensionless parameters
for the isotropic temperature case, we obtain the key relations

$\frac{v_{\parallel}}{v_{ti}}=\frac{\omega}{k_{\parallel}v_{ti}}=\overline{\omega}\beta_{i}^{-1/2},\quad\quad\quad\quad\frac{v_{\parallel}}{v_{te}}=\frac{\omega}{k_{\parallel}v_{te}}=\overline{\omega}\left(\frac{\tau}{\beta_{i}\mu}\right)^{1/2},$
(24)

where $\overline{\omega}\equiv\omega/k_{\parallel}v_{A}$. For Alfvén waves,
this normalized frequency $\overline{\omega}$ is well approximated by (Howes
et al., 2014)

$\overline{\omega}(k_{\perp}\rho_{i},\beta_{i},\tau)\equiv\frac{\omega}{k_{\parallel}v_{A}}=\sqrt{1+\frac{(k_{\perp}\rho_{i})^{2}}{\beta_{i}+2/(1+1/\tau)}},$
(25)

an expression that is valid in the limits that the ion plasma beta is not very
low $\beta_{i}\gg\mu^{-1}$ and that the wave frequency remains smaller than
the ion cyclotron frequency $\omega/\Omega_{i}\ll 1$, equivalent to the
condition $k_{\parallel}d_{i}\ll 1$, or
$k_{\parallel}\rho_{i}\ll\beta_{i}^{1/2}$.

Figure 6: Diagram demonstrating how Landau resonant damping qualitatively
depends on where the parallel phase velocity $\omega/k_{\parallel}$ of the
wave falls within the velocity distribution of particle species $s$, with weak
damping for $\omega/k_{\parallel}\ll v_{ts}$ or $\omega/k_{\parallel}\gg
v_{ts}$, and strong damping for $\omega/k_{\parallel}\sim v_{ts}$.

Whether damping by the Landau-resonant mechanisms is weak or strong primarily
depends on where the parallel wave phase velocity falls within the equilibrium
particle velocity distribution, as quantified by the expressions in (24). The
application of the weak-growth-rate approximation to the kinetic theory for
Landau damping (Krall & Trivelpiece, 1973) indicates that the linear damping
rate is proportional to the slope of the distribution at the resonant
velocity, as depicted in Fig. 6. For phase velocities that fall deep in the
core of the velocity distribution $\omega/k_{\parallel}\ll v_{ts}$ (blue
dotted) where the slope of $f_{s}(v_{\parallel})$ is flat, the damping is weak
because nearly equal numbers of particles with
$v_{\parallel}<\omega/k_{\parallel}$ are accelerated by the wave as particles
with $v_{\parallel}>\omega/k_{\parallel}$ are decelerated by the wave, leading
to little net energy exchange. Similarly, for phase velocities far out in the
tail of the velocity distribution $\omega/k_{\parallel}\gg v_{ts}$ (black
dashed), there are few particles to resonate with the waves and the slope of
$f_{s}(v_{\parallel})$ is also flat, so the damping is weak. Only if the
parallel phase velocity falls in the steep region of the velocity distribution
at $\omega/k_{\parallel}\sim v_{ts}$ (red solid) is there a significant net
energy transfer to the particles: more particles with
$v_{\parallel}<\omega/k_{\parallel}$ are accelerated by the wave than those
particles with $v_{\parallel}>\omega/k_{\parallel}$ are decelerated by the
wave, so that particles gain net energy at the expense of the wave.

The ultimate goal is to characterize the turbulent plasma heating—
specifically $Q_{i}/Q_{e}$, $Q_{\perp,i}/Q_{i}$,and $Q_{\perp,e}/Q_{e}$—due to
Landau damping and transit-time damping as a function of the fundamental
dimensionless parameters of turbulence. To accomplish this goal, it is
necessary to integrate the damping rates over the full range of the turbulent
cascade, as discussed in §5.2 and illustrated in Figures 4 and 5. This
integration over all scales effectively eliminates the dependence of the
damping mechanisms on the components of the wavevector $k_{\perp}\rho_{i}$ and
$k_{\parallel}\rho_{i}$, leaving only the desired dependence on the plasma and
turbulence parameters enumerated in Tables 2 and 3.

With the scale dependence on $k_{\perp}\rho_{i}$ eliminated, the ion Landau
damping (iLD) and ion transit-time damping (iTTD) rates depend primarily on
$\beta_{i}$, and the electron Landau damping (eLD) and electron transit-time
damping (eTTD) rates depend on $\beta_{i}$, $\tau$, and $\mu$, as shown by
(24). Note that the dependence of eLD and eTTD on $\tau$ and $\mu$ primarily
arises from those parameters controlling the scaling of the ion to electron
Larmor radius scales, $\rho_{i}/\rho_{e}\sim(\tau\mu)^{1/2}$.

It is worthwhile pointing out that the physics of collisionless damping via
the Landau resonance is fully captured in the gyrokinetic approximation, with
the linear wave frequencies and damping rates both linearly proportional to
$k_{\parallel}$, as shown in (23). Thus, the ratio of the damping-to-cascade
rate888Here we have adopted the condition of criticial balance for strong
plasma turbulence $\omega_{nl}\simeq\omega$ (Goldreich & Sridhar, 1995; Mallet
et al., 2015)., $\gamma/\omega_{nl}\simeq\gamma/\omega$, is formally
independent of $k_{\parallel}$, leaving the expected parameter dependence on
the plasma parameters $\beta_{i}$ and $\tau$, along with the mass ratio $\mu$.

Note that dependence of the linear gyrokinetic dispersion relation on
$(k_{\perp}\rho_{i},\beta_{i},\tau)$ in (23) is consistent with the results of
our analysis of the parameter dependence of the Landau-resonant damping
mechanisms. The linear gyrokinetic dispersion relation also determines the
amplitudes and phases of the key electromagnetic fields $E_{\parallel}$ and
$\delta B_{\parallel}$ that mediate the energization by Landau damping and
transit-time damping. As done in the _H10_ turbulent heating model, the linear
dispersion relation can simply be used to calculate directly the linear
collisionless damping rates as a function of
$(k_{\perp}\rho_{i},\beta_{i},\tau)$, but these linear damping rate
predictions may be modified under strongly turbulent plasma conditions.
However, field-particle correlation analyses of numerical simulations (Klein
et al., 2017, 2020; Huang et al., 2024) and spacecraft observations (Chen et
al., 2019; Afshari et al., 2021) of plasma turbulence have clearly shown that
Landau damping and transit-time damping do indeed play a role in the damping
of strong plasma turbulence.

One final step is to estimate how the ion-to-electron heating ratio
$Q_{i}/Q_{e}$ depends on the remaining key parameters $\beta_{i}$, $\tau$, and
$\mu$. As illustrated in Fig. 5(c) and Fig. 7, the ion damping rate
$\gamma_{i}$ peaks at scales $k_{\perp}\rho_{i}\sim 1$, and this Landau-
resonant damping onto the ions is essentially independent of the ratio ion-to-
electron temperature ratio. The electron damping rate $\gamma_{e}$, on the
other hand, increases monotonically as the perpendicular wavevector increases
to the electron Larmor radius scale $k_{\perp}\rho_{e}\rightarrow 1$,
eventually yielding a normalized collisionless damping rate that is
sufficiently strong to terminate the turbulent cascade with
$\gamma_{e}/\omega\rightarrow 1$, as shown in Fig. 7. Thus, any energy in the
turbulent cascade that passes beyond the ion kinetic scales at
$k_{\perp}\rho_{i}\sim 1$ will ultimately be transferred to the electrons.
Since the ion damping rate depends dominantly on $\beta_{i}$ with very weak
dependence on $\tau$, yielding $Q_{i}(\beta_{i})$ to lowest order, this means
that the ratio $Q_{i}/Q_{e}$ is effectively determined by the ion plasma beta
alone given a total turbulent energization rate that eventually must balance
the total cascade rate $Q\sim\varepsilon_{0}$, since it can be shown that
$Q_{e}/Q_{i}=\varepsilon_{0}/Q_{i}(\beta_{i})-1$. This dominant dependence of
the Landau-resonant collisionless damping mechanisms on $\beta_{i}$ is
consistent with existing turbulent heating models (Howes, 2010; Kawazura et
al., 2019).

One final unexpected physical behavior of Landau damping and transit-time
damping is also illustrated in Fig. 7. For isotropic Maxwellian velocity
distributions of the ions and electrons, the imaginary component of the
complex frequency, $\omega_{c}\equiv\omega-i\gamma$, determined by the linear
dispersion relation is always negative (solid black curve), corresponding to a
damping of the wave with $\gamma>0$. However, for this particular plasma
parameter choice of $\beta_{i}=3$, $T_{i}/T_{e}=1$, and $m_{i}/m_{e}=1836$,
the individual energy transfer rates to the ions through Landau damping
mediated by $E_{\parallel}$ are actually negative over the range
$k_{\perp}\rho_{i}\lesssim 2$ (red dotted), indicating a net transfer of
energy from the ions to the wave fields when averaged over the full period
(or, alternatively, over $2\pi$ in phase) of the Alfvén wave. But the
contribution from transit-time damping mediated by $E_{\perp}$ (Huang et al.,
2024; Howes et al., 2024) is positive and much larger (green long dashed),
leading to a net damping of the wave ($\gamma_{i}>0$, thin blue) by the sum of
the Landau-resonant interactions of ions with the Alfvén wave. Whether the
wave-period-averaged energy transfer to ions is positive or negative depends
on the relative phases of $E_{\parallel}$ and $E_{\perp}$ to the self-
consistently determined ion current density $\mathbf{j}_{i}$. For the case of
$\beta_{i}=3$ pictured in Fig. 7, at scales $k_{\perp}\rho_{i}\lesssim 2$ (red
dotted), this relative phase $\delta\phi$ between $E_{\parallel}$ and
$j_{\parallel,i}$ must fall in the range $\pi/2\leq\delta\phi\leq 3\pi/2$
(corresponding to wave growth at the expense of ion energy), while the
relative phase between $E_{\perp}$ and $j_{\perp,i}$ must fall in the range
$-\pi/2\leq\delta\phi\leq\pi/2$ (corresponding to wave damping and ion energy
gain). Because each individual mechanism can possibly lead to negative energy
transfer while the sum always leads to damping, it suggests that perhaps a
separation of the particle energization rates by Landau damping and transit-
time damping is not particularly important in unraveling the turbulent heating
of the plasma (especially since both of these mechansims ultimately lead to
energization of the parallel degree of freedom of the particles), so that only
the summed effect of the Landau-resonant collisionless damping rates is
important.

Figure 7: Plot of the normalized damping rate $\gamma/\omega$ (black dashed)
for the Alfvén and kinetic Alfvén wave over the perpendicular wavenumber range
$10^{-2}\leq k_{\perp}\rho_{i}\leq 10^{2}$ for a plasma with $\beta_{i}=3$,
$T_{i}/T_{e}=1$, and $m_{i}/m_{e}=1836$, where the complex eigenfrequency is
given by $\omega_{c}=\omega-i\gamma$, so that $\gamma>0$ corresponds to wave
damping and $\gamma<0$ to wave growth. The contributions from the ion damping
rate $\gamma_{i}$ (thin blue) and electron damping rate $\gamma_{e}$ (thin
magenta) are separately plotted, with this decomposition failing at strong
damping rates with $\gamma/\omega\gtrsim 0.5$, occurring for
$k_{\perp}\rho_{i}\gtrsim 40$. Ion damping due to Landau damping
$\gamma_{i,LD}>0$ (red dashed) and transit-time damping $\gamma_{i,TTD}>0$
(green long dashed) are separated, showing a region of net energy transfer
_from_ ions _to_ the wave for ion Landau damping with $\gamma_{i,LD}>0$ (red
dotted) at $k_{\perp}\rho_{i}\leq 2$.

Summarizing our results for Landau-resonant collisionless damping of the
turbulent cascade, the resulting ion-to-electron heating is predicted to be
dominantly a function of the ion plasma beta, $Q_{i}/Q_{e}(\beta_{i})$, with a
weak dependence on the ion-to-electron temperature ratio $\tau$ and mass ratio
$\mu$. Both Landau damping, mediated by the parallel component of the electric
field $E_{\parallel}$, and transit-time damping, mediated by the perpendicular
component of the electric field $E_{\perp}$, lead to energization of the
degrees of freedom parallel to the magnetic field, so $Q_{\perp,i}/Q_{i}=0$
and $Q_{\perp,e}/Q_{e}=0$.

###  5.5 Cyclotron Damping

For the fundamental ($n=1$) cyclotron resonance, the resonance condition in
(20) simplifies to $(\omega-k_{\parallel}v_{\parallel})/\omega_{s}=1$. The
numerator involves the fluctuation frequency modified by the Doppler shift of
the particles motion through spatial variations parallel to the magnetic
field. When that Doppler-shifted frequency matches the ion cyclotron
frequency, the perpendicular compnent of the electric field $E_{\perp}$ can
resonantly interact with the particles to energize the particles in the
directions perpendicular to the magnetic field, yielding
$Q_{\perp,s}/Q_{s}=1$. To estimate the conditions required for cyclotron
resonant damping, it is often sufficient to require the condition
$\omega/\omega_{s}\rightarrow 1$.

In typical space and astrophysical plasmas, under the conditions derived below
the turbulent cascade can reach the ion cyclotron frequency
$\omega/\omega_{i}\rightarrow 1$, but since the electron cyclotron frequency
is higher by the mass ratio, $\Omega_{e}/\Omega_{i}=m_{i}/m_{e}$, there is
typically little energy remaining in the cascade at the very high frequencies
needed to transfer is significant amount of energy to the electron through
their cyclotron resonance. Thus, here we consider only the physics of the ion
cyclotron resonance as a possible physical mechanism for the dissioation of
plasma turbulence.

Evaluating the conditions needed to obtain $\omega/\omega_{i}\rightarrow 1$ in
terms of the dimensionless parameters of turbulence in the isotropic
temperature case, we find

$\frac{\omega}{\Omega_{i}}=\overline{\omega}k_{\parallel}d_{i}=\overline{\omega}\frac{k_{\parallel}\rho_{i}}{\beta_{i}^{1/2}},$
(26)

where $\overline{\omega}$ is defined by (25). Note that the middle expression
indicates that the most natural unit to normalize the length scale parallel to
the magnetic field is the ion inertial length
$d_{i}=c/\omega_{pi}=v_{A}/\Omega_{i}$, but to maintain consistency with the
length normalizations in Tables 2 and 3, we substitute
$d_{i}=\rho_{i}\beta_{i}^{-1/2}$. Utilizing the general scaling of the
parallel wavenumber in terms of the isotropic driving wavenumber,

$k_{\parallel}\rho_{i}=(k_{0}\rho_{i})^{(1+\alpha)/(3+\alpha)}(k_{\perp}\rho_{i})^{2/(3+\alpha)},$
(27)

to determine the frequency achieved when the cascade reaches
$k_{\perp}\rho_{i}=1$, we obtain the condition

$\frac{\omega}{\Omega_{i}}=\frac{\overline{\omega}}{\beta_{i}^{1/2}}(k_{0}\rho_{i})^{(1+\alpha)/(3+\alpha)}.$
(28)

This indicates clearly that ion cyclotron damping is more likely to occur for
small ion plasma beta $\beta_{i}\ll 1$ or a smaller turbulent inertial range
with a value of $k_{0}\rho_{i}$ this is not too small. For the GS95 scaling
with $\alpha=0$, we obtain a condition for the onset of ion cyclotron damping
of $k_{0}\rho_{i}=\beta_{i}^{3/2}/\overline{\omega}^{3}$; for the B06 scaling
with $\alpha=1$, we obtain
$k_{0}\rho_{i}=\beta_{i}^{/2}/\overline{\omega}^{2}$. For simplicity, note
that at the small-scale end of the inertial range with $k_{\perp}\rho_{i}\sim
1$, the normalized wave frequency $\overline{\omega}\sim 1$.

In summary, ion cyclotron damping is predicted to play a role in the damping
of plasma turbulence for conditions with a smaller inertial range and a low
plasma beta values $\beta_{i}\ll 1$. The nature of the resulting plasma
heating is strictly ion energization $Q_{i}/Q=1$ in the perpendicular degrees
of freedom $Q_{\perp,i}/Q_{i}=1$.

Note that recent theoretical and numerical study of turbulence with a
sufficiently large amount of imbalance $Z_{0}^{+}/Z_{0}^{-}\gg 1$, such as may
be encountered in the inner heliosphere, has uncovered a new phenomenon
denoted the _helicity barrier_ , which can lead to the generation of
sufficiently high-frequency fluctuations, specifically ion cyclotron waves, to
lead to iCD energization of the plasma particles (Squire et al., 2022, 2023).

###  5.6 Stochastic Heating

The stochastic heating of ions has long been proposed as a means to energize
ions in turbulent plasmas (Johnson & Cheng, 2001; Chen et al., 2001; White et
al., 2002; Voitenko & Goossens, 2004; Bourouaine et al., 2008), with more
recent work providing significant theoretical and numerical evidence for ion
stochastic heating as a significant mechanism for turbulent damping in space
and astrophysical plasmas (Chandran et al., 2010; Chandran, 2010; Chandran et
al., 2011; Bourouaine & Chandran, 2013; Chandran et al., 2013; Klein &