3.1. Linear eigenmodes in low-$\hat {s}$ configurations

Typically, as $\beta$ increases, the normalized ITG growth rates $\gamma$ in units of $c_{s} / a$, where $c_{s}$ is the sound speed, steadily decrease (Pu & Migliuolo Reference Pu and Migliuolo1985; Dong, Guzdar & Lee Reference Dong, Guzdar and Lee1987; Kim, Horton & Dong Reference Kim, Horton and Dong1993). This is known as finite-$\beta$ ITG stabilization, or linear electromagnetic stabilization. The mechanism by which this stabilization occurs is the coupling between the ITG mode and the shear Alfvén wave. This efficiently transfers energy out of the ITG mode, reducing its growth. Alterations of the ion Landau resonance due to finite-equilibrium-$\beta$ $\boldsymbol {\nabla } B$ modifications from the bending of perturbed field lines can also play a role in the stabilization (Jarmén, Anderson & Malinov Reference Jarmén, Anderson and Malinov2015). This reduction in growth rate will tend to reduce the nonlinear heat flux, but additional physics can also have an impact on the dynamics.

Another mechanism that affects ITG transport levels as $\beta$ increases is nonlinear electromagnetic stabilization (Pueschel et al. Reference Pueschel, Kammerer and Jenko2008; Pueschel & Jenko Reference Pueschel and Jenko2010; Whelan, Pueschel & Terry Reference Whelan, Pueschel and Terry2018; Whelan et al. Reference Whelan, Pueschel, Terry, Citrin, McKinney, Guttenfelder and Doerk2019). As $\beta$ increases, one observes a reduction in nonlinear heat flux that is greater than the reduction in linear growth rates or quasilinear fluxes. As has been shown in Whelan et al. (Reference Whelan, Pueschel and Terry2018, Reference Whelan, Pueschel, Terry, Citrin, McKinney, Guttenfelder and Doerk2019), the difference between the (quasi-)linear and nonlinear reduction stems from changes in the efficiency of zonal-flow-mediated energy transfer to stable modes. This effect has been shown to be active in a range of configurations, both of tokamak and of stellarator type. It has also been demonstrated that the inclusion of a self-consistent Shafranov shift can eliminate nonlinear electromagnetic stabilization for different, specific parameter regimes (Ishizawa et al. Reference Ishizawa, Urano, Nakamura, Maeyama and Watanabe2019).

The physical mechanism underlying KBM destabilization is the following. As with most microinstabilities, KBMs are destabilized when the driving force associated with the pressure gradient is sufficiently strong in the bad curvature region to overcome the stabilizing force from the magnetic field (Cheng Reference Cheng1982). More specific to KBMs, the modifications to magnetic drifts by $\beta$-induced magnetic fluctuations result in destabilizing effects when a kinetic treatment is applied to the governing equations. This is due to a resonance between the KBM and thermal ions, giving the KBM an additional free energy source.

For the cases studied here, the dominant eigenmode for $\beta < \beta _{\rm crit}^{\rm KBM}$ is ITG, whose growth rate decreases with $\beta$. For $\beta \geq \beta _{\rm crit}^{\rm KBM}$, KBM has a larger growth rate than the ITG mode with the growth rate increasing rapidly with $\beta$. This critical value $\beta _{\rm crit}^{\rm KBM}$ is generally of the order of 1 % in standard tokamaks and is comparable to the ideal MHD ballooning $\beta$ limit for physically relevant gradient values (see Pueschel et al. Reference Pueschel, Kammerer and Jenko2008; Pueschel & Jenko Reference Pueschel and Jenko2010; Ishizawa et al. Reference Ishizawa, Maeyama, Watanabe, Sugama and Nakajima2013). However, $\beta _{\rm crit}^{\rm KBM}$ is not always close to $\beta _{\rm crit}^{\rm MHD}$ and is, in general, a complicated function of the magnetic geometry.

