The gyrokinetic dispersion relation of microtearing modes in collisionless toroidal plasmas

We solve the linearized gyrokinetic equation, quasineutrality condition and Ampere's law to obtain the dispersion relation of microtearing modes (MTMs) in collisionless low-$\beta$ toroidal plasmas. Consistent with past studies, we find that MTMs are driven unstable by the electron temperature gradient and that this instability drive is mediated by magnetic drifts. The dispersion relation that we derive can be evaluated numerically very quickly and may prove useful for devising strategies to mitigate MTM instability in fusion devices.


Introduction
The free energy stored in magnetically confined, toroidal plasmas gives rise to a number of plasma instabilities.As fluctuations excited by these instabilities grow, they develop into turbulence, which rapidly transports thermal energy from the plasma core to the plasma edge.The resulting heat loss is a major obstacle to developing a commercially viable fusion reactor, and finding ways to reduce turbulent transport is one of the primary goals of current fusion research.An important step towards achieving this goal is to determine the linear stability thresholds of the relevant plasma modes.
At sufficiently small values of β (the ratio of plasma pressure to magnetic pressure), it is difficult for plasma fluctuations to perturb the magnetic field, and the dominant instabilities, such as the ion-and electron-temperature-gradient modes, are electrostatic (see, e.g.Cowley et al. 1991;Dorland et al. 2000).However, as β increases, electromagnetic instabilities, such as the microtearing mode (MTM) and kinetic ballooning mode (KBM), eventually become the main drivers of turbulence.Such electromagnetic instabilities are of particular relevance to spherical tokamaks, in which β is typically several times larger than in conventional tokamaks (see, e.g., Giacomin et al. 2023;Kennedy et al. 2023, and references therein).The purpose of this paper is to derive the gyrokinetic MTM dispersion relation in the collisionless limit, which is relevant to the hot plasmas in the cores of existing and planned fusion devices.
We have organized the remainder of this paper as follows.In Sections 1.1 through 1.5, we highlight selected results from the literature and preview the main steps in our derivation of the MTM dispersion relation, which follows in detail in Section 2. In Section 3, we present several numerical examples, and in Section 4 we discuss our principal findings and conclude.

Magnetic drift waves
A simple but useful reference point for understanding the MTM is the isobaric magnetic drift wave in a plasma in which the equilibrium magnetic field B is uniform and static and neither the equilibrium electron pressure p e nor the fluctuating quantities vary along B. If we neglect electron inertia, then we can write the component of the electron momentum equation along the total magnetic field B + δB as or, equivalently, where e is the proton charge, c is the speed of light, δA = b • δA, b = B/B, δA is the perturbation to the vector potential, δB ⊥ = (∇δA ) × b is the component of δB perpendicular to B, and v * e = −cB × ∇p e /(en 0 B 2 ) is the electron diamagnetic drift velocity.This equation is equivalent to equation (D20) of Adkins et al. (2022) in the limit that λ ≫ d e , where λ is the perpendicular wavelength and d e is the electron skin depth.Equation (1.2) describes magnetic drift waves, in which δA is advected at velocity v * e .

