2. Near-axis expansion

The properties of a magnetohydrodynamic (MHD) equilibrium can be described in the vicinity of the magnetic axis by performing a Taylor expansion in the parameter $\epsilon = a/R$, the inverse aspect ratio, where $a$ and $R$ are the minor and major radius, respectively. This procedure, using Boozer coordinates ($\psi,\theta,\varphi$), was shown in Garren & Boozer (Reference Garren and Boozer1991) to allow for the numerical construction of quasi-symmetric configurations, and was extended to consider QI configurations at first order in Plunk et al. (Reference Plunk, Landreman and Helander2019). We restrict our study to the case of first-order, QI, stellarator-symmetric equilibria as described in Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022).

This method allows us to directly construct QI equilibria at first order by prescribing the shape of a magnetic axis $x_0$, the intensity of the magnetic field on the axis $B_0$, the number of field periods $N$, two toroidal geometric functions $\alpha (\varphi )$ and $d(\varphi )$, which are later discussed in more detail, and a distance from the axis at which we desire to construct the equilibrium. In this formalism, the magnetic field intensity is given by

(2.1)\begin{equation} B(\epsilon, \theta, \varphi) \approx B_{0}(\varphi) + \epsilon B_{1}(\theta,\varphi)= B_{0}(\varphi) \left(1 + \epsilon d(\varphi) \cos[\theta - \alpha(\varphi) ] \right) \end{equation}

and the plasma boundary is described by

(2.2)\begin{equation} \boldsymbol{x} \approx \boldsymbol{r}_{0} + \epsilon \left( X_{1}\boldsymbol{n}^{s} + Y_{1}\boldsymbol{b}^{s} \right), \end{equation}

with

(2.3)\begin{gather} X_{1} = \frac{d(\varphi)}{\kappa^{s}} \cos{[\theta - \alpha (\varphi) ]}, \end{gather}

(2.4)\begin{gather}Y_{1} = \frac{2\kappa^{s}}{B_{0}(\varphi)d(\varphi)} \left( \sin{[\theta - \alpha (\varphi)]} + \sigma(\varphi) \cos{[\theta - \alpha (\varphi) ]} \right). \end{gather}

Here $\kappa$ and $\tau$ are the curvature and torsion of the magnetic axis and the quantity $\sigma (\varphi )$ can be found by solving

(2.5)\begin{equation} \sigma{'} + (\iotaslash -\alpha{'})\left(\sigma^{2}+ 1 +\frac{B_{0}^{2}{\bar{d}}^{4}}{4} \right)-G_{0}\bar{d}^{2}\left(\tau + I_{2}/2\right) = 0 \end{equation}

self-consistently with the rotational transform on axis $\iotaslash$. Here, $\bar {d}(\varphi ) = d(\varphi ) /\kappa ^{s}(\varphi )$, primes denote derivatives with respect to the toroidal Boozer angle $\varphi$ and $G_{0}$ and $I_2$ are related to the first non-zero terms of the poloidal and toroidal current functions on axis; $I_2$ can be assumed to be zero for QI configurations.

The function $d(\varphi )$ needs to be specified as an input for the construction, but the choice is not entirely free; a set of conditions need to be fulfilled to be consistent with omnigenity and stellarator symmetry. It must vanish at all extrema of the on-axis magnetic field $B_0$:

(2.6)\begin{equation} d(\varphi_{\mathrm{min}})=d(\varphi_{\mathrm{max}})=0. \end{equation}

Given how the definition of the boundary depends directly on the quantity $\bar {d}$ ((2.3) and (2.4)), the curvature must necessarily have zeros at the same toroidal locations and of the same order as $d$, to ensure a smooth plasma boundary. Additionally $d$ must be an odd function about $\varphi _{\mathrm {min}}$, which implies this function will always have derivatives of odd order equal to zero at such points. Then $\kappa$ must also be zero up to odd order. The simplest case corresponds to both functions having first-order zeros at the minima of $B_0$.

Finally, the last ingredient for constructing a near-axis QI equilibrium is the function $\alpha (\varphi )$. Just like $d$, it must be specified in a way consistent with omnigenity, albeit as shown in Plunk et al. (Reference Plunk, Landreman and Helander2019) this condition being not compatible with a periodic plasma boundary. For this reason, omnigenity must be broken around the maxima of $B_0$. The method we use is that described in § 4 of Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022), with

(2.7)\begin{equation} \alpha(\varphi) = \iotaslash (\varphi-\varphi_{\mathrm{min}}^{i})+{\rm \pi} \left(2 m i + \frac{1}{2}\right) + {\rm \pi}\left( m - \iotaslash /N \right) \left(\frac{\varphi-\varphi_{\mathrm{min} }^{i}}{{\rm \pi}/N} \right)^{2k+1}, \end{equation}

where $m$ is the axis helicity, to which the next section is devoted, and $k$ is an integer that controls the size of the region in which the omnigenity conditions are fulfilled; increasing $k$ results in a larger toroidal domain in which $\alpha$ is omnigenous but at the cost of a sharp behaviour around the maxima to maintain a periodic behaviour. The superscript $i$ is an integer number labelling the minima of $B_0$; given how we are just considering a minimum per field period, this is equivalent to the field-period label.

3. Helicity

The geometric properties of a three-dimensional space curve $\boldsymbol {r}(\ell )$ can be described by the Frenet–Serret vectors $(\boldsymbol {t},\boldsymbol {n},\boldsymbol {b})$. These orthogonal unitary vectors, named tangent, normal and binormal, respectively, are defined by

(3.1)\begin{equation} \left. \begin{gathered} \boldsymbol{t}(\ell) = \frac{\boldsymbol{r}{'}(\ell)}{\|\boldsymbol{r}{'}(\ell)\|},\\ \boldsymbol{n}(\ell) = \frac{\boldsymbol{t}{'}(\ell)}{\|\boldsymbol{t}{'}(\ell)\|},\\ \boldsymbol{b}(\ell) = \boldsymbol{t}(\ell) \times \boldsymbol{n}(\ell), \end{gathered} \right\} \end{equation}

and are related to the curvature and torsion ($\kappa, \tau$) of the curve through the Frenet–Serret formulas:

(3.2)\begin{equation} \left. \begin{gathered} \frac{{\rm d}\boldsymbol{t}}{{\rm d}\ell} = \|\boldsymbol{r}{'}(\ell)\| \kappa \boldsymbol{n} ,\\ \frac{{\rm d}\boldsymbol{n}}{{\rm d}\ell} ={-} \|\boldsymbol{r}{'}(\ell)\|\left[ \kappa\boldsymbol{t} + \tau\boldsymbol{b} \right],\\ \frac{{\rm d}\boldsymbol{b}}{{\rm d}\ell}={-} \|\boldsymbol{r}{'}(\ell)\|\tau \boldsymbol{n} . \end{gathered} \right\} \end{equation}

This frame plays an important role in the construction of exactly omnigenous stellarators at first order around the magnetic axis as described in Landreman & Sengupta (Reference Landreman and Sengupta2018) and Plunk et al. (Reference Plunk, Landreman and Helander2019), given that it is used to describe the plasma boundary as seen in (2.2).

One interesting property of a magnetic axis described using this apparatus is the so-called helicity, $M$, the number of times the normal vector $\boldsymbol {n(\ell )}$ rotates poloidally around the axis after one toroidal turn. In this work we often refer to the per-field-period helicity $m=M/N$, for example, in (2.7), the definition of the $\alpha$ function. In the case of quasi-symmetric configurations, the helicity effectively divides the space of solutions as shown in Landreman & Sengupta (Reference Landreman and Sengupta2018) and Rodriguez et al. (Reference Rodriguez, Sengupta and Bhattacharjee2022): a value of $M=0$ corresponds to axisymmetric configurations while an integer value of $M\neq 0$ corresponds to quasi-helically symmetric solutions.

For QI configurations, which are the main focus of this work, calculating helicity requires extra considerations. As shown in the previous section, points of zero curvature are required in magnetic axes consistent with quasi-isodynamicity in the near-axis formalism. But we can see from the first equation in (3.2) that $\kappa = 0$ results in an undefined normal vector. In the presence of points of first-order zero curvature,Footnote ¹ the Frenet frame undergoes a jump-like rotation which makes the frame, specifically the normal and binormal vectors, non-continuous (Aicardi Reference Aicardi2000), as can be seen in figure 1(a). This discontinuity can be alleviated by using the signed curvature $\kappa ^{s}= s \kappa$, where $s = \pm 1$ changes at each point of odd-order zeros of curvature. In this new frame (Carroll, Köse & Sterling Reference Carroll, Köse and Sterling2013), the sign of the traditional Frenet normal and binormal vectors $\boldsymbol {n}(\ell ),\boldsymbol {t}(\ell )$ also changes when $\kappa = 0$ resulting in a continuous frame as shown in figure 1(b). Helicity for QI configurations is hence calculated using this modified Frenet frame to avoid any discontinuities in the plasma boundary shape.

Figure 1. (a) Unsigned Frenet frame and (b) signed Frenet frame for a two-field-period axis. Side view and top view are shown. Normal $\hat {n}$ and binormal $\hat {b}$ vectors are shown in red and blue, respectively. We can see that the discontinuities in the frame are alleviated when using the signed curvature, $\kappa ^s$.

As in the case of quasi-symmetry, helicity also effectively divides the space of QI near-axis solutions as can be seen in figure 2. Here, the helicity value for 1156 one-field-period, near-axis configurations is shown, each for a different axis shape while keeping the rest of the near-axis parameters constant. All these configurations have the same axis radial component

(3.3)\begin{equation} R(\phi) = 1 - 0.2 \cos (2\phi), \end{equation}

which ensures two zeros of curvature are present at the minima and maxima of the magnetic field on axis. The vertical component of the axis shape is different for each configurations and chosen as

(3.4)\begin{equation} z(\phi) = z_s(1)\sin(\phi) + z_s(2)\sin(2\phi). \end{equation}

Each square in figure 2 corresponds to a configuration with $z_s(1)$ ranging from 0 to 0.8 and $z_s(2)$ from $-$0.425 to 0.4, constructed using the near-axis method described in § 2. The aspect ratio is set to $A=20$, the function $\alpha (\varphi )$ is described by (2.7) with $k=2$, the magnetic field on axis is $B_0 = 1+0.15\cos (\varphi )$ and $d = 0.73 \kappa ^s$.

Figure 2. Space of one-field-period near-axis QI configurations and their helicity values. Each square corresponds to a near-axis QI configuration and its helicity value, green for $m=1$ and blue $m=0$. The solid lines show the transition between helicity regions and correspond to second-order zeros of curvature at $\varphi _{\mathrm {min}}$ (upper line) and $\varphi _{\mathrm {max}}$ (lower line).

The solid black lines in figure 2, bounding each helicity region, can be obtained analytically and correspond to the values of $z_s(1){-}z_s(2)$ that result in second-order zeros of curvature, i.e. ${\textrm {d}\kappa }/{\textrm {d}\phi }|_{\varphi }=0$, at $\varphi _{\mathrm {min}}$ (upper line) and $\varphi _{\mathrm {max}}$ (lower line). The analytical derivation is shown in Appendix A.

It is clear from figure 2 that the space of QI solutions is divided into subregions by the helicity of the magnetic axis. But perhaps more surprising is the fact that this division is also present when calculating the effective ripple of each configuration in this space, as shown in figure 3. Here, each of the plasma boundary shapes of the configurations in figure 2 is used to solve for a MHD equilibrium using the VMEC code (Hirshman & Whitson Reference Hirshman and Whitson1983), and $\epsilon _{\mathrm {eff}}$, the effective ripple, is calculated as described in Drevlak et al. (Reference Drevlak, Heyn, Kalyuzhnyj, Kasilov, Kernbichler, Monticello, Nemov, Nührenberg and Reiman2003). The maximum value of this quantity is shown colour-coded in the figure. The white squares correspond to configurations for which VMEC could not find an equilibrium consistent with the provided plasma boundary. Once again, the division between regions of helicity is observed, and consistently lower values of $\epsilon _{\mathrm {eff}}$ and rotational transform $\iotaslash$ are found for zero-helicity configurations. The intensity of $B$ on the boundary for a configuration from each region is also shown, but constructed at a smaller aspect ratio $A=10$, to exemplify these differences. On the left-hand side, one corresponding to an $m=1$ axis, with a maximum $\epsilon _{\mathrm {eff}}=1.1\,\%$ and rotational transform $\iotaslash = 0.4351$; while the configuration with $m=0$ has a substantially lower $\epsilon _{\mathrm {eff}}=0.35\,\%$ and $\iotaslash = 0.0587$.

Figure 3. Space of one-field-period, $A=20$, near-axis QI configurations. Parameters $z_s(1)$ and $z_s(2)$ correspond to the parameters of the cylindrical Fourier representation of the axis (3.4). Each square is a near-axis configuration constructed with the parameters described in § 3. Colour coded is the maximum value of the effective ripple ${\epsilon _{\mathrm {eff}}}{}$. Black solid lines correspond to the transition between regions with helicity 0 and 1. One configuration of each region is shown as an example of the properties in that space; the last closed surface and the intensity of the magnetic field on the boundary are shown. Zero-helicity configurations have lower ${\epsilon _{\mathrm {eff}}}{}$ but also lower rotational transform than their $m=1$ counterparts.

Configurations with $m=1$ and remarkable low values of $\epsilon _{\mathrm {eff}}$ can be found as shown in Jorge et al. (Reference Jorge, Plunk, Drevlak, Landreman, Lobsien, Camacho Mata and Helander2022) for one field period, but it requires a significant effort in choosing the initial parameters. The lower values of effective ripple observed in the zero-helicity configurations suggest that less effort might be necessary in this region of the parameter space.

Helicity seems to play an important role in the near-axis QI space. Although the integer-helicity space seems to correlate low rotational transform and low $\epsilon _{\mathrm {eff}}$, this occurs because helicity is directly related to the integrated torsion of the magnetic axis, which is an important mechanism for rotational transform generation as discussed in Helander (Reference Helander2014). A more in-depth discussion about the relation between torsion and rotational transform in the presence of a discontinuous Frenet frame, without using the signed curvature correction, can be found in Pfefferlé et al. (Reference Pfefferlé, Gunderson, Hudson and Noakes2018). Optimized examples of QI configurations with low $\epsilon _{\mathrm {eff}}$ and higher per-field-period rotational transform exist. The exploration of such examples may shed light light on novel ways of generating near-axis equilibria with higher $N$ and small $\epsilon _\textrm {eff}$.

4. Half-helicity

In the previous section, we restricted ourselves to the case of integer helicity, implying that the only case with less than complete rotation of the normal is the zero-helicity case. However, it is also possible to construct curves with fractional values of per-field-period helicity. In fact, such behaviour is exhibited by the magnetic axis of QIPC (Subbotin et al. Reference Subbotin, Mikhailov, Shafranov, Isaev, Nührenberg, Nührenberg, Zille, Nemov, Kasilov and Kalyuzhnyj2006), a QI stellarator found through numerical optimization. This configuration has poloidally closed contours of the magnetic field close to the axis and an effective ripple profile that increases radially, making it more similar to a near-axis QI configuration than W7-X.

Although QIPC's magnetic axis does not have points of exactly zero curvature at the extrema of $B_0$, the curvature does reach relatively small values, and its signed Frenet frame behaves in a similar way to near-axis configurations, having a flip at points where the magnetic field on axis has its minimum. Figure 4 shows this axis together with the normal and binormal vectors associated with its signed Frenet frame, the change in sign being applied at the position of the minima of $B_0$. As we can see from the figure, the normal vector does not complete a whole rotation around the axis each field period;instead the normal vector seems to complete exactly a half-rotation, ending at the opposite direction to that which it started at the beginning of the period. We note that all QI designs found in the past by numerical optimization share this half-rotation feature ($m=1/2$). However, this case has not been considered when constructing near-axis QI configurations in the literature. As we show, it allows us to access a region in the solution space where good confinement can be achieved with a relatively high number of field periods.

Figure 4. Unsigned and signed Frenet frame for the magnetic axis of QIPC, a six-field-period QI optimized configuration (Subbotin et al. Reference Subbotin, Mikhailov, Shafranov, Isaev, Nührenberg, Nührenberg, Zille, Nemov, Kasilov and Kalyuzhnyj2006). Green circles indicate the toroidal locations of the maxima of $B_0$, marking the beginning/end of a magnetic field period. Green diamonds indicate the minima of $B_0$. Normal and binormal Frenet vectors are shown in red and blue lines, respectively. A discontinuity in the frame at minima of $B_0$ is evident in the signed frame from the direction change of the normal vector. This discontinuity is alleviated when using the signed frame and is transferred to the maxima of $B_0$.

Let us now discuss the particularities of constructing a magnetic axis with $m=1/2$. We know a magnetic axis consistent with omnigenity needs to have points of zero curvature at all minima and maxima of the magnetic field (Plunk et al. Reference Plunk, Landreman and Helander2019), which requires calculating helicity in the signed Frenet frame, since the axis helicity is otherwise not well defined. Additionally, for a curve to have half-helicity, we require the normal vector to perform half a rotation per field period. This situation is not consistent with first-order zeros of curvature at both $\varphi _{\mathrm {min}}$ and $\varphi _{\mathrm {max}}$, at least not for the case of analytic stellarator-symmetric axis curves, which can be seen by re-examining the expansion near zero-curvature points, as follows.

As shown in § 3 of Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022) and discussed after (2.6), omnigenity conditions imposed in the NAE for stellarator-symmetric configurations require odd-order zeros of $\kappa$ at the minima of $B_0$. If the axis is described as a Fourier series on the cylindrical $R$ and $z$ coordinates,