Before presenting analyses of nonlinear simulations, the behaviour and scaling of linear instability will be elucidated, for both dominant and subdominant eigenmodes. Typically, ITG and KBM growth rates decrease and increase, respectively, as $\beta$ increases. Figure 1 highlights this behaviour in HSX at normalized binormal wavenumber $k_y \rho _{s} = 0.6$. Henceforth, normalized wavenumbers are denoted $k_y \rho _{s} \rightarrow k_y$. The normalized gradients, $\beta$ and the temperature ratio used in these calculations and throughout this work are $a / L_{T{i}} = 3$ and $a / L_{T{e}} = a / L_n = 1$, $\beta = 0.48\,\%$ and $T_{i} / T_{e} = 1$, respectively, unless otherwise stated. Note the factor-of-three reduction in the ITG growth rate as $\beta$ increases before the KBM becomes dominant. This significant reduction in ITG growth rate is not observed for all values of $k_y$, as evidenced by figure 2, where at $k_y = 0.1$ there is no significant reduction in $\gamma$ before KBMs become dominant for either the $N_{\rm pol} = 1$ or the $N_{\rm pol} = 4$ case. This is partially due to the fact that, in the $k_y = 0.1$ case, the KBM becomes dominant at a small critical normalized plasma pressure $\beta _{\rm crit}^{\rm KBM} \approx 0.18\,\%$ versus $\beta _{\rm crit}^{\rm KBM} \approx 2.2\,\%$ for the $k_y = 0.6$ case. Also note that $N_{\rm pol} = 4$ is required to achieve convergence, as elongated eigenmodes can artificially self-reinforce via the parallel boundary condition.

Figure 1. The growth rate and real frequency for the eigenmode with the largest growth rate at $k_y = 0.6$ over a range of $\beta$ values for HSX with $a / L_{T{i}} = 3$ and $a / L_{T{e}} = a / L_n = T_{i} / T_{e} = 1$. Note the stabilization of ITG growth rates (shown as triangles) as $\beta$ increases and approaches $\beta _{\rm crit}^{\rm KBM} \approx 2.2\,\%$, where there is a clear discontinuity in real frequency, which highlights the change in dominant mode branch from ITG to KBM (shown as diamonds).

Figure 2. Growth rate and real frequency for $k_y = 0.1$ of the HSX configuration with $a / L_{T{i}} = 3$ and $a / L_{T{e}} = a / L_n = T_{i} / T_{e} = 1$ for both $N_{\rm pol} = 1$ (blue diamonds) and $N_{\rm pol} = 4$ (black triangles). Note the lack of ITG stabilization relative to the $k_y = 0.6$ case and the significant difference between the two cases, an indication that $N_{\rm pol} = 1$ is insufficient for convergence. Also, $\beta _{\rm crit}^{\rm KBM}$ is much lower for this case.

The scaling of $\beta _{\rm crit}^{\rm KBM}$ with $k_y$ is presented in figure 3 for both HSX (blue squares) and NCSX (black triangles), a quasi-axisymmetric stellarator configuration optimized for operation at $\beta \approx 4.2\,\%$ (Neilson, Zarnstorff & Lyon Reference Neilson, Zarnstorff and Lyon2002) with $\hat {s} \approx -0.5$ at $s_0 = 0.5$. For HSX, $\beta _{\rm crit}^{\rm KBM}$ is small relative to the ideal MHD ballooning limit $\beta _{\rm crit}^{\rm MHD}$, and considerably smaller than the high-$\hat {s}$ NCSX analogue where $\beta _{\rm crit}^{\rm KBM}$ is comparable to $\beta _{\rm crit}^{\rm MHD}$. The normalized gradients and temperature ratio used in the NCSX calculations are the same as were used for HSX. The ideal ballooning limit $\beta _{\rm crit}^{\rm MHD}$ is also shown as a horizontal dashed line for each configuration. The NCSX curve smoothly approaches the MHD limit as $k_y$ decreases until $k_y < 0.1$ where there is a slight uptick in $\beta _{\rm crit}^{\rm KBM}$. The HSX curve is well below the MHD limit in the transport-relevant $k_y < 0.25$ range, and dips to values $\approx 0.1\,\%$ that are an order of magnitude smaller than $\beta _{\rm crit}^{\rm MHD}$. This suggests that conditions in HSX are more conducive to KBM excitation due to ion kinetic physics (see § 4) than they are in NCSX. As will be shown, this is associated with the low average magnetic shear in HSX.

