3.2.1 Instability in two dimensions

The dispersion relation for the unstable 2-D ($k_\parallel = 0$) modes is (3.2). These modes were studied carefully in Ivanov et al. (Reference Ivanov, Schekochihin, Dorland, Field and Parra2020); let us recap some important points.

The 2-D modes exist at large perpendicular scales, viz., $k_{\perp } < \min \{k_{\perp, \max,\textrm {FLR}}, k_{\perp, \max, \chi }\}$, where the collisionless and collisional cutoffs are given by

(3.7a,b)\begin{equation} k_{{\perp}, \max,\text{FLR}}^2 = \frac{1+2\sqrt{\kappa_T}}{\kappa_T}, \quad k_{{\perp}, \max, \chi}^2 = \sqrt{\frac{\kappa_T}{a \chi^2}}, \end{equation}

respectively. As shown in Ivanov et al. (Reference Ivanov, Schekochihin, Dorland, Field and Parra2020), the Dimits threshold in two dimensions satisfies $\kappa _T \sim \chi$. Here, however, we shall be interested in the strongly driven limit of $\kappa _T \gg \chi$ and $\kappa _T \gg ~1$, for which a saturated state exists only in three dimensions. In this limit, (3.7a) tells us that the cITG modes exist at (and below) wavenumbers

(3.8)\begin{equation} k_{{\perp}} \sim k_{{\perp}, \max,\text{FLR}} \sim \kappa_T^{{-}1/4} \ll 1. \end{equation}

Solving (3.2) shows that these modes also satisfy

(3.9)\begin{equation} {\text{Re}{\left(\omega_{\boldsymbol{k}}\right)}} \sim {\text{Im}{\left(\omega_{\boldsymbol{k}}\right)}} \sim \kappa_T^{1/4}. \end{equation}

3.2.2 $k_\parallel \neq 0$ corrections

Let us now see how $k_\parallel \neq 0$ affects the strongly driven modes at the curvature-driven scales $k_{\perp } \sim \kappa _T^{-1/4} \ll 1$. At these large perpendicular scales, the effects of collisions are negligible, so we may set $\chi = 0$. Note that the scaling $k_{\perp } \sim \kappa _T^{-1/4} \ll 1$ implies $\kappa _T k_{\perp }^2 \sim \sqrt {\kappa _T} \gg 1$, in which case the dispersion (3.1) becomes

(3.10)\begin{equation} \omega_{\boldsymbol{k}}\left[\omega_{\boldsymbol{k}} \left(\omega_{\boldsymbol{k}} + \kappa_Tk_{{\perp}}^2k_y\right) + \kappa_T k_y^2\right] = k_\parallel^2 \left(\omega_{\boldsymbol{k}} + \kappa_T k_y\right), \end{equation}

where the dispersion relation for the curvature-driven 2-D modes is the expression in the square brackets on the left-hand side. Using the results in § 3.2.1, we can estimate that for these modes, the left-hand and right-hand sides of (3.10) satisfy

(3.11a,b)\begin{equation} \omega_{\boldsymbol{k}}\left[\omega_{\boldsymbol{k}} \left(\omega_{\boldsymbol{k}} + \kappa_Tk_{{\perp}}^2k_y\right) + \kappa_T k_y^2\right] \sim \kappa_T^{3/4}, \quad k_\parallel^2 \left(\omega_{\boldsymbol{k}} + \kappa_T k_y\right) \sim k_\parallel^2 \kappa_T^{3/4}. \end{equation}

We can then conclude that for $k_\parallel \ll 1$, the solutions are essentially 2-D, i.e., the dispersion relation (3.10) is well-approximated by (3.2), whereas $k_\parallel \gg 1$ is expected to introduce qualitative changes to the modes. Let us now investigate the $k_\parallel \gg 1$ sITG instability.

3.3.1 Large-scale sITG instability: $k_{\perp } \ll 1$ modes

To lowest order in $k_{\perp } \ll 1$, (3.12) simplifies to

(3.15)\begin{equation} \hat{\omega}_{\boldsymbol{k}}^3 - \hat{k}_{{\parallel}}^2 \hat{\omega}_{\boldsymbol{k}} - \hat{k}_{{\parallel}}^2 = 0. \end{equation}

This is the well-known sITG dispersion relation without FLR effects and in the absence of a density gradient (Cowley et al. Reference Cowley, Kulsrud and Sudan1991). In this limit, the instability boundaries (3.14) become

(3.16a,b)\begin{equation} \hat{k}_{{\parallel},-} = \frac{2}{3\sqrt{3}}k_{{\perp}}^3 + O\left(k_{{\perp}}^5\right), \quad \hat{k}_{{\parallel},+} = \frac{3\sqrt{3}}{2} + \frac{\sqrt{3}}{4}k_{{\perp}}^2 + O\left(k_{{\perp}}^4\right). \end{equation}

For small $\hat {k}_{\parallel }$, the linearly unstable solution of (3.15) is $\hat {\omega }_{\boldsymbol {k}} \approx |\hat {k}_{\parallel }^{2/3}| (-1 + {i}\sqrt {3})/2$. Thus, the linear growth rate for small $\hat {k}_{\parallel }$, or $k_\parallel \ll \kappa _T k_y$, is

(3.17)\begin{equation} {\text{Im}{\left(\omega_{\boldsymbol{k}}\right)}} \approx \frac{\sqrt{3}}{2} \left(\kappa_T k_yk_\parallel^2 \right)^{1/3}. \end{equation}

This is the most widely recognised expression for the sITG growth rate at long wavelengths, however, it is not the fastest-growing mode at $k_y \ll 1$. From (3.16a,b), we know that (3.15) has ${\textrm {Im}{(\hat {\omega }_{\boldsymbol {k}})}} > 0$ solutions up to $\hat {k}_{\parallel } = O(1)$, or $k_\parallel \sim \kappa _T k_y$. The growth rate of these modes evidently satisfies ${\textrm {Im}{(\hat {\omega }_{\boldsymbol {k}})}} = O(1)$, or ${\textrm {Im}{(\omega _{\boldsymbol {k}})}} \sim \kappa _T k_y$.

We are now able to confirm that neglecting the magnetic drift in deriving (3.12), and hence (3.15), was appropriate. As we saw in § 3.2.1, the strongly driven ($\kappa _T \gg 1$) 2-D cITG modes satisfy $k_y \sim \kappa _T^{-1/4}$. At these wavenumbers, the sITG modes exist at scales $k_\parallel \sim \kappa _T^{3/4} \gg 1$ (as expected and assumed) and have a growth rate ${\textrm {Im}{(\omega _{\boldsymbol {k}})}} \sim \kappa _T k_y \sim \kappa _T^{3/4}$, which is asymptotically larger than the growth rate ${\textrm {Im}{(\omega _{\boldsymbol {k}})}}\sim \kappa _T^{1/4}$ of the cITG modes.

3.3.2 Small-scale sITG instability: $k_{\perp } \gg 1$ modes

Expanding (3.13) for $k_{\perp } \gg 1$ and using $\hat {\gamma }_{\boldsymbol {k}}=O\left (k_{\perp }^{-1}\right )$, we find

(3.18)\begin{equation} \hat{k}_{{\parallel},\pm} = k_{{\perp}} (1 \pm 2\hat{\gamma}_{\boldsymbol{k}}) + O\left(k_{{\perp}}^{{-}1}\right). \end{equation}

Therefore, at small perpendicular scales, the sITG is localised at $\hat {k}_{\parallel } = \pm k_{\perp }$, or, equivalently, at

(3.19)\begin{equation} k_\parallel \approx \pm \kappa_T k_y k_{{\perp}}. \end{equation}

In terms of the dimensional wavenumbers, (3.19) tells us that this instability is localised at $k_\parallel L_B/2 \approx \pm \kappa _T k_yk_{\perp } \rho _s^2$, or, equivalently, $k_\parallel L_T \approx \pm k_yk_{\perp }\rho _i^2$. For $\hat {\gamma }_{\boldsymbol {k}}^2 = 1/2k_{\perp }^2$, (3.18) is

(3.20)\begin{equation} \hat{k}_{{\parallel},\pm} = k_{{\perp}} \pm \sqrt{2} + O\left(k_{{\perp}}^{{-}1}\right), \end{equation}

which implies, for the unstable modes,

(3.21)\begin{equation} \left\lvert \hat{k}_{{\parallel}} - k_{{\perp}} \right\rvert = \left\lvert \frac{k_\parallel}{\kappa_T k_y} - k_{{\perp}} \right\rvert < \sqrt{2}. \end{equation}

We can also find the $k_y$ width of the region of instability at fixed $k_\parallel$. Substituting $k_y = k_\parallel / \kappa _T k_{\perp } + \delta k_y$ into (3.21) and expanding for $\delta k_y \ll k_y$, we find

(3.22)\begin{equation} |\delta k_y| < \sqrt{2}\frac{k_{{\perp}} k_y}{|k_{{\perp}}^2 + k_y^2|} \leq \frac{\sqrt{2}}{2}. \end{equation}

To find the growth rate, consider the $k_{\perp } \gg 1$ limit of (3.12) and set $\hat {k}_{\parallel } = \pm k _{\perp } + \delta \hat {k}_{\parallel }$ and $\hat {\omega }_{\boldsymbol {k}} = -1 + \delta \hat {\omega }_{\boldsymbol {k}}$, where $\delta \hat {\omega }_{\boldsymbol {k}} \sim \delta \hat {k}_{\parallel } / k_{\perp } \sim O(k_{\perp }^{-1}) \ll 1$. Keeping terms of order up to $O(k_{\perp }^{-2})$, (3.12) becomes

