2.1. Distribution function

Velocity space can be separated into the trapped–barely passing region and the freely passing region. We denote the distribution function in the freely passing region as $g^p$ and in the trapped–barely passing region as $g^{t,bp}$ . It will be shown later that the distribution functions in these two regions are governed by different equations. The integration over the respective regions is defined as

(2.21) \begin{equation} \int _{V_{tbp}}\mathrm{d}^3v\: \equiv \lim _{W\rightarrow \infty }\int _0^\infty \mathrm{d}\mu \: \int _{-W}^W\mathrm{d}w\:2\pi B \end{equation}

in the trapped–barely passing region and

(2.22) \begin{equation} \int _{V_{p}}\mathrm{d}^3v\: \equiv \lim _{\delta \rightarrow 0}\int _0^\infty \mathrm{d}\mu \: \int _{-\infty }^{-u-\delta }\mathrm{d}v_\parallel \: 2\pi B+\lim _{\delta \rightarrow 0}\int _0^\infty \mathrm{d}\mu \: \int _{-u+\delta }^{\infty }\mathrm{d}v_\parallel \: 2\pi B, \end{equation}

where $\delta \gt 0$ , in the freely passing region.

The distribution functions have to match asymptotically, i.e.

(2.23) \begin{equation} g^p(v_\parallel \rightarrow -u^+)=g^{t,bp}(w\rightarrow \infty )\quad \text{and}\quad g^p(v_\parallel \rightarrow -u^-)=g^{t,bp}(w\rightarrow -\infty ). \end{equation}

Here, $-u^+$ ( $-u^-$ ) indicates that $-u$ is approached from the positive (negative) side. The distribution function $g^p$ is discontinuous across the trapped–barely passing region. This jump $\Delta g^p$ matches the jump in the trapped–barely passing distribution function

(2.24) \begin{equation} \Delta g^p\equiv g^p(v_\parallel \rightarrow -u^+)-g^p(v_\parallel \rightarrow -u^-)=g^{t,bp}(w\rightarrow \infty )-g^{t,bp}(w\rightarrow -\infty ). \end{equation}

With these definitions, we can also introduce the linearised collision operator in the freely passing region

(2.25) \begin{align} C_p^{(l)}[g]= \boldsymbol{\nabla} _v\boldsymbol{\cdot }\Bigg [f_M {\unicode{x1D648}} \boldsymbol{\cdot }\boldsymbol{\nabla} _v\left (\frac {g^p}{f_M}\right )-\lambda f_M\int _{V_{tbp}}\mathrm{d}^3v' f'_M \boldsymbol{\nabla} _\omega \boldsymbol{\nabla} _\omega \omega \boldsymbol{\cdot }\boldsymbol{\nabla} _v'\left (\frac {g^{t,bp'}}{f'_M}\right )\nonumber \\ -\lambda f_M\int _{V_{p}}\mathrm{d}^3v' f'_M \boldsymbol{\nabla} _\omega \boldsymbol{\nabla} _\omega \omega \boldsymbol{\cdot }\boldsymbol{\nabla} _v'\left (\frac {g^{p'}}{f'_M}\right )\Bigg ] \equiv \boldsymbol{\nabla} _v\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p, \end{align}

and the linearised collision operator in the trapped–barely passing region

(2.26) \begin{align} C_{t,bp}^{(l)}[g]= \boldsymbol{\nabla} _v\boldsymbol{\cdot }\Bigg [f_M {\unicode{x1D648}} \boldsymbol{\cdot }\boldsymbol{\nabla} _v\left (\frac {g^{t,bp}}{f_M}\right )-\lambda f_M\int _{V_{tbp}}\mathrm{d}^3v' f'_M \boldsymbol{\nabla} _\omega \boldsymbol{\nabla} _\omega \omega \boldsymbol{\cdot }\boldsymbol{\nabla} _v'\left (\frac {g^{t,bp'}}{f'_M}\right )\nonumber \\ -\lambda f_M\int _{V_{p}}\mathrm{d}^3v' f'_M \boldsymbol{\nabla} _\omega \boldsymbol{\nabla} _\omega \omega \boldsymbol{\cdot }\boldsymbol{\nabla} _v'\left (\frac {g^{p'}}{f'_M}\right )\Bigg ] \equiv \boldsymbol{\nabla} _v\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp}. \end{align}

We expand the distribution function $g$ in the size of the trapped region $\epsilon ^{1/2}$ , i.e.

(2.27) \begin{align} g=g_0+g_1+\ldots , \end{align}

where $g_1\sim \epsilon ^{1/2}g_0$ and $g_0\sim \epsilon ^{1/2}f_M$ . This expansion is valid in both the trapped–barely passing and the freely passing regions.

In the freely passing region, $g^p$ is independent of $\theta$ to lowest order in $\sqrt {\epsilon }$ , as shown by Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023). The $\theta$ -dependent part was calculated by Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023) in their equation (4.53) and yields

(2.28) \begin{align} & g^p-\langle g^p\rangle _\psi = -\frac {I}{\varOmega }\frac {r}{R}\frac {\big(v_\parallel ^2+\mu B\big) \cos {\theta }-ZeR\phi _\theta /mr}{v_\parallel +u}\Bigg [\frac {\partial }{\partial \psi }\ln {p} +\frac {m(v_\parallel - V_\parallel )}{T}\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}\right )\nonumber\\& \!\!\quad + \left (\frac {m(v_\parallel -V_\parallel )^2}{2T}+\frac {m\mu B}{T}-\frac {5}{2}\right )\frac {\partial }{\partial \psi }\ln {T}\Bigg ]f_M -\frac {r}{R}\cos {\theta }\frac {m}{T}\left [v_\parallel (v_\parallel -V_\parallel )+\mu B \right ]\!f_M\sim \epsilon f_M. \end{align}

Here, we introduced the flux surface average

(2.29) \begin{equation} \langle \ldots \rangle _\psi =\int \frac {\mathrm{d}\theta }{2\pi }(\ldots ). \end{equation}

The expression (2.28) can also be found by noting that the terms proportional to $\partial g/\partial v_\parallel$ and $\partial g/\partial \psi$ in (2.9) are negligible in the freely passing region, where $\partial /\partial v_\parallel \sim 1/v_t$ . The collision operator is negligible to lowest order for collisionalities $\nu _\ast \ll \epsilon ^{-3/2}$ . The $\theta$ dependence in the trapped–barely passing region has to match the freely passing solution according to (2.23). Thus, any $\theta$ dependence in $g_0^{t,bp}$ has to decay as $w\rightarrow \pm \infty$ and cannot contribute to the jump $\Delta g^p$ .

In the trapped–barely passing region, we can use the fact that $w\simeq 0$ to lowest order. The linearised collision operator (2.26) to lowest order in $\epsilon$ is

(2.30) \begin{equation} C_{t,bp}^{(l)}[g]=\mathsf{M}_\parallel \frac {\partial ^2 {g^{t,bp}_0}}{\partial {w}^{2}}, \end{equation}

where we defined

(2.31) \begin{equation} \mathsf{M}_\parallel \equiv \boldsymbol{\hat {b}}\boldsymbol{\cdot }{\unicode{x1D648}}\boldsymbol{\cdot }\boldsymbol{\hat {b}}\Big \rvert _{w=0}= \frac {\nu _\perp }{2}\mu B+\frac {\nu _\parallel }{2}(u+V_\parallel )^2. \end{equation}

Here, $\nu _\perp$ and $\nu _\parallel$ are evaluated at $w=0$ . We find that (2.18) is, to lowest order in $\sqrt {\epsilon }$ and for $\nu _\ast \sim 1$ ,

(2.32) \begin{align} &\frac {\partial }{\partial \theta }\left (\frac {w}{qR}g_0^{t,bp}\right )-S\frac {\partial }{\partial w}\left [\mathcal{P}(\theta )\frac {r}{R}\frac {1}{qR}g_0^{t,bp}\right ]\nonumber \\ &\qquad\qquad\qquad\qquad\qquad\qquad =\frac {I}{\varOmega }\mathcal{P}(\theta )\frac {r}{R}\frac {1}{qR}\mathcal{D}f_M(w=0)+\mathsf{M}_\parallel \frac {\partial ^2{g^{t,bp}_0}}{\partial {w}^2}, \end{align}

where we introduced the poloidal dependence

(2.33) \begin{equation} \mathcal{P}(\theta )=(u^2+\mu B)\sin \theta +\frac {ZeR}{mr}\frac {\partial \phi _\theta }{\partial \theta } \end{equation}

and the derivative term

(2.34) \begin{equation} \mathcal{D}\equiv \frac {\partial }{\partial \psi }\ln {p} -\frac {m(u+ V_\parallel )}{T}\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}\right ) +\left (\frac {m(u+V_\parallel )^2}{2T}+\frac {m\mu B}{T}-\frac {5}{2}\right )\frac {\partial }{\partial \psi }\ln {T}. \end{equation}

