Hostname: page-component-84b7d79bbc-l82ql Total loading time: 0 Render date: 2024-07-29T23:21:28.307Z Has data issue: false hasContentIssue false

Local gyrokinetic collisional theory of the ion-temperature gradient mode

Published online by Cambridge University Press:  13 May 2022

B.J. Frei*
Affiliation:
École Polytechnique Fédérale de Lausanne (EPFL), Swiss Plasma Center (SPC), CH-1015 Lausanne, Switzerland
A.C.D. Hoffmann
Affiliation:
École Polytechnique Fédérale de Lausanne (EPFL), Swiss Plasma Center (SPC), CH-1015 Lausanne, Switzerland
P. Ricci
Affiliation:
École Polytechnique Fédérale de Lausanne (EPFL), Swiss Plasma Center (SPC), CH-1015 Lausanne, Switzerland
*
Email address for correspondence: baptiste.frei@epfl.ch
Rights & Permissions [Opens in a new window]

Abstract

We present a study of the linear properties of the ion-temperature gradient (ITG) modes with collisions modelled for the first time by the linearized gyrokinetic (GK) Coulomb collision operator (Frei et al., J. Plasma Phys., vol. 87, issue 5, 2021, 905870501) in the local limit. The study is based on a Hermite–Laguerre polynomial expansion of the perturbed ion distribution function applied to the linearized GK Boltzmann equation, yielding a hierarchy of coupled equations for the expansion coefficients, referred to as gyromoments. We explore the collisionless and high-collisional limits of the gyromoment hierarchy analytically. Parameter scans revealing the dependence of the ITG growth rate on the collisionality modelled using the GK Coulomb operator are reported, showing strong damping at small scales as the collisionality increases and, therefore, the need for a steeper gradient for the ITG onset at high collisionality to overcome the finite Larmor radius (FLR) collisional stabilization. The predictions on the ITG growth rate by the GK Coulomb operator are compared with other collision operator models, such as the Sugama, the Dougherty, as well as the momentum-conserving pitch-angle scattering and the Hirshman–Sigmar–Clarke collision operators derived for the first time in terms of gyromoments. The importance of FLR terms in the collision operators is pointed out by the appearance of a short wavelength ITG branch when collisional FLR terms are neglected, this branch being completely suppressed by FLR collisional effects. Energy diffusion is shown to be important at high collisionality and at small scale lengths. Among the GK collision operators considered in this work, the GK Sugama collision operator yields, in general, the smallest deviation compared with the GK Coulomb collision operator, while the largest deviations are found with the GK Dougherty operator. Convergence studies of the gyromoment method are reported and show that the drifts associated with the gradient and curvature of the magnetic field increase the required number of gyromoments at low collisionality. Nevertheless, the low number of gyromoments necessary for convergence at high collisionality constitutes an attractive numerical and analytical feature of the gyromoment approach to study the plasma dynamics in the boundary of fusion devices.

Type
Research Article
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
Copyright © The Author(s), 2022. Published by Cambridge University Press

1. Introduction

Belonging to the class of instabilities that develop on the ion gyroscale, the ion-temperature gradient (ITG) modes, driven unstable by the parallel plasma compression (slab ITG (sITG)) or by the presence of magnetic drifts (toroidal ITG), are widely recognized as the main candidate to explain the experimental observations of anomalous ion heat turbulent transport in the tokamak core (Garbet et al. Reference Garbet, Mantica, Angioni, Asp, Baranov, Bourdelle, Budny, Crisanti, Cordey and Garzotti2004). Because of its important role and since the plasma collision frequency is considerably smaller than the typical ITG mode frequency in the core region such that collisions can be neglected, the gyrokinetic (GK) collisionless theory of the ITG modes is well documented in the literature (see, e.g. Hahm & Tang Reference Hahm and Tang1989; Romanelli Reference Romanelli1989; Romanelli & Briguglio Reference Romanelli and Briguglio1990; Chen, Briguglio & Romanelli Reference Chen, Briguglio and Romanelli1991). On the other hand and because of the temperature drop, the effects of collisions on the ITG modes are expected to be important in the tokamak boundary, the region that encompasses the edge and the scrape-off layer (SOL). In fact, and more generally, collisions are known to modify the linear and nonlinear properties of microinstabilities. For instance, collisions contribute to smearing out small scale structures in velocity space (Watanabe & Sugama Reference Watanabe and Sugama2004; Ajay, Brunner & Ball Reference Ajay, Brunner and Ball2021), to affect their linear properties (Manas et al. Reference Manas, Camenen, Benkadda, Hornsby and Peeters2015; Belli & Candy Reference Belli and Candy2017), to suppress short-wavelength structures (Barnes et al. Reference Barnes, Abel, Dorland, Ernst, Hammett, Ricci, Rogers, Schekochihin and Tatsuno2009) and can, ultimately, affect their saturation mechanism via the damping of the zonal flow (Pan, Ernst & Crandall Reference Pan, Ernst and Crandall2020; Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). Since the boundary region sets the heat exhaust and the thermal load on the vessel walls and determines the overall confinement capabilities of the device, it is of primary importance to assess the role of the ITG modes in this key region that is characterized by a wide range of plasma collisionality. For instance, the collisionality dependence of the ITG mode is important in the prediction of, e.g. the transition from resistive ballooning to ITG driven turbulence regime in the boundary region (Zeiler et al. Reference Zeiler, Biskamp, Drake and Rogers1998; Dorf & Dorr Reference Dorf and Dorr2022).

In the past, the effects of finite collisionality on the ITG linear properties have been introduced in kinetic models using simplified collision operators such as, e.g. an energy-conserving Krook and a pitch-angle operator model (Romanelli & Briguglio Reference Romanelli and Briguglio1990; Romanelli, Chen & Briguglio Reference Romanelli, Chen and Briguglio1991; Pusztai et al. Reference Pusztai, Fülöp, Candy and Hastie2009). More generally, the effect of collisions on microinstabilities (including the trapped electron and ITG modes) have been documented in previous studies using GK codes (Catto & Ernst Reference Catto and Ernst2009; Belli & Candy Reference Belli and Candy2017; Pan & Ernst Reference Pan and Ernst2019; Crandall et al. Reference Crandall, Jarema, Doerk, Pan, Merlo, Görler, Navarro, Told, Maurer and Jenko2020; Pan et al. Reference Pan, Ernst and Crandall2020) showing discrepancies produced by the choice of collision operator models when compared with the linearized Fokker–Planck collision operator (Rosenbluth, MacDonald & Judd Reference Rosenbluth, MacDonald and Judd1957; Hazeltine & Meiss Reference Hazeltine and Meiss2003), that we refer to as the Coulomb collision operator in the present study. However, studies based on the GK Coulomb collision operator remain very limited (Pan et al. Reference Pan, Ernst and Crandall2020), due to the very recent availability of such an operator (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). Additionally, we lack a comparison between this operator and simplified collision operator models for the ITG case.

The limited number of collisional GK investigations stems from the frequent use of drift-reduced Braginskii fluid models (Zeiler, Drake & Rogers Reference Zeiler, Drake and Rogers1997) for the simulation of the boundary region (Dudson et al. Reference Dudson, Umansky, Xu, Snyder and Wilson2009; Halpern et al. Reference Halpern, Ricci, Jolliet, Loizu, Morales, Mosetto, Musil, Riva, Tran and Wersal2016; Tamain et al. Reference Tamain, Bufferand, Ciraolo, Colin, Galassi, Ghendrih, Schwander and Serre2016; Paruta et al. Reference Paruta, Ricci, Riva, Wersal, Beadle and Frei2018; Zhu, Francisquez & Rogers Reference Zhu, Francisquez and Rogers2018; Stegmeir et al. Reference Stegmeir, Ross, Body, Francisquez, Zholobenko, Coster, Maj, Manz, Jenko and Rogers2019; Giacomin, Stenger & Ricci Reference Giacomin, Stenger and Ricci2020). These fluid models are based on the assumptions that $k_\parallel \lambda _{{\rm mfp}} \ll 1$, i.e. that the particle mean free path $\lambda _{{\rm mpf}}$ is much shorter than the typical parallel length scales, described by the parallel wavenumber $k_\parallel$, and that $k_\perp \rho _i \ll 1$, i.e. that the modes develop a perpendicular wavenumber $k_\perp$ associated with a scale length much larger than the ion gyroradius $\rho _i$. Therefore, while the drift-reduced Braginskii models retain finite ITG effects (Hallatschek & Zeiler Reference Hallatschek and Zeiler2000; Mosetto et al. Reference Mosetto, Halpern, Jolliet, Loizu and Ricci2015), their applicability to the study of ITG modes is questionable since, for instance, these modes develop on a $k_\perp \rho _i \sim 1$ scale. In addition, while the $k_\parallel \lambda _{{\rm mfp}} \ll 1$ assumption is well established at the separatrix and in the SOL (with temperature $T \lesssim 50$ eV and particle density $n \sim 10^{18}$ m$^{-3}$ yielding typically $k_\parallel \lambda _{{\rm mfp}} \lesssim 10^{-2}$ in, e.g. JET, ITER and smaller tokamaks such as TCV (Frassinetti et al. Reference Frassinetti, Dunne, Sheikh, Saarelma, Roach, Stefanikova, Maggi, Horvath, Pamela and De La Luna2019)), the mean free path is larger inside the separatrix, where the temperature increases. Hence, for typical H-mode pedestal values ($T \sim 10^{2}\text{--}10^{3}$ eV and $n \sim 10^{20}$ m$^{-3}$ (Leyland et al. Reference Leyland, Beurskens, Frassinetti, Giroud, Saarelma, Snyder, Flanagan, Jachmich, Kempenaars and Lomas2014; Schneider et al. Reference Schneider, Angioni, Frassinetti, Horvath, Maslov, Auriemma, Cavedon, Challis, Delabie and Dunne2021)), an estimate of $k_\parallel \lambda _{{\rm mfp}} \sim 1$ can be obtained. Thus, one expects a level of collisionality sufficient to justify the use of the drift-reduced Braginskii fluid models only in the SOL and, possibly, in the outermost region of the edge in typical tokamak conditions.

Gyrofluid models have been developed by considering a finite number of velocity space moments of the GK equation to extend the validity of Braginskii-like fluid models (Hammett, Dorland & Perkins Reference Hammett, Dorland and Perkins1992; Dorland & Hammett Reference Dorland and Hammett1993; Beer & Hammett Reference Beer and Hammett1996; Snyder & Hammett Reference Snyder and Hammett2001). Collisionless kinetic effects are introduced in these models by designing ad hoc closures for the high-order moments that mimic the linear kinetic response in the absence of collisions. While usually neglected, collisional effects are introduced in, for example, the gyrofluid model proposed by Snyder & Hammett (Reference Snyder and Hammett2001) by considering a particle, momentum and energy-conserving BGK-like operator for ion collisions and a pitch-angle scattering operator for the collisions of the electrons. However, BGK-like operators and the absence of energy diffusion in the pitch-angle scattering operator can potentially lead to significant deviations compared with the Coulomb collision operator.

In this work, we aim to present a study of the linear properties of the ITG mode within a GK framework able to capture accurately the collisionless kinetic physics and, at the same time, the collisional effects described through the GK Coulomb collision operator for the first time. This is done by leveraging previous work that led to a formulation of the GK model for the boundary region based on a Hermite–Laguerre expansion of the gyrocentre distribution function (Frei, Jorge & Ricci Reference Frei, Jorge and Ricci2020), and the development of a GK Coulomb collision operator derived using the same technique expansion, as presented in Jorge, Frei & Ricci (Reference Jorge, Frei and Ricci2019a). The present study represents the first physical application of this model, which has the potential of providing an efficient description of plasma turbulence in the boundary region. In fact, the Hermite–Laguerre expansion technique shows its ability to treat the ITG mode, numerically and analytically, for a wide range of plasma collisionality and with sophisticated collision operators (Jorge, Ricci & Loureiro Reference Jorge, Ricci and Loureiro2018; Jorge et al. Reference Jorge, Ricci, Brunner, Gamba, Konovets, Loureiro, Perrone and Teixeira2019b). More precisely, the ITG model we use in this work is based on the linearized (also known as $\delta f$) and local limits of the gyromoment formulation of the GK model coupled with a linearized GK Coulomb collision operator, recently derived and numerically implemented in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). The Hermite–Laguerre expansion allows us to express the self-consistent GK model in terms of moments of the gyrocentre distribution function that we refer to as gyromoments.

By using the gyromoment approach, we first retrieve the collisionless and the high-collisional limits of the ITG mode. In particular, the ability of the gyromoment method to be reduced to simpler fluid models at high collisionality is analytically demonstrated and numerically tested. Then, using the linearized GK Coulomb collision operator, wide parameter scans of the ITG growth rate are performed. It is demonstrated that a steeper temperature gradient is required for the ITG onset at high collisionality (although higher than typical SOL experimental values) to overcome the finite Lramor radius (FLR) collisional stabilization effects. We also compare, for the first time, the GK Coulomb collision operator with the drift-kinetic (DK) Coulomb collision operator, a GK/DK momentum-conserving pitch-angle scattering operator (Helander & Sigmar Reference Helander and Sigmar2002), the zeroth-order DK Hirshman–Sigmar–Clarke (HSC) collision operator (Hirshman & Sigmar Reference Hirshman and Sigmar1976a) and the DK/GK Dougherty collision operator (Dougherty Reference Dougherty1964), in addition to the DK/GK Sugama collision operator implemented in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). Among the considered collision operator models, the GK Sugama operator yields the smallest deviation with respect to the GK Coulomb collision, while the largest deviations are found when the GK Dougherty collision operator is considered. We show the importance of energy diffusion and FLR terms in the collisional damping of ITG at small scales. Also, the importance of retaining FLR terms in the collision operator models is demonstrated by revealing that a short wavelength ITG (SWITG) mode at high collisionality and steep ITGs can be destabilized when the DK limit of the collision operators is used. We remark that the assessment of the potential deviations between collision operator models is of primary importance, since simplified and DK collision operators are often used in turbulent GK codes (Bernard et al. Reference Bernard, Shi, Gentle, Hakim, Hammett, Stoltzfus-Dueck and Taylor2019; Francisquez et al. Reference Francisquez, Juno, Hakim, Hammett and Ernst2021; Ulbl, Michels & Jenko Reference Ulbl, Michels and Jenko2021). Finally, we investigate the convergence properties of the gyromoment approach as a function of relevant physical parameters in different collisionality regimes, analysing the Hermite–Laguerre spectrum of the ITG distribution function and showing that magnetic drifts broaden the gyromoment spectrum significantly in the absence of collisions. The capability of reducing the number of gyromoments with increasing collisionality is illustrated by the present study, a key feature of the present approach to simulate the boundary region.

The structure of the paper is the following. In § 2 we introduce the gyromoment hierarchy equation. Then, in § 3, we derive the collisionless limit of the ITG dispersion relation by using the gyromoment hierarchy equation, and § 4 focuses on the high-collisional limits. In § 5, we evaluate the ITG growth rate over a wide range of parameters with collisions modelled by the GK Coulomb collision operator. Section 6 reports on the comparison of collision operator models with the GK Coulomb collision operator. Finally, we analyse the Hermite–Laguerre spectrum and report convergence studies in §§ 7 and 8, respectively. The conclusions follow in § 9. The Hermite–Laguerre expansions of the collision operators used in this work are reported in appendices A–C, while analytical results are detailed in Appendix D.

2. Gyrokinetic model

Focusing on a shearless configuration and adopting a local GK approach where quantities do not vary along the magnetic field lines (Kadomtsev & Shafranov Reference Kadomtsev and Shafranov2012), we consider the linearized electrostatic GK Boltzmann equation expressed in gyrocentre phase-space coordinates $(\boldsymbol {R}, \mu, v_\parallel,\theta )$ and assume that the electrons are adiabatic. Here, $\boldsymbol {R} = \boldsymbol {r} - \boldsymbol \rho$ is the ion gyrocentre position (with $\boldsymbol {r}$ being the particle position and $\boldsymbol \rho (\boldsymbol {R}, \mu, \theta ) = \boldsymbol {b} \times \boldsymbol v/\varOmega$ the ion gyroradius, with $\varOmega = q B / m$ the ion gyrofrequency, $q = Z e$ the ion charge with $Z =1$ and $\boldsymbol {b} = \boldsymbol {B} / B$), $\mu = m v^{2}/[2 B(\boldsymbol {R})]$ is the magnetic moment, $v_\parallel = \boldsymbol b \boldsymbol {\cdot } \boldsymbol v$ is the component of the velocity parallel to the magnetic field and $\theta$ is the gyroangle. Since only the evolution of the ion distribution is considered, we drop the ion species label for simplicity. We denote the equilibrium gradient scale lengths of ion density $N$, ion temperature $T$ and magnetic field strength $B$, by $L_N = | \boldsymbol {\nabla }_\perp \ln N |^{-1}$, $L_{T} = | \boldsymbol {\nabla }_\perp \ln T |^{-1}$ and $L_B = | \boldsymbol {\nabla }_\perp \ln B |^{-1}$, respectively. We assume that the gyroaveraged ion gyrocentre distribution function, $F = F (\boldsymbol {R}, \mu, v_\parallel )$, is a perturbed Maxwellian, that is $F = F_{M} + {g}$, where ${g}$ is the perturbed part (${g}/ F_{M} \ll 1$) and $F_{M}$ a local Maxwellian distribution function defined by $F_{M} = N(\boldsymbol {R}) {\rm e}^{- s_\parallel ^{2} - x} / ({\rm \pi} ^{3/2} v_{T}(\boldsymbol {R})^{3})$ with $s_\parallel = v_\parallel / v_{T}(\boldsymbol {R})$, $x= \mu B(\boldsymbol {R}) / T(\boldsymbol {R})$ and $v_T(\boldsymbol {R})^{2} = 2 T(\boldsymbol {R})/m$ the ion thermal velocity. The linearized electrostatic GK Boltzmann equation describing the time evolution of a single Fourier harmonic ${g}(\boldsymbol k,t) = \int \,{\rm d} \boldsymbol {R} {g}(\boldsymbol {R} ,t) {\rm e}^{- {\rm i} \boldsymbol R \boldsymbol {\cdot } \boldsymbol k}$ is then given by

(2.1)\begin{gather} \frac{\partial}{\partial t } {g} + \left({\rm i} \omega_{B} + {\rm i} k_\parallel v_\parallel \right) {g} - {\rm i} \omega_{*}^{T} {\rm J}_0(b\sqrt{x}) \frac{e \phi}{T_e} F_{M} + {\rm i} \omega_{\boldsymbol{\nabla} B} {\rm J}_0(b\sqrt{x}) \frac{q \phi }{T} F_{M}\nonumber\\ + \frac{q}{T} {\rm i} k_\parallel v_\parallel {\rm J}_0(b\sqrt{x}) \phi F_{M} = \mathcal{C}, \end{gather}

where $e$ is the electron charge and $\tau = T / T_e$. In (2.1), we express the wavevector as $\boldsymbol k = k_\parallel \boldsymbol {b} + \boldsymbol k_\perp$, and introduce the velocity-dependent resonant frequencies $\omega _{\boldsymbol {\nabla } B} = v_{T}^{2} ( \boldsymbol b \times \boldsymbol {\nabla } B /B) \boldsymbol {\cdot } \boldsymbol k /(2 \varOmega )(x + 2 s_\parallel ^{2})$ and $\omega _{*}^{T} = \omega _*[1 + \eta ( x + s_\parallel ^{2} - 3/2)]$, with the ion diamagnetic frequency $\omega _* = T_e (\boldsymbol b \times \boldsymbol {\nabla } \ln N ) \boldsymbol {\cdot } \boldsymbol k /(eB) > 0$. We also define the normalized temperature gradient $\eta = L_N / L_{T}$, a measure of the relative strength between the density and temperature gradients. The FLR effects are taken into account through the zeroth-order Bessel function, ${\rm J}_0$, with argument $k_\perp v_\perp / \varOmega = b \sqrt {x}$ being $b = k_\perp v_{T} /\varOmega$. Collisional effects are introduced in (2.1) through the GK linearized ion collision operator defined by

(2.2)\begin{equation} \mathcal{C} = \left\langle {C} \right\rangle_{\boldsymbol{R}} = \frac{1}{2 {\rm \pi}}\int_0^{2 {\rm \pi}} \,{\rm d} \theta C, \end{equation}

where $C$ is the linearized ion collision operator acting on the ion perturbed particle distribution function, and $\left \langle {\dots } \right \rangle _{\boldsymbol {R}} =\int _0^{2 {\rm \pi}} \,{\rm d} \theta \dots /(2 {\rm \pi})$ denotes the gyroaverage operator evaluated at constant $\boldsymbol {R}$. As the effects of collisions between ions and electrons develop over a time scale that is longer by a factor $O(\sqrt {m_i/ m_e})$ than the ones due to ion–ion collisions, we neglect the former, i.e. we assume $C = C_{ii}$ in (2.2). The linearized GK Boltzmann equation is closed by the GK quasineutrality condition that determines the Fourier component of the electrostatic potential $\phi (\boldsymbol k, t) = \int \,{\rm d} \boldsymbol {r} \phi (\boldsymbol {r},t) {\rm e}^{{\rm i} \boldsymbol k \boldsymbol {\cdot } \boldsymbol {r}}$, that is

(2.3)\begin{equation} \left[ \frac{ {\rm e}^{2}}{T_e} + \frac{ q^{2} }{T} \left( 1 - \varGamma_0(a) \right) \right] \phi = q \frac{\left\langle {\delta n} \right\rangle}{N}, \end{equation}

where $\left \langle {\delta n} \right \rangle = 2 {\rm \pi}\int \,{\rm d} v_\parallel \,{\rm d} \mu B / m {\rm J}_0(b) {g}$ is the perturbed gyroaveraged ion gyrocentre density, $a= b^{2} / 2$ and $\varGamma _0(a) = {\rm I}_0(a) {\rm e}^{- a}$, with ${\rm I}_0(a)$ the modified Bessel function. The first term on the left-hand side of (2.3) is the adiabatic electron response assumed to be given by $n_e \simeq e \phi / T_e$, while the second term constitutes the ion polarization term. We remark that the contribution from the flux-surface-average of $\phi$ in the electron adiabatic response is neglected in (2.3) due to the linear and local limit adopted in the present study (Dorland & Hammett Reference Dorland and Hammett1993).

We now normalize (2.1) and (2.3). The physical quantities $t$, $v_\parallel$, $k_\perp$, $k_\parallel$, $\phi$ and $q$, are normalized, respectively, to $L_N/c_{s}$, $c_s$, $1/ \rho _s$, $1/ L_N$, $T_{e}/e$ and $e$, with $c_{s} = \sqrt {T_{e}/m}$ the ion sound speed and $\rho _s = c_s / \varOmega$ the ion sound Larmor radius. The ion–ion collision frequency, $\nu = 4 N q^{4} \ln \varLambda / [ 3 \sqrt {{\rm \pi} } \sqrt {m} T^{3/2}]$ that enters in the ion–ion collision operator $\mathcal {C}$, is normalized to $c_s / (N L_N)$ and the gyrocentre density to $N$. The same dimensionless units are used in the rest of the present paper. The normalized GK Boltzmann equation, (2.1), then becomes

(2.4)\begin{gather} \frac{\partial }{\partial t} {g} + {\rm i} ( \omega_{\boldsymbol{\nabla} B} + k_\parallel \sqrt{2 \tau} s_{{\parallel} }) {g} - {\rm i} \omega_{*}^{T} {\rm J}_0(b \sqrt{x}) \phi F_{M} + \frac{q}{\tau}{\rm i} \omega_{\boldsymbol{\nabla} B} {\rm J}_0(b \sqrt{x}) \phi F_{M} \nonumber\\ + \frac{\sqrt{2}}{\sqrt{\tau}} {\rm i} k_\parallel s_\parallel {\rm J}_0(b\sqrt{x}) \phi F_{M} = \mathcal{C}, \end{gather}

where $\omega _{\boldsymbol {\nabla } B} = \tau \omega _{B} ( x + 2 s_\parallel ^{2}) /q$ (with $\omega _{ B} = k_\perp R_B$, $\tau = T_i / T_e$ and $R_B = L_N / L_B$ the normalized magnetic gradient), $\omega _*^{T} = \omega _{*} [1 + \eta (x + s_\parallel ^{2} - 3/2)]$ (with $\omega _* = k_\perp$ the normalized ion diamagnetic frequency). The argument of the Bessel function, $b \sqrt {x}$, is normalized to $b = k_\perp \sqrt {2 \tau }$. The self-consistent electrostatic potential, $\phi$, is obtained by the normalized linearized GK quasineutrality,

(2.5)\begin{equation} \left[ 1 + \frac{q^{2}}{\tau}\left( 1 - \varGamma_0(a) \right) \right] \phi = q \left\langle {\delta n} \right\rangle. \end{equation}

Equations (2.4) and (2.5) constitute a closed system that evolves self-consistently ${g}$ and $\phi$, in the presence of the temperature and magnetic equilibrium gradients, once the collision operator is specified (see § 2.2).

2.1. gyromoment hierarchy equation