(3.23)\begin{equation} \left(\delta \hat{\omega}_{\boldsymbol{k}} \pm \frac{\delta \hat{k}_{{\parallel}}}{k_{{\perp}}} \right) \delta \hat{\omega}_{\boldsymbol{k}} + \hat{\gamma}_{\boldsymbol{k}}^2 \approx 0 \implies \delta \hat{\omega}_{\boldsymbol{k}} \approx{-}\frac{\delta \hat{k}_{{\parallel}}}{2k_{{\perp}}} \pm \sqrt{\frac{\delta \hat{k}_{{\parallel}}^2}{4k_{{\perp}}^2} - \hat{\gamma}_{\boldsymbol{k}}^2}. \end{equation}

Thus, in agreement with (3.21), the instability exists only for $|\delta \hat {k}_{\parallel }| < 2k_{\perp }\hat {\gamma }_{\boldsymbol {k}}$, and its growth rate is

(3.24)\begin{equation} {\text{Im}{\left(\hat{\omega}_{\boldsymbol{k}}\right)}} \approx \sqrt{\hat{\gamma}_{\boldsymbol{k}}^2 - \frac{\delta \hat{k}_{{\parallel}}^2}{4k_{{\perp}}^2}}. \end{equation}

The maximum growth rate is then achieved for $\delta \hat {k}_{\parallel } = 0$, i.e., at $k_\parallel =\pm \kappa _T k_yk_{\perp }$, and is given by ${\textrm {Im}{\left (\hat {\omega }_{\boldsymbol {k}}\right )}} \approx \hat {\gamma }_{\boldsymbol {k}}$. Since $\hat {\gamma }_{\boldsymbol {k}}^2 = 1/2k_{\perp }^2$, this is

(3.25)\begin{equation} {\text{Im}{\left(\omega_{\boldsymbol{k}}\right)}} \approx \frac{\kappa_T k_y}{\sqrt{2} k_{{\perp}}} \leq \frac{\kappa_T}{\sqrt{2}}. \end{equation}

The characteristics of the small-scale sITG instability are summarised in figure 3.

Note that for $\kappa _T \gg 1$, the linear growth rate (3.25) of the small-scale sITG modes scales as $O(\kappa _T)$ and, therefore, dominates both the curvature-driven modes (§ 3.2.1) and the large-scale slab modes (§ 3.3.1). This time-scale separation will prove critical for the saturation of strongly driven turbulence (see § 4.2.2).

Finally, an important feature of the $k_{\perp } \gg 1$ sITG modes is the approximate relation $T \approx - \varphi$, or equivalently, $p / \varphi \ll 1$, where $p = \varphi + T$ is the perturbed pressure. Indeed, using (2.12) and (3.23), we find

(3.26)\begin{equation} \frac{T_{\boldsymbol{k}}}{\varphi_{\boldsymbol{k}}} = \frac{\kappa_T k_y}{\omega_{\boldsymbol{k}}} = \frac{1}{\hat{\omega}_{\boldsymbol{k}}} ={-}1 + \frac{\delta \hat{k}_{{\parallel}}}{2 k_{{\perp}}} - {i}\sqrt{\hat{\gamma}_{\boldsymbol{k}}^2 - \frac{\delta \hat{k}_{{\parallel}}^2}{4k_{{\perp}}^2}} + O\left(k_{{\perp}}^{{-}2}\right) \end{equation}

for the modes with ${\textrm {Im}{\left (\hat {\omega }_{\boldsymbol {k}}\right )}} > 0$. Thus, these modes generally have $\,p_{\boldsymbol {k}} / \varphi _{\boldsymbol {k}} \sim O(k_{\perp }^{-1}) \ll 1$, while the most unstable of them ($\delta \hat {k}_{\parallel }=0$) satisfy $\,p_{\boldsymbol {k}} / \varphi _{\boldsymbol {k}} \sim O(k_{\perp }^{-2})$ and ${\textrm {Re}{\left (T_{\boldsymbol {k}} / \varphi _{\boldsymbol {k}}\right )}} = -1 + O(k_{\perp }^{-2})$. This relationship between $T_{\boldsymbol {k}}$ and $\varphi _{\boldsymbol {k}}$ will allow us to identify the sITG modes in the saturated state, and will prove useful in understanding their role in maintaining the Dimits state (see §§ 4.2.1 and 4.2.5).

4.1.1 The 2-D picture

Recall that the 2-D Dimits transition is a sharp transition from a finite-amplitude saturated state with strong ZFs to a ‘blow-up’ state dominated by ever-growing streamers (Ivanov et al. Reference Ivanov, Schekochihin, Dorland, Field and Parra2020). The key to understanding this is the equation for the zonal electrostatic potential

(4.1)\begin{equation} \partial_t \overline{\varphi} + \Pi_\varphi + \Pi_T + \Pi_\chi = 0, \end{equation}

where

(4.2)\begin{equation} \Pi_\varphi \equiv{-} \overline{(\partial_x \varphi) (\partial_y \varphi)}, \quad \Pi_T \equiv{-} \overline{(\partial_x \varphi) (\partial_y T)}, \quad \Pi_\chi \equiv{-} \chi \partial_x^2 \left(a \overline{\varphi} - b \overline{T}\right) \end{equation}

are the Reynolds, diamagnetic, and diffusive stresses, respectively. Equation (4.1) describes how the ZFs are generated or eroded by turbulence (via the Reynolds and diamagnetic stresses, depending on their sign) and damped by collisional viscosity. We then consider a region of nearly constant zonal shear (a ‘shear zone’) of radial width $d$ and find that the integral of the total turbulent stress $\Pi _t = \Pi _\varphi + \Pi _T$ over such a region can be written as

(4.3)\begin{equation} \frac{1}{d}\int \,{{\rm d} x} \Pi_t ={-}\sum_{\boldsymbol{k}} k_xk_y |\varphi_{\boldsymbol{k}}|^2 \left[1 + {\text{Re}{\left(\frac{T_{\boldsymbol{k}}}{\varphi_{\boldsymbol{k}}}\right)}} \right]. \end{equation}

Thus, the effect of the mode with wavenumber $\boldsymbol {k}$ on the ZFs depends on the ratio ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$. Namely, ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} < -1$ implies that the mode will destabilise the ZFs, while ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} > -1$ means that the mode will reinforce the ZFs. This observation is based on the fact that sheared (by the ZFs) turbulence is ‘tilted’ and the sign of $k_x k_y$ is correlated with the sign of the zonal shear in each shear zone. In Ivanov et al. (Reference Ivanov, Schekochihin, Dorland, Field and Parra2020), we derived a simple estimate for the Dimits threshold at large $\kappa _T$ that was based on applying these ideas to the linear modes of the 2-D system. More generally, we argued that the Dimits transition occurred at the threshold of a nonlinear version of the secondary instability – when sheared by ZFs, turbulence either reinforced these flows and thus a Dimits state was maintained (the Reynolds stress won), or it failed to do so (the diamagnetic stress won) and saturation had to be reached via a different route that did not rely on zonal shear. In the 2-D case, no such alternative route for finite-amplitude saturation existed. This description of how a ZF-dominated state was maintained was demonstrated to be accurate by calculating the turbulent viscosity

(4.4)\begin{equation} \nu_t \equiv{-}\frac{ \langle \int_0^{L_x}\,{{\rm d} x} \ \Pi_t S \rangle_{\varDelta t}}{ \langle \int_0^{L_x}\, {{\rm d} x} \ S^2 \rangle_{\varDelta t}} \end{equation}

in numerical simulations with an imposed static ZF profile; here $\langle \dots \rangle _{\varDelta t}$ is a saturated-state time average and $S \equiv \partial _x^2 \overline {\varphi }$ is the zonal shear. Essentially, $\nu _t$ is a measure of the correlation between the turbulent stress $\Pi _t$ and the zonal shear $S$. We found that $\nu _t < 0$ on the Dimits side of the threshold, indicating that sheared turbulence was feeding the ZFs, which were shearing it. Accordingly, we also found $\nu _t > 0$ beyond the threshold, implying that the turbulent stress was actively suppressing the ZFs.

4.1.2 The influence of $L_{\parallel }$ on the Dimits state

Taking the limit $L_{\parallel } \to 0$ effectively restricts our model equations (2.11)–(2.13) to two dimensions, and thus their saturated Dimits state converges to that of the 2-D model. In figure 6(a,b), we show what happens to the turbulent heat flux $Q$ with increasing $L_{\parallel }$ for two cases: far below the 2-D Dimits threshold ($\kappa _T = 0.36$, $\chi = 0.1$), where turbulent bursts dominate the 2-D state; and closer to it ($\kappa _T = 0.8$, $\chi = 0.1$), where the bursts start to overlap in time. As expected, if $L_{\parallel }$ is small enough, we recover the 2-D results. As $L_{\parallel }$ increases, $Q$ converges in a monotonic way to a definite 3-D value that is smaller than the 2-D heat flux. Figure 6(c–f) show that, for larger values of $L_{\parallel }$, the turbulent bursts become more frequent, but shorter in duration and lower in amplitude. There are two effects responsible for this: parallel localisation of turbulence; and the development of the ‘parasitic’ small-scale sITG modes.

Figure 6. (a) Dependence of the time-averaged heat flux $Q$ on the parallel size of the box $L_{\parallel }$ for $\chi = 0.1$, $L_x = L_y = 80$, and $\kappa _T=0.36$. The orange dotted line shows the time-averaged heat flux for the 2-D state ($L_{\parallel } = 0$). (b) Same as (a), but with $\kappa _T = 0.8$. (c–f) Time evolution of the heat flux $Q$ for $\kappa _T = 0.8$, $\chi = 0.1$, $L_x = 80$, $L_y = 80$, and four different values of $L_{\parallel }$ (notated on each panel). As $L_{\parallel }$ increases, the turbulent bursts become more frequent and less violent, and the time-averaged $Q$ drops.