In the large $w$ limit, (2.32) reduces to

(2.35) \begin{equation} \frac {\partial }{\partial \theta }\left (\frac {w}{qR}g_0^{t,bp}\right ) =\frac {I}{\varOmega }\mathcal{P}(\theta )\frac {r}{R}\frac {1}{qR}\mathcal{D}f_M(w=0). \end{equation}

From this, we find that the $\theta$ -dependent piece of $g^{t,bp}_0$ decays as $\sim 1/w$ for $w\rightarrow \infty$ as expected. Furthermore, we recover the terms in expression (2.28) that diverge as $1/w$ for small $w$ . The other terms in (2.28) match the large $w$ limit of the next-order correction $g_1^{t,bp}$ . Equation (2.32) is sufficient to determine neoclassical transport quantities and we will not solve for $g_1^{t,bp}$ in this paper.

2.2. Moment equations

We are interested in calculating particle, parallel momentum and energy transport in the limit of $\epsilon \ll 1$ and $\nu _\ast \sim 1$ . The transport relations are associated with moments of the drift kinetic equation. The integration has to treat freely passing and trapped–barely passing particles separately. First, we discuss the freely passing region.

The transit average of (2.9) appears to have many contributions by terms proportional to $f_M$ and terms proportional to derivatives of $g$ . However, the $\theta$ dependence of $f_M$ and $g^p$ is of order O $(\epsilon )$ ; see (2.28). Thus, the transit average of (2.9) in the freely passing region to lowest order in $\epsilon$ reduces to

(2.36) \begin{equation} \big\langle C_p^{(l)}[g]\big\rangle _\tau \simeq -\langle \varSigma \rangle _\tau . \end{equation}

We define the transit average for passing particles as

(2.37) \begin{align} \langle \ldots \rangle _\tau =\frac {1}{\tau } \int _0^{2\pi }\mathrm{d}\theta \:\frac { qR}{v_\parallel +u}(\ldots ), \end{align}

where

(2.38) \begin{align} \tau =\int _0^{2\pi }\mathrm{d}\theta \:\frac {qR}{v_\parallel +u}. \end{align}

Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023) introduced a new set of variables that is based on conserved quantities, which proved (2.36) to be exact if the new variables are held fixed in the transit average; see their equation (4.21). The difference between the new variables and $\lbrace v_\parallel ,\psi \rbrace$ is small in $\epsilon$ for freely passing particles and negligible for the purpose of this paper. We derive transport equations of particles, momentum and energy by taking moments of the drift kinetic equation (2.36).

The integration of (2.36) over the freely passing region to find the particle transport is

(2.39) \begin{equation} \int _{V_p}\mathrm{d}^3v\: \big\langle C_p^{(l)}[g]\big\rangle _\tau =-\int _{V_p}\mathrm{d}^3v\:\langle \varSigma \rangle _\tau . \end{equation}

The integration over the freely passing region is defined as

(2.40) \begin{equation} \int _{V_p}\mathrm{d}^3v\:(\ldots )=\int _{-\infty }^{-u^-}\mathrm{d}v_\parallel \:\int _0^\infty \mathrm{d}\mu \:2\pi B(\ldots )+\int _{-u^+}^{\infty }\mathrm{d}v_\parallel \:\int _0^\infty \mathrm{d}\mu \:2\pi B(\ldots ), \end{equation}

where the limits $-u^-$ and $-u^+$ indicate that $g^p$ is not continuous at $-u$ and we need to distinguish between the two sides of $-u$ .

Because the change in parallel velocity along the orbit is small for freely passing particles, it is possible to replace transit averages by flux surface averages, $\langle \ldots \rangle _\tau \simeq \langle \ldots \rangle _\psi$ . The integration over $\mu$ in (2.39) cancels out the $\mu$ derivative in the divergence of the collision operator in (2.13) and (2.39) becomes

(2.41) \begin{equation} \bigg \langle\! \int _{V_p}\mathrm{d}^3v\:\frac {\partial }{\partial v_\parallel }\big( \boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p\big)\bigg \rangle _\psi \simeq -\bigg \langle\! \int _{V_p} \mathrm{d}^3 v\:\varSigma \bigg \rangle _\psi . \end{equation}

The integration over $v_\parallel$ gives a jump contribution,

(2.42) \begin{equation} \bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \big[\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p \big] \bigg \rangle _\psi \simeq \bigg \langle\! \int _{V_p} \mathrm{d}^3 v\:\varSigma \bigg \rangle _\psi . \end{equation}

In the freely passing region, the jump of a function $\Delta \mathcal{F}^p$ across the trapped–barely passing region is defined as

(2.43) \begin{equation} \Delta \mathcal{F}^p\equiv \mathcal{F}^p(v_\parallel \rightarrow -u^+)-\mathcal{F}^p(v_\parallel \rightarrow -u^-). \end{equation}

By integrating the drift kinetic equation over the freely passing region, we connected the particle transport and the associated particle flux to a jump across the trapped–barely passing region.

Similarly, we can take the parallel momentum moment of (2.36), i.e.

(2.44) \begin{equation} \int _{V_p}\mathrm{d}^3v\: mv_\parallel \big\langle C_p^{(l)}[g]\big\rangle _\tau =-\int _{V_p}\mathrm{d}^3v\:mv_\parallel \langle \varSigma \rangle _\tau \end{equation}

to find that

(2.45) \begin{equation} \bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \big [mv_\parallel \boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p\big ]\, \bigg \rangle _\psi +\bigg \langle\! \int _{V_p}\mathrm{d}^3v\: m\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p\bigg \rangle _\psi \simeq \bigg \langle\! \int _{V_p} \mathrm{d}^3 v\:mv_\parallel \varSigma \bigg \rangle _\psi . \end{equation}

For the energy moment of (2.36),

(2.46) \begin{equation} \int _{V_p}\mathrm{d}^3v\: \frac {mv^2}{2} \big\langle C_p^{(l)}[g]\big\rangle _\tau =-\int _{V_p}\mathrm{d}^3v\:\frac {mv^2}{2} \langle \varSigma \rangle _\tau , \end{equation}

we find that

(2.47) \begin{align} \bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \left [ \frac {mv^2}{2}\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p\right ] \bigg \rangle _\psi +\bigg \langle\! \int _{V_p}\mathrm{d}^3v\:m\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p\bigg \rangle _\psi \simeq \bigg \langle\! \int _{V_p} \mathrm{d}^3 v\:\frac {mv^2}{2}\varSigma \bigg \rangle _\psi .\nonumber\\[2pt] \end{align}

Particle, parallel momentum and energy conservation (2.42), (2.45) and (2.47) all depend on jump contributions. We need to determine the jumps from the trapped–barely passing region.

For the particle transport, we can integrate the drift kinetic equation in (2.18) over trapped–barely passing velocity space and flux surface average. The integration over the collision operator gives the jump contribution

(2.48) \begin{align} \bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \big[\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp} \big]\bigg \rangle _\psi &=-\,\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\varSigma \bigg \rangle _\psi \nonumber \\ & \quad +\, \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \mathcal{L}[g^{t,bp}]\bigg \rangle _\psi +\, \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \mathcal{L}[f_M]\bigg \rangle _\psi . \end{align}

The jump across the trapped–barely passing region in said region is defined as

(2.49) \begin{equation} \Delta \mathcal{F}^{t,bp}=\mathcal{F}^{t,bp}(w\rightarrow \infty )-\mathcal{F}^{t,bp}(w\rightarrow -\infty ). \end{equation}

The jump of a quantity $\mathcal{F}$ has to be the same in the region of overlap between the regions of freely passing and trapped–barely passing particles, according to (2.24), and thus,

(2.50) \begin{equation} \bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \big[\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp} \big]\bigg \rangle _\psi = \bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \big[\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p \big]\bigg \rangle _\psi . \end{equation}

The jump (2.48) can be substituted into the particle transport equation (2.42) to find that

(2.51) \begin{equation} \bigg \langle\! \int \mathrm{d}^3v\:\varSigma \bigg \rangle _\psi = \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \mathcal{L}[g^{t,bp}]\bigg \rangle _\psi +\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \mathcal{L}[f_M]\bigg \rangle _\psi . \end{equation}

The jump in parallel momentum follows from multiplying (2.18) by $mv_\parallel$ and integrating over the trapped–barely passing particle velocity space. We find that

(2.52) \begin{align} &\bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \big[mv_\parallel \boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp} \big]\bigg \rangle _\psi =\bigg \langle\! \int _{V_{tbp}}\mathrm{d}^3v\: m\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp}\bigg \rangle _\psi -\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:mv_\parallel \varSigma \bigg \rangle _\psi \nonumber \\ &\qquad \quad + \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: m(w-u)\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi +\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: m(w-u)\mathcal{L}[f_M]\bigg \rangle _\psi , \end{align}

