# The maximum-$J$ property in quasi-isodynamic stellarators

E. Rodríguez  P. Helander  and A. G. Goodman  Max Planck Institute for Plasma
Physics, 17491 Greifswald, Germany

###### Abstract

Some stellarators tend to benefit from favourable average magnetic curvature
for trapped particles when the plasma pressure is sufficiently high. This so-
called maximum-$J$-property has several positive implications such as good
fast-particle confinement, magnetohydrodynamic stability, and suppression of
certain trapped-particle instabilities. This property cannot be attained in
quasisymmetric stellarators, in which deeply trapped particles experience
average bad curvature and therefore precess in the diamagnetic direction close
to the magnetic axis. However, quasi-isodynamic stellarators offer greater
flexibility and allow the average curvature to be favourable and the
precession to be reversed. We find that it is possible to design such
stellarators so that the maximum-$J$ condition is satisfied for the great
majority of all particles, even when the plasma pressure vanishes. The
qualitative properties of such a stellarator field can be derived analytically
by examining the most deeply and the most shallowly trapped particles,
although some small fraction of the latter will inevitably not behave as
desired. However, through numerical optimisation, we construct a vacuum field
in which 99.6% of all trapped particles satisfy the maximum-$J$ condition.

##  1 Introduction

A basic requirement for magnetic fusion devices is the ability to confine
collisionless particle orbits, a property which not all stellarators possess.
Although the existence of nested magnetic surface is sufficient to confine
circulating (untrapped) particle orbits to lowest order in the smallness of
the gyroradius (Helander, 2014), particles that are trapped in local minima
along the magnetic field tend to drift out of the plasma unless the field is
carefully tailored to avoid this phenomenon.

If we write the magnetic field in Clebsch form (D’haeseleer et al., 2012),
$\mathbf{B}=\nabla\psi\times\nabla\alpha$, then the coordinates
$(\psi,\alpha)$ are constant along $\bf B$ and can thus be used to label field
lines. In a field tracing out toroidal magnetic surfaces, we take $\psi$ to be
the magnetic flux enclosed by such a surface divided by $2\pi$, and
$\alpha\in[0,2\pi]$ thus labels the different field lines on each surface. In
Boozer coordinates $\alpha=\theta-\iota\varphi$, where $\theta$ denotes the
poloidal angle, $\varphi$ the toroidal angle, and $\iota(\psi)$ the rotational
transform (Boozer, 1981).

As a trapped particle drifts across the field, the second adiabatic invariant
(Hastie et al., 1967)

$\mathcal{J}_{\parallel}=\int mv_{\|}\mathrm{d}\ell$

is conserved to high accuracy if the gyroradius is small. Here, $m$ denotes
the mass of the particle in question, $v_{\|}=v\sqrt{1-\lambda B}$ its speed
along the magnetic field line (if the electric field is sufficiently small),
and $\ell$ the arc length in this direction. The ratio of the magnetic moment
to the kinetic energy is denoted by $\lambda=v_{\perp}^{2}/(v^{2}B)$, where
$B$ is the field strength, and the integration is carried out between two
consecutive bounce points where $v_{\|}$ vanishes. $\mathcal{J}_{\|}$ is thus
a function of the variables $(\psi,\alpha,v,\lambda)$.

The net cross-field magnetic drift of a magnetically trapped particle of
charge $q$ travelling between two consecutive bounce points is given by the
Hamiltonian equations (Kadomtsev & Pogutse, 1967; Helander, 2014)

$\Delta\psi=\frac{1}{q}\frac{\partial\mathcal{J}_{\|}}{\partial\alpha},$ (1a)
$\Delta\alpha=-\frac{1}{q}\frac{\partial\mathcal{J}_{\|}}{\partial\psi}.$ (1b)

In particular, the particle remains in the vicinity of one flux surface if
$\p\mathcal{J}_{\|}/\p\alpha=0$. Magnetic fields in which this property holds
for all orbits are called omnigenous (Hall & McNamara, 1975; Cary &
Shasharina, 1997). The direction and magnitude with which particles precess
within the surface is governed by $\p\mathcal{J}_{\parallel}/\partial\psi$. It
is convenient to define a precession frequency $\omega_{\alpha}$, by dividing
$\Delta\alpha$ by the bounce time $\Delta t=\int\mathrm{d}\ell/v_{\parallel}$.
Fields in which this quantity is negative for electrons (that is,
$q\omega_{\alpha}>0$) and thus opposed to the diamagnetic drift
($q\omega_{\star}\sim p^{\prime}(\psi)/n<0$), so-called
maximum-$J$-configurations, have long been known to possess favourable
stability properties for trapped-particle modes (Rosenbluth, 1968; Proll et
al., 2012; Helander et al., 2013).

The maximum-$J$ property is additionally correlated with, but not identical to
(Helander, 2014), the existence of a magnetic well (Greene, 1997), which is
beneficial for MHD-stability. This circumstance has to do with its relation to
magnetic curvature. If the magnetic curvature vector is decomposed into
components in the direction of $\nabla\psi$ and $\nabla\alpha$,

$\mbox{\boldmath$\kappa$}={\bf b}\cdot\nabla{\bf
b}=\kappa_{\psi}\nabla\psi+\kappa_{\alpha}\nabla\alpha,$

then a positive value of $\kappa_{\psi}$ is referred to as favourable (or
“good”) curvature. Since the guiding-centre curvature drift is in the
direction ${\bf b}\times\mbox{\boldmath$\kappa$}$ and

$\kappa_{\psi}=\frac{({\bf
b}\times\mbox{\boldmath$\kappa$})\cdot\nabla\alpha}{B},$

good curvature is clearly correlated with the maximum-$J$ property. Indeed,
the latter can be interpreted as the condition that a certain average of the
normal curvature should be favourable for all trapped-particle orbits.

The trapped particles in an omnigenous stellarator can precess in the
poloidal, toroidal or helical directions, and if the precession is poloidal,
the field is called quasi-isodynamic (QI) (Helander & Nührenberg, 2009;
Nührenberg, 2010). This feature is determined by the topology of the contours
of constant $|\mathbf{B}|$. The aim of the present paper is to explore
conditions under which a QI field may possess the maximum-$J$-property. It has
long been known that this can be the case at finite plasma $\beta$ (thermal
pressure divided by magnetic pressure), very recently achieved at low values
(Sánchez et al., 2023; Goodman, 2024) and, as we shall see, it is also
possible (to a very good approximation) when $\beta$ vanishes. This is in
stark contrast to quasisymmetric magnetic fields (including tokamaks), where
the maximum-$J$-property is unattainable, at least in a region near the
magnetic axis. If the maximum-$J$ condition is satisfied for some, but not
all, particle orbits, then ${\partial\mathcal{J}_{\|}}/{\partial\psi}$ must
vanish for certain orbits and these will usually undergo large radial
excursions, resuling in super-banana (Velasco et al., 2021) or banana-drift
convective (Paul et al., 2022) losses. It is clearly preferable that all
orbits satisfy the maximum-$J$ condtion.

Quasi-isodynamic stellarators are sometimes referred to as a system of linked-
mirrors (Boozer, 1998). Following this analogy, the maximum-$J$-condition can
be seen to be closely related to the so-called minimum-$B$ property in
magnetic mirrors, which was proposed in the early 1960’s (Taylor, 1963). To
see why, consider particles trapped in the vicinity of a local minimum,
$B_{\rm min}(\psi)$. Since

$\frac{\partial\mathcal{J}_{\parallel}}{\partial\psi}=-\frac{mv}{2}\int\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,l}\frac{\lambda dl}{\sqrt{1-\lambda B}},$ (2)

it follows that the maximum-$J$-property requires

$\left(\frac{\partial B}{\partial\psi}\right)_{\alpha,l}>0$

at the field minimum. In other words, $B_{\rm min}$ increases with radius and
the magnetic field strength assumes a global minimum on the magnetic axis of a
maximum-$J$ stellarator (Helander, 2014). In the vicinity of this minimum, the
surfaces of constant field strength are ellipsoids aligned with the magnetic
field.

In the following sections, we explore these ideas further and derive
conditions for attaining the maximum-$J$-condition, focusing particularly on
deeply trapped particles close to the magnetic axis.

##  2 Trapped-particle precession

A derivation of the general expression for the precession frequency of trapped
particles in an arbitrary stellarator is presented in Appendix A. The result
is expressed as an integral along the magnetic field between two consecutive
bounce points, which in Boozer coordinates becomes

$\omega_{\alpha}=\frac{mv^{2}}{q}\int\left[\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,\varphi}\frac{1-\lambda
B/2}{B^{2}\sqrt{1-\lambda B}}+\frac{\mu_{0}p^{\prime}(\psi)}{\langle
B^{2}\rangle}\frac{\sqrt{1-\lambda
B}}{B}\right]d\varphi\bigg{/}\int\frac{d\varphi}{B\sqrt{1-\lambda B}}.$ (3)

For simplicity, we have assumed that the net toroidal plasma current enclosed
by the magnetic surface under consideration vanishes, $I(\psi)=0$.

In Equation (3) there are two separate terms contributing to the precession:
the radial derivative of $B=|\mathbf{B}|$ and the pressure gradient
$p^{\prime}$, which is negative for typical pressure profiles.111In the vacuum
limit, the curvature drift and the grad-B drift are equal, and thus only the
gradient of $B$ appears in the expression for $\omega_{\alpha}$, but if
$\beta>0$ this is generally not the case. The difference between the two
drifts is then proportional to the pressure gradient, which is responsible for
the additional pressure term in Eq. (3). However, it would be wrong to
conclude that a non-zero pressure gradient only affects the precession
frequency through the explicit pressure term in Eq. (3). This naïve
interpretation would imply that a negative pressure gradient tends to make
$q\omega_{\alpha}$ more negative and thus oppose the maximum-$J$ condition.
However, the radial derivative of $|\mathbf{B}|$ depends on the pressure
profile in an equilibrium field. If the central plasma $\beta$ is increased,
the magnetic field strength usually drops in the centre of the plasma to hold
the plasma in place, and $\partial_{\psi}B$ thus grows. The net result turns
out, in practice, to make the maximum-$J$ condition much easier to attain at
high plasma beta compared with the case of a vacuum field. It is thus
important to address the question of how $\partial_{\psi}B$ depends on the
pressure, which we shall consider later.

In general, $\omega_{\alpha}$ is a function, for each trapping well, of the
pitch-angle parameter $\lambda$, the flux-surface label $\psi$, and the field-
line label $\alpha$. The dependence on $\alpha$ is present only when the field
is not omnigeneous. Then the second adiabatic invariant depends on $\alpha$,
and, as a result, the precession frequency of trapped particles varies as they
drift from field line to field line. Numerical examples are shown in Figure 1,
where the normalised precession frequency

$\hat{\omega}_{\alpha}=\frac{4q\psi_{a}\omega_{\alpha}}{mv^{2}}$ (4)

is displayed for three different numerically optimised QI fields (Goodman et
al., 2023). Here $\psi_{a}$ denotes the value of $\psi$ at the plasma boundary
and $\hat{\omega}_{\alpha}$ has been plotted as a function of the trapping
parameter (Roach et al., 1995),

$k^{2}=\frac{B_{\rm min}^{-1}-\lambda}{B_{\rm min}^{-1}-B_{\rm max}^{-1}}.$
(5)

As discussed below, the magnitude of $\hat{\omega}_{\alpha}$ can be
interpreted as the precession frequency (for thermal particles) divided by the
diamagnetic frequency (within a factor of order unity). The grey lines
correspond to field lines on the same flux surface with different values of
$\alpha$. The fact that these lines are different from each other is an
indication of departure from exact omnigeneity. For the most part of the
discussion to follow, we shall assume, unless otherwise stated, that the field
is omnigeneous. We may then choose any field line to compute the bounce-
average (2) of $\partial_{\psi}B$, or alternatively its Boozer form in (3). In
a vacuum field with $p^{\prime}(\psi)=0$, these integrals are weighted
averages of $\partial_{\psi}B$, where the weight is positive for all trapped
particles (i.e., all $\lambda$). It thus follows that the so-called
minimum-$B$ condition, $\partial_{\psi}B>0$, is a sufficient condition for
maximum-$J$ behaviour.

Figure 1: Example of precession in highly quasi-isodynamic fields. Normalised
precession frequency $\hat{\omega}_{\alpha}$, Eq. (4), as a function of the
trapping parameter $k^{2}$, Eq. (5), for three recently optimised quasi-
isodynamic configurations with different number of field periods $N$ (Goodman
et al., 2023), in which omnigeneity is attained to a high degree. For each
configuration, $\hat{\omega}_{\alpha}$ is displayed at two values of the
normalised radial flux coordinate $\varrho$ on a number of different field
lines (gray curves), with the black curve representing $\alpha=0$. A positive
value of $\hat{\omega}_{\alpha}$ corresponds to the maximum-$J$-property. Note
that deeply trapped particles ($k^{2}\rightarrow 0$) tend to behave worse than
barely trapped ones, and the intermediate population can exhibit non-monotonic
dependence on $k^{2}$.

###  2.1 Necessary conditions: special points

It is however not necessary to make $\partial_{\psi}B>0$ everywhere in order
to attain the maximum-$J$ condition. Even if $\partial_{\psi}B<0$ at certain
positions along the field line, the average precession can still be
favourable. The condition $\partial_{\psi}B>0$ must nevertheless be satisfied
at all local minima and maxima of $|\mathbf{B}|$ along the field line. At
these points, deeply and barely trapped particles spend nearly all their time,
other points along the orbit not contributing much to the average. It is
therefore necessary, both at the bottom and top of every magnetic well, that
$\partial_{\psi}B|_{\alpha,\varphi}>0$, even if $p^{\prime}\neq 0$.

This result may be formally derived as the appropriate limits of Eq. (3), but
a simpler argument has already been given in Section 1: since the time
required for a trapped particle to move from one bounce point to the next is
equal to

$\Delta t=\int\frac{dl}{v_{\|}}=\int\frac{dl}{v\sqrt{1-\lambda B}},$

the precession frequency for deeply trapped particles, for which $1-\lambda
B\ll 1$, becomes

$\omega_{\alpha}=\frac{\Delta\alpha}{\Delta
t}=\frac{mv^{2}}{2qB}\left(\frac{\partial B}{\partial\psi}\right)_{\alpha,l}.$
(6)

Since $\partial B/\partial\ell=0$ at the point of minimum field strength, the
derivative $\partial B/\partial\psi$ can equally well be computed at constant
Boozer angles $(\theta,\varphi)$ for an omnigeneous stellarator. A similar
argument holds for barely trapped particles. A more detailed discussion can be
found in Appendix A, including a proof that the explicit pressure term does
not contribute, Eq. (38).

Thus, if the minimum-$B$-condition is satisfied at all local maxima and minima
of $B$ along the field line, i.e. if

$B^{\prime}_{\mathrm{max}}(\psi)>0\qquad\mbox{and}\qquad
B^{\prime}_{\mathrm{min}}(\psi)>0,$ (7)

then deeply and barely trapped particles are guaranteed to precess in the
desired direction ($q\omega_{\alpha}>0$). It would however be wrong to infer
that all particles trapped at intermediate depths in the magnetic wells then
precess in the same direction. Such a conclusion cannot be drawn by solely
considering the radial gradient of $B$ at its extrema along the field. The
conditions in Eq. (7) are necessary for maximum-$J$, but not sufficient.
Enforcing them as part of an optimisation effort may be helpful, but does not
guarantee the maximum-$J$ behaviour everywhere.

This inference is however possible under additional assumptions about
$|\mathbf{B}|$. Obviously, if $\partial_{\psi}B$ were bounded from below by
either of the local $|\mathbf{B}|$ gradients, then Eq. (7) would become a
necessary and sufficient criterion. One scenario in which this holds true is
when the shape of the magnetic well along field lines behaves ‘rigidly’. That
is, in going from one surface to another,
$|\mathbf{B}|=B_{0}(\psi)+\eta(\psi)f(\ell)$ along the field line. The case of
a simple sinusoidal well in QI stellarators has recently been analysed and
proposed as a potentially useful optimisation target (Velasco et al.,
2023).222The flat-mirror’ criterion presented in Velasco et al. (2023)
guarantees Eq. (7) regardless of the model assumed for $|\mathbf{B}|$. In the
case $|B|=B_{00}+B_{M}\cos\varphi$ considered by the authors, this criterion
implies maximum-$J$ behaviour of all trapped particles, but not in general. In
general, though, the behaviour of $\partial_{\psi}B$, and thus
$\omega_{\alpha}$, can be rather complicated. Figure 1 shows a recently
published example of omnigeneity-optimised QI stellarators, where the
behaviour at $k=(0,1)$ (deeply and barely trapped particles respectively) are
not always good indicators for intermediate particles. Nevertheless, given the
simplicity and necessary nature of the condition in Eq. (7), we shall use it
as a representative feature of the maximum-$J$-condition, and extend it to
more global considerations.

###  2.2 Quasi-isodynamic fields

The treatment so far has considered little information about the magnetic
field other than properties that follow directly from the requirements of the
MHD equilibrium and omnigeneity. In the interest of describing the particulars
of QI fields, though, we introduce to the discussion features unique to this
class of stellarators. To do so, and to further simplify the discussion, we
focus our attention to the vicinity of the magnetic axis.

By expansion in the minor-radius coordinate $r=\sqrt{2\psi/\bar{B}}$, where
$\bar{B}$ is some reference magnetic-field strength, which we choose to be
equal to that at the minimum, to second order the magnetic field magnitude may
be written in the following form,

$B=B_{0}(\varphi)+rB_{1}(\alpha,\varphi)+r^{2}B_{2}(\alpha,\varphi),$

where, as shown in Appendix A, the coefficients satisfy the following relation
in an exactly QI stellarator-symmetric field (Rodríguez & Plunk, 2023):