Parallel localisation is inevitable because the turbulent nonzonal modes cannot propagate information infinitely quickly along the field lines. As we increase $L_{\parallel }$ away from 0, we see elongated nonzonal modes that eventually lose the ability to stay coherent along the field lines if $L_{\parallel }$ is large enough. Figure 7 shows that the typical Dimits-state ferdinons are not true 2-D structures and develop a finite parallel extent if the parallel size of the box allows it. This is in contrast with the ZFs, which do stay perfectly coherent along the entire domain regardless of $L_{\parallel }$. This puts the ZFs at an advantage because the turbulent stresses in (4.1) are parallel averages and so a turbulent burst that is localised to a fraction $\varDelta L _{\parallel }/L_{\parallel }$ of the parallel extent of the box has its turbulent stress diminished by a factor of $\varDelta L_{\parallel }/L_{\parallel }$. As we increase $L_{\parallel }$, every such localised burst provides a smaller restoring ‘kick’ to the ZFs and so it takes less time for the ZFs to decay to a level that permits the development of a new burst. The turbulent heat flux $Q$ is also a spatial average of the turbulent fields and it too is diminished for a localised burst. Thus, we expect smaller, more frequent bursts, and this is precisely what is observed. Note that the ability of ZFs to communicate infinitely fast along the field lines is a consequence of the asymptotic limit of small mass ratio and the modified adiabatic electron response (2.1), which is itself due to the assumed infinitely fast parallel electron streaming. Therefore, the inclusion of kinetic electron effects in the equations would lead to qualitative changes for a large enough $L_{\parallel }$. Naturally, this is outside the scope of this work, but is certainly an important consideration for real devices.

Figure 7. Snapshots of the 3-D temperature perturbations associated with a ferdinon. The plots in each row are cross-sections in different planes at the same $t$ taken from simulations that have the same $\kappa _T = 0.36$, $\chi =0.1$, $L_x = L_y = 80$, but (a) $L_{\parallel } = 32$, (b) $L_{\parallel } = 64$, and (c) $L_{\parallel } = 256$. The black dashed lines visualise the intersections of the cross-sectional planes. As we increase $L_{\parallel }$, turbulence loses the ability to stay coherent along the parallel extent of the box and the bursts become localised in $z$.

Secondly, we find small-scale sITG modes that feed off the perpendicular temperature gradients associated with the ferdinons. The presence of this 3-D small-scale ‘parasitic’ instability can be detected via the parallel velocity ${u}$ because the latter is only involved in the 3-D sITG modes and not in the 2-D cITG modes. Figure 8 shows an example of a ferdinon that is ‘infected’ with such small-scale sITG instabilities. As we shall discuss in § 4.2, the small-scale instability leads to an effective increase in thermal diffusion, and thus an increase in the effective damping at large scales that reduces the large-scale temperature perturbations. This additional damping likely contributes to the reduced $Q$ of the 3-D Dimits state. It also enables saturation at finite amplitude when the Dimits state is broken (§ 4.3).