where we have used the fact that $v_\parallel =w-u$ . It is important at this point to keep the distinction between $v_\parallel$ and $-u$ in the trapped–barely passing region. The parallel momentum equation (2.45) with the jump (2.52) becomes

(2.53) \begin{align} &\bigg \langle\! \int \mathrm{d}^3v\:mv_\parallel \varSigma \bigg \rangle _\psi = \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: m(w-u)\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi +\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: m(w-u)\mathcal{L}[f_M]\bigg \rangle _\psi . \end{align}

Here, we used the fact that the collision operator conserves momentum, i.e.

(2.54) \begin{align} \int \mathrm{d}^3v\: v_\parallel \boldsymbol{\nabla} _v\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}=-\int \mathrm{d}^3v \:\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}=-\int _{V_{p}}\mathrm{d}^3v \:\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p-\int _{V_{tbp}}\mathrm{d}^3v\: \boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp}=0. \end{align}

To find the jump in the energy, we can multiply the drift kinetic equation in (2.18) by $mv^2/2$ and integrate over the trapped–barely passing velocity space:

(2.55) \begin{align} &\bigg \langle\! \int \mathrm{d}\mu \:2\pi B\Delta \left [\frac {mv^2}{2}\boldsymbol{\hat {b}}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp} \right ]\bigg \rangle _\psi = \, \bigg \langle\! \int _{V_{tbp}}\mathrm{d}^3v\: m\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp} \bigg \rangle _\psi -\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\frac {mv^2}{2}\varSigma \bigg \rangle _\psi \nonumber \\ &\qquad \qquad \quad + \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \frac {mv^2}{2}\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi +\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \frac {mv^2}{2}\mathcal{L}[f_M]\bigg \rangle _\psi . \end{align}

The energy moment (2.47) gives

(2.56) \begin{align} \bigg \langle\! \int \mathrm{d}^3v\:\frac {mv^2}{2}\varSigma \bigg \rangle _\psi &= \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \left (\frac {m(w-u)^2}{2}+m\mu B\right )\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi \nonumber \\&\qquad \quad +\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\left (\frac {m(w-u)^2}{2}+m\mu B\right ) \mathcal{L}[f_M]\bigg \rangle _\psi . \end{align}

Here, we used the energy conservation property of the collision operator:

(2.57) \begin{align} &\bigg \langle\! \int \mathrm{d}^3v\:\frac {mv^2}{2}\boldsymbol{\nabla} _v\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}\bigg \rangle _\psi =-\,\bigg \langle\! \int \mathrm{d}^3v\:m\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}\bigg \rangle _\psi \nonumber \\ &\qquad \qquad =-\,\bigg \langle\! \int _{V_{p}} \mathrm{d}^3v\:m\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^p\bigg \rangle _\psi -\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:m\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{M}_{\text{col}}^{t,bp}\bigg \rangle _\psi =0. \end{align}

2.3. Particle, momentum and energy transport

The particle transport expression (2.51) requires evaluating $\langle \int _{V_{tbp}}\mathcal{L}[g^{t,bp}]\mathrm{d}^3v\rangle _\psi$ and $\langle \int _{V_{tbp}}\mathcal{L}[f_M]\mathrm{d}^3v\rangle _\psi$ . The integral $\langle \int _{V_{tbp}}\mathcal{L}[f_M]\mathrm{d}^3v\rangle _\psi$ vanishes due to the integral over $\theta$ . To lowest order, using expression (2.32), the integral $\langle \int _{V_{tbp}}\mathcal{L}[g^{t,bp}]\mathrm{d}^3v\rangle _\psi$ gives

(2.58) \begin{equation} \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi \simeq -\, \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\frac {S}{qR}\frac {r}{R}\frac {\partial }{\partial w} \big[\mathcal{P}(\theta )g^{t,bp}_{0}\big]\bigg \rangle _\psi . \end{equation}

The integral over the $\theta$ -independent piece of $g^{t,bp}_0$ vanishes when averaged over the flux surface. The $\theta$ -dependent piece of $g^{t,bp}_0$ decays for $w\rightarrow \pm \infty$ in order to match the $\theta$ dependence in the freely passing region as argued below (2.28) and shown in (2.35). Thus, the integration over the $\theta$ -dependent piece of $g^{t,bp}_0$ vanishes too.

To next order in $\sqrt {\epsilon }$ of (2.51), using expression (2.18), we find that

(2.59) \begin{align} &\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi \simeq -\,\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\frac {S}{qR}\frac {r}{R}\frac {\partial }{\partial w} \big[\mathcal{P}(\theta )g^{t,bp}_{1}\big]\bigg \rangle _\psi \nonumber \\[4pt] &\quad + \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: \frac {\partial u}{\partial \theta }\frac {1}{qR}\frac {\partial }{\partial w}\big(w g^{t,bp}_{0}\big) \bigg \rangle _\psi -\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\frac {\partial }{\partial \psi }\left [\frac {I}{\varOmega }\frac {1}{qR}\frac {r}{R}\mathcal{P}(\theta )g^{t,bp}_{0}\right ]\bigg \rangle _\psi \nonumber \\[4pt] &\quad +\bigg \langle \left (1+\frac {2I}{\varOmega }\frac {\partial u}{\partial \psi }\right )\int _{V_{tbp}} \mathrm{d}^3v\:\frac {\partial }{\partial w} \left (\frac {wu}{qR}\frac {r}{R}\sin \theta g^{t,bp}_{0}\right )\bigg \rangle _\psi . \end{align}

We show in Appendix A that the $\theta$ -dependent piece of $g_1^{t,bp}$ decays for $w\rightarrow \pm \infty$ and, consequently, the first term in (2.59) vanishes. The $\theta$ dependence of $g^{t,bp}_{0}$ is such that $wg^{t,bp}_0\rightarrow K$ for $w\rightarrow \pm \infty$ , with $K$ being the same constant for $w\rightarrow +\infty$ and $w\rightarrow -\infty$ . Thus, $\langle wg^{t,bp}_{0}\sin \theta \rangle _\psi$ does not jump across the trapped–barely passing region, and the second and fourth terms vanish when integrated over $w$ . The third term is the only term that does not vanish. Thus, we find that

(2.60) \begin{align} \bigg \langle\! \int \mathrm{d}^3v\:\varSigma \bigg \rangle _\psi =\frac {\partial \varGamma }{\partial \psi }, \end{align}

where

(2.61) \begin{align} \varGamma =-\,\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\frac {I}{\varOmega }\frac {1}{qR}\frac {r}{R}\mathcal{P}(\theta )g^{t,bp}_{0}\bigg \rangle _\psi \end{align}

is the ion neoclassical particle flux. We need to calculate $g^{t,bp}_{0}$ to evaluate this integral.

Using the same manipulations that we used for the particle conservation equation, we find that the parallel momentum moment of $\mathcal{L}[g^{t,bp}]$ in the trapped–barely passing region is

(2.62) \begin{align} \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: m(w-u)\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi \simeq \bigg \langle mu\int _{V_{tbp}} \mathrm{d}^3v\:\frac {\partial }{\partial \psi }\left [\frac {I}{\varOmega }\frac {1}{qR}\frac {r}{R}\mathcal{P}(\theta )g^{t,bp}_{0}\right ]\bigg \rangle _\psi \nonumber \\[4pt] -\,\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:mw\frac {S}{qR}\frac {r}{R}\frac {\partial }{\partial w} \big[\mathcal{P}(\theta )g^{t,bp}_{0}\big]\bigg \rangle _\psi . \end{align}

Note that the second term is the result of keeping the small difference between $v_\parallel$ and $-u$ in (2.45). In (2.58) the particle moment of this term did not give a contribution, but here we need to keep it. We can integrate the second term in (2.62) by parts to relate it to the particle flux (2.61):

(2.63) \begin{equation} \bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\: m(w-u)\mathcal{L}[g^{t,bp}]\bigg \rangle _\psi \simeq -mu\frac {\partial \varGamma }{\partial \psi } -mS\frac {\varOmega }{I}\varGamma . \end{equation}

There was no jump term when we integrated by parts because $\langle w g_0^{t,bp}\mathcal{P}(\theta )\rangle _\psi$ does not jump across the trapped–barely passing region for the same reason as argued below (2.59). Note that we have neglected the $\theta$ dependence of $u$ , $S$ and $\varOmega$ , which is small in $\epsilon$ . Going forward, we can drop all $\theta$ dependence in $u$ , $S$ and $\varOmega$ . We can use (2.63) in the parallel momentum relation (2.53) and find the parallel momentum conservation equation

(2.64) \begin{equation} \bigg \langle\! \int \mathrm{d}^3v\:mv_\parallel \varSigma \bigg \rangle _\psi =\frac {\partial }{\partial \psi }\left (-mu\varGamma \right )-m\frac {\varOmega }{I}\varGamma . \end{equation}