To approach the solution of the GK model, given by (2.4) and (2.5), we use a Hermite–Laguerre moment expansion of the ion distribution function. This allows us to reduce the dimensionality of the linearized GK equation, (2.1), by projecting it onto velocity-space Hermite–Laguerre basis polynomials. More precisely, we decompose the perturbed gyrocentre distribution function as (Jorge, Ricci & Loureiro Reference Jorge, Ricci and Loureiro2017; Mandell, Dorland & Landreman Reference Mandell, Dorland and Landreman2018; Jorge et al. Reference Jorge, Frei and Ricci2019a; Frei et al. Reference Frei, Jorge and Ricci2020)

(2.6)\begin{equation} {g} = \sum_{p = 0}^{\infty} \sum_{j = 0}^{\infty} N^{pj} \frac{H_p(s_{{\parallel}}) L_j(x)}{\sqrt{2^{p} p!}} F_{M}, \end{equation}

where we define the Hermite–Laguerre velocity moments of ${g}$, as follows:

(2.7)\begin{equation} N^{pj} = \frac{1}{N} \int \,{\rm d} \mu \,{\rm d} {v_{{\parallel}}} \,{\rm d} \theta\frac{B}{m} {g} \frac{H_p(s_\parallel) L_j(x)}{\sqrt{2^{p} p!}}, \end{equation}

with $N = \int \,{\rm d} \mu \,{\rm d} {v_{\parallel }} \,{\rm d} \theta B F_M / m$. The Hermite–Laguerre coefficients $N^{pj}$ in (2.6) are referred to as the gyromoments. In (2.6), we introduce the Hermite and Laguerre polynomials, $H_p$ and $L_j$, via their Rodrigues’ formulae (Gradshteyn & Ryzhik Reference Gradshteyn and Ryzhik2014), that is

(2.8)\begin{equation} H_p(y) = ({-}1)^{p} {\rm e}^{y^{2}} \frac{{\rm d}^{p}}{{{\rm d} y}^{p}} ( {\rm e}^{- y^{2}} ) \end{equation}

and

(2.9)\begin{equation} L_j(y) = \frac{{\rm e}^{y}}{j!} \frac{{\rm d}^{j} }{{{\rm d} y}^{j}} ({\rm e}^{- y } y^{j}) \end{equation}

and their orthogonality relations, that is

(2.10)\begin{equation} \int_{- \infty}^{\infty} {{\rm d} y} H_p(y) H_{p'}(y) {\rm e}^{- y^{2}} = 2^{p} p! \sqrt{\rm \pi} \delta_{p}^{p'} \end{equation}

and

(2.11)\begin{equation} \int_0^{\infty} {{\rm d} y} L_j(y) L_{j'}(y) {\rm e}^{{-}y} = \delta_j^{j'}, \end{equation}

respectively. We now project (2.4) and (2.5) onto the Hermite–Laguerre polynomial basis. For this purpose, we note that the Bessel function ${\rm J}_m$ appearing in both (2.4) and (2.5) and in the GK collision operators (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021) can be conveniently projected by expanding the Bessel functions, ${\rm J}_m$, onto associated Laguerre polynomials, $L^{m}_n(x) = (-1)^{m} \,{\rm d}^{m} L_{n + m}(x) / {{\rm d} x}^{m}$,

(2.12)\begin{equation} {\rm J}_m(b \sqrt{x}) = \left(\frac{b \sqrt{x}}{2}\right)^{m}\sum_{n=0}^{\infty} \frac{n!\mathcal{K}_{n}}{(n + m)!} L^{m}_n(x), \end{equation}

where we introduce the $n$th-order kernel function

(2.13)\begin{equation} \mathcal{K}_{n} = \frac{1}{n!}\left(\frac{b}{2}\right)^{2n } {\rm e}^{- b^{2} /4}, \end{equation}

with argument $b = \sqrt {2 \tau } k_\perp$ (Frei et al. Reference Frei, Jorge and Ricci2020). The projections of (2.4) and (2.5) onto the Hermite–Laguerre basis yields the linearized gyromoment hierarchy equation that states the time evolution of the gyromoments, $N^{pj}$, that is

(2.14)\begin{gather} \frac{\partial}{\partial t} N^{pj} + {\rm i} k_\parallel \sqrt{\tau} \sum_{p'} \mathcal{H}^{{\parallel} p}_{p'} N^{p'j} + {\rm i} \frac{\tau}{q}\omega_{ B } \sum_{p'j'} ( \mathcal{H}_{B p'j'}^{{\parallel} pj} + \mathcal{H}_{B p'j'}^{{\perp} pj} ) N^{p'j'} \nonumber\\ = (\mathcal{F}^{pj} - \mathcal{F}_{B}^{{\parallel} pj}) \phi + \mathcal{C}^{pj} , \end{gather}

where we introduce the phase-mixing operators, associated with the parallel and perpendicular drifts,

(2.15a)\begin{gather} \mathcal{H}^{{\parallel} p}_{p'} = \sqrt{p +1 } \delta_{p'}^{p+1}+ \sqrt{p} \delta_{p'}^{p-1}, \end{gather}
(2.15b)\begin{gather}\mathcal{H}_{B p'j'}^{{\parallel} pj} = \delta_{j'}^{j} (\sqrt{(p+1)(p+2)} \delta_{p'}^{p+2} + \sqrt{2p+1} \delta_{p'}^{p} + \sqrt{p(p-1)} \delta_{p'}^{p-2}), \end{gather}
(2.15c)\begin{gather}\mathcal{H}_{B p'j'}^{{\perp} pj} = \delta_{p'}^{p} ( (2j+1)\delta_{j'}^{j} - j \delta_{j'}^{j-1} - (j+1) \delta_{j'}^{j+1}), \end{gather}

and the operators that provide the instability drive and are proportional to the equilibrium gradients

(2.16a)\begin{gather} \mathcal{F}^{pj} = {\rm i} \omega_* \left[ \mathcal{K}_{j} \delta_p^{0} + \eta \left( \mathcal{K}_{j} \frac{\delta_p^{2}}{\sqrt{2}} + \delta_p^{0} \left(2 j\mathcal{K}_{j} - j \mathcal{K}_{j-1} - (j+1)\mathcal{K}_{j+1} \right) \right) \right] , \end{gather}
(2.16b)\begin{gather}\mathcal{F}_{B}^{{\parallel} pj} =\frac{\delta_p^{1}}{\sqrt{\tau} } {\rm i} k_\parallel \mathcal{K}_{j} + {\rm i} \omega_{ B} \delta_p^{0} \left[ 2 ( j +1) \mathcal{K}_{j} - j \mathcal{K}_{j-1} - (j+1)\mathcal{K}_{j+1} \right]+ {\rm i} \omega_B \mathcal{K}_{j} \sqrt{2} \delta_p^{2} . \end{gather}

In (2.14), the Hermite–Laguerre projection of the collision operator, $\mathcal {C}$, is defined by

(2.17)\begin{equation} \mathcal{C}^{pj} = \int \,{\rm d} \theta \,{\rm d} \mu \,{\rm d} v_\parallel \frac{B}{m} \frac{H_p(s_\parallel) L_j(x)}{ \sqrt{2^{p} p!}} \mathcal{C}. \end{equation}

In this work, different collision operator models for $\mathcal {C}$ are considered and detailed in § 2.2. The linearized gyromoment hierarchy equation is coupled with the self-consistent GK quasineutrality condition that, projected onto the Hermite–Laguerre basis, is given by

(2.18)\begin{equation} \left[ 1 + \frac{q^{2}}{\tau}\left( 1 - \sum_{n=0}^{\infty} \mathcal{K}_{n}^{2} \right)\right] \phi = q \sum_{n=0}^{\infty} \mathcal{K}_{n} N^{0n}. \end{equation}

Equations (2.14) and (2.18) constitute an infinite set of linear coupled fluid equations for the gyromoments $N^{pj}$, and depend only on the Fourier mode wavevector $\boldsymbol k$ and time $t$. In order to solve numerically (2.14) together with (2.18), we apply a simple closure by truncation, i.e. we solve (2.14) and (2.18) up to $(p,j) \leq (P,J)$ and set $N^{pj} = 0$ for all $(p,j) > (P,J)$ given $0 < P, J < \infty$. When applying the closure by truncation to the gyromoment hierarchy equation, the infinite sums appearing in the GK quasineutrality, (2.18), are consistently truncated up to $J$, that is

(2.19)\begin{equation} \left[ 1 + \frac{q^{2}}{\tau}\left( 1 - \sum_{n=0}^{J} \mathcal{K}_{n}^{2} \right)\right] \phi = q \sum_{n=0}^{J} \mathcal{K}_{n} N^{0n}. \end{equation}

Finally, we remark that a high-collisional closure of the gyromoment hierarchy is derived in § 4.1.

A gyromoment hierarchy equation, similar to the one given in (2.14), is obtained by Mandell et al. (Reference Mandell, Dorland and Landreman2018). In their formulation, the probabilist's Hermite polynomials, defined by $\hat {H}_p(y) = 2^{- p/2} H_p( y / \sqrt {2})$ and satisfying the orthogonality relation given in (2.10) but weighted by ${\rm e}^{-y^{2} / 2}$, are used, while collisional effects are modelled by the GK Dougherty collision operator model because of the simplicity of its expression on the Hermite–Laguerre basis (see Appendix C). Since the collisional effects are the main focus of the present study, we consider here more realistic collision operator models.

2.2. Gyrokinetic linearized collision operator models

While numerous collision operator models are available in the literature (see, e.g. Hirshman & Sigmar Reference Hirshman and Sigmar1976a; Abel et al. Reference Abel, Barnes, Cowley, Dorland and Schekochihin2008; Sugama, Watanabe & Nunami Reference Sugama, Watanabe and Nunami2009; Li & Ernst Reference Li and Ernst2011; Madsen Reference Madsen2013b; Sugama et al. Reference Sugama, Matsuoka, Satake, Nunami and Watanabe2019), the main focus of this work is a linearized GK Coulomb collision operator, and we leverage its formulation derived and implemented in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). This advanced linearized GK collision operator is obtained from the nonlinear full-F GK Coulomb operator derived in Jorge et al. (Reference Jorge, Frei and Ricci2019a). To highlight the importance of FLR terms in the Coulomb collision operator, the DK limit of the GK Coulomb collision operator is also considered. In addition, collision operators found in the literature and widely used in numerical codes are also considered for comparison. In particular, we consider the GK and DK Sugama collision operators (Sugama et al. Reference Sugama, Watanabe and Nunami2009), the GK and DK momentum-conserving pitch-angle scattering operators (Helander & Sigmar Reference Helander and Sigmar2002), the zeroth-order DK HSC collision operator (Hirshman & Sigmar Reference Hirshman and Sigmar1976a) and the GK and DK Dougherty collision operators (Dougherty Reference Dougherty1964). The gyromoment expansion of the linearized GK/DK Coulomb and the GK/DK Sugama collision operators were reported in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021) where the linearized DK Coulomb and GK/DK Sugama operators are successfully benchmarked against the GK continuum code GENE (Jenko et al. Reference Jenko, Dorland, Kotschenreuther and Rogers2000). The Hermite–Laguerre expansion of the pitch-angle, HSC and Dougherty collision operators are reported in appendices A–C, respectively.

Before turning to the investigation of the impact of collisions on the ITG mode, here we recall some important conservation properties of linearized like-species collision operators. These properties are useful to derive the high-collisional limit of the gyromoment hierarchy (see § 4.1) and to verify the numerical implementation of the collision operators. A model of a linearized collision operator for ion–ion collisions, $C$, is obtained by linearizing a nonlinear collision operator, that we denote by $C^{NL}$. Since collisions act at the particle position $\boldsymbol {r}$ (rather than at the gyrocentre position $\boldsymbol {R}$), the nonlinear collision operator is most often derived as acting on the full particle distribution function and is expressed in the particle phase space $(\boldsymbol {r}, \boldsymbol v)$. While the functional form of $C^{NL}$ depends on the choice of collision model, the linearized collision operator, $C$, can always be expressed as

(2.20)\begin{equation} C = C^{T} + C^{F}, \end{equation}

where we introduce the test component of the linearized collision operator, $C^{T} = C^{NL}(f,f_{M})$, and the field component, $C^{F} = C^{NL}(f_{M},f)$, with $f_M = f_M(\boldsymbol r, \boldsymbol v)$ the particle Maxwellian distribution function and $f = f (\boldsymbol {r},\boldsymbol v)$ a small perturbation. Conservation properties of the collision operator constrain the test and field components. In fact, while particle conservation is satisfied by both components, such that

(2.21)\begin{equation} \int \,{\rm d} \boldsymbol{v} C^{T} = \int \,{\rm d} \boldsymbol{v} C^{F} = 0, \end{equation}

the momentum and energy conservations require that

(2.22)\begin{gather} \int \,{\rm d} \boldsymbol{v} \boldsymbol v C^{T} ={-} \int \,{\rm d} \boldsymbol{v}\boldsymbol v C^{F} , \end{gather}
(2.23)\begin{gather}\int \,{\rm d} \boldsymbol{v} v^{2} C^{T} ={-} \int \,{\rm d} \boldsymbol{v} v^{2} C^{F} , \end{gather}

respectively. From the conservation properties of the linearized collision operator $C$, (2.21)(2.23), constraints on the coefficients $\mathcal {C}^{pj}$ can be derived in the zero gyroradius limit, when the gyrocentre and particle position correspond. In fact, using (2.17) combined with (2.21)(2.23) yields the relations

(2.24a–c)\begin{equation} \mathcal{C}^{00} = 0, \quad \mathcal{C}^{10}= 0, \quad \mathcal{C}^{20} = \sqrt{2} \mathcal{C}^{01}, \end{equation}

which express the conservation of gyrocentres, momentum and energy. The constraints in (2.24ac) are satisfied by all linearized collision operator models considered in the present work in the long wavelength limit and are fulfilled numerically to machine precision.

3. Collisionless limit

In this section, we derive the collisionless limit of the gyromoment hierarchy equation given in (2.14), showing that the ITG collisionless results (see, e.g. Hahm & Tang Reference Hahm and Tang1989; Romanelli Reference Romanelli1989; Romanelli & Briguglio Reference Romanelli and Briguglio1990; Chen et al. Reference Chen, Briguglio and Romanelli1991) can be retrieved as an asymptotic limit of the gyromoment hierarchy. We first derive the collisionless ITG dispersion relation in § 3.1 by solving directly the GK equation, (2.4). Then, we obtain closed semianalytical expressions of the collisionless linear gyromoment response in § 3.2. The expressions of the gyromoments we obtain are valid at arbitrary order in the normalized magnetic frequency, $\omega _{\boldsymbol {\nabla } B}/\omega$, and are written in terms of linear combinations of definite integrals of velocity-independent hypergeometric functions, which reduce to linear combinations of derivatives of plasma dispersion functions in the slab limit. Finally, in § 3.3 and thanks to the expressions of the gyromoments developed in § 3.2, we show the equivalence between the gyromoment and GK dispersion relations by considering the infinite gyromoment limit, focusing on the slab branch. Ultimately, the results of the present section (in particular the expressions of the collisionless gyromoment response) allow us to investigate the collisionless gyromoments spectrum and to progress in the analysis of the convergence properties of the gyromoment approach, discussed in §§ 7 and 8.

3.1. Collisionless GK dispersion relation

We first derive the ITG dispersion relation by solving directly the GK equation, (2.4). Fourier transforming (2.1) in time with $\mathcal {C} =0$, such that $\partial _t \to - {\rm i} \omega$, and introducing $\omega = \omega _r + {\rm i} \gamma$, with $\omega _r$ and $\gamma$ the real mode frequency and the growth rate, respectively, the normalized perturbed ion gyrocentre distribution function, ${g}/ \phi$, can be expressed as

(3.1)\begin{equation} \frac{{g}}{\phi}= \sum_{l=1}^{3} \hat{{g}}_{l} \end{equation}

with

(3.2a)\begin{gather} \hat{{g}}_{1 } ={-} \frac{q}{\tau} {\rm J}_0(b \sqrt{x}) F_{M}, \end{gather}
(3.2b)\begin{gather}\hat{{g}}_{2 } = \frac{q}{\tau} \frac{ \omega }{ \omega - \omega_{\boldsymbol{\nabla} B}- z_\parallel s_\parallel } {\rm J}_0(b \sqrt{x}) F_{M}, \end{gather}
(3.2c)\begin{gather}\hat{{g}}_{3} ={-}\frac{\omega_{T}^{*}}{\omega - \omega_{\boldsymbol{\nabla} B}- z_\parallel s_\parallel }{\rm J}_0(b \sqrt{x}) F_{M}, \end{gather}

where we introduce $z_\parallel = \sqrt {2 \tau } k_\parallel$. The $l=1$ contribution to the sum in (3.1) is the adiabatic response, while the $l=2$ and $l=3$ terms are associated with the Landau damping and equilibrium gradient drive. We remark the presence of the velocity-dependent magnetic frequency $\omega _B$ in the resonant denominator of the $l=2$ and $l=3$ contributions.

The collisionless GK dispersion relation is deduced from the linearized GK quasineutrality condition given in (2.5),

(3.3)\begin{equation} 1 + \frac{q^{2}}{\tau} \left( 1 - \varGamma_0(a) \right) - q \hat{n} =0, \end{equation}

where $\hat {n} = \sum _{l=1}^{3} \hat {n}_l$ is the density perturbation with $\hat {n}_l = \int \,{\rm d} \boldsymbol {v} {\rm J}_0 \hat {{g}}_l$ and $\hat {{g}}_l$ given in (3.2). By considering an unstable mode, i.e. $\gamma > 0$, the velocity integrals appearing in $\hat {n}_l$ can be performed explicitly by writing the resonant term, $1/(\omega - \omega _{\boldsymbol {\nabla } B} - z_\parallel s_\parallel )$, in integral form (Beer & Hammett Reference Beer and Hammett1996)

(3.4)\begin{equation} \frac{1}{ \omega - \omega_{\boldsymbol{\nabla} B} - z_\parallel s_\parallel} ={-} {\rm i} \int_0^{\infty} \,{\rm d} \sigma \exp({ {\rm i} \sigma (\omega - \omega_{\boldsymbol{\nabla} B} - z_\parallel s_\parallel)}). \end{equation}

We remark that (3.4) allows us to evaluate the velocity integrals, appearing in $\hat {n}_l$, exactly and to retain the effects of magnetic drifts in $\hat {N}_{l}^{pj}$ at arbitrary order in $\omega _B / \omega$. Using (3.4), one derives

(3.5a)\begin{align} \hat{n}_{1} & ={-} \frac{q}{\tau} \varGamma_0(a), \end{align}
(3.5b)\begin{align} \hat{n}_{2} & ={-} \frac{ {\rm i}q}{\tau} \omega \int_0^{\infty} \,{\rm d} \sigma {\rm e}^{{\rm i} \sigma \omega} {\rm I}_\perp(\sigma) {\rm I}_\parallel(\sigma), \end{align}
(3.5d)\begin{align} \hat{n}_{3} & = {\rm i} k_\perp \int_0^{\infty} \,{\rm d} \sigma {\rm e}^{{\rm i} \sigma \omega} \left[\vphantom{\frac{3}{2}} {\rm I}_\parallel (\sigma) {\rm I}_\perp (\sigma) \right. \nonumber\\ & \quad \left.+ \eta \left( {\rm I}_{{\parallel} 2} (\sigma) {\rm I}_\perp (\sigma) + {\rm I}_{{\parallel} } (\sigma) {\rm I}_{{\perp} 1} (\sigma) - \frac{3}{2} {\rm I}_\parallel (\sigma) {\rm I}_\perp (\sigma) \right) \right], \end{align}

where we introduce

(3.6a)\begin{align} {\rm I}_\perp(\sigma)& = \frac{1}{1 + {\rm i}\alpha \sigma} {\rm I}_0\left(\frac{a }{1 + {\rm i}\alpha \sigma} \right) \exp({-a /(1 + {\rm i} \alpha\sigma)}), \end{align}
(3.6b)\begin{align} {\rm I}_\parallel(\sigma) & = \frac{1}{\sqrt{1 + 2{\rm i} \alpha \sigma}} \exp({ - z_\parallel^{2} \sigma^{2} /[4 (1 + 2 {\rm i} \alpha \sigma)]}), \end{align}
\begin{align} {\rm I}_{{\perp} 1}(\sigma) & = \frac{\exp({- a / (1 + {\rm i} \alpha \sigma)})}{2 (1 + {\rm i} \alpha \sigma)^{3}} \nonumber\\ & \quad \times \left[ (2 (1 + {\rm i} \alpha \sigma) - 2 a) {\rm I}_0 \left( \frac{a}{ (1 + {\rm i} \alpha \sigma)} \right)+ 2 a {\rm I}_1 \left( \frac{a}{ (1 + {\rm i} \alpha \sigma)} \right) \right], \end{align}
(3.6c)\begin{align} {\rm I}_{{\parallel} 2}(\sigma) & = \frac{(2(1 + 2 {\rm i} \sigma \alpha ) - \sigma^{2} z_\parallel^{2}) }{4 (1 + 2 {\rm i} \sigma \alpha)^{5/2}} \exp({- z_\parallel^{2} \sigma^{2} / ( 4 (1 + 2 {\rm i} \sigma \alpha))}), \end{align}

with $\alpha = q \omega _B / \tau$. Equation (3.3), with the definitions in (3.5) and the velocity integrals in (3.6) constitute the collisionless GK dispersion relation for the unstable modes $\gamma > 0$. We note that the exponential factors ${\rm e}^{{\rm i} \sigma \omega }$, appearing in (3.5), ensure the convergence of the integrals at $\sigma \to \infty$ for $\gamma > 0$. On the other hand, the inclusion of stable modes ($\gamma < 0$) in (3.3) can be obtained by expressing the velocity integrals in $\hat {n}$ using generalized plasma dispersion functions (Gültekin & Gürcan Reference Gültekin and Gürcan2019) and efficient numerical algorithms can be used to solve the dispersion relation in the entire complex plane (Xie et al. Reference Xie, Li, Lu, Ou and Li2017; Gültekin & Gürcan Reference Gültekin and Gürcan2018). However, here we focus on the upper half of the complex plane since one of the advantages of the transformation (3.4) is that it yields one-dimensional integrals that can be easily evaluated numerically.

We finally note that the collisionless GK dispersion relation can be simplified in the sITG case. Neglecting $\omega _B$ in (3.2), the velocity integrals appearing in $\hat {n}$ and given in (3.5), can be expressed in terms of the plasma dispersion function. Then, the GK dispersion relation, (3.3), becomes

(3.7)\begin{align} & 1 + \frac{1}{\tau} + \frac{1}{\tau} \xi Z(\xi) \varGamma_0(a) - \frac{k_\perp}{k_\parallel} \frac{1}{\sqrt{2 \tau}} \left[ \varGamma_0(a) Z(\xi) \right. \nonumber\\ & \qquad \left. + \eta \left( \varGamma_0(a) \xi \left( 1 + \xi Z(\xi)\right) - \frac{\varGamma_0(a)}{2}Z(\xi) + a Z(\xi)\left(\varGamma_1(a) - \varGamma_0(a) \right) \right)\right] =0, \end{align}

where $\xi = \omega \sqrt {2 \tau }/ k_\parallel$ is the normalized mode phase velocity, $Z(\xi ) \!=\! \int _{-\infty }^{\infty } \,{\rm d} s_\parallel {\rm e}^{- s_\parallel ^{2}}/(s_\parallel \!-\! \xi ) /\sqrt {{\rm \pi} }$ is the plasma dispersion function and $\varGamma _1(a) = {\rm I}_1(a) {\rm e}^{- a}$.

3.2. Collisionless linear gyromoment response

We now obtain the collisionless expressions of the gyromoments by projecting (3.1) onto the Hermite–Laguerre basis. This yields

(3.8)\begin{equation} \frac{N^{pj}}{\phi} = \sum_{l=1}^{3} \hat{N}_{l}^{pj}, \end{equation}

where we introduce

(3.9)\begin{equation} \hat{N}_{l}^{pj} = 2 {\rm \pi}\int \,{\rm d} v_\parallel \,{\rm d} \mu \frac{B}{m} \frac{H_p(s_\parallel ) L_j(x) \hat{{g}}_l }{ \sqrt{2^{p} p !}}, \end{equation}

with $\hat {{g}}_l$ defined in (3.2). Performing the velocity integral in (3.9) using (3.4) for an unstable mode with $\gamma >0$ and identities involving hypergeometric functions and Hermite–Laguerre polynomials (described in Appendix D), we obtain

(3.10a)\begin{align} \hat{N}_1^{pj} & ={-} \frac{q}{\tau} \mathcal{K}_{j} \delta_p^{0}, \end{align}
(3.10b)\begin{align} \hat{N}_{2}^{pj}& ={-} \frac{{\rm i} q \omega}{\tau} \sum_{n=0}^{\infty} \frac{ \mathcal{K}_{n} }{\sqrt{2^{p} p!} }\mathcal{I}_{n}^{pj} , \end{align}
(3.10c)\begin{align} \hat{N}_{3}^{pj}& = \frac{{\rm i} k_\perp}{\sqrt{2^{p} p!}}\sum_{n=0}^{\infty}\mathcal{K}_{n} \left\{\vphantom{\frac{3}{2}} \mathcal{I}_n^{pj} + \eta \left[\vphantom{\frac{3}{2}} (2j+1) \mathcal{I}_n^{pj} - j \mathcal{I}_n^{pj-1}- (j+1) \mathcal{I}_n^{pj+1} \right.\right. \nonumber\\ & \quad \left. \left. + \frac{\mathcal{I}_{n}^{p+2j}}{4} + \left(p + \frac{1}{2} \right) \mathcal{I}_n^{pj} + p(p-1) \mathcal{I}_n^{p-2j}- \frac{3}{2}\mathcal{I}_n^{pj}\right] \right\}, \end{align}