Ballooning transformation, quasimodes, and mode rational surfaces
Throughout the rest of this paper, we consider axisymmetric toroidal equilibria and focus on individual Fourier modes ∝ exp (inζ − iωt) with infinitesimal amplitudes, where ζ is the toroidal angle, n is the toroidal mode number, and ω is the frequency.In order to enforce rapid spatial variation perpendicular to B, slow variation along B, and periodicity in the poloidal angle θ, we set n ≫ 1 and employ the ballooning transformation (Connor et al. 1978(Connor et al. , 1979;;Tang et al. 1980): where u is a vector whose components are the various fluctuating quantities, ψ is the poloidal flux, and k(ψ) is a function that is discussed in the run-up to (2.9).The triad (α, ψ, θ) is a Clebsch coordinate system, in which α(ψ, θ, ζ) (defined in (2.4)) and ψ are constant along magnetic-field lines, while the poloidal angle θ serves to measure position along B. Although the (position-space) mode u(ψ, θ, ζ) is periodic in θ, the (ballooningspace) 'quasimode' û is not.Instead, û(ψ, θ) → 0 as θ → ±∞ to ensure that the sum in (1.3) converges.As discussed further in Section 2.2, the very broad θ envelope of the MTM's electrostatic potential eigenfunction δ Φ implies that δΦ in position space is peaked around mode rational surfaces on which nq(ψ) is an integer, where q(ψ) is the safety factor defined in (2.5) (Cowley et al. 1991;Hardman et al. 2023) equilibrium field line perturbed field line θ θ Figure 1.The blue lines are segments of an equilibrium magnetic-field line that traces out a mode rational surface in a hypothetical spherical tokamak.Left: the black dashed curve shows, in an exaggerated fashion, how a segment of this field line might be perturbed by an MTM.Right: the red line highlights one of the blue field-line segments, and the black dotted line is a nearby equilibrium magnetic-field line at slightly larger ψ.We assume dq/dψ > 0, where q is defined in (2.5), so that the black dotted field line rotates through a smaller θ interval than the solid red line as the two lines traverse the same interval of toroidal angle.This magnetic shear rotates the phase fronts of the MTMs, causing them to draw closer together in the ∇ψ direction as one follows the red field-line segment from the lower right-hand side of the figure to the upper left-hand side, as illustrated schematically by the blue-and-gray-striped squares.

Tearing parity
MTMs involve δA perturbations that behave like the magnetic drift waves described in Section 1.1, propagating at a velocity ≃ v * e (see figure 3).As we discuss in greater detail in Section 2.2, a defining feature of the MTM is 'tearing parity,' which means that δ Â has a non-vanishing line integral along the magnetic field (Hatch 2010;Dickinson et al. 2011;Ishizawa et al. 2015;Patel et al. 2022).This in turn implies that, as one follows a perturbed magnetic-field line at a mode rational surface, the field line wanders secularly in the ψ direction (Hardman et al. 2023), as illustrated in the left half of figure 1. Magneticfield lines perturbed by MTMs thus create channels for electrons to transport heat down the temperature gradient, enabling MTMs to tap into the free energy stored in the electron temperature profile (Drake et al. 1980;Guttenfelder et al. 2012b).In contrast, in KBMs, a perturbed magnetic-field line at a mode rational surface returns to its initial equilibrium magnetic flux surface after each poloidal revolution about the plasma (see Section 2.2).This essential difference between the MTM and KBM is why the MTM (in contrast to the KBM) is driven by the electron temperature gradient (and the rapid transport of heat along perturbed magnetic-field lines by electrons) and not by the density gradient (see Section 2.10 and, e.g., Hazeltine et al. 1975;Drake & Lee 1977;Hassam 1980;Applegate et al. 2007;Predebon & Sattin 2013;Guttenfelder et al. 2012a;Zocco et al. 2015;Hamed et al. 2019;Geng et al. 2020;Patel et al. 2022;Hardman et al. 2023;Yagyu & Numata 2023;Giacomin et al. 2023).

Characteristic scales and quasimode eigenfunctions
The characteristic MTM binormal wavenumber satisfies k ∧ ρ e ≪ 1, where ρ e is the electron gyroradius.As we discuss further in Section 2, the δ Â fluctuation of an MTM is approximately localized to a θ interval of order unity (e.g., Applegate et al. 2007;Hamed et al. 2019;Hardman et al. 2023).Because the MTM has tearing parity, this localized δ Â fluctuation creates parallel current density δ ĵ via two very powerful mechanisms: the rapid streaming of passing electrons along perturbed magnetic-field lines that connect different equilibrium magnetic flux surfaces, and the parallel inductive electric field −c −1 (∂/∂t)δ Â .In Section 4, we label the current densities created by these two mechanisms δ ĵδB ψ and δ ĵδE , respectively, and give mathematical expressions for each.Because of the rapid motion of electrons, the non-Boltzmann part of the perturbation to the passing-electron gyrokinetic distribution function, denoted by ĥe,passing , created by these two current-generation mechanisms persists out to great distances along the magnetic field (i.e., out to |θ| ≫ 1) and generates, via the quasineutrality condition, a perturbation to the electrostatic potential δ Φ that likewise extends to |θ| ≫ 1 (Hardman et al. 2022(Hardman et al. , 2023)).
The width of the δ Φ eigenfunction in θ is ultimately limited by several factors.Magnetic shear (i.e., non-zero dq/dψ) endows δ Φ and ĥe,passing at large |θ| with spatial structure at scales ≪ k −1 ∧ in the ∇ψ direction, as illustrated in the right panel of figure 1.In addition, at large |θ|, passing electrons at the same position but different velocities that are moving towards larger |θ| will have previously interacted with electromagnetic fluctuations at smaller |θ| that were at substantially different phases, which causes ĥe,passing at sufficiently large |θ| to become a rapidly varying function of velocity.As we discuss further in Appendix B, the rapid variation of ĥe,passing (in space and velocity) causes δ Φ to decay (via gyroaveraging and phase mixing) at |θ| (k ∧ ρ e ) −1 (cf., Hardman et al. 2022Hardman et al. , 2023)).† As the MTM δ Φ eigenfunction extends out to |θ| ∼ (k ∧ ρ e ) −1 in ballooning space, the δΦ fluctuations in position space have a characteristic scale (k ∧ θ) −1 ∼ ρ e in the ∇ψ direction.