Finally, using again the same manipulations that we used for the particle and parallel momentum conservation equation, we find that the energy transport (2.56) is described by the equation

(2.65) \begin{align} \bigg \langle\! \int \mathrm{d}^3v\:\frac {mv^2}{2}\varSigma \bigg \rangle _\psi &=-\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\left (\frac {mu^2}{2}+m\mu B\right )\frac {\partial }{\partial \psi }\left [\frac {I}{\varOmega }\frac {1}{qR}\frac {r}{R}\mathcal{P}(\theta )g^{t,bp}_{0}\right ]\bigg \rangle _\psi \nonumber \\ & \quad -\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:muw\frac {S}{qR}\frac {r}{R}\frac {\partial }{\partial w}\big[\mathcal{P}(\theta )g^{t,bp}_{0}\big]\bigg \rangle _\psi . \end{align}

The second term on the right-hand side is due to the small difference of $v_\parallel$ and $-u$ in the trapped–barely passing region. We can integrate the second term by parts, because $\langle w \mathcal{P}(\theta )g_0^{t,bp}\rangle _\psi$ does not jump across the trapped–barely passing region, and relate it to the particle flux:

(2.66) \begin{align} \bigg \langle\! \int \mathrm{d}^3v\:\frac {mv^2}{2}\varSigma \bigg \rangle _\psi & =-\frac {\partial }{\partial \psi }\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\left (\frac {mu^2}{2}+m\mu B\right )\left [\frac {I}{\varOmega }\frac {1}{qR}\frac {r}{R}\mathcal{P}(\theta )g^{t,bp}_{0}\right ]\bigg \rangle _\psi \nonumber \\& \quad -mu\frac {\partial u}{\partial \psi }\varGamma +muS\frac {\varOmega }{I}\varGamma . \end{align}

The energy conservation (2.66) can be written as

(2.67) \begin{equation} \bigg \langle\! \int \mathrm{d}^3v\:\frac {mv^2}{2}\varSigma \bigg \rangle _\psi =\frac {\partial Q}{\partial \psi } +mu\frac {\varOmega }{I}\varGamma , \end{equation}

with the energy flux

(2.68) \begin{equation} Q= -\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:\left (\frac {mu^2}{2}+m\mu B\right )\frac {I}{\varOmega }\frac {1}{qR}\frac {r}{R}\mathcal{P}(\theta )g^{t,bp}_{0}\bigg \rangle _\psi . \end{equation}

2.4. Electron neoclassical transport

The electron neoclassical transport relations can be derived following the exact same steps as for the ions. The main differences are that the square root of the mass ratio $\sqrt {m_e/m} \ll \epsilon \ll 1$ , where $m_e$ is the electron mass, introduces a second small parameter and that electron–ion collisions must be kept in addition to electron–electron collisions. We showed in Trinczek et al. (Reference Trinczek, Parra and Catto2025, Equation (2.23)) that this leads to a collision operator where we can replace the term $\mathcal{\sf {M}}_\parallel$ by

(2.69) \begin{equation} \mathcal{\sf {M}}_{\parallel e}=\frac {3}{2}\sqrt {\frac {\pi }{2}}\frac {T_e}{m_e}\frac {\nu _{ee}}{x_e}\left [\varTheta (x_e)-\varPsi (x_e)+Z\right ]\!. \end{equation}

Here, $T_e$ is the electron temperature, $\nu _{ee}=4\sqrt {\pi }e^4n_e\log \varLambda /(3T_e^{3/2}m_e^{1/2})$ , $n_e$ is the electron density and $x_e=v/v_{te}$ , where $v_{te}=\sqrt {2T_e/m_e}$ is the electron thermal speed. Note that the definition of $\nu _{ee}$ is consistent with the definition by Trinczek et al. (Reference Trinczek, Parra and Catto2025) (see text above (2.1) in that paper) but differs by a factor of $\sqrt {2}$ from the standard Braginskii definition. Dropping terms that are small in mass ratio and using the fact that the electron charge is $-e$ , we find the electron neoclassical particle flux, i.e.

(2.70) \begin{equation} \varGamma _e=-\,\bigg \langle\! \int _{V_{tbp}}\mathrm{d}^3v\:\frac {I}{\varOmega _e}\frac {1}{qR}\frac {r}{R}\mathcal{P}_e(\theta )g^{t,bp}_{e0}\bigg \rangle _\psi , \end{equation}

where $\varOmega _e=-eB/m_ec$ is the electron Larmor frequency and

(2.71) \begin{equation} \mathcal{P}_e(\theta )=\mu B\sin \theta -\frac {eR}{m_er}\frac {\partial \phi _\theta }{\partial \theta }. \end{equation}

The electron neoclassical energy flux is

(2.72) \begin{equation} Q_e= -\,\bigg \langle\! \int _{V_{tbp}} \mathrm{d}^3v\:m_e\mu B\frac {I}{\varOmega _e}\frac {1}{qR}\frac {r}{R}\mathcal{P}_{e}(\theta )g^{t,bp}_{e0}\bigg \rangle _\psi . \end{equation}

The bootstrap current can also be calculated from the electron distribution function. We follow the derivation for the bootstrap current in § 2.3 of Trinczek et al. (Reference Trinczek, Parra and Catto2025), but generalise it for $\nu _\ast \sim 1$ . The bootstrap current is defined as

(2.73) \begin{equation} j_\parallel ^B\equiv Ze\int \mathrm{d}^3v\: v_\parallel f - e\int \mathrm{d}^3v \: v_\parallel f_e=ZenV_\parallel -en_eV_{\parallel } -e\int \mathrm{d}^3v\: v_\parallel g_e. \end{equation}

Quasineutrality reduces this expression to

(2.74) \begin{equation} j_\parallel ^B=-e\int \mathrm{d}^3v \: v_\parallel g_{e}. \end{equation}

Using the same approach as Trinczek et al. (Reference Trinczek, Parra and Catto2025, § 2.3), we can use the Spitzer–Härm function in (2.74). The Spitzer–Härm function is defined to satisfy the property

(2.75) \begin{align} v_\parallel f_{Me}=C_e[f_{e,SH}] \end{align}

and is of the form

(2.76) \begin{align} f_{e,SH}=\frac {v_\parallel }{\sqrt {2}\nu _{ee}}f_{Me}A_{SH}\big(x_e^2\big), \end{align}

where

(2.77) \begin{align} A_{SH}\big(x_e^2\big)=\sum _i a_i L_i^{3/2}\big(x_e^2\big). \end{align}

Here, $L_i^{3/2}\big(x_e^2\big)$ are Laguerre polynomials and $a_i$ are coefficients determined by (2.75). The first three coefficients for $Z=1$ are $a_0=-1.975$ , $a_1=0.558$ , $a_2=0.015$ . We use (2.75) in (2.74) to write the flux surface averaged bootstrap current as

(2.78) \begin{equation} \big\langle j_\parallel ^B\big\rangle _\psi =-e\Bigg \langle \int \mathrm{d}^3v\: \frac {g_{e}}{f_{Me}}C_e[f_{e,SH}]\Bigg \rangle _\psi =-e\Bigg \langle \int \mathrm{d}^3v\: \frac {f_{e,SH}}{f_{Me}}C_e[g_{e}]\Bigg \rangle _\psi , \end{equation}

where we have used the self-adjointness of the collision operator. The bootstrap current can be expressed as a moment of the collision operator,

(2.79) \begin{equation} \big\langle j_\parallel ^B\big\rangle _\psi \simeq -\Bigg \langle \int \mathrm{d}^3v\: \frac {e}{\sqrt {2}\nu _{ee}} v_\parallel A_{SH} C_e[g_e]\Bigg \rangle _\psi . \end{equation}

At this point, we can employ (2.20),

(2.80) \begin{align} \big\langle j_\parallel ^B\big\rangle _\psi &\simeq -\Bigg \langle \int _{V_{tbp}}\mathrm{d}^3v\: \frac {e}{\sqrt {2}\nu _{ee}} v_\parallel A_{SH} \mathcal{L}[g_e^{t,bp}]\Bigg \rangle _\psi \nonumber \\ & \simeq \Bigg \langle \int _{V_{tbp}}\mathrm{d}^3v\: \frac {e}{\sqrt {2}\nu _{ee}} w A_{SH} \frac {1}{qR}\frac {r}{R}\frac {\partial }{\partial w} \big[\mathcal{P}_e(\theta )g^{t,bp}_{e0}\big]\Bigg \rangle _\psi . \end{align}

Upon integration by parts, where we use the fact that $\langle w\mathcal{P}_e(\theta )g_{e0}^{t,bp}\rangle _\psi$ does not jump across the trapped–barely passing region as argued below (2.59) and $x_e^2\simeq m_e\mu B/T_e$ , the bootstrap current is given by