(4.1)\begin{gather} R(\phi) = \sum_{n=0}^{n_\mathrm{max}} R_c(n) \cos(n N \phi), \end{gather}

(4.2)\begin{gather}z(\phi) = \sum_{n=1}^{n_\mathrm{max}} z_s(n) \sin(n N \phi). \end{gather}

Then having a point of zero curvature at first order requires $({\textrm {d}^{2}x}/{\textrm {d}\phi ^2})\cdot \hat {R} = 0$ at such a point, as shown in Appendix A (appendix II in Camacho Mata et al. Reference Camacho Mata, Plunk and Jorge2022), which given the definition of $\boldsymbol {n}$ (3.1), is equivalent to the condition

(4.3)\begin{equation} \boldsymbol{n}\boldsymbol{\cdot} \hat{R}=0, \end{equation}

at $\varphi = \varphi _{\mathrm {min}}$. The case $m=1/2$ corresponds to the normal vector performing a ${\rm \pi}$-rotation in the $\hat {R}$–$\hat {z}$ plane per field period. Hence, after half a period, at $\varphi = \varphi _{\mathrm {max}}$, a ${{\rm \pi} }/{2}$-rotation would require the normal vector to point in the radial direction:

(4.4)\begin{equation} \boldsymbol{n} = \hat{R}. \end{equation}