where $\mathcal {I}_{n}^{pj}$ is given in (D8) as a finite sum of terms that involve velocity-independent definite integrals of hypergeometric functions with complex argument. The expressions in (3.10) are semianalytical because they still contain a real definite integrals that need to be computed numerically. We remark that generalized plasma dispersion relations can also be used to express the collisionless gyromoments $\hat {N}_{l}^{pj}$, given in (3.10), extending the formulae in Appendix D to the lower half of the complex plane of $\omega$. We use the collisionless expressions of the gyromoments in § 7 to investigate their magnitude in the presence of kinetic effects, such as Landau damping and magnetic gradient drift resonance.

In the case of the sITG branch, the velocity integrals containing the resonant term, $1 /(\omega - z_\parallel s_\parallel )$ (being $\omega _B = 0$), can be performed without using (3.4). In fact, in this case, the collisionless expression of the gyromoments becomes

(3.11)\begin{align} \frac{N^{pj}}{\phi}& ={-}\frac{q}{\tau} \mathcal{K}_{j} \delta_{p}^{0} + \frac{1}{\sqrt{2 \tau}}\frac{k_\perp}{k_\parallel} \frac{1 }{\sqrt{2^{p} p!}} \left( \mathcal{K}_{j} ({-}1)^{p} Z^{(p)} (\xi) \right. \nonumber\\ & \quad \left. + \eta \left[\!\mathcal{K}_{j} \left( \frac{({-}1)^{p+2}}{4} Z ^{(p+2)}(\xi) + (p + 2j) ({-}1)^{p} Z ^{(p)}(\xi) + p(p-1) ({-}1)^{p-2} Z ^{(p-2)}(\xi)\!\right) \right. \right.\nonumber\\ & \quad \left. \left. - \left[ j \mathcal{K}_{j-1} + (j+1) \mathcal{K}_{j+1} \right] ({-}1)^{p} Z ^{(p)}(\xi) \right) \right] - \frac{q}{\tau} \frac{ \xi \mathcal{K}_{j}}{\sqrt{2^{p} p!}} ({-}1)^{p} Z^{(p)}(\xi), \end{align}

where $Z^{(p)} (\xi )= {\rm d}^{p} Z(\xi )/ {\rm d} \xi ^{p}$ is the $p$th-order derivative of the plasma dispersion relation, which can be defined in terms of the Hermite polynomials, $H_p$, by $Z^{(p)} (\xi ) = (-1)^{p} \int _{-\infty }^{\infty } \,{\rm d} s_\parallel H_p(s_\parallel ) {\rm e}^{- s_\parallel ^{2}}/(s_\parallel - \xi )/\sqrt {{\rm \pi} }$.

As an aside, we also note that (3.8) can be related to the kinetic response of the monomial velocity-space moments, $M_{j,k} = v_T^{j + 2k}\int \,{\rm d} \boldsymbol {v} {g}s_\parallel ^{j} x^{k} / 2^{k}$, used in Beer & Hammett (Reference Beer and Hammett1996) to derive general closure expressions for toroidal gyrofluid equations, in the $k_\perp = 0$ limit, such as

(3.12)\begin{equation} N^{pj} = \sum_{p_1=0}^{\lfloor p/2 \rfloor} \sum_{j_1=0}^{j} \frac{({-}1)^{p_1} 2^{p-2p_1}}{p_1!(p-2p_1)!} \frac{({-}1)^{j_1} j !}{(j - j_1)!(j_1!)^{2}} \frac{2^{j_1}}{v_T^{2j_1 + p - 2p_1}} \frac{M_{p-2p_1,j_1}}{\sqrt{2^{p} p!}}. \end{equation}

Equation (3.12) provides a direct link between the gyromoment hierarchy equation described here and the moment equations derived in Beer & Hammett (Reference Beer and Hammett1996) and associated closure expressions.

3.3. gyromoment dispersion relation and infinite gyromoment limit

The ITG dispersion relation expressed in terms of the gyromoments $N^{pj}$ (with the closure up to $(p,j) \leq (P,J)$) is obtained from the truncated GK quasineutrality condition (2.19), yielding

(3.13)\begin{equation} 1 + \frac{q^{2}}{\tau}\left( 1 - \sum_{j=0}^{J} \mathcal{K}_{j}^{2} \right) - \sum_{j=0}^{J} \mathcal{K}_{j} \sum_{l=1}^{3 }\hat{N}^{0j}_l = 0, \end{equation}

where $\hat {N}_l^{0j}$ is given by (3.10) with $p=0$.

Focusing on the sITG branch, we now demonstrate that, in the $P \to \infty$ and $J \to \infty$ limits, the ITG gyromoment dispersion relation, (3.13), is equivalent to the GK dispersion relation given in (3.7). In order to solve the gyromoment dispersion relation (3.13), we use the expressions for the collisionless gyromoments, $N^{pj} / \phi$ given in (3.11) in the sITG case. Setting $p =0$ in (3.11) and plugging the resulting expression for $N^{p0} / \phi$ into (3.13) yields

(3.14)\begin{align} & 1 + \frac{1}{\tau} + \frac{1}{\tau} \xi Z(\xi) \hat{\varGamma_0}(a) - \frac{k_\perp}{k_\parallel} \frac{1}{\sqrt{2 \tau}} \left[\vphantom{\frac{\hat{\varGamma_0}(a,J)}{2}} \hat{\varGamma_0}(a,J) Z(\xi) + \eta\left(\vphantom{\frac{\hat{\varGamma_0}(a,J)}{2}}\hat{\varGamma_0}(a,J) \xi \left( 1 + \xi Z(\xi)\right) \right. \right. \nonumber\\ & \quad \left. \left. - \frac{\hat{\varGamma_0}(a,J)}{2}Z(\xi) +Z(\xi) \hat{\varPi}_1(a,J) \right)\right] =0, \end{align}

where we introduce the quantities

(3.15a)\begin{gather} \hat{\varGamma_0}(a,J)= \sum_{j=0}^{J} \mathcal{K}_{j}^{2}, \end{gather}
(3.15b)\begin{gather}\hat{\varPi}_1(a,J) = \sum_{j=0}^{J} \mathcal{K}_{j} \left(2j \mathcal{K}_{j} - j \mathcal{K}_{j-1} -(j+1)\mathcal{K}_{j+1} \right), \end{gather}

which are associated with FLR effects. We remark that the sITG dispersion relation, given in (3.14) obtained using (3.11) and (3.13), agrees with the one found in Mandell et al. (Reference Mandell, Dorland and Landreman2018). By comparing (3.14) with (3.7), one observes that they are equivalent in the infinite moment limits provided that $\hat {\varGamma _0}(a,\infty ) = \varGamma _0$ and $\hat {\varPi }_1(a,\infty ) = a (\varGamma _1 - \varGamma _0)$. This can be verified explicitly by considering the $J \to \infty$ limit in (3.15) and using the definition of the kernel function $\mathcal {K}_{j}$ in (2.13). In fact, we derive

(3.16)\begin{align} \hat{\varGamma_0}(a,\infty) & = \sum_{j=0}^{\infty} \frac{1}{(j!)^{2}} \left( \frac{b}{2}\right)^{4j} {\rm e}^{- b^{2} /2} \nonumber\\ & = {\rm e}^{- a} {\rm I}_0\left( a\right) = \varGamma_0(a) \end{align}

and

(3.17)\begin{align} \hat{\varPi}_1(a,\infty) & = 2 {\rm e}^{ - b^{2} /2}\sum_{j=0}^{\infty} \left[ \frac{1 }{j!(j+1)!} \left( \frac{b}{2} \right)^{4(j+1)} - \frac{1}{(j!)^{2}} \left( \frac{b}{2}\right)^{4j + 2}\right] \nonumber\\ & = a {\rm e}^{- a} \left( {\rm I}_1(a) - {\rm I}_0(a)\right) = a \left( \varGamma_1(a)- \varGamma_0(a)\right). \end{align}

Figure 1 displays the quantities, $\hat {\varGamma _0}(a,J)$ and $\hat {\varPi }_1(a,J)$ defined in (3.15), for different values of $J$, as well as their $J \to \infty$ limits. It can be observed that $\hat {\varGamma _0}(a,J)$ and $\hat {\varPi }_1(a,J)$ converge to $\varGamma _0$ and $a(\varGamma _1 - \varGamma _0)$ when $J\to \infty$. We notice that convergence is faster for low values of $a = b^{2}/2$, i.e. at long wavelength. This is in agreement with the fact that the truncation of the sums in (3.16) and (3.17) to include only $j \leq J$ terms yields an error $O( b^{4J})$ (Mandell et al. Reference Mandell, Dorland and Landreman2018).

Figure 1. Plots of $\hat {\varGamma _0}(a,J)$ (a) and $\hat {\varPi }_1(a,J)$ (b) as a function of $a = b^{2}/2$ for increasing values of $J$. The solid black lines represent their asymptotic limits, $\varGamma _0$ and $- a (\varGamma _1 - \varGamma _0)$, respectively.

The impact on the ITG mode due to approximating $\hat {\varGamma _0}(a,\infty )$ and $\hat {\varPi }_1(a,\infty )$ with $\hat {\varGamma _0}(a,J)$ and $\hat {\varPi }_1(a,J)$ can be illustrated by solving the ITG dispersion relation, (3.14), for the mode complex frequency $\omega = \omega _r + {\rm i} \gamma$, as a function of $J$, and comparing its prediction with the one from the GK dispersion relation, given in (3.7). Figure 2 shows the ITG growth rate $\gamma$ and the corresponding real frequency $\omega _r$ solution of the gyromoment sITG dispersion relation given in (3.14), plotted as a function of the perpendicular wavenumber $k_\perp$ for increasing values of $J$. For comparison, the collisionless GK solution, obtained from (3.7), is shown. Both $\gamma$ and $\omega _r$, obtained by the dispersion relation expressed in terms of gyromoments, (3.14), approach the one predicted by (3.7) as $J \to \infty$. Convergence is achieved with small $J$ at long wavelengths (e.g. the growth rate peak at $k_\perp \simeq 0.6$ is well reproduced with $J \sim 5$). On the other hand, for small perpendicular wavelengths and as a rule of thumb, $J \gtrsim k_\perp ^{2}$ must be retained in order to retrieve accurately the GK collisionless prediction (e.g. the $k_\perp \simeq 5$ case is well reproduced with $J \gtrsim 25$). We remark that, if $J$ is not sufficiently large, the truncation of the quantities, used in (3.15), can yield significant deviations in the growth rate $\gamma$ and frequency $\omega _r$. The presence of collisions is expected to affect the convergence estimate reducing $J$ even for short wavelength modes, as numerically demonstrated in § 8.

Figure 2. Collisionless sITG growth rate, $\gamma$ (a), and frequency, $\omega _r$ (b), as a function of the perpendicular wavenumber $k_\perp$ for different values of $J$. The coloured lines are the solution of the gyromoment hierarchy dispersion relation, given in (3.14), while the black lines are solutions of (3.7), i.e. the $J = \infty$ solution. Here, $\eta =3$ and $k_\parallel = 0.1$.

4. High-collisional limit

We now investigate the high-collisional limit of the gyromoment hierarchy, (2.14), focusing on two reduced models that retain only a finite number of gyromoments. In the high-collisionality limit, high-order gyromoments are damped by the presence of collisions and the ITG mode can be accurately described by considering only the evolution of the lowest-order gyromoments. Indeed, the perturbed gyrocentre distribution function is expected not to significantly deviate from a perturbed Maxwellian at high collisionalities. We consider, first, a reduced model based on six gyromoments (6GM), obtained by using a closure by truncation of the gyromoment hierarchy. Second, we consider a reduced model based on four gyromoments (4GM), rigorously derived as an asymptotic limit of the 6GM using the Chapman–Enskog closure to express the highest-order gyromoments as a function of the lower-order ones in the limit of $\epsilon \sim k_\parallel \lambda _{{\rm mfp}} \ll 1$. While the Chapman–Enskog procedure allows us to reduce the number of gyromoments, closures for the FLR terms associated with the $\mathcal {K}_{n}$ kernel functions and inherent to gyrofluid models (Brizard Reference Brizard1992; Dorland & Hammett Reference Dorland and Hammett1993; Beer & Hammett Reference Beer and Hammett1996; Waltz et al. Reference Waltz, Staebler, Dorland, Hammett, Kotschenreuther and Konings1997; Madsen Reference Madsen2013a; Held, Wiesenberger & Kendl Reference Held, Wiesenberger and Kendl2020), are still an open issue in the literature. The closure issue appears explicitly in the Fourier exponential form of the kernel function defined in (2.13). Appendix E is dedicated to exploring possible FLR closures by using a Padé approximation technique applied to the 4GM model. To make analytical progress, we focus here on long wavelength modes. Consistently, we model collisions using the DK Coulomb operator (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). This has the advantage of avoiding infinite sums with FLR terms that need to be truncated depending on the value of the perpendicular wavenumber, $k_\perp$.

The present section is structured as follows. In § 4.1, we derive the 6GM and 4GM reduced models. The predictions are compared with the ones of the full gyromoment hierarchy and are also verified against the GK code GENE. In § 4.2, we explore the high-collisional limit by deriving a general algebraic collisional dispersion relation from which we discuss the limits of the slab and toroidal branches of the ITG mode.

4.1. 6GM and 4GM reduced models

We start by deriving the 6GM model. We introduce the lowest-order normalized fluid gyromoments, i.e. the gyrocentre density $N = N^{00}$, the parallel velocity $V = N^{10}/\sqrt {2}$, the parallel and perpendicular temperatures, $T_\parallel = \sqrt {2} N^{20} + N$ and $T_\perp = N - N^{01}$ and, finally, the parallel and perpendicular heat fluxes, associated with the deviations of ${g}$ from the Maxwellian distribution function, $Q_\parallel = N^{30}$ and $Q_\perp = N^{11}$. The evolution equations of these gyromoments are obtained by setting $(p,j) = (0,0)$, $(1,0)$, $(2,0)$, $(0,1)$, $(3,0)$ and $(1,1)$ in (2.14), respectively, and by neglecting all higher-order gyromoments. This yields

(4.1a)\begin{align} & \frac{\partial}{\partial t } N + {\rm i} k_\parallel \sqrt{2 \tau} V + {\rm i} \tau \omega_{B} \left( T_\parallel{-} T_\perp \right) + {\rm i} \varGamma_1 \phi = 0 \end{align}
(4.1b)\begin{align} & \sqrt{2} \frac{\partial}{\partial t }V + \frac{{\rm i} k_\parallel}{\sqrt{\tau}} \mathcal{K}_{0} \phi + {\rm i} k_\parallel \sqrt{\tau} T_\parallel{+} {\rm i} \tau \omega_B \left( 4 \sqrt{2} V- Q_\perp{+} \sqrt{6} Q_\parallel \right)= 0,\end{align}
(4.1c)\begin{align} & \frac{\partial}{\partial t } T_\parallel{+} {\rm i} k_\parallel \sqrt{6\tau} Q_\parallel{+} {\rm i} \tau \omega_B \left( 5 T_\parallel{-} 4 N + T_\perp \right) \nonumber\\ & \qquad + {\rm i} (\varGamma_2 - \varGamma_1) \phi ={-} \nu \alpha \left( T_\parallel{-} T_\perp \right),\end{align}
(4.1d)\begin{align} & \frac{\partial}{\partial t }T_\perp{+} {\rm i} k_\parallel \sqrt{\tau} ( \sqrt{2 } V+Q_\perp ) - {\rm i} \tau \omega_B \left( 3 N - T_\parallel{-}3 T_\perp \right) \nonumber\\ & \qquad - {\rm i} (\varGamma_{3}- \varGamma_1) \phi = \nu \frac{\alpha}{2}\left( T_\parallel{-} T_\perp\right), \end{align}
(4.1e)\begin{align} & \frac{\partial}{\partial t } Q_\parallel{+} {\rm i} k_\parallel \sqrt{\frac{3 \tau}{2}} (T_\parallel{-} N) + {\rm i} \tau \omega_B \left( 2 \sqrt{3} V + 8Q_\parallel\right)={-} \nu \left( \frac{8}{5}\sqrt{\frac{2}{\rm \pi}} Q_\parallel{+} \frac{8}{5 \sqrt{3 {\rm \pi}}} Q_\perp \right), \end{align}
(4.1f)\begin{align} & \frac{\partial}{\partial t }Q_\perp{+} \frac{{\rm i} k_\parallel }{\sqrt{\tau}} \mathcal{K}_{1} \phi + {\rm i} k_\parallel \sqrt{\tau} ( N - T_\perp)\nonumber\\ & \qquad - {\rm i} \tau \omega_B \left( \sqrt{2} V - 6 Q_\perp \right)={-} \nu \left( \frac{28}{15}\sqrt{\frac{2}{\rm \pi}} Q_\perp{+} \frac{8}{5 \sqrt{3 {\rm \pi}}} Q_\parallel \right), \end{align}

where we introduce the numerical coefficients $\alpha = 16/15 \sqrt {2 / {\rm \pi}}$. We remark that the collisional terms in (4.1) are obtained by using the lowest-order gyromoments of the DK Coulomb collision operator (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021) and the energy constraints expressed in (2.24ac). The quantities appearing in the evolution equations, (4.1), are defined by

(4.2a)\begin{gather} \varGamma_1 = \mathcal{K}_{0} \left( 2 \omega_B - \omega_* \right) + \mathcal{K}_{1} \left( \omega_* \eta - \omega_B\right) , \end{gather}
(4.2b)\begin{gather}\varGamma_2 = \left( 2 \omega_B- \omega_* \eta \right) \mathcal{K}_{0}, \end{gather}
(4.2c)\begin{gather}\varGamma_3 = \mathcal{K}_{0} \left( \omega_* \eta - \omega_B \right) + \mathcal{K}_{1} \left( 4 \omega_B- \omega_* ( 1 + 2 \eta)\right) + 2 \mathcal{K}_{2} \left( \eta \omega_* - \omega_B \right). \end{gather}

Additionally, the electrostatic potential, $\phi$, is determined by the self-consistent quasineutrality condition obtained from (2.18), i.e.

(4.3)\begin{equation} \left[ 1 + \frac{q^{2}}{\tau} \left( 1 - \sum_{j=0}^{1} \mathcal{K}_{j}^{2} \right) \right]\phi = q(\mathcal{K}_{0} + \mathcal{K}_{1} ) N - \mathcal{K}_{1}T_\perp, \end{equation}

where the infinite sum appearing in the polarization term in (2.18) is truncated to $j = 1$ consistently with the fact that the FLR terms are represented up to $\mathcal {K}_{1}$ in the right-hand side of the same equation. The system in (4.1) with (4.3) constitutes a closed set of fluid-like equations that we refer to as the 6GM model.

To derive the 4GM model, further reducing the number of gyromoments appearing in the 6GM model, we apply the Chapman–Enskog asymptotic closure scheme (Chapman & Cowling Reference Chapman and Cowling1941; Jorge et al. Reference Jorge, Ricci and Loureiro2017) to the 6GM model. We introduce the dimensionless small parameter $\epsilon \sim \lambda _{{\rm mfp}} k_\parallel \sim \sqrt {\tau } k_\parallel / \nu \ll 1$, which quantifies the importance of the non-Maxwellian part of ${g}$. Thus, we expect that the Chapman–Enskog closure and the reduced models derived in this section to be valid when typically $\nu \gg \sqrt {\tau } k_\parallel$. The small parameter $\epsilon$ allows us to express the heat fluxes $Q_\parallel$ and $Q_\perp$ as a function of the lowest-order gyromoments, by applying the ordering $Q_\parallel \sim Q_\perp \sim \epsilon N$ (with $T_\parallel \sim T_\perp \sim V \sim N$) and $\partial _t \sim \gamma \sim \epsilon \nu$. We neglect the terms proportional to $\omega _B$, which are smaller than the terms proportional to $k_\parallel T_\parallel$ and $k_\parallel T_\perp$ since $\omega _B / k_\parallel \lesssim 1$ in the boundary region. With these assumptions, (4.1c) and (4.1d) can be solved at the leading order in $\epsilon$, yielding

(4.4a)\begin{gather} Q_{{\parallel} } = {\rm i} k_\parallel \frac{\chi_{{\perp}}^{{\parallel}} }{\tau} \mathcal{K}_{1} \phi + {\rm i} k_\parallel \chi_{{\perp}}^{{\parallel}} ( N - T_\perp) - {\rm i} k_\parallel \chi_{{\parallel}}^{{\parallel}} (T_\parallel{-} N ), \end{gather}
(4.4b)\begin{gather}Q_\perp{=}{-} {\rm i} k_\parallel \frac{\chi_\perp^{{\perp}} }{\tau} \mathcal{K}_{1} \phi - {\rm i} k_\parallel \chi_\perp^{{\perp}} ( N - T_\perp) + {\rm i} k_\parallel \chi^{{\perp}}_{{\parallel}} (T_\parallel{-} N ), \end{gather}

where we introduce the thermal conductivities,

(4.5a–d)\begin{align} \chi_{{\perp}}^{{\parallel}} = \frac{5}{16}\sqrt{\frac{\rm \pi}{3}} \frac{\sqrt{\tau}}{\nu}, \quad \chi_{{\parallel}}^{{\parallel}} = \frac{35}{32} \sqrt{\frac{2 {\rm \pi}}{3}}\frac{\sqrt{\tau}}{\nu}, \quad \chi_\perp^{{\perp}} = \frac{5\sqrt{2{\rm \pi}}} {16}\frac{\sqrt{\tau}}{\nu}, \quad \chi_\parallel^{{\perp}} = \frac{5\sqrt{\rm \pi}}{16} \frac{\sqrt{\tau}}{\nu}. \end{align}

We remark that $\chi _{\parallel }^{\parallel } > \chi _{\perp }^{\parallel }$ and that $\chi _\perp ^{\perp } > \chi _\parallel ^{\perp }$. Equation (4.1), with the closure relations expressed in (4.4), yields the 4GM for the lowest-order gyromoments fluid quantities $N$, $V$, $T_\parallel$ and $T_\perp$ and provides the rigorous $\nu \gg \sqrt {\tau } k_\parallel$ asymptotic limit of the gyromoment hierarchy. We remark that an analytical procedure to close the gyromoment hierarchy similar to the one considered here for the ITG mode, can be applied in the more general cases that include variation of equilibrium quantities along the magnetic field line and electron dynamics.

To study the limits of validity of the 6GM and 4GM models, we compare in figure 3 the estimate of the ITG growth rate $\gamma$ and real frequency $\omega _r$ with the results of the gyromoment hierarchy using the DK Coulomb collision operator (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021) and the GK collisionless dispersion relation given in (3.3). To compute $\gamma$ and $\omega _r$ from the gyromoment hierarchy equation (see (2.14)) and from the 6GM and 4GM models (see (4.1) and (4.4)), an initial-value calculation as well as a direct eigenvalue solver are used, and their results are compared for consistency. The same numerical methods are used to evaluate all the results reported in the present paper. For verification purposes, the results of GENE simulations using the same operator (Jenko et al. Reference Jenko, Dorland, Kotschenreuther and Rogers2000) are also shown. Here, we consider $(P,J) = (18,6)$ with a closure by truncation (convergence studies are reported in § 8) and $k_\perp = 0.5$, which corresponds approximately to the fastest growing ITG mode. The results are shown as a function of the collisionality $\nu$, for two different normalized temperature gradients, $\eta =3$ and $\eta =5$. First, an excellent agreement is observed between the gyromoment approach and the GENE simulations at arbitrary collisionality $\nu$, using the same DK Coulomb collision operator, confirming the preliminary study reported in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). Second, both GENE and the gyromoment ITG growth rates and real frequencies agree well with the collisionless limit (i.e. $\nu \ll 1$) and with the high-collisional predictions (i.e. $\nu \gg \sqrt {\tau } k_\parallel \sim 0.1$), obtained from the 6GM and 4GM models. Overall, the 6GM provides a better approximation to the low-collisionality solution than the 4GM, while they both agree in the high-collisional limits. It is remarkable that the 6GM retrieves surprisingly well the ITG growth rate $\gamma$ when $\eta =5$, but fails to predict correctly the mode frequency $\omega _r$. Finally, we remark that the 4GM yields an unphysical $\omega _r < 0$, when the collisionality is below $\nu \lesssim 10^{-2}$ with $\eta =3$. This is mainly due to the fact that, in the low-collisionality limit, the $1/\nu$ dependence of the thermal conductivities given in (4.5ad) produces arbitrarily large heat fluxes. The good agreement observed in figure 3 illustrates the multi-fidelity nature of the gyromoment approach that allows the derivation of reduced models able to address both the collisionless and high-collisional limits of the GK equation, an attractive feature in the assessment of the plasma dynamics in the boundary region.

Figure 3. The sITG growth rate $\gamma$ (a) and real frequency $\omega _r$ (b) as a function of collisionality $\nu$ obtained using the gyromoment approach with the DK Coulomb collision operator (solid blue line), GENE with the same operator (red markers), the collisionless GK dispersion relation (dotted black) and the 4GM and 6GM high-collisional limits (red dashed–dotted and red dotted lines, respectively). The solution of the dispersion relation given by (4.6), derived in the case of $Q_{\parallel, \perp } = 0$, is also plotted for comparison by the solid thin black lines. Here, we consider $k_\parallel = 0.1$, $k_\perp = 0.5$ and $\tau = 1$.