(2.81) \begin{equation} \big\langle j_\parallel ^B\big\rangle _\psi =-\bigg \langle\! \int _{V_{tbp}}\mathrm{d}^3v\: \frac {e}{\sqrt {2}\nu _{ee}}A_{SH}\frac {1}{qR}\frac {r}{R}\mathcal{P}_e(\theta )g_{e0}^{t,bp}\bigg \rangle _\psi . \end{equation}

All results up to this point are valid to lowest order in $\sqrt {\epsilon }$ and $\nu _\ast \sim 1$ for asymptotically small $\epsilon$ . When $\epsilon$ is a small but finite number, this approach breaks down when $\nu _\ast \gtrsim \epsilon ^{-3/2}$ because a distinction between trapped–barely passing and freely passing particles is no longer possible. To calculate the particle flux and the energy flux, it is necessary to distinguish between the different collisionality regimes. The calculation for the banana regime, $\nu _\ast \ll 1$ , was carried out in Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023). We present the derivation for the plateau regime, $\epsilon ^{-3/2}\gg \nu _\ast \gg 1$ , in the next section.

3.1. Ion transport in the plateau regime

In order to determine $g^l_0$ in the plateau regime, we can use the fact that the second term on the left-hand side of (2.32), the term proportional to $\partial g_0^l/\partial w$ , is smaller than other terms that contain $g^l_{0}$ by $1/\nu _\ast ^{2/3}\ll 1$ , due to (3.3). Thus, (2.32) simplifies to

(3.4) \begin{equation} \frac {\partial }{\partial \theta }\left (\frac {w}{qR}g^l_{0}\right )-\mathsf{M}_\parallel \frac {\partial ^2{g^l_{0}}}{\partial {w}^2}=\frac {I}{\varOmega }\mathcal{P}(\theta )\frac {r}{R}\frac {1}{qR} \mathcal{D}f_M. \end{equation}

It follows from this equation that the lowest order correction to the distribution function in the plateau regime must be of the size $g^l_{0}/f_M\sim \epsilon ^{1/2}/\nu _\ast ^{1/3}\ll 1$ . Equation (3.4) is analytically solvable if the $\theta$ dependence of $u$ and $\varOmega$ is neglected as small in $\epsilon$ and the form of the potential is taken to be

(3.5) \begin{equation} \phi _\theta =\phi _c(\psi )\cos \theta +\phi _s(\psi )\sin \theta . \end{equation}

This choice will be justified for concentric circular flux surfaces in the calculation of the electric potential in § 3.3. Introducing the amplitudes

(3.6) \begin{align} A_S\equiv u^2+\mu B-\frac {Ze\phi _c R}{mr}, && A_C\equiv \frac {Ze\phi _sR}{m r}, \end{align}

we can write

(3.7) \begin{equation} \mathcal{P}(\theta )=A_S\sin \theta +A_C\cos \theta \end{equation}

in (3.4). The full derivation of the solution to (3.4) is in Appendix B, which follows Su & Oberman (Reference Su and Oberman1968). The solution is

(3.8) \begin{equation} g^l_{0}=\left [A_SF_s(\xi ,\theta )+A_CF_c(\xi ,\theta )\right ]\frac {I}{\varOmega } \frac {\mathcal{D}}{v_{\textrm {ref}}}\frac {r}{R} f_M(w=0), \end{equation}

where $\xi =w/v_{\textrm {ref}}$ , $v_{\textrm {ref}}=(qR\mathsf{M}_\parallel )^{1/3}$ ,

(3.9) \begin{equation} F_s(\xi ,\theta )=\int _0^\infty \mathrm{d}p \:\exp \left (-\frac {p^3}{3}\right )\sin (\theta -p\xi ) \end{equation}

and

(3.10) \begin{equation} F_c(\xi ,\theta )=\int _0^\infty \mathrm{d}p\:\exp \left (-\frac {p^3}{3}\right )\cos (\theta -p\xi ). \end{equation}

Importantly, $\langle g^l_{0}\rangle _\tau =0$ and $g^l_{0}$ matches (2.28) for $\xi \rightarrow \infty$ , a result demonstrated in Appendix B. We also need the integral

(3.11) \begin{equation} \lim _{W\rightarrow \infty } \int _{-W}^{W}\text{d}w\:g^l_{0}=\frac {{\pi}I}{\varOmega }\frac {r}{R}\mathcal{D}\left (A_C\cos \theta +A_S\sin \theta \right )f_M(w=0) \end{equation}

to calculate the particle flux $\varGamma$ in (2.61). The derivation of (3.11) is presented in Appendix B.

We can use (3.11) in (2.61) to find the ion neoclassical particle flux, i.e.

(3.12) \begin{equation} \varGamma =-\frac {BI^2\pi ^2}{qR\varOmega ^2}\left (\frac {r}{R}\right )^2\int \text{d}\mu \left (A_S^2+A_C^2\right )\mathcal{D}f_M(w=0), \end{equation}

where we also carried out the flux surface average. The integration over $\mu$ gives

(3.13) \begin{align} \varGamma & =-\frac {\sqrt {\pi }}{4}\frac {n}{qR}\frac {I^2}{\varOmega ^2}\left (\frac {r}{R}\right )^2\left (\frac {2T}{m}\right )^{3/2}\exp \left (-\frac {m(u+V_\parallel )^2}{2T}\right )\nonumber \\[4pt] & \quad \times \Bigg \lbrace \frac {1}{2}\left [\left (\frac {mu^2}{T}-\frac {ZeR\phi _c}{Tr}\right )^2+\left (\frac {ZeR\phi _s}{Tr}\right )^2\right ]\mathcal{D}_{-3/2} +\left (\frac {mu^2}{T}-\frac {ZeR}{Tr}\phi _c\right )\mathcal{D}_{-1/2}+\mathcal{D}_{1/2}\Bigg \rbrace .\nonumber\\[2pt] \end{align}

Here, we have introduced the notation,

(3.14) \begin{equation} \mathcal{D}_l\equiv \frac {\partial }{\partial \psi }\ln p-\frac {m(u+V_\parallel )}{T}\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}\right )+\left [\frac {m(u+V_\parallel )^2}{2T}+l\right ]\frac {\partial }{\partial \psi }\ln T, \end{equation}

where $l$ is a rational number.

The energy flux $Q$ in (2.68) is

(3.15) \begin{equation} Q=-\frac {BTI^2\pi ^2}{qR\varOmega ^2}\left (\frac {r}{R}\right )^2\int \text{d}\mu \:\left (\frac {mu^2}{2T}+\frac {m \mu B}{T}\right )\left (A_S^2+A_C^2\right )\mathcal{D}f_M(w=0), \end{equation}

which can be calculated to give

(3.16) \begin{align} Q & = \frac {mu^2}{2}\varGamma _i-\frac {3\sqrt {\pi }}{4}\frac {nT}{qR}\frac {I^2}{\varOmega ^2 }\left (\frac {r}{R}\right )^2 \left (\frac {2T}{m}\right )^{3/2}\exp \left (-\frac {m(u+V_\parallel )^2}{2T}\right )\nonumber \\& \quad\times \Bigg \lbrace \frac {1}{6}\left [\left (\frac {mu^2}{T}-\frac {Ze\phi _c R}{Tr}\right )^2 + \left (\frac {Ze\phi _s R}{Tr}\right )^2\right ]\mathcal{D}_{-1/2}\nonumber \\ & \qquad\qquad +\frac {2}{3}\left (\frac {mu^2}{T}-\frac {Ze\phi _cR}{Tr}\right )\mathcal{D}_{1/2}+\mathcal{D}_{3/2}\Bigg \rbrace . \end{align}

In the limit of weak gradients, the poloidal variation vanishes, and the neoclassical ion particle flux (3.13) and energy flux (3.16) agree with the results of weak gradient neoclassical theory, as shown in Appendix C.

We compare the results in (3.13) and (3.16) to previous work on strong gradient neoclassical transport in the plateau regime by Pusztai & Catto (Reference Pusztai and Catto2010) and Seol & Shaing (Reference Seol and Shaing2012) in Appendix D, and find that our results disagree due to inconsistencies in previous work.

3.2. Electron transport in the plateau regime

The electron distribution function $g^l_{e0}$ can be derived following the exact same steps as for the ion distribution function. The result is

(3.17) \begin{equation} g^l_{e0}=\left [A_{Se}F_s(\xi _e,\theta )+A_{Ce}F_c(\xi _e,\theta )\right ]\frac {I}{\varOmega _e} \frac {\mathcal{D}_e}{v_{\textrm {ref,e}}}\frac {r}{R} f_{Me}(v_\parallel =0), \end{equation}