$\displaystyle B_{1}(\alpha,\varphi)=-d(\varphi)\sin\alpha,$ (8a)
$\displaystyle B_{2}(\alpha,\varphi)=B_{20}(\varphi)-B_{2s}(\varphi)\sin
2\alpha-\frac{\partial}{\partial\varphi}\left(\frac{B_{0}^{2}d^{2}}{4B_{0}^{\prime}}\right)\cos
2\alpha.$ (8b)

Here $B_{0}(\varphi)$ and $B_{20}(\varphi)$ are even in $\varphi$ whereas
$d(\varphi)$ and $B_{2s}(\varphi)$ are odd, if we choose the angle $\varphi$
to vanish at the point where $B_{0}(\varphi)$ attains its minimum. As a
result, the first-order term $B_{1}$ vanishes at the minimum, and the
precession frequency Eq. (6) for the most deeply trapped particles, being
proportional to $\partial B/\partial\psi$, is determined by the second-order
term $B_{2}$. It becomes

$\omega_{\alpha}=\frac{mv^{2}}{qB_{0}^{2}}\left[B_{20}-\left(\frac{B_{0}^{2}d^{2}}{4B_{0}^{\prime}}\right)^{\prime}\right]+\omega_{\alpha}^{\mathrm{non-
QI}},$ (9)

where everything is evaluated at the minimum, $\varphi=0$, the primes denote a
derivative with respect to the toroidal angle $\varphi$, and
$\omega_{\alpha}^{\mathrm{non-QI}}$ represents the contribution from
deviations from omnigeneity, which will be discussed later. To achieve
$q\omega_{\alpha}>0$, the magnetic field must be carefully tailored in such a
way that for deeply trapped particles the quantity in the square brackets of
Eq. (9) is positive. As we shall see, this is eminently possible but requires
careful consideration of the terms within the square brackets.

As shown explicitly in Appendix B, Eq. (52), departures from the condition of
omnigeneity lead to additional terms in the expression for precession, here
denoted by $\omega_{\alpha}^{\mathrm{non-QI}}$. These contributions may be
interpreted as non-intrinsic contributions due to departures from QI. Breaking
omnigeneity at first order in $r$ leads to a contribution that scales as
$\omega_{\alpha,-1}^{\mathrm{non-QI}}\propto\cos\alpha/r$, Eq. (51), while
deviations at second-order lead to $\omega_{\alpha,0}^{\mathrm{non-
QI}}\propto\cos 2\alpha$, Eq. (52). Because of the dependence on $\alpha$,
there always exists a field line on which the non-omnigeneous contribution is
$q\omega_{\alpha}^{\mathrm{non-QI}}<0$, i.e., detrimental to the maximum-$J$
condition. The deviations that arise at first order are particularly worrying
near the magnetic axis due to the $1/r$ scaling. We explain the origin of this
behaviour in the following subsection. From the forms above we may
nevertheless conclude that omnigeneity is necessary for the maximum-$J$
property.

###  2.3 Comparison with quasisymmetric fields

Trapped-particle precession in a QI field is markedly different from that in
quasisymmetric (QS) configurations, the other important class of optimised
stellarators. Particle orbits in the latter are similar to those in tokamaks,
and the maximum-$J$ property cannot be satisfied throughout the volume for all
particle classes. In fact, close to the magnetic axis, the magnetic field
strength in any quasisymmetric field (characterised by the helicity of the
symmetry $N$) is

$B(r,\chi=\theta-N\varphi)=B_{0}(1-r\eta\cos\chi)$ (10)

and the precession frequency of trapped particles becomes (Rodríguez &
Mackenbach, 2023; Kadomtsev & Pogutse, 1967; Helander & Sigmar, 2005)

$\omega_{\alpha}^{\mathrm{QS}}\approx-\frac{mv^{2}}{2q}\frac{\eta}{rB_{0}}\left(2\frac{E(k)}{K(k)}-1\right),$
(11)

where $E$ and $K$ denote elliptic integrals of the first and second kind
(Olver et al., 2020, Sec. 19). Their argument is defined by
$k=\sin(\chi_{b}/2)$, where $\chi_{b}$ is the bounce point in the well, so
that $k=0$ refers to deeply trapped particles and $k=1$ to the trapped-passing
boundary. This definition matches Eq. (5) to leading order. We present a plot
of the behaviour in quasisymmetric configurations in Figure 2.

Figure 2: Precession in zero-$\beta$ tokamaks and quasisymmetric stellarators.
The left plot shows the predicted dependence of the leading-order precession
frequency in a tokamak or quasisymmetric stellarator as a function of the
trapping parameter, $k$, Eq. (11). The right plots show the precession
frequency on the boundary of actual QS vacuum configurations (Landreman &
Paul, 2022), normalised as in Eq. (4) (see Rodríguez & Mackenbach (2023) for
more details and discussion). As in Figure 1, the plots on the right show the
variation with field line label $\alpha$, which is difficult to discern due to
the high degree of quasisymmetry in these configurations.

To make the comparison to the QI case more explicit, we turn to the general
expression for the precession frequency close to the magnetic axis given in
Eq. (47) and consider the limit in which the magnetic-field strength only
varies slightly along the magnetic axis. Then trapped particles have
$1-\lambda B\ll 1$ and Eq. (47) reduces to

$\omega_{\alpha}=\frac{2mv^{2}}{q\bar{B}}\int\frac{f(\varphi)d\varphi}{\sqrt{1-\lambda
B_{0}}}\bigg{/}\int\frac{d\varphi}{\sqrt{1-\lambda B_{0}}},$ (12)

where

$f(\varphi)=\frac{B_{20}(\varphi)}{B_{0}}-\frac{1}{B_{0}}\frac{d}{d\varphi}\left(\frac{B_{0}^{2}d^{2}}{4B_{0}^{\prime}}\right).$
(13)

The precession frequency is thus proportional to a bounce-average of the
function $f(\varphi)$ between the two turning points. This expression reveals
several important differences between the precession frequency in QI and QS
stellarators (including tokamaks):

  1. 1.

In a QS field the particle precession becomes ‘infinite’ ($\sim 1/r$) as the
magnetic axis is approached, whereas it remains finite in QI fields. This
$1/r$ behaviour comes from a finite poloidal component of the curvature drift
as the axis is approached. In a QI configuration, the curvature of the
magnetic field vanishes at the minimum on the axis, eliminating this
behaviour. In the QS case, deeply trapped particles reside in the bad-
curvature region, and the normal curvature is even about the minimum of $B$
(see Figure 10). This difference in behaviour follows directly from the
difference in the topology of $|\mathbf{B}|$ contours, and its implications on
the order-$r$ correction to $B_{0}$, Eqs. (8) and (10). Formally, this
difference in parity explains the leading-order cancellation of the precession
in QI fields and the elimination of the $1/r$ contribution. The cancellation
ceases to be exact, though, whenever the field deviates from omnigeneity. In
that event, there will always exist some $r_{a}$ such that, for $r<r_{a}$, the
field ceases to satisfy the maximum-$J$ condition (see Fig. 6 for example). It
is known that omnigeneity must be broken at first order near the tops of the
magnetic well, which unavoidably leads to a small, but finite, $r_{a}$ below
which trapped particles near the magnetic well tops are not maximum-$J$. We
discuss this issue further in Appendix A.4, but shall otherwise make the
assumption of exactly omnigenity, so that this term may be neglected.

Figure 3: Curvature about the minmum in a QS and QI field. The plots show a
comparison of the curvature ($\mathbf{B}\times\kappa\cdot\nabla\alpha$) in (a)
a QA (the precise QA in Landreman & Paul (2022)) and (b) a QI (the $N=2$
configuration in Goodman et al. (2023)) configurations. The left plots show a
3D rendition of the boundary of the configurations, in which the colours
display the quantity $\mathbf{B}\times\kappa\cdot\nabla\alpha$, with blue/red
representing negative/positive values (good/bad curvatures) respectively. The
plots to the right show the magnetic field $|\mathbf{B}|$ and the curvature
along the field line, using the cylindrical coordinate $\phi$ and the PEST
poloidal angle $\vartheta$ as coordinates along the fieldline. The scatter
plot indicates the position of the minimum. These two examples illustrate the
qualitative difference in the curvature parity between a QI and a QS field.

  2. 2.

In QS stellarators, the fraction of trapped particles decreases as one
approaches the magnetic axis, as it is the first-order poloidal variation of
the field strength that defines the trapping well. In QI stellarators, this is
not the case and there remains a finite trapped population on axis. This
difference is important close to the axis, where there is a region with so-
called potato-orbits (Helander & Sigmar, 2005; Rodríguez & Mackenbach, 2023)
in QS stellarators where the thin-orbit approximation to
$\mathcal{J}_{\parallel}$ breaks down.

  3. 3.

In QS, it is always the case that there are trapped particles precessing in
opposite directions. Deeply trapped particles precess like in a minimum-$J$
field ($q\omega_{\alpha}<0$), while barely trapped ones do so in the opposite
way. This is generally not the case in QI, where the function $f$, Eq. (13),
has more freedom than in the QS case. This stems from differences in the
distribution of good and bad curvature about the minimum of $|\mathbf{B}|$.

These conclusions hold rigorously on the magnetic axis and, by continuity,
also in its vicinity. The differences in the behaviour of curvature can be
seen to hold approximately true in practice, as is shown through two examples
in Figure 3. In addition, in the QS case it has been shown that the near-axis
description remains instructive as a model of QS configurations beyond its
asymptotic regime (see Figure 2).

The most important observation is that, in the context of QI configurations,
maximum-$J$ behaviour appears possible even in a vacuum field. Within the
near-axis framework, whether this is the case depends on the sign of the
bounce-average of the function $f(\varphi)$, which is analysed below.

###  2.4 Shallowly and deeply trapped particles in QI

The special character of the bottom and top of the magnetic well provides us
with a simple way to assess the possibility of maximum-$J$ behaviour in QI
stellarators.

Satisfying the maximum-$J$-condition for the most shallowly trapped particles
tends to be easy in practice. These particles spend most of their time close
to the turning point where the field strength reaches its maximum along the
field line, $B(\psi,\alpha,\varphi)=B_{\rm max}(\psi)$, which is independent
of the field-line label $\alpha$ in a QI stellarator. Moreover, the maximum is
located at a constant value of the toroidal Boozer angle,
$\varphi=\varphi_{\rm max}$, which is independent of $\alpha$ (Cary &
Shasharina, 1997; Landreman & Catto, 2012; Helander, 2014). At the maximum,
the Taylor expansion of the field strength in the coordinates $\alpha$ and
$\varphi$ is thus of the form

$B(\psi,\alpha,\varphi)=B_{\rm
max}(\psi)+\frac{B_{\varphi\varphi}}{2}(\varphi-\varphi_{\rm
max})^{2}+\cdots,$ (14)

since, at the toroidal angle of the maximum, the following derivatives all
vanish:

$\frac{\partial B}{\partial\alpha}=\frac{\partial
B}{\partial\varphi}=\frac{\partial^{2}B}{\partial\alpha\partial\varphi}=\frac{\partial^{2}B}{\partial\alpha^{2}}=0.$
(15)

Therefore
$\nabla^{2}B=\nabla^{2}B_{\mathrm{max}}+B_{\varphi\varphi}|\nabla\varphi|^{2}+\dots$,
where the second term is negative at the top of the well. Now, in a vacuum
field (zero plasma pressure), the quantity
$B^{2}=B_{x}^{2}+B_{y}^{2}+B_{z}^{2}$ is a subharmonic function since

$\nabla^{2}B^{2}=\nabla\cdot\nabla(B_{x}^{2}+B_{y}^{2}+B_{z}^{2})=2(|\nabla
B_{x}|^{2}+|\nabla B_{y}|^{2}+|\nabla B_{z}|^{2})\geq 0,$

where we have used $\nabla^{2}B_{x}=0$ etc. for each Cartesian component of
$\bf B$ (Solov’ev & Shafranov, 1970). Because of the maximum principle for
subharmonic functions, the maximum of $B^{2}$ over any closed, bounded domain
cannot be attained in its interior unless $B^{2}$ is constant (Evans, 2022,
Thrm. 1, Ch. 6.4.1). It thus follows that the function $B_{\rm max}(\psi)$ in
Eq. (14) cannot have a local maximum at any value of $\psi$. Therefore,
$B_{\rm max}^{\prime}(\psi)\geq 0$ for all $\psi$ and, as a result of Eq. (6),
shallowly trapped particles must satisfy the maximum-$J$-condition,
$q\omega_{\alpha}>0$, throughout the plasma.333Of course, $B^{\prime}(\psi)$
could also vanish without $B(\psi)$ acquiring a local maximum, but we do not
consider such singular cases that are not robust to arbitrarily small
perturbutions. It should be clear that this argument holds for any omnigenous
magnetic configutation, including a quasisymmetric stellarator.

The situation is very different at the opposite end of the trapped population.
The most deeply trapped particles at the bottom of the well tend to precess in
the unfavourable direction. Indeed, these particles frequently (though not
always, see Figure 1), have the smallest value of $q\omega_{\alpha}$. The
relevant question here is then whether this tendency be reversed to make the
system acquire the maximum-$J$ property at the bottom of the trapping well.

In order to understand the maximum-$J$ condition on deeply trapped particles,
we resort to Eq. (9), which provides an explicit expression for
$\omega_{\alpha}$ in terms of $B_{20}$, $B_{0}$ and $d$. The expressions at
this point may be further simplified by realising that the curvature of the
magnetic axis vanishes at the point of minimum $|\mathbf{B}|$, i.e.
$\kappa(\varphi)=0$ when $\varphi=0$. As a result, the behaviour of
$|\mathbf{B}|$ is significantly constrained. The field around this point thus
resembles a straight magnetic mirror, and thus we can gain some intuition
about QI fields by analysing a straight magnetic mirror, which is done in
Appendix C.

With these remarks in mind, we proceed to consider the governing set of
equations in the near-axis framework. Formally, we explore the properties of a
magnetic field in MHD equilibrium that satisfies the solenoidal condition, and
possess flux surfaces at $\varphi=0$. The analysis is presented in Appendix B,
assuming perfect QI at first order for simplicity, and leads to an explicit
form for the average radial derivative of the magnetic pressure ($B_{20}$),
which is primary ingredient of Eq. (9), explicitly constructed in Eqs. (59) in
Appendix B. When the condition of omnigeneity at second order is also imposed
(that is, the choice of functions is made so that QI is satisfied at the
bottom of the well at second order), the function $f$ in Eq. (12) can be
written in terms in the following form

$f_{0}=f_{p}+f_{B_{0}^{\prime\prime}}+f_{I_{2}}+f_{\tau_{0}^{2}}+f_{\mathrm{QI}},$

with

$\displaystyle f_{p}=-\frac{\mu_{0}p_{2}}{B_{0}^{2}},$ (16a) $\displaystyle
f_{B_{0}^{\prime\prime}}=-\frac{1}{2(l^{\prime})^{2}}\frac{\sqrt{\bar{\alpha}}}{\bar{\alpha}+1}\frac{B_{0}^{\prime\prime}}{B_{0}},$
(16b) $\displaystyle
f_{\tau_{0}^{2}}=\frac{\sqrt{\bar{\alpha}}}{1+\bar{\alpha}}\tau_{0}^{2},$
(16c) $\displaystyle
f_{I_{2}}=-\frac{2\sqrt{\bar{\alpha}}}{1+\bar{\alpha}}\frac{I_{2}}{B_{0}}\tau_{0},$
(16d) $\displaystyle
f_{\mathrm{QI}}=-\frac{1}{4B_{0}}\frac{2}{1+\bar{\alpha}}\left.\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)^{\prime}\right|_{\varphi=0},$
(16e)

where $\bar{d}=d/\kappa$ and $\bar{\alpha}=\bar{d}^{4}B_{0}^{2}/\bar{B}^{2}$
is directly related to the elongation (along the direction of curvature) of
the cross-section at the toroidal position $\varphi=0$. There are thus five
different terms contributing to $f$, which we proceed to analyse individually.

####  2.4.1 Role of pressure

As expected, increasing the pressure gradient supported by the field,
$|p_{2}|=(\bar{B}/2)|\mathrm{d}p/\mathrm{d}\psi|$, increases $f_{0}$. That is,
it makes deeply trapped particles more likely to satisfy the maximum-$J$
condition $q\omega_{\alpha}>0$. This is the well-known effect of finite
$\beta$ improving the maximum-$J$ property of QI fields, which is precisely
the opposite to the naïve interpretation of the role of the pressure from the
explicit term in Eq. (3), which vanishes at the extremal points. This effect
of $p^{\prime}$ is the same at both the bottom and top of the well, thus
bringing both deeply and barely trapped particles towards
$q\omega_{\alpha}>0$. To a large extent, then, we would expect the rest of the
trapped population to do likewise, and although we cannot prove it, it is
convenient to think of this effect as an overall upshift of
$q\omega_{\alpha}$. This diamagnetic behaviour turns out to be the same in the
axisymmetric/quasisymmetric limit (Rodríguez & Mackenbach, 2023, Eq.
(3.6a))(Rosenbluth & Sloan, 1971).

Whether a non-zero pressure gradient is needed to ensure $f_{0}>0$ depends on
the relative size of the other terms in Eqs. (16). At a minimum, we need to
overcome the detrimental effect of being located at the minimum of the
magnetic field, $f_{B_{0}^{\prime\prime}}$. The plasma beta,
$\beta_{0}=2\mu_{0}p_{0}/B_{0}^{2}$, that neutralises this term,
$|f_{p}|=|f_{B_{0}^{\prime\prime}}|$, we define to be $\beta_{J}$ and provides
an estimate of the critical $\beta$ that leads to maximum-$J$ behaviour. If we
consider a simple quadratic pressure profile for a field with minor radius
(defined as in the near-axis) $r=a$, then
$\mathrm{d}p/\mathrm{d}\psi\approx-2p_{0}/B_{0}a^{2}$, and thus $f_{0}>0$ when