4.2. High-collisional dispersion relation

We now derive an algebraic dispersion relation using the 4GM equations and discuss the limits of the slab and toroidal branch of the ITG mode. As collisional effects enter through the parallel and perpendicular thermal conductivities such that $\chi ^{\parallel }_{\perp } \sim \chi ^{\parallel }_{\parallel } \sim \chi ^{\perp }_{\parallel } \sim \chi ^{\perp }_{\perp } \sim 1/\nu$ (see (4.5ad)), we assume that the heat fluxes become negligible when $\nu \gg k_\parallel$, i.e. $Q_{\parallel } \simeq 0$ and $Q_\perp \simeq 0$ (see (4.4)). Hence, a dispersion relation can be derived from the 4GM equations, (4.1), yielding

(4.6)\begin{equation} a \omega^{4} + b \omega^{3} + c \omega^{2} + d \omega + e = 0, \end{equation}

where

(4.7a)\begin{align} a & = \varGamma_\phi , \end{align}
(4.7b)\begin{align} b & = \frac{3 }{2} {\rm i} \alpha \nu \varGamma_\phi - \varGamma_1 \mathcal{K}_{0} - \varGamma_3 \mathcal{K}_{1}- 16 \varGamma_\phi \tau \omega_{B}, \end{align}
(4.7c)\begin{align} c & = \frac{ {\rm i} \alpha \nu}{2} \left( \mathcal{K}_{1}\left( \varGamma_2 - 2 \varGamma_3\right) - 32 \tau \varGamma_\phi \omega_B \varGamma_1 \mathcal{K}_{0} \right) \nonumber\\ & \quad - 3 k_\parallel^{2} \tau \varGamma_\phi - k_\parallel^{2} \mathcal{K}_{0}^{2} + \tau \omega_B \left( \mathcal{K}_{0} \left( 14 \varGamma_1 - \varGamma_2 + \varGamma_3 \right) + \mathcal{K}_{1} \left(\varGamma_1 + 12 \varGamma_3 \right) \right), \end{align}
(4.7d)\begin{align} d& = {\rm i} \alpha \nu \left[ k_\parallel^{2} \left( \mathcal{K}_{0} \mathcal{K}_{1} - \frac{5}{2} \tau \varGamma_\phi - \frac{3}{2} \mathcal{K}_{0}^{2} \right)\right.\nonumber\\ & \quad \left.+ \tau \omega_B \left( \varGamma_1 \left(13 \mathcal{K}_{0} + 2 \mathcal{K}_{1} \right) - (\varGamma_2 - 2 \varGamma_3)(\mathcal{K}_{0} + 3 \mathcal{K}_{1}) \right)\vphantom{\left( \mathcal{K}_{0} \mathcal{K}_{1} - \frac{5}{2} \tau \varGamma_\phi - \frac{3}{2} \mathcal{K}_{0}^{2} \right)} \right] \nonumber\\ & \quad + \tau k_\parallel^{2} \left[ \omega _{B } \left(18 \tau \varGamma _{\phi }+\mathcal{K}_{0} \left(8 \mathcal{K}_{0} +\mathcal{K}_{1} \right)\right)+ \left(2 \varGamma _1-\varGamma _2\right) \mathcal{K}_{0} +3 \varGamma _3 \mathcal{K}_{1} \right], \end{align}
(4.7e)\begin{align} e & = {\rm i} \alpha \nu \left[ 5 k_\parallel^{2} \tau \omega_B \left( \tau \varGamma _{\phi } + \mathcal{K}_{0}^{2} \right) +k_\parallel^{2} \tau \left( \varGamma _1- \frac{\varGamma _2}{2}+ \varGamma _3\right) \left(\mathcal{K}_{0}+ \mathcal{K}_{1}\right] \right) \nonumber \\ &+ k_\parallel^{2} \tau^{2} \omega_B \left(\mathcal{K}_{1} \left(\varGamma _2-6 \varGamma _3 -2 \varGamma _1 \right) -2 \mathcal{K}_{0} \left(4 \varGamma _1-2 \varGamma _2+\varGamma _3\right) \right), \end{align}

being terms proportional to $\omega _B^{2}$ neglected ($\omega _B \ll 1$) and having defined $\varGamma _\phi = 1 + q^{2} (1 - \mathcal {K}_{0}^{2} - \mathcal {K}_{1}^{2}) / \tau$. In (4.6), collisional effects are represented through the terms proportional to ${\rm i}\alpha \nu$. The solution of (4.6) is plotted in figure 3 and approaches the 6GM and 4GM models when $\nu \gg k_\parallel$, such that $Q_\parallel = Q_\perp \simeq 0$.

We now consider the slab and toroidal limits of (4.6) separately in detail. An estimate of the sITG peak can be derived from (4.6) as a function of $k_\parallel$, $\nu$ and $\eta$. Besides neglecting the terms proportional to $\omega _B$ and FLR effects (i.e. imposing $\mathcal {K}_{0} \simeq 1$, $\mathcal {K}_{1} \simeq 0$, $\mathcal {K}_{2} \simeq 0$, such that $\varGamma _\phi \simeq 1$), we assume that the sITG mode is far from the marginal stability limit $\eta \sim \nu \gg 1$. Then, the dispersion relation in (4.6) reduces to

(4.8)\begin{align} \hat{\omega}_*^{4} + \frac{3}{2} {\rm i} \alpha \hat{\nu}_* \hat{\omega}_*^{3} - \hat{k}_{{\parallel} *}^{2} \left( 1 + 3 \tau\right) \hat{\omega}_*^{2} + \hat{k}_{{\parallel} *}^{2} \left[ \tau - \frac{3}{2} {\rm i} \alpha \hat{\nu}_* \left( 1 + \frac{5}{2} \tau\right)\right] \hat{\omega}_* + \frac{3}{2} {\rm i} \alpha \hat{k}_{{\parallel} *} \hat{\nu}_* \tau =0, \end{align}

where we introduce $\hat {\omega }_* = \omega / (\eta \omega _*)$, $\hat {\nu }_* = \nu / (\eta \omega _*)$ and $\hat {k}_{\parallel *} = \hat {k}_\parallel / (\eta \omega _*)$. Given the normalization of (4.8), the maximum growth rate is expected to be proportional to $\gamma \simeq g(\tau, \hat {\nu }_*) \eta \omega _*$ and occurs at $k_\parallel \simeq f(\tau, \hat {\nu }_*) \eta \omega _*$. In figure 4, we show the maximum ITG growth rate, $\hat {\gamma }_* = \gamma / (\eta \omega _*)$ (i.e. the imaginary part of $\hat {\omega }_*$) obtained by solving (4.8) as a function of $\hat {k}_{\parallel *}$ and $\hat {\nu }_*$. The effect of increasing $\hat {\nu }_*$ is primarily to shift the growth rate peak to larger values of $\hat {k}_\parallel ^{*}$ for all values of $\tau$. Additionally, while $g(\tau,\nu )$ increases with $\tau$ for all $\nu$, it is found that $f(\tau, \hat {\nu }_*)$ is an increasing function of $\tau$ at large values of $\hat {\nu }_*$.

Figure 4. Estimates of the normalized sITG peak, $\hat {\gamma }_*$, obtained from (4.8), for increasing values of $\hat {\nu }_* = 0.015$ (from (a) to (b)), as a function of $\hat {k}_{\parallel *}$ and the temperature ratio $\tau$. The local maxima of the growth rate is indicated by the dashed black lines. The growth rates $\hat {\gamma }_*$ are normalized to their maximal values.

We remark that the collisionality dependence of the functions $f$ and $g$ is not present in previous fluid ITG models based on the drift-reduced fluid equations (Mosetto et al. Reference Mosetto, Halpern, Jolliet, Loizu and Ricci2015). This dependence arises here because the 4GM model (and more generally the gyromoment hierarchy equation in (2.14)) allows for a different evolution of the perpendicular and parallel temperature fluctuations assuming, in general, $T_\parallel \neq T_\perp$. On the other hand, the drift-reduced Braginskii fluid model in, e.g. Zeiler et al. (Reference Zeiler, Drake and Rogers1997), assumes that $T_\parallel = T_\perp$. In fact, the $\nu$ dependence in (4.8) enters in the dispersion relation, given in (4.8), through terms proportional to $\nu \alpha (T_\parallel - T_\perp )$ which arise from the gyromoment $\mathcal {C}^{20}$ and $\mathcal {C}^{01}$ components of the DK Coulomb collision operator.

To investigate the effect of anistotropic temperature fluctuations, we derive the evolution equation of the temperature, $T = (T_\parallel + 2 T_\perp ) /3 - N$. Using (4.1c) and (4.1d) and the continuity equation, (4.1a), this can be expressed as

(4.9)\begin{equation} \frac{\partial }{\partial t} T - {\rm i} k_\parallel \frac{ \sqrt{2 \tau}}{3}V + \frac{{\rm i}}{3} \left( \varGamma_2 - 2 \varGamma_3 - 2 \varGamma_1\right) \phi =0. \end{equation}

From (4.9) with (4.1a) and (4.1b), and considering the limits used to derive (4.8), we obtain the estimate of the maximum sITG growth rate,

(4.10)\begin{equation} \hat{\omega}_*^{3} - \hat{k}_{{\parallel} *}^{2} \left( 1 + \frac{2}{3}\tau \right) \hat{\omega}_* + \hat{k}_{{\parallel}*}^{2} \tau = 0. \end{equation}

From (4.10), one finds that the peak of the sITG growth rate is proportional to $\gamma \simeq g(\tau ) \eta \omega _*$ and occurs at $k_{\parallel } \simeq f(\tau ) \eta \omega _*$. It is found that $g(\tau )$ is an increasing function of $\tau$, while $f(\tau )$ decreases with $\tau$. We remark that (4.10) agrees with Mosetto et al. (Reference Mosetto, Halpern, Jolliet, Loizu and Ricci2015) by replacing the numerical factor $2/3$ by $5/3$, a difference due to the additional degree of freedom assumed in the 4GM model.

We now turn to the toroidal ITG mode, driven unstable by the magnetic drifts and by the presence of a finite temperature gradient. We derive a simple dispersion relation for the toroidal ITG from (4.6) indicative of its properties. Since, in contrast to the sITG, the toroidal ITG subsists at $k_\parallel = 0$, we set $k_\parallel =0$ thus removing the coupling with sound waves, and we adopt the long perpendicular wavelength limit, keeping terms up to $O(k_\perp ^{2})$, such that $\mathcal {K}_{0} \simeq 1 - \tau k_\perp ^{2}/2$, $\mathcal {K}_{1} \simeq \tau k_\perp ^{2}/2$, and $\mathcal {K}_{2} \simeq 0$. Since $\omega _B \ll \omega _*$ in the boundary region, we obtain

(4.11)\begin{equation} \omega^{3} + b_t \omega^{2} + c_t \omega + d_t = 0, \end{equation}

where

(4.12a)\begin{align} b_t & = \frac{3}{2} {\rm i} \alpha \nu + \omega_* - 2 \omega_B \left( 1 + 8 \tau \right) + k_\perp^{2} \left[ 2 \omega_B \left( 1 + \frac{3}{2} \tau \right) - \omega_* \left( 1 + \tau + \eta \tau \right) \right], \end{align}
(4.12b)\begin{align} c_t & = \frac{3}{2 } {\rm i} \alpha \nu \left[ \omega_* - 2 \omega_B - \frac{32}{3} \tau \omega_B + k_\perp^{2} \left( 2 \omega_B \left( 1 + \frac{19}{12} \tau \right) -\omega_* (1 + \tau + \tau \eta) \right) \right] \nonumber\\ & \quad + \omega_B \omega_* \left[ 2 \tau ( \eta -7) + k_\perp^{2} \left(10 \eta \tau \left(\tau - \frac{1}{5}\right) +\tau\left(13 \tau +14\right) \right) \right], \end{align}
(4.12c)\begin{align} d_t & ={\rm i} \alpha \nu \omega _{B } \omega _* \tau \left[ \left(3 \eta -13 \right) + k_\perp^{2} \left(3 \eta \left(2 \tau-1 \right) +11 \tau +13 \right)\right]. \end{align}

An estimate of the value of $k_\perp$ that yields the largest growing mode can be derived from (4.11) by neglecting terms proportional to $\nu$, for simplicity. Then, the dispersion relation in (4.11) becomes a second-order equation for $\omega$ presenting the largest solution $\omega$ when the coefficient $b_t$ is minimized. Hence, by minimizing $b_t$ and considering $\omega _B \ll \omega _*$, the growth rate is maximized when $k_\perp ^{2} \simeq 1 / [ 1 + \tau (1 + \eta )]$. Figure 5 shows the toroidal ITG growth rate obtained by solving (4.11), as a function of the perpendicular wavenumber $k_\perp$ and temperature ratio $\tau$. We observe that the decrease of the perpendicular wavenumber of the fastest growing mode with $\tau$ is present also at finite collisionality, $\nu$. We also remark that collisions tend to destabilize long wavelength modes.

Figure 5. Toroidal ITG growth rate as a function of the perpendicular wavenumber, $k_\perp$ and the temperature ratio $\tau$ obtained from (4.11), for a low (a) and a high (b) value of collisionality. The local maxima of the growth rate are indicated by the dashed black lines. The growth rates are normalized to the maximum value. Here, $\eta = 7$ and $R_B = 0.1$.

5. ITG mode with GK coulomb collision operator

We now perform parameter scans to explore the dependence of the ITG linear growth rate on the perpendicular and parallel wavenumbers, on the level of collisionality and on the temperature and magnetic gradient strengths. Collisional effects are modelled with the GK Coulomb collision operator (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). As a reference value, we consider $k_\parallel \sim 0.1$. This is justified by the strong ballooning character of the ITG mode in the boundary region where $k_\parallel$ is of the order of $L_N / q_0 R_0$ (with $q_0$ the safety factor and $R_0$ the major radius of the tokamak devices) and typical H-mode edge values of $q_0 \sim 4$ and $L_N / R_0 \lesssim 1/50$. In addition, we note that $1 \lesssim \tau \lesssim 4$ in the boundary. Therefore, we consider a temperature ratio of $\tau = 1$ for simplicity, and note that the effect of increasing $\tau$ can be inferred from figures 4 and 5 for the slab and toroidal cases, respectively. We remark, that given the value of $k_\parallel$ and $\tau$ used in this section, the high-collisional regime is reached when $\nu \gg \sqrt {\tau } k_\parallel =0.1$ in agreement with § 4.1. For the present numerical calculations, we fix $(P,J) = (18,6)$. As we show in § 8, this choice ensures good convergence of the numerical results in the range of parameters explored in this section.

Figure 6 displays the ITG growth rate in the $(k_\parallel,k_\perp )$ parameter space when $\eta =3$ for increasing collisionality $\nu$ (from figure 6a to 6c and from figure 6d to 6f) and magnetic gradient $R_B$ (slab in figure 6ac and toroidal in figure 6ac). For comparison purposes, we show the stability boundaries obtained numerically from the gyromoment hierarchy (dotted white) and from the 6GM and 4GM reduced models (red dotted and dashed–dotted). First, we observe that, as collisionality increases, FLR collisional stabilization effects become important and eventually suppress the ITG modes when $k_\perp \gtrsim 1$ and $\nu \gtrsim 0.5$. In fact, we notice that, at intermediate and high collisionalities, the gyromoment hierarchy stability boundary follows closely the ones predicted by the 6GM and 4GM models for long parallel ($k_\parallel \lesssim 0.15$) and perpendicular ($k_\perp \lesssim 0.2$) wavelengths where the FLR collisional effects are negligible and kinetic effects are small. However, while the ITG mode is completely suppressed for $k_\perp \gtrsim 0.5$ when $\nu > 0.5$ according to the full gyromoment calculation, the 6GM and 4GM stability boundaries extend to $k_\perp \gtrsim 1$. This highlights the importance of FLR terms in the Coulomb collision operator model. In fact, despite the lack of collisionless FLR terms in the 6GM and 4GM models, these models agree with the gyromoment hierarchy at high collisionality when collisional effects are modelled by the DK Coulomb even when $k_\perp \sim 2$ (see, e.g. figure 12). We remark that the complete stabilization of the ITG mode by FLR collisional damping requires a level of collisionality (i.e. $\nu \gtrsim 0.5$) that is even higher than the typical SOL experimental values, i.e. $T \gtrsim 10$ eV, $N \sim 10^{19}$ m$^{-3}$ and $L_N \sim 5$ cm yields a simple estimate of $\nu \lesssim 0.3$. Second, at a given collisionality, the toroidal ITG mode extends to shorter perpendicular wavelength with respect to the slab case. At the same time, the peak of the growth rate of the toroidal ITG (which occurs near $k_\perp \sim 0.5$) is shifted at longer wavelengths in the parallel direction compared with the slab branch, consistently with the fact that the toroidal mode subsists in the $k_\parallel = 0$ limit. In all cases, the stabilization of the growth rate at large $k_\parallel$ is caused by Landau damping. Third, the peak growth rate of the toroidal ITG is approximately twice as large as the slab one at all collisionalities.

Figure 6. The ITG linear growth rate, $\gamma$, in the $(k_\perp, k_\parallel )$ parameter space, for increasing collisionality $\nu$ (from (a) to (c) and from ( d) to ( f)) and magnetic gradient strength $R_B$ (from (ac) to (df)) at a temperature gradient strength of $\eta =3$. The stability boundary of the gyromoment hierarchy (dotted white line) and 6GM and 4GM (dotted and dot–dashed lines, respectively) are shown for comparison. Here, we fix $(P,J) = (18,6)$.

We now assess the effects of collisions on the critical value of the temperature gradient above which the ITG develops. For this purpose, we compute the ITG growth rate, maximized over $k_\parallel$ and $k_\perp$, as a function of $\nu$ and $\eta$. This allows us to evaluate the stability boundary shown in figure 7(a) for the slab and in figure 7(b) for the toroidal ITG modes and compared them with the predictions of the 6GM and 4GM models. The slab and toroidal ITG modes are destabilized by increasing $\eta$ at all collisionality and are damped (and/or suppressed) by FLR collisional effects. We remark that the gyromoment stability boundary agrees with the collisionless one as obtained from (3.3). The sITG stability boundary is highly sensitive to the level of collisionality, in particular for $\nu \gtrsim 0.1$. In fact, the critical temperature gradient for the sITG mode onset is increased from $\eta \sim 1$, when $\nu \lesssim 1$, to $\eta \sim 3$, when $\nu \gtrsim 1$, in order to overcome the strong FLR collisional stabilization effects. A similar observation can be made when using other GK collision operators, but not shown (see § 6). However, we remark that the increase of $\eta$ for the mode onset due to FLR collisional damping is expected to be much weaker in experiments because of lower values of $\nu$. On the other hand, the toroidal ITG mode stability boundary is not significantly affected by collisions and remains near the collisionless threshold. The reason of this difference can be inferred from figure 6 which shows a peak of the growth rate near $k_\perp \simeq 0.3$ when $\nu = 5$ in the toroidal case (figure 6df), while the mode is completely suppressed in the slab case (figure 6ac). Finally, we notice that the 6GM and 4GM stability boundaries converge to a constant value near $\eta \simeq 0.5$ in the high-collisional limit. At low collisionality, $\nu \lesssim 0.1$, the ITG mode is completely suppressed according to the 4GM model because of kinetic effects associated with the parallel streaming.

Figure 7. The ITG growth rate maximized over $k_\parallel$ and $k_\perp$, denoted by $\max _{k_\parallel, k_\perp }(\gamma )$, as a function of the collisionality $\nu$ and normalized temperature gradient $\eta$ for the slab (a) and toroidal (b) ITG modes. The stability boundaries (which separates the unstable and stable modes) predicted by the gyromoment approach using the GK Coulomb collision operator and the 6GM and 4GM models are plotted by the white dotted, red dotted and red dotted–dashed lines, respectively. The values of the growth rates in the low collisionality limit are in agreement with (3.3).

6. Comparisons between collision operator models

In this section, we investigate the differences between the collision operator models introduced in § 2.2. In particular, we compare them as a function of collisionality and perpendicular wavenumber and measure their relative deviation with respect to the GK Coulomb collision operator in § 6.1. We consider different collisionality regimes and temperature gradient strengths. We demonstrate that energy diffusion and FLR effects are important in the collisional damping of ITG modes, and find that the deviations, when compared with the GK Coulomb collision operator, decrease with temperature gradient strength at all collisionalities. One of the main results of this comparison is that the smallest relative deviation (${\lesssim }15\,\%$) is produced by the GK Sugama operator, while the largest (${\lesssim }50\,\%$) is obtained from the GK Dougherty operator. The importance of FLR effects is highlighted in § 6.2 by investigating a SWITG typically peaking near $k_\perp \sim 1.5$ that can be destabilized if collisional FLR terms are neglected. The present results are particularly important since DK collision operators are often used in the continuum GK codes to simulate the plasma boundary. In this section, we truncate the gyromoment hierarchy, (2.14), at $(P,J) = (18,6)$ for all cases. We have checked that convergence is achieved in all cases, in agreement with the convergence studies performed in § 8.

6.1. Comparisons between the GK coulomb and collision operators models

We focus on the differences between the GK Coulomb collision operator, the GK/DK Sugama collision operator (Sugama et al. Reference Sugama, Watanabe and Nunami2009), the GK/DK momentum-conserving pitch-angle scattering (pitch) operator (Helander & Sigmar Reference Helander and Sigmar2002), the zeroth-order DK HSC collision operator (Hirshman & Sigmar Reference Hirshman and Sigmar1976a) and the GK/DK Dougherty collision operator (Dougherty Reference Dougherty1964). The gyromoment expansion of the Sugama collision operator is detailed in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021), while the gyromoment expansions of the other collision operators are reported in appendices A–C.

Figure 8 shows the slab (figure 8ac) and the toroidal (figure 8df) ITG linear growth rate at a relative temperature gradient strength of $\eta =3$ from low (figure 8a,d) to high (figure 8cf) collisionality as a function of the perpendicular wavenumber $k_\perp$ using the GK collision operators. The collisionless and high-collisional limits are plotted for comparison. It is remarkable that the ITG growth rate is sensitive to both the collision operator models and to the presence of FLR related terms. In fact, the GK Dougherty collision operator produces the results that mostly differ with respect to the other GK operator models, yielding the strongest ITG damping. The effect of energy diffusion can be observed by comparing the GK pitch-angle and GK Sugama collision operators. While the former systematically overestimates the ITG growth rate (compared with the GK Coulomb operator) due to the absence of energy diffusion, the latter damps stronger the ITG mode with growth rates smaller (but closer) than the GK Coulomb. Finally, we remark that the GK Coulomb operator predictions approach the ones of the 6GM and 4GM models in the limit of long perpendicular wavelengths for all collisionalities, since the FLR collisional damping becomes negligible when $k_\perp \ll 1$.

Figure 8. Comparisons between GK (solid lines) collision operator models in the case of the slab (ac) and toroidal (df) ITG growth rate as a function of the perpendicular wavenumber $k_\perp$, from the low (a,d) to high (cf) collisionality regime. The collisionless (dashed black) and 6GM and 4GM (red dotted and dashed–dotted lines, respectively) growth rates are plotted for comparisons. Here, $\eta = 3$, $k_\parallel = 0.1$.

By neglecting the FLR terms, the DK collision operators lead to a larger ITG growth rate than the GK operators as shown in figure 9. Contrary to the latter, the DK collision operators can eventually destabilize a SWITG peaking near $k_\perp \simeq 1.5$. The presence of this SWITG can be inferred from figure 9 from the increase of the ITG growth rate occurring for $k_\perp \gtrsim 1$ when $\nu \gtrsim 0.1$ (see § 6.2). This behaviour is particularly pronounced in the pitch angle and HSC collision operators at short wavelengths, and it is related to the absence of energy diffusion in these operators (Ricci, Rogers & Dorland Reference Ricci, Rogers and Dorland2006; Barnes et al. Reference Barnes, Abel, Dorland, Ernst, Hammett, Ricci, Rogers, Schekochihin and Tatsuno2009), which is retained in the Coulomb, Sugama and Dougherty collision operators. We also notice that the absence of energy diffusion in the pitch angle and HSC collision operators affects the high-collisional limit. While the deviations between the DK Coulomb, Sugama and Dougherty collision operators are reduced as the collisionality increases and, ultimately, becomes negligible when $\nu \gg 1$, converging to the 6GM and 4GM solutions, a difference subsists in the pitch-angle and HSC collision operators. Finally, among the DK collision operators considered, the DK Sugama approaches better the DK Coulomb collision operator.

Figure 9. Same as figure 8 using the DK collision operator models.

The deviations of the collision operator models with respect to the GK Coulomb operator can be quantified by computing their signed relative difference of the growth rate, $\sigma _C = (\gamma _C - \gamma ) / \gamma _C$ where $\gamma _C$ is the growth rate using the GK Coulomb and $\gamma$ is the growth rate obtained using the other collision operator models considered in figures 8 and 9. The results are shown in figure 10 where we plot $\sigma _C$ as a function of collisionality, $\nu$, and of temperature gradient strength, $\eta$, for the toroidal branch ($R_B = 0.1$). The GK Sugama collision operator provides the results that best approach the GK Coulomb collision operator. Among the other GK operators considered, the largest deviation is produced by the GK Dougherty where the relative difference is larger than $30\,\%$ when $\nu \gtrsim 0.5$. Surprisingly, while the GK pitch-angle operator systematically overestimates the ITG growth rate (as observed from $\sigma _C< 0$), it yields an error smaller than (or of the order of) $30\,\%$ in almost all the parameter space considered. In particular, the relative signed difference of the GK pitch-angle collision operator is similar in absolute value to the one of the GK Sugama at high $\eta$ and all collisionalities, but with opposite sign. We remark also that the negative sign of $\sigma _C$ for the DK collision operators is the result of ignoring the FLR collisional damping. For all the considered operators, the relative deviation with GK Coulomb operator decreases at steeper temperature gradients while it increases with collisionality.