As shown above and explained in greater detail in Appendix A, odd-order zeros of curvature are generally inconsistent with a normal vector having a component in the radial direction, but consistent with even-order zeros. We can intuitively visualize each increase in the order of the zeros of curvature as introducing a flip on the Frenet frame, hence the need for different behaviour of $\kappa$ at the minima and maxima of $B_0$ to preserve a discontinuity of the frame that is not corrected by the sign change described in the previous section. We can then construct half-helicity closed curves, just by requiring a first-order zero of curvature at $\varphi = \varphi _{\mathrm {min}}$ and a second-order zero at $\varphi = \varphi _{\mathrm {max}}$. We note that these arguments rely on continuity of the derivatives of the curve, which can safely be assumed for curves constructed by Fourier series, but need not be the case more generally, as with curves found by solving Frenet–Serret equations directly. We also remark that, when assuming stellarator symmetry as is done in this work, the only helicities possible are either integer or half-helicity, since this is the only fractional helicity allowing for points of zero curvature and continuity of the plasma boundary. Different fractional helicities might be possible for certain non-stellarator-symmetric configurations, but this issue needs to be explored further.

Following the same procedure described in Appendix A, we perform a local expansion of a stellarator-symmetric magnetic axis described by (4.1) and (4.2), and impose conditions on the curvature of such an axis, which we can transform into conditions for the Fourier coefficients $R_{c}(n)$ and $z_{s}(n)$. The local condition for having first-order zeros of curvature is

(4.5)\begin{equation} R_{2} = R_{0}, \end{equation}

where $R_n$ and $z_n$ denote the $n$th derivatives with respect to $\phi$. For a second-order zero of $\kappa$, we additionally need to fulfil

(4.6)\begin{equation} z_{3} = 2z_{1}. \end{equation}

We can apply the previous conditions to a truncated Fourier representation in order to obtain conditions on its coefficients that can be used to construct curves with the curvature properties required to obtain a half-helicity behaviour. We can consider, for example, the family of curves

(4.7)\begin{equation} \left. \begin{gathered} R = 1 + R_c(1) \cos(N \phi) + R_c(2) \cos(2 N \phi) + R_c(3) \cos(3 N \phi),\\ z = z_s(1) \sin(N \phi) + z_s(2) \sin(2 N \phi) + z_s(3) \sin(3 N \phi).\end{gathered} \right\} \end{equation}

Applying condition (4.5) at $\phi =2{\rm \pi} /N$ and both the previous condition and (4.6) to $\phi =0$ returns the following restrictions on the Fourier coefficients necessary to construct curves with half-helicity per field period:

(4.8)\begin{gather} R_c(1) ={-}\frac{R_c(3)+9N^2R_c(3)}{1+N^2}, \end{gather}

(4.9)\begin{gather}R_c(2) ={-}\frac{1}{4 N^2+1}, \end{gather}

(4.10)\begin{gather}z_s(1) ={-}\frac{4z_s(2)+8N^2z_s(2)+6z_s(3)+27N^2z_s(3) }{2+N^2}. \end{gather}

An axis constructed using this method is shown in figure 5, where the parameters were chosen as $N=2$, $R_c(3) = 5 {\times } 10^{-3}$, $z_s(2) = 3 {\times } 10^{-2}$ and $z_s(3) = 1.5 {\times } 10^{-3}$. We can see the rotation of the signed normal and binormal vectors has the expected behaviour, with one jump per field period.

Figure 5. Signed Frenet frame for a half-helicity two-field-period curve. Note the discontinuity that remains at the location where the magnetic field strength is taken to be maximal, where the normal vector points in the radial direction.

5. Constructing half-helicity near-axis stellarators

Now we can proceed to use this axis to construct a near-axis configuration, as described in § 2. A natural question to ask is whether the discontinuity in the signed Frenet frame at maxima of the magnetic field, required for half-helicity magnetic axes, has a detrimental effect on the continuity of the plasma boundary. We confirm that the formalism developed in Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022) is valid for fractional helicity and results in continuous plasma boundaries.

We know from the previous section that the signed Frenet frame, specifically the normal and binormal vectors are discontinuous at maxima of $B_{0}$ for half-helicity configurations. Given how the plasma boundary of a NAE equilibrium is constructed using (2.2),

(5.1)\begin{equation} \boldsymbol{x} \approx \boldsymbol{r}_{0} + \epsilon \left( X_{1}\boldsymbol{n}^{s} + Y_{1}\boldsymbol{b}^{s} \right), \end{equation}

it would be natural to assume the discontinuity of $\boldsymbol {n}^{s}$ and $\boldsymbol {b}^{s}$ translates into a discontinuity of the boundary shape, but this is only true if $X_1$ and $Y_1$ are continuous. Let us analyse the behaviour of this first-order correction to the plasma boundary:

(5.2)\begin{gather} X_{1} = \frac{d(\varphi)}{\kappa^{s}} \cos{[\theta - \alpha (\varphi) ]}, \end{gather}

(5.3)\begin{gather}Y_{1} = \frac{2\kappa^{s}}{B_{0}(\varphi)d(\varphi)} \left( \sin{[\theta - \alpha (\varphi)]} + \sigma(\varphi) \cos{[\theta - \alpha (\varphi) ]} \right). \end{gather}