Figure 8. Snapshots of perturbed nonzonal (a) temperature, (b) potential, (c) pressure, and (d) parallel velocity at fixed $z$ in the Dimits state with parameters $\kappa _T = 0.8$, $\chi = 0.1$, $L_x = L_y = 80$, $L_{\parallel } = 1$, parallel hyperviscosity $\nu = 2.4\times 10^{-8}$, and Fourier-space resolution $(n_x, n_y, n_z) = (171, 171, 21)$. The colour scale is relative to the maximum absolute amplitude in each panel (given in the panel's title). Small-scale sITG modes driven by the gradients of the ferdinon are evident in panel (d).

4.2.1 Numerical evidence

Let us now address the small-scale sITG instability. This instability exists only in the 3-D model and is the most important distinction between it and its 2-D counterpart. It is the presence of this instability that enables, in three dimensions and with finite $L_{\parallel }$, the existence of a strongly turbulent saturated state, i.e., one in which there are no strong, coherent ZFs (the zonal profiles of such a state are shown in figure 4d–f). We find that the most distinctive feature of this state is the concentration of pressure perturbations at perpendicular scales that are much larger than the typical (small) scales for the perturbations in $\varphi$ and $T$ (or, to be more precise, the absence of pressure perturbations in the small-scale structure present in $\varphi$ and $T$); this is manifest in figure 9.

Figure 9. Snapshots of perturbed nonzonal (a) temperature, (b) potential, (c) pressure and (d) parallel velocity at fixed $z$ in the strongly turbulent state with parameters $\kappa _T = 3$, $\chi = 0.05$, $L_x = L_y = 80$, $L_{\parallel } = 1$, parallel hyperviscosity $\nu = 1.5\times 10^{-10}$, and Fourier-space resolution $(n_x, n_y, n_z) = (285, 285, 83)$. The colour scale is relative to the maximum absolute amplitude in each panel (given in the panel's title). Time-averaged spectra from the same simulation are shown in figure 10.

Figure 10. Time-averaged spectra (a) $W_{\boldsymbol {k}}$ and (b) $I_{\boldsymbol {k}}$, defined by (2.17) and (2.18), respectively, in the strongly turbulent state with parameters $\kappa _T = 3$, $\chi = 0.05$, $L_x = 80$, $L_y = 80$, and $L_{\parallel } = 1$. The solid black lines demarcate the region of linear instability for $k_x = 0$, and the red dashed line is $k_\parallel = \kappa _T k_{\perp }^2$, where the collisionless modes with largest growth rate reside (see § 3.3.2). We can see that the largest contributions to the two conserved quantities are offset from the region of linear instability. The dotted black line denotes the peak $k_{\perp, I}(k_\parallel )$ of $I_{\boldsymbol {k}}$ at fixed $k_\parallel$. Zonal profiles and cross-sectional snapshots from the same simulation are shown in figures 4 and 9, respectively. The spectra of the saturated state with the same parameters, but with $\kappa _T$ set to $0$ for all $k_\parallel \neq 0$ modes, are given in (c,d). Turning off the equilibrium gradient for the 3-D modes does not alter the spectra noticeably. Snapshots from this modified simulation are shown in figure 11.

In § 3.3.2, we showed that the smallness of the pressure perturbations (compared with the perturbations of the electrostatic potential and temperature) was characteristic of the small-scale ($k_{\perp } \gg 1$) sITG instability, see (3.26). However, the small-scale structure that we see in figure 9 is not produced by the equilibrium-driven instability. In fact, the equilibrium-driven sITG instability is inconsequential in the saturated state. To show this, we ran artificially modified simulations where $\kappa _T$ was set to $0$ for all modes with $k_\parallel \neq 0$ (this is straightforward to do in our spectral code). This removed the equilibrium-driven linear instability from all 3-D ($k_\parallel \neq 0$) modes. Examining the spectra of the two conserved quantities $W_{\boldsymbol {k}}$ and $I_{\boldsymbol {k}}$ (see § 2.2), we see that turning off the equilibrium temperature gradient for the 3-D modes has no noticeable effect on the structure of turbulence (see figure 10). As figure 11 shows, the modified simulations are also visually indistinguishable from the unmodified ones shown in figure 9. The total heat flux $Q$ changes by approximately 20%–30%, likely due to the loss of radial-symmetry breaking for the 3-D modes, which are now free to transport heat in either direction equally, so, on average, they have zero radial heat flux. The nonlinear interactions between the 2-D ($k_\parallel = 0$) modes cannot produce the 3-D modes that we see in the modified simulations. Therefore, these 3-D modes must be produced by a ‘parasitic’ sITG instability of the 2-D fields (into which energy is injected by the equilibrium gradient).

Figure 11. Same as figure 9, but for a modified simulation, i.e., with $\kappa _T$ set to zero for the $k_\parallel \neq 0$ modes. Visually, the saturated state is identical to that shown in figure 9.

Furthermore, the spectra of $W_{\boldsymbol {k}}$ and $I_{\boldsymbol {k}}$ measured in regular simulations are inconsistent with the region of linear instability of the dispersion relation (3.1). Namely, figure 10(a,b) show that $W_{\boldsymbol {k}}$ and $I_{\boldsymbol {k}}$ of the linearly unstable modes of (3.1) are orders of magnitude smaller than the corresponding spectral peaks of the two conserved quantities. We can quantify this by using the turbulent spectra to determine the ‘dominant’ perpendicular scale as a function of the parallel scale $k_\parallel$. As figure 10(a,b) show, this is the scale at which $I_{\boldsymbol {k}}$ peaks and the dependence of $W_{\boldsymbol {k}}$ on $k_{\perp }$ changes from flat to steeply declining. To extract this scale, we define $k_{\perp, I} (k_\parallel )$ as the $k_{\perp }$ that maximises $I_{\boldsymbol {k}}$ at a fixed $k_\parallel$. Figure 12(a) shows that $k_{\perp, I} (k_\parallel )$ lies outside of the region of linear instability of (3.1). Thus, the 3-D structure of the saturated state is not produced by the linear sITG instability driven by the equilibrium gradient.

Figure 12. (a) Comparison of the location of the spectral peak $k_{\perp, I} (k_\parallel )$ of $I_{\boldsymbol {k}}$ (blue line) and the location of the peak of the growth rate of the collisionless linear instability driven by the equilibrium gradient (orange dashed line), given by $k_\parallel = \kappa _T k_{\perp }^2$. The black curve circumscribes the region of linear instability, i.e., all ${\textrm {Im}{\left (\omega _{\boldsymbol {k}}\right )}} \geq 0$ solutions to (3.1) are inside it and outside of it, all solutions satisfy ${\textrm {Im}{\left (\omega _{\boldsymbol {k}}\right )}} \leq 0$. (b) Comparison of the equilibrium temperature gradient $\kappa _T$ and the ‘effective’ temperature gradient $\kappa _T^\textrm {eff}(k_\parallel ) \equiv k_\parallel / k_{\perp, I}^2$. The data is from the same simulation as shown in figure 9. The spectra of this simulation are given in figure 10. This rough estimate of $\kappa _T^\textrm {eff}$ being approximately 5–10 times larger than $\kappa _T$ is consistent with the calculated growth rate of the parasitic small-scale instability (see figure 13).

In § 3.3.2, we showed that the equilibrium-driven sITG instability is localised at $k_\parallel \approx \kappa _T k_{\perp }^2$. A similar relationship holds for $k_{\perp, I} (k_\parallel )$, viz., $k_\parallel \approx \kappa _T^\textrm {eff} k_{\perp, I}^2$, where $\kappa _T^\textrm {eff}$ can be thought of as an effective temperature gradient. Figure 12(b) shows that this $\kappa _T^\textrm {eff}$ is several times larger than the equilibrium temperature gradient. As we shall see shortly, $\kappa _T^\textrm {eff}$ is actually the gradient of the large-scale 2-D temperature perturbations.

4.2.2 Scale-separated equations for cITG and sITG modes

The numerical analysis above leads us to believe that the 3-D structure of the saturated state is a consequence of an instability driven not by the equilibrium gradient $\kappa _T$, but rather by the gradients of the 2-D perturbations. Let us attack on the analytical front. As we discussed at the start of § 4, the 3-D sITG modes are naturally scale-separated from the 2-D cITG modes both in wavenumber and in frequency. We introduce the parallel average

(4.5)\begin{equation} \left\langle f \right\rangle_\parallel{\equiv} \int \frac{dz'}{L_z} f(z'). \end{equation}

This average allows us to split (2.11)–(2.13) into separate equations for the slow 2-D modes governed by the cITG instability at large perpendicular scales ($k_{\perp } \ll 1$), and for the fast 3-D sITG modes, which live at small perpendicular scales ($k_{\perp } \gg 1$).Footnote ⁴ We define the small-scale perturbations as $\widetilde {f} \equiv f - \left \langle f \right \rangle _\parallel$. The large-scale equations are then

(4.6)\begin{gather} \partial_t \left( \left\langle \varphi' \right\rangle_\parallel{-} \nabla_{{\perp}}^2 \left\langle \varphi \right\rangle_\parallel\right) - \partial_y \left( \left\langle \varphi \right\rangle_\parallel{+} \left\langle T \right\rangle_\parallel \right) + \kappa_T \partial_y \nabla_{{\perp}}^2 \left\langle \varphi \right\rangle_\parallel \nonumber\\ +\left\{ \left\langle \varphi \right\rangle_\parallel, \left\langle \varphi' \right\rangle_\parallel{-} \nabla_{{\perp}}^2 \left\langle \varphi \right\rangle_\parallel \right\} + \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \left\{ \boldsymbol{\nabla}_\perp \left\langle \varphi \right\rangle_\parallel, \left\langle T \right\rangle_\parallel \right\} + \chi \nabla_\perp^4 \left(a\left\langle \varphi \right\rangle_\parallel{-} b\left\langle T \right\rangle_\parallel\right) \nonumber\\ ={-}\left\langle \left\{ \widetilde{\varphi}, \widetilde{\varphi}' - \nabla_{{\perp}}^2 \widetilde{\varphi} \right\} - \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \left\{ \boldsymbol{\nabla}_\perp \widetilde{\varphi}, \widetilde{T} \right\} \right\rangle_\parallel, \end{gather}

(4.7)\begin{gather}\partial_t\left\langle T \right\rangle_\parallel{+}\kappa_T\partial_y \left\langle \varphi \right\rangle_\parallel{+}\left\{ \left\langle \varphi \right\rangle_\parallel, \left\langle T \right\rangle_\parallel \right\} - \chi \nabla_{{\perp}}^2 \left\langle T \right\rangle_\parallel = {-}\!\left\langle \left\{ \widetilde{\varphi}, \widetilde{T} \right\} \right\rangle_\parallel, \end{gather}

(4.8)\begin{gather}\partial_t \left\langle {u} \right\rangle_\parallel{+} \left\{ \left\langle \varphi \right\rangle_\parallel, \left\langle {u} \right\rangle_\parallel \right\} - s \chi \nabla_{{\perp}}^2 \left\langle {u} \right\rangle_\parallel = {-}\!\left\langle \left\{ \widetilde{\varphi}, \widetilde{{u}} \right\} \right\rangle_\parallel. \end{gather}

The right-hand sides of (4.6)–(4.8) represent the influence of the 3-D sITG modes on the large-scale fields. Subtracting (4.6)–(4.8) from (2.11)–(2.13), we find the small-scale equations:

(4.9)\begin{gather} \partial_t \left( \widetilde{\varphi}' - \nabla_{{\perp}}^2 \widetilde{\varphi}\right) + \partial_\parallel \widetilde{{u}} - \partial_y (\widetilde{\varphi} + \widetilde{T}) + \kappa_T \partial_y \nabla_{{\perp}}^2 \widetilde{\varphi} + \widetilde{\left\{ \widetilde{\varphi}, \widetilde{\varphi}' - \nabla_{{\perp}}^2 \widetilde{\varphi} \right\} } + \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \widetilde{\Big\{ \boldsymbol{\nabla}_\perp \widetilde{\varphi}, \widetilde{T} \Big\}} \nonumber\\ + \chi \nabla_\perp^4 \left(a\widetilde{\varphi} - b\widetilde{T}\right) ={-}\left[\Big\{ \left\langle \varphi \right\rangle_\parallel, \widetilde{\varphi}' - \nabla_{{\perp}}^2 \widetilde{\varphi} \Big\} + \Big\{ \widetilde{\varphi}, \left\langle \varphi' \right\rangle_\parallel{-} \nabla_{{\perp}}^2 \left\langle \varphi \right\rangle_\parallel \Big\} \right.\nonumber\\ \left.+ \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \Big\{ \boldsymbol{\nabla}_\perp \left\langle \varphi \right\rangle_\parallel, \widetilde{T} \Big\} + \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \Big\{ \boldsymbol{\nabla}_\perp \widetilde{\varphi}, \left\langle T \right\rangle_\parallel \Big\}\right], \end{gather}

(4.10)\begin{gather}\partial_t\widetilde{T} +\kappa_T\partial_y \widetilde{\varphi} +\widetilde{\Big\{ \widetilde{\varphi}, \widetilde{T} \Big\}} - \chi \nabla_{{\perp}}^2 \widetilde{T} ={-}\left[\Big\{ \left\langle \varphi \right\rangle_\parallel, \widetilde{T} \Big\} + \Big\{ \widetilde{\varphi}, \left\langle T \right\rangle_\parallel \Big\}\right], \end{gather}

(4.11)\begin{gather}\partial_t \widetilde{{u}} + \partial_\parallel \left(\widetilde{\varphi} + \widetilde{T}\right) + \widetilde{\left\{ \widetilde{\varphi}, \hat{u} \right\}} - s \chi \nabla_{{\perp}}^2 \hat{u} ={-}\left[\Big\{ \left\langle \varphi \right\rangle_\parallel, \widetilde{{u}} \Big\} + \Big\{ \widetilde{\varphi}, \left\langle {u} \right\rangle_\parallel \Big\}\right]. \end{gather}

In order to simplify the following analysis, we shall assume both temporal and spatial scale separation between (4.6)–(4.8) and (4.9)–(4.11), i.e., that the large-scale fields are constant in time in (4.9)–(4.11) and that the spatial and temporal scales of (4.9)–(4.11) are short compared with the respective scales of (4.6)–(4.8). In particular, we shall assume that the perpendicular scales of the 2-D modes are sufficiently large for the derivatives of their gradients to be ignored. This assumption turns out to be equivalent to $k_{\perp } q_\perp \ll 1$, where $k_{\perp }$ and $q_\perp$ are the typical perpendicular wavenumbers of the 2-D and parasitic modes, respectively. According to (3.8) and (3.33), the linearly unstable modes satisfy $k_{\perp } \sim \kappa _T^{-1/4}$ and $q_\perp \sim (\kappa _T/\chi )^{1/3}$ in the limit $\kappa _T \gg \chi$. Therefore, the condition $k_{\perp } q_\perp \ll 1$ is equivalent to $\chi \gg \kappa _T^{1/4}$. Additionally, recall that the Dimits threshold in two dimensions is found at $\kappa _T \sim \chi$ (Ivanov et al. Reference Ivanov, Schekochihin, Dorland, Field and Parra2020) and that, as we showed in § 4.1, the 2-D Dimits regime is qualitatively unchanged when we include 3-D effects. Thus, for the remainder of this section, we shall consider the limit

(4.12)\begin{equation} \kappa_T \gg \chi \gg \kappa_T^{1/4} \gg 1, \end{equation}

which puts us beyond the 2-D Dimits transition (i.e., in two dimensions, such a state blows up). Importantly, we limit our analysis to $\kappa _T/\chi \ll 830$, in which case the $\chi$ITG instability can safely be neglected (see Appendix C). The limit (4.12) then allows us to simplify (4.6)–(4.8) and (4.9)–(4.11) significantly, and thus to describe the interplay between 2-D and parasitic modes analytically. These analytical results agree with our simulations, even though the latter do not strictly conform to (4.12).

4.2.3 Parasitic sITG instability

First, we investigate the small-scale sITG instability in the presence of large-scale 2-D modes. Linearising (4.9)–(4.11) in the limit (4.12), we obtain

(4.13)\begin{gather} \left(\partial_t + \left\langle \boldsymbol{V_{E}} \right\rangle_\parallel \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \right) \left( \widetilde{\varphi}' - \nabla_{{\perp}}^2 \widetilde{\varphi}\right) + \partial_\parallel \widetilde{{u}} - \partial_y (\widetilde{\varphi} + \widetilde{T}) \nonumber\\ \quad +\, \boldsymbol{\kappa}_T \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \nabla_{{\perp}}^2 \widetilde{\varphi} + \boldsymbol{\kappa}_n \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \widetilde{\varphi} ={-} \chi \nabla_\perp^4 (a\widetilde{\varphi} - b\widetilde{T}), \end{gather}

(4.14)\begin{gather}\left(\partial_t + \left\langle \boldsymbol{V_{E}} \right\rangle_\parallel \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \right)\widetilde{T} +\boldsymbol{\kappa}_T \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \widetilde{\varphi} = \chi \nabla_{{\perp}}^2 \widetilde{T}, \end{gather}

(4.15)\begin{gather}\left(\partial_t + \left\langle \boldsymbol{V_{E}} \right\rangle_\parallel \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \right) \widetilde{{u}} + \partial_\parallel \left(\widetilde{\varphi} + \widetilde{T}\right) = s \nabla_{{\perp}}^2 \widetilde{{u}}, \end{gather}

where the ‘local-equilibrium’ quantities

(4.16)\begin{equation} \left\langle \boldsymbol{V_{E}} \right\rangle_\parallel \equiv \hat{\boldsymbol{z}} \times \boldsymbol{\nabla}_\perp \left\langle \varphi \right\rangle_\parallel, \quad \boldsymbol{\kappa}_n \equiv{-} \hat{\boldsymbol{z}} \times \boldsymbol{\nabla}_\perp \left\langle \varphi' \right\rangle_\parallel, \quad \boldsymbol{\kappa}_T \equiv \kappa_T \hat{\boldsymbol{y}} - \hat{\boldsymbol{z}} \times \boldsymbol{\nabla}_\perp \left\langle T \right\rangle_\parallel \end{equation}

are the $\boldsymbol {E}\times \boldsymbol {B}$ advecting flow, the local density gradient, and the total local temperature gradient (large-scale perturbation plus equilibrium), respectively. We assume that $|\boldsymbol {\kappa }_T|\sim |\boldsymbol {\kappa }_n|$ (see § 4.2.4). Note that only the nonzonal electrostatic potential $\varphi '$ gives rise to a density perturbation – this is a consequence of the modified adiabatic electron response (2.1). Note also that we have ignored the large-scale 2-D parallel flow $\left \langle {u} \right \rangle _\parallel$. Since $\left \langle {u} \right \rangle _\parallel$ is not involved in any linear instability, the only way it could be driven is via the small-scale response, viz., the right-hand side of (4.8). In Appendix E, we show that a small initial $\left \langle {u} \right \rangle _\parallel$ will decay under the influence of growing small-scale modes. Accordingly, in our numerical simulations, we find that $\left \langle u \right \rangle _\parallel$ is many orders of magnitude smaller than the other two 2-D fields and is irrelevant for the saturated state.

Ignoring collisions (i.e., setting $\chi = 0$) and taking the gradients of the large-scale fields to be constant over the small scales at which (4.13)–(4.15) hold, we can investigate the small-scale linear instability in a way analogous to what we did in § 3. In particular, we shall focus on the $k_\parallel \sim k_{\perp }^2 \gg 1$ regime analysed in § 3.3.2. We look for Doppler-shifted solutions to (4.13)–(4.15) of the form $\widetilde {\varphi }_{\boldsymbol {k}}, \widetilde {T}_{\boldsymbol {k}}, \widetilde {{u}}_{\boldsymbol {k}} \propto \exp {[-{i}(\omega _{\boldsymbol {k}} + \left \langle \boldsymbol {V_{E}} \right \rangle _\parallel\ \boldsymbol{\cdot }\ \boldsymbol {k} ) t + {i}\boldsymbol {k} \boldsymbol{\cdot} \boldsymbol {r}]}$. Note that we ignore the shear in the $\boldsymbol {E}\times \boldsymbol {B}$ flow $\left \langle \boldsymbol {V_{E}} \right \rangle _\parallel$. We also ignore the magnetic-drift term $-\partial _y(\widetilde {\varphi } + \widetilde {T})$ in (4.13) because it is subdominant for the sITG modes with $k_{\perp } \gg 1$. The resulting dispersion relation for these modes is

(4.17)\begin{equation} \left(\omega_{\boldsymbol{k}}^2 - \frac{k_\parallel^2}{1 + k_{{\perp}}^2}\right) \left( \omega_{\boldsymbol{k}} + \boldsymbol{\kappa}_T \boldsymbol{\cdot} \boldsymbol{k} \right) = \frac{\omega_{\boldsymbol{k}}^2}{1 + k_{{\perp}}^2} \left[ \left(\boldsymbol{\kappa}_n + \boldsymbol{\kappa}_T\right) \boldsymbol{\cdot} \boldsymbol{k} \right]. \end{equation}

Since (4.13)–(4.15) describe real fields, (4.17) must be invariant under $\boldsymbol {k} \mapsto -\boldsymbol {k}$ and $\omega _{\boldsymbol {k}} \mapsto -\omega _{\boldsymbol {k}}^*$. We may, therefore, assume that $\boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k} > 0$ without loss of generality. Repeating the arguments of § 3.3.2, we define $\omega _{\boldsymbol {k}} \equiv \hat {\omega }_{\boldsymbol {k}}\boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k}$ and $k_\parallel \equiv \hat {k}_{\parallel } \boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k}$. Then (4.17) turns out to be formally the same as our old dispersion relation (3.12), but now with

