1 Introduction
Energetic particles, such as those utilized to heat plasma in magnetic confinement experiments, are subject to being redistributed by turbulent transport. The net flux of these non-Maxwellian energetic particles between flux surfaces is a strong function of energy (Hauff et al. Reference Hauff, Pueschel, Dannert and Jenko2009), and this causes the lower-energy (though still not thermalized) part of the velocity distribution to be ‘carved out’, sometimes leaving a bump-on-tail feature (Wilkie et al. Reference Wilkie, Abel, Landreman and Dorland2016). Here, a simplified model which describes this effect is introduced, with a particular focus on alpha particles produced by the deuterium–tritium (DT) fusion reaction. It is a generalization of the analytic slowing-down distribution of Gaffey (Gaffey Reference Gaffey1976; Helander & Sigmar Reference Helander and Sigmar2002) (which will repeatedly be referred to as the ‘classical’ slowing-down distribution) that includes an additional term mimicking the velocity dependence of microturbulent transport. A related analytic distribution was also previously derived by Anderson, Batistoni & Anderson (Reference Anderson, Batistoni and Anderson1991) and Sigmar, Gormley & Kamelander (Reference Sigmar, Gormley and Kamelander1993), wherein anomalous transport was treated with a constant diffusion coefficient, rather than one with velocity dependence.
Alpha particles are born isotropically from reacting deuterium and tritium nuclei at an energy of
$E_{\unicode[STIX]{x1D6FC}}=m_{\unicode[STIX]{x1D6FC}}v_{\unicode[STIX]{x1D6FC}}^{2}/2=3.52~\text{MeV}$
. The particle density, mass, temperature and thermal speed of a species
$s$
are denoted
$n_{s}$
,
$m_{s}$
,
$T_{s}$
and
$v_{ts}=\sqrt{2T_{s}/m_{s}}$
, respectively. Plasma heating is caused by the friction of these particles against the bulk ions/electrons (which have characteristic thermal speeds
$v_{ti/e}=\sqrt{2T_{i/e}/m_{i/e}}$
, respectively). As they slow down and their magnetic orbits become smaller, alpha particles become more and more susceptible to ion-scale microturbulence. As they transport outward, they continue slowing down on outer flux surfaces. This, combined with the fact that fewer alpha particles are produced near the cooler edge of the plasma, means that the alpha particle energy leaving the plasma will be relatively small (Kurki-Suonio et al.
Reference Kurki-Suonio, Asunta, Hirvijoki, Koskela, Snicker, Hauff, Jenko, Poli and Sipilä2011). Even so, alpha particles will be redistributed radially, changing the global heating profile. Alphas then eventually become thermalized in equilibrium with the bulk ions, creating a population of helium ‘ash’ with densityFootnote
1
$n_{\text{ash}}>n_{\unicode[STIX]{x1D6FC}}$
and temperature
$T_{\text{ash}}=T_{i}$
. This population can then be treated as a Maxwellian impurity.
All of the aforementioned effects are captured by the kinetic transport equation from the low-collisionality ordering of gyrokinetics for the equilibrium distribution
$f_{\unicode[STIX]{x1D6FC}}$
:
$$\begin{eqnarray}\displaystyle & & \displaystyle \frac{\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}t}-\frac{1}{V^{\prime }}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}r}\left(V^{\prime }D_{rr}\frac{\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}r}+V^{\prime }D_{rv}\frac{\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}v}\right)-\frac{1}{v^{2}}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}v}\left(v^{2}D_{vr}\frac{\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}r}+v^{2}D_{vv}\frac{\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}v}\right)\nonumber\\ \displaystyle & & \displaystyle \quad =C[f_{\unicode[STIX]{x1D6FC}}]+S_{\unicode[STIX]{x1D6FC}}.\end{eqnarray}$$
The collision operator
$C$
includes collisions of trace particles against the Maxwellian bulk, and
$S_{\unicode[STIX]{x1D6FC}}(r,v)$
is the source of alpha particles from fusion reactions. Flux surfaces are labelled by the radial coordinate
$r$
, defined as the surface half-width at the height of the magnetic axis, and enclose a volume
$V(r)$
, and primes denote differentiation with respect to radius. The radius of the last closed flux surface is the plasma minor radius
$a$
. The various turbulent diffusion coefficients
$D_{xx}$
are functions of radius and speed
$v$
, and these can be calculated from nonlinear gyrokinetic simulations. Equation (1.1) is valid for either trace impurities or for energetic particles (Wilkie et al.
Reference Wilkie, Iantchenko, Abel, Highcock and Pusztai2018). Both of these conditions are typically satisfied for alpha particles (Wilkie et al.
Reference Wilkie, Abel, Highcock and Dorland2015).
The main body of this work derives (§ 2), benchmarks and discusses (§ 3) a particularly simple analytic distribution using a local loss term that models the approximate effect of energy-dependent transport. The approximate velocity dependence used is derived from first principles in appendix A. Appendix B relaxes the isotropic approximation and leads to a velocity dependence with the same energy scaling, but a different pitch-angle scaling than Hauff et al. (Reference Hauff, Pueschel, Dannert and Jenko2009). An analytic model distribution with pitch-angle dependence is thereby also derived in appendix B. Appendix C relaxes the local-sink approximation to solve analytically for the radial transport in a plasma that is uniform except for the energetic particle source. These generalizations are combined to obtain an analytic distribution in three-dimensional phase space
$f_{SD,3D}(r,v,\unicode[STIX]{x1D709})$
in (C 8).
2 Model distribution function
An approximate local solution to (1.1) will be found after a series of assumptions. In this section, it will be assumed that the energetic particles are isotropic in velocity space, but this is relaxed in the more general treatment of appendix B. Also, in the limit where
$v_{ti}\ll v\ll v_{te}$
, the collision operator acting on an isotropic alpha distribution reduces to the simplified form
$C[f_{\unicode[STIX]{x1D6FC}}]\approx \unicode[STIX]{x1D70F}_{s}^{-1}v^{-2}(\unicode[STIX]{x2202}/\unicode[STIX]{x2202}v)(v_{c}^{3}+v^{3})f_{\unicode[STIX]{x1D6FC}}$
, where the critical speed is defined as
with the sum over all non-trace ionic species with charge number
$Z_{i}$
. The slowing-down time is:
where
$e$
is the elementary charge and
$\ln \unicode[STIX]{x1D6EC}_{\unicode[STIX]{x1D6FC}e}$
is the Coulomb logarithm for collisions of alpha particles against electrons.
Previous studies (Wilkie et al.
Reference Wilkie, Pusztai, Abel, Dorland and Fülöp2017) found that the
$D_{rv}$
,
$D_{vr}$
and
$D_{vv}$
coefficients play a sub-dominant role compared to radial diffusion (
$D_{rr}$
) at high energy, so this will be the only transport term kept from (1.1). To see that this is justified, consider the two source terms in the gyrokinetic equation:
$\unicode[STIX]{x1D714}e\langle \unicode[STIX]{x1D719}\rangle _{\boldsymbol{R}}\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}{\mathcal{E}}$
and
$(c/B)k_{y}\langle \unicode[STIX]{x1D719}\rangle _{\boldsymbol{R}}\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}r$
for the case of ion-scale turbulence (
$\unicode[STIX]{x1D714}\sim v_{ti}/a$
is the frequency,
$k_{y}\sim \unicode[STIX]{x1D70C}_{i}^{-1}$
is a perpendicular wavenumber,
$\unicode[STIX]{x1D70C}_{i}$
is the thermal ion Larmor radius,
${\mathcal{E}}$
is energy,
$B$
is the magnetic field magnitude,
$c$
is the vacuum speed of light, and
$\langle \cdot \rangle _{\boldsymbol{R}}$
denotes the gyroaverage). The former term contributes to
$D_{rv}$
and the latter to
$D_{rr}$
. Therefore, the latter will dominate when
$T_{i}\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}{\mathcal{E}}\ll a\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}r$
, as is the case for energetic particles. Similar arguments apply to
$D_{vr}$
and
$D_{vv}$
when comparing their definitions to
$D_{rr}$
and
$D_{rv}$
(Wilkie et al.
Reference Wilkie, Abel, Landreman and Dorland2016). See appendix A for more details.
An analytic model is derived in appendix C which includes the global radial diffusion operator, but also includes other generous assumptions about the radial plasma profile. In this section, a local approximation will be made for the radial diffusion term, where the radial derivatives are parametrized by a length scale such that:
This is done to capture the approximate effect of transport on the local velocity distribution without needing to solve the radial transport equation as was done in Wilkie et al. (Reference Wilkie, Abel, Landreman and Dorland2016). In this framework, transport acts as a velocity-dependent sink of alpha particles. In reality, transported particles move to other flux surfaces, where they continue to slow down. Therefore, this approximation represents a ‘worst-case scenario’. This excludes cases where the energetic particle density might be greater at some radii due to transport from inner flux surfaces. In such cases, a local model is inadequate and full radial transport solution of (1.1) is required (Wilkie et al. Reference Wilkie, Abel, Landreman and Dorland2016).
Since one does not have the luxury of the full radial and velocity dependence of
$f_{\unicode[STIX]{x1D6FC}}$
, an a priori proxy is desired for the purposes of defining
$L_{\unicode[STIX]{x1D6FC}}$
. For the classical slowing-down distribution, the energetic particle density
$n_{\unicode[STIX]{x1D6FC}}$
is proportional to the source
$S_{0_{\unicode[STIX]{x1D6FC}}}$
, the number of alpha particles produced per unit volume per unit time. The latter is therefore chosen to represent the radial derivative of the energetic particle distribution so that
$\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}r\approx -f_{\unicode[STIX]{x1D6FC}}/L_{\unicode[STIX]{x1D6FC}}$
. This defines the radial scale length:
where
$D_{\unicode[STIX]{x1D6FC}}=D_{rr}(v=v_{\unicode[STIX]{x1D6FC}})$
is typically small, less than about. This form of
$L_{\unicode[STIX]{x1D6FC}}$
is kept in the main body of this work for generality, but it does have the drawback of choosing a relevant energy to represent the radial distribution of energetic particles. For an analytic model that relaxes this assumption, see appendix C. The neoclassical diffusion of slowing-down alphas due to collisional pitch-angle scattering (Catto Reference Catto1987) is expected to even smaller than turbulent transport, though this may be enhanced by a lack of perfect toroidal symmetry near the plasma edge (Kurki-Suonio et al.
Reference Kurki-Suonio, Asunta, Hirvijoki, Koskela, Snicker, Hauff, Jenko, Poli and Sipilä2011).
The model energy dependence of the diffusion coefficient used for this model is borrowed from Hauff et al. (Reference Hauff, Pueschel, Dannert and Jenko2009):
See appendix A for details on how such a scaling is found. In appendix B, it is shown that this is good approximation to the more general case where the pitch-angle dependence of
$D_{rr}$
is taken into account. This is because the range of pitch angles for which
$D_{rr}$
departs from this
$v^{-3}$
scaling is very narrow (Pueschel et al.
Reference Pueschel, Jenko, Schneller, Hauff, Günter and Tardini2012).
One benefit of the locality approximation is that it reduces the effect of transport to a single dimensionless parameter, which shall be defined as:
The appearance of the factor
$v_{\unicode[STIX]{x1D6FC}}^{3}/v_{c}^{3}$
is curious, but is motivated by the way in which
$b$
appears in the solution that follows. This dimensionless parameter can be simplified further by scaling the diffusion coefficient similarly to Pueschel et al. (Reference Pueschel, Jenko, Schneller, Hauff, Günter and Tardini2012). There, the diffusion coefficient at thermal energies was deemed proportional to the effective thermal diffusivity
$\unicode[STIX]{x1D712}_{\text{eff}}$
. If it is supposed that
$D_{rr}(v=v_{c})\approx \unicode[STIX]{x1D712}_{\text{eff}}$
, then one can approximate
$b\approx \unicode[STIX]{x1D712}_{\text{eff}}\unicode[STIX]{x1D70F}_{s}/L_{\unicode[STIX]{x1D6FC}}^{2}$
.
With the aforementioned approximations, equation (1.1) reduces to:
in steady state. In this model, the transport term diverges as
$v\rightarrow 0$
. More realistically, one could define a cutoff speed below which the diffusion coefficient is constant – this also results in an analytic solution. However, due to the other approximations made, equation (2.7) is invalid at low energy, where the distribution is dominated by thermalized helium ash. Therefore, a bare
$D_{rr}\propto v^{-3}$
scaling, without an ad hoc cutoff, is used without significant reservation. This is done with the knowledge that, like the classical slowing-down distribution, it is not valid for
$v\lesssim v_{ti}$
. Unfortunately, this precludes a useful analogous generalization for heated minority ions (Pusztai et al.
Reference Pusztai, Wilkie, Kazakov and Fülöp2016) since it is much more sensitive to the diffusion coefficient at low energy, for which a reduced model has yet to be developed. However, the generalization to injected auxiliary fast ions (e.g. neutral beams) is straightforward, even if their injection energy is not as high as for alpha particles.
The solution to (2.7) is:
Now the choice of including
$v_{\unicode[STIX]{x1D6FC}}^{3}/v_{c}^{3}$
in the definition of
$b$
, equation (2.6), becomes clear: this factor appears in the exponent of this distribution, and is thereby important in quantifying the effect of turbulence on the energetic particle distribution. Even for fusion-produced alpha particles, the source
$S_{\unicode[STIX]{x1D6FC}}$
has a significant spread in energy proportional to the geometric mean of the fuel temperature and the alpha particle birth energy (Brysk Reference Brysk1973). However, in order to obtain an analytic solution, it will be convenient to use a mono-energetic source such that
$S_{\unicode[STIX]{x1D6FC}}\approx (S_{0}/4\unicode[STIX]{x03C0}v_{\unicode[STIX]{x1D6FC}}^{2})\unicode[STIX]{x1D6FF}(v-v_{\unicode[STIX]{x1D6FC}})$
. The following analytic model distribution is thereby found:
where
$H$
is the Heaviside step function. Equation (2.9) reduces to the classical slowing-down distribution when transport is negligible:
$b=0$
. This is a generalization of the analytic distribution found by Sigmar et al. (Reference Sigmar, Gormley and Kamelander1993) to the case of energy-dependent diffusion as given by (2.5). Representative examples of this distribution at different values of
$b$
are shown in figure 1.
Comparing the classical Gaffey slowing-down distribution with the analytic transport-modified version – (2.9) – for different values of
$b$
. Lighter shades corresponding to higher relative turbulence intensity. The shaded area represents the Maxwellian ash population for
$n_{\text{ash}}=10n_{\unicode[STIX]{x1D6FC}}$
and
$T_{\text{ash}}=E_{\unicode[STIX]{x1D6FC}}/350$
. For these cases,
$v_{c}=0.3v_{\unicode[STIX]{x1D6FC}}$
.