The quantity $\bar {d} = d(\varphi )/\kappa ^{s}$ is, by definition, constructed to be symmetric and continuous at every point, in particular at $\varphi _{\mathrm {min}}$ and $\varphi _{\mathrm {max}}$ as discussed in § 2 and evident from (2.6) (Camacho Mata et al. Reference Camacho Mata, Plunk and Jorge2022). The magnetic field on axis $B_0(\varphi )$ is always continuous, also at points where the signed Frenet frame is discontinuous. Hence, any possible discontinuity in $X_1$ and $Y_1$ is going to depend on the behaviour of $\alpha (\varphi )$, given by (2.7), reproduced here:

(5.4)\begin{equation} \alpha(\varphi) = \iotaslash (\varphi-\varphi_{\mathrm{min}}^{i})+{\rm \pi} \left(2 m i + \frac{1}{2}\right) + {\rm \pi}\left( m - \iotaslash /N \right) \left(\frac{\varphi-\varphi_{\mathrm{min} }^{i}}{{\rm \pi}/N} \right)^{2k+1}. \end{equation}

Note the first term ensures that the magnetic perturbation is approximately in phase at two bounce points (omnigenity), while the final term is designed to break this criterion away from minima of the field, while ensuring the boundary retains the correct periodicity; see Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022) for more details. For half-helicity curves, $m=1/2$, and let us consider the first field period, which means $i=1$ and the minimum of $B$ is located at $\varphi _{\mathrm {min}}^{i}={\rm \pi} /N$. Now let us calculate the value that the $\alpha$-function takes at two consecutive maxima of $B$, i.e. the beginning and end of a field period, corresponding to $\varphi ^{1}_{\mathrm {max}}=0$ and $\varphi ^{2}_{\mathrm {max}}=2{\rm \pi} /N$:

(5.5)\begin{gather} \alpha(\varphi^{1}_{\mathrm{max}}) = \iotaslash (-{\rm \pi}/N)+ {\rm \pi}\left(1 + \frac{1}{2}\right) + {\rm \pi}\left( \frac{1}{2} - \iotaslash /N \right) \left(\frac{-{\rm \pi}/N}{{\rm \pi}/N} \right)^{2k+1}= {\rm \pi}, \end{gather}

(5.6)\begin{gather}\alpha(\varphi^{2}_{\mathrm{max}}) = \iotaslash (2{\rm \pi}/N-{\rm \pi}/N)+ {\rm \pi}\left(1 + \frac{1}{2}\right) + {\rm \pi}\left( \frac{1}{2} - \iotaslash /N \right) \left(\frac{{\rm \pi}/N}{{\rm \pi}/N} \right)^{2k+1}= 2{\rm \pi}. \end{gather}

We can see that for $m=1/2$ the function increases by ${\rm \pi}$ per field period, contrary to the $2{\rm \pi}$-periodicity necessary in the integer-helicity curves. Now let us see the impact this has on the functions $X_1$ and $Y_1$ describing the plasma boundary, (2.3) and (2.4):

(5.7)\begin{gather} X_{1}(\theta,\varphi=\varphi^{1}_{\mathrm{max}}) = \frac{d(\varphi)}{\kappa^{s}} \cos{[\theta - {\rm \pi}]} ={-} \frac{d(\varphi)}{\kappa^{s}} \cos{[\theta]}, \end{gather}

(5.8)\begin{gather}X_{1}(\theta,\varphi=\varphi^{2}_{\mathrm{max}}) = \frac{d(\varphi)}{\kappa^{s}} \cos{[\theta - 2{\rm \pi} ]} = \frac{d(\varphi)}{\kappa^{s}} \cos{[\theta]} ={-} X_{1}(\theta,\varphi=\varphi^{1}_{\mathrm{max}}),\end{gather}

(5.9)\begin{gather} Y_{1}(\theta,\varphi=\varphi^{1}_{\mathrm{max}}) = \frac{2\kappa^{s}}{B_{0}(\varphi)d(\varphi)} \left( \sin{[\theta - {\rm \pi}]} + \sigma(\varphi) \cos{[\theta - {\rm \pi}]} \right)\nonumber\\ ={-}\frac{2\kappa^{s}}{B_{0}(\varphi)d(\varphi)} \left( \sin{[\theta]} + \sigma(\varphi) \cos{[\theta]} \right),\end{gather}

(5.10)\begin{gather} Y_{1}(\theta,\varphi=\varphi^{2}_{\mathrm{max}}) = \frac{2\kappa^{s}}{B_{0}(\varphi)d(\varphi)} \left( \sin{[\theta - 2{\rm \pi}]} + \sigma(\varphi) \cos{[\theta - 2{\rm \pi}]} \right) \nonumber\\ ={-}Y_{1}(\theta,\varphi=\varphi^{1}_{\mathrm{max}}). \end{gather}

As shown above, the boundary components $X_1$ and $Y_1$ are discontinuous at maxima of $B_{0}$. This discontinuity is, necessarily, of the same type as that of the Frenet frame; continuity of the coordinate mapping (2.2) is thus achieved despite a discontinuous signed Frenet frame. Now we can use the formalism developed in Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022) to construct a magnetic configuration that is QI at first order, using a half-helicity magnetic axis.

5.1. A two-field-period half-helicity near-axis stellarator

Let us now proceed to show an example of a two-field-period configuration constructed around a half-helicity magnetic axis:

(5.11)\begin{gather} R = 1 - 0.074\cos{(2\phi)} - \tfrac{1}{17} \cos{(4\phi)} + 0.01 \cos{(6\phi)} , \end{gather}

(5.12)\begin{gather} z ={-}0.474 \sin{(2\phi)} + 0.06 \sin{(4\phi)} + 6\times10^{{-}3}\sin{(6\phi)}. \end{gather}

This axis was constructed to have $m=1/2$ as described in § 4, and the values of the near-axis parameters were chosen such that the corresponding VMEC equilibrium has small $\epsilon _{\mathrm {eff}}$ on its boundary, as discussed below. The shape of this magnetic axis, together with its signed Frenet frame, and curvature and torsion profiles are shown in figure 6. We can observe the expected curvature behaviour for half-helicity: a first-order zero at half-period and second-order zeros at the end/beginning of the period.