Figure 3. The $\beta _{\rm crit}^{\rm KBM}$ spectrum of NCSX (black triangles) and HSX (blue squares) with $a / L_{T{i}} = 3$ and $a / L_{T{e}} = a / L_n = T_{i} / T_{e} = 1$. Horizontal dashed lines correspond to $\beta _{\rm crit}^{\rm MHD}$. Note the low $\beta _{\rm crit}^{\rm KBM}$ relative to $\beta _{\rm crit}^{\rm MHD}$ for HSX compared to NCSX.

Figure 4 shows the dominant growth rates $\gamma$ and real frequencies $\omega$ as functions of $k_y$ at a constant $\beta = 0.48\,\%$. It is worth noting that KBM is dominant over a certain range of $k_y$ values, namely $k_y \in [0.1,0.2]$, while ITG is dominant over the remainder of the wavelength range. This means that nonlinear calculations of HSX at this value of $\beta$ will potentially exhibit concurrent ITG and KBM drive characteristics, a phenomenon that is also observed in nonlinear simulations of tokamaks if $\beta \geq \beta _{\rm crit}^{\rm KBM}$ (Pueschel et al. Reference Pueschel, Kammerer and Jenko2008; Pueschel & Jenko Reference Pueschel and Jenko2010).

Figure 4. Growth rate and real frequency spectra for the HSX configuration (both with an equilibrium $\beta = 0$ (black) and equilibrium $\beta = 0.48\,\%$ (blue)) with $a / L_{T{i}} = 3$, $\beta = 0.48\,\%$ and $a / L_{T{e}} = a / L_n = T_{i} / T_{e} = 1$. The KBMs (ITG modes) are denoted by diamonds (triangles). Note the transition from ITG to KBM (at roughly $k_y \approx 0.1$) and back to ITG (at roughly $k_y \approx 0.2$), highlighted by the discontinuity in real frequency. Also note the similarity in the two spectra when using a finite-$\beta$ magnetic equilibrium versus using the vacuum case.

Lastly, the subdominant mode spectrum of HSX at $k_y = 0.2$ with $\beta = 0.5\,\%$, shown in figure 5, indicates that there is not just a single KBM that is destabilized as $\beta$ increases, but rather there are two families of KBMs, one of which is a set of modes that are centred at the outboard midplane ($\theta _p=0$) and another that consists of pairs of sibling modes that peak away from the outboard midplane ($\theta _p \neq 0$). It should be noted that some of the centred KBMs exhibit tearing parity (blue and green curves in figure 5b), i.e. the $\varPhi$ eigenfunctions are odd in ballooning angle. The low background magnetic shear in the present scenario enables destabilization of these modes, termed tearing-parity KBMs, which had previously been conjectured not to exist based on a study involving high-shear equilibria (Pueschel et al. Reference Pueschel, Hatch, Ernst, Guttenfelder, Terry, Citrin and Connor2019). It is worth noting that it is unlikely that the tearing-parity KBMs discussed here are micro tearing modes because of the unique ITG dependence of the modes, the fact that the modes propagate in the ion diamagnetic direction and the relatively small electron thermal transport nonlinearly, three qualities that are not characteristic of micro tearing modes. For each mode that peaks away from the outboard midplane, a sibling mode exists with the same growth rate and frequency with an eigenmode structure that is near-perfectly mirrored with respect to $\theta _p = 0$. The physics implication of this is that HSX is qualitatively different from typical high-$\hat {s}$ tokamak cases, where only a single KBM is destabilized. Eigenvalue calculations with different gradients, $a / L_{T{i}} = 4$ for example (with $a / L_n$ and $a / L_{T{e}}$ kept fixed), also exhibit the two branches of KBMs. However, the exact location of minimal $\beta _{\rm crit}^{\rm KBM}$ shifts in $k_y$ space relative to the $a / L_{T{i}} = 3$ case. The region of $k_y$ space in which the outboard-midplane-peaked KBMs are dominant shifts to smaller $k_y$ as $\omega _{T{i}}$ increases. It should also be noted that the two branches of KBMs are robustly present when ensuring numerical convergence. This also complicates analysis of the nonlinear energy dynamics, as there may be a number of unstable KBMs involved in nonlinear energy transfer.