3 Benchmarking and discussion
There are several curious features of the distribution described by (2.9) which are nevertheless backed up by numerical simulations and physical insight.
Firstly,
$F_{SD,\text{mod}}(v)$
has a local maximum at
$v_{\ast }=v_{\unicode[STIX]{x1D6FC}}(D_{\unicode[STIX]{x1D6FC}}\unicode[STIX]{x1D70F}_{s}/3L_{\unicode[STIX]{x1D6FC}}^{2})^{1/3}$
. Therefore, the distribution exhibits a ‘bump-on-tail’ similar to what was previously observed (Wilkie et al.
Reference Wilkie, Abel, Landreman and Dorland2016, Reference Wilkie, Pusztai, Abel, Dorland and Fülöp2017). In order for this inversion to be observable, this maximum must be well separated from the helium ash population. Therefore, a bump-on-tail caused by turbulent transport is possible only if:
Using
$T_{i}\sim 10~\text{keV}$
,
$\unicode[STIX]{x1D70F}_{s}\sim 1~\text{s}$
and
$L_{\unicode[STIX]{x1D6FC}}\sim 0.4~\text{m}$
, the presence of a bump requires
$D_{\unicode[STIX]{x1D6FC}}\gg 7\times 10^{-5}~\text{m}^{2}~\text{s}^{-1}$
. Estimates of
$D_{\unicode[STIX]{x1D6FC}}$
from the literature range from approximately
$6\times 10^{-4}$
to
$2\times 10^{-2}~\text{m}^{2}~\text{s}^{-1}$
using gyrokinetic simulations (Hauff et al.
Reference Hauff, Pueschel, Dannert and Jenko2009; Zhang et al.
Reference Zhang, Decyk, Holod, Xiao, Lin and Chen2010; Pueschel et al.
Reference Pueschel, Jenko, Schneller, Hauff, Günter and Tardini2012; Wilkie et al.
Reference Wilkie, Abel, Highcock and Dorland2015), making this condition marginally or well-satisfied depending on the scenario. The presence of a local maximum in velocity space has the potential to destabilize waves, such as Alfvén eigenmodes, via inverse Landau damping. A flattening or inversion in velocity space is destabilizing, but it accompanied by a flattening of the radial profile, which is stabilizing. Without a detailed analysis, it is not immediately clear what the net effect is.
Also, notice that, at the birth speed
$v=v_{\unicode[STIX]{x1D6FC}}$
, there is no turbulent correction to the classical slowing-down distribution. No matter how strong the turbulence is, it will not affect newly born alpha particles because it cannot possibly compete with a delta function at its peak. As the alpha particle slows down, it becomes doubly affected because it is further away in velocity space from the source in addition to turbulent transport becoming stronger. Since most of the electron heating occurs near the injection speed, this remains immune to turbulent transport of alphas. What might be affected is the heating of ions at lower energy.
The model distribution was benchmarked against a numerical solution of (1.1). The test case is an ITER-like ELMy H-mode DT scenario that is unstable to electrostatic ion-temperature gradient turbulence beyond
$r=0.5a$
. It can be found on the public tokamak profile database (Roach et al.
Reference Roach, Walters, Budny, Imbeaux, Fredian, Greenwald, Stillerman, Alexander, Carlsson and Cary2008) as ITER scenario 10010100. The diffusion coefficient was calculated directly from first-principles nonlinear gyrokinetic simulations (and thus did not diverge as
$v\rightarrow 0$
). Further, the full test-particle collision operator was used to capture physically consistent thermalization of helium ash. The numerical solution of the global phase space transport (1.1) was found with the t3core solver. For more details about the transport of alpha particles in this case, see Wilkie et al. (Reference Wilkie, Abel, Landreman and Dorland2016).
The local distribution functions at two different radii are shown in figure 2. Despite the generous assumptions built into the analytic model, it qualitatively captures the basic features of the numerical solution including: well-confined high-energy alphas, a depletion of the distribution at more moderate energies, and the bump-on-tail feature. At these radii (
$r=0.5a$
and
$0.6a$
), the values of
$b$
were found to be 5.3 and 2.6, respectively.
Comparing the local velocity distribution function of alpha particles for two different radii of the ITER-like scenario of (Wilkie et al. Reference Wilkie, Abel, Landreman and Dorland2016). Shown are the Gaffey slowing-down distribution (cyan dashed line), equation (2.9) (red line), and a numerical phase-space transport solution (black line). Each of the analytic distributions are shown with a population of Maxwellian ash added, equal to the ash density calculated with t3core.