Figure 10. Signed relative difference, $\sigma _C$ (defined in the main text), of the growth rate obtained using the collision operator models, used in figures 8 and 9, with respect to the GK Coulomb operator as a function of collisionality $\nu$ and temperature gradient strength $\eta$. The colourbars are saturated at a maximal relative deviation of $\sigma _C = \pm 0.4$. Here, $k_\parallel = 0.1$, $k_\perp = 0.5$ with $R_B = 0.1$ are considered, corresponding approximately to the peak ITG growth rate.

The gyromoment approach allows for the investigation of the ITG eigenvalue spectrum at finite collisionality. An example is shown in figure 11 for the toroidal case when $\nu = 0.1$ using the DK (figure 11a) and GK (figure 11b) pitch angle, Dougherty, Sugama and Coulomb collision operators near to the peak growth rate occurring at $k_\perp = 0.6$. First, focusing on the spectrum of the DK collision operators, we note that, while the DK Sugama reproduces qualitatively well the eigenvalues of the DK Coulomb, the DK Dougherty collision operator displays an eigenvalue spectrum with strongly damped subdominant modes characterized by $\gamma < 0$. Similar differences between the DK Dougherty and the DK Coulomb collision operators in the eigenvalue spectrum is reported in the study of electron plasma waves at arbitrary collisionality by Jorge et al. (Reference Jorge, Ricci, Brunner, Gamba, Konovets, Loureiro, Perrone and Teixeira2019b). In addition, despite the absence of energy diffusion, the DK pitch-angle operator features similar eigenvalues than the DK Coulomb. Focusing now on the spectrum of the GK collision operators, we remark that a large difference exists between the GK Dougherty operator and the other GK collision models. It is noticeable that the GK Coulomb collision operator produces subdominant modes that are less damped than the ones predicted by the GK Sugama and pitch-angle operators but, nevertheless, yields a similar eigenvalue spectrum. We remark that the eigenvalues have all a positive real frequency, $\omega _r > 0$, corresponding to the ion diamagnetic direction.

Figure 11. Toroidal ITG eigenvalue spectrum when $\nu = 0.1$ and $k_\perp = 0.6$ using the DK (a) and GK (b) collision operator models. Here, $R_B = 0.1$ and $\eta = 3$.

6.2. Destabilization of the SWITG mode

An initial analysis of the effects of FLR terms in the collision operators can be carried out by observing, first, the discrepancies between the GK Coulomb and the 6GM and 4GM stability boundaries in the $\nu =1$ case (see figure 6) and, second, from the comparison between GK and DK collision operators (see figures 8 and 9). A clear deviation between GK and DK Coulomb operators occurs at sufficiently steep temperature gradients, $\eta \gtrsim 3$, and intermediate collisionality, i.e. when $\nu \gtrsim 0.1$. For a more detailed analysis, we plot in figure 12 the ITG growth rate evaluated with the DK Coulomb in the slab case for two different values of temperature gradient strength.

Figure 12. Slab SWITG growth rate (unstable when the DK Coulomb model is used) when $\eta =3$ (ac) and $\eta =5$ (df). The stability boundaries of the gyromoment hierarchy (dotted white line), 6GM and 4GM (red dotted and dot–dashed lines, respectively) are plotted for comparison. The main ITG branch is identified near $k_\perp \simeq 0.5$ and the SWITG near $k_\perp \simeq 2$, appearing as the collisionality increases when $\nu \gtrsim 0.5$. The colourbar is saturated at the maximum of $\gamma$ when $\nu = 0.05$.

The comparison between figures 6 and 12 reveals that a short wavelength branch of the ITG, referred to as the SWITG mode and typically peaking at $k_\perp \sim 2$, is present when $\nu \gtrsim 0.5$ if the DK Coulomb operator is used. We remark that the SWITG is also captured by the 6GM and 4GM models. This mode has been identified in the collisionless limit as a continuous extension of the main branch of the ITG mode in the $k_\perp > 1$ region. In fact, the collisionless ITG growth rate reveals a typical ‘double-humped’ behaviour when plotted as a function of $k_\perp$ (Smolyakov, Yagi & Kishimoto Reference Smolyakov, Yagi and Kishimoto2002; Gao et al. Reference Gao, Sanuki, Itoh and Dong2003, Reference Gao, Sanuki, Itoh and Dong2005). Additionally, figure 12 shows the presence of a stability region that increases with collisionality, isolating the conventional ITG and the SWITG modes. We note that, in the toroidal case, the SWITG is stabilized (or even suppressed) by increasing magnetic gradients, a similar feature observed also in the collisionless limit (Gao et al. Reference Gao, Sanuki, Itoh and Dong2005).

The SWITG mode is also predicted by all the other DK collision operators considered in this work, but it is suppressed by FLR collision damping in the GK collision operators. This questions the applicably of DK collision operator models to correctly predict the turbulent transport level in ITG driven turbulence in the boundary region, since the SWITG can produce an additional (yet spurious) contribution in addition to the main ITG mode peaking at larger $k_\perp$ (Gao et al. Reference Gao, Sanuki, Itoh and Dong2005). However, further investigations are required to go beyond the local approximation used in the present work.

7. Gyromoment spectrum

We now study the ITG gyromoment spectrum. In particular, we numerically evaluate the collisionless gyromoment spectrum using the semianalytical expressions obtained in § 3.2. This reveals that magnetic drift resonance effects broaden the gyromoment spectrum compared with the slab case. The effect of collisions is also considered by deriving asymptotic expressions of the spectrum in the limit of a large number of gyromoments. We first focus on the slab limit and, then, consider the purely toroidal ($k_\parallel = 0$) case. The results presented in this section are used to analyse the convergence properties in § 8.

Considering the sITG mode, we focus on the parameters $k_\parallel = 0.1$, $\eta =5$ and $k_\perp = 0.5$ (i.e. near the growth rate peak). We first note that the solution of the GK dispersion relation, given in (3.3), yields the solution $\omega = 0.1754 + 0.1504 {\rm i}$. Given $\omega$, the normalized collisionless gyromoments spectrum, $N^{pj} / \phi$, can be evaluated using (3.11) (or equivalently (3.3) with $R_B = 0$). This is shown in figure 13 as a function of $(p,j)$ on a logarithmic scale. As observed, the collisionless gyromoment spectrum decreases rapidly with $j$, while the decrease of the magnitude of the gyromoments is slower along the $p$ direction, as highlighted in figure 13(b) that shows $N^{p0} / \phi$ and $N^{0j} / \phi$. The finite amplitude of the gyromoments $N^{pj} / \phi$ at $j \neq 0$ (and small $p$) is due to the presence of FLR effects, yielding terms proportional to $\mathcal {K}_{j}$ in the collisionless gyromoment response (see (3.11)) and, to a smaller extent, to the presence of temperature gradients, yielding terms proportional to $\eta \mathcal {K}_{j \pm 1}$ in (3.11). As a consequence, the gyromoment amplitude rapidly decreases with $j$ because of the behaviour of the FLR kernel $\mathcal {K}_{j}$ for large $j$, such that $\mathcal {K}_{j} \sim (b/2)^{2j}/j!$. We remark that the decrease of modulus of $\hat {N}^{pj}$ with $j$ is slower at large $b$ and follows, in general, $1/j!$. On the other hand, the gyromoments with $p \neq 0$ arise because of the parallel resonance terms proportional to the $p$th derivative of the plasma dispersion functions, $Z^{(p)}(\xi ) / \sqrt {2^{p} p!}$, as illustrated by figure 13. The extend of the spectrum in figure 13 demonstrates that the sITG mode requires a larger number of Hermite than Laguerre gyromoments, i.e. $P > J$, to be resolved.

Figure 13. Modulus of the collisionless gyromoments spectrum $N^{pj} / \phi$ (normalized to the maximum value) in the case of sITG, as a function of $p$ and $j$ plotted on a logarithmic scale and artificially saturated at $10^{-4}$. The modulus of the spectrum $N^{p0} / \phi$ and $N^{0j} / \phi$ are shown in panel (b) (solid blue and red lines, respectively). Here, $k_\parallel = 0.1$, $\eta =5$ and $k_\perp = 0.5$.

We now consider the collisionless gyromoment spectrum in the presence of magnetic gradients. As in the slab case, we focus on $k_\perp = 0.5$ and $\eta = 5$, but we choose $R_B = 0.5$ and $k_\parallel = 0$. We remark that setting $k_\parallel =0$ allows us to investigate the kinetic effects associated with the ${\rm i} \omega _{\boldsymbol {\nabla } B}$ term in (2.4) independently of Landau damping. The solution of the GK dispersion relation (3.3), for the complex mode frequency, yields $\omega = 0.7972 + 0.4364 {\rm i}$. The expressions given in (3.10) are used to evaluate the gyromoment spectrum. More precisely, the quantities $\mathcal {I}_n^{pj}$ are evaluated numerically using (D10) and the infinite sums in (3.10) are truncated when machine precision is reached. Similarly to figure 13, we plot the toroidal collisionless gyromoment spectrum, $N^{pj} / \phi$, on a logarithmic scale in figure 14. Compared with the slab case, the gyromoment spectrum is significantly broader. This is attributed to the velocity-dependence of the magnetic drift resonant coefficients, proportional to $1/[\omega - \omega _{\boldsymbol {\nabla } B}]$, that drives high-order gyromoments. In fact, while the dependence on $j$ of the gyromoments stems mostly from terms proportional to $\mathcal {K}_{j}$ in the slab case, the definite integrals of Gauss hypergeometric functions, appearing in $\mathcal {I}_{n}^{pj}$ (see (D10)), extend the gyromoment spectrum, producing a considerably slower decay than $1/j!$. We also notice that, since $k_\parallel =0$, the gyromoment spectrum vanishes for odd values of $p$, where $\mathcal {I}_n^{pj} = 0$ (see (D10)). Finite $k_\parallel$ yields non-vanishing odd gyromoments in $p$ due to the term proportional to $(1 + (-1)^{p-2p_1 +1})$ in (D8).

Figure 14. Modulus of the collisionless gyromoments spectrum $N^{pj} / \phi$ (normalized to the maximum value) in the case of purely toroidal ($k_\parallel = 0$) ITG with $R_B = 0.5$, as a function of $p$ and $j$ plotted on a logarithmic scale and artificially saturated at $10^{-4}$. The moduli of the spectrum of $N^{p0} / \phi$ and $N^{0j} / \phi$ are shown in panel (b) (solid blue and red lines, respectively). Here, $k_\perp = 0.5$, $\eta =5$.

Despite the fact that the gyromoment spectrum is significantly broader with respect to the slab limit when the magnetic drifts are important (large $R_B$), we numerically demonstrate the convergence of the gyromoment representation of the perturbed distribution function, ${g} /\phi$. Given the complex frequency $\omega$, a truncated distribution function ${g}^{PJ}$ can be defined from the Hermite–Laguerre expansion in (2.6) according to

(7.1)\begin{equation} {g}^{PJ} = \sum_{p=0}^{P} \sum_{j=0}^{J} \frac{\hat{N}^{pj}}{\phi} \frac{H_p(s_{{\parallel}}) L_j(x)}{\sqrt{2^{p} p!}} F_{M}. \end{equation}

From ${g}^{PJ}$ and ${g} /\phi$ obtained from (3.1), the error norm $e(P,J)$ can introduced by

(7.2)\begin{equation} e(P,J) = \int_{- \infty}^{\infty} \,{\rm d} s_\parallel \int_0^{\infty} x \left| {g}^{PJ} - \frac{{g}}{\phi} \right|^{2}, \end{equation}

and can be computed as a function of $(P,J)$. The results are shown in figure 15 for the slab and toroidal cases. Convergence is observed in both cases. However, while the magnitude of $e(P,J)$ becomes negligibly small when $(P,J) \gtrsim (20,10)$ (corresponding to $e(20,10) \lesssim 10^{-3}$) in the slab case, a larger number of gyromoments is necessary to achieve a similar error, i.e. $(P,J ) \gtrsim (50,20)$.

Figure 15. Error norm $e(P,J)$, (7.2), on a logarithmic scale for the slab (a) and toroidal (b) ITG cases. The same parameters as in figures 13 and 14 are used, respectively. The colourbars are saturated artificially at $\log _{10} e(P,J) = -4$.

We now turn to the collisional effects, which provide a natural cutoff of the gyromoment hierarchy. We aim to derive an ad hoc proxy of the scaling of $(P,J)$ at which the gyromoment hierarchy can be truncated as a function of the collisionality. In particular, we consider an equation for the generalized gyromoment amplitude $G ^{pj} = \sqrt {p} (g^{pj})^{2} /2$ (where $g^{pj} = (-1)^{j} i^{ p +j } \text {sgn}(k_\parallel )^{p} N^{pj}$ is introduced to remove the oscillatory behaviour of $N^{pj}$) using the gyromoment hierarchy equation, (2.14) in the limit $p,j \gg 1$ (Zocco & Schekochihin Reference Zocco and Schekochihin2011; Loureiro, Schekochihin & Zocco Reference Loureiro, Schekochihin and Zocco2013; Jorge et al. Reference Jorge, Ricci and Loureiro2018). In this limit, the amplitude $g^{pj}$ can be considered as a continuous complex function of $p$ and $j$. This allows us to approximate $g^{p \pm 1j} \simeq g^{pj} \pm \partial _p g^{pj}$ (and $g^{p j \pm 1} \simeq g^{pj} \pm \partial _j g^{pj}$). Assuming that an unstable mode is considered, $\partial _t G^{pj} = 2 \gamma G^{pj}$, we derive the equation for $G^{pj}$, up to first order in $1/p$ and $1/j$, i.e.

(7.3)\begin{align} |k_\parallel |\sqrt{\tau} \partial_{s} G^{pj} + \frac{ \tau}{q} \omega_{B} ( 2 j +1) \partial_j G^{pj}+ \frac{2 \tau}{q} \omega_B \left( 1 + {\rm i} \left( 2j +1\right)\right) G^{pj} = 2 \nu (f_{p,j} + b^{\alpha} h_{p,j} - \gamma ) G^{pj}, \end{align}

with $s = \sqrt {p}$. In deriving (7.3), collisional effects are modelled using the ansatz $\mathcal {C}^{pi} \simeq \nu [ f_{p,j} + b^{\alpha } h_{p,j}] N^{pj}$, where the functions $f_{p,j} = f(p,j)$ and $h_{p,j} = h(p,j)$, and the constant $\alpha$ depends on the choice of the collision operator. We remark that the coupling along the Hermite direction associated with the curvature drifts vanishes at the leading order in $1/p$ in (7.3). We solve (7.3) for the slab and toroidal $k_\parallel = 0$ ITG cases using the GK Dougherty and GK Coulomb collision operators. While the functions $f_{p,j}$ and $h_{p,j}$ are evaluated numerically for the GK Coulomb collision operator, the sparse gyromoment expansion of the GK Dougherty collision operator, given in Appendix C, allows us to solve (7.3) analytically.

Using the gyromoment expansion of the GK Dougherty collision operator, we solve first (7.3) in the case of the sITG mode. With $f_{p,j} = - p - 2j$, $h_{p,j} = -1/2$ and $\alpha = 2$ (see (C3)), the solution of (7.3) is

(7.4)\begin{equation} \left| N^{pj}\right| \sim \frac{1}{p^{1/4}} \exp\left\{ - \left( \frac{p}{p_\gamma}\right)^{1/2} - \left( \frac{p}{p_{j \perp}}\right)^{1/2}- \left(\frac{p}{p_{D\nu}} \right)^{3/2} \right\}, \end{equation}

where we introduce $p_\gamma = \tau / \gamma _\parallel ^{2}$ (with $\gamma _\parallel = \gamma / k_\parallel$), $p_{j\perp } = \tau /[ \nu _\parallel (2 j + b^{2}/2)]^{2}$ ($\nu _\parallel = \nu / k_\parallel$) and $p_{D \nu } = (9 \tau / \nu _\parallel ^{2})^{1/3}$. The quantities $p_\gamma$, $p_{j \perp }$ and $p_{D \nu }$ introduce different decay scale lengths in the gyromoment spectrum. They are associated with the growth of mode ($p_\gamma$), with the $j$ dependence of $\mathcal {C}^{pj}$ ($p_{j \perp }$) and with the $p$ dependence of $\mathcal {C}^{pj}$ ($p_{D \nu }$). The $1 / \nu _\parallel$ dependence of $p_{D \nu }$ indicates that the exponential decay of (7.4) is slower at high $k_\parallel$ and faster at high $\nu$. At high $k_\parallel$, higher-order $p$ gyromoments are excited due to parallel resonance effects. At high $\nu$, the Chapman–Enskog closure procedure ensures that $\nu _\parallel \sim \sqrt {\tau } / \epsilon \gg 1$, yielding $p_{D\nu } \sim \nu _\parallel ^{-2/3} \sim \epsilon ^{2/3} / \tau ^{1/3}$ such that that the ITG dynamics can be described by the lowest-order gyromoments at high collisionality when $\epsilon \ll 1$ (see § 4.1). Moreover, when $\nu \gg \gamma$ (such that $p_{D\nu }$ precedes $p_\gamma$) and in the DK limit (such that the $b^{2}$ term in $p_{j \perp }$ is neglected), the collisional damping provides a truncation at $P \gtrsim p_{D \nu }$. We remark that FLR terms in the collision operator have a small effect on the scaling of (7.4) at finite collisionality. In fact, near the ITG peak growth rate (i.e. $k_\perp \sim 0.5$), one typically has $j \gtrsim b^{2}/4$. Finally, we remark that the $j$ dependence in (7.3) becomes important when the second exponent becomes comparable in magnitude to the last one. Despite that $p_{j \perp } / p_{ D\nu } \sim \epsilon ^{4/3} / j^{2} \ll 1$ at large collisionality and at large $j$, the different exponents $1/2$ and $3/2$ in (7.4) imply that this happens when

(7.5)\begin{equation} \frac{p_{D \nu}^{3/2}}{p (p_{j \perp})^{1/2}} \sim 1, \end{equation}

which is satisfied when $j \sim p / 10$. This shows that the decay of the gyromoments becomes faster with $j$ because of the increase of the second exponent in magnitude.

We now derive the solution of (7.3) for the toroidal $k_\parallel = 0$ ITG case, still considering the GK Dougherty collision operator. It yields

(7.6)\begin{equation} \left| N^{pj}\right| \sim \frac{1}{p^{1/4}} \exp\left\{ - \nu_B j- \frac{1}{2}\ln\left( 1 + 2j \right) \left[ \gamma_B + 1 + \nu_B \left( p - 1+ b^{2}/2 \right) \right]\right\}, \end{equation}

where we introduce $\gamma _B = q \gamma /( \tau \omega _B)$ and $\nu _B = q \nu / (\tau \omega _B)$. The $\nu _B$ dependence highlights the broadening of the gyromoment spectrum associated with the presence of magnetic gradient drifts (see figure 14), an effect reduced at high collisionality. Similar to the slab case, FLR terms do not affect the scaling of (7.4), since the linear term dominates over the second term at large $j$ values, and $p-1 \gg b^{2} /2\sim k_\perp ^{2}$ for ITG modes.

Because the Coulomb collision operator has a non-trivial block matrix structure when projected onto the Hermite–Laguerre basis (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021), an approximate expression of $\mathcal {C}^{pj}$ needs to be determined numerically. Based on the results derived from the GK Dougherty collision operator, we assume that $\mathcal {C}^{pj} \simeq \nu f_{p,j} N^{pj}$, thus neglecting the FLR effects in the spectrum scaling. Fitting numerically the values of $\mathcal {C}^{pj}$ using the GK Coulomb collision operator yields $f_{p,j} \simeq - A j - B \sqrt {p}$ ($A \simeq 0.3$ and $B \simeq 0.8$). Similar values (within a $10\,\%$ accuracy) are typically found in the case of the DK Coulomb and show that FLR terms have a small effect on the ITG spectrum scaling when $k_\perp \sim 1$. We remark that the $\sqrt {p}$ dependence of $f_{p,j}$ is consistent with previous DK moment hierarchy using the DK Coulomb collision operator (Jorge et al. Reference Jorge, Ricci and Loureiro2018). Given $f_{p,j}$, we again solve (7.3) for the case of the slab and toroidal $k_\parallel =0$ ITG modes yielding, respectively,

(7.7)\begin{equation} \left| N^{pj} \right| \sim \frac{1}{p^{1/4}} \exp\left\{ - \left( \frac{p}{p_\gamma}\right)^{1/2} - \left( \frac{p}{p_j} \right)^{1/2} - \frac{p}{p_\nu} \right\}, \end{equation}

with $p_j = \tau / \nu _\parallel A j^{2}$ and $p_\nu = 2 \sqrt {\tau } / (B \nu _\parallel )$ and

(7.8)\begin{equation} \left| N^{pj} \right| \sim \frac{1}{p^{1/4}} \exp\left\{ - \frac{ A}{2} \nu_B j - \frac{1}{2}\ln\left( 1 + 2j \right) \left[ \gamma_B + 1 + \nu_B \left( B \sqrt{p} - A /2 \right) \right] \right\}. \end{equation}

The main differences between the results of the GK Coulomb collision operator, given in (7.7) and (7.8) and the results of the GK Dougherty operator, in (7.4) and (7.6), are that the scaling in $p$ predicted by the GK Coulomb collision operator is weaker in the slab case (compared with a $3/2$ exponent in (7.4)), and that the numerical factor in the linear term in $j$ in the toroidal case, $A / 2 \simeq 0.15$, is smaller.

We compare the scaling of the gyromoment spectra obtained using the GK Dougherty and GK Coulomb collision operators in figure 16 for the slab and toroidal ITG modes for different values of $\nu _\parallel$ (at $k_\parallel = 0.1$) and $\nu _B$ (at $R_B = 0.1$). The differences between the $p$ exponent and the numerical factors between GK Dougherty and GK Coulomb, do not significantly impact the gyromoment spectra that feature the same qualitative behaviour. The differences between the GK Dougherty and GK Coulomb observed when $\nu _\parallel = 1$ and $\nu _B =2$ stem from the difference in the exponents and numerical factors in the spectrum. Finally, we remark that the asymptotic dependence of the expansion coefficients, i.e. $N^{pj}$, on their indices $p$ and $j$ (see, e.g. (7.7) and (7.8)) shows that the gyromoment method features a spectral convergence rate that increases with collisionality (Boyd Reference Boyd2001).

Figure 16. Normalized slab (a) and toroidal (b) gyromoment spectra as a function of $p$ and $j$, respectively, obtained using the GK Coulomb collision operator from (7.7) and (7.8) (solid lines) and using the GK Dougherty collision operator from (7.4) and (7.6) (dashed lines). The slab and toroidal growth rates are estimated at $\gamma \simeq 0.1$ and $\gamma \simeq 0.2$ (see figure 8), respectively, near the ITG peaks ($k_\perp \simeq 0.5$).

8. Convergence studies

Following the analysis of the gyromoment spectrum presented in § 7, we perform convergence studies of the ITG growth rate as a function of collisionality in § 8.1, perpendicular wavenumber in § 8.2 and magnetic gradient in § 8.3 for various numbers of gyromoments. Collisions are modelled using the GK and DK Coulomb collision operators (Frei et al. Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021), since similar convergence properties are observed with the other collision operator models considered in this work. We show that the convergence of the gyromoment approach improves with collisionality.

8.1. Collisionality convergence

We first consider the convergence of the ITG linear growth mode as a function of the collisionality, $\nu$, using the GK and DK Coulomb collision operators. The collisionless and high-collisional limits are obtained by solving the collisionless GK dispersion relation, (3.3), and the 6GM and 4GM models (see § 4.1). The results are shown in figure 17(a) for the slab and figure 17(b) for the toroidal ITG modes. Different numbers of gyromoments, $(P,J)$, are used. The collisionless GK solution is retrieved for $(P,J) \simeq (18,6)$ with both the GK (dashed) and DK (solid) Coulomb collision operators. In the slab case, convergence is achieved with $P \lesssim 8$, when $\nu \simeq 0.1$, and with $P \gtrsim 18$, when $\nu \lesssim 10^{-3}$. This is in qualitative agreement with the estimate of $P$ that can be obtained from (7.7). In fact, solving the gyromoment hierarchy up to $P \gtrsim p_\nu$ yields $p_\nu \sim 3$ for $\nu = 0.1$ and $p_\nu \sim 20$ for $\nu = 10^{-3}$. Also consistent with the observations made in § 7 is the fact that only a few numbers of the Laguerre gyromoments are necessary to resolve the sITG, as illustrated by the overlap between the $(10,4)$ and $(10,6)$ lines in figure 17. This is because the sITG mode is mainly driven by the parallel dynamics, resolved by the Hermite part of the gyromoment spectrum (see figure 13), and the relatively small value of $k_\perp$ considered in this case (see figure 18). This remains also true for the toroidal ITG case considered in figure 17 because of the relatively small value of $R_B$ used. In fact, larger values of $R_B$ deteriorates the convergence as discussed below. While the GK Coulomb yields a strong FLR collisional damping as $\nu$ increases, the gyromoment hierarchy with DK Coulomb agrees with the 6GM and 4GM when $\nu \gg 1$. Nevertheless, we notice that, because of finite magnetic drifts (see § 3), agreements with the 6GM and the 4GM models occur at higher collisionality in the toroidal case. In fact, in the slab case, the high-collisional limits are retrieved when $\nu \gtrsim 1$, while for the toroidal case ($R_B = 0.1$) when $\nu \gtrsim 10^{2}$.