$\beta_{0}>\beta_{J}\equiv
R_{M}\frac{\sqrt{\bar{\alpha}}}{1+\bar{\alpha}}\left(\frac{a}{L_{B}}\right)^{2},$
(17)

and $R_{M}/L_{B}^{2}=\partial_{\ell\ell}\ln B_{0}$,
$R_{M}=(B_{\mathrm{max}}-B_{\mathrm{min}})/B_{\mathrm{min}}$ is the mirror
ratio and $L_{B}$ the characteristic parallel length-scale of the bottom of
the well. For typical field parameters (a moderately elongated cross-section
with $\bar{\alpha}\sim 1$, a well of the scale of the major radius, and
$R_{M}\sim 0.2$), we get critical plasma betas around (or below) a percent.
This suggests that maximum-$J$ behaviour is readily achievable at a moderate
plasma $\beta$.

####  2.4.2 Difficulties at the minimum

The special character of the bottom of the well is captured by the
$f_{B_{0}^{\prime\prime}}$ term, which is negative (since
$B_{0}^{\prime\prime}>0$) and thus opposes maximum-$J$ behaviour. The
situation is however reversed at the maximum of $B$, where
$B_{0}^{\prime\prime}<0$ and $f_{B_{0}^{\prime\prime}}>0$, implying an
‘intrinsic’ tendency for maximum-$J$ behaviour at this location. In fact, in a
vacuum field, $f_{0}>0$, in agreement with the general proof presented in the
previous sections.

The detrimental contribution of $f_{B_{0}^{\prime\prime}}$ at the bottom of
the well can be mitigated by reducing $B_{0}^{\prime\prime}$ to a minimum,
i.e. by flattening the bottom of the well (increasing $L_{B}$) or reducing the
mirror ratio $R_{M}$. The elongation of the flux-surfaces may also be tweaked
to reduce $f_{B_{0}^{\prime\prime}}$, which is maximal for circular cross-
sections ($\bar{\alpha}=1$). However, note that the other geometric
contributions to $f_{0}$, Eqs. (16) also depend on the elongation, and thus in
relative terms this shaping may not be as effective although it does reduce
$\beta_{J}$.

Figure 4: Schematic depiction of trapping wells that favour maximum-$J$
behavior. Certain magnetic field features make them more prone to maximum-$J$
behaviour than others. In particular, wide and flat magnetic trapping wells
(left) at straighter sections (right) of the stellarator favour maximum-$J$.
The diagram on the right is a schematic top-down view on an $N=2$ stellarator,
with the shading denoting $|\mathbf{B}|$.

Flattening the bottom of the well comes with potential drawbacks in the form
of sensitivity to deviations from omnigeneity. Having an extended region of
small $\partial B_{0}/\partial l$ makes the configuration more susceptible to
error fields, as secondary shallow trapping wells may be created by small
perturbations. Such wells can be seen in the near-axis description of
$|\mathbf{B}|$, where the requirement of omnigeneity (and more particularly,
that of poloidally closed $B$-contours) limits the behaviour of $B_{1}$
depending on how flat $B_{0}$ is. As analysed in detail in Rodríguez & Plunk
(2023), for a field near the minimum described by $B_{1}\sim\varphi^{v}$ and
$B_{0}^{\prime}\sim\varphi^{u-1}$, if $B_{0}$ is too flat, $u\geq 2v$ (except
$u=2,~{}v=1$), this will introduce defects in $|\mathbf{B}|$ that lead to
losses of deeply trapped particles. In fields where this situation is avoided,
that is, for $u<2v$, the quantity $f_{\mathrm{QI}}$ vanishes, Eq. (16e). Only
in the special case of a first-order curvature zero ($v=1$) and a quadratic
well ($u=2$) is this contribution finite. Because $d^{2}/B_{0}^{\prime}$
vanishes at $\varphi=0$ and is positive for $\varphi>0$, its derivative must
be greater than or equal to zero. Hence, $f_{\mathrm{QI}}\leq 0$ at the bottom
of the well; that is, its contribution is detrimental. This makes the
‘standard’ first-order curvature zero and quadratic-well field (Camacho Mata
et al., 2022) particularly unfavourable for maximum-$J$ behaviour. Conversely,
making the section of the field where minimum $|\mathbf{B}|$ is located as
straight as possible (larger $u$ and $v$) should be beneficial (see Fig. 4).

####  2.4.3 Role of local shear, twist and elongation

The contribution of the torsion of the axis to $f_{0}$ is always beneficial
since $f_{\tau_{0}^{2}}>0$. In fact, the larger the torsion, the larger
$f_{0}$, and thus, the closer the behaviour of deeply trapped particles will
be to the maximum-$J$ requirement. The role of torsion may be surprising, but
is understandable from the perspective of a straight magnetic mirror. At the
bottom of a straight magnetic mirror, attaining maximum-$J$ is only possible
if magnetic field lines are locally twisted (see Appendix C and Eqs.
(86a)-(86b)), i.e., if they experience some form of left-right asymmetry, and
in this way possess non-zero local magnetic shear. In the context of our QI
stellarator, the appearance of the $\tau_{0}^{2}$ term is simply a statement
of the necessity of this twist about the magnetic axis, which is the geometric
meaning of torsion.

The comparison with a straight magnetic mirror is also helpful for
understanding the role of flux-surface elongation, which in a mirror needs to
increase away from the bottom of the magnetic well (see Figure 11). In the
context of a stellarator-symmetric, QI stellarator, the elongation of the
cross-section at the minimum is given by $\bar{d}^{2}$, as previously
mentioned. Thus, a change in elongation would be expected to involve a term
proportional to $\bar{d}^{\prime\prime}$. 444Identifying
$\bar{d}^{\prime\prime}$ with the change in elongation at the minimum is not
quite correct, as the elongation can vary through other means (Landreman &
Sengupta, 2018; Rodríguez, 2023) even when $\bar{d}$ does not change. In fact,
one may show that at the minimum the change in elongation is affected by
$\bar{d}^{\prime\prime}$, $B_{0}^{\prime\prime}$ and $\tau^{2}$. It can be
rigorously shown that in the omnigeneous case increasing the torsion and
$B_{0}^{\prime\prime}$ (when $\bar{\alpha}\neq 1$ for the latter) always
increase the elongation. Anyhow, the interpretation of
$\bar{d}^{\prime\prime}$ as an addedd change to elongation is correct and
illustrating. There is no term in $f_{0}$ which involves
$\bar{d}^{\prime\prime}$, but this is a result of the condition of omnigeneity
at second order. In fact, upon relaxing the latter, $B_{20}$ at the bottom of
the well does depend on $\bar{d}^{\prime\prime}$, Eq. (60), and it is only
through the omnigeneity condition of Eq. (64) that this explicit dependence
can be eliminated. From the omnigeneity condition it follows that
$\bar{d}\bar{d}^{\prime\prime}(1+1/\bar{\alpha})=-(\tau_{0}\bar{d}l^{\prime})^{2}(3+1/\bar{\alpha})+\dots$,
meaning that increasing the torsion to favour $f_{\tau_{0}^{2}}$, requires one
to modify $\bar{d}^{\prime\prime}$ accordingly. The plasma cross-section must
become elongated in the binormal direction away from the bottom of the well,
which is a consequence purely of the omnigeneity condition. Looking at the
contribution of the $\bar{d}^{\prime\prime}$ term to the equilibrium equation
of $B_{20}$, Eq. (60), we see that
$\bar{d}\bar{d}^{\prime\prime}(1-1/\bar{\alpha})>0$ promotes maximum-$J$
behaviour there. Shaping of elongation is once again vital, now for attaining
maximum-$J$. For this shaping to be synergistic between the omnigeneity
condition and maximum-$J$, we need a binormally elongated cross-section
$\bar{\alpha}<1$ (i.e., binormal elongation at the minimum, $|\bar{d}|<1$). In
practice, we seek shapes like the exagerated schematic in Figure 5, features
that appear to be common in many QI configurations (Camacho Mata et al., 2022;
Jorge et al., 2022).

Figure 5: Prototypical shape of flux surfaces of an omnigeneous, maximum-$J$
field near the point of minimum field strength. 3D illustration of features
generally expected from maximum-$J$, QI magnetic fields. Note that field lines
twist and that the elongation grows with increasing distance from the field-
strength minimum.

Given the central role played by the local field-line twist, a non-zero plasma
current can be either beneficial or detrimental depending on its alignment
with the torsion of the axis. As can be inferred from its contribution to the
rotational transform (C. Mercier & N. Luc, 1974), a negative current acts
constructively with a positive torsion (and vice-versa). The term $f_{I_{2}}$
reflects precisely this fact. In QI configurations, toroidal currents tend to
be small and this contribution may then be disregarded.

All in all, we learn from this analysis that in a vacuum magnetic field there
is only one way of possibly attaining maximum-$J$ behaviour for deeply-trapped
particles, which is to have a large torsion at the point of minimum $B$, in
the sense that $\tau_{0}L_{B}>\sqrt{R_{M}/2}$, and thereby also a strong
growth of binormal elongation. In other words, the binormal vector must rotate
significantly within the magnetic well. From this analysis, it follows that it
is possible to make deeply trapped particles precess in the favourable
maximum-$J$ direction, unlike the case in quasi-symmetric or tokamak fields.
This result is robust even if we acknowledge the impossibility of exact QI at
first order in the near-axis expansion, because the violation of omnigeneity
is only necessary near the maximum of $|\mathbf{B}|$.

The construction above poses some important practical difficulties. First of
all, the field will generally tend to develop large shaping, both because the
maximum-$J$ criterion requires significant torsion and omnigeneity also an
increase in elongation (see the depiction in Figure 5). Comparing the various
contributions to $f_{0}$ in Eq. (16) reveals that the shaping needs to be
specific and substantial to endow the deeply trapped particles with the
maximum-$J$ property. If the shaping is constrained by practical
considerations, as is usually the case in stellarator optimisation studies,
the opportunity for attaining the maximum-$J$ property is correspondingly
limited. How much shaping is sufficient varies from case to case, and the
competition between the various field and geometric quantities (we shall
recall that the shaping arguments above are a simplified local view near the
minimum). The second thing to bear in mind is that our considerations here are
limited to deeply (and barely) trapped particles. Intermediate orbits could
behave differently, which is something that needs to be checked by computing
$\omega_{\alpha}$ for all $k$.

####  2.4.4 Beyond stellarator-symmetry

The realisation that a certain level of asymmetry about the bottom of the
magnetic well provides the means to enhance the precession of deeply trapped
particles opens the door to several possibilities. In the case of a straight
mirror, one needs to break left-right symmetry in order to achieve an
omnigeneous, maximum-$J$ field (see Appendix C and Catto & Hazeltine (1981)).
In the case of a QI stellarator it is then natural to ask the question whether
breaking stellarator symmetry can be exploited to further improve the
behaviour of the deeply trapped particles.

The procedure followed for the discussion above can be extended to the non-
stellarator-symmetric case. The details of the derivation are presented in
Appendix B. As a result of this extension, the expression for $f_{0}$ acquires
a number of additional terms, and the stellarator-symmetric contributions are
also modified. The additions are however limited if we specialise, for
simplicity, to the case in which the elliptical cross-section at the bottom of
the magnetic well remains up-down symmetric in the Frenet-Serret frame, i.e.,
if we choose $\sigma(\varphi=0)=0$. In that case, there is only a single
additional term that arises from stellarator symmetry, which is proportional
to $(\bar{d}^{\prime})^{2}$,

$f_{(\bar{d}^{\prime})^{2}}=\frac{(\bar{d}^{\prime}/\ell^{\prime})^{2}}{1+\bar{\alpha}}.$
(18)

It is clear that $f_{(\bar{d}^{\prime})^{2}}\geq 0$, and thus one can exploit
stellarator-symmetry breaking to improve the maximum-$J$ behaviour of deeply
trapped particles. In this particular case, $f_{(\bar{d}^{\prime})^{2}}>0$ and
thus any amount of symmetry breaking will help the precession of deeply
trapped particles. When the level of asymmetry is such that it changes
$\bar{d}$ significantly within the well, $L_{\bar{d}}<L_{B_{0}}$, where
$L_{\bar{d}}^{-1}\sim\partial_{\ell}\ln\bar{d}$, this is capable of
overturning the unfavourable contribution from $B_{0}^{\prime\prime}$ at the
bottom of the well. This may be regarded as evidence that breaking stellarator
symmetry might be beneficial.

##  3 Some examples

We showed in the previous section that it is possible to achieve maximum-$J$
behaviour for deeply trapped particles in QI stellarators but not in
quasisymmetric ones. To attain this goal, the shaping of the magnetic field
must however be carefully tailored and, in particular, significant shaping of
the magnetic field is necessary to promote the correct precession behaviour.
In this section, we construct concrete examples of QI fields and analyse their
maximum-$J$ properties.

###  3.1 Near-axis constructions

We begin by considering configurations found through the near-axis expansion.
To this end, we use the QI-specific developments of Plunk et al. (2019) and
Rodríguez & Plunk (2023) as well as the general equilibrium framework of
Landreman & Sengupta (2019). To diagnose the resulting fields, we use the
normalised precession frequency $\hat{\omega}_{\alpha}$, Eq. (4), where we
express $\psi_{a}$, the value of the flux at the boundary, in terms of an
effective aspect ratio.

To this end, we denote the length of the magnetic axis by

$2\pi R=\oint dl,$ (19)

where the integral is taken along the axis once around the torus. To lowest
order in the distance from the magnetic axis, the volume enclosed by a flux
surface $\psi=\psi_{a}$ is equal to (Helander, 2014)

$V=2\pi\psi_{a}\oint\frac{dl}{B_{0}},$ (20)

where $B_{0}$ is the magnetic strength on axis. It is natural to define an
“average” minor radius $a$ by setting $V=2\pi R\cdot\pi a^{2}$, and a logical
definition of the aspect ratio is then

$A=\frac{R}{a}=\frac{1}{\pi}\sqrt{\left(\oint
dl\right)^{3}\bigg{/}\left(8\psi_{a}\oint\frac{dl}{B}\right)},$ (21)

which may be expressed in terms of averages over the Boozer toroidal angle
$\varphi$,

$\psi_{a}=\frac{1}{2}\left(\frac{G_{0}}{A}\right)^{2}\frac{\overline{1/B_{0}}^{3}}{\overline{1/B_{0}^{2}}},$
(22)

where the overline represents a toroidal average
$\overline{(\dots)}=\int_{0}^{2\pi}(\dots)\mathrm{d}\varphi/2\pi$. We shall
use this form of the toroidal flux in the presentation of numerical results to
follow (choosing the representative value of $A\sim 10$).

To illustrate the behaviour of the precession in some near-axis constructions,
we first consider configurations recently constructed by Camacho Mata et al.
(2022), which were designed to be quasi-isodynamic to first order in the
distance from the magnetic axis. This work emphasised the reduction of shaping
and neoclassical transport losses when a global equilibrium was constructed
using the near-axis field. Because only first-order considerations were taken
into account, there is in principle not a unique precession frequency
characterising these configurations since a certain degree of freedom exists
at second order to complete the construction. Nevertheless, there is a
‘natural’ second-order extension of these configurations, which we refer to as
the minimal-shaping construction, namely, the one that makes the $X_{2c}$ and
$X_{2s}$ modulations in the near-axis vanish. The resulting field should be
representative of the first-order construction, especially if one considers
the construction of a global solution using the first-order fields. We
construct the field following the equilibrium equations in Landreman &
Sengupta (2019) (see Appendix B), using the code pyQIC (Jorge et al., 2023).
With such a second-order field in place, we may calculate the precession
frequency, which is plotted in Fig. 6 for the $N=2$ and $N=3$ configurations
of Camacho Mata et al. (2022), where the necessary integrals were computed
following Eq. (51) and (52). Some additional details are included in Table 1.

Figure 6: Precession of trapped electrons for two near-axis QI examples. 3D
rendition of flux surfaces at $r=0.1$ (with the colour map representing the
magnetic field strength) and normalised precession frequency
$\hat{\omega}_{\alpha}$ at two different radii for near-axis fields with $N=2$
and $N=3$ of Camacho Mata et al. (2022). The precession is computed using the
analytical expressions derived in this paper. The near-axis fields have been
constructed to second order, taking a ‘minimal shaping’ construction
$X_{2c}=0=X_{2s}$. The gray curves denote the variation of precession between
different field lines (with the black curve corresponding to the average over
$\alpha$), reflecting the non-omnigeneous nature of the fields. Two different
origins of the variability are apparent: a roughly $r$ independent variation
from the second order contribution, and a variation proportional to $1/r$ due
to breaking omnigeneity at first order.

Similarly to the global QI-optimised stellarators in Fig. 1, these
configurations do not have the maximum-$J$ property. Most of the trapped
particles precess in the diamagnetic direction, but some particles with
turning points close to the magnetic-field maxima behave as expected in a
maximum-$J$ field. This behaviour is similar to that in
axisymmetric/quasisymmetric fields: deeply trapped particles tend to co-
precess with the diamagnetic frequency, while barely trapped ones do the
opposite. The configurations importantly exhibit a significant field-line-to-
field-line variation (see the different grey lines), as a result of deviations
from omnigeneity. This should not come as a surprise given that the second-
order construction is not QI. But in addition, as noted in the construction of
Camacho Mata et al. (2022), the tops of the wells deviate from QI already at
first order. The result is a variability of $\hat{\omega}_{\alpha}$ that
diverges as $r\rightarrow 0$ and primarily affects shallowly trapped
particles, see Fig. 6. Note, however, that this first-order effect is only
noticeable very close to the axis. This is testimony to the quality of the QI
optimisation performed by Camacho Mata et al. (2022).