Figure 5. The subdominant spectrum (a) for HSX with $k_y = 0.2$ and $\beta = 0.5\,\%$ consisting of KBMs (diamonds) and ITG modes (triangles), and the associated electrostatic potential $\varPhi$ eigenmode structures (b) for various KBMs. Note the two distinct families of KBMs: one with $\varPhi$ symmetric about $\theta _p = 0$ and one with off-centre peaking. For each mode of the latter family, there is a sibling mode that is mirrored across the $\theta _p = 0$ axis with identical $\gamma$ and $\omega$; hence they are indistinguishable in (a). Also note that some of the centred KBMs (blue, green) have tearing parity, i.e. are odd functions of ballooning angle.

3.2. Nonlinear characteristics of the ITG–KBM system

Previous studies of KBM instability and turbulence have found critical $\beta$ values of the order of one or a few per cent for experimentally relevant gradient values and magnetic geometries (Pueschel et al. Reference Pueschel, Kammerer and Jenko2008; Pueschel & Jenko Reference Pueschel and Jenko2010; Ishizawa et al. Reference Ishizawa, Maeyama, Watanabe, Sugama and Nakajima2013). These analyses also show that for $\beta < \beta _{\rm crit}^{\rm KBM}$, the resulting turbulence is solely ITG or trapped-electron mode and for $\beta > \beta _{\rm crit}^{\rm KBM}$, signatures of multiple modes are present and observable in the simulations concurrently. Contrary to these cases, when $\beta > \beta _{\rm crit}^{\rm KBM}$, nonlinear simulations of the present HSX case do not achieve saturation when the minimum binormal wavenumber $k_y^{\rm min}$ of the system is unstable to KBMs. In such simulations, streamer modes span the length of the periodic radial domain, regardless of how large the radial box is, and self-reinforce via the radial periodic boundary condition of the flux tube. However, if one chooses $k_y^{\rm min}$ so that it is stable to KBMs, saturation can be achieved in electromagnetic flux-tube simulations, where a significant reduction in transport is observed relative to the near-electrostatic ($\beta = 0.05\,\%$) case. The time traces of heat flux (both ion electrostatic and electron electromagnetic) and particle flux for $k_y^{\rm min} = 0.025$ are shown in figure 6. Despite the KBM being an electromagnetic mode, below the threshold it does not drive significant electromagnetic flutter heat transport given the choice of small electron temperature gradient $a / L_{T{e}} = 1$. While the condition that the smallest finite $k_y$ must be KBM-stable is empirical, one possible underlying cause is the ability of the system to terminate an inverse cascade by means of ITG (stable) eigenmodes. The precise nature of the phenomenon is left for future study.

Figure 6. Heat (both ion electrostatic (black) and electron electromagnetic (blue)) and particle flux (red) time traces for HSX with $k_y^{\rm min} = 0.025$ and $\beta = 0.48\,\%$. There is a significant (factor of $\approx 5$) reduction in $Q_{i}^{\rm es}$ relative to the near-electrostatic ($\beta \approx 0.05\,\%$) case (McKinney et al. Reference McKinney, Pueschel, Faber, Hegna, Talmadge, Anderson, Mynick and Xanthopoulos2019). There is negligible, slightly negative magnetic flutter transport, a result of the low electron temperature gradient $a / L_{T{e}} = 1$.

The heat flux spectrum given in figure 7 shows that the ion electrostatic heat flux $Q_{i}^{\rm es}(k_y)$ peaks in the $k_y$ range that is dominated by KBMs, consistent with the fact that $\beta$ is roughly three times as large as $\beta _{\rm crit}^{\rm KBM}$ for $k_y = 0.1$. Figure 7 also shows that KBM transport is the dominant physics mechanism nonlinearly as long as KBMs are destabilized linearly, even if KBM growth rates are subdominant, where linear destabilization is shown in figure 4. As discussed above, $k_y^{\rm min}$ is sufficiently small that no KBM instability occurs at that mode; this is consistent with findings in high-$\hat {s}$ scenarios where nonlinear saturation is possible for $\beta$ values up to $\beta _{\rm crit}^{\rm KBM}$ as $k_y \rightarrow 0$ (Pueschel & Jenko Reference Pueschel and Jenko2010). The nonlinear electrostatic potential $\varPhi$ spectrum is depicted in figure 8, showing that the simulation is not dominated by a zonal flow, as the $k_y=0.1$ non-zonal component of $\varPhi$ is a factor of nearly two larger than the zonal component and the sum of the non-zonal contributions is significantly larger than the zonal contribution. This is an indication that zonal flows may play at most a minor role in saturation. This view is refined, however, based on nonlinear energy transfer analysis.