The MTM dispersion relation
Our derivation of the MTM dispersion relation in Section 2 consists of four steps to determine the four unknowns ĥe , δ Â , δ Φ, and ω, where ĥe is the non-Boltzmann part of the perturbed electron distribution function for both passing and trapped electrons.First, we integrate the gyrokinetic equation (2.25) to solve for ĥe in terms of δ Â , δ Φ, and ω.Second, we take ∂/∂θ of 1/B times the parallel component of Ampere's law to show that δ Â is localized at θ ∼ O(1).This localization implies that further appearances of δ Â are always inside its line integral along the magnetic field, ∞ −∞ JB δ Â dθ, where J = (B •∇θ) −1 is the Jacobian of the (α, ψ, θ) coordinate system.As long as it is non-zero (as it is for a tearing-parity mode), this integral can be factored out of the remaining equations as an overall normalization constant.Third, we evaluate the parallel component of Ampere's law at θ = 0 to obtain a single equation for the two remaining unknowns, ω and δ Φ.Finally, we use the quasineutrality condition to solve for δ Φ in terms of ω and plug this value back into the parallel component of Ampere's law at θ = 0.
Our analysis is greatly simplified by the two-scale nature of the problem.In particular, the contribution of δ Φ to the parallel current at θ = 0, denoted by δ ĵδΦ , is dominated by the δ Φ fluctuations at θ ∼ (k ∧ ρ e ) −1 .As a consequence, when we use the quasineutrality condition to solve for δ Φ in step 4 of the programme described in previous paragraph, we can, to leading order, restrict our attention to values of |θ| that are sufficiently large that: (1) the non-Boltzmann part of the perturbed ion distribution function ĥi can be neglected because of gyroaveraging and phase mixing, and (2) δ Â enters the quasineutrality condition only via the quantity ∞ −∞ JB δ Â dθ, as already mentioned.We note that δ ĵδΦ does not arise from the parallel electric field associated with δΦ, whose effects are included in the Boltzmann response, which is an even function of the parallel velocity and hence does not generate parallel current.Instead, δ ĵδΦ arises from electron energization or de-energization caused by the partial time derivative of the electrostatic potential energy −e∂δ Φ/∂t and from δ Φ causing electrons to E × B-drift across the equilibrium flux surfaces.†

Derivation of the MTM dispersion relation
We consider a gyrokinetic model of a single-ion-species plasma whose equilibrium state is axisymmetric.A comprehensive derivation of the equations describing this system is reviewed by Abel et al. (2013), who included plasma rotation, which we neglect for simplicity.In this model, the equilibrium distribution function of species s (with s = i for ions and s = e for electrons) is a Maxwellian, and the number density n 0 and temperature T s are flux functions: where is the thermal speed of species s, and m s is the mass of a particle of species s.
The equilibrium magnetic field B can be written in two equivalent ways: the standard form for axisymmetric equilibria, and the Clebsch form (Kruskal & Kulsrud 1958) Here, ζ is the toroidal angle, I(ψ) is the axial current divided by 2π, is the safety factor, and The θ integrals in (2.5) and (2.6) are evaluated at constant ψ.Unlike α, ν is a singlevalued, periodic function of θ.As mentioned in Section 1.2, in Clebsch coordinates (α, ψ, θ), α and ψ serve to label the magnetic-field lines, and θ determines the position along a magnetic-field line.