Although the precise form of the curves in Fig. 6 depends on how the near-axis
magnetic field has been completed at second order, the behaviour of deeply
trapped particles is independent of this detail, as we learnt in previous
sections. For these particles, it suffices to compute the various terms that
make up $\omega_{\alpha}$ in Eq. (9) at the minimum of $B_{0}$. For this
calculation, though, one should not use the form of $f_{0}$ in Eq. (16) but
instead relax the assumption of QI at second order. Setting the explicit non-
QI contributions aside for now, we must use the expression for $B_{20}$ needed
for Eq. (9). The necessary expressions are given in Appendix B, Eqs. (60),
where we use the notation $\mathcal{P}_{i}$ to denote the contribution of a
quantity $i$ to $B_{20}$, in analogy to $f_{i}$ in $f_{0}$. These terms can be
used to assess the near-axis construction at first order as in Table 1. As
$\bar{d}^{\prime\prime}=0$ (and thus
$\mathcal{P}_{\bar{d}\bar{d}^{\prime\prime}}=0$) and $\bar{d}<1$ in these
examples, it is a priori clear that it is impossible for the deeply-trapped
particles to have the maximum-$J$ property in these configurations. A finite
plasma $\beta$ would be necessary to attain this property. Much like the
critical $\beta_{J}$ in Eq. (17), we can define
$\beta^{\star}=a^{2}(\sum\mathcal{P}_{i})/(2(\ell^{\prime})^{2})$, where $a$
is the value of $r$ at the plasma edge, which we take to be at roughly $a\sim
R/10$. The quantity $\beta^{\star}$ represents the plasma $\beta$ necessary to
make the vacuum QI configuration reverse the behaviour of deeply trapped
particles. The key features for the equilibria in Figure 6 are collected in
Table 1. We have at this point neglected the contribution of
$\omega_{\alpha}^{\mathrm{non-QI}}$ to Eq. (9), and since the configurations
listed in Table 1 are not exactly omnigenous, it is actually necessary to
increase $\beta$ further in order to overcome the $\alpha$-dependence. The
requisite $\beta$ can be estimated from Eqs. (53) and (62), and we may thus
define $\beta^{\star}_{\mathrm{tot}}$ as the total $\beta$ needed to make the
‘least-maximum-$J$’ field line be so.

| $\bar{\alpha}$ | $\mathcal{P}_{B_{0}^{\prime\prime}}$ | $\mathcal{P}_{\tau_{0}^{2}}$ | $\mathcal{P}_{\bar{d}\bar{d}^{\prime\prime}}$ | $\mathcal{P}_{\mathrm{QI}}$ | $\beta^{\star}$ | $\beta^{\star}_{\mathrm{tot}}$
---|---|---|---|---|---|---|---
N2 | 0.21 | -1.83 | -10.48 | 0 | -8.82 | 5.5% | 15.5%
N3 | 0.16 | -10.01 | 0 | 0 | -17.14 | 4.3% | 4.7%
Table 1: Critical $\beta$ and geometric contributions. Geometric parameters
and critical $\beta$ for the $N=2$ and $N=3$ near-axis QI examples in Camacho
Mata et al. (2022). Here $\beta^{\star}$ is a measure of the required plasma
beta to prevent deeply trapped particles from precessing in the diamagnetic
direction in an idealised omnigeneous field ($\beta_{\mathrm{tot}}^{\star}$ if
the non-omnigeneous nature of the field is considered).

The analytical value of $\beta^{\star}$ can be seen to be correct in the
example of Figure 7, where the plasma $\beta$ of the configuration in Fig. 6a
increases. Here, the pressure gradient has been varied whilst the shape of the
magnetic axis and the ellipticity of the flux surfaces in its vicinity is kept
fixed to first order. Note that different authors and contexts mean different
things by ‘increasing $\beta$’ (i.e., different features of the equilibrium
are kept constant), making direct comparison difficult. Of course, reaching
this $\beta$ is necessary but not sufficient for maximum-$J$ of the whole
trapped population. This is especially true in the examples in Fig. 6, where
it is not the most deeply trapped particles that have the largest precession
frequency. From the preceeding analysis, we expect a non-zero plasma $\beta$
to introduce an overall upshift of the precession $q\omega_{\alpha}$, which is
indeed seen numerically.

Figure 7: Reversal of trapped-particle precession with increasing $\beta$.
Normalised precession frequency $\hat{\omega}_{\alpha}$ for three different
values of plasma $\beta$ for the ‘minimally-shaped’ near-axis field $N=2$ from
Camacho Mata et al. (2022) at $r=0.01$. The black lines correspond to the
$\alpha=\pi/4$ field line, and gray curves reflect the variation in precession
frequency due to non-omnigeneity. The analytical estimate for the normalised
pressure at which deeply trapped particles reverse their precession is
$\beta^{\star}=5.5\%$ and corresponds to the middle set of curves. The legend
on the right shows how the plasma cross section at the radii $r=0.05,~{}0.1$
and toroidal angle $\phi=0$ changes with $\beta$.

A central conclusion of our analysis thus far is the possibility of making
deeply trapped particles acquire maximum-$J$ behaviour without the need of a
nonzero plasma pressure. Accordingly, we now attempt to construct such a field
first through optimisation within the near-axis framework. We construct an
optimisation measure for maximum-$J$ by summing $\omega_{\alpha}^{2}$ (which
we have learnt how to compute) over values of $k$ that satisfy
$q\omega_{\alpha}<0$, call it $g_{\omega_{\alpha}}$. We are also interested in
imposing the condition of QI, especially at second order. At the maxima and
minima of the field strength we learnt in this work that we must satisfy Eq.
(64), from which we may construct an additional cost function, call it
$\check{g}_{\mathrm{QI}}$. Under the assumption of satisfying this condition,
one can show that it is possible to choose the near-axis construction at
second order in such a way that it guarantees the correct second order QI
behaviour elsewhere. We give the most essential elements of this in Appendix
D, but leave a full exploration to a later publication. Note that this way of
completing the solution is formally correct, but in practice (i.e., when
taking into account shaping, QI breaking in buffer regions, etcetera) it may
not be the best choice. For our proof-of-principle, though, it should suffice.
With this, then, we construct our near-axis cost function simply as the
weighted sum of the negative $q\omega_{\alpha}$ and the QI condition, Eq.
(64), at the $|\mathbf{B}|$ extrema,
$g=g_{\omega_{\alpha}}+\check{g}_{\mathrm{QI}}$.

With the cost function thus defined, we must explicitly state which our
minimal degrees-of-freedom are. We shall allow only the axis shape and the
function $\bar{d}(\varphi)$ to vary, while keeping the field strength
$B_{0}(\varphi)$ fixed. The idea is not to find a practical field, which would
require limiting the shaping and other additional practical features (and
unlocking other degrees of freedom such as $B_{0}$ and the order of curvature
zeroes). We are simply aiming at the construction of a proof-of-principle
field that exhibits maximum-$J$ behaviour in vacuum, especially for the deeply
trapped particles. Other details of the optimisation are left to Appendix D.
In Figure 8 we present the resulting optimised configuration for $N=1$. The
shaping is forbiddingly large, owing to the fact that no attention was paid to
limit it. Flux surfaces are extremely shaped and limit the physical radius of
the configuration. However, we may from this approach formally construct an
asymptotic $B$ at any $r$ (as a model if one wishes, $r_{c}\approx 0.002$
(Landreman, 2021)), and represent the geometry to first order to exhibit the
features resulting from the optimisation.

The optimised field is one in which the deeply trapped particles exhibit
maximum-$J$ behaviour (see Fig. 8c).555There exists some noise very close to
the minimum owing to the appearence of ‘puddles’ in $|\mathbf{B}|$ (Rodríguez
& Plunk, 2023). This is the result of the simple axis shapes considered for
this case. The noise should however disappear as $r\rightarrow 0$. To achieve
this, the optimiser has found a field with the features identified in the
preceding section. Specifically, the flux surfaces are twisted in the
necessary manner and elongated away from the region of low magnetic field
strength. Although maximum-$J$ behaviour has thus been achieved for deeply
trapped particles, it is clear that the field is is not exactly omnigenous
especially at larger $k$. This is mainly due to the first order deviations
from omnigeneity which the optimisation target did not include (clearly seen
by the growth of the variation going from the bottom to the top plot of Fig.
8). Reducing the ‘buffer region’ (the region where omnigeneity is violated at
first order) would reduce number of trapped-particles that do not satisfy the
maximum-$J$ condition.

Figure 8: Proof-of-principle of a maximum-$J$-optimised near-axis field. The
plots show a near-axis field optimised so as to exhibit quasi-isodynamicity
and maximum-$J$ behaviour, especially near the minimum of the magneic field on
each flux surface. The shaping of the configuration is very large, and thus
the resulting field impracticable, but shows that the optimisation criteria
can be met. a) A 3D rendition of the field for $r=0.02$, using only the 1st
order description. (Second-order contributions are large and obscure the
visualisation.) b) First-order cross sections near the minimum of
$|\mathbf{B}|$ (in $(R,Z)$ coordinates and between $\phi=3\pi/4,~{}5\pi/4$).
The dotted curve represents the position of the magnetic axis. c) Precession
$\hat{\omega}_{\alpha}$ computed numerically using $|\mathbf{B}|$ from the
near-axis-expansion at two different radii on a number of field lines (gray
curves, with the black representing the average). The increase in variability
at low $r$ is due to the buffer region in which omnigeneity is broken. For
this example, this contribution was not minimised. Close to the magnetic axis,
the field satisfies the maximum-$J$ criterion for almost every orbit, except
those trapped in secondary minima (Rodríguez & Plunk, 2023).

###  3.2 Traditionally optimised equilibrium

Motivated by the theoretical possibility of a vacuum, QI, maximum-$J$ field,
and proceeding beyond the near-axis expansion, we now turn to more
traditional, global, stellarator optimisation. To this aim, we employ three
target functions: $g_{QI}$ ensures that the field is QI,
$g_{B_{\mathrm{min}}}$ ensures that the most deeply-trapped particles satisfy
the maximum-$J$ criterion, and $g_{J}$ ensures that the other particles also
do so. The total target function that we minimise is thus

$g=g_{QI}+g_{B_{\textrm{min}}}+g_{J}.$ (23)

Both $g_{QI}$ and $g_{J}$ are complicated targets, which are explained in
greater detail by Goodman (2024), and so is the starting point for the
optimisation. Broadly speaking, $g_{QI}$ penalises the difference between
$J$-contours on a flux surface by computing the difference in the second
adiabatic invariant and a closely-related, artificially-constructed, perfectly
QI flux surface. From the calculation of $\mathcal{J}_{\|}$, we also evaluate
the term $g_{J}$, which imposes $\partial_{\psi}J<0$ in this constructed
field. Thus, as $g_{QI}$ and $g_{J}$ decrease, the field becomes more QI and
more maximum-$J$.

The term $g_{B_{\textrm{min}}}$ is a more straightforward target function.
Using the fact that the maximum-$J$ condition corresponds to minimum-$B$ (for
the most deeply- and shallowly-trapped particles), we simply designed this
function to encourage the flux-surface’s minimum field strength,
$B_{\textrm{min}}(s)$, to have a positive derivative, where $s=(r/a)^{2}$. To
do this, for every consecutive pair of flux surface, the target calculates
$B_{\textrm{min}}$ on $s_{0}$ and $s_{1}=s_{0}+\delta s$, and then the
fractional difference between the two

$\delta_{B}(s_{0})=\frac{1}{\delta
s}\frac{B_{\textrm{min}}(s_{1})-B_{\textrm{min}}(s_{0})}{B_{\textrm{min}}(s_{1})+B_{\textrm{min}}(s_{0})}.$
(24)

We can thus define the target as
$f_{B_{\textrm{min}}}=\max(0.01-\delta_{B},0)^{2}$, where the value of 0.01
has been chosen as an arbitrary positive number here.

We initialise optimisation from a near-axis construction with a flatter bottom
of the magnetic well (which we know from the work above that should favour
maximum-$J$), using techniques to be presented by Plunk (2024). The resulting
optimised configuration is indeed (mostly) maximum-$J$ and QI. This can be
seen in the second adiabatic invariant contours of Figure 10, which show that
$\mathcal{J}$ decreases from the centre of the $(s,\alpha)$ plot radially
outward, and the contours are approximately circular. It is evident though
that for trapped particles sufficiently close to the bottom of the wells (see
left-most plot), there is a fraction of the population that does not behave in
a maximum-$J$ fashion. This is signalled here by a relative shift of the
approximately circular contours respect to the centre. As a result, there is a
portion of field lines that are minimum-$J$ (see also the top plot of Figure
9c). In the spirit of a more quantitative measure, it would be convenient to
come up with a single scalar that indicates what ‘fraction’ of the
configuration truly behaves in a maximum-$J$ fashion. For a Maxwellian
distribution function, the fraction of the trapped-particle population whose
rotation satisfies $q\omega_{\alpha}>0$ is equal to

$f_{\mathrm{max-J}}=\frac{1}{\mathcal{N}}\int_{0}^{1}\varrho\mathrm{d}\varrho\int_{0}^{2\pi}\mathrm{d}\alpha\int_{1/B_{\mathrm{max}}}^{1/B_{\mathrm{min}}}\Theta[q\omega_{\alpha}(\lambda,\alpha,\varrho)]\hat{\tau}_{b}(\lambda,\alpha,\varrho)\mathrm{d}\lambda,$
(25)

where $\Theta$ is the Heaviside step function, $\varrho=\sqrt{\psi/\psi_{a}}$,

$\mathcal{N}=\int_{0}^{1}\varrho\mathrm{d}\varrho\int_{0}^{2\pi}\mathrm{d}\alpha\int_{1/B_{\mathrm{max}}}^{1/B_{\mathrm{min}}}\hat{\tau}_{b}(\lambda,\alpha,\varrho)\mathrm{d}\lambda,$
(26)

and

$\hat{\tau}_{b}=\int_{\ell_{L}}^{\ell_{R}}\frac{\mathrm{d}\ell}{\sqrt{1-\lambda
B}}$ (27)

is proportional to the bounce time of trapped particles. This way each flux
surface is considered equally important (the Maxwellian distribution does not
change from surface to surface). With this definition, for the optimised
configuration above, $f_{\mathrm{max-J}}\approx 0.996$. That is, a $99.6\%$ of
the volume is maximum-$J$; i.e., for all intents and purposes, the
configuration is maximum-$J$. The omnigeneous nature of the configuration can
be quantified by the effective ripple (Nemov et al., 1999), which measures the
$1/\nu$ neoclassical transport and in the present case remains in the range
$\epsilon_{\mathrm{eff}}\sim 0.15\%-0.45\%$ across the volume.

The resulting configuration also has a significant so-called “vacuum magnetic
well” of about 3.3% (as defined in Landreman & Paul (2022)), a property that
is important for MHD stability (Landreman, 2021). A field is said to have a
vacuum magnetic well if $V^{\prime\prime}(\psi)<0$ (where $V$ is the volume
enclosed by a flux surface). Generally, unless special care is taken,
optimised configurations tend to have a “magnetic hill”, $i.e.$
$V^{\prime\prime}(\psi)>0$, so it is notable that this configuration has such
a substantial vacuum well.

The condition for a magnetic well means that, on a suitable average, the
magnetic field strength should increase with minor radius. This can be seen
from the circumstance that the volume enclosed by a flux surface $\psi$ is

$V(\psi)=\int_{0}^{\psi}d\psi^{\prime}\int_{0}^{2\pi}d\theta\int_{0}^{2\pi}\mathcal{J}d\varphi,$
(28)

where the Jacobian in Boozer coordinates is $\mathcal{J}=(G+\iota I)/B^{2}$.
Since the flux-surface average is defined as

$\left\langle\cdots\right\rangle=\frac{1}{V^{\prime}(\psi)}\int_{0}^{2\pi}d\theta\int_{0}^{2\pi}(\cdots)\mathcal{J}d\varphi,$
(29)

it follows that $\langle B^{2}\rangle=4\pi^{2}(G+\iota I)/V^{\prime}$ and
consequently

$\frac{V^{\prime\prime}}{V^{\prime}}=-\frac{d\ln\langle
B^{2}\rangle}{d\psi}+\frac{d}{d\psi}\ln(G+\iota I),$ (30)

where the last term vanishes in a vacuum field, where $G^{\prime}(psi)$ and
$I\psi)$ vanish. The magnetic well is sometimes defined without this term
(Freidberg, 2014). Thus we see, as argued in Section 1, that increasing the
radial derivative of $B$ is directly related to MHD stability. Indeed,
(Helander, 2014, p. 24) showed that there is a mathematical relation between
maximum-$J$ and magnetic-well criteria, but they are not identical. In the
limit of large mirror ratio, where most particles are magnetically trapped,
then the maximum-$J$ property rigorously implies a magnetic well. It is
therefore not surprising that the particular configuration optimised here
possesses a vacuum magnetic well.