Figure 7. Heat flux spectra for HSX with $k_y^{\rm min} = 0.025$ and $\beta = 0.48\,\%$ (black) associated with the time trace shown in figure 6 and with $k_y^{\rm min} = 0.025$ and $\beta = 0.24\,\%$ (blue) for comparison. Note that the peaks of the spectra are almost fully contained in the binormal wavelength range over which KBM is destabilized, which is denoted by the shaded region. This shows that KBM transport dominates nonlinearly when KBMs are linearly destabilized, even if KBM growth rates are subdominant to ITG growth rates.

Figure 8. The nonlinear $\varPhi$ spectrum for HSX with $k_y^{\rm min} = 0.025$ and $\beta = 0.48\,\%$. The spectrum is truncated at $k_y \rho _{s} = 0.5$ since $|\varPhi |^2$ amplitudes are negligible above this threshold. Note that the zonal component is much weaker than the integrated non-zonal amplitudes.

Nonlinear frequencies are shown in figure 9. Interestingly, there is no discontinuity in the nonlinear frequency spectrum of the HSX configuration, a feature that was present linearly at the ITG–KBM transition points in figure 4. Both dominant and subdominant linear KBM frequencies are overlaid in figure 9 to facilitate comparison with the nonlinear data. Clearly, throughout the range where linear data are shown, nonlinear frequency signatures match the values associated with the dominant outboard-centred KBMs. An exact frequency match is not expected to occur, as turbulence may result in a $k_y$-dependent nonlinear frequency shift and broadening. This constitutes evidence that KBMs do indeed play an important role in the nonlinear turbulent state.

Figure 9. The nonlinear frequency spectrum for HSX with $k_y^{\rm min} = 0.025$ and $\beta = 0.48\,\%$. The colour scale has arbitrary units and is linear and normalized at each $k_y$ separately. Linear frequencies of the most unstable (gold crosses), maximal-frequency (green triangles), minimal-frequency (red squares) and most unstable outboard-centred (black–white dashed diamonds) KBMs are also included for comparison. Note the agreement between the dominant linear outboard-centred KBM frequencies and the nonlinear signal.

Analysis of the nonlinear energy transfer to a given ($k_x$,$k_y$) point in Fourier space due to interaction with ($k^\prime _x$,$k^\prime _y$) and ($k_x-k^\prime _x,k_y-k^\prime _y$) also suggests that modes in the KBM-dominated $k_y$ range play an important role in the energy transfer dynamics of the turbulence and therefore in the dynamics that leads to saturation. Figure 10 shows nonlinear energy transfer at a given ($k_x=0.06$, $k_y=0.1$), denoted by the tip of the black arrow, to and from various ($k^\prime _x$, $k^\prime _y$). Regions of blue correspond to locations which cause energy input into ($k_x=0.06$, $k_y=0.1$) while red regions correspond to locations which draw energy from ($k_x=0.06$, $k_y=0.1$). The dominant method by which ($k_x=0.06$, $k_y=0.1$) receives energy is zonal energy transfer, indicated by the blue clouds at both $k^\prime _y=0$ and $k^\prime _y = 0.1$. There is also significant non-zonal energy transfer via ($k^\prime _x=0.12$, $k^\prime _y=-0.025$). It is important to note both that the largest energy sinks are due to non-zonal transfer and that the largest energy sinks are even larger than the largest energy inputs from zonal transfer. Lastly, the zonal energy transfer can occur due to both the zonal flow and the zonal field, i.e. the zonal component of the fluctuating magnetic field.

Figure 10. Nonlinear energy transfer functions indicate locations which give (blue) and receive (red) energy to and from ($k_x = 0.06$, $k_y = 0.1$), denoted by the tip of the black arrow. There is significant zonal transfer from ($k_x = 0.06$, $k_y = 0.1$) to the blue clouds near $k^\prime _y=0$ and $k^\prime _y = 0.1$. Significant non-zonal energy transfer is also observed at $k^\prime _y = -0.025$ and $k^\prime _y = 0.125$.