The alpha particle density as a function of radius for the test case is shown in figure 3. There, the analytic solution is compared to the numerical transport solution and to the slowing-down distribution without radial transport. Even when the classical density differs from the numerical one due to transport, equation (2.9) captures this difference fairly well. This motivates using the model as a heuristic to answer the question: ‘can one expect significant energetic particle transport from microturbulence?’ An affirmative answer indicates that further study is warranted, but if
$b<1$
, the classical slowing-down distribution is reasonably accurate.
Comparing the radial density profile of energetic alpha particles as predicted by the Gaffey slowing-down distribution (cyan dashed line), equation (2.9) (red line), and a numerical phase-space transport solution (black line).

4 Conclusion
This work introduced and briefly discussed a simplified model for slowing-down alpha particles in the presence of microturbulence. Although heuristic, this model distribution captures the approximate effect of turbulence through a single dimensionless parameter, which has a stronger impact than one might expect when the birth speed
$v_{\unicode[STIX]{x1D6FC}}$
significantly exceeds the critical speed
$v_{c}$
. When relaxing several assumptions made for this simple model, more robust analytic distributions are found (equations (B 5), (C 5) and (C 8)) which maintain the same basic character of (2.9). Although more detailed simulations are generally required to robustly predict the impact of microturbulence on plasma heating, Alfvén eigenmode drive and other effects of by energetic particles, this model provides physical insight into the expected impact of microturbulence on these phenomena.
Acknowledgements
The author would like to thank I. Pusztai, I. Abel and T. Fülöp for helpful comments in preparing this manuscript. The author was supported by the EUROfusion Researcher Grant AWP18-ERG-VR. The views and opinions expressed herein do not necessarily reflect those of the European Commission.
Appendix A. Energy scaling of the diffusion coefficient
In this section, the fluctuating distribution of energetic particles is analysed to motivate the energy scaling of the diffusion coefficient used in this work (Hauff et al.
Reference Hauff, Pueschel, Dannert and Jenko2009). Only electrostatic fluctuations will be considered here. After finding a solution to the gyrokinetic equation in the high-energy limit as was done in Wilkie et al. (Reference Wilkie, Iantchenko, Abel, Highcock and Pusztai2018), the radial flux is averaged over the spatial and temporal scale of the fluctuations. The diffusion coefficient is what remains after factoring out
$\unicode[STIX]{x2202}F_{0\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}r$
. In the meantime, it will become apparent why radial diffusion is dominant over the other turbulent terms of (1.1) in this limit.
The collisionless gyrokinetic equation for non-Maxwellian fast particles is (Frieman & Chen Reference Frieman and Chen1982; Wilkie Reference Wilkie2015; Abel, Wilkie & Schekochihin Reference Abel, Wilkie and Schekochihin2018):
where it is expressed in the coordinates
${\mathcal{E}}=Z_{\unicode[STIX]{x1D6FC}}e\unicode[STIX]{x1D719}+m_{\unicode[STIX]{x1D6FC}}v^{2}/2$
, and
$\unicode[STIX]{x1D707}$
, which is the exactly conserved magnetic moment. With this definition,
$\unicode[STIX]{x2202}F_{0\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}\unicode[STIX]{x1D707}$
does not appear in (A 1). It does appear, along with the adiabatic contribution, in the equations that determine the fluctuating fields, which will not be used here. The distribution is decomposed into an equilibrium part
$F_{0\unicode[STIX]{x1D6FC}}$
and a fluctuating partFootnote
2
$h_{\unicode[STIX]{x1D6FC}}$
such that
$f_{\unicode[STIX]{x1D6FC}}=F_{0\unicode[STIX]{x1D6FC}}+h_{\unicode[STIX]{x1D6FC}}$
. The magnetic field is given in the Clebsch representation as
$\boldsymbol{B}=B\boldsymbol{b}=\unicode[STIX]{x1D735}\unicode[STIX]{x1D6FC}\times \unicode[STIX]{x1D735}\unicode[STIX]{x1D713}$
, where
$\unicode[STIX]{x1D713}$
is the poloidal flux in the flux surface that it labels and
$\unicode[STIX]{x1D6FC}$
(not to be confused with the subscript denoting energetic particles) labels the magnetic field line within said flux surface. The velocity parallel to the magnetic field is
$v_{\Vert }=\unicode[STIX]{x1D70E}_{\Vert }\sqrt{(2/m_{\unicode[STIX]{x1D6FC}})({\mathcal{E}}-\unicode[STIX]{x1D707}B(\unicode[STIX]{x1D703}))}$
, where
$\unicode[STIX]{x1D70E}=\pm 1$
, and
$\unicode[STIX]{x1D703}$
is the ballooning-extended poloidal angle, used as a coordinate along the field line. The fluctuating static potential is
$\unicode[STIX]{x1D719}$
. The gyroaverage at constant gyrocentre
$\boldsymbol{R}$
is denoted
$\langle \rangle _{\boldsymbol{R}}$
, and the spatial gradients in (A 1) are taken with respect to
$\boldsymbol{R}$
. The drift associated with the fluctuating fields is given by
$\boldsymbol{v}_{\unicode[STIX]{x1D719}}\equiv (c/B)\boldsymbol{b}\times \unicode[STIX]{x1D735}\langle \unicode[STIX]{x1D719}\rangle _{\boldsymbol{R}}$
, and
$\boldsymbol{v}_{D}$
is the magnetic drift due to the curvature and gradient in
$\boldsymbol{B}$
.
Consider that, in ion-scale turbulence, the bulk ions set the spatial and temporal scales for the variation of
$\unicode[STIX]{x1D719}$
and
$h_{\unicode[STIX]{x1D6FC}}$
. In the limit where
$v\gg v_{ti}$
, the magnetic drift, parallel streaming and equilibrium gradient are the largest terms in (A 1). The ‘nonlinear’
$\boldsymbol{v}_{\unicode[STIX]{x1D719}}$
term is thereby small. It is convenient and conventional to express the fluctuating quantities in an eikonal representation so that
$\unicode[STIX]{x1D719}=\hat{\unicode[STIX]{x1D719}}(\unicode[STIX]{x1D703})\exp [\text{i}S(\unicode[STIX]{x1D713},\unicode[STIX]{x1D6FC})]$
, where
$S$
contains the small-scale variation perpendicular to the magnetic field and
$\hat{\unicode[STIX]{x1D719}}$
is the relatively long-scale parallel variation (Antonsen & Lane Reference Antonsen and Lane1980). In this representation, the energetic limit of (A 1) reduces to:
Even near the source, where
$\unicode[STIX]{x2202}F_{0\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}{\mathcal{E}}\sim F_{0\unicode[STIX]{x1D6FC}}/T_{i}$
, the first term on the right-hand side of (A 1) (which ultimately contributes to
$D_{rv}$
) is smaller than the first term (which leads to
$D_{rr}$
) on the right-hand side of (A 2) by virtue of the large spatial gradients in the source of energetic particles (and therefore in
$F_{0\unicode[STIX]{x1D6FC}}$
itself). In any case, at energies this high, collisions acting on the equilibrium dominate both
$D_{rr}$
and
$D_{rv}$
. Away from the source, where turbulence is observed to have an effect,
$\unicode[STIX]{x2202}F_{0\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}{\mathcal{E}}\sim F_{0\unicode[STIX]{x1D6FC}}/E_{\unicode[STIX]{x1D6FC}}$
, and this term is smaller still. Similar arguments apply to the turbulent heating terms (
$D_{vr}$
and
$D_{vv}$
), only more so. Furthermore, if the
$\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D719}\rangle _{\boldsymbol{R}}/\unicode[STIX]{x2202}t$
term were kept, its contribution to the time-averaged radial flux vanishes in the energetic (thereby, linear) limit. It is for these reasons that radial diffusion is dominant over other terms in (1.1).
In finding an integrating factor to solve (A 2), it will be convenient to define:
where
$\unicode[STIX]{x1D714}_{D}\equiv \boldsymbol{v}_{D}\boldsymbol{\cdot }\unicode[STIX]{x1D735}S$
and
$\unicode[STIX]{x1D714}_{D}(\unicode[STIX]{x1D703}_{0})=0$
. Primes denote quantities taken at
$\unicode[STIX]{x1D703}^{\prime }$
rather than
$\unicode[STIX]{x1D703}$
. A saddle point is located at
$\unicode[STIX]{x1D703}_{0}$
such that
$[\text{d}z/\text{d}\unicode[STIX]{x1D703}]_{\unicode[STIX]{x1D703}_{0}}=0$
. With this definition, the solution to (A 2) is (Kim, Horton & Dong Reference Kim, Horton and Dong1993; Wilkie et al.
Reference Wilkie, Iantchenko, Abel, Highcock and Pusztai2018):
where
$\text{J}_{0}$
is the zeroth-order Bessel function of the first kind and
$\unicode[STIX]{x1D6FA}_{\unicode[STIX]{x1D6FC}}\equiv Z_{\unicode[STIX]{x1D6FC}}eB/m_{\unicode[STIX]{x1D6FC}}c$
is the gyrofrequency. The lower limit of integration in (A 4) is either
$\pm \infty$
for passing particles depending on the sign of
$v_{\Vert }$
, or
$\pm \unicode[STIX]{x1D703}_{b}$
if the particles are trapped. Expanding
$z$
in a Taylor series about
$\unicode[STIX]{x1D703}_{0}$
leaves a Gaussian which is integrated over the contour of steepest descent. Take advantage of the fact that
$z\gg 1$
to obtain:
where
$[\,]_{0}$
denotes quantities taken at
$\unicode[STIX]{x1D703}_{0}$
, and
$\unicode[STIX]{x1D701}_{0}$
is defined as:
The approximate solution of (A 5) excludes deeply trapped particles; namely those that bounce at a poloidal location
$|\unicode[STIX]{x1D703}_{b}|\lesssim |\unicode[STIX]{x1D703}_{0}|$
.
The radial flux is kept as a function of velocity in the kinetic phase-space transport equation, and is defined as:
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D6FC}}\equiv \left\langle \mathop{\sum }_{\unicode[STIX]{x1D70E}_{\Vert }}\int _{-\infty }^{\infty }h_{\unicode[STIX]{x1D6FC}}\boldsymbol{v}_{\unicode[STIX]{x1D719}}\boldsymbol{\cdot }\unicode[STIX]{x1D735}r\,\text{d}\unicode[STIX]{x1D703}\right\rangle _{\bot ,t}, & & \displaystyle\end{eqnarray}$$
where
$\langle \,\rangle _{\bot }$
is a perpendicular spatial average over the flux tube and a time average over the fluctuation time scale (Abel et al.
Reference Abel, Plunk, Wang, Barnes, Cowley, Dorland and Schekochihin2013). In performing the
$\unicode[STIX]{x1D703}$
integral, apply the same large-
$z$
approximation that led to (A 5). Also note that
$h_{\unicode[STIX]{x1D6FC}}$
is proportional to
$\unicode[STIX]{x2202}F_{0\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}\unicode[STIX]{x1D713}$
, allowing one to write
$\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D6FC}}=-D_{rr}\unicode[STIX]{x2202}F_{0\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}r$
for
$\unicode[STIX]{x1D713}=\unicode[STIX]{x1D713}(r)$
.Footnote
3
This gives:
$$\begin{eqnarray}\displaystyle D_{rr}\approx \left\langle \!2\unicode[STIX]{x03C0}c^{2\!}\left(\frac{\text{d}r}{\text{d}\unicode[STIX]{x1D713}}\right)^{\!2}\!\left(\frac{\unicode[STIX]{x2202}S}{\unicode[STIX]{x2202}\unicode[STIX]{x1D6FC}}\right)^{\!2}\text{J}_{0}^{2}\left(\frac{|\unicode[STIX]{x1D735}S|v_{\bot }}{\unicode[STIX]{x1D6FA}_{\unicode[STIX]{x1D6FC}}(\unicode[STIX]{x1D703}_{0})}\right)\!|\unicode[STIX]{x1D719}(\unicode[STIX]{x1D703}_{0})|^{2}\displaystyle \left|(v_{\Vert }\boldsymbol{b}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703})_{\unicode[STIX]{x1D703}_{0}}\left.\!\frac{\text{d}}{\text{d}\unicode[STIX]{x1D703}}\right|_{\unicode[STIX]{x1D703}_{0}}\frac{\unicode[STIX]{x1D714}_{D}}{v_{\Vert }\boldsymbol{b}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}}\right|^{-1}\right\rangle _{\!\bot ,t}\!. & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$
The origin of the
$v^{-3}$
scaling of Hauff et al. (Reference Hauff, Pueschel, Dannert and Jenko2009) is now apparent by counting powers of
$v$
, noting that
$\text{J}_{0}(x\gg 1)\sim x^{-1/2}$
and
$\unicode[STIX]{x1D714}_{D}\propto v^{2}$
. It was with this scaling that the isotropic slowing-down distribution function was found in (2.9).
Appendix B. Pitch-angle dependence and anisotropic distribution
In this section, the pitch-angle dependence of the approximate diffusion coefficient in (A 8) is kept to obtain a slowing-down distribution that includes pitch-angle-dependent transport. Thus an anisotropic equilibrium distribution for weakly collisional energetic particles is derived. Doing so will require making further assumptions about the turbulent spectrum that should have minimal impact at high energy. By using the value of the diffusion coefficient at
$\unicode[STIX]{x1D709}=0$
, it will be found that (2.9) remains an adequate approximation to the anisotropic distribution up to
$|\unicode[STIX]{x1D709}|\approx 0.9$
.
To facilitate calculating the perpendicular spatial average in (A 8), make the ansatz that the turbulent spectrum is on the ion Larmor radius scale and has the following form:
where
$k_{y}\equiv a^{-1}\unicode[STIX]{x2202}S/\unicode[STIX]{x2202}\unicode[STIX]{x1D6FC}$
and
$k_{\text{max}}\sim \unicode[STIX]{x1D70C}_{i}^{-1}$
is the wavenumber at the peak of the turbulent spectrum. The exact form of the spectrum does not significantly impact the functional dependence of the result at high energy (especially when
$D_{\unicode[STIX]{x1D6FC}}$
is assumed given as it is here), but provides a meaningful cutoff that will prove necessary. Integrating over all
$k_{y}$
smooths the oscillations in the Bessel function and note that
$\unicode[STIX]{x1D714}_{D}\propto k_{y}(v_{\Vert }^{2}+v_{\bot }^{2}/2)$
. Equation (A 8) is thereby well approximated by:
$$\begin{eqnarray}\displaystyle D_{rr}\approx D_{\unicode[STIX]{x1D6FC}}\frac{v_{\unicode[STIX]{x1D6FC}}^{3}}{v^{3}}\frac{1}{1+\unicode[STIX]{x1D709}^{2}}\frac{1}{\displaystyle \sqrt{1-\unicode[STIX]{x1D709}^{2}}+\frac{v_{ti}}{v}\frac{1}{2\sqrt{\unicode[STIX]{x03C0}}k_{\text{max}}\unicode[STIX]{x1D70C}_{i}}}, & & \displaystyle\end{eqnarray}$$
where
$D_{\unicode[STIX]{x1D6FC}}$
is the diffusion coefficient at
$\unicode[STIX]{x1D709}=0$
and
$v=v_{\unicode[STIX]{x1D6FC}}\gg v_{ti}$
. The rightmost factor comes from the Padé approximation
$\int _{0}^{\infty }k^{3}\text{e}^{-k^{2}}\text{J}_{0}^{2}(ak)\,\text{d}k\approx (4\sqrt{\unicode[STIX]{x03C0}}a+2)^{-1}$
. This approximate form contains all the velocity dependence of the diffusion coefficient, and is shown in figure 4. Equation (B 2) recovers the energy dependence of both electrostatic limits considered in Hauff et al. (Reference Hauff, Pueschel, Dannert and Jenko2009), but avoids singularities in pitch angle. The pitch-angle dependence in (B 2) is approximately constant except for a narrow region near
$\unicode[STIX]{x1D709}=\pm 1$
. There, the diffusion coefficient becomes large and the energy dependence changes to
$v^{-2}$
(Pueschel et al.
Reference Pueschel, Jenko, Schneller, Hauff, Günter and Tardini2012). This represents the part of phase space where the energetic particle’s Larmor radius becomes comparable to
$\unicode[STIX]{x1D70C}_{i}$
, despite its large energy, by virtue of its relatively small perpendicular velocity. In this regime, the details of the turbulence do matter and the ansatz made in (B 1) becomes unreliable. In any case, one can nevertheless expect a significant increase in the diffusion coefficient when
$v\approx |v_{\Vert }|$
.
Consider the low-collisionality transport equation with pitch-angle scattering included as well as the
$\unicode[STIX]{x1D709}$
dependence of transport given by (B 2):
$$\begin{eqnarray}\displaystyle & & \displaystyle b\frac{v_{c}^{3}}{v^{3}}\left[(1+\unicode[STIX]{x1D709}^{2})\left(\sqrt{1-\unicode[STIX]{x1D709}^{2}}+\frac{v_{ti}}{2v\sqrt{\unicode[STIX]{x03C0}}k_{\text{max}}\unicode[STIX]{x1D70C}_{i}}\right)\right]^{-1}f_{\unicode[STIX]{x1D6FC}}\nonumber\\ \displaystyle & & \displaystyle \quad =\frac{1}{2}\frac{v_{ti}^{3}}{v^{3}}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D709}}(1-\unicode[STIX]{x1D709}^{2})\frac{\unicode[STIX]{x2202}f_{\unicode[STIX]{x1D6FC}}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D709}}+\frac{1}{v^{2}}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}v}[(v_{c}^{3}+v^{3})f_{\unicode[STIX]{x1D6FC}}]+\frac{S_{0}\unicode[STIX]{x1D70F}_{s}}{4\unicode[STIX]{x03C0}v_{\unicode[STIX]{x1D6FC}}^{2}}\unicode[STIX]{x1D6FF}(v-v_{\unicode[STIX]{x1D6FC}}),\end{eqnarray}$$
where the
$v\gg v_{ti}$
limit was used for the deflection collision frequency on the Lorentz operator. Note that, unless
$\unicode[STIX]{x2202}\ln f_{\unicode[STIX]{x1D6FC}}/\unicode[STIX]{x2202}\unicode[STIX]{x1D709}\gtrsim O(v_{c}^{3}/v_{ti}^{3})$
, the pitch-angle scattering term is small. In this case, a complicated, but analytic
$f_{\unicode[STIX]{x1D6FC}}$
can be found with the help of the integrating factor
$(v_{c}^{3}+v^{3})\exp (g)$
, where
$$\begin{eqnarray}\displaystyle g(v,\unicode[STIX]{x1D709}) & = & \displaystyle \frac{-b}{(1+\unicode[STIX]{x1D709}^{2})}\frac{\unicode[STIX]{x03C0}}{24\unicode[STIX]{x03C0}^{2}v_{c}^{3}(1-\unicode[STIX]{x1D709})^{3/2}-3\sqrt{\unicode[STIX]{x03C0}}v_{ti}^{3}}\nonumber\\ \displaystyle & & \displaystyle \times \left\{2\sqrt{3}v_{c}v_{ti}(2\sqrt{\unicode[STIX]{x03C0}}v_{c}\sqrt{1-\unicode[STIX]{x1D709}^{2}}-v_{ti})\tan ^{-1}\left(\frac{2v-v_{c}}{\sqrt{3}v_{c}}\right)\right.\nonumber\\ \displaystyle & & \displaystyle +\,v_{c}v_{ti}(2\sqrt{\unicode[STIX]{x03C0}}v_{c}\sqrt{1-\unicode[STIX]{x1D709}^{2}}+v_{ti})\ln \left[\frac{v^{2}-vv_{c}+v_{c}^{2}}{(v_{c}+v)^{2}}\right]\nonumber\\ \displaystyle & & \displaystyle \left.-\,8\unicode[STIX]{x03C0}v_{c}^{3}(1-\unicode[STIX]{x1D709}^{2})\ln \left[\frac{v_{c}^{3}+v^{3}}{(2\sqrt{\unicode[STIX]{x03C0}}v\sqrt{1-\unicode[STIX]{x1D709}^{2}}+v_{ti})^{3}}\right]\right\}.\end{eqnarray}$$
The approximate distribution function with the full pitch-angle dependence of (B 2) is therefore:
This distribution is shown for select parameters in figure 5, and the pitch-angle dependence is shown in figure 6. With this form, the assumption of negligible pitch-angle scattering breaks down in the region
$|\unicode[STIX]{x1D709}|\gtrsim 0.99$
, where pitch-angle derivatives become sufficiently large compared to the distribution. Also in this region, the turbulent spectrum used to derive (B 5) becomes unreliable because the transport will depend on the details of the turbulence, which can only be found from direct simulation.
Pitch-angle dependence of the generalized slowing-down distribution, equation (B 5), for several values of turbulent intensity (parametrized by
$b$
). Here,
$v_{c}=0.3v_{\unicode[STIX]{x1D6FC}}$
and
$v_{ti}=0.15k_{\text{max}}\unicode[STIX]{x1D70C}_{i}v_{\unicode[STIX]{x1D6FC}}$
.