Figure 9: Example of a maximum-$J$-optimised global equilibrium. The plots
show a magnetic field optimised by minimising Eq. (23) so as to attain quasi-
isdynamicity and maximum-$J$ behaviour in a vacuum magnetic field. The
configuration has an aspect ratio of $A_{\mathrm{VMEC}}\sim 7$, three field
periods, and exhibits strong shaping, which was not constrained in this proof-
of-principle example. a) 3D rendition of the outermost surface of the field,
where the colour map represents $|\mathbf{B}|$. b) Detail of cross-sections
near the core (at $\varrho=0.1$, and for reference to indicate the large
shaping of surfaces $\varrho=0.5$ as broken contours) and about the minimum of
$|\mathbf{B}|$ (in $(R,Z)$ coordinates and between $\phi=3\pi/4,~{}5\pi/4$),
showing features of twist and shape studied analytically. The dotted curve
represents the position of the magnetic axis. c) Precession frequency at two
radii on a number of field lines (gray curves, with the black $\alpha=\pi/2$).
The increase in variability at low $r$ is a consequence of the breaking of
omnigeneity. The overwhelming majority of all trapped particles satisfy the
maximum-$J$ criterion, $\hat{\omega}_{\alpha}>0$. Figure 10: Contours of the
second adiabatic invariant $\mathcal{J}_{\parallel}$ in polar coordinates
$(s,\alpha)$. Contours showing the second adiabatic invariant as a function of
the polar coordinates $(s,\alpha)$ for three different trapped particle
classes. The values of $\lambda$ used are
$\lambda=[(B_{\mathrm{max}}-B_{\mathrm{min}})\lambda_{n}+B_{min}]^{-1}$ where
$B_{\mathrm{max}}$ and $B_{\mathrm{min}}$ denote the maximum and minimum field
strength on the flux surface in question. An ideal omnigeneous field would
have concentric circular contours. The wiggles in the contours is indicative
of QI-breaking, which is particularly prominent close to the trapped-passing
boundary (small $\lambda$). An ideally maximum-$J$ field would show a
monotonic decrease of $\mathcal{J}_{\parallel}$ along any ray emanating from
the origin.

##  4 Discussion and conclusions

Before concluding the paper, we comment on the connection between the
maximum-$J$ condition and the suppression of trapped-particle instabilities
(Proll et al., 2012; Helander et al., 2013) as well as the general reduction
of curvature-driven instabilities and turbulence (Proll et al., 2022). These
phenomena are sensitive to the relative size of the drift frequency
$\omega_{d}=k_{\alpha}\omega_{\alpha}$ and the diamagnetic frequency

$\omega_{\star}=\frac{k_{\alpha}T}{q}\frac{\ln n}{d\psi},$ (31)

where $T$ denotes the plasma temperature in units of energy, $n$ the density,
and the wave vector has been taken to be ${\bf k}=k_{\alpha}\nabla\alpha$. It
is thus of interest to assess the magnitude of the ratio

$\frac{\omega_{d}}{\omega_{\ast}}=\frac{mv^{2}}{4T}\frac{\hat{\omega}_{\alpha}}{\psi_{a}d\ln
n/d\psi},$ (32)

which for a density profile $n=n_{0}(1-\psi/\psi_{a})$ becomes
$\omega_{d}/\omega_{\ast}=-\hat{\omega}_{\alpha}/2$ for particles moving at
the thermal speed, $v=\sqrt{2T/m}$. In the QI cases studied numerically in the
present paper, $\hat{\omega}_{\alpha}$ is of order $1/10$, making
$\omega_{d}/\omega_{\ast}\ll 1$ for most of the particles in a thermal
distribution. This circumstance is stabilising to trapped-electron modes even
if the absence of the maximum-$J$ property. Relatively few electrons are able
to satisfy the resonance condition in the limit
$\omega_{d}/\omega_{\ast}\rightarrow 0$ (Connor et al., 2006), and the
nonlinear energy available to trapped electrons for driving instabilities is
small (Mackenbach et al., 2022, 2023a). The latter was recently studied and
labelled as the ‘strong’ regime (Rodríguez & Mackenbach, 2023), in which the
available energy is only supplied by a small fraction of the trapped
population $\sim\omega_{\alpha}/\omega_{\star}$. From a practical point of
view, it is true that as long as $|\hat{\omega}_{\alpha}|\ll 1$ the magnetic
field may not need to strictly satisfy the maximum-$J$ condition in order to
substantially improve trapped-electron-mode stability. This condition is
nevertheless helpful for MHD stability and fast-particle confinement, since
particles with $\omega_{\alpha}=0$ are prone collisionless losses (Velasco et
al., 2021; Paul et al., 2022). In any case, whenever the magnetic field of any
concrete stellarator is optimised, a compromise will need to find between the
maximum-$J$ condition and other desirable properties.

It is fortunate, though not merely a happy coincidence, that the maximum-$J$
property is related to both MHD stability, trapped-electron-mode suppression,
and good fast-ion confinement: three very beneficial properties that are
otherwise considered separate and independent of each other. Quasi-symmetric
stellarators struggle to possess the maximum-$J$ property, but we show that
quasi-isodynamic ones are more amenable to it, and can be designed
accordingly. Physically, this property is realised if all trapped particles
experience average favourable magnetic curvature, and the precession frequency
of these particles then has the opposite sign from the diamagnetic frequency.
The behaviour of the most deeply and the most shallowly trapped particle
orbits can be examined analytically. It is found that the former usually
precess in the wrong (i.e. diamagnetic) direction. This behaviour may be
modified by appropriate twisting and flattening of the field near the bottom
of the trapping well. Barely trapped particles always precess as desired,
provided the magnetic field is perfectly omnigenous. Nevertheless, omnigenity
always fails for some of these particles, making them, too, precess in the
wrong direction (as explicitly shown). As we have seen, despite this formal
limitation, it is possible to design a vacuum QI field in which the vast
majority of all particles satisfy the maximum-$J$ condition.

## Data availability

The data that support the findings of this study are openly available at the
Zenodo repository with DOI/URL 10.5281/zenodo.10203417.

## Acknowledgements

Many of the numerical utilities to handle VMEC equilibria (Hirshman & Whitson,
1983) are based on the stellarator package simsopt (Landreman et al., 2021).
We gratefully acknowledge fruitful discussion with J. L. Velasco and R. J. J.
Mackenbach. E. R. was supported by a grant by Alexander-von-Humboldt-Stiftung,
Bonn, Germany, through a postdoctoral research fellowship, and A. G. supported
by a grant from the Simons Foundation (560651, PH). This work has been carried
out within the framework of the EUROfusion Consortium, funded by the European
Union via the Euratom Research and Training Programme (Grant Agreement No
101052200 — EUROfusion). Views and opinions expressed are however those of the
authors only and do not necessarily reflect those of the European Union or the
European Commission. Neither the European Union nor the European Commission
can be held responsible for them.

##  Appendix A Precession frequency

###  A.1 General expressions

In accordance with the discussion in Section 1, we first turn our attention to
the derivatives of the adiabatic invariant $\mathcal{J}_{\parallel}$. In
Boozer coordinates, the magnetic field can be written as (Boozer, 1981;
Helander, 2014)

${\bf
B}=G(\psi)\nabla\varphi+I(\psi)\nabla\theta+B_{\psi}(\psi,\theta,\varphi)\nabla\psi,$
(33)

the arc length element is $d\ell=(G+\iota I)d\varphi/B$, so the adiabatic
invariant becomes

$\mathcal{J}_{\parallel}=mv(G+\iota I)\int\frac{\sqrt{1-\lambda
B}}{B}d\varphi,$

and its derivatives are

$\left.\frac{\partial\mathcal{J}_{\parallel}}{\partial\alpha}\right|_{\psi,\lambda}=-mv(G+\iota
I)\int\left(\frac{\partial
B}{\partial\alpha}\right)_{\psi,\varphi}F(\lambda,B)d\varphi,$ (34)
$\left.\frac{\partial\mathcal{J}_{\parallel}}{\partial\psi}\right|_{\alpha,\lambda}=mv\left[(G+\iota
I)^{\prime}\int\frac{\sqrt{1-\lambda B}}{B}d\varphi-(G+\iota
I)\int\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,\varphi}F(\lambda,B)d\varphi\right],$ (35)

with

$F(\lambda,B)=\frac{1-\lambda B/2}{B^{2}\sqrt{1-\lambda B}},$

and the contribution from the boundary terms vanishing at the end points of
the bounce integral. Equation (35) can be simplified by using the MHD force
balance equation. In particular, using Eq. (33) to find the plasma current,
computing $\mathbf{j}\times\mathbf{B}$ using the Boozer Jacobian and flux-
surface averaging the $\nabla\psi$ component of the equation (Kruskal &
Kulsrud, 1958; Helander, 2014)

$G^{\prime}+\iota I^{\prime}=-\mu_{0}p^{\prime}\frac{G+\iota I}{\langle
B^{2}\rangle},$

where $p(\psi)$ denotes the plasma pressure and angular brackets indicate the
usual flux-surface average. Hence

$\left.\frac{\partial\mathcal{J}_{\parallel}}{\partial\psi}\right|_{\alpha,\lambda}=mv(G+\iota
I)\left[\left(\frac{\iota^{\prime}I}{G+\iota
I}-\frac{\mu_{0}p^{\prime}}{\langle
B^{2}\rangle}\right)\int\frac{\sqrt{1-\lambda
B}}{B}d\varphi-\int\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,\varphi}F(\lambda,B)d\varphi\right].$ (36)

The bounce time is

$\Delta t=\int\frac{d\ell}{v_{\|}}=\frac{G+\iota
I}{v}\int\frac{d\varphi}{B\sqrt{1-\lambda B}},$

and the precession frequency $\omega_{\alpha}=\Delta\alpha/\Delta t$, using
Eq. (1b), thus becomes

$\omega_{\alpha}=\frac{mv^{2}}{q}\left[\int\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,\varphi}F(\lambda,B)\;d\varphi-\left(\frac{\iota^{\prime}I}{G+\iota
I}-\frac{\mu_{0}p^{\prime}}{\langle
B^{2}\rangle}\right)\int\frac{\sqrt{1-\lambda
B}}{B}d\varphi\right]\bigg{/}\int\frac{d\varphi}{B\sqrt{1-\lambda B}}.$ (37)

The same expressions may be reached by considering the bounce integral of the
poloidal and radial drifts directly.

###  A.2 Precession for deeply-trapped particles

The expression in Eq. (37) is valid in any stellarator with flux surfaces, but
becomes much simpler for deeply trapped particles. In this case $1-\lambda
B\ll 1$ in the full integration domain, which also shrinks as $\lambda^{-1}$
approaches the minimum of $|\mathbf{B}|$. This makes the first term in the
square brackets dominate over the second one, and $F(\lambda,B)\approx
1/(2B^{2}\sqrt{1-\lambda B})$ independently of the shape of the magnetic well.
To show this, one can approximate the bottom of the well as $B\approx
B_{0}(1+B^{(n)}\varphi^{n}/n!)$ and evaluate for a deeply trapped particle
with bounce point $\varphi_{b}$,

$\displaystyle\int_{0}^{\varphi_{b}}\frac{\sqrt{1-\lambda
B}}{B}\mathrm{d}\varphi\approx\frac{\sqrt{\pi}}{2B_{0}}\frac{\Gamma(1+1/n)}{\Gamma(3/2+1/n)}\left(\lambda
B_{0}\frac{B^{(n)}}{n!}\right)^{-1/n}(1-\lambda B_{0})^{1/2+1/n},$ (38a)
$\displaystyle\int_{0}^{\varphi_{b}}\frac{1}{\sqrt{1-\lambda
B}}\mathrm{d}\varphi\approx\sqrt{\pi}\frac{\Gamma(1+1/n)}{\Gamma(1/2+1/n)}\left(\lambda
B_{0}\frac{B^{(n)}}{n!}\right)^{-1/n}(1-\lambda B_{0})^{-1/2+1/n}.$ (38b)

The second term in the square bracket of Eq. (37) thus yields a contribution
that is $\sim(1-\lambda B_{0})$ smaller than the first, and thus vanishes for
deeply trapped particles regardless of the shape of the well.

As a result, at the bottom of the magnetic well (i.e., at the point where
${\bf B}\cdot\nabla B=0$) the precession frequency becomes simply

$\omega_{\alpha}=\frac{mv^{2}}{2qB}\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,\varphi}.$ (39)

This result regarding the local role of $\partial_{\psi}B$ also applies to the
mostly shallowly trapped particles, which spend most of the time close to the
turning point where the magnetic field strength reaches its maximum on the
field line in question, $\lambda B_{\rm max}\simeq 1$.

The partial derivative of $B$ keeping $(\alpha,\varphi)$ constant can be
expressed in terms of common Boozer derivatives at fixed $(\theta,\varphi)$ by
using the chain rule,

$\partial_{\psi}B|_{\alpha,\varphi}=\partial_{\psi}B|_{\theta,\varphi}+\iota^{\prime}\varphi\partial_{\theta}B|_{\psi,\varphi}.$
(40)

In this expression, the second term must vanish at the bottom of the magnetic
well if the field is QI (in fact pseudosymmetry is sufficient (Mikhailov et
al., 2002; Skovoroda, 2005)), as can be seen from the following argument. By
definition, $(\iota\partial_{\theta}+\partial_{\varphi})B=0$ at the bottom of
the magnetic well. Moreover, in a QI field, the level contours of $B=|{\bf
B}|$ are poloidally closed, which implies that
$(c\partial_{\theta}+\partial_{\varphi})B=0$ for some
$c(\psi,\theta,\varphi)\neq 0$ that cannot equal $\iota$, at least not at all
minima of $B$. At the bottom of the well, then, the directional derivative of
$B$ vanishes in two independent directions tangential to the magnetic surface,
which implies that it vanishes in any such direction including
$\partial_{\theta}B=0$. As a result, Eq. (40) becomes

$\partial_{\psi}B|_{\alpha,\varphi}=\partial_{\psi}B|_{\theta,\varphi}.$ (41)

A similar argument also applies to maxima of the field strength.

###  A.3 Near-axis expansion

We now consider the region close to the magnetic axis, expanding in the minor-
radius variable $r=\sqrt{2\psi/\bar{B}}$, where $\bar{B}$ is some reference
magnetic-field strength. To second order, this expansion is (Garren & Boozer,
1991; Rodríguez & Plunk, 2023)

$B=B_{0}(\varphi)+rB_{1}(\alpha,\varphi)+r^{2}B_{2}(\alpha,\varphi),$

with

$\displaystyle
B_{1}(\alpha,\varphi)=B_{0}(\varphi)d(\varphi)\cos\left[\alpha+\nu(\varphi)\right],$
(42a) $\displaystyle B_{2}(\alpha,\varphi)=B_{20}(\varphi)+B_{2c}(\varphi)\cos
2\left[\alpha+\nu(\varphi)\right]+B_{2s}(\varphi)\sin
2\left[\alpha+\nu(\varphi)\right].$ (42b)

To check for omnigeneity, let us write the integral in Eq. (34) for
$\partial_{\alpha}\mathcal{J}_{\parallel}$,

$\int\frac{\partial
B}{\partial\alpha}F(\lambda,B)d\varphi\approx\int\left(r\frac{\partial
B_{1}}{\partial\alpha}+r^{2}\frac{\partial
B_{2}}{\partial\alpha}\right)\left[F(\lambda,B_{0})+rB_{1}\frac{\partial
F}{\partial B}\right]d\varphi,$ (43)

where $\partial F/\partial B=(\partial
F(\lambda,B_{0}))/\partial\varphi)/B_{0}^{\prime}(\varphi)$ and the integral
is meant between bouncing points $\varphi_{R,L}$. Expanding to second order in
$r$ and integrating by parts666Some care is called for when computing
perturbations of bounce integrals since the limits of integration (bounce
points) are themselves perturbed and one must guarantee that the argument in
$\sqrt{1-\lambda B}$ is positive. This can be nevertheless done carefully to
arrive at the same final QI conditions.

$\int F(\lambda,B_{0})\left[r\frac{\partial
B_{1}}{\partial\alpha}+r^{2}\frac{\partial
B_{2}}{\partial\alpha}-r^{2}\frac{\partial}{\partial\varphi}\left(\frac{B_{1}}{B_{0}^{\prime}}\frac{\partial
B_{1}}{\partial\alpha}\right)\right]d\varphi\\\ =\int
F(\lambda,B_{0})\left[-rdB_{0}\sin(\alpha+\nu)+2r^{2}\left(B_{2s}\cos
2(\alpha+\nu)-B_{2c}\sin 2(\alpha+\nu)\right)\right.\\\
\left.+r^{2}\frac{\partial}{\partial\varphi}\left(\frac{B_{0}^{2}d^{2}}{2B_{0}^{\prime}}\sin
2(\alpha+\nu)\right)\right]d\varphi,$ (44)

where the integral is taken between bounce points of $B_{0}(\varphi)$. If the
field is to be omnigenous, $\partial_{\alpha}\mathcal{J}_{\parallel}=0$, this
quantity must vanish for all $\lambda\in[B_{\rm max}^{-1},B_{\rm min}^{-1}]$
and all $\alpha\in[0,2\pi]$. For a stellarator-symmetric field,
$B(r,\alpha,\varphi)=B(r,-\alpha,-\varphi)$, we have the relations

$\displaystyle B_{0}(\varphi)=B_{0}(-\varphi),$ (45a) $\displaystyle
d(\varphi)\cos[\alpha+\nu(\varphi)]=d(-\varphi)\cos[-\alpha+\nu(-\varphi)],$
(45b) which are satisfied if $d$ is chosen to be odd in $\varphi$, and
$\nu=(n+1/2)\pi$. As a result of these constraints, the first-order
contribution to Eq. (44) in $r$ vanishes. That at 2nd-order will also vanish
iff $\displaystyle B_{2s}(\varphi)=-B_{2s}(-\varphi),$ (45c) $\displaystyle
B_{2c}(\varphi)=\frac{\partial}{\partial\varphi}\left(\frac{B_{0}^{2}d^{2}}{4B_{0}^{\prime}}\right).$
(45d)

These relations form the basic requirements for omnigenity close to the
magnetic axis and were recently derived by Rodríguez & Plunk (2023) by a
different method that avoids the explicit integrals. We refer to that paper
for further details and discussion.

Employing the same near-axis explansion to (36) gives in analogy to Eq. (44)
for a general stellarator,