This result appears to contradict the earlier finding in figure 8 that the zonal-flow amplitudes are low. However, when coupling is sufficiently resonant, zonal flows can very efficiently mediate energy transfer even if they are only excited to low amplitudes.

Further analysis of nonlinear energy transfer yields insight into the wavelength ranges which dominate energy transfer between the turbulence and the zonal ($k_y=0$) modes. Figure 11 shows normalized nonlinear energy transfer corresponding to $k_y=0$ modes summed over both $k_x$ and $k^\prime _x$ versus $k^\prime _y$. The peak of the data occurs at $k^\prime _y = 0.125$, well within the KBM-dominated $k_y$ range. A significant portion of the energy transfer to and from zonal modes is facilitated by modes in the KBM-dominated $k_y$ range, constituting additional evidence that KBMs are important nonlinearly. One can also construct a quantity to gauge the relative importance of a given $k_y$ in the overall nonlinear energy transfer dynamics, as highlighted by figure 12, a plot of averaged nonlinear energy transfer functions versus $k_y$. This quantity is the sum (over $k_x \neq 0$) of the root mean squares (of each $k^\prime _x$–$k^\prime _y$ plane) of the nonlinear energy transfer functions corresponding to a given ($k_x$,$k_y$). Both the zonal modes and KBMs contribute to the overall nonlinear energy transfer dynamics, with KBMs contributing roughly twice as much if one integrates over the entire KBM-dominated $k_y$ range.

Figure 11. The sum (over $k_x \neq 0$) of the root mean squares (over coupled $k_x^\prime$) of time-averaged nonlinear energy transfer functions for zonal ($k_y = 0$) modes as a function of $k^\prime _y$. The data are normalized to the value of the point at $k^\prime _y = 0.125$. Note that the peak is in the KBM-dominated $k^\prime _y$ range, evidence that KBMs play an important role in zonal dynamics.

Figure 12. The sum (over $k_x \neq 0$) of the root mean squares (over coupled $k_x^\prime$, $k_y^\prime$) of time-averaged nonlinear energy transfer functions, which are normalized by the value of the same quantity at $k_y = 0.1$.

Lastly, figure 13 shows both the electrostatic ion heat and particle fluxes as $\beta$ increases. Simulations with $k_y^{\rm min} = 0.025$ for $\beta = 0.75\,\%$ and $1\,\%$ grow without bounds and therefore do not have associated data points. This is consistent with the requirement that $k_y^{\rm min}$ be stable to KBMs for saturation to occur, as $\beta _{\rm crit}^{\rm KBM} \approx 0.6\,\%$ for $k_y = 0.025$, as shown in figure 3. Note the decrease in heat flux as $\beta$ increases until $\beta$ is sufficiently close to $\beta _{\rm crit}^{\rm KBM}(k_y^{\rm min})$, a result that is consistent with ITG nonlinear finite-$\beta$ stabilization. The uptick in $Q_{i}^{\rm es}$ as $\beta$ increases from $0.48\,\%$ to $0.55\,\%$ is an expression of the KBMs, which, unlike the ITG modes, become more virulent as $\beta$ is increased.

Figure 13. Normalized ion heat flux $Q_{i}^{\rm es}$ (black diamonds) and particle flux $\varGamma ^{\rm es}$ (red squares) as a function of $\beta$. Observe the significant reduction of transport for $\beta \gtrsim 0.2\,\%$ relative to $\beta \approx 0.05\,\%$ until $\beta$ approaches the $\beta _{\rm crit}^{\rm KBM} \approx 0.6\,\%$ threshold for $k_y = 0.025$. Above $\beta = 0.6\,\%$, simulations no longer achieve a saturated state. The vertical dashed black line indicates $\beta _{\rm crit}^{\rm KBM}=0.6\,\%$ at $k_y = 0.025$ and the vertical dotted black line indicates $\beta _{\rm crit}^{\rm KBM} = 0.18\,\%$ at $k_y = 0.1$. The electron electromagnetic heat flux is negligible (normalized $Q_{e}^{\rm em} \approx -0.2$ for $\beta = 0.48\,\%$) and therefore not included in this analysis.