Figure 6. (a) Two-field-period half-helicity magnetic axis described by expressions (5.10) and (5.11). The normal and binormal vectors are shown. The discontinuity of the signed Frenet frame is evident at the beginning/end of the period. (b) Curvature (solid) and torsion (dashed) per-field-period profiles of the axis. Note the first-order zero of curvature at $\phi = {\rm \pi}/ N$ and second-order zeros at $\phi = 0, 2{\rm \pi} /N$, necessary to achieve values of helicity $m=1/2$.

The intensity of the magnetic field on this axis must also be specified and was chosen as

(5.13)\begin{equation} B_{0} (\varphi) = 1 + 0.165\cos{(2\varphi)}, \end{equation}

which, as required by the near-axis formalism, has extrema at points where the curvature of the magnetic axis is zero. Another parameter playing an important role in the quality of the solutions, specifically in their elongation and maximum values of $\epsilon _{\mathrm {eff}}$, is the function $d(\varphi )$. It was chosen to be proportional to the curvature of the axis as $d=d_{\kappa }\kappa ^s$, with

(5.14)\begin{equation} d_{\kappa}(\varphi) = \sqrt{0.4/B_{0}(\varphi)}. \end{equation}

This particular form of $d(\varphi )$, inversely proportional to $B_0$, is useful for controlling elongation, and here helps flatten the elongation profile, i.e. reduce its variation in $\varphi$, resulting in smaller maximum elongations than using a constant $\bar {d}$ choice.Footnote ² It is important to remember the elongation dependence on $\bar {d}$ is complicated and includes other quantities of the near-axis construction, making it difficult to predict which combination of parameters will result in a desired elongation behaviour. The freedom in $\bar {d}$, however, actually allows more control of the elongation profile, as compared with quasi-symmetric configurations where it must be constant; this is demonstrated in § 7.

The $\alpha$ function is defined as in (2.7), with $k=2$ chosen for this configuration. A near-axis solution is constructed using the aforementioned parameters at an aspect ratio $A=\sqrt {{\overline {B_{0}}}/{2}}({R_c(0)}/{\epsilon })=10$, where $\overline {B_0}$ is the average value of $B_0$. The boundary is converted to cylindrical coordinates using the method described in § 4.1 of Landreman et al. (Reference Landreman, Sengupta and Plunk2019), which results in elliptical cross-sections as shown in figure 7. This plasma boundary is then used to calculate a MHD equilibrium using the VMEC code. The intensity of the magnetic field on the last-closed flux surface is depicted in figure 8, which shows, even at this small aspect ratio, the behaviour expected for a QI configuration, with poloidally closing contours of $|B|$. These contours in Boozer coordinates are shown in figure 9. As expected, the contours look ‘more QI’ the closer we look to the axis, in the sense of having a greater degree of poloidal closure, as well as fewer island defects.

Figure 7. (a) Cross-sections of the two-field-period configuration shown in figure 8 at different toroidal values. (b) Elongation toroidal profile for the same configuration.

Figure 8. Magnetic field intensity in the plasma boundary for a two-field-period half-helicity configuration. (a) Side view and (b) top view.

Figure 9. Contours of magnetic field intensity for the configuration in § 5.1 at $s=0.1$ (a), $s=0.5$ (b) and $s=1$ (c). The poloidally closed contours of $B$, characteristic of QI configurations, degrade with distance from the axis.

The resulting configuration has maximum elongation under $5.5$ and a dependence on $\phi$ as shown in figure 7. Elongation here is taken at $\phi = \textrm {constant}$ in real space, and calculated as the ratio between the maximum and minimum distances from the axis on the resulting cross-section. One of the interesting properties of half-helicity configurations is the higher values of rotational transform they can achieve compared with zero-helicity cases due to the fact that the axis curves have higher values of integrated torsion. Here, a rotational transform spanning from $\iotaslash = 0.351$ to $0.367$ is found for the VMEC equilibria as seen in figure 10.

Figure 10. Rotational transform profile for the two-field-period half-helicity configuration.

The effective ripple, $\epsilon _{\mathrm {eff}}$, is a proxy often used to estimate neoclassical transport in the $1/\nu$ regime in stellarators. We have calculated it as described in Drevlak et al. (Reference Drevlak, Heyn, Kalyuzhnyj, Kasilov, Kernbichler, Monticello, Nemov, Nührenberg and Reiman2003) for 16 radial points. A maximum value of $\epsilon _{\mathrm {eff}}=1\,\%$ is found, as shown in figure 18. The values of the effective ripple are similar to those of the standard configuration of W7-X, with the important distinction that this configuration consists entirely of elliptical cross-sections and no boundary optimization was required. Half-helicity configurations seem to provide the high rotational transform benefits from $m=1$ configurations while preserving the low values of neoclassical transport found in $m=0$ equilibria. In the next section we show that we can take advantage of this property to construct $N>3$ equilibria with low effective ripple.

6. The half-helicity space

Now that we have adjusted the near-axis method to deal with the possibility of half-helicity axes, we proceed to identify in which region of the space of solutions we can find such configurations. Let us analyse the space of five-field-period configurations using the same procedure as that in § 3.

Once again, we hold all parameters fixed, apart from those of the magnetic axis, and we parametrize the axis in cylindrical coordinates as

(6.1)\begin{gather} R(\phi) = 1 - \frac{1}{101} \cos (5\phi), \end{gather}

(6.2)\begin{gather} z(\phi) = z_s(1)\sin(5\phi) + z_s(2)\sin(10\phi). \end{gather}

The aspect ratio is chosen as $A=20$, the function $\alpha (\varphi )$ is described by (2.7) with $k=2$, the magnetic field on axis is $B_0 = 1+0.15\cos (\varphi )$ and $d = 0.73 \kappa ^s$. Each colour-coded square in figure 11 is an integer-helicity axis configuration, with values of $z_s(1)$ ranging from 0 to 0.18 and $z_s(2)$ from $-$0.425 to 0.425, constructed using the near-axis method described in § 2.