$\int F(\lambda,B)\left(\frac{\partial
B}{\partial\psi}\right)_{\alpha,\varphi}\mathrm{d}\varphi=\\\
\frac{2}{\bar{B}}\left[\frac{1}{2r}\int
F(\lambda,B_{0})B_{1}\mathrm{d}\varphi+\int
F(\lambda,B_{0})\left[B_{2}-\frac{\partial}{\partial\varphi}\left(\frac{B_{1}^{2}}{2B_{0}^{\prime}}\right)\right]\mathrm{d}\varphi\right],$
(46)

where the integrals are, once again, taken between bounce points of
$B_{0}(\varphi)$. Substituting the omnigeneous forms of the magnetic field
functions into this integral, and noting that the enclosed toroidal current
$I\sim r^{2}$ (in the absence of an infinite current density on axis), one
obtains $\omega_{\alpha}=\omega_{\alpha,\mathrm{vac}}+\omega_{\alpha,p}$ with

$\omega_{\alpha,\mathrm{vac}}=\frac{2mv^{2}}{q\bar{B}}\frac{\int_{\varphi_{L}}^{\varphi_{R}}\frac{1-\lambda
B_{0}/2}{B_{0}\sqrt{1-\lambda
B_{0}}}\left[\frac{B_{20}}{B_{0}}-\frac{1}{4B_{0}}\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)^{\prime}\right]\mathrm{d}\varphi}{\int_{\varphi_{L}}^{\varphi_{R}}\frac{\mathrm{d}\varphi}{B_{0}\sqrt{1-\lambda
B_{0}}}},$ (47)
$\omega_{\alpha,p}=\frac{2mv^{2}}{q\bar{B}}\frac{\mu_{0}p_{2}}{\langle
B_{0}^{2}\rangle}\frac{\int_{\varphi_{L}}^{\varphi_{R}}\frac{\sqrt{1-\lambda
B_{0}}}{B_{0}}\mathrm{d}\varphi}{\int_{\varphi_{L}}^{\varphi_{R}}\frac{\mathrm{d}\varphi}{B_{0}\sqrt{1-\lambda
B_{0}}}},$ (48)

where $p(r)=p_{0}+p_{2}r^{2}+\cdots$ and ${\varphi_{L}}$ and $\varphi_{R}$
denote the toroidal locations of the bounce point to the left (L) and right
(R) of the bottom of the well (for stellarator symmetry
$\varphi_{L}=-\varphi_{R}$). Although we have used the notation
$\omega_{\alpha,\mathrm{vac}}$ and $\omega_{\alpha,p}$, one should not forget
that the former includes pressure effects through $B_{20}$. The flux-average
of $B_{0}^{2}$ appearing in this expression is equal to

$\langle
B_{0}^{2}\rangle=2\pi\bigg{/}\int_{0}^{2\pi}\frac{\mathrm{d}\varphi}{B_{0}^{2}}.$
(49)

###  A.4 Omnigeneity breaking

Because the expression in Eq. (37) is valid for any stellarator, it provides
an opportunity to learn how deviations from exact omnigenity affect the
precession frequency $\omega_{\alpha}$. We consider this question from the
near-axis perspective.

####  A.4.1 First order breaking

Let us start with the QI requirement at first order in $r$. As we saw
following Eq. (44), the condition of omnigeneity requires the functions $d$
and $\nu$ to be odd and even, respectively, in $\varphi$ (see Plunk & Helander
(2018); Camacho Mata et al. (2022); Rodríguez & Plunk (2023)). This symmetry
is crucial in a quasi-isodynamic stellarator, as it causes the radial drift on
either side of the magnetic well to cancel to leading order (as discussed in
the main text). When the condition of omnigeneity is satisfied exactly, the
$1/r$ contribution to $\omega_{\alpha}$ in Eq. (46) vanishes. This balance is
delicate, however, as any small deviation can lead to a large precession,
restoring the $O(1/r)$ dependence. This behaviour can be understood as the
result of a finite poloidal drift occurring over an increasingly smaller flux
surface as the magnetic axis is approached. The contribution may be written as

$\omega_{\alpha,-1}^{\mathrm{non-
QI}}=\frac{mv^{2}}{q\bar{B}}\frac{1}{r}\left[{\int_{\varphi_{L}}^{\varphi_{R}}\frac{\mathrm{d}\varphi}{B_{0}\sqrt{1-\lambda
B_{0}}}}\right]^{-1}\int_{\varphi_{L}}^{\varphi_{R}}F(\lambda,B_{0})B_{1}(\alpha,\varphi)\mathrm{d}\varphi$
(50)

following directly from Eq. (46).

This term is in principle always present since omnigeneity can never be
realised exactly except in trivial cases (Cary & Shasharina, 1997; Landreman &
Catto, 2012; Plunk et al., 2019; Rodríguez & Plunk, 2023). In any non-
axisymmetric torus, the requirements of omnigeneity clash with the required
periodicity of the solution in $\theta$ and $\varphi$. Formally, periodicity
requires $\nu(-\pi)=\nu(\pi)+2\pi(\iota-N)$ in our near-axis description,
which is incompatible with its even symmetry. This condition is necessarily
broken at the points of maximum field strength, which we may write
$\nu=\pi/2+\tilde{\nu}$ and $\tilde{\nu}$ is odd in $\varphi$.Barely trapped
particles will always exerience this lack of omnigeneity.

Using this notation, we may rewrite the $O(1/r)$ contribution to the
precession,

$\omega_{\alpha,-1}^{\mathrm{non-
QI}}=-\frac{mv^{2}}{q\bar{B}}\frac{\cos\alpha}{r}\left[{\int_{\varphi_{L}}^{\varphi_{R}}\frac{\mathrm{d}\varphi}{B_{0}\sqrt{1-\lambda
B_{0}}}}\right]^{-1}\int_{\varphi_{L}}^{\varphi_{R}}F(\lambda,B_{0})B_{0}d\sin\tilde{\nu}\mathrm{d}\varphi,$
(51)

to which contributions arise from any particles that venture into the region
where omnigeneity is broken, i.e. $\sin\tilde{\nu}\neq 0$. For these orbits,
the non-omnigeneous part of the precession will always dominate for
sufficiently small $r$. Because of the $\cos\alpha$ dependence of the
precession, these particles will precess in either poloidal direction
depending on the field line considered. The immediate implication of this
behaviour is that in a realistic magnetic field (one that allows for the
necessary breaking of QI), the maximum-$J$ property is technically not
possible close to the magnetic axis, independently of plasma pressure.

This somewhat bleak result appears to rule out maximum-$J$ behaviour in QI
configurations, much as in quasisymmetric stellarators or tokamaks. It should
however be noted that the factors that determine the size of the $1/r$ terms
are naturally small. At the top of the magnetic well, where the conditions of
exact omnigeneity need to be relaxed, it is by construction the case that the
function $d$ vanishes. That is, the magnetic field perturbation $B_{1}$ is
small in the vicinity of the maximum, as is necessary for the consistency with
a locally straight magnetic field (and for the contour of $B_{\mathrm{max}}$
to be independent of the poloidal angle). By shrinking the region in which
omnigeneity is broken, the number of particles for which the precession
frequency diverges as $1/r$ can be made arbitrarily small777There is however a
price to be paid in terms of geometric complexity, see Camacho Mata et al.
(2022).. Only extremely close to the axis is the divergent behaviour of the
precession frequency significant.

The size of this region can be gauged by estimating the minor radius $r_{a}$
below which the $1/r$ behaviour, Eq. (51), dominates over the second-order
contributions, Eq. (47). By way of example, we define and compute the critical
radius $r_{a}=\sqrt{\sum_{k}w_{\alpha,-1}^{2}/\sum_{k}w_{\alpha,0}^{2}}$ (with
the $-1$ and 0 indices referring to the respective $r$ scalings) for the near-
axis fields in Fig. 6. The amplification of the precession and its variation
from field line to field line is especially noticeable in Fig. 6b for which
$r_{a}\sim 0.003$, driven by the presence of ‘buffer’ regions (Plunk et al.,
2019; Camacho Mata et al., 2022). It may perhaps appear counterintuitive that
close to the axis, where the field is most omngninous, the variation between
field lines (an indicator of departure from omnigeneity) increases. However,
the bounce-averaged radial drift $\omega_{\psi}$ does decrease linearly in $r$
as expected from Eq. (44)), indicating that omngeneity does become perfect in
the limit $r\rightarrow 0$. The existence of a region very close to the axis
in which the maximum-$J$ condition cannot be satisfied is nevertheless
unavoidable, but this region can in principle be made arbitrarily small.

####  A.4.2 Second order breaking

From Fig. 6 it is clear that there are non-QI contributions even if the
leading first order field is precisely QI. In the exactly omnigenous limit,
the contribution to $\omega_{\alpha}$ from second-order terms is
$\alpha$-independent, but one regains some $\alpha$ dependence when
omnigeneity is broken. Taking into account the symmetry conditions on the
various functions that constitute the near-axis description, we may write the
non-QI contribution to the precession as

$\omega_{\alpha,0}^{\mathrm{non-
QI}}=-\frac{2mv^{2}}{q\bar{B}}\left[{\int_{\varphi_{L}}^{\varphi_{R}}\frac{\mathrm{d}\varphi}{B_{0}\sqrt{1-\lambda
B_{0}}}}\right]^{-1}\cos 2\alpha\times\\\
\times\int_{\varphi_{L}}^{\varphi_{R}}F(\lambda,B_{0})\left\\{\left[B_{2c}-\frac{1}{4}\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)^{\prime}\right]\cos
2\tilde{\nu}+\left[B_{2s}+\frac{B_{0}^{2}d^{2}}{2B_{0}^{\prime}}\tilde{\nu}^{\prime}\right]\sin
2\tilde{\nu}\right\\}\mathrm{d}\varphi.$ (52)

Two limits of this expression are of particular interest. The first limit
concerns deeply trapped particles, which are relatively easily well confined
and may therefore be omnigenous to second order. Then $\tilde{\nu}=0$ and

$\omega_{\alpha,0}^{\mathrm{non-
QI}}=-\frac{2mv^{2}}{qB_{0}^{2}}\left[B_{2c}-\frac{1}{4}\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)^{\prime}\right]\cos
2\alpha.$ (53)

The second limit of interest is realised when $B_{2c}$ has its QI form, Eq.
(45d). Then the first term in square brackets of the second integral vanishes,
leaving a contribution that vanishes unless omnigeneity is broken at first
order. The QI contribution $\omega_{\alpha,\mathrm{vac}}$ in Eq. (47) remains
unchanged.

###  A.5 Stellarator-symmetry breaking

The expressions above were derived under the assumption of stellarator
symmetry. This is a convenient simplifying assumption, under which functions
acquire a definite parity in $\varphi$, which simplifies bounce integrals. If
$B_{0}$ is not symmetric in $\varphi$, the expressions need to be revisited
although most of the considerations for symmetric fields continue to hold.

To deal with the stellarator asymmetric case, we employ the notation of
Rodríguez & Plunk (2023) and denote the ‘bounce mapping’ by
$\eta_{0}(\varphi)$. This is the function that maps a point on one side of the
magnetic well defined by $B_{0}$ to the other side of the well, i.e.,
$B_{0}(\varphi)=B_{0}[\eta_{0}(\varphi)]$. Symmetry conditions on the
functions $B_{0}$, $d$ etc. then apply with respect to the mapping $\eta_{0}$
(Rodríguez & Plunk, 2023).

With this in mind, we proceed from Eq. (46) without making any assumption of
symmetry.888It can be rigurously shown that such an expression holds even upon
breaking of stellarator symmetry to the order shown. The first term,
corresponding to the $1/r$ order, vanished due to the QI condition and the odd
symmetry of $B_{1}$. It can be shown that this vanishing holds even when
stellarator symmetry is abandoned,

$\int_{\varphi_{L}}^{\varphi_{R}}F(\lambda,B_{0})B_{1}(\varphi)\mathrm{d}\varphi=\int_{\varphi_{R}}^{\varphi_{L}}F(\lambda,B_{0})\underbrace{B_{1}[\eta_{0}(\bar{\varphi})]\eta_{0}^{\prime}(\bar{\varphi})}_{\stackrel{{\scriptstyle\mathrm{QI}}}{{=}}B_{1}(\bar{\varphi})}\mathrm{d}\bar{\varphi}=0,$
(54)

where we changed variables $\varphi=\eta_{0}(\bar{\varphi})$ for the first
equality, applied the first-order QI condition in a generally non-symmetric QI
field [Eq. (19) in Rodríguez & Plunk (2023)], and noted that the second
integral is exactly the negative of the first one. This proves that the QI
property is sufficient to make the integral vanish, without any stellarator
symmetry requirement being necessary. This should not come as a surprise.

The QI condition at second order may be succinctly written as an equation for
the $\alpha$-dependent part of $B_{2}$ as (Eq. (C12) in Rodríguez & Plunk
(2023)),

$B_{2}(\theta,\varphi)-\bar{B}_{2}(\varphi)-\partial_{\varphi}\left(\frac{B_{1}^{2}}{2B_{0}^{\prime}}\right)=\eta^{\prime}_{0}\left[B_{2}(\theta,\varphi)-\bar{B}_{2}(\varphi)-\partial_{\varphi}\left(\frac{B_{1}^{2}}{2B_{0}^{\prime}}\right)\right]_{b},$
(55)

where the subscript $b$ denotes the evaluation of the expression inside the
brackets at $\eta_{0}$ and the barred $\bar{B}_{2}$ represents the poloidal
average of $B_{2}$. Following the same substitution trick as for Eq. (54), the
$\theta$ dependent part of the order $O(1)$ contribution to $\omega_{\alpha}$
can be shown to vanish, only leaving the $\theta$ (or $\alpha$) dependent
piece. The conclusionis that the expression for the omnigeneous part of
$\omega_{\alpha,0}$ is the same whether or not the field is stellarator
symmetric, with the difference that in the asymmetric case usually
$\varphi_{L}\neq-\varphi_{R}$. In addition, symmetry breaking will also affect
the value of $B_{20}$.

##  Appendix B Radial derivative of $B$ for deeply trapped particles in QI
configurations

As we have seen in the main text, the average radial derivative of $B$ at the
minimum along the field, $B_{20}=\bar{B}\partial_{\psi}B/2|_{\psi=0}$, is key
in determining the tangential drift of trapped particles. This term is related
to the magnetic well (30), which plays an important role for MHD stability
(Freidberg, 2014; Landreman & Jorge, 2020; Rodríguez, 2023).

Regardless of whether the field is QI, the term $B_{20}$ must have a form that
is consistent with the equilibrium and solenoidal properties of the magnetic
field. That is, it must be consistent with $\nabla\cdot\mathbf{B}=0$ and
$\mathbf{j}\times\mathbf{B}=\nabla p$. Within the asymptotic framework of the
near-axis description of the field, the function $B_{20}$ is only partially
constrained by the properties of the magnetic field on the axis and to first
order. There is in general some freedom in its choice. However, at ‘straight’
sections where the curvature of the magnetic axis vanishes, the constraints of
the equations become particularly stringent leaving no freedom. $B_{20}$
becomes uniquely determined by lower-order choices and the plasma pressure. In
this Appendix we find determine its value and explore formal consequences of
it.

For this purpose, we focus on the minimum of $|\mathbf{B}|$ along field, where
deeply trapped particles reside and $\kappa=0$. We are interested in relating
the magnetic field magnitude $B_{20}$ to lower-order quantities in the near-
axis expansion and to the flux-surface geometry. To do so, and as part of what
is known as the inverse-coordinate approach to the near-axis expansion (Garren
& Boozer, 1991), we describe flux surfaces in terms of Boozer coordinates (our
independent set of coordinates) respect to a signed Frenet-Serret frame (Plunk
et al., 2019) $\\{\hat{b},\hat{\kappa},\hat{\tau}\\}$ as
$\mathbf{x}-\mathbf{r}_{0}=X\hat{\kappa}+Y\hat{\tau}+Z\hat{b}$. Here
$\mathbf{r}_{0}$ represents the magnetic axis and $\\{X,Y,Z\\}$ are functions
describing the flux-surface shape. The near-axis expansion framework is
designed to uncover the connection between the $|\mathbf{B}|$ and these
features. The relevant equations and constructions can be found in Landreman &
Sengupta (2019), to which we refer for further details. We shall invoke
results from this work, reproducing the necessary expressions. In what
follows, the symbol $f_{ij}$ represents the expansion coefficients of the
function $f$, where the subscript $i$ corresponds to the power of $r$ in the
expansion, and $j$ equals $0,~{}c$ or $s$ depending on the $\theta$ harmonic
it represents (constant, cosine or sine). At higher order generalised notation
can be devised (Rodriguez & Bhattacharjee, 2021).

We first consider the shape of flux surfaces in the direction normal to the
magnetic axis and begin with the function $X_{20}$, which describes a rigid
shift of the flux surfaces relative to each other (akin to the Shafranov shift
(Rodríguez, 2023; Landreman, 2021)). This function can be read off from Eq.
(A34) of Landreman & Sengupta (2019),

$X_{20}=\frac{1}{\kappa\ell^{\prime}}\left\\{Z_{20}^{\prime}-\frac{1}{\ell^{\prime}}\left[-\frac{G_{0}^{2}B_{20}}{B_{0}^{3}}+\frac{3G_{0}^{2}(B_{1c}^{2}+B_{1s}^{2})}{4B_{0}^{4}}+\right.\right.\\\
\left.\left.+\frac{G_{0}(G_{2}+\iota_{0}I_{2})}{B_{0}^{2}}-\frac{X_{1c}^{2}+X_{1s}^{2}}{4}(\kappa\ell^{\prime})^{2}-\frac{q_{c}^{2}+q_{s}^{2}+r_{c}^{2}+r_{s}^{2}}{4}\right]\right\\},$
(56)