(4.18)\begin{equation} \hat{\gamma}_{\boldsymbol{k}}^2 = \frac{(\boldsymbol{\kappa}_n+\boldsymbol{\kappa}_T) \boldsymbol{\cdot} \boldsymbol{k}}{2k_{{\perp}}^2\boldsymbol{\kappa}_T \boldsymbol{\cdot} \boldsymbol{k}}. \end{equation}

Thus, the results of § 3.3.2 carry over to the parasitic instability described by (4.17). In particular, the sITG instability exists if $\hat {\gamma }_{\boldsymbol {k}}^2 > 0$, i.e., if $(\boldsymbol {\kappa }_n+\boldsymbol {\kappa }_T) \boldsymbol{\cdot} \boldsymbol {k}$ and $\boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k}$ have the same sign, and is localised to $k_\parallel \approx \pm \boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k} k_{\perp }$. Its growth rate is given by

(4.19)\begin{equation} {\text{Im}{\left(\omega_{\boldsymbol{k}}\right)}} \approx \text{Re} \sqrt{\frac{\boldsymbol{\kappa}_T \boldsymbol{\cdot} \hat{\boldsymbol{k}} \left(\boldsymbol{\kappa}_n + \boldsymbol{\kappa}_T\right) \boldsymbol{\cdot} \hat{\boldsymbol{k}}}{2}}, \end{equation}

where $\hat {\boldsymbol {k}} = \boldsymbol {k} / k_{\perp }$. As expected, this is the same as (3.25) if $\boldsymbol {\kappa }_n = 0$ and $\boldsymbol {\kappa }_T = \kappa _T \hat {\boldsymbol {y}}$. In figure 13(b), we show the maximum growth rate obtained from the numerical solution of the full (with collisionality and magnetic curvature turned back on) dispersion relation (3.1) with the addition of the local temperature and density gradients of the large-scale fields. As expected from the numerical analysis in § 4.2.1, the small-scale instability driven by the large-scale gradients is significantly (${\approx }5$ times in this case) stronger than the equilibrium-driven instability. This is consistent with the estimate of the effective temperature gradient $\kappa _T^\textrm {eff}$ for the sITG instability that we showed in figure 12(b).

Figure 13. (a) Snapshot of the 2-D temperature perturbation $\left \langle T \right \rangle _\parallel$ in the $(x, y)$ plane. The data is taken from the same $\kappa _T = 3$, $\chi = 0.05$ simulation that we showed in figure 9. The 2-D temperature perturbations lack the small-scale structure that was seen in figure 9(a), confirming that the parallel average (4.5) removes small-scale perpendicular structure. (b) Small-scale growth rate in the $(x, y)$ plane. This plot is obtained by finding the maximum growth rate of the full (including collisionality and magnetic curvature) dispersion relation (3.1) with the addition of the local temperature and density gradients of the large-scale fields at every point. For this simulation, $\kappa _T = 3$, and so the largest collisionless growth rate, given by (3.25), is $\kappa _T / \sqrt {2} \approx 2.1$. It is thus evident that the influence of the gradients of the large-scale fields dominates over that of the equilibrium gradient $\kappa _T$ by a factor of $5$. The ‘effective’ $\kappa _T^\textrm {eff}$ that we estimated for the same simulation in figure 12(b) is, indeed, a factor of 5–10 larger that the equilibrium gradient $\kappa _T$.