Figure 11. Space of five-field-period near-axis configurations. Squares correspond to integer-helicity values ($m=0,1$) and circles to half-helicity ($m=1/2$). Colours indicate the maximum effective ripple. An example of each type of helicity is shown on the right. Note that the half-helicity configurations attain values of ${\epsilon _{\mathrm {eff}}}{}$ similar to those of integer configurations but have higher values of rotational transform.

We also include half-helicity axes for which $R_c(3)$ and $z_s(3)$ in expression (4.8) are zero. Thirty-seven values of $z_s(1)$ are considered, ranging from 0 to 0.18. In order to fulfil the half-helicity constraints on the axis, $z_s(2)$ is calculated following expression (4.8):

(6.3)\begin{equation} z_s(2) = \frac{2+N^2}{4+8N^2} z_1. \end{equation}

The configurations for these axes are constructed as in the previous section: the only difference from the integer-helicity cases is the need to set $m=1/2$ when calculating the function $\alpha (\varphi )$. All other initial parameters are as described earlier in this section. Given how $z_s(1)$ depends on the chosen values for $z_s(2)$ when restricting ourselves to two parameter curves, the half-helicity space has one dimension fewer than the equivalent integer-helicity space. This is always the case, independently of the number of Fourier coefficients included in the axis representation.

In figure 11, the circles correspond to the half-helicity configurations and are distributed on the line where the conditions for second-order zeros of curvature at the maxima of $B_0$ are fulfilled. The line in the upper quadrant is for those axes in which $\kappa$ is zero at second order at the minima of $B_0$. These are not consistent with near-axis QI equilibria; hence the lack of configurations in said region.

A configuration from the zero-helicity section is shown above a half-helicity configuration in figure 11. Both are the configurations with lowest ${\epsilon _{\mathrm {eff}}}{}$ in the boundary, in each space. We can see that both have similar boundary effective ripple values (1.99 % for $m=0$ and 2 % for $m=1/2$) but the rotational transform is substantially larger for the half-helicity case ($\iotaslash _{m=0}=0.23$ and $\iotaslash _{m=1/2}=0.85$). This behaviour is consistently observed between all configurations in this space. The half-helicity space, apart from its lower dimensionality, also has higher values of rotational transform, closer to the per-field-period values encountered in traditional optimized stellarators. Higher values of rotational transform are desirable as they are related to improved confinement (Ascasíbar et al. Reference Ascasíbar, López-Bruna, Castejon, Vargas, Tribaldos, Maassberg, Beidler, Brakel, Dinklage and Geiger2008) and given the dependence of ideal MHD $\beta$-limit on $\iotaslash$ (Miyamoto Reference Miyamoto2005; Loizu et al. Reference Loizu, Hudson, Nührenberg, Geiger and Helander2017).

7. Five-field-period configuration

Up until now it has proven difficult to find configurations with integer helicity, high numbers of field periods and reasonably low aspect ratio, which are accurate to first order, using VMEC as described previously. This is true even when attempting to optimize within the space of near-axis solutions. We now show that this is not the case when searching in the lower-dimensional space of half-helicity magnetic axes in the NAE.

A magnetic axis shape is chosen by varying the Fourier coefficients $R_c(3), z_s(2)$ and $z_s(3)$ of expression (4.7) and imposing the half-helicity conditions in (4.8). An axis is selected such that low values of effective ripple are achieved in the boundary. The one used for the construction of the following example is described by

(7.1)\begin{gather} R = 1 + \tfrac{0.2260}{26}\cos{(5\phi)} - \tfrac{1}{101} \cos{(10\phi)} - 1\times10^{{-}3} \cos{(15\phi)} , \end{gather}

(7.2)\begin{gather} z ={-}\tfrac{3.1986}{27} \sin{(5\phi)} +7\times10^{{-}3} \sin{(10\phi)} + 2.6\times10^{{-}3}\sin{(15\phi)}. \end{gather}

The curvature and torsion per field period are shown in figure 12. For this high number of field periods, the magnetic field on axis and the function $d(\varphi )$ need to be carefully chosen to reduce $\epsilon _{\mathrm {eff}}$ to acceptable levels. These quantities are parametrized as

(7.3)\begin{gather} B_{0} (\varphi) = 1 + B0_1 \cos{(5\varphi)} + B0_2 \cos{(10\varphi)}, \end{gather}

(7.4)\begin{gather} d(\varphi) = \sqrt{d_{\kappa}/B_{0}(\varphi))}\kappa^s + d_{\kappa c} \kappa^s \cos{(5\phi)} + d_{\kappa s} \kappa \sin{(5\phi)}. \end{gather}

Figure 12. Curvature (solid) and torsion (dashed) profiles of the five-field-period magnetic axis described by (7.1) and (7.2).

Adding more parameters to these expressions allows for a finer profile control. These parameters are varied for a set axis shape until a desired configuration is found. For this example $d_{\kappa }=0.28$, $d_{\kappa c}=-0.065$, $d_{\kappa s}=0.04$, $B0_1 = 0.12$ and $B0_2 = -0.002$. The resulting shapes of $B_0(\varphi )$ and $\bar {d}(\varphi )$ are shown in figure 13. Here we can observe a non-analytic behaviour present at $\varphi ={\rm \pi} /N$. This arises due to the sine term in expression (7.4) and can be set to zero but including it resulted in lower values of $\epsilon _{\mathrm {eff}}$, and, most importantly, it does not result in any sharp behaviour on the plasma boundary. As in the previous section, the function $\alpha$ is defined through (2.7) with $k=2$. Following the NAE formalism, we construct a solution at a set distance from the axis, in this case $r=1/12$ corresponding to an aspect ratio $A=12$. However, using the method of § 4.1 of Landreman et al. (Reference Landreman, Sengupta and Plunk2019) was not enough to obtain equilibria with good confinement; instead the method of § 4.2 of the same publication was used. This method, albeit more complicated, results in boundary shapes identical at first order to the boundary in the $X$–$Y$ space, with the caveat that cross-sections in the $R$–$z$ plane are no longer elliptical, as can be seen in figure 14. Having a boundary defined by elliptical cross-sections requires the use of less poloidal modes when using equilibrium codes like VMEC and might be useful for plasma boundary optimization, where initial points with low modes are preferred, although not necessary (Landreman Reference Landreman2022). 