Pitch-angle dependence of the generalized slowing-down distribution, equation (B 5), at
$v=v_{c}=0.3v_{\unicode[STIX]{x1D6FC}}$
and
$v_{ti}=0.15k_{\text{max}}\unicode[STIX]{x1D70C}_{i}v_{\unicode[STIX]{x1D6FC}}$
. Dotted lines represent the value of the isotropic distribution, equation (2.9) using the diffusion coefficient at
$\unicode[STIX]{x1D709}=0$
.

Appendix C. Global model
In this section, an analytic distribution with radial dependence is derived similarly, but more generally, than that of Anderson et al. (Reference Anderson, Batistoni and Anderson1991), which will aid in interpreting the length scale
$L_{\unicode[STIX]{x1D6FC}}$
, and thereby the model distribution itself, more rigorously. This is possible for a fictitious radial profile in the large aspect ratio limit where
$D_{rr}$
,
$\unicode[STIX]{x1D70F}_{s}$
, and
$v_{c}$
are independent of radius, while the energetic particle source (which is typically more strongly peaked than the other parameters) retains full radial dependence. This establishes that the velocity dependence of the model maintains its general form even when consistently accounting for radial transport. Equation (C 8) provides a turbulence-modified analytic distribution with radial, energy and pitch-angle dependence.
Expand the radial dependence of the source in a Fourier–Bessel series:
where
$\unicode[STIX]{x1D70E}_{1}$
is the velocity-integrated source at the magnetic axis and
$a_{0,j}$
is the
$j$
th positive root of the zero-order Bessel function
$\text{J}_{0}$
. This diagonalizes the radial diffusion operator in the circular large aspect ratio limit where
$V^{\prime }(r)\propto r$
. Expand
$f_{\unicode[STIX]{x1D6FC}}(r,v)=\sum f_{\unicode[STIX]{x1D6FC},j}(v)\text{J}_{0}(a_{0,j}r/a)$
and use the same approximations of § 2 for radial diffusion, the collision operator and source to obtain:
After multiplying by
$r\text{J}_{0}(a_{0,k}r/a)$
and integrating from
$r=0$
to
$a$
, the approximate equation for each Bessel mode
$k$
is:
For each
$j$
, this can be solved the same way as (2.7) to obtain:
where:
and
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D70E}_{j}\equiv \frac{\displaystyle 2\int _{0}^{a}S_{0}(r)\text{J}_{0}(a_{0,j}r/a)r\,\text{d}r}{a^{2}\text{J}_{1}^{2}(a_{0,j})}. & & \displaystyle\end{eqnarray}$$
This is a radially global distribution approximately valid where the plasma profile (particularly,
$D_{rr}$
,
$\unicode[STIX]{x1D70F}_{s}$
and
$v_{c}$
) is uniform, but the energetic particle source is not. It has a similar velocity dependence to the local distribution, equation (2.9). If only the first term is used to approximate the source, then
$S_{\unicode[STIX]{x1D6FC}}\propto \text{J}_{0}(a_{0,1}r/a)$
, and the correspondence is direct with a more rigorous definition of
$L_{\unicode[STIX]{x1D6FC}}=a/a_{0,1}$
. Otherwise, the global distribution will be a sum of distributions like
$F_{SD,\text{mod}}$
, with each term having its own value of
$b$
.
If the pitch-angle and radial dependence are retained simultaneously and one further assumes that
$v_{ti}$
is uniform in
$r$
, an analytic distribution is obtainable in three-dimensional phase space:
where
$g_{j}$
is given by (B 4), except
$b$
replaced by
$b_{j}$
from (C 6).