Note that if $\boldsymbol {\kappa }_n \neq 0$, (4.19) implies that modes with $\boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k} \left (\boldsymbol {\kappa }_n + \boldsymbol {\kappa }_T\right ) \boldsymbol{\cdot} \boldsymbol {k} < 0$ are linearly stable. Is there a $\boldsymbol {\kappa }_n$ that quenches the sITG instability for all $\boldsymbol {k}$? Suppose $\boldsymbol {\kappa }_n\nparallel \boldsymbol {\kappa }_T$; then we can choose $\hat {\boldsymbol {k}} \boldsymbol{\cdot} \boldsymbol {\kappa }_n=0$, but $\hat {\boldsymbol {k}} \boldsymbol{\cdot} \boldsymbol {\kappa }_T\neq 0$. By (4.19), any such $\hat {\boldsymbol {k}}$ is an unstable mode. Therefore, to stabilise all modes, we require $\boldsymbol {\kappa }_n\parallel \boldsymbol {\kappa }_T$. In this case, it is evident that $\boldsymbol {\kappa }_n \boldsymbol{\cdot} \hat {\boldsymbol {k}} = (\boldsymbol {\kappa }_n \boldsymbol{\cdot} \boldsymbol {\kappa }_T)(\boldsymbol {\kappa }_T \boldsymbol{\cdot} \hat {\boldsymbol {k}})/|\boldsymbol {\kappa }_T|^2$. Therefore, in order to quench the sITG instability for all $\boldsymbol {k}$, we need

(4.20)\begin{equation} \boldsymbol{\kappa}_n \parallel \boldsymbol{\kappa}_T, \quad \frac{\boldsymbol{\kappa}_n \boldsymbol{\cdot} \boldsymbol{\kappa}_T}{|\boldsymbol{\kappa}_T|^2} \leq{-}1. \end{equation}

We shall now show that the effect of the growing small-scale modes on the large-scale 2-D fields can be expressed as an enhanced thermal diffusivity for the latter.

4.2.4 Anomalous heat flux due to parasitic sITG modes

We expect that the growth of small-scale sITG modes, which are driven by the gradients associated with the large-scale fluctuations, will check the growth of the amplitudes of the driving large-scale fields. This is an intuitive consequence of the conservation laws (2.14) and (2.16). As the parasitic instability is driven by the nonlinear terms that conserve $W=\sum _{\boldsymbol {k}} W_{\boldsymbol {k}}$ and $I =\sum _{\boldsymbol {k}} I_{\boldsymbol {k}}$, an excitation of parasitic small-scale modes should show up as a sink in the large-scale equations. Let us now calculate explicitly the influence of small-scale sITG modes on the large-scale modes and show that this is indeed true. This influence is represented by the terms of the form $\left \langle \left \{ ., . \right \} \right \rangle _\parallel$ on the right-hand sides of (4.6)–(4.8).

First, consider the temperature equation (4.7). The relevant term is

(4.21)\begin{equation} -\Big\langle \Big\{ \widetilde{\varphi}, \widetilde{T} \Big\} \Big\rangle_\parallel = {-}\boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \Big\langle \Big[ (\hat{\boldsymbol{z}}\times\boldsymbol{\nabla}_\perp\widetilde{\varphi}) \widetilde{T} \Big] \Big\rangle_\parallel{\equiv} -\boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \langle\widetilde{\boldsymbol{Q}}\rangle_\parallel, \end{equation}

where $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ is the turbulent heat flux associated with the small-scale modes. Let us compute it quasilinearly (i.e., assuming that $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ is determined by the most unstable small-scale modes), assuming scale separation. As stated in § 4.2.2, we imagine that the small-scale equations are solved in an infinitesimal (compared with the large scales) box, and thus the parallel average is equivalent to an average over such a small-scale box,

(4.22)\begin{equation} \langle\widetilde{\boldsymbol{Q}}\rangle_\parallel = {-}\sum_{\boldsymbol{q}} {i} \hat{\boldsymbol{z}}\times\boldsymbol{q}\widetilde{\varphi}_{-\boldsymbol{q}} \widetilde{T}_{\boldsymbol{q}} \approx{-}\sum_{\boldsymbol{q}} \hat{\boldsymbol{z}}\times\hat{\boldsymbol{q}} \sqrt{\frac{(\boldsymbol{\kappa}_T + \boldsymbol{\kappa}_n) \boldsymbol{\cdot} \hat{\boldsymbol{q}}}{2\boldsymbol{\kappa}_T \boldsymbol{\cdot} \hat{\boldsymbol{q}}}} |\widetilde{\varphi}_{\boldsymbol{q}}|^2, \end{equation}

where $\hat {\boldsymbol {q}} = \boldsymbol {q} / q_\perp$ and we have assumed that the sum is dominated by the wavenumbers $\boldsymbol {q}$ corresponding to the largest linear growth rate of the parasitic sITG instability, and so have replaced $\widetilde {T}_{\boldsymbol {q}}/\widetilde {\varphi }_{\boldsymbol {q}}$ with the collisionless expression (3.26) for the modes with $k_\parallel = \boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k}k_{\perp }$ that maximise this growth rate. Note that the small-scale fields $\widetilde {T}_{\boldsymbol {q}}$ and $\widetilde {\varphi }_{\boldsymbol {q}}$, and thus $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ itself, depend implicitly on the position variable of the large-scale equations (4.6)–(4.8).

In order to verify that $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ does indeed damp the large-scale temperature perturbations $\left \langle T \right \rangle _\parallel$, we multiply (4.7) by $\left \langle T \right \rangle _\parallel$ and integrate over space to find

(4.23)\begin{align} & \partial_t \int \,{\rm d}^3\boldsymbol{r} \frac{1}{2} \left\langle T \right\rangle_\parallel^2 + \text{linear terms} = \int \,{\rm d}^3\boldsymbol{r} \langle\widetilde{\boldsymbol{Q}}\rangle_\parallel \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \left\langle T \right\rangle_\parallel \nonumber\\ & \quad\approx{-} \int \,{\rm d}^3\boldsymbol{r} \sum_{\boldsymbol{q}} (\hat{\boldsymbol{z}}\times\hat{\boldsymbol{q}}) \boldsymbol{\cdot} \boldsymbol{\nabla}_\perp \left\langle T \right\rangle_\parallel \sqrt{\frac{(\boldsymbol{\kappa}_T + \boldsymbol{\kappa}_n) \boldsymbol{\cdot} \hat{\boldsymbol{q}}}{2\boldsymbol{\kappa}_T \boldsymbol{\cdot} \hat{\boldsymbol{q}}}} |\widetilde{\varphi}_{\boldsymbol{q}}|^2 \nonumber\\ & \quad={-} \int \,{\rm d}^3\boldsymbol{r} \sum_{\boldsymbol{q}} \sqrt{\frac{\boldsymbol{\kappa}_T \boldsymbol{\cdot} \hat{\boldsymbol{q}}(\boldsymbol{\kappa}_T + \boldsymbol{\kappa}_n) \boldsymbol{\cdot} \hat{\boldsymbol{q}}}{2}} |\widetilde{\varphi}_{\boldsymbol{q}}|^2 ={-} \int \,{\rm d}^3\boldsymbol{r} \sum_{\boldsymbol{q}} {\text{Im}{\left(\omega_{\boldsymbol{k}}\right)}} |\widetilde{\varphi}_{\boldsymbol{q}}|^2, \end{align}

where ${\textrm {Im}{\left (\omega _{\boldsymbol {k}}\right )}}$ is the sITG growth rate (4.19). Thus, the linearly unstable small-scale modes have a sign-definite effect on $\left \langle T \right \rangle _\parallel$: they provide additional dissipation.

The heat flux (4.22) depends on $\left \langle T \right \rangle _\parallel$ in a nontrivial way. Let us quantify its influence on $\left \langle T \right \rangle _\parallel$ by working out its direction as a function of $\boldsymbol {\kappa }_T$. Let us assume that $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ is dominated by the fastest-growing sITG modes, and let their wavevector direction be $\hat {\boldsymbol {q}}_\textrm {max}$, so $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ is parallel to $\hat {\boldsymbol {z}}\times \hat {\boldsymbol {q}}_\textrm {max}$. In figure 14, we illustrate the influence on $\boldsymbol {\kappa }_T$ of the contribution to $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ from the most unstable small-scale modes. As expected, we find that the turbulent heat flux due to the small-scale modes pushes the large-scale gradient $\boldsymbol {\kappa }_T$ towards the linearly stable configuration (4.20).

Figure 14. This plot shows the direction in which the heat flux (4.22) of the most unstable small-scale mode ($\hat {\boldsymbol {q}} = \hat {\boldsymbol {q}}_{\textrm {max}}$) pushes the temperature gradient $\boldsymbol {\kappa }_T = -\hat {\boldsymbol {z}}\times \boldsymbol {\nabla }_\perp \left \langle T \right \rangle _\parallel$. We have chosen a coordinate system in which the large-scale density gradient is $\boldsymbol {\kappa }_n = (0, 1)$, denoted by the green arrow. The red line shows the values of $\boldsymbol {\kappa }_T$ for which the sITG instability has zero growth rate, according to (4.20). The black arrows represent the direction of $-\hat {\boldsymbol {z}}\times \langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$. We see that $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$ pushes the large-scale temperature gradient $\boldsymbol {\kappa }_T$ towards the linearly stable region.

Now consider (4.6), the evolution equation for $\left \langle \varphi \right \rangle _\parallel$. The relevant nonlinear terms are