where all the expressions on the right hand side except $B_{20},~{}G_{2}$ and
$I_{2}$ are first-order quantities.

Since $\kappa$ appears in the denominator of this expression, to avoid an
unphysical diverging shift of flux surfaces at the point where $\kappa=0$,
$B_{20}$ is highly constrained at this point. If the curvature has a zero of
order $v$, then we expect the first $v-1$ derivatives of $B_{20}$ to be
determined by Eq. (56). Let us momentarily focus on $B_{20}(\varphi=0)$ for
$\kappa(\varphi=0)=0$. We shall assume stellarator symmetry and the
configuration to be QI to first order, so that $\bar{d}=d/\kappa$ is even in
$\varphi$, the (signed) curvature is odd, the torsion even, $B_{0}$ even, and
$\sigma$ odd ($\sigma=(B_{1s}Y_{1s}+B_{1c}Y_{1c})/\bar{B}\kappa$ is a measure
of up-down asymmetry). With this in mind, using the expression for $Z_{20}$
and $q_{i}$ from Eqs. (A24), (A27) and (A37-40) of Landreman & Sengupta
(2019), we find999All functions of $\varphi$ are meant to be evaluated at
$\varphi=0$, but we avoid writing this explicitly for simplicity of notation.

$\frac{B_{20}}{B_{0}}=-\frac{\mu_{0}p_{2}}{B_{0}^{2}}+\frac{1}{4(\ell^{\prime})^{2}}\left[\bar{d}\bar{d}^{\prime\prime}\left(1-\frac{\bar{B}^{2}}{B_{0}^{2}\bar{d}^{4}}\right)+\left(\frac{\bar{B}}{B_{0}\bar{d}}\right)^{2}(\sigma^{\prime})^{2}-\frac{\bar{B}^{2}B_{0}^{\prime\prime}}{B_{0}^{3}\bar{d}^{2}}\right]-\\\
-\frac{1}{4}\left[\left(\frac{\bar{B}\tau_{0}}{B_{0}\bar{d}}\right)^{2}+\left(\frac{\bar{B}\sigma^{\prime}}{B_{0}\ell^{\prime}\bar{d}}+\bar{d}\tau_{0}\right)\right].$
(57)

To further simplify this expression, we eliminate $\sigma^{\prime}$ in favour
of the properties of the axis and first order quantities. Using the near-axis
Riccati $\sigma$-equation (Landreman & Sengupta, 2019, Eq. (A21)), assuming an
exact QI magnetic field at first order $B_{1c}=-B_{0}d\sin\iota_{0}\varphi$
and $B_{1s}=B_{0}d\cos\iota_{0}\varphi$ (Plunk et al., 2019; Rodríguez &
Plunk, 2023),

$\sigma^{\prime}=2\left(\frac{I_{2}}{\bar{B}}-\tau_{0}\right)\frac{G_{0}\bar{d}^{2}}{\bar{B}}.$
(58)

With this and the relation $\ell^{\prime}=G_{0}/B_{0}$ (assuming $G_{0}>0$)
(Landreman & Sengupta, 2019, Eq. (A20)) we obtain

$\frac{B_{20}}{B_{0}}=-\frac{\mu_{0}p_{2}}{B_{0}^{2}}+\frac{1}{4(\ell^{\prime})^{2}}\sum_{i}\mathcal{P}_{i},$
(59)

where

$\displaystyle\mathcal{P}_{\bar{d}\bar{d}^{\prime\prime}}=\bar{d}\bar{d}^{\prime\prime}\left(1-\frac{\bar{B}^{2}}{B_{0}^{2}\bar{d}^{4}}\right),$
(60a)
$\displaystyle\mathcal{P}_{B_{0}^{\prime\prime}}=-\frac{B_{0}^{\prime\prime}}{B_{0}}\frac{\bar{B}^{2}}{B_{0}^{2}\bar{d}^{2}},$
(60b)
$\displaystyle\mathcal{P}_{\tau_{0}^{2}}=(\bar{d}\tau_{0}\ell^{\prime})^{2}\left(3-\frac{\bar{B}^{2}}{B_{0}^{2}\bar{d}^{4}}\right),$
(60c)
$\displaystyle\mathcal{P}_{I_{2}}=-4(\ell^{\prime})^{2}\bar{d}^{2}\frac{\tau_{0}I_{2}}{\bar{B}},$
(60d)

which are the forms used in the main text. It will be convenient to choose our
reference field $\bar{B}$ to equal $B$ at the bottom of the well, $B_{0}$, but
we keep the expressions here general.

An analogous approach to avoid divergence of $X_{2c}$ yields the consistent
value of $B_{2c}$ at the minimum of the magnetic field. In this case, we need
Eq. (A36) from Landreman & Sengupta (2019),

$X_{2c}=\frac{1}{\kappa\ell^{\prime}}\left\\{Z_{2c}^{\prime}+2\bar{\iota}_{0}Z_{2s}-\frac{1}{\ell^{\prime}}\left[-\frac{G_{0}^{2}B_{2c}}{B_{0}^{3}}+\frac{3G_{0}^{2}(B_{1c}^{2}-B_{1s}^{2})}{4B_{0}^{4}}-\right.\right.\\\
\left.\left.-\frac{X_{1c}^{2}-X_{1s}^{2}}{4}(\kappa\ell^{\prime})^{2}-\frac{q_{c}^{2}-q_{s}^{2}+r_{c}^{2}-r_{s}^{2}}{4}\right]\right\\}.$
(61)

Then, as before, we may write for $B_{2c}$,

$\frac{B_{2c}}{B_{0}}=-\frac{1}{4(\ell^{\prime})^{2}}\left[\bar{d}\bar{d}^{\prime\prime}\left(1+\frac{\bar{B}^{2}}{\bar{d}^{4}B_{0}^{2}}\right)+\left(\frac{\bar{B}}{\bar{d}B_{0}}\right)^{2}\frac{B_{0}^{\prime\prime}}{B_{0}}-4(\ell^{\prime})^{2}\bar{d}^{2}\frac{\tau_{0}I_{2}}{\bar{B}}+\right.\\\
\left.+(\bar{d}\tau_{0}\ell^{\prime})^{2}\left(3+\frac{\bar{B}^{2}}{\bar{d}^{4}B_{0}^{2}}\right)\right],$
(62)

with all quantities evaluated at the point $\varphi=0$. For a simple zero, the
equivalent for $B_{2s}$ yields that $B_{2s}(\varphi=0)=0$. This is satisfied
by the parity conditions derived from stellarator symmetry, and thus bring no
additional information to the table. The expressions here presented may be
found derived in a Mathematica script in the Zenodo repository associated with
this paper. We re-emphasise that the form of these expressions assume exact QI
at first order in the near-axis expansion. For more general forms of the
expressions see the comments in Appendix D.

###  B.1 Implications of the QI condition

We showed in Appendix A and learnt from Rodríguez & Plunk (2023) that, in an
omnigeneous field, $B_{2c}$ must have a very particular form in a stellarator
symmetric field. At the same time, from the conditions of equilibrium, we have
also found in Eq. (62) that the value of $B_{2c}$ is constrained at the bottom
of the magnetic well. Thus, in principle, there are two different constraints
that must be simultaneously satisfied in a QI field.

To further look into this issue, we note that the $B_{2c}$ employed in this
section is somewhat different from the QI form of $B_{2c}$ employed elsewhere
in the text. In fact, $B_{2c}$ here (see Landreman & Sengupta (2019)) is
defined respect to the helical angle $\chi=\theta-N\varphi$, where $N$ is
related to the self-linking number of the signed Frenet frame of the axis
(Rodriguez et al., 2022). In contrast, the QI version of $B_{2c}$ (which we
may denote by $B_{2c}^{\mathrm{QI}}$) is defined in Eqs. (42) respect to the
angle $\alpha+\nu$. Therefore,

$\displaystyle B_{2c}=B_{2c}^{\mathrm{QI}}\cos
2(\bar{\iota}\varphi-\nu)-B_{2s}^{\mathrm{QI}}\sin 2(\bar{\iota}\varphi-\nu),$
(63a) $\displaystyle B_{2s}=B_{2s}^{\mathrm{QI}}\cos
2(\bar{\iota}\varphi-\nu)+B_{2c}^{\mathrm{QI}}\sin 2(\bar{\iota}\varphi-\nu),$
(63b)

where $\bar{\iota}=\iota-N$. Now, for a stellarator symmetric QI field at the
bottom of the well, $\nu=\pi/2$ and $\varphi=0$, meaning that
$B_{2c}=-B_{2c}^{\mathrm{QI}}$. Due to the sign flip, the resulting constraint
from Eq. (45d) becomes

$\bar{d}\bar{d}^{\prime\prime}\left(1+\frac{\bar{B}^{2}}{\bar{d}^{4}B_{0}^{2}}\right)+\left(\frac{\bar{B}}{\bar{d}B_{0}}\right)^{2}\frac{B_{0}^{\prime\prime}}{B_{0}}-4(\ell^{\prime})^{2}\bar{d}^{2}\frac{\tau_{0}I_{2}}{\bar{B}}+\\\
+(\bar{d}\tau_{0}\ell^{\prime})^{2}\left(3+\frac{\bar{B}^{2}}{\bar{d}^{4}B_{0}^{2}}\right)=\frac{(\ell^{\prime})^{2}}{B_{0}}\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)_{\varphi=0}^{\prime}.$
(64)

The condition of omnigeneity at the bottom of the trapping well thus
translates into an additional constraint on the zeroth and first order
components that make up the magnetic field. In the limit of a vacuum field
with no secondary minima of $|\mathbf{B}|$, it follows that
$\bar{d}\bar{d}^{\prime\prime}<0$: geometrically, flux surfaces must become
more elongated in the binormal direction away from the minimum. This is a
geometric consequence of the QI conditions.

With this constraint at hand, we may eliminate $\bar{d}^{\prime\prime}$
altogether from the expression we have for $B_{20}$ in Eq. (59). Defining
$f_{0}=B_{20}/B_{0}-(B_{0}^{2}d^{2}/B_{0}^{\prime})^{\prime}/4B_{0}=f_{p}+f_{B_{0}^{\prime\prime}}+f_{I_{2}}+f_{\tau_{0}^{2}}+f_{\mathrm{QI}}$,
the resulting expression is

$\displaystyle f_{p}=-\frac{\mu_{0}p_{2}}{B_{0}^{2}},$ (65a) $\displaystyle
f_{B_{0}^{\prime\prime}}=-\frac{1}{2(l^{\prime})^{2}}\frac{\sqrt{\bar{\alpha}}}{1+\bar{\alpha}}\frac{\bar{B}}{B_{0}}\frac{B_{0}^{\prime\prime}}{B_{0}},$
(65b) $\displaystyle
f_{\tau_{0}^{2}}=\frac{\sqrt{\bar{\alpha}}}{1+\bar{\alpha}}\frac{\bar{B}}{B_{0}}\tau_{0}^{2},$
(65c) $\displaystyle
f_{I_{2}}=-\frac{2\sqrt{\bar{\alpha}}}{1+\bar{\alpha}}\frac{I_{2}}{B_{0}}\tau_{0},$
(65d) $\displaystyle
f_{\mathrm{QI}}=-\frac{1}{2B_{0}}\frac{1}{1+\bar{\alpha}}\left.\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)^{\prime}\right|_{\varphi=0},$
(65e)

where $\bar{\alpha}=\bar{d}^{4}B_{0}^{2}/\bar{B}^{2}$.

###  B.2 Stellarator symmetry breaking

The derivation above may be extended to non-stellarator-symmetric fields,
which leads to additional terms in various equations. For simplicity, we shall
assume once againa that the field is exactly omnigeneous to first order and
that $\nu=\pi/2$. The latter condition is not necessary as one could add an
even function $\tilde{\nu}$ while still satisfying the conditions of
omnigeneity (Rodríguez & Plunk, 2023). In the interest of brevity, however, we
only consider this simplest case.

As a result of the lack of symmetry, the expression for the magnetic well
$B_{20}$ changes. Proceeding in an equivalent way to that above, it can be
shown that, in the expression (59), the coefficients of $B_{20}/B_{0}$ become
(assuming $G_{0}>0$)

$\displaystyle\mathcal{P}_{\bar{d}\bar{d}^{\prime\prime}}$
$\displaystyle=\bar{d}\bar{d}^{\prime\prime}\left(1-\frac{\bar{B}^{2}(1+\sigma^{2})}{B_{0}^{2}\bar{d}^{4}}\right)$
$\displaystyle\mathcal{P}_{\bar{d}^{\prime}}$
$\displaystyle=4\bar{d}^{\prime}\ell^{\prime}\frac{\bar{B}\sigma\tau}{B_{0}\bar{d}}$
$\displaystyle\mathcal{P}_{B_{0}^{\prime\prime}}$
$\displaystyle=-\frac{B_{0}^{\prime\prime}}{B_{0}}\frac{\bar{B}^{2}(1+\sigma^{2})}{B_{0}^{2}\bar{d}^{2}}$
$\displaystyle\mathcal{P}_{(\bar{d}^{\prime})^{2}}$
$\displaystyle=2(\bar{d}^{\prime})^{2}\frac{\bar{B}^{2}(1+\sigma^{2})}{B_{0}^{2}\bar{d}^{4}}$
$\displaystyle\mathcal{P}_{\tau_{0}^{2}}$
$\displaystyle=(\bar{d}\tau_{0}\ell^{\prime})^{2}\left(3-\frac{\bar{B}^{2}(1+\sigma^{2})}{B_{0}^{2}\bar{d}^{4}}\right)$
$\displaystyle\mathcal{P}_{\tau_{0}^{\prime}}$
$\displaystyle=-2\tau_{0}^{\prime}\ell^{\prime}\frac{\bar{B}\sigma}{B_{0}}$
$\displaystyle\mathcal{P}_{I_{2}}$
$\displaystyle=-4(\ell^{\prime})^{2}\bar{d}^{2}\frac{\tau_{0}I_{2}}{\bar{B}}.$
(66)

The main difference to the stellarator symmetric case is the modification of
terms through the breaking of up-down symmetry, $\sigma$, and the appearence
of three new terms. If the field is up-down symmetric at the minimum of
$|\mathbf{B}|$ ($\sigma=0$), the only new term compared with the stellarator-
symmetric case is that proportional to $(\bar{d}^{\prime})^{2}$. This term is
always positive, thus indicating a potential benefit of breaking stellarator-
symmetry to increase the depth of the radial well.

Proceeding similarly with $B_{2c}$, we find

$\frac{B_{2c}}{B_{0}}=-\frac{1}{4(\ell^{\prime})^{2}}\left[\bar{d}\bar{d}^{\prime\prime}\left(1+\frac{\bar{B}^{2}(1-\sigma^{2})}{\bar{d}^{4}B_{0}^{2}}\right)+\left(\frac{\bar{B}}{\bar{d}B_{0}}\right)^{2}\frac{B_{0}^{\prime\prime}}{B_{0}}(1-\sigma^{2})-4(\ell^{\prime})^{2}\bar{d}^{2}\frac{\tau_{0}I_{2}}{\bar{B}}+\right.\\\
\left.+(\bar{d}\tau_{0}\ell^{\prime})^{2}\left(3+\frac{\bar{B}^{2}(1-\sigma^{2})}{\bar{d}^{4}B_{0}^{2}}\right)+4\bar{d}^{\prime}\ell^{\prime}\frac{\bar{B}}{B_{0}\bar{d}}\sigma\tau_{0}-2(\bar{d}^{\prime})^{2}\frac{\bar{B}^{2}}{B_{0}^{2}\bar{d}^{4}}(1-\sigma^{2})-2\tau_{0}^{\prime}\ell^{\prime}\sigma\frac{\bar{B}}{B_{0}}\right]$
(67)

Now, when stellarator symmetry is broken, the condition of QI at second order
changes, and the self-consistency condition for $B_{2c}(0)$ will change
accordingly. That is, we must check the implications of Eq. (55), drawn from
Rodríguez & Plunk (2023) (see Eq. (26b) there). At the bottom of the well it
is always true that $\eta_{0}^{\prime}\approx-1$, by virtue of having a
minimum (and thus in a local Taylor expansion in $\varphi$, a leading even
power of $\varphi$). Because $d$ must vanish at this position rapidly enough
to avoid the loss of confinement of deeply trapped particles (see the
discussion on pseudosymmetry and the behaviour at the bottom of the well in
Rodríguez & Plunk (2023)), $B_{2c}(0)$ must be the same as in the stellarator
symmetric case, at least under the assumption made on the form of $\nu$.

This way, we form the equivalent to Eq. (65) when stellarator symmetry is
broken, which becomes