Figure 17. Convergence of the ITG linear growth rate, $\gamma$, as a function of the collisionality $\nu$ for different values of $(P,J)$, using the GK (dashed lines) and DK (solid lines) Coulomb operators, in the case of slab (a) and toroidal ((b) with $R_B = 0.1$) ITG branches. The collisionless limit (dotted black lines) and the high-collisional limits, 6GM and 4GM models (red dotted and dash–dot lines, respectively) are shown for comparison. Similar plots are obtained for the other GK and DK operator models. The parameters are $k_\perp = 0.5$, $k_\parallel = 0.1$, $\eta = 3$.

Figure 18. Slab ITG growth rate $\gamma$ as a function of $k_\perp$ for different values of $J$, in the low (a), intermediate (b) and high (c) collisionality regimes, with the GK (dashed lines) and DK (solid lines) Coulomb collision operators. Similar plots are obtained for the other GK and DK operator models.

8.2. Perpendicular wavenumber convergence

We now investigate the convergence associated with the FLR effects due to the gyroaverage of the electrostatic potential, i.e. with the number of Laguerre gyromoments. The coupling between the Laguerre gyromoments driven by magnetic gradients are neglected by focusing on the sITG. We scan the growth rate $\gamma$ as a function of the perpendicular wavenumber, $k_\perp$, for different values of $J$ and collisionality $\nu$. The number of the Hermite gyromoment is fixed at $P = 18$. The results are shown in figure 18. We observe that convergence occurs at lower $J$ at high collisionalities. Nevertheless, the value of $J$ necessary for convergence increases with $k_\perp$ at all collisionalities. This is in agreement with the unfavourable scaling of the FLR kernel with $b$, i.e. $\mathcal {K}_{j} \sim (b/2)^{2j}/j!$ (see § 7). We also notice that the ITG mode is strongly damped by the GK Coulomb collision operator that departs rapidly from the DK Coulomb as $k_\perp$ and $\nu$ increases, in agreement when comparing figure 8 with figure 9. The same qualitative observations hold in the toroidal case.

8.3. Magnetic gradient convergence

Finally, we consider the convergence of the gyromoment hierarchy as a function of the normalized magnetic gradient, $R_B$. As discussed in § 7, the presence of finite magnetic drifts broadens significantly the collisionless gyromoment spectrum. However, even at low collisionality, the gyromoment approach converges correctly to the collisionless solution $\gamma _{\text {GK}}$ obtained from (3.3). This is shown in figure 19 where the ratio of the ITG growth rate obtained using the gyromoment hierarchy to $\gamma _{\text {GK}}$ is plotted as a function of $P$ and $J$. It is remarkable that, consistently with § 7, a larger number of gyromoments is needed as $R_B$ increases.

Figure 19. Ratio between the ITG growth rate $\gamma$ obtained using the gyromoment hierarchy and the collisionless solution $\gamma _{\textrm {GK}}$ as a function of $P$ with $J = 10$ (a) and as a function of $J$ with $P=32$ (b) for two different values of $R_B$. Here, $k_\perp =0.5$, $k_\parallel = 0.1$ and $\nu = 0.001$.

Collisions can compete against the strong kinetic drive of the magnetic drifts by damping higher-order gyromoments (see figure 16). To investigate the convergence of the gyromoment approach at finite $R_B$ in the presence of collisions, we focus on the toroidal ITG with $k_\parallel =0$ and, hence, neglecting the coupling between the Hermite gyromoments associated with the parallel streaming. We scan the toroidal ITG linear growth rate $\gamma$ as a function of $R_B$ and increase the collisionality. The results are reported in figure 20. We observe that, at all collisionalities, the convergence of the gyromoment approach deteriorates as $R_B$ increases. Nevertheless and consistently with § 7, the number of gyromoment to achieve convergence is reduced as the collisionality increases, demonstrating the competition between the kinetic effects driven by the magnetic gradient drifts and the collisional effects. We also remark that the 4GM limit is retrieved at higher collisionality when $R_B$ increases.

Figure 20. Toroidal ($k_\parallel = 0$) ITG growth rate $\gamma$ as a function of the normalized magnetic gradient, $R_B$, for different $(P,J)$ and increasing collisionality (from (a) to (c)). Here, the GK (dashed lines) and DK (solid lines) Coulomb collision operators are shown, with the collisionless (black dotted) and high-collisional 4GM limit (red dotted–dashed lines).

9. Conclusion

In this paper, the properties of the ITG mode are studied using the GK Coulomb collision operator for the first time. The investigation is based on a Hermite–Laguerre moment expansion of the perturbed gyrocentre distribution function, which we refer to as the gyromoment approach. By projecting the GK Boltzmann equation onto the Hermite–Laguerre basis, a gyromoment hierarchy equation is deduced, which retains the effects of like-species collisions between ions thanks to the numerical implementation of the GK Coulomb collision operator described in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). Using the gyromoment hierarchy equation, the collisionless and high-collisional limits are derived analytically. In the collisionless limit, we find that the magnetic drift resonance effects significantly broaden the gyromoment spectrum compared with the slab case. In the high-collisional limit, the gyromoment hierarchy is reduced to a fluid model retaining only a finite number of gyromoments where the Chapman–Enskog asymptotic closure scheme is used. Using these analytical results, we show that the gyromoment hierarchy can retrieve the collisionless limit, where kinetic features are essential, and, at the same time, the high-collisional limit using a reduced number of gyromoments.

Numerical experiments are performed to study the ITG linear growth rate using the GK Coulomb collision operator and investigate the importance of the FLR collisional terms. We compare the GK Coulomb collision operator with other collision operator models such as the Sugama (Sugama et al. Reference Sugama, Watanabe and Nunami2009), the momentum-conserving pitch-angle scattering operator (Helander & Sigmar Reference Helander and Sigmar2002), the zeroth-order DK HSC collision operator (Hirshman & Sigmar Reference Hirshman and Sigmar1976a) and the Dougherty collision operator (Dougherty Reference Dougherty1964). The gyromoment expansions of the pitch angle and HSC collision operators are also derived for the first time here. We find that the ITG mode is strongly damped (or even suppressed) as the collisionality increases and that a steeper temperature gradient for the mode onset is necessary to overcome the FLR collisional stabilization when using the GK Coulomb collision operator. We reveal the importance of FLR terms in the Coulomb collision operator by demonstrating that neglecting these FLR terms destabilizes a short wavelength branch of the ITG mode, peaking near $k_\perp \rho _s \sim 1.5$. These observations are also found when using the other collision operator models considered in this work. In addition to these findings, the main outcome of the systematic comparison between the collision operators considered in this work is that the GK Sugama slightly underestimates the ITG growth rate compared to the GK Coulomb operator and that the GK pitch-angle scattering operator systematically yields a larger ITG linear growth because of the absence of energy diffusion in the latter. The largest deviations with respect to the GK Coulomb collision operator are observed when the GK Dougherty operator is used. In general, we observe that energy diffusion in the collision operators is important at high collisionality and in the damping of ITG at small scales. Finally, the numerical efficiency of the gyromoment approach is highlighted by performing convergence studies and show that the number of the gyromoment can be reduced as the collisionality increases. In particular, the ability of the presented method to naturally reduce its dimensionality (i.e. the number of gyromoments) without compromising the accuracy of kinetic effects and the details of advanced collision operators makes it an ideal framework to develop high-fidelity simulations of the boundary region, a possible advantage over present continuum GK codes.

Complementary to previous works (see, e.g. Jorge et al. Reference Jorge, Ricci and Loureiro2017, Reference Jorge, Ricci and Loureiro2018; Frei et al. Reference Frei, Jorge and Ricci2020, Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021), the present study of the ITG mode at arbitrary collisionality shows the ability of the gyromoment approach to capture the relevant physics in the kinetic and collisional regimes, with a fidelity cost that depends on the collisionality. The present study can be extended to more realistic conditions, including variation along the parallel direction, electron dynamics, electromagnetic effects and nonlinear physics. In particular, nonlinear simulations using the gyromoment approach are required to further investigate and to extrapolate the effects of collision operator models, FLR terms and convergence reported in the present work to, for example, ITG driven turbulence. For instance, the role of FLR collisional damping by the GK collision operators of the ITG (see, e.g. figure 11) is a mechanism that is expected to affect the level of the saturated turbulent fluxes. Additionally, the choice of the collision operator models (e.g. Coulomb, Dougherty and pitch-angle operators) used in the simulations can play an important role in the nonlinear saturated state, especially in the high-collisional regime where their analytical details can no longer be ignored. Nonlinear simulations are being carried out and will be subject to a future publication. Finally, the $\delta f$ formulation adopted in this work might not be valid for the study of the turbulent dynamics in the boundary region, where the large amplitude fluctuations require a full-F gyromoment hierarchy equation, a subject of future work.

Acknowledgements

The authors acknowledge helpful discussions with S. Brunner, P. Donnel and M. Held. The simulations presented herein were carried out on the CINECA Marconi supercomputer under the TSVVT421 project.

Editor A. Schekochihin thanks the referees for their advice in evaluating this article.

Funding

This work has been carried out within the framework of the EUROfusion Consortium, funded by the European Union via the Euratom Research and Training Programme (grant agreement no. 101052200 – EUROfusion). Views and opinions expressed are, however, those of the author(s) only and do not necessarily reflect those of the European Union or the European Commission. Neither the European Union nor the European Commission can be held responsible for them. This work was supported in part by the Swiss National Science Foundation.

Declaration of interests

The authors report no conflict of interest.

Appendix A. Pitch-angle scattering collision operator

The linearized pitch-angle scattering operator model is the sum of the pitch-angle scattering operator for the test component coupled and an ad hoc field component such that the momentum conservation is satisfied. This operator neglects energy diffusion, contrary to the Sugama and Coulomb collision operators, but conserves particle, momentum and energy. The test component of the GK pitch-angle scattering operator is given by (Helander & Sigmar Reference Helander and Sigmar2002)

(A1)\begin{equation} \mathcal{C}^{T} = \nu^{D}(v) \left\{- \mathcal{L}^{2} h - \frac{b^{2}}{2 v_T^{2}} (2 v_\parallel^{2} + v_\perp^{2}) h \right\}, \end{equation}

and the field component is defined by

(A2)\begin{equation} \mathcal{C}^{F} = 2 F_M\nu^{D}(v) \left( v_\perp {\rm J}_1 \bar{U}_\perp^{D}[h] + v_\parallel {\rm J}_0 \bar{U}_\parallel^{D}[h] \right), \end{equation}

where we introduce

(A3a,b)\begin{equation} \bar{U}_\perp^{D}[h] = \frac{3\int \,{\rm d} \boldsymbol{v} v_\perp {\rm J}_1 \nu^{D}(v) h}{\int \,{\rm d} \boldsymbol{v} v^{2} \nu^{D}(v) F_M}, \quad \bar{U}_\parallel^{D}[h] = \frac{3 \int \,{\rm d} \boldsymbol{v} v_\parallel {\rm J}_0 \nu^{D}(v) h}{\int \,{\rm d} \boldsymbol{v} v^{2} \nu^{D}(v) F_M} \end{equation}