(4.24)\begin{align} & \left\langle \left\{ \widetilde{\varphi}, \widetilde{\varphi}' - \nabla_{{\perp}}^2 \widetilde{\varphi} \right\} + \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \left\{ \boldsymbol{\nabla}_\perp \widetilde{\varphi}, \widetilde{T} \right\} \right\rangle_\parallel = \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \left\langle \left\{ \boldsymbol{\nabla}_\perp \widetilde{\varphi}, \widetilde{p} \right\} \right\rangle_\parallel \nonumber\\ & \quad= \boldsymbol{\nabla}_\perp \boldsymbol{\nabla}_\perp \boldsymbol{:} \sum_{\boldsymbol{q}} (\hat{\boldsymbol{z}}\times\boldsymbol{q})\boldsymbol{q} \left(1 + {\text{Re}{\frac{\widetilde{T}_{\boldsymbol{q}}}{\widetilde{\varphi}_{\boldsymbol{q}}}}}\right) |\widetilde{\varphi}_{\boldsymbol{q}}|^2 \equiv \boldsymbol{\nabla}_\perp \boldsymbol{\nabla}_\perp \boldsymbol{:} \langle\widetilde{\Pi}\rangle_\parallel. \end{align}

The collisionless calculations of § 3.3.2 are straightforward to generalise for the collisionless parasitic small-scale instability. They yield the same relation for ${\textrm {Re}{(\widetilde {T}_{\boldsymbol {q}} / \widetilde {\varphi }_{\boldsymbol {q}})}}$, viz., $1 + {\textrm {Re}{(\widetilde {T}_{\boldsymbol {q}} / \widetilde {\varphi }_{\boldsymbol {q}})}} = O(q_\perp ^{-2})$. However, as we shall discuss in § 4.2.5, the presence of nonzero $\chi$ alters this to $1 + {\textrm {Re}{(\widetilde {T}_{\boldsymbol {q}} / \widetilde {\varphi }_{\boldsymbol {q}})}} = O(q_\perp ^{-1})$. Assuming therefore that the dominant parasitic modes satisfy $1 + {\textrm {Re}{(\widetilde {T}_{\boldsymbol {q}}/\widetilde {\varphi }_{\boldsymbol {q}})}} \lesssim O(q_\perp ^{-1})$, we find $\boldsymbol {\nabla }_\perp \boldsymbol {\nabla }_\perp \boldsymbol {:} \langle \widetilde {\Pi }\rangle _\parallel \lesssim k_{\perp }^2 q_\perp |\hat {\varphi }|^2$. However, (4.22) implies $\boldsymbol {\nabla }_\perp \boldsymbol{\cdot} \langle \widetilde {\boldsymbol {Q}}\rangle _\parallel \sim k_{\perp } |\hat {\varphi }|^2$, and so

(4.25)\begin{equation} \frac{\boldsymbol{\nabla}_\perp \boldsymbol{\nabla}_\perp \boldsymbol{:} \langle\widetilde{\Pi}\rangle_\parallel}{\boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \langle\widetilde{\boldsymbol{Q}}\rangle_\parallel} \lesssim k_\perp q_\perp\sim O\left[\left(\frac{\kappa_T^{1/4}}{\chi}\right)^{1/3}\right] \ll 1, \end{equation}

in line with the assumption on scales formulated at the end of § 4.2.2. Therefore, assuming that $\left \langle \varphi ' \right \rangle _\parallel \sim \left \langle T \right \rangle _\parallel$Footnote ⁵ and that they evolve on the same time scale, we conclude that the main effect of the small-scale modes is to provide a feedback to the large-scale temperature in the form of the additional heat flux $\langle \widetilde {\boldsymbol {Q}}\rangle _\parallel$.

We can thus summarise the equations that govern the evolution of $\left \langle \varphi ' \right \rangle _\parallel$ and $\left \langle T \right \rangle _\parallel$ as

(4.26)\begin{gather} \partial_t \left\langle \varphi' \right\rangle_\parallel{-} \partial_y \left( \left\langle \varphi' \right\rangle_\parallel{+} \left\langle T' \right\rangle_\parallel \right) + \kappa_T \partial_y \nabla_{{\perp}}^2 \left\langle \varphi' \right\rangle_\parallel \nonumber\\ +\left\{ \left\langle \varphi \right\rangle_\parallel, \left\langle \varphi' \right\rangle_\parallel{-} \nabla_{{\perp}}^2 \left\langle \varphi \right\rangle_\parallel \right\} ' + \boldsymbol{\nabla}_\perp \boldsymbol{\cdot} \left\{ \boldsymbol{\nabla}_\perp \left\langle \varphi \right\rangle_\parallel, \left\langle T \right\rangle_\parallel \right\} ' + \chi \nabla_\perp^4 \left(a\left\langle \varphi' \right\rangle_\parallel{-} b\left\langle T' \right\rangle_\parallel\right) = 0, \end{gather}

(4.27)\begin{gather}\partial_t\left\langle T \right\rangle_\parallel{+}\kappa_T\partial_y \left\langle \varphi \right\rangle_\parallel{+}\left\{ \left\langle \varphi \right\rangle_\parallel, \left\langle T \right\rangle_\parallel \right\} - \chi \nabla_{{\perp}}^2 \left\langle T \right\rangle_\parallel = {-}\left\langle \left\{ \widetilde{\varphi}, \widetilde{T} \right\} \right\rangle_\parallel, \end{gather}

where the influence of small-scale fields appears only in the temperature equation (4.27). The system of (4.13)–(4.15) and (4.26)–(4.27) respects the conservation of the two conserved quantities described in § 2.2; this is shown in Appendix F.

The above reasoning does not apply to the ZFs. Indeed, the equation for $\overline {\varphi }$ is

(4.28)\begin{equation} \partial_t \overline{\varphi} - \overline{\partial_x \left\langle \varphi \right\rangle_\parallel \partial_y \left\langle \varphi + T \right\rangle_\parallel} - \partial_x^2(a\overline{\varphi} - b\overline{T}) = \overline{\partial_x \widetilde{\varphi} \partial_y (\widetilde{\varphi} + \widetilde{T})} = \overline{\widetilde{\Pi}_{xx}}. \end{equation}

This shows that the small-scale stress $\widetilde {\Pi }$ influences the zonal electrostatic potential $\overline {\varphi }$ more strongly (by a factor of $k_{\perp }^{-2}$) than it does the nonzonal $\left \langle \varphi ' \right \rangle _\parallel$. This is a consequence of the electron response (2.1) and the asymptotically smaller ‘inertia’ (i.e., the factor in front of the time derivative) $\propto k_{\perp }^2$ of the ZFs compared with the ‘inertia’ $\propto (1+k_{\perp }^2)$ of the nonzonal $\varphi '$. Thus, the right-hand side of (4.28) cannot be ignored. In fact, as we showed in § 4.1, the addition of 3-D effects, and hence of parasitic modes, has a profound impact on the stability of the Dimits-state ZFs, viz., the momentum flux $\widetilde {\Pi }_{xx}$ extends the Dimits state to higher temperature gradients than the 2-D system allows. Let us show why this is the case.

4.2.5 Turbulent stress due to parasitic sITG modes

In Ivanov et al. (Reference Ivanov, Schekochihin, Dorland, Field and Parra2020), we obtained a prediction for the critical gradient $\kappa _T^{c, \textrm {2D}} (\chi )$ above which a Dimits state with strong ZFs could not be sustained. This prediction was based on considerations of the ratio ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$ for the linear modes with largest growth rate. As explained in § 4.1.1, this ratio determines the balance of Reynolds and diamagnetic stresses for an individual Fourier mode: if ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} > -1$, then the Reynolds stress is larger and the mode favours a Dimits state, otherwise its diamagnetic stress is larger and the mode helps suppress the coherent ZFs needed for the Dimits state. In two dimensions, this ratio is sensitive to both the temperature gradient $\kappa _T$ and the collisionality $\chi$, and thus an appropriate balance between these two parameters is required in order to have ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} > -1$ for the dominant modes and thus to keep the system in the Dimits state. In particular, for $\kappa _T \gg 1$, the Dimits threshold is given by $\kappa _T / \chi = \textrm {const}$.

Let us adopt a similar approach for the fastest-growing small-scale sITG modes located at $k_\parallel \approx \boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k} k_{\perp }$. Equation (3.26) tells us that these modes satisfy ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} = -1 + O\left (k_{\perp }^{-2}\right )$ for $k_{\perp } \gg 1$. Therefore, to lowest order, the sITG modes are Dimits-marginal, i.e., their Reynolds and diamagnetic stresses balance out. This means that the lowest-order collisionless calculations of § 3.3.2 are insufficient for our needs. While we can extend these calculations to $O\left (k_{\perp }^{-2}\right )$, ignoring $\chi$ is, in fact, an unacceptable oversimplification. As we are about to see, nonzero $\chi$ provides a $O\left (k_{\perp }^{-1}\right )$ correction to ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$ and hence renders any collisionless higher-order corrections irrelevant.

As discussed in the beginning of § 4, the dominant collisionless sITG modes are found at $k_{\perp }^3 \sim \kappa _T / \chi$. It turns out that at those scales, it is collisional effects that determine ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$. The details of the relevant calculation can be found in Appendix D. We find that for $\chi$ ordered as $\kappa _T \sim \chi k_{\perp }^3$, the most unstable small-scale sITG modes are still located at $k_\parallel = \boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k}k_{\perp }$ and satisfy

(4.29)\begin{equation} \frac{T_{\boldsymbol{k}}}{\varphi_{\boldsymbol{k}}} ={-}1 - \sqrt{-\hat{\gamma}_{\boldsymbol{k}}^2 + \frac{{i}(a + b - 1)\chi k_{{\perp}}^2}{2\boldsymbol{\kappa}_T \boldsymbol{\cdot} \boldsymbol{k}}} + O\left(k_{{\perp}}^{{-}2}\right), \end{equation}