Figure 13. (a) Profiles of $\bar {e}$ for the two five-field-period configurations presented. The dotted line correspondsto the configuration with parameter $d(\varphi )$ optimized for low elongation. (b) Magnetic field on axis $B_0(\varphi )$ used for the construction of these configurations.

Figure 14. (a) Cross-sections of the five-field-period configuration shown in figure 15. (b) Elongation profiles for the two five-field-period configurations shown in § 7. The dotted line corresponds to the configuration with parameter $d(\varphi )$ optimized for low elongation.

Figure 15. Intensity of the magnetic field on the boundary of the five-field-period configuration. (a) Side view and (b) top view.

The intensity of the magnetic field on the boundary, as obtained from VMEC, is shown in figure 15. Despite the contours of $|B|$ looking far from the expected behaviour for a QI field, this configuration has good neoclassical confinement properties, as is evident from its effective ripple profile in figure 18. With ${\epsilon _{\mathrm {eff}}}{}$ $\approx 1.3\,\%$ it is comparable to that of W7-X and QIPC without the need of performing a boundary optimization. As is typical for near-axis configurations, the contours of $|B|$ resemble more straight poloidally closed contours the closer we look to the axis, as seen in figure 16.

Figure 16. Contours of magnetic field intensity for the configuration in § 7 at $s=0.1$ (a), $s=0.5$ (b) and $s=1$ (c). The poloidally closed contours of $B$, characteristic of QI configurations, degrade with distance from the axis.

The elongation of the plasma boundary cross-sections in the cylindrical angle $\phi$ is shown in figure 14, where a maximum value of $e=5.05$ is found. We can observe a clear difference in the elongation profile of this configuration and that of the two-field-period example shown in figure 7. The rotational transform radial profile for this configuration is shown in figure 17, where is clear the per-field-period rotational transform value achieved is close to those in existing optimized stellarators.

Figure 17. Rotational transform profile for the five-field-period configuration of § 7.

7.1. Optimizing $d$ for low elongation

The choice of the function $d(\varphi )$ plays an important role in the resulting elongation of the plasma boundary and needs to be chosen carefully, as most values will result in impractically large elongations. This function can be tuned to obtain a desired elongation profile, although the dependence on this and other quantities is not straightforward, as is clear from equation (5.1) of Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022):Footnote ³

(7.5)\begin{equation} e(\varphi) = \frac{B_{0}\bar{d}^{2}}{4} + \frac{1}{B_{0}\bar{d}^2}(1+\sigma^2) + \sqrt{\left( \frac{B_{0}\bar{d}^{2}}{4}\right)^2 + \frac{(\sigma^2-1)}{2} + \left(\frac{1}{B_{0}\bar{d}^2}(1+\sigma^2)\right)^{2}}. \end{equation}

For this reason, in the following configuration, $d(\varphi )$ is chosen through an optimization procedure to minimize elongation. All other near-axis parameters are kept the same as in the previous example. The function $d$ is parametrized as

(7.6)\begin{equation} d(\varphi) = \sqrt{d_{\kappa}/B_{0}(\varphi))}\kappa^s + \sum_{i=1}^4 d_{\kappa c}(i)\kappa^s \cos{(5i\phi)}, \end{equation}

and the coefficients $d_{\kappa c}(i)$ are optimized using an off-the-shelf optimization methodFootnote ⁴ to reduce the elongation variation, i.e. the difference between the maximum and minimum elongation, both outputs of the NAE construction. For this configuration the values resulting in the elongation shown with a dotted line in figure 14 are

(7.7)\begin{equation} d_{\kappa} = 0.4, \qquad d_{\kappa c} = [0.0591,-0.0282,-0.0157,-0.0038], \end{equation}

corresponding to $\bar {e}(\varphi )$ shown in figure 13. We can see that the variation in the elongation profile is greatly reduced compared with the previous configuration, which is also evident by looking at the toroidal cross-sections of the plasma boundary shown in figure 19, and comparing them with the previous configuration as done in figure 20. But even using four Fourier modes in the representation of $d(\varphi )$ is not enough to obtain a constant elongation profile, which is a particularly attractive case, allowed by the NAE QI theory, due to its especially smooth boundary behaviour. This fact, together with the complicated dependence of elongation, and the parameters required for the equilibrium construction (see (7.6)), as well as the non-analyticity of $\bar {d}(\varphi )$ necessary to achieve good confinement in the example in § 7 are clear indications of the unsuitability of $d(\varphi )$ as a natural input parameter for the NAE. A possible solution to this problem is specifying the elongation function as an input parameter instead. This will be described in a future publication.

Figure 18. Effective helical ripple ${\epsilon _{\mathrm {eff}}}{}$ for the three configurations shown in this work (solid lines) and two QI optimized configurations: W7-X (dotted line) and QIPC (dashed line).

Figure 19. (a) Cross-sections of the configuration shown in figure 21. (b) Rotational transform profile for the same configuration.

Figure 20. Cross-section comparisons between the two five-field-period configurations shown in § 7 at different toroidal locations. Dotted lines correspond to the configuration optimized for constant elongation.

The rotational transform profile is shown in figure 19. It crosses the value $\iotaslash =1$, which might result in the presence of magnetic islands degrading the nested flux surfaces necessary for confinement. No effort has been made in this work to avoid such regions and further optimization of this configuration with specific targets to limit $\iotaslash$ can be applied if necessary.

The intensity of the magnetic field on the boundary ($A=12$) for the configuration constructed in this section is shown in figure 21. We can see that the contours at maximum field do not close as expected but nonetheless $\epsilon _{\mathrm {eff}}$ remains under $1\,\%$ for half-flux and is comparable to that of QIPC in the remaining region as shown in figure 18. As expected, a reduction in elongation results in worse QI quality as discussed in Camacho Mata et al. (Reference Camacho Mata, Plunk and Jorge2022).

Figure 21. Magnetic field intensity on the boundary for the five-field-period low-elongation configuration described in § 7.1. (a) Side view and (b) top view.