Ballooning transformation
As mentioned in Section 1.2, we represent all fluctuating quantities as Fourier series in the toroidal angle ζ and focus on a single Fourier mode with toroidal mode number n. † These latter two effects are represented mathematically by the terms proportional to ωδ Φ and Ω * e(E)δ Φ, respectively, on the right-hand side of (2.25).
Because the spatial variation of MTMs in the plane perpendicular to B is much more rapid than their spatial variation along B, it would be natural to take all fluctuating quantities to be of the form f (ψ, θ) exp{in[α + g(ψ)]} with |n| ≫ 1, where f (ψ, θ) and g(ψ) are slowly varying functions.However, as pointed out by Connor et al. (1978), fluctuations of this form are unphysical when q is irrational, because they are not periodic in θ.This is problematic because q is irrational in essentially all of the plasma volume when q ′ (ψ) = 0.
To circumvent this difficulty, we follow Connor et al. (1978), Tang et al. (1980), and others by employing the ballooning transformation, where u is a vector whose components are the various fluctuating quantities, |n| ≫ 1, û is a slowly varying function of ψ and θ, and B • ∇S = 0.These last three conditions guarantee rapid spatial variation, but only in directions perpendicular to B. We require that û(ψ, θ) → 0 sufficiently rapidly as |θ| → ∞ for the sum in (2.7) to converge.As mentioned in Section 1.2, we refer to u as the 'mode' and û as the 'quasimode.'The ballooning transformation represents u as the sum of an infinite number of copies of ûe inS that are translated in θ by successive integer multiples of 2π, thereby ensuring that u is periodic in θ.
As we are taking u(ψ, θ, ζ) to be ∝ e inζ with no other ζ dependence, the condition B • ∇S = 0 implies that (Tang et al. 1980) where k(ψ) is some function of ψ alone.Thus, the eikonal form conjectured in the first paragraph of this section describes the rapid cross-field spatial variation of the summand in (2.7) rather than the spatial variation of u(ψ, θ, ζ) in its entirety.The function k(ψ) can in principle be determined through a global analysis, but here we carry out a local analysis about some flux surface ψ = ψ 0 , with k(ψ 0 ) a free parameter that is related to the ballooning angle (see, e.g., Hardman et al. 2022) (2.9) If we represent the linear eigenvalue problem that determines the MTM eigenfunctions and dispersion relation in the form where L is a linear operator whose coefficients are periodic in θ with period 2π, then the condition is sufficient for u to solve (2.10).In subsequent sections, we will solve (2.11) rather than (2.10) and simplify notation by writing û(ψ, θ) = û(θ) without explicitly referencing the slow dependence on ψ.
The perpendicular wavevector is (2.12) The component of k ⊥ in the ∇ψ direction can be written in the form (Hardman et al.

Mode rational surfaces and tearing parity
With the aid of (2.4) and (2.8), we can rewrite (2.7) in the form As discussed in Section 1, δ Φ retains a comparable magnitude as |θ| increases to values ≫ 1.If one were to treat δ Φ(θ) as approximately constant out to large values of |θ|, then the sum on the right-hand side of (2.15) would add to large values (exhibiting constructive interference of quasimodes) at mode rational surfaces on which nq(ψ) = m, where m is an integer (Cowley et al. 1991).For this reason, in position space, the electrostatic-potential fluctuations of MTMs are peaked on mode rational surfaces.Mode rational surfaces have an additional significance related to the perturbed magnetic-field lines.As mentioned in Section 1.3, MTMs, unlike KBMs, satisfy the tearing-parity condition (Hatch 2010;Dickinson et al. 2011;Ishizawa et al. 2015 is the Jacobian of both the (ζ, ψ, θ) and (α, ψ, θ) coordinate systems that was previously mentioned in Section 1.5.Equation (2.16) implies that perturbed magnetic-field lines at mode rational surfaces wander secularly towards either larger or smaller ψ (Hardman et al. 2023).To show this, we parameterize the perturbed magnetic-field line that passes through position (α 1 , ψ 1 , θ 1 ) using the Clebsch coordinate functions α(θ) = α 1 + δα(θ) and ψ(θ) = ψ 1 + δψ(θ), with δα(θ 1 ) = 0 and δψ(θ 1 ) = 0. We define l(θ) to be the distance along this perturbed magnetic-field line and s ⊥ (θ) to be the distance between this perturbed field line and the equilibrium flux surface ψ = ψ 1 .To leading order in the (infinitesimal) MTM amplitude, (2.17) where we have taken s ⊥ to increase in the direction of increasing ψ and l to increase in the direction of b, δ Bψ = δ B • ∇ψ/|∇ψ|, S 1 = α 1 + ψ1 k(ψ)dψ, q 1 = q(ψ 1 ), and the θ integral in (2.17) is carried out at α = α 1 and ψ = ψ 1 .To leading order in 1/n, δ Bψ = −ik ∧ δ Â .If we take ψ = ψ 1 to be a mode rational surface on which nq 1 is an integer, then, with the aid of (2.14), we can rewrite (2.17) as (2.18) Equation (2.16) implies that the right-hand side of (2.18) is non-zero.(In contrast, ∞ −∞ JB δ Â dθ vanishes for KBMs in the low and intermediate-frequency regimes (Tang et al. 1980).)Because the right-hand side of (2.18) is a function of α 1 and ψ 1 but not θ 1 , a perturbed magnetic-field line on a mode rational surface keeps wandering in the same direction in ψ each time it winds around the plasma in the poloidal direction.

Orderings
We assume that where B p is the poloidal magnetic field, a is the plasma minor radius, and R is the plasma major radius.We take the mode's frequency ω to satisfy where is the diamagnetic drift frequency of species s, Z s e is the charge of species s, e is the proton charge, and c is the speed of light.We also assume that where ρ e = v T e /|Ω e | is the electron gyroradius, and Ω s = Z s eB/(m s c) is the cyclotron frequency of species s.We note that, from (2.14), (2.20), and (2.22), (2.24)

Linearized gyrokinetic equation
In the limit of infinitesimal fluctuation amplitudes, the ballooning-space representation of the non-Boltzmann, gyrotropic part of the perturbed gyrokinetic distribution function ĥs satisfies the linearized gyrokinetic equation (Tang et al. 1980), is the magnetic drift frequency, is the guiding-centre drift velocity, J l denotes the Bessel function of the first kind of order l, α s = k ⊥ v ⊥ /Ω s (not to be confused with the Clebsch coordinate α), and η s = d ln T s /d ln n 0 .In (2.25), the partial derivative ∂/∂θ is taken at constant ψ, α, µ, and E, where µ = v 2 ⊥ /(2B), and v ⊥ is the velocity component perpendicular to B. In writing (2.25), we neglected a term involving the parallel component of the fluctuating magnetic field, which leads to only a small correction to the MTM dispersion relation when β e ≪ 1 (Applegate et al. 2007;Patel et al. 2022;Kennedy et al. 2023).

Passing electrons
To determine ĥe for passing electrons, we solve (2.25) subject to the boundary condition which, as noted in section 2.1, is required in order for the sum in (2.7) to converge.The unique solution for Im ω > 0 is given by (Frieman et al. 1980;Tang et al. 1980) Here and in the following, the ± sign indicates the sign of v , σ J = J/|J|, J is the Jacobian defined following (2.16), and is (−|v |/v times) the change in the MTM phase factor nS − ωt at the position of a passing electron as it propagates from θ = a to θ = b.In (2.30) and (2.32), the θ ′ and θ integrals are carried out at constant ψ, α, E and µ.In (2.30) and in the following, if a function of θ appears in an integral over θ ′ but the function's arguments are not listed, the function is to be evaluated at θ ′ rather than θ.The lower limit of integration in (2.30) is chosen to ensure that ĥe± → 0 as θ → ∓σ J ∞.The condition Im ω > 0 ensures that ĥe± → 0 as θ → ±σ J ∞ because δ Φ(θ) and δ Â (θ) also vanish as |θ| → ∞.
2.6.Leading-order parallel component of Ampere's law and its θ derivative at θ ∼ O(1) In ballooning space, the parallel component of Ampere's law is (Tang et al. 1980) (2.33) We only need to evaluate (2.33) at θ ∼ O(1).The ion contribution to the parallel current can be neglected as it is ∼ (m e /m i ) 1/2 times the passing-electron contribution, the ions being much slower than the electrons.The trapped-electron contribution to the parallel current can also be neglected, as we show in Appendix A. Using the solution (2.30) in (2.33) and dividing by B, we obtain where I b a = σ J I b a , B max is the maximum value of the magnetic field on the flux surface, and the upper limit of integration of the µ integral restricts the integral to the passing region of velocity space.The circled numbers in (2.34) are shorthand for the values of the terms underneath them after all multiplications and integrations have been carried out.In Appendix B, we will show that (2.35) at θ ∼ O(1), a set of relations that we will use in our derivation of (2.36).
As a first step towards solving (2.34), we consider its θ derivative.When ∂/∂θ acts upon the right-hand side of (2.34), the resulting quantity is the sum of three terms.The first results from taking the θ derivative of J 0 (α e (θ)), which is ∝ J 1 (α e (θ)); this term is a factor ∼ k ∧ ρ e smaller than the right-hand side of (2.34), because α e (θ) ∼ k ∧ ρ e when θ ∼ O(1).The second term results from taking the θ derivative of exp ±iI θ 0 ; it follows from (2.32) that this term is a factor ∼ k ∧ ρ e smaller than the right-hand side of (2.34) when θ ∼ O(1).The third term results from evaluating the integrands in 2a , 2b , 3a and 3b at the endpoints of the θ ′ integrations.The resulting δ Â terms vanish identically.The resulting δ Φ terms are a factor of ∼ k ∧ ρ e smaller than the sum of 2a and 2b because, as we will show in Appendix B, the θ ′ integrand in 2a has a similar magnitude throughout the interval 0 < |θ| (k ∧ ρ e ) −1 before decaying at larger |θ|, and likewise for the θ ′ integrand in 2b .The integrands in terms 2a and 2b are therefore of order ∼ k ∧ ρ e times the integrals.To summarize, at θ ∼ O(1), the θ derivative of the right-hand side of (2.34) is of order ∼ k ∧ ρ e times the right-hand side of (2.34).

Quasineutrality: determining δ Φ0 at |θ| ≫ 1
In ballooning space, the quasineutrality condition is (Tang et al. 1980) where the first term on the right-hand side of (2.44) is the Boltzmann response.Our goal in this section is to use (2.44) to determine δ Φ(θ) for any given value of ω so that we can evaluate the θ integral in (2.39).As shown in Appendix B, this integral is dominated by |θ| ∼ (k ∧ ρ e ) −1 .The contribution from ĥi to (2.44) at such large values of |θ| is much smaller than the ion Boltzmann term because of gyroaveraging: at θ ∼ (k ∧ ρ e ) −1 , α i ∼ (m i /m e ) 1/2 and J 0 (α i ) ≪ 1.With the aid of (2.30), we obtain (2.45) Upon substituting (2.45) into (2.44),neglecting ĥi , and making use of the results in Appendix A for the value of ĥe for trapped electrons, we find that where τ = T i /T e , j is the largest integer such that π(2j − 1) < θ (we assume that B attains its maximum value at θ = π), (2.49) , and . . .b and τ b are, respectively, the bounce average and bounce time defined in (A 7).Equation (2.46) can be solved numerically by discretizing the integral over θ ′ , which converts (2.46) into a matrix equation for δ Φ.We present examples of such solutions in figure 5. † We note that (2.39) and (2.46) can be combined into a single eigenvalue equation for δ Φ and ω by multiplying both equations by ψ ,∞ , using (2.39) to solve for ψ ,∞ , and substituting that value into (2.46) to obtain where (2.52)Although (2.51) could be used to determine ω, in this work we obtain the MTM dispersion relation from (2.39) and (2.46).

Maximum unstable binormal wavenumber
Given the scaling estimates in Appendix B, the first term (ω), second term (ω 0 ), and fourth term (containing δ Φ) in (2.39) are all comparable, and the ratio of the third term (iπ 1/2 v T e /L) to these other terms is k ∧ ρ e /β e .These estimates also follow from (2.35) upon noting that the iπ 1/2 v T e /L term in (2.39) comes from term 1 in (2.34), the ω and ω 0 terms come from the sum of terms 3a and 3b in (2.34), and the term containing δ Φ in (2.39) comes from the sum of terms 2a and 2b in (2.34).As the iπ 1/2 v T e /L term in (2.39) is always stabilizing, the MTM can only be unstable if (cf.Hardman et al. 2023) k ∧ ρ e β e .
(2.53) 2.9.The limit of long wavelength and cold ions In this section, we specialize our results to the limit in which (2.54) As W p (θ, θ ′ ) and W tr (θ, θ ′ ) are ∝ τ when τ ≡ T i /T e ≪ 1, the approximate magnitudes of W p (θ, θ ′ ) and W tr (θ, θ ′ ) in the cold-ion limit are the same as estimated in Appendix B, only multiplied by τ .To leading order in τ , (2.46) thus yields δ Φ = τ Γ. (2.55) Upon substituting (2.55) into (2.39) and defining ω r and γ to be the real and imaginary parts of ω, respectively, we find that, to leading order in τ , where the subscript 'ω → ω 0 ' means that ω is replaced with ω 0 in (2.32) and (2.41) when evaluating Γ .

The limit of long wavelength and zero temperature gradient
In this section, we assume that (2.58) When k ∧ ρ e ≪ β e , we can drop the iπ 1/2 v T e /L term in (2.39), and when η e = 0, Ω * e (E) = ω * e .Equation (2.39) then has the solution ω = ω * e , and for this solution ω − Ω * e (E), Γ (θ), W p (θ, θ ′ ), W tr (θ, θ ′ ), and δΦ(θ) all vanish.The stability of the MTM at k ∧ ρ e ≪ β e when η e = 0 implies that the MTM at k ∧ ρ e ≪ β e is driven by the electron temperature gradient and not by the density gradient.

Numerical examples
We consider a Miller-Mercier-Luc model (Mercier & Luc 1974;Miller et al. 1998) of a local Grad-Shafranov equilibrium with parameters taken from Table 2 of Patel et al. (2022), except for the electron collision frequency ν e , which we set equal to zero.These parameters were chosen to model a flux surface in the core of the proposed STEP spherical tokamak (Wilson et al. 2020).We plot the shape of this flux surface and the θ profiles of the poloidal and total magnetic-field strengths on this surface in figure 2. Throughout this section, we take θ to be the particular poloidal angle used by Miller et al. (1998).
The linear average of β e around the flux-surface contour in the poloidal plane in this equilibrium is β e,av = 0.125. (3.1) To obtain the MTM dispersion relation for this equilibrium, we solve (2.39) using Newton's method.At each step in Newton's method, we need to determine δ Φ in (2.39) for some assumed value of ω.To do so, we discretize in θ, which converts (2.46) into a matrix equation for δ Φ.When we solve this matrix equation, we first compute I θ 0 = σ J I θ 0 from (2.32) on a θ grid with 512 uniformly spaced grid points per 2π increment in θ.We evaluate all integrals over θ, E, and µ using the trapezoid rule.We evaluate I θ 0 and all other functions of E and µ on a grid of 64 uniformly spaced grid points along the velocity (v) axis ranging from 0.1v T e to 6v T e and, at each v, 64 evenly spaced grid points in µ within the passing region of velocity space and another 64 evenly spaced grid points in µ within the trapped region of velocity space.We evaluate Γ (θ), W p (θ, θ ′ ), W tr (θ, θ ′ ), and δ Φ on a coarser θ grid with only 32 points per 2π increment in θ.We adjust the total width of the θ grid as we vary n to ensure that δ Φ decays to small values before the edge of the grid is reached.In all the examples in this section, we set θ 0 = 0.
In figure 3, we plot the real and imaginary parts of ω, denoted by ω r and γ, respectively, for seven values of n: 25, 50, 100, 200, 400, 800, and 1060.At n < 800, ω r lies between the cold-ion MTM frequency ω 0 = ω * e (1 + η e /2) and the magnetic-drift-wave frequency that follows from (1.2) and (2.14), which is ω mdw = k ⊥ • v * e = ω * e (1 + η e ), where the electron diamagnetic drift velocity v * e is defined below (1.2).The MTM frequency differs from ω mdw because (1.2) is the fluid-theory requirement for avoiding infinite current in a perfectly conducting plasma with massless electrons, whereas (2.39) is the gyrokinetic statement of the parallel component of Ampere's law, which accounts for the finite values of the current and electron mass, the current produced by δ Φ, and how an electron's response to the fluctuating fields depends upon the electron's velocity.Nevertheless, the approximate equality ω ≃ ω mdw indicates that the MTM phase velocity is approximately v * e and that the MTM approximates the force balance that arises in a fluid-theory magnetic drift wave.Across these same n values (n < 800), the MTM growth rate is approximately 1/6 to 1/5 of ω r .However, as n increases above 800, k ∧ approaches  the approximate maximum unstable binormal wavenumber β e /ρ e given in (2.53), and γ drops sharply.
In figure 4, we plot δ Â from (2.37).We note that the θ profile of δ Â in (2.37) is independent of n, and thus this same profile applies to all of the data points in figure 3.
In figure 5, we plot δ Φ for n = 25, 50, 100, 200, and 400, with n increasing from the top row to the bottom row.As n increases, the width ∆θ of the δ Φ eigenfunction decreases, in agreement with the estimate in Section 2 that ∆θ ∼ (k ∧ ρ e ) −1 .Aside from the decrease in ∆θ, the qualitative shape of the δ Φ(θ) eigenfunction and the value of γ/ω r remain similar as n ranges from 25 to 400, even though this range of n corresponds to k ∧ ρ i values ranging from less than 1 to greater than 1, where ρ i is the ion gyroradius.This latter point is consistent with the analysis of Section 2, in which the ions play no particular role in the MTM other than via their Boltzmann response.
Although we have retained the trapped-electron ĥe terms in our analysis, they have only a modest effect on the MTM growth rate for the spherical-tokamak equilibrium that we investigated in Section 3.For example, these terms increase γ by ≃ 10% for the n = 50 data point in figure 3 and by ≃ 0.6% for the n = 400 data point.

Conclusion
In this paper, we have derived the collisionless gyrokinetic MTM dispersion relation, which is given by (2.39), supplemented by the quasineutrality condition, (2.46), which determines the function δ Φ(θ, ω).In agreement with past studies, we find that the MTM is driven unstable by the electron temperature gradient rather than by the density gradient (Section 2.10), and that the MTM instability mechanism requires the electrostatic potential fluctuation to be present: when the term containing δ Φ in (4.1) is neglected, the imaginary part of ω is strictly negative.The instability mechanism also depends on magnetic drifts in a way that is quantified by the I θ 0 terms in (2.39), (2.41) and (2.46) through (2.50).As discussed in Section 2.6, (4.1) is just the parallel component of Ampere's law evaluated at θ = 0. Upon multiplication by C 2 = σ J n 0 e 2 Bv T e ψ ,∞ /(π 1/2 T e B max ω), (4.1) becomes δ ĵδE + δ ĵδB ψ − δ ĵnet + δ ĵδΦ = 0, (4.2) where δ ĵδE = C 2 ω is the parallel current density at θ = 0 produced by the inductive parallel electric field ic −1 ω δ Â , δ ĵδB ψ = −C 2 ω 0 is the parallel current density at θ = 0 produced by δ Bψ = δ B • ∇ψ/|∇ψ| = −ik ∧ δ Â (i.e., by electrons streaming along perturbed magnetic-field lines that wander across the equilibrium flux surfaces), δ ĵnet is the net parallel current density at θ = 0 (which is given by k 2 ⊥ c δ Â (0)/4π, or, equivalently, −C 2 times the term proportional to 1/L in (4.1)), and δ ĵδΦ is C 2 times the term proportional to δ Φ in (4.1), or, equivalently, the parallel current density at θ = 0 produced by δ Φ, the nature of which is discussed at the end of Section 1.5.Figure 3 shows that ω r ≃ ω 0 and γ/ω r 1/5 for MTMs at k ∧ ρ e ≪ β e in the spherical-tokamak equilibrium considered in Section 3, which implies that the δ Φ term in (4.1) is significantly smaller than the ω and ω 0 terms, and hence that δ ĵδΦ is significantly smaller than the current densities δ ĵδE and δ ĵδB ψ that result from δ Â .The MTMs in Section 3 at k ∧ ρ e ≪ β e are thus basically magnetic drift waves (with ω ≃ ω 0 ≃ k ⊥ • v * e -see Section 1.1 and the discussion of figure 3 in Section 3) that are driven weakly unstable (γ ≪ ω r ) by the electron temperature gradient via δ Φ.
The analogy to the magnetic drift wave suggests the following heuristic way of thinking about how the MTM frequency is determined.In the case of the magnetic drift wave (Section 1.1), ω is fixed by requiring that the force from the inductive parallel electric field ic −1 ω δ Â cancel out the component of the electron pressure force parallel to the total magnetic field B + δB in order to avoid unphysically large currents.Analogously, in the case of the MTM, ω is determined by requiring that δ ĵδE offset any difference between δ ĵnet and δ ĵδB ψ + δ ĵδΦ .In an approximate sense (with corrections of order δ ĵδΦ /δ ĵδB ψ ), the real part of ω is determined by the condition that the C 2 ω r part of δ ĵδE cancel δ ĵδB ψ , which is always 90 • out of phase with δ ĵnet .This first condition is very similar to the force balance that arises in the magnetic drift wave and is the reason that the MTM propagates at approximately the electron diamagnetic drift velocity v * e .The MTM growth rate is then determined by the condition that the C 2 • iγ part of δ ĵδE match the difference between δ ĵnet and the part of δ ĵδΦ that is in phase with δ ĵnet .
There is some tension between the numerical examples presented in Section 3 and previous studies that suggest that γ → 0 as ν e → 0 when k ∧ ρ i < 1 (Applegate et al. 2007;Guttenfelder et al. 2012a;Patel et al. 2022), and further work is needed to clarify the extent of, and reasons for, this discrepancy.Other potentially fruitful directions for future research include applying the results of this paper to other tokamak equilibria and stellarator equilibria and generalizing the analysis to account for collisions.

Figure 2 .
Figure2.The left, middle, and right plots show, respectively, the shape of the flux surface about which the local Grad-Shafranov equilibrium is calculated in Section 3, the strength of the poloidal magnetic field Bp on this flux surface as a function of the poloidal angle θ, and the strength of the total magnetic field B on this flux surface as a function of θ.

Figure 3 .
Figure3.The real and imaginary parts of the MTM frequency (ωr and γ, respectively) as a function of the toroidal mode number n (top axis) and k∧ρe θ=0 (bottom axis), where k∧ is the binormal wavenumber defined in (2.14).Frequencies are given in units of vTe/a, where a is the plasma minor radius.The dotted line is a plot of ω0, which is defined in (2.40), and the dash-dotted line is a plot of the magnetic-drift-wave frequency ω mdw = ω * e(1 + ηe) that follows from (1.2).The vertical dashed line shows the approximate instability threshold k∧ρe βe from (2.53), where we have set βe equal to the value in (3.1).As in all the numerical examples in this paper, we have set the ballooning angle θ0 equal to zero.

Figure 4 .
Figure4.The parallel component of the fluctuating vector potential δ Â (θ) from (2.37) normalized to its value at θ = 0 when the ballooning angle θ0 is zero.The θ profile of δ Â (θ) that follows from (2.37) is independent of n.

Figure 5 .
Figure 5.The real (left) and imaginary (right) parts of the normalized fluctuating electrostatic potential δ Φ(θ)c/[vTeδ Â (0)] when θ0 = 0. From top to bottom, the five rows correspond to the toroidal mode numbers 25, 50, 100, 200, and 400, respectively, which can be converted to k∧ρe values by comparing the upper and lower horizontal axes in figure 3.