where we take the branch of the square root with a positive imaginary part. Since $\hat {\gamma }_{\boldsymbol {k}}^2 > 0$ and $\boldsymbol {\kappa }_T \boldsymbol{\cdot} \boldsymbol {k} > 0$ (as stipulated after (4.17)), we find that the sign of the real part of the square root in (4.29) is set by the sign of $a + b - 1$. Plugging in the numerical values $a = 9/40$ and $b = 67/160$, we find $a + b - 1 < 0$, hence the square root has a negative real part and ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} > -1$. Thus, collisionality always pushes the otherwise Dimits-marginal small-scale sITG modes to side with the Reynolds stress and reinforce the ZFs.Footnote ⁶ This is evident in figure 15.

Figure 15. (a) Linear growth rate and (b) the ratio ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$ of the most unstable ($k_x=0$) modes versus $k_\parallel$ and $k_y$ for $\kappa _T = 1$ and $\chi = 0.1$. The green dashed line is ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} = -1$. The black dashed line is the location of the largest collisionless growth rate $k_\parallel = \kappa _T k_y^2$. While the green and black lines would coincide to $O(k_{\perp }^{-2})$ for the collisionless modes, we see that the addition of collisions shifts the linearly unstable modes towards the Dimits-favourable ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}} > -1$ ratio.

The sensitivity of ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$ to the numerical factors $a$ and $b$ allows us to carry out a simple test of the above theory. We pick a simulation that is in the Dimits state in three dimensions, but above the 2-D Dimits threshold, i.e., has $\kappa _T > \kappa _T^{c, \textrm {2D}}$. We restart this simulation, but set $a = 1$ for all nonzonal modes. Linearly, this increases nonzonal viscosity and reduces growth rates, without affecting zonal physics. Naïvely, one might expect that with an increased damping of the turbulence, the Dimits state should become ‘stronger’. However, such reasoning does not take into account the structure of the 3-D modes and the change in the balance of Reynolds and diamagnetic stresses stemming from the change of the sign of $a + b - 1$. Indeed, in this numerical experiment, we discover that the Dimits regime is destroyed and strong turbulence sets in, just as the analysis above predicts. This is clear evidence that the most consequential role of collisionality for the Dimits regime of (2.11)–(2.13) is not to dissipate turbulent energy, but rather to regulate the turbulent stress via the ratio ${\textrm {Re}{\left (T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}}\right )}}$. This also suggests, for future analysis of the Dimits transition in different models of ITG turbulence, that the Dimits threshold may prove to be sensitive to the details of dissipation effects on the unstable modes, especially if, in the absence of collisions, these modes are Dimits-marginal, i.e., if they satisfy ${\textrm {Re}{(T_{\boldsymbol {k}}/\varphi _{\boldsymbol {k}})}} \approx -1$.

Let us also note that in the simple case of sITG modes in slab geometry, a more general calculation that includes kinetic effects is possible. In Appendices G and H, we derive the kinetic sITG dispersion relation and the kinetic equivalent of (4.1). Then, applying the ideas developed in Ivanov et al. (Reference Ivanov, Schekochihin, Dorland, Field and Parra2020) and in this work, we show that ZF-driving small-scale sITG modes are not limited to the cold-ion limit and could play a role outside of the realm of simple fluid approximations. This, of course, can be conclusively confirmed only by appropriate GK simulations.

4.3.1 Effect of parallel system size on the Dimits state

Figure 18(a) shows a typical example of the dependence of the saturated turbulent heat flux $Q$ on the parallel size of the box $L_{\parallel }$ for parameters $\kappa _T$ and $\chi$ that lie beyond the 2-D Dimits regime. For such parameters, the $L_{\parallel } = 0$ system does not reach finite-amplitude saturation. For $L_{\parallel }$ large enough, $Q$ is independent of $L_{\parallel }$, just as it was for parameters that were within the 2-D Dimits threshold (see § 4.1). As $L_{\parallel }$ is decreased, the ZFs break up and the system enters a strongly turbulent state. In figure 18(a), this happens for $L_{\parallel } < 1$. As $L_{\parallel }$ approaches zero, $Q$ starts to increase rapidly, signifying the approach to the 2-D state, where a blow up occurs.

Figure 18. Dependence of the saturated turbulent heat flux $Q$ on (a) the parallel size of the box $L_{\parallel }$ and (b) the largest parallel Fourier mode $k_{\parallel, \textrm {max}}$ that is included in the simulation.

Therefore, for each pair of values of $\kappa _T$ and $\chi$, there exists a critical $L_\parallel ^{c}$ such that the system is in the Dimits state for $L_{\parallel } > L_\parallel ^{c}$ and in the strongly turbulent regime for $L_{\parallel } < L_\parallel ^{c}$. It is clear that $L_\parallel ^{c} = 0$ if $\kappa _T <\kappa _T^{c, \textrm {2D}}$, i.e., if $\kappa _T$ and $\chi$ are such that the 2-D system is able to reach saturation. The dependence of $L_\parallel ^{c}$ on $\kappa _T$ and $\chi$ for $\kappa _T > \kappa _T^{c, \textrm {2D}}$ is not known at this point, due to the numerical cost of resolving simultaneously both the large $k_\parallel$ of the small-scale modes (see § 4.3.2) and the box-sized $k_\parallel \sim L_{\parallel }^{-1}$.

4.3.2 Effect of parallel resolution on the Dimits state

The scale separation between the large-scale cITG modes and the small-scale sITG modes increases the numerical cost of solving (2.11)–(2.13). When the parallel resolution, i.e., the largest $k_\parallel$ in the simulation, is too small, the Dimits state is destroyed numerically and the system is pushed into a strong-turbulence regime for parameters for which a Dimits state would have existed if given sufficient parallel resolution. This is shown in figure 18(b). Empirically, we have found that a good rule of thumb is ‘not to chop the leaves’ of the instability, i.e., to make sure that the wavenumbers that lie within the unstable ‘leaves’ at $k_\parallel \sim \kappa _T k_{\perp }^2$ (see figure 2) are fully included in the simulation.Footnote ⁷ This, however, rapidly increases the numerical cost of the simulations. Recall that according to (3.33), the collisionless sITG instability satisfies $k_{\perp } \sim (\kappa _T/\chi )^{1/3}$. Therefore, for a fixed $\chi$, the dimensional $k_\parallel$ of the unstable modes is given by

(4.30)\begin{equation} k_\parallel L_B \sim \kappa_T k_{{\perp}}^2 \sim \frac{\kappa_T^{5/3}}{\chi^{2/3}}. \end{equation}

The number of Fourier modes required to resolve a simulation properly then scales as $\kappa _T^{5/3}$, in addition to scaling linearly with $L_{\parallel }$. This quickly renders numerical efforts futile, even for a fluid code.

Of course, the infinitely extending ‘leaves’ of the instability in our 3-D model will, in reality, be ‘chopped off’ by phenomena that have been ordered out of our equations by (2.3). For example, (2.3) orders out the parallel thermal diffusion (Braginskii Reference Braginskii1965), but we can nonetheless estimate the dimensional $k_\parallel$ at which this effect will become important. This is the parallel scale at which the collisional heat conduction rate $v_{\textrm {th}i}^2k_\parallel ^2 / \nu _i$ becomes comparable to $\partial _t \sim c_s / L_B$ in (2.11)–(2.13), which happens at

(4.31)\begin{equation} k_\parallel L_B \sim \sqrt{\frac{L_B}{\sqrt{\tau}\lambda_\text{mfp}}}, \end{equation}

where $\lambda _\textrm {mfp} = v_{\textrm {th}i} / \nu _i$ is the mean-free path. In our ordering (2.3), $\lambda _\textrm {mfp} / L_B \sim \tau ^{3/2}$, so we find that the collisional heat conduction comes into play at $k_\parallel L_B \sim 1/\tau$. Formally, this is outside of the regime $k_\parallel L_B \sim 1$ assumed in (2.11)–(2.13), but physically, we conclude that the Dimits regime could be broken if the collisional cutoff (4.30) is superseded by the Braginskii scale (4.31), i.e., if

(4.32)\begin{equation} \frac{L_B}{L_T} \gtrsim \left(\frac{L_B}{\lambda_\text{mfp}}\right)^{7/10}\tau^{{-}11/20}, \end{equation}

where we used $\kappa _T \sim \tau L_B/L_T$ and $\chi \sim L_B\tau ^{3/2}/\lambda _\textrm {mfp}$. In a real fusion device, this condition will not be very difficult to reach, but, in fact, the more relevant mechanism for limiting the parallel wavenumber of the sITG instability is parallel streaming rather than collisional heat conduction. In Appendix G, we show that this too imposes a limit on the parallel wavenumber that is $O(\tau ^{-1})$ too large to be included in our ordering of $k_\parallel$. Namely, the sITG cutoff is given by

(4.33)\begin{equation} k_\parallel^{({c})} L_B = \frac{L_B}{2\sqrt{\rm \pi}(1+\tau)L_T}, \end{equation}

which supersedes the collisional cutoff (4.30) if

(4.34)\begin{equation} \frac{L_T}{\lambda_\text{mfp}} \lesssim \tau (1+\tau)^{3/2}. \end{equation}

Again, such a regime is entirely plausible for a real fusion device.

We conclude that in a more realistic physical regime than the one assumed in the derivation of our model equations (2.11)–(2.13), the behaviour (or even existence) of parasitic sITG modes may be influenced by parallel thermal diffusion or parallel streaming in a way that breaks the Dimits regime at large enough temperature gradients.