where $\xi _e=v_\parallel /v_{\textrm {ref,e}}$ , $v_{\textrm {ref,e}}=(qR\mathsf{M}_{\parallel e})^{1/3}$ ,

(3.18) \begin{align} A_{Se}=\mu B+\frac {e\phi _c R}{m_e r}, \qquad \qquad \qquad \qquad A_{Ce}=-\frac {e\phi _sR}{m_er}, \\[-12pt]\nonumber \end{align}

(3.19) \begin{align} \mathcal{D}_e=\frac {\partial }{\partial \psi }\ln {p_e} +\frac {m_e(u+ V_\parallel )}{T}\frac {\varOmega _e}{I} +\left (\frac {m_e\mu B}{T_e}-\frac {5}{2}\right )\frac {\partial }{\partial \psi }\ln {T_e} \end{align}

and $p_e=n_eT_e$ is the electron pressure. We can use the distribution function (3.17) in (2.70) and integrate to find that

(3.20) \begin{align} \varGamma _e=-\frac {\sqrt {\pi }}{4}\frac {n_e}{qR}\frac {I^2}{\varOmega _e^2}\left (\frac {r}{R}\right )^2\left (\frac {2T_e}{m_e}\right )^{3/2} &\Bigg \lbrace \frac {1}{2}\left [\left (\frac {e\phi _cR}{T_er}\right )^2+\left (\frac {e\phi _s R}{T_e r}\right )^2\right ]\mathcal{D}_{e,-3/2}\nonumber \\ & \quad \, +\frac {e\phi _c R}{T_e r}\mathcal{D}_{e,-1/2}+\mathcal{D}_{e,1/2}\Bigg \rbrace , \end{align}

where

(3.21) \begin{equation} \mathcal{D}_{e,l}\equiv \frac {\partial }{\partial \psi }\ln p_e+\frac {m_e(u+V_\parallel )}{T_e}\frac {\varOmega _e}{I}+l\frac {\partial }{\partial \psi }\ln T_e. \end{equation}

Similarly, the electron energy flux (2.72) yields

(3.22) \begin{align} Q_e=-\frac {3\sqrt {\pi }}{4}\frac {n_eT_e}{qR}\frac {I^2}{\varOmega _e^2}\left (\frac {r}{R}\right )^2\left (\frac {2T_e}{m_e}\right )^{3/2} &\Bigg \lbrace \frac {1}{6}\left [\left (\frac {e\phi _cR}{T_er}\right )^2 +\left (\frac {e\phi _s R}{T_e r}\right )^2\right ]\mathcal{D}_{e,-1/2}\nonumber \\ &\qquad \qquad +\frac {2}{3}\frac {e\phi _c R}{T_e r}\mathcal{D}_{e,1/2}+\mathcal{D}_{e,3/2}\Bigg \rbrace . \end{align}

In the limit of weak gradients, the neoclassical electron particle flux (3.20) and energy flux (3.22) agree with the results of weak gradient neoclassical theory, as shown in Appendix C.

We can insert the electron distribution function (3.17) in (2.81) to calculate the bootstrap current. We take the integral (2.81) by using properties of the Laguerre polynomials. The derivation is shown in detail in Appendix E. The final result for the bootstrap current reads

(3.23) \begin{align} \big\langle j_\parallel ^B\big\rangle _\psi &= -\frac {1}{8}\sqrt {\frac {\pi }{2}}\frac {en_e}{\nu _{ee}qR}\frac {I}{\varOmega _e}\left (\frac {r}{R}\right )^2\left (\frac {2T_e}{m_e}\right )^{3/2}\Bigg \lbrace \left [\left (\frac {e\phi _cR}{T_er}\right )^2+\left (\frac {e\phi _sR}{T_er}\right )^2\right ]\nonumber \\[4pt] &\quad \times \sum _{i=0} \frac {2a_i\varGamma (i+ {3}/{2})}{\sqrt {\pi }\varGamma (1+i)}\mathcal{D}_{e,\alpha _i} +2\frac {e\phi _cR}{T_e r}\left [a_0\mathcal{D}_{e,-({1}/{2})}+\sum _{i=1}\frac {a_i\varGamma (i+ {1}/{2})}{\sqrt {\pi }\varGamma (1+i)}\mathcal{D}_{e,\beta _i}\right ]\nonumber \\[4pt] &\quad +\left [2a_0\mathcal{D}_{e,1/2}-a_1\mathcal{D}_{e,13/2}-\sum _{i=2} \frac {a_i\varGamma (i-{1}/{2})}{\sqrt {\pi }\varGamma (1+i)} \mathcal{D}_{e,\gamma _i} \right ]\Bigg \rbrace , \end{align}

where

(3.24) \begin{align} \alpha _i=-\frac {3/2+5i}{1+2i},&&\beta _i=\frac {1-10i}{4i-2},&&\gamma _i=-\frac {(3+10i)(i-1/2)}{4}. \end{align}

Note that $\varGamma$ here is the Euler Gamma function and not the particle flux. The result agrees with the weak gradient limit in Appendix C and we compare it to previous results by Pusztai & Catto (Reference Pusztai and Catto2010) in Appendix D.

3.3. Electric potential

The calculation of the electric potential in the plateau regime follows the derivation in the banana regime in § 4.4 of Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023). We assume a Boltzmann response for the electrons such that quasineutrality reads

(3.25) \begin{equation} Zn_\theta =\frac {en_e}{T_e}\phi _\theta . \end{equation}

The poloidal density variation as derived in Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023, Equation (4.56)) is

(3.26) \begin{equation} n_\theta =\lim _{W\rightarrow \infty }\int _0^\infty \mathrm{d}\mu \int _{-W}^W\mathrm{d}w\: 2\pi B (g^l-\langle g^l\rangle _\psi )+\int _{V_p}\mathrm{d}^3v \:2\pi B(g^p-\langle g^p\rangle _\psi ) \equiv n_\theta ^l+n_\theta ^p, \end{equation}

where the integration is performed over both the collisional layer and the freely passing region of velocity space, and hence, one must calculate the contributions from $g^l$ and $g^p$ separately. The freely passing particle contribution is the same as in the banana regime. The distribution function $g^p$ diverges as $v_\parallel \rightarrow -u$ ; see (2.28). This divergence is compensated for by the contribution from the collisional layer. We comment on the validity of this integration in Appendix B.2. The freely passing particle contribution was calculated by Trinczek et al. (Reference Trinczek, Parra, Catto, Calvo and Landreman2023) in their equation (4.57) and yields

(3.27) \begin{align} n_\theta ^p & =- n\frac {Ir}{\varOmega R}\Bigg \lbrace \sqrt {\frac {2T}{m}}J \bigg [\left (\frac {mV_\parallel ^2}{T}\cos \theta +\cos \theta -\frac {Ze\phi _\theta R}{Tr}\right )\left (\frac {\partial }{\partial \psi }\ln p-\frac {3}{2}\frac {\partial }{\partial \psi }\ln T\right )\nonumber \\ &\quad +\cos \theta \frac {\partial }{\partial \psi }\ln T \bigg ] +\left [1-2\sqrt {\frac {m}{2T}}(V_\parallel +u)J\right ]\Bigg \lbrace (V_\parallel -u)\cos \theta \left (\frac {\partial }{\partial \psi }\ln p-\frac {3}{2}\frac {\partial }{\partial \psi }\ln T\right )\nonumber \\ &\quad -(V_\parallel +u)\Bigg [\left (\frac {mV_\parallel ^2}{2T}+\frac {1}{2}\right )\cos \theta -\frac {Ze\phi _\theta R}{2Tr}\Bigg ]\frac {\partial }{\partial \psi }\ln T+\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}\right )\nonumber \\ &\quad \times \left [\left (\frac {m u^2}{T}+1-\frac {m(V_\parallel +u)^2}{T}\right )\cos \theta -\frac {Ze\phi _\theta R}{Tr}\right ]\Bigg \rbrace \nonumber \\ &\quad +\left [1+2\frac {m(V_\parallel +u)^2}{2T}-4\left (\frac {m}{2T}\right )^{3/2}(V_\parallel +u)^3J\right ]\cos \theta \!\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}+\frac {V_\parallel -u}{2}\frac {\partial }{\partial \psi }\ln T\right )\!\Bigg \rbrace \nonumber \\ &\hspace{23.5pc}-2n\frac {r}{R}\cos \theta , \end{align}

where $J$ is defined as

(3.28) \begin{equation} J=\frac {\sqrt {\pi } }{2}\exp \left (-\frac {m(u+V_\parallel )^2}{2T}\right )\textrm {erfi}\left (\sqrt {\frac {m}{2T}}(u+V_\parallel )\right )\!. \end{equation}

Note that the dependence on $\cos \theta$ in this expression is coming from the magnetic field strength. We provide additional comments on the derivation of (3.27) and (3.28) in Appendix B.3.