with $h = {g} + q {\rm J}_0 \phi F_M / \tau$ the non-adiabatic part of the perturbed gyrocentre distribution function ${g}$. In (A1), $\mathcal {L}^{2}$ is the spherical angular operator defined by $v^{2} \boldsymbol {\nabla }_{\boldsymbol v}^{2} f= \partial _{v}(v^{2} \partial _v f) - \mathcal {L}^{2} f$ and $\nu ^{D}(v) = \nu [ \mathrm {erf}(s) -( \mathrm {erf}(s) - s \,\mathrm {erf}'(s)) / 2 s^{2}]/ s^{3}$ (with $s= v / v_T$) is the pitch-angle scattering (or deflection) velocity-dependent frequency.

Following Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021), we now project (A1) and (A2) onto the Hermite–Laguerre basis in terms of the non-adiabatic moments, $n^{pj}$, defined by

(A4)\begin{equation} n^{pj} = \frac{1}{N} \int \,{\rm d} \mu \,{\rm d} {v_{{\parallel}}} \,{\rm d} \theta\frac{B}{m} h \frac{H_p(s_{{\parallel}}) L_j(x)}{\sqrt{2^{p} p!}}, \end{equation}

and expressed in terms of the gyromoments, $N^{pj}$, by $n^{pj} = N^{pj} + q \mathcal {K}_{j} \phi \delta _p^{0} / \tau$. The gyromoment expansion of the test and field components, (A1) and (A2), are given by

(A5)\begin{equation} \mathcal{C}^{Tlk} = \mathcal{C}^{Dlk} + \bar{R}^{Dlk} \end{equation}

and

(A6)\begin{equation} \mathcal{C}^{Flk} = \bar{R}_\parallel^{Dlk} + \bar{R}_{{\perp}}^{Dlk}, \end{equation}

respectively. In (A5) and (A6), we introduce

(A7)\begin{equation} \mathcal{C}^{Dlk} ={-} \sum_{j=0}^{\infty} \sum_{p=0}^{l+2k} \sum_{h=0}^{k + \lfloor l/2\rfloor} \sum_{g=0}^{p+2j}\sum_{t=0} ^{j + \lfloor l/2\rfloor} \sum_{d=0}^{h} \frac{2^{p} (p!)^{2}\sqrt{2^{g} g!}}{(2p)!} \frac{ \left(T^{{-}1} \right)^{ph}_{lk} }{\sigma_p^{j} \sqrt{2^{l} l!}} \frac{ p(p+1)}{(2p +1)}T _{pj}^{gt}L_{hd}^{p} \bar{\nu}_{*}^{Dpjd}n^{gt}, \end{equation}

with

(A8)\begin{equation} \bar{\nu}_{*}^{Dpjd} = \frac{4 }{\sqrt{\rm \pi}} \sum_{j_1=0}^{j} L_{jj_1}^{p} \bar{\nu}^{Dp+d+j_1+1} \end{equation}

and

(A9a)\begin{gather} \bar{R}^{D lk} ={-} \frac{b^{2}}{\sqrt{\rm \pi}} \sum_{p,j} \sum_{g=0}^{l + 2 k} \sum_{h=0}^{k + \lfloor l/2\rfloor} \sum_{r=0}^{ p + 2j} \sum_{s=0}^{j + \lfloor p/2\rfloor} \sum_{s_1 =0}^{s}\sum_{h_1 =0}^{h} \frac{\left( T^{{-}1} \right)^{gh}_{lk}}{\sqrt{2^{l} l!}}\frac{\left( T^{{-}1} \right)^{rs}_{pj}}{\sqrt{2^{p} p!}} L^{r}_{ss_1} L^{g}_{hh_1} \mathcal{\alpha}_{{\perp}}^{grs_1 h_1} n^{pj}, \end{gather}
(A9b)\begin{gather}\bar{R}_\perp^{lk} = \bar{U}_\perp[h] \frac{4 }{\sqrt{\rm \pi}} \sum_{n=0}^{\infty} \sum_{s =0}^{n + k + 1} \sum_{h =0}^{s + \lfloor l/2\rfloor} \sum_{h_1=0}^{h} L_{hh_1}^{0} \frac{ d_{nks}^{1} b \mathcal{K}_{n}}{ (n+1)} \frac{\left( T^{{-}1} \right)_{ls}^{gh}}{\sqrt{2^{l} l!}} \bar{\nu}^{Dh_1 +1}, \end{gather}
(A9c)\begin{gather}\bar{R}_\parallel^{lk} = \bar{U}_\parallel[h] \frac{8 }{3 \sqrt{\rm \pi}}\sum_{n=0}^{\infty} \sum_{s =0}^{n+k} \sum_{g =0}^{l + 2s} \sum_{h =0}^{s + \lfloor l/2\rfloor} \sum_{h_1 =0}^{h} d_{nks}^{0} \mathcal{K}_{n} \frac{L_{hh_1}^{1} \left( T^{{-}1} \right)_{ls}^{gh}}{\sqrt{2^{l} l!}} \bar{\nu}^{Dh_1 + 2} [l > 0 \cup s > 0], \end{gather}

where $[\boldsymbol {\cdot }]$ denotes the Iverson bracket. In (A9), we define

(A10)\begin{equation} \mathcal{\alpha}_{{\perp}}^{grs_1 h_1} = \frac{2\delta_g^{r}}{(2g+1)} \bar{\nu}^{D g + 2 + s_1 + h_1 }+ d_{g}^{r} \bar{\nu}^{D(g+r)/2 + 2 + s_1 + h_1 } , \end{equation}

with the numerical coefficient

(A11)\begin{equation} d_g^{r} = \begin{cases} \dfrac{2 (r + 1)(r+2)}{(2 r+1)(2r +3)(2r+5)}, & \text{for } g = r+2\\ \dfrac{2(2r^{2} + 2r -1)}{(2r-1)(2r+1)(2r+3)}, & \text{for } g=r \\ \dfrac{2 r(r-1)}{(2r -3)(2r -1)(2r +1)}, & \text{for } g=r-2 \end{cases} \end{equation}

and

(A12a)\begin{gather} \bar{U}_\parallel[h] = \frac{1}{ \bar{\nu}^{D2}} \sum_{p,j} \sum_{n=0}^{\infty} \sum_{s =0}^{j + n }\sum_{g =0}^{p + 2s} \sum_{h =0}^{s + \lfloor p/2\rfloor} \sum_{h_1=0}^{h} d_{njs}^{0} \mathcal{K}_{n} \frac{ L_{hh_1}^{1} \left( T^{{-}1} \right)_{ps}^{gh} }{\sqrt{2^{p} p!}} \bar{\nu}^{Dh_1 +2} n^{pj} [p > 0 \cup s > 1], \end{gather}
(A12b)\begin{gather}\bar{U}_\perp[h] = \frac{3 }{2 \bar{\nu}^{D2}} \sum_{p,j} \sum_{n=0}^{\infty} \sum_{s =0}^{n + j + 1} \sum_{h =0}^{s + \lfloor p/2\rfloor} \sum_{h_1=0}^{h} L_{hh_1}^{0}\frac{b \mathcal{K}_{n}}{ (n+1)} \frac{ d_{njs}^{1} \left( T^{{-}1} \right)_{ps}^{0h}}{\sqrt{2^{p} p !}} \bar{\nu}^{Dh_1 +1 } n^{pj}. \end{gather}

In (A9) and (A12), we introduce the speed integrated deflection frequency by $\bar {\nu }^{Dk} = \int {\rm d s} {\rm e}^{- s^{2}} s^{2k} \nu ^{D}(v)$. Using the definitions of $\nu ^{D}(v)$, we derive

(A13)\begin{equation} \bar{\nu}^{Dk} = \frac{1}{2} \nu (2 E^{k-2} - E^{k-3} + {\rm e}^{k-2}), \end{equation}

where the definitions of $E^{k}$ and ${\rm e}^{k}$ can be found in Frei et al. (Reference Frei, Ball, Hoffmann, Jorge, Ricci and Stenger2021). We remark that the definitions of the numerical coefficients $( T^{-1} )_{ps}^{gh}$ and $T _{pj}^{gt}$, appearing in, e.g. (A7) and (A12), are reported in Jorge et al. (Reference Jorge, Ricci and Loureiro2017).

The DK limit of the pitch-angle scattering operator is obtained from (A5) and (A6) by taking the zero gyroradius limit, $b \to 0$. Noticing that $n^{pj} \simeq N^{pj}$, we derive the DK test and the field components given by

(A14)\begin{gather} \mathcal{C}_{}^{Tlk} = \mathcal{C}^{Dlk}, \end{gather}
(A15)\begin{gather} \mathcal{C}_{}^{Flk} = U_\parallel[{g}] \frac{8 }{3 \sqrt{\rm \pi}} \sum_{g =0}^{l + 2k} \sum_{h =0}^{k + \lfloor l/2\rfloor} \sum_{h_1 =0}^{h} \frac{L_{hh_1}^{1} \left( T^{{-}1} \right)_{lk}^{gh}}{\sqrt{2^{l} l!}} \bar{\nu}^{Dh_1 + 2} [ l > 0 \cup k > 0 ], \end{gather}

respectively, where $C^{Dlk}$ is given in (A7) and with

(A16)\begin{equation} U_\parallel^{D}[{g}] = \frac{1}{ \bar{\nu}^{D2}} \sum_{p,j} \sum_{g =0}^{p + 2j} \sum_{h =0}^{j + \lfloor p/2\rfloor} \sum_{h_1=0}^{h} N^{pj} \frac{ L_{hh_1}^{1} \left( T^{{-}1} \right)_{pj}^{gh} }{\sqrt{2^{p} p!}} \bar{\nu}^{Dh_1 +2} [ p > 0 \cup k > j ]. \end{equation}

Appendix B. HSC collision operator

We consider the gyromoment expansion of the HSC collision operator (Hirshman, Sigmar & Clarke Reference Hirshman, Sigmar and Clarke1976), also referred to as the zeroth-order Hirshman–Sigmar collision operator (Hirshman & Sigmar Reference Hirshman and Sigmar1976b). The HSC collision operator model has been widely used for neoclassical transport calculations (Hirshman et al. Reference Hirshman, Sigmar and Clarke1976; Belli & Candy Reference Belli and Candy2008). This operator is only considered in the DK limit. Additionally, it conserves particle, momentum and energy. While the test component of the HSC operator is the pitch-angle scattering operator, i.e.

(B1)\begin{equation} \mathcal{C}^{T} ={-} \nu^{D}(v) \mathcal{L}^{2} {g}, \end{equation}

the field component of the HSC operator is defined by

(B2)\begin{equation} \mathcal{C}^{F} = v_\parallel \left[ \nu^{S}(v) \frac{2r[{g}]}{v_T^{2}} + \frac{2 u[{g}]}{v^{2}} \left( \nu^{D}(v)- \nu^{S}(v) \right)\right]F_{M}, \end{equation}

where

(B3)\begin{gather} r[{g}] = \dfrac{3 \int \,{\rm d} \boldsymbol{v} v_\parallel \nu^{S}(v) {g}}{\int \,{\rm d} \boldsymbol{v} s^{2} \nu^{S}(v) F_{M}}, \end{gather}
(B4)\begin{gather}u[{g}] = \frac{3}{4 {\rm \pi}}\int \,{\rm d} \varOmega v_\parallel \frac{{g}}{F_{M}}, \end{gather}

with $\int \,{\rm d} \varOmega = \int _0^{2 {\rm \pi}} \,{\rm d} \theta \int _{-1}^{1} \,{\rm d}\xi$ the integral over the velocity-space solid angle, and $\nu ^{S} (v) = 2 \nu ( \mathrm {erf}(s) - s \,\mathrm {erf}'(s) ) /s^{3}$.

We project the HSC collision operator, given in (B2), onto the Hermite–Laguerre basis. As the Hermite–Laguerre projection of $- \nu ^{D}(v) \mathcal {L}^{2} g$ is given by (A7), we focus on the projection of (B2). Performing the Hermite–Laguerre projection yields

(B5)\begin{equation} \mathcal{C}^{Flk} = R^{Slk}_1 + R^{Dlk}, \end{equation}

where we introduce

(B6a)\begin{align} R_{1}^{Slk} & = \sum_{g = 0}^{l + 2k} \sum_{h = 0}^{k + \lfloor l/2 \rfloor} \sum_{h_1=0}^{h} L_{hh_1}^{1}\frac{\left( T ^{{-}1} \right)^{gh}_{lk}}{\sqrt{2^{l}l!}} \frac{8}{3 \sqrt{\rm \pi}} \frac{r[{g}]}{v_T} \bar{\nu}^{Sh_1 +2} [ l > 0 \cup k > 0], \end{align}
(B6b)\begin{align} R^{SDlk} & = \frac{8}{3 \sqrt{\rm \pi}} \sum_{p,j} \sum_{g=0}^{p + 2j} \sum_{h=0}^{j + \lfloor p/2\rfloor} \sum_{s=0}^{l + 2k} \sum_{t=0}^{k + \rfloor l /2 \lfloor} \sum_{h_1 =0}^{h} \sum_{t_1 = 0}^{t} \frac{\left( T ^{{-}1} \right)^{st}_{lk}}{\sqrt{2^{l} l! }} \frac{\left( T ^{{-}1} \right)^{gh}_{pj}}{\sqrt{2^{p} p! }} \nonumber\\ & \quad \times L_{tt_1}^{1} L_{hh_1}^{1} \left( \bar{\nu}^{D2 + t_1 + h_1} - \bar{\nu}^{S2 + t_1 + h_1}\right) N^{pj} [ p > 0 \cup j > 0] [ l > 0 \cup k > 0] \end{align}

and

(B7)\begin{equation} \frac{ r[{g}] }{v_T} = \frac{1}{\bar{\nu}^{S2}} \sum_{p,j} \sum_{s=0}^{p + 2j} \sum_{t=0}^{j + \lfloor p/2\rfloor} \sum_{t_1=0}^{t} \left( T ^{{-}1} \right)^{st}_{pj} L_{tt_1}^{1} \frac{N^{pj}}{\sqrt{2^{p} p!}} \bar{\nu}^{St_1 +2} [ p > 0 \cup j > 0 ]. \end{equation}

In, (B6a) and (B7), we introduce the speed integrated frequency by $\bar {\nu }^{Sk} = \int d s {\rm e}^{- s^{2}} s^{2k} \nu ^{S}(v)$. Using the definitions of $\nu ^{S}(v)$, we derive

(B8)\begin{equation} \bar{\nu}^{Sk} = 2 \nu (E^{k-2} - {\rm e}^{k-1}). \end{equation}

Appendix C. Dougherty collision operator

The linearized GK Dougherty collision operator is given by the test and field components (Dougherty Reference Dougherty1964)

(C1)\begin{equation} \mathcal{C}^{T} = \nu \left[ \frac{T}{m} \frac{\partial^{2}}{\partial v_\parallel^{2}} h + 2 \frac{T}{B} \frac{\partial}{\partial \mu} \left( \mu \frac{\partial}{\partial \mu} h \right)- \frac{1}{2} b^{2} h + 3 h+ v_\parallel \frac{\partial}{\partial v_\parallel } h + 2 \mu \frac{\partial}{\partial \mu} h \right] \end{equation}

and

(C2)\begin{equation} \mathcal{C}^{F} = \nu \left[ 2 {\rm J}_0 T[h]\left( x + s_{{\parallel} }^{2} - \frac{3}{2} \right) + 2 u_{{\parallel}}[h] s_\parallel {\rm J}_0(b) + {\rm J}_1(b) 2 u_\perp[h] \sqrt{x} \right] F_{M}, \end{equation}

respectively. Here, the particle fluid quantities are defined by $T[h] = ( T_{\parallel }[h]+ 2 T_{\perp }[h])/3 - n[h]$ where $T_{\parallel } [h]= \int \,{\rm d}^{3} \boldsymbol {v} m v_\parallel ^{2} {\rm J}_0 h / T$, $T_{\perp }[h] = \int \,{\rm d}^{3} \boldsymbol {v} {\rm J}_0 \mu B h/T$, $n[h] = \int \,{\rm d}^{3} \boldsymbol {v} {\rm J}_0 h/N$, and $u_\parallel [h] = \int \,{\rm d}^{3} \boldsymbol {v} {\rm J}_0 s_\parallel h$, $u_\perp [h] = \int \,{\rm d}^{3} \boldsymbol {v} \sqrt {x} {\rm J}_1 h$.

Projecting the GK Dougherty collision operator, (C1) and (C2), onto the Hermite–Laguerre basis yields,

(C3)\begin{equation} \mathcal{C}^{Tlk} ={-}\nu \left( 2k + l+ \frac{b^{2}}{2} \right) n^{lk } \end{equation}

and

(C4)\begin{align} \mathcal{C}^{Flk} & = \nu ( T[h] [2 \delta_{l}^{0} ( 2 k \mathcal{K}_{k}-(k+1) \mathcal{K}_{k+1} - k \mathcal{K}_{k-1}) + \mathcal{K}_{k} \sqrt{2} \delta_l^{2} ]\nonumber\\ & \quad + u_{{\parallel} }[h] \mathcal{K}_{k} \delta_l^{1} + u_{{\perp}}[h] b ( \mathcal{K}_{k} - \mathcal{K}_{k-1} )\delta_l^{0} \vphantom{\sqrt{2} \delta_l^{2}}) , \end{align}

where the particle fluid quantities are expressed as

(C5a)\begin{gather} n[h] = \sum_{j=0}^{\infty} \mathcal{K}_{j} n^{0j}, \end{gather}
(C5b)\begin{gather}T_{{\parallel} }[h]= \sum_{j=0}^{\infty} \mathcal{K}_{j} [\sqrt{2} n^{2j} + n^{0j}], \end{gather}
(C5c)\begin{gather}T_{{\perp} } [h] = \sum_{j=0}^{\infty} n^{0j} [(2 j+1) \mathcal{K}_{j} - (j+1) \mathcal{K}_{j+1} - j \mathcal{K}_{j-1}], \end{gather}
(C5d)\begin{gather}u_{{\parallel} }[h] = \sum_{j=0}^{\infty} \mathcal{K}_{j} n^{1j}, \end{gather}
(C5e)\begin{gather}u_{{\perp}} [h]= b \sum_{j=0}^{\infty} \frac{1}{2 } [\mathcal{K}_{j} - \mathcal{K}_{j-1}] n^{0j} . \end{gather}

We remark that the GK Dougherty collision operator, given in (C1) and (C2), is used in the previous Hermite–Laguerre GK formulation (Mandell et al. Reference Mandell, Dorland and Landreman2018), because of its sparse representation on this basis. In fact, $\mathcal {C}^{Tlk} = -\nu ( 2k + l+ b^{2} / 2 ) n^{lk}$ and $\mathcal {C}^{Flk}=0$ for $l > 2$.

The DK Dougherty collision operator is obtained in the zero gyroradius limit of (C3) and (C4), and is given by

(C6)\begin{equation} \mathcal{C}^{lk} = \nu \left[ - (2 k + l ) N^{lk} - ( \sqrt{2} N^{20} - 2 N^{01})\left( \frac{2}{3} \delta_l^{0} \delta_k^{1}- \frac{\sqrt{2}}{3} \delta_l^{2} \delta_k^{0} \right) + N^{10}\delta_l^{1} \delta_k^{0} \right]. \end{equation}

We remark that the DK Dougherty operator, given in (C6), is equivalent to the one used in Jorge et al. (Reference Jorge, Ricci and Loureiro2018).

Appendix D. Collisionless gyromoment expressions

This appendix reports on the derivations of the collisionless gyromoment expressions, defined in (3.8). We start from (3.1) that we multiply by the Hermite–Laguerre basis yielding

(D1)\begin{equation} \frac{N^{pj}}{\phi } = \sum_{l=1}^{3} \hat{N}_{l}^{pj}, \end{equation}

where we introduce the velocity integrals,

(D2a)\begin{gather} \hat{N}_{1}^{pj} ={-} \frac{q}{\tau}\int_0^{\infty} {{\rm d} x} \int_{- \infty}^{\infty} \frac{{\rm d} s_\parallel}{\sqrt{\rm \pi}} \frac{H_p(s_\parallel ) L_j(x)}{\sqrt{2^{p} p!}} {\rm J}_0(b\sqrt{x}) {\rm e}^{ - s_\parallel^{2} -x}, \end{gather}
(D2b)\begin{gather}\hat{N}_{2}^{pj} = \frac{q}{\tau} \int_0^{\infty} {{\rm d} x} \int_{- \infty}^{\infty} \frac{{\rm d} s_\parallel}{\sqrt{\rm \pi}} \frac{\omega {\rm J}_0(b \sqrt{x}) \exp({-x - s_\parallel^{2}})} { \omega - \omega_{\boldsymbol{\nabla} B}- \sqrt{2 \tau} k_\parallel s_\parallel } \frac{H_p(s_\parallel ) L_j(x)}{\sqrt{2^{p} p!}}, \end{gather}
(D2c)\begin{gather}\hat{N}_{3}^{pj} ={-} \int_0^{\infty} {{\rm d} x} \int_{- \infty}^{\infty} \frac{{\rm d} s_\parallel}{\sqrt{\rm \pi}} \omega_{T}^{*}\frac{ {\rm J}_0(b \sqrt{x}) {\rm e}^{{-}x - s_\parallel^{2}}} { \omega - \omega_{\boldsymbol{\nabla} B }- \sqrt{2 \tau} k_\parallel s_\parallel } \frac{H_p(s_\parallel ) L_j(x)}{\sqrt{2^{p} p!}}. \end{gather}

We now perform analytically the velocity integrals that contain the resonant term proportional to $1 / ( \omega - \omega _{\boldsymbol {\nabla } B }- \sqrt {2 \tau } k_\parallel s_\parallel )$ using the transformation given (3.4) valid for unstable modes. The adiabatic part of the collisionless gyromoment response, $\hat {N}_1^{pj}$, is deduced from the orthogonality relations, (2.10) and (2.12), such that

(D3)\begin{equation} \hat{N}_1^{pj} ={-} \frac{q}{\tau} \mathcal{K}_{j} \delta_p^{0}. \end{equation}

Expanding the Bessel function in terms of Laguerre polynomials using (2.12), the term $\hat {N}^{pj}_2$, defined in (D2b), can be written as

(D4)\begin{align} \hat{N}_{2}^{pj}& ={-}{\rm i} \frac{q}{\tau} \sum_{n=0}^{\infty} \mathcal{K}_{n} \int_0^{\infty} \,{\rm d} \sigma \omega {\rm e}^{{\rm i} \sigma \omega} \int_{- \infty}^{\infty} \frac{{\rm d} s_\parallel}{\sqrt{\rm \pi}} \exp({- (1 + {\rm i} 2 \sigma \alpha )s_\parallel^{2}- {\rm i} \sigma z_\parallel s_\parallel})\frac{H_p(s_\parallel ) }{\sqrt{2^{p} p!}} \nonumber\\ & \quad \times \int_0^{\infty} \,{\rm d} xL_n(x)L_j(x) \exp({-(1 + {\rm i} \sigma \alpha)x}), \end{align}

where $\alpha = \tau \omega _{ B}/q$. The $x$-integration can be performed analytically using the identity (Gradshteyn & Ryzhik Reference Gradshteyn and Ryzhik2014)

(D5)\begin{equation} \int_0^{\infty} {{\rm d} x} {\rm e}^{- \beta x} L_n( x) L_m( x) = \frac{(m + n ) !}{m! n!} \frac{(\beta - 1)^{n+m}}{\beta^{m + n + 1}} \mathcal{F}^{m}_n\left[\frac{\beta(\beta - 2)}{(\beta - 1)^{2}}\right], \end{equation}

valid for $Re \beta > 0$ and where $\mathcal {F}^{a}_b[z] = F[-a,-b;-a-b,z]$ is the Gauss hypergeometric function. We remark that, in the slab limit obtain by setting $\alpha =0$ in (D4) (and thus $\beta = 1$ in (D5)), (D5) reduces to the orthogonality relation between Laguerre polynomials, (2.11). The $s_\parallel$-integration in (D4) is evaluated by expanding $H_p(s_\parallel )$ such that

(D6)\begin{equation} H_p (s_\parallel) = \sum_{p_1=0}^{\lfloor p/2\rfloor} \frac{({-}1)^{p_1}2^{p - 2p_1}p! }{p_1!(p -2 p_1)!} s_\parallel^{p - 2p_1}. \end{equation}

Using the above identities, we derive

(D7)\begin{equation} \hat{N}_{2}^{pj} ={-} \frac{{\rm i} q \omega}{\tau} \sum_{n=0}^{\infty} \frac{ \mathcal{K}_{n} }{\sqrt{2^{p} p!} }\mathcal{I}_{n}^{pj} , \end{equation}

where we introduce the real definite integral

(D8)\begin{align} \mathcal{I}_{n}^{pj} & = \sum_{p_1=0}^{\lfloor p/2\rfloor} \frac{({-}1)^{p_1 }2^{p - 2p_1 -1}p! }{p_1!(p -2 p_1)!} \int_0^{\infty} \frac{{\rm d} \sigma}{\sqrt{\rm \pi}} {\rm e}^{{\rm i} \sigma \omega} (1 + 2 {\rm i} \sigma \alpha) ^{p_1-{p}/{2}-{1}/{2}} \nonumber\\ & \quad \times \frac{(n+j ) !}{n! j!} \frac{({\rm i} \alpha \sigma)^{n+j}}{(1 + {\rm i} \alpha \sigma)^{j + n + 1}}\mathcal{F}^{j}_n\left[\frac{(1 + {\rm i} \alpha \sigma)( {\rm i} \sigma \alpha-1 )}{( {\rm i} \alpha \sigma)^{2}}\right] \nonumber\\ & \quad \times \left\{\left(1+({-}1)^{p - 2p_1}\right) \left( \frac{p}{2} - p_1 - \frac{1}{2}\right)! \mathcal{M}^{{p/2 - p_1 +1/2}}_{1/2}\left[- \frac{z_\parallel^{2} \sigma^{2}}{4(1 + 2 \alpha {\rm i} \sigma )} \right] \right. \nonumber\\ & \quad \left. -2 {\rm i} \frac{z_\parallel \sigma \left(1+({-}1)^{p - 2 p_1+1}\right) }{2 \sqrt{(1 + 2 {\rm i} \sigma \alpha)}} \left(\frac{p}{2} -p_1 \right)! \mathcal{M}^{{p}/{2}-p_1+1}_{3/2} \left[- \frac{z_\parallel^{2} \sigma^{2}}{4(1 + 2 \alpha {\rm i} \sigma )} \right]\right\}, \end{align}

where $\mathcal {M}^{a}_b[z] = M(a;b;z)$ is the Kummer confluent hypergeometric function, which stems from the $s_\parallel$-integration.

Finally, using the recursive properties of Hermite and Laguerre polynomials, we deduce that

(D9)\begin{align} \hat{N}_{3}^{pj}& = \frac{{\rm i} k_\perp}{\sqrt{2^{p} p!}}\sum_{n=0}^{\infty}\mathcal{K}_{n} \left\{\vphantom{\frac{1}{2}} \mathcal{I}_n^{pj} + \eta \left[\vphantom{\frac{1}{2}} (2j+1) \mathcal{I}_n^{pj} - j \mathcal{I}_n^{pj-1}- (j+1) \mathcal{I}_n^{pj+1} \right.\right. \nonumber\\ & \quad \left. \left. + \frac{\mathcal{I}_{n}^{p+2j}}{4} + \left(p + \frac{1}{2}\right) \mathcal{I}_n^{pj} + p(p-1) \mathcal{I}_n^{p-2j}- \frac{3}{2}\mathcal{I}_n^{pj}\right] \right\}. \end{align}

With (D3), (D7) and (D9), the closed semianalytical expressions of the collisionless gyromoment response, (3.8), can be evaluated. We remark that a simpler expression for $\mathcal {I}_{n}^{pj}$ can be obtained in the purely toroidal limit, i.e. $k_\parallel =0$, neglecting the resonance effects due to the parallel streaming. Hence, from (D8), one obtains

(D10)\begin{align} \mathcal{I}_{n}^{pj} & = \sum_{p_1=0}^{\lfloor p/2\rfloor} \frac{({-}1)^{p_1 }2^{p - 2p_1}p! }{p_1!(p -2 p_1)!}\frac{(n+j ) !}{n! j! \sqrt{\rm \pi}} \left( \frac{p}{2} - p_1 - \frac{1}{2}\right)! \nonumber\\ & \quad \times \int_0^{\infty} \,{\rm d} \sigma {\rm e}^{{\rm i} \sigma \omega} (1 + 2 {\rm i} \sigma \alpha) ^{p_1-{p}/{2}-{1}/{2}} \frac{({\rm i} \alpha \sigma)^{n+j}}{(1 + {\rm i} \alpha \sigma)^{j + n + 1}} \mathcal{F}^{j}_n \left[ \frac{1 + \alpha^{2} \sigma^{2}}{ \alpha^{2} \sigma^{2}} \right], \end{align}

if $p$ is even and $0$ otherwise.

Appendix E. FLR closures

We explore possible FLR closures for the terms, proportional to $\mathcal {K}_{0}$, $\mathcal {K}_{1}$ and $\mathcal {K}_{2}$ appearing in the 6GM and 4GM using the Padé approximation technique. In order to deduce Padé approximation of the kernel functions, we rewrite the exact form of $\mathcal {K}_{n}$ defined in (2.13) as

(E1)\begin{equation} \mathcal{K}_{j} = \left.\frac{({-}1)^{j}}{j!} \frac{\partial^{j}}{\partial \beta^{j}} \varUpsilon( a \beta ) \right\rvert_{\beta =1}, \end{equation}

where we introduce the basic FLR operator $\varUpsilon ( a \beta ) = \sum _{l=0}^{\infty } (\beta a)^{l}/ (l !2^{l} ) = {\rm e}^{- \beta a/2}$. From (E1), approximation to $\mathcal {K}_{j}$ can be constructed by specifying the basic FLR operator $\varUpsilon$. We remark that in contrast to the choice $\varUpsilon (a \beta ) = {\rm e}^{- \beta a/2}$ introduced in Brizard (Reference Brizard1992), Dorland & Hammett (Reference Dorland and Hammett1993) proposed a modified basic FLR operator such as $\varUpsilon ( a \beta ) = \sqrt { \varGamma _0( a \beta )}$. The latter approximation is motivated by the fact that basic FLR operator, $\varUpsilon ( a \beta ) = {\rm e}^{- \beta a/2}$, yields a large damping at short wavelength with respect to Dorland & Hammett (Reference Dorland and Hammett1993) when applied to, for example, the collisionless ITG case. We focus on the basic FLR operator $\varUpsilon ( a \beta ) = {\rm e}^{- \beta a/2}$, as little difference with respect to $\sqrt { \varGamma _0(a \beta )}$ is expected in the high-collisional limits as collisions damp short wavelength modes (see figure 8).

We now aim to construct a Padé approximant, of order $(p,q)$, of the basic FLR operator, that we denote by ${}_q^{p}\check {\varUpsilon }(x) = R^{(p)}(x) / S^{(q)}(x)$ where $R^{(p)}(x)$ and $S^{(q)}(x)$ are polynomials in $x$ of order $p$ and $q$, respectively. The Padé approximant ${}_q^{p}\check {\varUpsilon }(x)$ satisfies $\varUpsilon ^{(p+q)}(0) = {}_q^{p}\check {\varUpsilon }^{(p+q)}(0)$. Once ${}_q^{p}\check {\varUpsilon }^{(p+q)}(0)$ is specified, the self-consistent Padé approximant, of order $(p,q)$, of the $j$th-order kernel function, denoted by ${}_q^{p}\check {\mathcal {K}_{j}} \simeq \mathcal {K}_{j}$, is expressed by

(E2)\begin{equation} {}_q^{p}\check{\mathcal{K}_{j}} = \left.\frac{({-}1)^{j}}{j!}\frac{\partial^{j}}{\partial \beta^{j}} {}_q^{p}\check{\varUpsilon}(\beta a) \right\rvert_{\beta =1}. \end{equation}

We consider the Padé approximant of order $(p,q) =(1,2)$ and $(1,4)$ for $j = 0,1,2$ to approximate $\mathcal {K}_{0}$, $\mathcal {K}_{1}$ and $\mathcal {K}_{2}$ appearing in the FLR terms contained in both the 6GM and 4GM models (see (4.2)). Equation (E2) yields the $(p,q) =(1,2)$ self-consistent Padé approximants

(E3a)\begin{gather} {}_2^{1}\check{\mathcal{K}_{0}} = \frac{1}{1 + b^{2}/4}, \end{gather}
(E3b)\begin{gather}{}_2^{1}\check{\mathcal{K}_{1}} = \frac{b^{2}}{4} \frac{1}{(1 + b^{2}/4)^{2}} , \end{gather}
(E3c)\begin{gather}{}_2^{1}\check{\mathcal{K}_{2}} = \frac{b^{4}}{16}\frac{1}{(1 + b^{2}/4)^{3}}, \end{gather}

and the $(p,q) =(1,4)$ Padé approximants

(E4a)\begin{gather} {}_4^{1}\check{\mathcal{K}_{0}} = \frac{1}{1 + b^{2}/4 +b^{4}/32 }, \end{gather}
(E4b)\begin{gather}{}_4^{1}\check{\mathcal{K}_{1}} =\frac{b^{2}/4 + b^{4}/16}{(1 + b^{2}/4 +b^{4}/32 )^{2}} , \end{gather}
(E4c)\begin{gather}{}_4^{1}\check{\mathcal{K}_{2}} = \frac{1}{2} \left[ \frac{2(b^{2}/4 + b^{4}/16)^{2}}{(1 + b^{2}/4 +b^{4}/32 )^{3}} - \frac{b^{4}}{16(1 + b^{2}/4 +b^{4}/32 )^{2}} \right]. \end{gather}

Figure 21(a) displays the lowest-order kernel functions, $\mathcal {K}_{j}$ for $j =0$, $1$, $2$, and the Padé approximants, given in (E3) and (E4), respectively. As observed, the $(p,q) = (1,2)$ Padé approximants decay slower in the larger $b$ limit, compared with the $(p,q) = (1,4)$ models, while they both agree with the kernel functions $\mathcal {K}_{j}$ in the small $b$ limit. Also, the $(p,q) = (1,4)$ case provides a better approximation near the $\mathcal {K}_{1,2}$ maxima. Using the Padé approximants, defined (E3) and (E4), allows us to compute the sITG growth rate using the 4GM. We consider the parameters $k_\parallel = 0.1$, $\eta = 5$ and $\nu = 1$ and show the results in figure 21(b) (similar results are obtained with the 6GM). Notice that, while the Padé approximants yield a good approximation when $k_\perp \lesssim 0.4$, the $(p,q) = (1,4)$ Padé approximant models accurately describe the linear growth $\gamma$ rate near the peak and behave qualitatively well when $k_\perp \gtrsim 1$. The second peak corresponds to the SWITG mode (see § 6.2). This simple application demonstrates the ability of the Padé approximant to accurately model FLR terms. Indeed, in conventional space, the kernel function, $\mathcal {K}_{n}$ introduced in (2.13), defines an infinite linear combination of differential operators in conventional space, i.e.

(E5)\begin{equation} \mathcal{K}_{j} = \sum_{l =0}^{\infty} \frac{({-}1)^{j}}{l!j!} \left(\frac{\rho^2}{4} \varDelta_\perp\right)^{j+l}, \end{equation}

where $\varDelta _\perp = \boldsymbol {\nabla } \boldsymbol {\cdot } \boldsymbol {\nabla }_\perp$ is the perpendicular Laplacian. From a numerical point of view, (E5) is not practical due to the presence of the infinite sum over the $l$ index. A naïve approach is by performing a truncation of the sum at a given order in $b$. But, it would yield an arbitrary loss of accuracy, in particular, for short wavelength mode, such as, e.g. ITG that peak at $k_\perp \sim 0.6 / \rho _s$ (in physical units). A second option is to approximate the functional form of (E5) using a Padé approximant as discussed in this appendix. Such approach can be generalized to describe FLR terms in the GK collision operators and in the gyromoment hierarchy for arbitrary number of gyromoments.

Figure 21. (a) Lowest-order kernel functions $\mathcal {K}_{j}$ (solid line) and the corresponding Padé approximants, ${}_q^{p}\check {\mathcal {K}_{j}}$ for $j =0,1,2$ when $(p,q) = (1,2)$ (dashed line) and $(p,q) = (1,4)$ (dotted line). (b) The ITG linear growth rate $\gamma$ as a function of the perpendicular wavenumber $k_\perp$, obtained by the 4GM with the exact kernel functions $\mathcal {K}_{j}$ (solid line) and with the Padé approximants ${}_q^{p}\check {\mathcal {K}_{j}}$, with $(p,q) = (1,2)$ (dashed line) and $(p,q ) = (1,4)$ (dotted line), respectively.

References

REFERENCES

Abel, I.G., Barnes, M., Cowley, S.C., Dorland, W. & Schekochihin, A.A. 2008 Linearized model Fokker-Planck collision operators for gyrokinetic simulations. I. Theory. Phys. Plasmas 15 (12), 122509.CrossRefGoogle Scholar
Ajay, C.J., Brunner, S. & Ball, J. 2021 Effect of collisions on non-adiabatic electron dynamics in ITG-driven microturbulence. Phys. Plasmas 28 (9), 092303.Google Scholar
Barnes, M., Abel, I.G., Dorland, W., Ernst, D.R., Hammett, G.W., Ricci, P., Rogers, B.N., Schekochihin, A.A. & Tatsuno, T. 2009 Linearized model Fokker-Planck collision operators for gyrokinetic simulations. II. Numerical implementation and tests. Phys. Plasmas 16 (7), 072107.CrossRefGoogle Scholar
Beer, M.A. & Hammett, G.W. 1996 Toroidal gyrofluid equations for simulations of tokamak turbulence. Phys. Plasmas 3 (11), 4046.CrossRefGoogle Scholar
Belli, E.A. & Candy, J. 2008 Kinetic calculation of neoclassical transport including self-consistent electron and impurity dynamics. Plasma Phys. Control. Fusion 50 (9), 095010.CrossRefGoogle Scholar
Belli, E.A. & Candy, J. 2017 Implications of advanced collision operators for gyrokinetic simulation. Plasma Phys. Control. Fusion 59 (4), 045005.CrossRefGoogle Scholar
Bernard, T.N., Shi, E.L., Gentle, K.W., Hakim, A., Hammett, G.W., Stoltzfus-Dueck, T. & Taylor, E.I. 2019 Gyrokinetic continuum simulations of plasma turbulence in the texas helimak. Phys. Plasmas 26 (4), 042301.CrossRefGoogle Scholar
Boyd, J.P. 2001 Chebyshev and Fourier Spectral Methods. Courier Corporation.Google Scholar
Brizard, A. 1992 Nonlinear gyrofluid description of turbulent magnetized plasmas. Phys. Fluids B: Plasma Phys. 4 (5), 1213.CrossRefGoogle Scholar
Catto, P.J. & Ernst, D.R. 2009 Alternate form of model like particle collision operator. Plasma Phys. Control. Fusion 51 (6), 062001.CrossRefGoogle Scholar
Chapman, S. & Cowling, T.G. 1941 The velocity of diffusion in a mixed gas; the second approximation. Proc. R. Soc. Lond. A 179 (977), 159.Google Scholar
Chen, L., Briguglio, S. & Romanelli, F. 1991 The long-wavelength limit of the ion-temperature gradient mode in tokamak plasmas. Phys. Fluids B: Plasma Phys. 3 (3), 611.CrossRefGoogle Scholar
Crandall, P., Jarema, D., Doerk, H., Pan, Q., Merlo, G., Görler, T., Navarro, A.B., Told, D., Maurer, M. & Jenko, F. 2020 Multi-species collisions for delta-f gyrokinetic simulations: implementation and verification with gene. Comput. Phys. Commun. 255, 107360.CrossRefGoogle Scholar
Dorf, M. & Dorr, M. 2022 Modelling of electrostatic ion-scale turbulence in divertor tokamaks with the gyrokinetic code cogent. Contrib. Plasma Phys., e202100162.Google Scholar
Dorland, W. & Hammett, G.W. 1993 Gyrofluid turbulence models with kinetic effects. Phys. Fluids B: Plasma Phys. 5 (3), 812.CrossRefGoogle Scholar
Dougherty, J.P. 1964 Model Fokker–Planck equation for a plasma and its solution. Phys. Fluids 7 (11), 1788.CrossRefGoogle Scholar
Dudson, B.D., Umansky, M.V., Xu, X.Q., Snyder, P.B. & Wilson, H.R. 2009 Bout++: a framework for parallel plasma fluid simulations. Comput. Phys. Commun. 180 (9), 1467.CrossRefGoogle Scholar
Francisquez, M., Juno, J., Hakim, A., Hammett, G.W. & Ernst, D.R. 2021 Improved multispecies Dougherty collisions. https://doi.org/10.48550/arXiv.2109.10381.CrossRefGoogle Scholar
Frassinetti, L., Dunne, M.G., Sheikh, U., Saarelma, S., Roach, C.M., Stefanikova, E., Maggi, C., Horvath, L., Pamela, S., De La Luna, E., et al. 2019 Role of the pedestal position on the pedestal performance in AUG, JET-ILW and TCV and implications for ITER. Nucl. Fusion 59 (7), 076038.CrossRefGoogle Scholar
Frei, B.J., Ball, J., Hoffmann, A.C.D., Jorge, R., Ricci, P. & Stenger, L. 2021 Development of advanced linearized gyrokinetic collision operators using a moment approach. J. Plasma Phys. 87 (5), 905870501.CrossRefGoogle Scholar
Frei, B.J., Jorge, R. & Ricci, P. 2020 A gyrokinetic model for the plasma periphery of tokamak devices. J. Plasma Phys. 86 (2), 905860205.CrossRefGoogle Scholar
Gao, Z., Sanuki, H., Itoh, K. & Dong, J.Q. 2003 Temperature gradient driven short wavelength modes in sheared slab plasmas. Phys. Plasmas 10 (7), 2831.CrossRefGoogle Scholar
Gao, Z., Sanuki, H., Itoh, K. & Dong, J.Q. 2005 Short wavelength ion temperature gradient instability in toroidal plasmas. Phys. Plasmas 12 (2), 022502.CrossRefGoogle Scholar
Garbet, X., Mantica, P., Angioni, C., Asp, E., Baranov, Y., Bourdelle, C., Budny, R., Crisanti, F., Cordey, G., Garzotti, L., et al. 2004 Physics of transport in tokamaks. Plasma Phys. Control. Fusion 46 (12B), B557.CrossRefGoogle Scholar
Giacomin, M., Stenger, L.N. & Ricci, P. 2020 Turbulence and flows in the plasma boundary of snowflake magnetic configurations. Nucl. Fusion 60 (2), 024001.CrossRefGoogle Scholar
Gradshteyn, I.S. & Ryzhik, I.M. 2014 Table of Integrals, Series, and Products. Academic Press.Google Scholar
Gültekin, Ö. & Gürcan, Ö.D. 2018 Stable and unstable roots of ion temperature gradient driven mode using curvature modified plasma dispersion functions. Plasma Phys. Control. Fusion 60 (2), 025021.CrossRefGoogle Scholar
Gültekin, Ö. & Gürcan, Ö.D. 2019 Generalized curvature modified plasma dispersion functions and Dupree renormalization of toroidal ITG. Plasma Phys. Control. Fusion 62 (2), 025018.CrossRefGoogle Scholar
Hahm, T.S. & Tang, W.M. 1989 Properties of ion temperature gradient drift instabilities in h-mode plasmas. Phys. Fluids B: Plasma Phys. 1 (6), 1185.CrossRefGoogle Scholar
Hallatschek, K. & Zeiler, A. 2000 Nonlocal simulation of the transition from ballooning to ion temperature gradient mode turbulence in the tokamak edge. Phys. Plasmas 7 (6), 2554.CrossRefGoogle Scholar
Halpern, F., Ricci, P., Jolliet, S., Loizu, J., Morales, J., Mosetto, A., Musil, F., Riva, F., Tran, T.-H. & Wersal, C. 2016 The GBS code for tokamak scrape-off layer simulations. J. Comput. Phys. 315, 388.CrossRefGoogle Scholar
Hammett, G.W., Dorland, W. & Perkins, F.W. 1992 Fluid models of phase mixing, landau damping, and nonlinear gyrokinetic dynamics. Phys. Fluids B: Plasma Phys. 4 (7), 2052.CrossRefGoogle Scholar
Hazeltine, R.D. & Meiss, J.D. 2003 Plasma Confinement. Courier Corporation.Google Scholar
Helander, P. & Sigmar, D.J. 2002 Collisional Transport in Magnetized Plasmas. Cambridge University Press.Google Scholar
Held, M., Wiesenberger, M. & Kendl, A. 2020 Padé-based arbitrary wavelength polarization closures for full-f gyro-kinetic and fluid models. Nucl. Fusion 60 (6), 066014.CrossRefGoogle Scholar
Hirshman, S.P. & Sigmar, D.J. 1976 a Approximate Fokker–Planck collision operator for transport theory applications. Phys. Fluids 19 (10), 1532.CrossRefGoogle Scholar
Hirshman, S.P. & Sigmar, D.J. 1976 b Approximate Fokker-Planck collision operator for transport theory applications. Phys. Fluids 19 (10), 1532.CrossRefGoogle Scholar
Hirshman, S.P., Sigmar, D.J. & Clarke, J.F. 1976 Neoclassical transport theory of a multispecies plasma in the low collision frequency regime. Phys. Fluids 19 (5), 656.CrossRefGoogle Scholar
Jenko, F., Dorland, W., Kotschenreuther, M. & Rogers, B.N. 2000 Electron temperature gradient driven turbulence. Phys. Plasmas 7 (5), 1904.CrossRefGoogle Scholar
Jorge, R., Frei, B.J. & Ricci, P. 2019 a Nonlinear gyrokinetic coulomb collision operator. J. Plasma Phys. 85 (6), 905850604.CrossRefGoogle Scholar
Jorge, R., Ricci, P., Brunner, S., Gamba, S., Konovets, V., Loureiro, N.F., Perrone, L.M. & Teixeira, N. 2019 b Linear theory of electron-plasma waves at arbitrary collisionality. J. Plasma Phys. 85 (2), 905850211.CrossRefGoogle Scholar
Jorge, R., Ricci, P. & Loureiro, N.F. 2017 A drift-kinetic analytical model for scrape-off layer plasma dynamics at arbitrary collisionality. J. Plasma Phys. 83 (6), 905830606.CrossRefGoogle Scholar
Jorge, R., Ricci, P. & Loureiro, N.F. 2018 Theory of the drift-wave instability at arbitrary collisionality. Phys. Rev. Lett. 121 (16), 165001.CrossRefGoogle ScholarPubMed
Kadomtsev, B.B. & Shafranov, V.D. 2012 Reviews of Plasma Physics, vol. 21. Springer Science & Business Media.Google Scholar
Leyland, M.J., Beurskens, M.N.A., Frassinetti, L., Giroud, C., Saarelma, S., Snyder, P.B., Flanagan, J., Jachmich, S., Kempenaars, M., Lomas, P., et al. 2014 The h-mode pedestal structure and its role on confinement in JET with a carbon and metal wall. Nucl. Fusion 55 (1), 013019.CrossRefGoogle Scholar
Li, B. & Ernst, D.R. 2011 Gyrokinetic Fokker–Planck collision operator. Phys. Rev. Lett. 106 (19), 195002.CrossRefGoogle ScholarPubMed
Loureiro, N.F., Schekochihin, A.A. & Zocco, A. 2013 Fast collisionless reconnection and electron heating in strongly magnetized plasmas. Phys. Rev. Lett. 111 (2), 025002.CrossRefGoogle ScholarPubMed
Madsen, J. 2013 a Full-f gyrofluid model. Phys. Plasmas 20 (7), 072301.CrossRefGoogle Scholar
Madsen, J. 2013 b Gyrokinetic linearized landau collision operator. Phys. Rev. E 87 (1), 011101.CrossRefGoogle ScholarPubMed
Manas, P., Camenen, Y., Benkadda, S., Hornsby, W.A. & Peeters, A.G. 2015 Enhanced stabilisation of trapped electron modes by collisional energy scattering in tokamaks. Phys. Plasmas 22 (6), 062302.CrossRefGoogle Scholar
Mandell, N.R., Dorland, W. & Landreman, M. 2018 Laguerre–Hermite pseudo-spectral velocity formulation of gyrokinetics. J. Plasma Phys. 84 (1), 905840108.CrossRefGoogle Scholar
Mosetto, A., Halpern, F.D., Jolliet, S., Loizu, J. & Ricci, P. 2015 Finite ion temperature effects on scrape-off layer turbulence. Phys. Plasmas 22 (1), 012308.CrossRefGoogle Scholar
Pan, Q. & Ernst, D.R. 2019 Gyrokinetic Landau collision operator in conservative form. Phys. Rev. E 99 (2), 023201.CrossRefGoogle ScholarPubMed
Pan, Q., Ernst, D.R. & Crandall, P. 2020 First implementation of gyrokinetic exact linearized landau collision operator and comparison with models. Phys. Plasmas 27 (4), 042307.CrossRefGoogle Scholar
Paruta, P., Ricci, P., Riva, F., Wersal, C., Beadle, C. & Frei, B. 2018 Simulation of plasma turbulence in the periphery of diverted tokamak by using the GBS code. Phys. Plasmas 25 (11), 112301.CrossRefGoogle Scholar
Pusztai, I., Fülöp, T., Candy, J. & Hastie, R.J. 2009 Collisional model of quasilinear transport driven by toroidal electrostatic ion temperature gradient modes. Phys. Plasmas 16 (7), 072305.CrossRefGoogle Scholar
Ricci, P., Rogers, B.N. & Dorland, W. 2006 Small-scale turbulence in a closed-field-line geometry. Phys. Rev. Lett. 97 (24), 245001.CrossRefGoogle Scholar
Romanelli, F. 1989 Ion temperature-gradient-driven modes and anomalous ion transport in tokamaks. Phys. Fluids B: Plasma Phys. 1 (5), 1018.CrossRefGoogle Scholar
Romanelli, F. & Briguglio, S. 1990 Toroidal semicollisional microinstabilities and anomalous electron and ion transport. Phys. Fluids B: Plasma Phys. 2 (4), 754.CrossRefGoogle Scholar
Romanelli, F., Chen, L. & Briguglio, S. 1991 Kinetic theory of the ion-temperature-gradient-driven mode in the long wavelength limit. Phys. Fluids B: Plasma Phys. 3 (9), 2496.CrossRefGoogle Scholar
Rosenbluth, M.N., MacDonald, W.M. & Judd, D.L. 1957 Fokker–Planck equation for an inverse-square force. Phys. Rev. 107 (1), 1.CrossRefGoogle Scholar
Schneider, P.A., Angioni, C., Frassinetti, L., Horvath, L., Maslov, M., Auriemma, F., Cavedon, M., Challis, C.D., Delabie, E., Dunne, M.G., et al. 2021 The dependence of confinement on the isotope mass in the core and the edge of AUG and JET-ILW H-mode plasmas. Nucl. Fusion 62 (2), 026014.CrossRefGoogle Scholar
Smolyakov, A.I., Yagi, M. & Kishimoto, Y. 2002 Short wavelength temperature gradient driven modes in tokamak plasmas. Phys. Rev. Lett. 89 (12), 125005.CrossRefGoogle ScholarPubMed
Snyder, P.B. & Hammett, G.W. 2001 A landau fluid model for electromagnetic plasma microturbulence. Phys. Plasmas 8 (7), 3199.CrossRefGoogle Scholar
Stegmeir, A., Ross, A., Body, T., Francisquez, M., Zholobenko, W., Coster, D., Maj, O., Manz, P., Jenko, F., Rogers, B.N., et al. 2019 Global turbulence simulations of the tokamak edge region with grillix. Phys. Plasmas 26 (5), 052517.CrossRefGoogle Scholar
Sugama, H., Matsuoka, S., Satake, S., Nunami, M. & Watanabe, T.-H. 2019 Improved linearized model collision operator for the highly collisional regime. Phys. Plasmas 26 (10), 102108.CrossRefGoogle Scholar
Sugama, H., Watanabe, T.-H. & Nunami, M. 2009 Linearized model collision operators for multiple ion species plasmas and gyrokinetic entropy balance equations. Phys. Plasmas 16 (11), 112503.CrossRefGoogle Scholar
Tamain, P., Bufferand, H., Ciraolo, G., Colin, C., Galassi, D., Ghendrih, P., Schwander, F. & Serre, E. 2016 The TOKAM3X code for edge turbulence fluid simulations of tokamak plasmas in versatile magnetic geometries. J. Comput. Phys. 321, 606.CrossRefGoogle Scholar
Ulbl, P., Michels, D. & Jenko, F. 2021 Implementation and verification of a conservative, multi-species, gyro-averaged, full-f, lenard-bernstein/dougherty collision operator in the gyrokinetic code gene-x. Contrib. Plasma Phys., e202100180.Google Scholar
Waltz, R.E., Staebler, G.M., Dorland, W., Hammett, G.W., Kotschenreuther, M. & Konings, J.A. 1997 A gyro-Landau-fluid transport model. Phys. Plasmas 4 (7), 2482.CrossRefGoogle Scholar
Watanabe, T.-H. & Sugama, H. 2004 Kinetic simulation of steady states of ion temperature gradient driven turbulence with weak collisionality. Phys. Plasmas 11 (4), 1476.CrossRefGoogle Scholar
Xie, H.S., Li, Y.Y., Lu, Z.X., Ou, W.K. & Li, B. 2017 Comparisons and applications of four independent numerical approaches for linear gyrokinetic drift modes. Phys. Plasmas 24 (7), 072106.CrossRefGoogle Scholar
Zeiler, A., Biskamp, D., Drake, J.F. & Rogers, B.N. 1998 Transition from resistive ballooning to $\eta$ i driven turbulence in tokamaks. Phys. Plasmas 5 (7), 2654.CrossRefGoogle Scholar
Zeiler, A., Drake, J.F. & Rogers, B. 1997 Nonlinear reduced Braginskii equations with ion thermal dynamics in toroidal plasma. Phys. Plasmas 4 (6), 2134.CrossRefGoogle Scholar
Zhu, B., Francisquez, M. & Rogers, B.N. 2018 GDB: a global 3D two-fluid model of plasma turbulence and transport in the tokamak edge. Comput. Phys. Commun. 232, 46.CrossRefGoogle Scholar
Zocco, A. & Schekochihin, A.A. 2011 Reduced fluid-kinetic equations for low-frequency dynamics, magnetic reconnection, and electron heating in low-beta plasmas. Phys. Plasmas 18 (10), 102309.CrossRefGoogle Scholar
Figure 0

Figure 1. Plots of $\hat {\varGamma _0}(a,J)$ (a) and $\hat {\varPi }_1(a,J)$ (b) as a function of $a = b^{2}/2$ for increasing values of $J$. The solid black lines represent their asymptotic limits, $\varGamma _0$ and $- a (\varGamma _1 - \varGamma _0)$, respectively.

Figure 1

Figure 2. Collisionless sITG growth rate, $\gamma$ (a), and frequency, $\omega _r$ (b), as a function of the perpendicular wavenumber $k_\perp$ for different values of $J$. The coloured lines are the solution of the gyromoment hierarchy dispersion relation, given in (3.14), while the black lines are solutions of (3.7), i.e. the $J = \infty$ solution. Here, $\eta =3$ and $k_\parallel = 0.1$.

Figure 2

Figure 3. The sITG growth rate $\gamma$ (a) and real frequency $\omega _r$ (b) as a function of collisionality $\nu$ obtained using the gyromoment approach with the DK Coulomb collision operator (solid blue line), GENE with the same operator (red markers), the collisionless GK dispersion relation (dotted black) and the 4GM and 6GM high-collisional limits (red dashed–dotted and red dotted lines, respectively). The solution of the dispersion relation given by (4.6), derived in the case of $Q_{\parallel, \perp } = 0$, is also plotted for comparison by the solid thin black lines. Here, we consider $k_\parallel = 0.1$, $k_\perp = 0.5$ and $\tau = 1$.

Figure 3

Figure 4. Estimates of the normalized sITG peak, $\hat {\gamma }_*$, obtained from (4.8), for increasing values of $\hat {\nu }_* = 0.015$ (from (a) to (b)), as a function of $\hat {k}_{\parallel *}$ and the temperature ratio $\tau$. The local maxima of the growth rate is indicated by the dashed black lines. The growth rates $\hat {\gamma }_*$ are normalized to their maximal values.

Figure 4

Figure 5. Toroidal ITG growth rate as a function of the perpendicular wavenumber, $k_\perp$ and the temperature ratio $\tau$ obtained from (4.11), for a low (a) and a high (b) value of collisionality. The local maxima of the growth rate are indicated by the dashed black lines. The growth rates are normalized to the maximum value. Here, $\eta = 7$ and $R_B = 0.1$.

Figure 5

Figure 6. The ITG linear growth rate, $\gamma$, in the $(k_\perp, k_\parallel )$ parameter space, for increasing collisionality $\nu$ (from (a) to (c) and from ( d) to ( f)) and magnetic gradient strength $R_B$ (from (ac) to (df)) at a temperature gradient strength of $\eta =3$. The stability boundary of the gyromoment hierarchy (dotted white line) and 6GM and 4GM (dotted and dot–dashed lines, respectively) are shown for comparison. Here, we fix $(P,J) = (18,6)$.

Figure 6

Figure 7. The ITG growth rate maximized over $k_\parallel$ and $k_\perp$, denoted by $\max _{k_\parallel, k_\perp }(\gamma )$, as a function of the collisionality $\nu$ and normalized temperature gradient $\eta$ for the slab (a) and toroidal (b) ITG modes. The stability boundaries (which separates the unstable and stable modes) predicted by the gyromoment approach using the GK Coulomb collision operator and the 6GM and 4GM models are plotted by the white dotted, red dotted and red dotted–dashed lines, respectively. The values of the growth rates in the low collisionality limit are in agreement with (3.3).

Figure 7

Figure 8. Comparisons between GK (solid lines) collision operator models in the case of the slab (ac) and toroidal (df) ITG growth rate as a function of the perpendicular wavenumber $k_\perp$, from the low (a,d) to high (cf) collisionality regime. The collisionless (dashed black) and 6GM and 4GM (red dotted and dashed–dotted lines, respectively) growth rates are plotted for comparisons. Here, $\eta = 3$, $k_\parallel = 0.1$.

Figure 8

Figure 9. Same as figure 8 using the DK collision operator models.

Figure 9

Figure 10. Signed relative difference, $\sigma _C$ (defined in the main text), of the growth rate obtained using the collision operator models, used in figures 8 and 9, with respect to the GK Coulomb operator as a function of collisionality $\nu$ and temperature gradient strength $\eta$. The colourbars are saturated at a maximal relative deviation of $\sigma _C = \pm 0.4$. Here, $k_\parallel = 0.1$, $k_\perp = 0.5$ with $R_B = 0.1$ are considered, corresponding approximately to the peak ITG growth rate.

Figure 10

Figure 11. Toroidal ITG eigenvalue spectrum when $\nu = 0.1$ and $k_\perp = 0.6$ using the DK (a) and GK (b) collision operator models. Here, $R_B = 0.1$ and $\eta = 3$.

Figure 11

Figure 12. Slab SWITG growth rate (unstable when the DK Coulomb model is used) when $\eta =3$ (ac) and $\eta =5$ (df). The stability boundaries of the gyromoment hierarchy (dotted white line), 6GM and 4GM (red dotted and dot–dashed lines, respectively) are plotted for comparison. The main ITG branch is identified near $k_\perp \simeq 0.5$ and the SWITG near $k_\perp \simeq 2$, appearing as the collisionality increases when $\nu \gtrsim 0.5$. The colourbar is saturated at the maximum of $\gamma$ when $\nu = 0.05$.

Figure 12

Figure 13. Modulus of the collisionless gyromoments spectrum $N^{pj} / \phi$ (normalized to the maximum value) in the case of sITG, as a function of $p$ and $j$ plotted on a logarithmic scale and artificially saturated at $10^{-4}$. The modulus of the spectrum $N^{p0} / \phi$ and $N^{0j} / \phi$ are shown in panel (b) (solid blue and red lines, respectively). Here, $k_\parallel = 0.1$, $\eta =5$ and $k_\perp = 0.5$.

Figure 13

Figure 14. Modulus of the collisionless gyromoments spectrum $N^{pj} / \phi$ (normalized to the maximum value) in the case of purely toroidal ($k_\parallel = 0$) ITG with $R_B = 0.5$, as a function of $p$ and $j$ plotted on a logarithmic scale and artificially saturated at $10^{-4}$. The moduli of the spectrum of $N^{p0} / \phi$ and $N^{0j} / \phi$ are shown in panel (b) (solid blue and red lines, respectively). Here, $k_\perp = 0.5$, $\eta =5$.

Figure 14

Figure 15. Error norm $e(P,J)$, (7.2), on a logarithmic scale for the slab (a) and toroidal (b) ITG cases. The same parameters as in figures 13 and 14 are used, respectively. The colourbars are saturated artificially at $\log _{10} e(P,J) = -4$.

Figure 15

Figure 16. Normalized slab (a) and toroidal (b) gyromoment spectra as a function of $p$ and $j$, respectively, obtained using the GK Coulomb collision operator from (7.7) and (7.8) (solid lines) and using the GK Dougherty collision operator from (7.4) and (7.6) (dashed lines). The slab and toroidal growth rates are estimated at $\gamma \simeq 0.1$ and $\gamma \simeq 0.2$ (see figure 8), respectively, near the ITG peaks ($k_\perp \simeq 0.5$).

Figure 16

Figure 17. Convergence of the ITG linear growth rate, $\gamma$, as a function of the collisionality $\nu$ for different values of $(P,J)$, using the GK (dashed lines) and DK (solid lines) Coulomb operators, in the case of slab (a) and toroidal ((b) with $R_B = 0.1$) ITG branches. The collisionless limit (dotted black lines) and the high-collisional limits, 6GM and 4GM models (red dotted and dash–dot lines, respectively) are shown for comparison. Similar plots are obtained for the other GK and DK operator models. The parameters are $k_\perp = 0.5$, $k_\parallel = 0.1$, $\eta = 3$.

Figure 17

Figure 18. Slab ITG growth rate $\gamma$ as a function of $k_\perp$ for different values of $J$, in the low (a), intermediate (b) and high (c) collisionality regimes, with the GK (dashed lines) and DK (solid lines) Coulomb collision operators. Similar plots are obtained for the other GK and DK operator models.

Figure 18

Figure 19. Ratio between the ITG growth rate $\gamma$ obtained using the gyromoment hierarchy and the collisionless solution $\gamma _{\textrm {GK}}$ as a function of $P$ with $J = 10$ (a) and as a function of $J$ with $P=32$ (b) for two different values of $R_B$. Here, $k_\perp =0.5$, $k_\parallel = 0.1$ and $\nu = 0.001$.

Figure 19

Figure 20. Toroidal ($k_\parallel = 0$) ITG growth rate $\gamma$ as a function of the normalized magnetic gradient, $R_B$, for different $(P,J)$ and increasing collisionality (from (a) to (c)). Here, the GK (dashed lines) and DK (solid lines) Coulomb collision operators are shown, with the collisionless (black dotted) and high-collisional 4GM limit (red dotted–dashed lines).

Figure 20

Figure 21. (a) Lowest-order kernel functions $\mathcal {K}_{j}$ (solid line) and the corresponding Padé approximants, ${}_q^{p}\check {\mathcal {K}_{j}}$ for $j =0,1,2$ when $(p,q) = (1,2)$ (dashed line) and $(p,q) = (1,4)$ (dotted line). (b) The ITG linear growth rate $\gamma$ as a function of the perpendicular wavenumber $k_\perp$, obtained by the 4GM with the exact kernel functions $\mathcal {K}_{j}$ (solid line) and with the Padé approximants ${}_q^{p}\check {\mathcal {K}_{j}}$, with $(p,q) = (1,2)$ (dashed line) and $(p,q ) = (1,4)$ (dotted line), respectively.