$\displaystyle f_{p}$ $\displaystyle=-\frac{\mu_{0}p_{2}}{B_{0}^{2}},$
$\displaystyle f_{\bar{d}^{\prime}}=$ $\displaystyle
2\frac{\bar{B}}{B_{0}\bar{d}\ell^{\prime}}\frac{\sigma\tau_{0}}{1+\bar{\alpha}-\sigma^{2}}\bar{d}^{\prime},$
$\displaystyle f_{B_{0}^{\prime\prime}}$
$\displaystyle=-\frac{1}{2(l^{\prime})^{2}}\frac{\sqrt{\bar{\alpha}}}{\bar{\alpha}+1-\sigma^{2}}\frac{\bar{B}}{B_{0}}\frac{B_{0}^{\prime\prime}}{B_{0}},$
$\displaystyle f_{(\bar{d}^{\prime})^{2}}=$
$\displaystyle\frac{1}{1+\bar{\alpha}-\sigma^{2}}\left(\frac{\bar{d}^{\prime}}{\ell^{\prime}}\right)^{2},$
$\displaystyle f_{\tau_{0}^{2}}$
$\displaystyle=\frac{\sqrt{\bar{\alpha}}}{1+\bar{\alpha}-\sigma^{2}}\frac{\bar{B}}{B_{0}}\tau_{0}^{2},$
$\displaystyle f_{\tau_{0}^{\prime}}=$
$\displaystyle-\frac{\tau_{0}^{\prime}}{\ell^{\prime}}\frac{\bar{B}}{B_{0}}\frac{\sigma}{1+\bar{\alpha}-\sigma^{2}},$
$\displaystyle f_{I_{2}}$
$\displaystyle=-\frac{2\sqrt{\bar{\alpha}}}{1+\bar{\alpha}-\sigma^{2}}\frac{I_{2}}{B_{0}}\tau_{0},$
$\displaystyle f_{\mathrm{QI}}=$
$\displaystyle-\frac{1}{2B_{0}}\frac{1}{1+\bar{\alpha}-\sigma^{2}}\left.\left(\frac{B_{0}^{2}d^{2}}{B_{0}^{\prime}}\right)^{\prime}\right|_{\varphi=0}.$
(68)

Some of these expressions are used in the main text.

##  Appendix C The magnetic mirror problem

There is an obvious analogy between magnetic mirrors and quasi-isodynamic
stellarators. Indeed, the latter are sometimes referred to as linked mirrrors.
This analogy is especially accurate near the maxima and minima of
$|\mathbf{B}|$ along field lines, which tend to be located in places where the
curvature is small (and vanishes on the magnetic axis) and where trapped
particles thus have little sense of the toroidicity of the field. This analogy
provides us with a practical tool to assess the behaviour of the QI field from
the perspective of a straight magnetic mirror.

###  C.1 Near-axis expansion procedure

The near-axis approximation (commonly referred to as the paraxial
approximation in the context of magnetic mirrors) is common in the study of
magnetic mirrors (Furth & Rosenbluth, 1964; Kadomtsev & Pogutse, 1967; Catto
et al., 1981; Savenko, 2006). In this appendix, we adopt this approximation
within the inverse-coordinate approach used for stellarators. We do so in an
attempt to simplify the comparison and to clarify the quadrupole ansatz
traditionally used in the analysis of magnetic mirrors.

We focus on the case of a vacuum field and describe the field lines as curves
of constant $\psi$ and $\alpha$, whose position is given by
$\mathbf{x}(\psi,\alpha,z)=X\hat{\mathbf{x}}+Y\hat{\mathbf{y}}+(z+Z)\hat{\mathbf{z}}$
in Cartesian coordinates. The $z$-axis can be thought of as the ‘magnetic
axis’, about which we will perform the expansion, and the function $Z$ can
without loss of generality be chosen to vanish.

In order to find the governing equations, we write the vacuum condition of the
magnetic field and its Clebsch representation as (D’haeseleer et al., 2012)

$\nabla\Phi=\nabla\psi\times\nabla\alpha,$ (69)

or in the inverse-coordinate representation $\\{\psi,\alpha,z\\}$,

$\partial_{z}\mathbf{x}=\partial_{\psi}\Phi\partial_{\alpha}\mathbf{x}\times\partial_{z}\mathbf{x}+\partial_{\alpha}\Phi\partial_{z}\mathbf{x}\times\partial_{\psi}\mathbf{x}+\partial_{z}\Phi\partial_{\psi}\mathbf{x}\times\partial_{\alpha}\mathbf{x}.$
(70)

To relate it to the magnetic field magnitude, we then write what we shall
refer to as the Jacobian equation,

$\left(\frac{\partial_{z}\Phi}{B}\right)^{2}=|\partial_{z}\mathbf{x}|^{2},$
(71)

which follows from taking the scalar product of Eq. (69) with its right-hand-
side, using the dual relation
$\nabla\psi\times\nabla\alpha=\mathcal{J}^{-1}\partial_{z}\mathbf{x}=\mathbf{B}$.

As usual, we need to expand these equations by expressing all functions as
Taylor-Fourier expansions in the pseudo-radial coordinate
$r=\sqrt{2\psi/\bar{B}}$ and $\alpha$. The latter is a field-line label, and
because we are in a straight, non-toroidal system, we may take it as our
poloidal angle. (In the absence of toroidal flux surfaces, there is no need to
introduce a toroidal or a helical angle.) The relevant functions in the
expansion are $B,~{}X$ and $Y$ as well as the scalar potential $\Phi$. We
shall use the same index notation as in the regular stellarator near-axis
expansion, so that, for instance,

$X(\psi,\alpha,z)=r\left[X_{1s}(z)\sin\alpha+X_{1c}(z)\cos\alpha\right]+\\\
+r^{2}\left[X_{20}(z)+X_{2s}(z)\sin 2\alpha+X_{2c}(z)\cos
2\alpha\right]+\cdots$ (72)

Here, the functions $X_{1c}(z)$ and $X_{1s}(z)$ are considered to be first-
order inputs to the problem alongside $B_{0}(z)$, the magnetic field strength
on axis.

####  C.1.1 Zeroth order (71)

To leading order $O(r^{0})$ of Eq. (71),

$\Phi_{0}=\int B_{0}(z)\mathrm{d}z,$ (73)

where $B_{0}$ is taken to be positive.

####  C.1.2 Leading order (70)

To leading order, $O(r^{-1})$, there are three different equations we obtain
from the three projections along the Cartesian basis of Eq. (70). From the
$\hat{z}$-projection we have

$X_{1c}Y_{1s}-X_{1s}Y_{1c}=\frac{\bar{B}}{B_{0}(z)},$ (74)

which is a statement of magnetic flux conservation along the flux tube
labelled by $r$. This is of course fully analogous to the situation in a
stellarator, as it is a property inherent to any flux tube.

For a non-vanishing magnetic field $B_{0}\neq 0$, the other two components of
Eq. (70) to $O(r)$ can be combined linearly using $X_{1c}$ and $X_{1s}$. The
result is that $\Phi_{1,c}=\Phi_{1,s}=0$. There is, in other words, no first-
order contribution to the vacuum magnetic potential.

Next, the first-order $O(r)$ form of Eq. (71) yields the natural conclusion
that $B_{1s}=0$ and $B_{1c}=0$, so that there is no $\alpha$-variation of the
magnetic field strength to first order. This is an expected result, as we well
know that, in the absence of a pressure gradient, the perpendicular gradient
of $|\mathbf{B}|$ is proportional to the curvature vector. In the vicinity of
a straight axis, there is therefore no first-order variation of the magnetic
field strength. This circumstance sets the straight mirror apart from a
typical QI stellarator, in which the toroidal nature of the field forces the
field to exhibit such a variation in most places along the magnetic axis. The
leading form of the field is then that at next order, which must be
quadrupole-like.

####  C.1.3 Zeroth order (70)

Let us now turn to the next order of Eq. (70), which is $O(r^{0})$. At this
order we have a total of six independent equations, as each of the projections
has $\cos\alpha$ and $\sin\alpha$ terms.

From the $\hat{z}$ projection, and constructing appropriate linear
combinations of the equations, we may write the following equations for
$Y_{2c}$ and $Y_{2s}$,

$Y_{2c}=\frac{X_{1c}^{2}-X_{1s}^{2}}{X_{1c}^{2}+X_{1s}^{2}}Y_{20}+\frac{X_{1c}Y_{1s}-X_{1s}Y_{1c}}{X_{1c}^{2}+X_{1s}^{2}}X_{2s}+\frac{X_{1c}Y_{1c}+X_{1s}Y_{1s}}{X_{1c}^{2}+X_{1s}^{2}}X_{2c}+\\\
+\frac{X_{1s}Y_{1s}-X_{1c}Y_{1c}}{X_{1c}^{2}+X_{1s}^{2}}X_{20},$ (75a)
$Y_{2c}=\frac{2X_{1c}X_{1s}}{X_{1c}^{2}+X_{1s}^{2}}Y_{20}+\frac{X_{1c}Y_{1c}+X_{1s}Y_{1s}}{X_{1c}^{2}+X_{1s}^{2}}X_{2s}+\frac{X_{1s}Y_{1c}-X_{1c}Y_{1s}}{X_{1c}^{2}+X_{1s}^{2}}X_{2c}-\\\
-\frac{X_{1s}Y_{1c}+X_{1c}Y_{1s}}{X_{1c}^{2}+X_{1s}^{2}}X_{20},$ (75b)

which are thus explicitly related to the second-order quantities $X_{2}$ and
$Y_{20}$. This circumstance is reminiscent of the situation in a stellarator
(Eqs. (A32)-(A33) in Landreman & Sengupta (2019)). In the interest of brevity,
we shall not proceed to find expressions for these coefficients. Full
expressions may be found in the Zenodo repository associated with this paper.

The remaining components of Eq. (70) can also be combined to yield expressions
for $\Phi_{2}$ (i.e., three equations) and a consistency condition on first-
order quantities. The former may be written as follows, using Eq. (74),

$\displaystyle\Phi_{20}=$
$\displaystyle\frac{B_{0}}{4}\frac{(X_{1c}^{2}+X_{1s}^{2})(Y_{1s}X_{1c}^{\prime}+Y_{1c}X_{1s}^{\prime})+(Y_{1c}^{2}+Y_{1s}^{2})(X_{1s}Y_{1c}^{\prime}+X_{1c}Y_{1s}^{\prime})}{X_{1s}Y_{1c}+X_{1c}Y_{1s}},$
(76a) $\displaystyle\Phi_{2c}=$
$\displaystyle\frac{B_{0}}{4}\frac{(X_{1c}^{2}-X_{1s}^{2})(Y_{1s}X_{1c}^{\prime}+Y_{1c}X_{1s}^{\prime})+(Y_{1c}^{2}-Y_{1s}^{2})(X_{1s}Y_{1c}^{\prime}+X_{1c}Y_{1s}^{\prime})}{X_{1s}Y_{1c}+X_{1c}Y_{1s}},$
(76b) $\displaystyle\Phi_{2s}=$
$\displaystyle\frac{B_{0}}{4}\frac{X_{1c}Y_{1c}[(X_{1s}^{2})^{\prime}+(Y_{1s}^{2})^{\prime}])+X_{1s}Y_{1s}[(X_{1c}^{2})^{\prime}+(Y_{1c}^{2})^{\prime}])}{X_{1s}Y_{1c}+X_{1c}Y_{1s}}+\frac{\bar{B}}{2}\frac{X_{1c}X_{1s}^{\prime}-X_{1s}X_{1c}^{\prime}}{X_{1s}Y_{1c}+X_{1c}Y_{1s}},$
(76c)

and the latter as

$(X_{1c}^{2}-X_{1s}^{2})(X_{1s}X_{1c}^{\prime}-X_{1c}X_{1s}^{\prime}+Y_{1s}Y_{1c}^{\prime}-Y_{1c}Y_{1s}^{\prime})=0.$
(77)

Because we have adopted an approach in which $X_{1s}$ and $X_{1c}$ are inputs
to the near-axis construction, the second bracket must vanish. This is an
ordinary differential equation for $Y_{1}$, which together with Eq. (74) give
both $Y_{1s}$ and $Y_{1c}$ in terms of $X_{1}$ and $B_{0}$. Introducing the
definition $Y_{1c}=Y_{1s}\sigma$, this equation can be cast, with Eq. (74),
into the form of a Riccati equation,

$\sigma^{\prime}=-(X_{1c}-X_{1s}\sigma)^{2}\left(\frac{B_{0}}{\bar{B}}\right)^{2}\left(X_{1s}X_{1c}^{\prime}-X_{1c}X_{1s}^{\prime}\right),$
(78)

while Eq. (74) becomes

$Y_{1s}=\frac{\bar{B}}{B_{0}(X_{1c}-\sigma X_{1s})}.$ (79)

####  C.1.4 Second order (71)

Having closed forms for the various terms in $\Phi_{2}$, we are in a position
to evaluate the magnetic field $B_{2}$ by considering the $O(r^{2})$ expansion
of Eq. (71). From each of its harmonics the following relations can be read
off

$\displaystyle\frac{B_{20}}{B_{0}}=$
$\displaystyle\frac{\Phi_{20}^{\prime}}{B_{0}}-\frac{1}{4}\left[(X_{1c}^{\prime})^{2}+(X_{1s}^{\prime})^{2}+(Y_{1c}^{\prime})^{2}+(Y_{1s}^{\prime})^{2}\right],$
(80a) $\displaystyle\frac{B_{2c}}{B_{0}}=$
$\displaystyle\frac{\Phi_{2c}^{\prime}}{B_{0}}-\frac{1}{4}\left[(X_{1c}^{\prime})^{2}-(X_{1s}^{\prime})^{2}+(Y_{1c}^{\prime})^{2}-(Y_{1s}^{\prime})^{2}\right],$
(80b) $\displaystyle\frac{B_{2s}}{B_{0}}=$
$\displaystyle\frac{\Phi_{2s}^{\prime}}{B_{0}}-\frac{1}{2}\left(X_{1c}^{\prime}X_{1s}^{\prime}+Y_{1c}^{\prime}Y_{1s}^{\prime}\right),$
(80c)

which uniquely describe the magnetic field magnitude to second order in terms
of first-order shaping. Note that the quadrupole nature of $|\mathbf{B}|$
arises naturally through the near-axis construction.

###  C.2 Magnetic field at the bottom of the well

We now have all the tools we need to understand the behaviour of the field at
the bottom of the magnetic well, which we define to be at $z=0$. We choose to
align our coordinate system with the leading elliptical shape at $z=0$ and
define the poloidal label $\alpha$ in such a way that $\sigma(0)=0$ and
$X_{1s}(0)=0$, without any loss of generality.101010We could have kept
$\sigma(0)$ non-zero in order to mirror the stellarator case, but in the
interest of reducing the algebra we choose it to vanish. The reason why in the
stellarator case this choice is not general is that the basis used is the
Frenet-Serret frame, which we are not free to redefine. By definition,
$B_{0}^{\prime}(0)=0$ and we shall choose our reference field as
$B_{0}(0)=\bar{B}$.

With these definitions and choices, we may then evaluate the magnetic field at
the bottom of the magnetic well for a general straight magnetic mirror. The
result is

$\displaystyle
B_{20}(z=0)=\frac{B_{0}}{4}\left\\{-\frac{B_{0}^{\prime\prime}}{B_{0}X_{1c}^{2}}+2\left[\left(\frac{X_{1c}^{\prime}}{X_{1c}^{2}}\right)^{2}+(X_{1s}^{\prime})^{2}\right]+X_{1c}^{\prime\prime}X_{1c}\left(1-\frac{1}{X_{1c}^{4}}\right)\right\\},$
(81a) $\displaystyle
B_{2c}(z=0)=\frac{B_{0}}{4}\left\\{\frac{B_{0}^{\prime\prime}}{B_{0}X_{1c}^{2}}-2\left[\left(\frac{X_{1c}^{\prime}}{X_{1c}^{2}}\right)^{2}+(X_{1s}^{\prime})^{2}\right]+X_{1c}^{\prime\prime}X_{1c}\left(1+\frac{1}{X_{1c}^{4}}\right)\right\\},$
(81b) $\displaystyle B_{2s}(z=0)=\frac{B_{0}}{2}X_{1c}X_{1s}^{\prime\prime},$
(81c)

where all quantities are evaluated at $z=0$.

###  C.3 Maximum-$J$ in a magnetic mirror

The analysis of the maximum-$J$ condition in a mirror closely follows that in
a QI stellarator and proceeds from the general expression for the precession
frequency (37). In the context of a straight mirror, relaxing the conditions
of omnigeneity and stellarator symmetry, but instead noting that there is no
first order variation in the field magnitude, $B_{1}$, we find in vacuum

$\omega_{\alpha}=\frac{2mv^{2}}{q\bar{B}}\frac{\int_{z_{L}}^{z_{R}}\frac{1-\lambda
B_{0}/2}{B_{0}\sqrt{1-\lambda
B_{0}}}B_{2}(z,\alpha)\mathrm{d}z}{\int_{z_{L}}^{z_{R}}\frac{\mathrm{d}z}{\sqrt{1-\lambda
B_{0}}}},$ (82)

where $B_{2}=B_{20}+B_{2c}\cos 2\alpha+B_{2s}\sin 2\alpha$, $z_{L,R}$ are the
bouncing points for given $\lambda$ and the magnetic field functions are given
in Eqs. (80). We may normalise $\omega_{\alpha}$ following the prescription in
Eq. (4), for which we need to define an edge flux $\psi_{a}$. Following the
same procedure as for the near-axis stellarator, defining the aspect ratio as
$A=L/a$, where $L$ is the length of the mirror,

$\psi_{a}=\frac{1}{2}\frac{L^{3}}{A^{2}}\left(\int\frac{\mathrm{d}z}{B_{0}}\right)^{-1}.$
(83)

At the bottom of the well, the expression for $\omega_{\alpha}$ for deeply
trapped particles reduces to a form analogous to Eq. (9),

$\omega_{\alpha}\approx\frac{mv^{2}}{qB_{0}^{2}}\left(B_{20}+B_{2c}\cos
2\alpha+B_{2s}\sin 2\alpha\right).$ (84)

For maximum-$J$ behaviour, i.e. $q\omega_{\alpha}>0$ for all $\alpha$, it is
thus necessary that

$B_{20}^{2}>B_{2c}^{2}+B_{2s}^{2},$ (85)

in addition to $B_{20}>0$. Because the mirror has no first-order variation of
$|\mathbf{B}|$, the maximum-$J$ and minimum-$B$ properties are equivalent,
i.e. both require $\partial_{\psi}B>0$ to leading order in the expansion in
$r$. This minimum-$B$ property is well known to endow the magnetic mirror with
MHD stability at sufficiently low $\beta$ (Berkowitz et al., 1958; Taylor,