In the banana regime the contribution from trapped and barely passing particles to the density vanishes because the $\theta$ -dependent part of the distribution function is odd in $w$ . In the plateau regime the contribution from the particles in the collisional layer does not vanish to lowest order. The term $\int \text{d}\mu \int \text{d}w\:2\pi B g^l_{0}$ needs to be evaluated, whereas $\langle g^l_{0}\rangle _\tau$ vanishes. The velocity integration of $g^l_{0}$ can be performed using the solution for $g^l_{0}$ , (3.8), the relation (3.11) and assuming (3.5):

(3.29) \begin{align} n^l_\theta \simeq& \lim _{W\rightarrow \infty }\int _0^\infty \text{d}\mu \int _{-W}^W\text{d}w \: 2\pi B g^l_{0}= n\sqrt {\frac {m}{T}}\sqrt {\frac {\pi }{2}}\exp \left (-\frac {m(u+V_\parallel )^2}{2T}\right )\nonumber \\[4pt] &\times \Bigg \lbrace \frac {I}{\varOmega }\Bigg [\frac {Ze}{m}\phi _s\cos \theta +\left (u^2\frac {r}{R}-\frac {Ze}{m}\phi _c\right )\sin \theta \Bigg ] \mathcal{D}_{-3/2}+\frac {ITr}{\varOmega m R}\sin \theta \mathcal{D}_{-1/2}\Bigg \rbrace .\nonumber\\[3pt] \end{align}

Combining (3.27) and (3.29) in (3.25) yields an equation for the poloidally varying component of the electric potential $\phi _\theta =\phi _c(\psi )\cos \theta +\phi _s(\psi )\sin \theta$ . The terms proportional to $\cos \theta$ give

(3.30) \begin{align} \Bigg \lbrace \frac {en_e}{T_e} & -\frac {Z^2enI}{T\varOmega }\left[\sqrt {\frac {2T}{m}}J\left (\frac {\partial }{\partial \psi }\ln p-\frac {3}{2}\frac {\partial }{\partial \psi } \ln T\right )+\left [1-2\sqrt {\frac {m}{2T}}(V_\parallel +u) J\right ]\left(\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}\right.\right.\nonumber \\ &\quad \left.\left.-\frac {(V_\parallel +u)}{2}\frac {\partial }{\partial \psi }\ln T\right)\right] \Bigg \rbrace \phi _c=-Zn\frac {Ir}{\varOmega R}\Bigg \lbrace \sqrt {\frac {2T}{m}}J\bigg [\left (\frac {mV_\parallel ^2}{T}+1\right )\nonumber \\ &\quad\times \left (\frac {\partial }{\partial \psi }\ln p-\frac {3}{2}\frac {\partial }{\partial \psi } \ln T\right )+\frac {\partial }{\partial \psi }\ln T\bigg ] +\left [1-2\sqrt {\frac {m}{2T}}(V_\parallel +u) J\right ]\nonumber \\ &\quad \times \bigg [(V_\parallel -u)\left (\frac {\partial }{\partial \psi }\ln p-\frac {3}{2}\frac {\partial }{\partial \psi }\ln T\right )+\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}\right )\left (\frac {mu^2}{T}+1-\frac {m(V_\parallel +u)^2}{T}\right )\nonumber \\ &\quad-\frac {V_\parallel +u}{2}\left (\frac {mV_\parallel ^2}{T}+1\right )\frac {\partial }{\partial \psi }\ln T\bigg ]\nonumber \\ &\quad +\left [1+\frac {m}{T}(V_\parallel +u)^2-4\left (\frac {m}{2T}\right )^{3/2}\left (V_\parallel +u\right )^3 J\right ]\left (\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}+\frac {V_\parallel -u}{2}\frac {\partial }{\partial \psi }\ln T\right )\Bigg \rbrace \nonumber \\ &\quad+Zn\sqrt {\frac {m}{T}}\sqrt {\frac {\pi }{2}}\exp \left (-\frac {m(u+V_\parallel )^2}{2T}\right )\frac {Ze}{m}\phi _s \mathcal{D}_{-3/2}-2Zn\frac {r}{R}, \end{align}

and all terms proportional to $\sin \theta$ give

(3.31) \begin{align} \Bigg \lbrace \frac {en_e}{T_e} & -\frac {Z^2enI}{T\varOmega }\bigg [\sqrt {\frac {2T}{m}}J\left (\frac {\partial }{\partial \psi }\ln p-\frac {3}{2}\frac {\partial }{\partial \psi } \ln T\right )+\left [1-2\sqrt {\frac {m}{2T}}(V_\parallel +u)J\right ]\nonumber \\[3pt]& \times \left(\frac {\partial V_\parallel }{\partial \psi }-\frac {\varOmega }{I}-\frac {V_\parallel +u}{2}\frac {\partial }{\partial \psi }\ln T\right)\bigg ] \Bigg \rbrace \phi _s=Zn\sqrt {\frac {m}{T}}\sqrt {\frac {\pi }{2}}\exp \left (-\frac {m(u+V_\parallel )^2}{2T}\right )\nonumber \\[3pt] & \times \Bigg [ \frac {I}{\varOmega }\left (u^2\frac {r}{R}-\frac {Ze}{m}\phi _c\right )\mathcal{D}_{-3/2} +\frac {TIr}{m\varOmega R}\mathcal{D}_{-1/2}\Bigg ] . \end{align}

Equations (3.30) and (3.31) form a linear system that can be solved for $\phi _c$ and $\phi _s$ . This result is consistent with our choice of $\phi _\theta$ in (3.5).

The usual neoclassical results for the electric potential can be retrieved by taking the limit $u\ll v_{t}$ , $V_\parallel \ll v_{t}$ and weak gradients. In this weak gradient limit, (3.30) and (3.31) reduce to

(3.32) \begin{align} \phi _c=0 \end{align}

and

(3.33) \begin{align} \left ( \frac {en_e}{T_e}+\frac {Z^2en}{T}\right )\phi _s =Zn\sqrt {\frac {T}{m}}\sqrt {\frac {\pi }{2}} \frac {Ir}{\varOmega R}\left [\frac {\partial }{\partial \psi }\ln p+\frac {m(u+V_\parallel )}{T}\frac {\varOmega }{I}-\frac {1}{2}\frac {\partial }{\partial \psi }\ln T\right ]\!.\nonumber\\[3pt] \end{align}

The vanishing parallel friction force in the weak gradient limit gives a relation for the mean parallel flow (C3), which can be used to simplify (3.33) to

(3.34) \begin{equation} \left ( \frac {en_e}{T_e}+\frac {Z^2en}{T}\right )\phi _s =-Zn\sqrt {\frac {T}{m}}\sqrt {\frac {\pi }{2}} \frac {Ir}{\varOmega R}\frac {\partial }{\partial \psi }\ln T. \end{equation}

A non-vanishing up–down asymmetry is consistent with the results of Hinton & Rosenbluth (Reference Hinton and Rosenbluth1973).

4.1. Radial force balance

Radial force balance as defined in (4.12) gives the radial electric field directly from the pressure gradient. The radial electric field for the input profiles in figure 1 is shown in figure 2. We observe a deep radial electric field well, as observed in pedestals (McDermott et al. Reference McDermott2009; Viezzer et al. Reference Viezzer2013). Radial force balance is independent of the mean parallel flow, so the radial electric field is the same for any choice of $\alpha$ in (4.16).

Figure 2.

The radial electric field as determined by the radial force balance (4.12).

Figure 3.

The ion neoclassical particle flux is mostly positive for $\alpha =-0.25$ and negative for $\alpha =0.59$ . The ion neoclassical energy flux exceeds the weak gradient energy flux in the strong gradient region.

Figure 4.

On the left is the poloidal variation of the electric potential $\bar {\phi }_\theta =\bar {\phi }_c\cos \theta +\bar {\phi }_s\sin \theta$ for $\alpha =0.59$ , and on the right for $\alpha =-0.25$ . Both cases show in–out and up–down asymmetry. The corresponding amplitudes $\bar \phi _c$ and $\bar \phi _s$ are shown as a function of $\bar \psi$ below the respective two-dimensional plots.

Figure 4. Long description

Panel A: A contour plot shows the poloidal variation of the electric potential for a specific case. The x-axis and y-axis range from -10 to 10. The color scale on the right indicates the electric potential values, ranging from -0.4 to 0.4. The plot exhibits in-out and up-down asymmetry. Below this contour plot, a line graph displays the amplitudes as a function of a variable. The x-axis ranges from 6 to 10, and the y-axis ranges from -0.4 to 0.4. Two lines, one in cyan and one in magenta, represent different amplitudes. Panel B: Another contour plot shows the poloidal variation of the electric potential for a different case. The x-axis and y-axis range from -10 to 10. The color scale on the right indicates the electric potential values, ranging from -0.4 to 0.4. This plot also exhibits in-out and up-down asymmetry. Below this contour plot, a line graph displays the amplitudes as a function of a variable. The x-axis ranges from 6 to 10, and the y-axis ranges from -0.4 to 0.4. Two lines, one in cyan and one in magenta, represent different amplitudes.

The ion neoclassical particle and energy fluxes are shown in figure 3. The ion neoclassical particle flux is mostly positive for $\alpha =-0.25$ and negative for $\alpha =0.59$ . Radial force balance thus implies that turbulence must exist in the system such that ambipolarity is achieved through a balance between the neoclassical particle flux (which we calculate) and the turbulent particle flux (not calculated in this paper); see Trinczek & Parra (Reference Trinczek and Parra2026) for a more extensive discussion.

We can compare the energy flux for both values of $\alpha$ with the weak gradient result (C4). The weak gradient energy flux (C4) does not explicitly depend on the mean parallel flow because the flow is known in the weak gradient limit; see (C3). This is not the case in the strong gradient neoclassical expression (4.8), where the dependence on $V_\parallel$ is kept explicitly. The strong gradient neoclassical energy flux exceeds the weak gradient theory result for both values of $\alpha$ . The maximum energy flux exceeds the maximum energy flux according to weak gradient theory for $\alpha =-0.25$ by $\bar {Q}^{FB}_{\textit{max}}/\bar {Q}^{wg}_{\textit{max}}\simeq 1.27$ and by $\bar {Q}^{FB}_{\textit{max}}/\bar {Q}^{wg}_{\textit{max}}=3.41$ for $\alpha =0.59$ . This is a significant increase of the neoclassical energy flux when strong gradient effects are kept. Clearly, the energy fluxes within strong gradient neoclassical theory depend strongly on the mean parallel flow. It is also worth noting that Seol & Shaing (Reference Seol and Shaing2012) claimed that strong gradient effects will always reduce the energy flux in comparison to weak gradient neoclassical theory. Our example case is in clear contradiction to this result. Further discussion of this discrepancy is provided in Appendix D.

We can use (3.5), (3.30) and (3.31) to obtain the poloidal variation of the electric potential for force balance for the two choices of $\alpha$ . In figure 4 we show $\bar {\phi }$ as a function of $\bar {\psi }$ and $\theta$ for $\bar {\psi }\in [6,10]$ , as well as the respective amplitudes $\bar \phi _c$ and $\bar \phi _s$ . Note that figure 4 is not true to scale in radius – the edge is enlarged relative to the rest of the poloidal plane. There is stronger in–out asymmetry for $\alpha =-0.25$ and stronger up–down asymmetry for $\alpha =0.59$ .

The bootstrap current is shown in figure 5. Again, we can compare our modification to the limit of weak gradient neoclassical theory (C10). For $\alpha =-0.25$ , the changes are not significant. However, for $\alpha =0.59$ , the maximum bootstrap current calculated with strong gradient neoclassical theory is smaller than the peak bootstrap current predicted by weak gradient neoclassical theory by $\bar {j}^{B,NA}_{\max}/\bar {j}^{B,wg}_{\max}=0.89$ . In this case, keeping strong gradient effects results in a lower bootstrap current than weak gradient neoclassical theory predicts.

The electron neoclassical particle flux, which is also plotted in figure 5, shows little to no deviation from weak gradient neoclassical theory for both $\alpha =-0.25$ and $\alpha =0.59$ . The neoclassical electron particle flux is just one contribution to the total electron flux and cannot balance the neoclassical ion particle flux by itself due to the smallness of $\varGamma _e$ in the mass ratio. The turbulent electron particle flux has to be of the same size as the ion neoclassical particle flux to satisfy ambipolarity.

Figure 5.

The bootstrap current is almost identical to weak gradient neoclassical theory for $\alpha =-0.25$ but smaller than weak gradient neoclassical predictions for $\alpha =0.59$ . The electron neoclassical particle flux is closer to the weak gradient limit for $\alpha =0.59$ and larger for $\alpha =-0.25$ .

Figure 6.

The radial electric field is shown as determined by neoclassical ambipolarity (4.14) for $\alpha =0.59$ and $\alpha =-0.25$ .

4.2. Neoclassical ambipolarity

In the absence of a parallel momentum source, neoclassical ambipolarity (4.14) can be used to determine the radial electric field by solving $\bar \varGamma =0$ , with $\bar \varGamma$ given by (4.6), for $\bar u$ for a specified $\bar V$ . This approach was not possible in § 4.1 because $\bar \varGamma$ was not known and, hence, force balance was used to determine the radial electric field. Solving (4.6) for $\bar \varGamma =0$ gives the radial electric field profiles for $\alpha =-0.25$ and $\alpha =0.59$ in figure 6. Here, the radial electric field depends on the choice of the mean parallel flow and we get a larger radial electric field for $\alpha =0.59$ . We can compare these profiles with the radial force balance in figure 2 and find that neoclassical ambipolarity predicts a weaker radial electric field for $\alpha =-0.25$ and a stronger one for $\alpha =0.59$ .

By definition, the ion neoclassical particle flux vanishes to lowest order for neoclassical ambipolarity. The ion neoclassical energy flux is shown in figure 7. Again, the energy flux with strong gradient effects exceeds the weak gradient prediction for both choices of the mean parallel flow. For $\alpha =0.59$ , the energy flux exceeds the weak gradient theory prediction by $\bar {Q}^{NA}_{\textit{max}}/\bar {Q}^{wg}_{\textit{max}}=1.61$ . For $\alpha =-0.25$ , the energy flux is very similar to the one we determined using force balance and $\bar {Q}^{NA}_{\textit{max}}/\bar {Q}^{wg}_{\textit{max}}\simeq 1.25$ . This shows that the mean parallel flow has a strong impact on the deviation of fluxes from weak gradient neoclassical theory in the case of neoclassical ambipolarity as well.

Figure 7.

The ion energy flux exceeds the weak gradient neoclassical prediction for $\alpha =0.59$ and $\alpha =-0.25$ .

The poloidal variation of the electric potential for neoclassical ambipolarity is shown in figure 8. For $\alpha =0.59$ , the variation is more in–out asymmetric whereas we observe stronger up–down asymmetry for $\alpha =-0.25$ . Overall, the peak amplitudes of the variation is weaker than with radial force balance. Noticeably, the amplitude for up–down asymmetry is positive for neoclassical ambipolarity and negative for radial force balance for $\alpha =0.59$ . For $\alpha =-0.25$ , the poloidal variation is qualitatively fairly similar for radial force balance and neoclassical ambipolarity.

Figure 8.

On the left is the poloidal variation of the electric potential $\bar {\phi }_\theta =\bar {\phi }_c\cos \theta +\bar {\phi }_s\sin \theta$ for $\alpha =0.59$ . The amplitudes $\bar \phi _c$ and $\bar \phi _s$ are given as functions of $\bar \psi$ in the lower figures. On the right is the poloidal variation for $\alpha =-0.25$ that shows stronger up–down asymmetry and, hence, a larger absolute value of $\bar \phi _s$ .

Figure 8. Long description

Panel A: The top left graph is a heat map showing the poloidal variation of the electric potential. The x-axis and y-axis both range from -10 to 10. The color bar on the right indicates the values of the electric potential, ranging from -0.4 to 0.4. The bottom left graph is a line graph showing the amplitudes as functions of psi bar. The x-axis ranges from 6 to 10, and the y-axis ranges from -0.4 to 0.4. Two lines are plotted: one in cyan representing phi bar c and another in magenta representing phi bar s. Panel B: The top right graph is a heat map showing the poloidal variation of the electric potential for a different condition. The x-axis and y-axis both range from -10 to 10. The color bar on the right indicates the values of the electric potential, ranging from -0.4 to 0.4. The bottom right graph is a line graph showing the amplitudes as functions of psi bar. The x-axis ranges from 6 to 10, and the y-axis ranges from -0.4 to 0.4. Two lines are plotted: one in cyan representing phi bar c and another in magenta representing phi bar s.

The bootstrap current is almost identical to weak gradient neoclassical theory for $\alpha =-0.25$ , as can be seen in figure 9. This is similar to the case where force balance was used. However, for $\alpha =0.59$ , we observe an increased bootstrap current of $\bar {j}^{B,NA}_{\textit{max}}/\bar {j}^{B,wg}_{\textit{max}}=1.25$ . This is the opposite to what we observed in radial force balance, where a reduced bootstrap current was found for $\alpha =0.59$ .

The electron neoclassical particle flux strongly exceeds the weak gradient prediction for $\alpha =0.59$ .

Figure 9.

The bootstrap current exceeds weak gradient neoclassical predictions for $\alpha =0.59$ and is indistinguishable from weak gradient theory for $\alpha =-0.25$ . The neoclassical electron flux exceeds weak gradient results in both cases, but more so for $\alpha =0.59$ .