1 Introduction
The problem of convective stability has a long history, one that has enjoyed over a century’s worth of attention by the fluids community in the study of atmospheric dynamics, oceanography, geophysics and stellar atmospheres (‘Strahlungsgleichgewicht’, to use the term coined by Schwarzschild (Reference Schwarzschild1906)). The study of convective stability in magnetised plasmas, on the other hand, is just over 50 years old, with rigorous astrophysical inquiry emerging only in the last decade and a half. Without this context, it may come as a surprise to find a contemporary paper on something as pedestrian as the linear Vlasov theory of a magnetised, thermally stratified atmosphere. But, as we demonstrate here, such surprise is not warranted – there are new and novel findings to be had. This is despite a large body of work in the magnetic-confinement-fusion literature on instability and anomalous transport in magnetised, thermally stratified plasmas (e.g. Rudakov & Sagdeev Reference Rudakov and Sagdeev1961; Mikhailovskii Reference Mikhailovskii1962; Coppi, Rosenbluth & Sagdeev Reference Coppi, Rosenbluth and Sagdeev1967; Cowley, Kulsrud & Sudan Reference Cowley, Kulsrud and Sudan1991; Horton Reference Horton1999; Dimits et al. Reference Dimits, Bateman, Beer, Cohen, Dorland, Hammett, Kim, Kinsey, Kotschenreuther and Kritz2000; Dorland et al. Reference Dorland, Jenko, Kotschenreuther and Rogers2000), as well as two encyclopædic texts by Mikhailovskii (Reference Mikhailovskii1974, Reference Mikhailovskii1992) on electromagnetic instabilities in inhomogeneous plasmas.
The culprit is that prior work on the stability of magnetised, inhomogeneous plasmas has focused almost entirely on low- $\unicode[STIX]{x1D6FD}$ plasmas, that is, those whose thermal energy is strongly subdominant to the energy of the magnetic field that threads them. This is only natural. Energetically dominant magnetic fields are routinely employed to confine terrestrial plasmas, both for experiments in basic plasma physics and in efforts to produce a self-sustaining and controlled nuclear fusion reaction, and it is largely in these contexts that studies of plasma stability have been carried out. But most astrophysical plasmas are high- $\unicode[STIX]{x1D6FD}$ , with the magnetic field energetically subdominant not only to the thermal energy but, often, to the bulk kinetic energy as well. While it may be tempting to ignore the field altogether in this case due to its dynamical weakness, one must exercise great care. Since magnetic fields serve as conduits of heat and momentum transport in magnetised plasmas, it matters where and how effectively the field is carried around by the plasma.
This paper adds to the study of magnetised, thermally stratified, high- $\unicode[STIX]{x1D6FD}$ plasmas that are commonplace in astrophysical environments such as galaxy clusters and certain classes of black-hole accretion flows. Here, we take a kinetic Vlasov (i.e. collisionless) approach to study linear stability, both complementing and extending recent work on the topic that alternatively used a magnetofluid (viz., Braginskii) treatment. In the next section (§ 2), we provide a brief review of what is known from these previous fluid analyses about the linear convective stability of a collisional magnetised plasma. We also highlight some additional physics that is brought into play by adopting a kinetic approach, including Landau resonances, collisionless wave damping and finite-Larmor-radius effects. The calculation proper begins in § 3, in which the basic equations are given (§ 3.1), the equilibrium state is constructed (§ 3.2) and the equilibrium particle trajectories are elucidated (§§ 3.3 and 3.4). In § 4 we derive the linear theory governing small-amplitude perturbations about this equilibrium and specialize it for a quasi-Maxwell–Boltzmann equilibrium distribution function. Certain technical aspects of the derivation are relegated to appendices A and B. The heart of the paper is § 5, in which the linear theory is used to assess convective stability across a wide range of perturbation wavelengths and frequencies. The physical interpretation of instability in each parameter regime is discussed in detail. Connection is made to drift-wave instabilities driven by ion- and electron-temperature gradients that are well known in the magnetic-confinement-fusion community. A supplementary appendix is provided in which the nonlinear gyrokinetic theory of a thermally stratified atmosphere is systematically derived (appendix C). Finally, a summary of our findings is given in § 6.
2 Prerequisites
2.1 Convective stability of a collisional magnetised plasma
Before embarking on a detailed kinetic analysis, it seems prudent to briefly review what is known about the convective stability of a collisional magnetised plasma, that is, one in which the collisional mean free path $\unicode[STIX]{x1D706}_{mfp}$ is intermediate between the Larmor radius $\unicode[STIX]{x1D70C}_{s}$ of plasma species $s$ and the scales of interest on which convection sets in (e.g. the thermal-pressure scale height $H$ ). Even in this case, the physics is surprisingly different than that which applies to hydrodynamic (or even ideal magnetohydrodynamic (MHD)) convection, for which the entropy gradient is the discriminating quantity for stability (Schwarzschild Reference Schwarzschild1906).
The culprit is that frozen-in magnetic-field lines serve as conduits along which heat and/or (parallel) momentum are transported from one tethered fluid element to the next, the rate of transport increasing the more aligned the magnetic field becomes with the temperature and/or velocity gradients. This is a consequence of the small Larmor radius: the magnetic field interferes strongly with motions across itself, and so collisional transport occurs most readily along magnetic lines of force (Braginskii Reference Braginskii1965). For sufficiently weak collisions or on sufficiently small length scales such that $k_{\Vert }^{2}\unicode[STIX]{x1D706}_{mfp}H\gg 1$ , where $k_{\Vert }$ is the wavenumber along the magnetic field, field lines tend towards isotherms and/or isotachs.
This has profound implications for the convective stability of such plasmas (Balbus Reference Balbus2000, Reference Balbus2001). Upward and downward fluid displacements that carry magnetic-field lines with them tend not to be isentropic but rather isothermal, since they remain magnetically and thus thermodynamically connected to conditions at their native altitude. For atmospheres in which the temperature increases in the direction of gravity, this implies that an upwardly (downwardly) displaced fluid element is always warmer (cooler) than the surroundings it is passing through. As it rises (sinks), the frozen-in field lines become ever more parallel to the temperature gradient, with the rate of entropy transfer increasing with separation. The result is instability, the temperature gradient taking precedence over the entropy gradient.
To offer a concrete example of this process, let us consider a non-rotating, vertically ( $z$ ) stratified atmosphere in hydrostatic equilibrium, threaded by a uniform, horizontal ( $x$ ) magnetic field, subthermal in its strength. For simplicity, we take the gravitational acceleration $\boldsymbol{g}=-g\hat{\boldsymbol{e}}_{z}$ to be constant and assume that sound waves propagate fast enough to ensure near incompressibility (i.e. the Boussinesq approximation). We also neglect any dynamical effect of the magnetic field; its only role is to channel the flow of entropy between magnetically tethered fluid elements. Under these conditions, a vertical displacement $\unicode[STIX]{x1D709}_{z}$ of a fluid element excites a linear perturbation in the temperature $\unicode[STIX]{x1D6FF}T$ according to
where $\unicode[STIX]{x1D714}_{cond}$ is the (field-aligned, or ‘parallel’) conduction frequency and $N$ is the Brunt–Väisälä frequency (Väisälä Reference Väisälä1925; Brunt Reference Brunt1927). (This follows from the energy equation.) When conduction is negligible ( $\unicode[STIX]{x1D714}_{cond}\rightarrow 0$ ), we find $\unicode[STIX]{x1D6FF}T/T=-(N^{2}/g)\unicode[STIX]{x1D709}_{z}$ . In an atmosphere with an upwardly increasing entropy profile ( $N^{2}>0$ ), an upwardly displaced fluid element ( $\unicode[STIX]{x1D709}_{z}>0$ ) adiabatically cools and falls back down towards its equilibrium position. Oscillations ensue with frequency $N$ , viz., $\unicode[STIX]{x2202}^{2}\unicode[STIX]{x1D709}_{z}/\unicode[STIX]{x2202}t^{2}=-N^{2}\unicode[STIX]{x1D709}_{z}$ . In contrast, when conduction is rapid ( $\unicode[STIX]{x1D714}_{cond}\gg \unicode[STIX]{x2202}/\unicode[STIX]{x2202}t$ ), the perturbed temperature satisfies $\unicode[STIX]{x0394}T/T\doteq \unicode[STIX]{x1D6FF}T/T+\unicode[STIX]{x1D709}_{z}\,\text{d}\ln T/\text{d}z=0$ ; i.e. the Lagrangian change in a fluid element’s temperature vanishes. Isentropic displacements are precluded, and the assessment of stability demands a comparison between the (unchanged) temperature of the fluid element and that of its new surroundings. For $\text{d}\ln T/\text{d}z<0$ , an upwardly (downwardly) displaced fluid element is too hot (cool), and the displacement grows exponentially: the $z$ -component of the force equation becomes $\unicode[STIX]{x2202}^{2}\unicode[STIX]{x1D709}_{z}/\unicode[STIX]{x2202}t^{2}=g|\text{d}\ln T/\text{d}z|\,\unicode[STIX]{x1D709}_{z}$ . This is the simplest version of the magnetothermal instability (MTI).
The dynamics illustrated here is complicated further by the field-aligned collisional transport of parallel momentum, which performs the dual role of viscously damping those motions which change the magnetic-field strength (at a rate $\unicode[STIX]{x1D714}_{visc}$ ) and of coupling Alfvén and slow modes so that slow-mode perturbations excite a buoyantly unstable Alfvénic response (Kunz Reference Kunz2011). Restoring the magnetic tension and the associated Alfvén frequency $\unicode[STIX]{x1D714}_{A}$ in our presentation, vertical displacements with spatial dependence $\exp (\text{i}k_{x}x+\text{i}k_{y}y+\text{i}k_{z}z)$ in the limit of rapid field-aligned conduction satisfy
where $k_{\bot }^{2}\doteq k_{y}^{2}+k_{z}^{2}$ is the square of the wavenumber perpendicular to the mean magnetic field and ${\mathcal{K}}\doteq k_{x}^{2}+k_{y}^{2}$ (cf. equation (38) of Kunz (Reference Kunz2011)). The first term in parentheses represents Alfvén waves that are polarized with the perturbed magnetic field oriented along the $y$ -axis. In the absence of parallel viscosity, they are unaffected by buoyancy. The second term in (2.2) represents slow modes, which are viscously damped at the rate $\unicode[STIX]{x1D714}_{visc}(k_{\bot }^{2}/k^{2})$ and subject to the MTI when $\text{d}\ln T/\text{d}z<0$ . For $k_{y}\neq 0$ , the right-hand side of (2.2) is non-zero, and parallel viscosity couples the Alfvén and slow modes; slow-mode perturbations induce an Alfvénic response. When the rate of transfer of parallel momentum $\unicode[STIX]{x1D714}_{visc}$ is much faster than the dynamical time scale, the slow mode is rapidly damped; in this limit, (2.2) becomes
If the temperature increases in the direction of gravity ( $\text{d}\ln T/\text{d}z<0$ ), slow-mode perturbations excite a buoyantly unstable Alfvénic response (‘Alfvénic MTI’, to borrow the term used by Kunz (Reference Kunz2011)), an effect otherwise absent in an inviscid fluid. By damping motions along field lines, parallel viscosity effectively reorients magnetic-field fluctuations to be nearly perpendicular to the background magnetic field (via flux freezing). These modes therefore display characteristics of both slow and Alfvén modes: they have density and temperature perturbations, and therefore are subject to buoyancy forces, but their velocity and magnetic-field perturbations are predominantly polarized across the mean field. As long as field-aligned conduction is rapid, regardless of whether or not parallel viscosity is effective, the maximum growth rate is $\sqrt{-g\,\text{d}\ln T/\text{d}z}$ .
This concludes our recapitulation of convective instability in a collisional magnetised plasma. We now train our focus on collisionless plasmas.Footnote 1
2.2 Collisionless physics in stratified atmospheres
The main difference between a collisional magnetised plasma comprising a stratified atmosphere and its collisionless counterpart is the ability of particles in the latter case to resonate with emergent waves in the system. For low-frequency, long-wavelength fluctuations – those with frequencies $\unicode[STIX]{x1D714}$ much smaller than the Larmor frequency $\unicode[STIX]{x1D6FA}_{s}$ and wavenumbers $k_{\Vert }$ much larger than the inverse Larmor radius $\unicode[STIX]{x1D70C}_{s}^{-1}$ ( $s$ denotes the species index) – the primary resonance is the Landau resonance, $v_{\Vert }=\unicode[STIX]{x1D714}/k_{\Vert }$ . In the presence of an electric field, an inhomogeneous magnetic field and a gravitational field, the parallel velocity $v_{\Vert }$ of a particle with mass $m_{s}$ , charge $q_{s}$ and magnetic moment $\unicode[STIX]{x1D707}_{s}$ evolves according to
where the subscript $\Vert$ on each field denotes the vector component along the local magnetic-field direction $\hat{\boldsymbol{b}}$ . Landau-resonant particles interacting with the parallel electric field $E_{\Vert }$ result in collisionless damping of the electric-field fluctuations (Landau Reference Landau1946). A similar effect is caused by the second term in (2.4), which goes by the name transit-time or Barnes (Reference Barnes1966) damping. Particles that are almost at rest with respect to the slow magnetosonic wave are subject to the action of the mirror force associated with the magnetic compressions in the wave. Since, for a monotonically decreasing distribution function ( $\unicode[STIX]{x2202}f/\unicode[STIX]{x2202}v_{\Vert }<0$ ), there are more particles with $v_{\Vert }<\unicode[STIX]{x1D714}/k_{\Vert }$ than with $v_{\Vert }>\unicode[STIX]{x1D714}/k_{\Vert }$ , the energy exchange between resonant particles and the wave leads to a net gain (loss) of energy by the particles (wave). This is true even for non-oscillatory disturbances (which will constitute the majority of the instabilities investigating in this paper): particles with $v_{\Vert }\sim 0$ are resonant with the zero-frequency wave and extract energy through betatron acceleration. Put differently, the only way to maintain perpendicular pressure balance for such a disturbance is to increase the energy of the resonant particles at the expense of the wave energy. The fewer particles there are at small parallel velocity to be resonant, the higher the damping rate (i.e. energy transfer) needs to be to maintain pressure balance. In both cases – oscillatory waves and aperiodic modes – the result is wave damping. The final term in (2.4) represents the parallel acceleration of particles once the magnetic field acquires a component along the direction of gravity: $g_{\Vert }=\boldsymbol{g}\boldsymbol{\cdot }\hat{\boldsymbol{b}}=-gb_{z}$ . In that case, the plasma streams downwards along the perturbed magnetic-field lines, with Landau-resonant particles transferring the free energy stored in the vertical magnetic-field fluctuations to the particle distribution function. As with Landau and Barnes damping, the result is collisionless damping of electromagnetic fluctuations.
There are two further differences between the physics investigated in this paper and that which was the focus of previous (fluid) studies. First, the ions and electrons need not remain in thermal equilibrium as the plasma is perturbed. In standard fluid (collisional) treatments of the MTI, the electrons are responsible for ensuring isothermal displacements, while the ions, with their larger mass, are responsible for the buoyant response. The two are connected via thermal equilibration. In the collisionless case, this connection is lost. Second, we account for finite-Larmor-radius (FLR) effects. The consideration of FLR effects does not, of course, require specialization to a collisionless plasma, but our kinetic treatment does afford a systematic and rigorous way of deriving their impact on plasma stability. In particular, we will show that fluctuations with $k_{\bot }\unicode[STIX]{x1D70C}_{i}\sim k_{\Vert }H\gg 1$ are unstable to a kinetic-Alfvén version of the MTI, in which the ions are convected across the field lines as the electrons remain isothermal along them.
These prerequisites fulfilled, our analysis begins.
3 Vlasov description of a magnetised thermally stratified atmosphere
3.1 Basic equations
We start with the kinetic Vlasov equation
which governs the phase-space evolution of the distribution function of species $s$ , $f_{s}=f_{s}(t,\boldsymbol{r},\boldsymbol{v})$ . Our notation is standard: $q_{s}$ is the charge and $m_{s}$ is the mass of species $s$ , $\boldsymbol{E}$ is the electric field, $\boldsymbol{B}$ is the magnetic field and $\boldsymbol{g}$ is the gravitational acceleration. We take the latter to lie along the $-\hat{\boldsymbol{e}}_{z}$ direction of an $(x,y,z)$ Cartesian coordinate system, $\boldsymbol{g}=-g\hat{\boldsymbol{e}}_{z}$ with $g=\text{const.}>0$ . In the non-relativistic limit, the electromagnetic fields satisfy the plasma quasi-neutrality constraint (which follows from the Poisson equation to lowest order in $k^{2}\unicode[STIX]{x1D706}_{D}^{2}$ , where $\unicode[STIX]{x1D706}_{D}$ is the Debye length),
and the pre-Maxwell version of Ampère’s law,
The magnetic field evolves according to Faraday’s law of induction,
while satisfying the solenoidality constraint $\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{B}=0$ . We neglect particle–particle collisions in this paper; the linear theory of a weakly collisional, magnetised, thermally stratified atmosphere, with which we make frequent contact, is presented in Kunz (Reference Kunz2011).
For future reference, we define the pressure tensor of species $s$ ,
as well as its mean velocity,
The total thermal pressure $P_{s}=n_{s}T_{s}\doteq (1/3)\,\text{tr}\,\unicode[STIX]{x1D64B}_{s}$ , where the temperature $T_{s}$ is defined in situ. In many of the limiting cases we investigate, the pressure tensor is diagonal in a coordinate system defined by the directions perpendicular ( $\bot$ ) and parallel ( $\Vert$ ) to the magnetic-field direction $\hat{\boldsymbol{b}}\doteq \boldsymbol{B}/B$ :
where $\unicode[STIX]{x1D644}$ is the unit dyadic, and
are the thermal pressures perpendicular and parallel to the magnetic-field direction, respectively, of species $s$ ; $T_{\bot s}$ ( $T_{\Vert s}$ ) is the perpendicular (parallel) temperature. The total thermal pressure $P_{s}=(2/3)P_{\bot s}+(1/3)P_{\Vert s}$ .
We will sometimes express the electric and magnetic fields in terms of scalar and vector potentials:
where $\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{A}=0$ (the Coulomb gauge). Working with potentials instead of fields is particularly useful when taking the gyrokinetic limit (§ 5.3, appendix C).
3.2 Equilibrium state
We seek equilibrium solutions to (3.1)–(3.4) in the presence of a uniform, stationary background magnetic field oriented perpendicular to gravity, $\boldsymbol{B}_{0}=B_{0}\hat{\boldsymbol{e}}_{x}$ , whose magnitude $B_{0}$ is such that the Larmor radius of species $s$ , $\unicode[STIX]{x1D70C}_{s}\doteq v_{th\bot s}/\unicode[STIX]{x1D6FA}_{s}$ , where $v_{th\bot s}\doteq (2T_{\bot 0s}/m_{s})^{1/2}$ is the (perpendicular) thermal speed and $\unicode[STIX]{x1D6FA}_{s}=q_{s}B_{0}/m_{s}c$ is the Larmor frequency, is much smaller than the thermal-pressure scale height $H$ of the equilibrium distribution function $f_{0s}$ , i.e. the plasma is ‘magnetised’. In a collisionless plasma, this scale separation demands that $\boldsymbol{B}_{0}\boldsymbol{\cdot }\unicode[STIX]{x1D735}f_{0s}=0$ , or $\unicode[STIX]{x2202}f_{0s}/\unicode[STIX]{x2202}x=0$ . Physically, this is due to the fast streaming along the field lines of particles setting up thermal equilibrium (see § C.2.2 for a proof). The cross-field drift of particles is instead constrained by the Larmor motion, and gradients in the field-perpendicular directions are thus allowed. In this paper, we assume $\unicode[STIX]{x2202}f_{0s}/\unicode[STIX]{x2202}y=0$ , and hence consider equilibrium distribution functions that are dependent only upon the vertical ( $z$ ) coordinate.
Under these conditions, the first moment of (3.1) becomes
where $\unicode[STIX]{x1D617}_{0s}^{(zz)}$ is the $zz$ -component of the equilibrium pressure tensor (cf. (3.5)). Summing this equation over species and using quasi-neutrality (3.2), we find
a statement of hydrostatic equilibrium. That an equilibrium electrostatic field $\boldsymbol{E}_{0}=E_{0}\hat{\boldsymbol{e}}_{z}$ is in general required may be seen by multiplying (3.11) by $q_{s}/m_{s}$ , summing over species, and using quasi-neutrality (3.2) to find
For this to vanish in a hydrogenic plasma with constant ion-to-electron temperature ratio, the thermal speed must be species independent – a rather limiting condition. Put differently, unless the ion and electron pressure gradients differ by the ion-to-electron mass ratio – an unlikely situation – an electric field must be set-up to maintain quasi-neutrality as the gravitational force is greater on the ions than the electrons.
Up to this point, we have made no assumptions about $f_{0s}$ other than its lack of $y$ dependence. We now investigate what is allowed. A necessary and sufficient condition for a stationary solution to the Vlasov equation (3.1) is that $f_{0s}$ be a function only of the integrals of the motion. These may be identified by considering the Lagrangian
where $\boldsymbol{A}_{0}=-B_{0}z\hat{\boldsymbol{e}}_{y}$ is the equilibrium vector potential and we have introduced the total (i.e. electrostatic plus gravitational) potential $\unicode[STIX]{x1D6F7}_{s}$ , defined via
the final equality following from (3.11). Since the Lagrangian (3.14) is independent of the $y$ coordinate, the $y$ -component of the canonical momentum $\boldsymbol{{\wp}}=\unicode[STIX]{x2202}{\mathcal{L}}/\unicode[STIX]{x2202}\boldsymbol{v}$ is conserved: ${\wp}_{y}=m_{s}v_{y}-m_{s}\unicode[STIX]{x1D6FA}_{s}z=\text{const}$ . This in turn implies that
Because the Lagrangian is time independent, the Hamiltonian
is also conserved. Expressing ${\mathcal{H}}$ in terms of the particle velocity $\boldsymbol{v}=\unicode[STIX]{x2202}{\mathcal{H}}/\unicode[STIX]{x2202}\boldsymbol{{\wp}}$ yields the total (i.e. kinetic plus potential) particle energy
which is an invariant of the motion. Finally, the Lagrangian is independent of $x$ , and so ${\wp}_{x}=m_{s}v_{x}$ is conserved as well. The equilibrium distribution function for species $s$ may therefore be written
In this work, we specialize to isotropic equilibrium distribution functions, i.e. $F_{0s}=F_{0s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ . In magnetised, high- $\unicode[STIX]{x1D6FD}$ plasmas, such as the intracluster medium, velocity-space anisotropy in the equilibrium distribution function can drive myriad linear instabilities at Larmor scales, all of which we seek to exclude from (at least this first pass at) the calculation.
Before proceeding any further, we forewarn the reader that, throughout the manuscript, the upper case $F_{0s}$ refers to the equilibrium distribution function written as a function of the integrals of motion ( ${\mathcal{Z}}_{s},{\mathcal{E}}_{s}$ ), while the lower case $f_{0s}$ denotes the equilibrium distribution function written as a function of the standard phase-space coordinates ( $z,\boldsymbol{v}$ ). Special attention must be paid to what variable is held fixed under partial differentiation.
3.3 Guiding-centre drifts and $\unicode[STIX]{x1D707}_{s}$ conservation
A particular feature of our equilibrium is that it supports guiding-centre drifts, namely
where we have used (3.15) to obtain the second equality (note that $v_{ds}>0$ ). The velocities of individual particles of species $s$ can thus be decomposed in terms of the parallel velocity $v_{\Vert }$ , the perpendicular velocity $v_{\bot }$ relative to the guiding-centre drift $v_{ds}$ and the gyrophase angle $\unicode[STIX]{x1D717}$ :
Despite the presence of guiding-centre drifts in the equilibrium state, the plasma itself has no net momentum. This is because the equilibrium pressure gradients drive diamagnetic drifts $\propto \unicode[STIX]{x1D735}P_{s}\boldsymbol{\times }\hat{\boldsymbol{b}}_{0}$ , which precisely cancel the guiding-centre drifts when the plasma is in hydrostatic equilibrium. This is shown explicitly in § 4.3.
Expressed in terms of the velocity-space variables $(v_{\Vert },v_{\bot },\unicode[STIX]{x1D717})$ , the integrals of the motion are
onto the $z$ axis. Differentiating (3.22), it is easy to show that the magnetic moment of a particle of species $s$ , $\unicode[STIX]{x1D707}_{s}\doteq m_{s}v_{\bot }^{2}/2B_{0}$ , is also conserved:
Were we to be considering anisotropic equilibrium distribution functions, $F_{0s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s},v_{x})$ could equivalently be replaced with $F_{0s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s},\unicode[STIX]{x1D707}_{s})$ .
3.4 Particle trajectories
Having established that the particles drift across magnetic-field lines while conserving ${\mathcal{E}}_{s}$ , ${\mathcal{Z}}_{s}$ and $\unicode[STIX]{x1D707}_{s}$ , we now determine their full phase-space trajectories. The equations of particle motion may be obtained from the characteristics of the steady-state Vlasov equation or, equivalently, from the Hamiltonian (3.17). For some set of phase-space coordinates $[\boldsymbol{r}^{\prime }(t^{\prime }),\boldsymbol{v}^{\prime }(t^{\prime })]$ , they are
Given some data $[\boldsymbol{r}(t),\boldsymbol{v}(t)]$ at time $t=t^{\prime }+\unicode[STIX]{x1D70F}$ , (3.25) may be straightforwardly integrated to yield the phase-space trajectories
4 Linear theory
In this section, we consider the linear evolution of small perturbations to the equilibrium state. By considering the response of the distribution function to the fluctuating electromagnetic fields along the phase-space characteristics (3.26), an integral relation is derived relating the current density to the fluctuating electric field (§ 4.1). This defines the conductivity tensor of the plasma and hence, in conjunction with Ampère’s law (3.3), the dispersion properties of the linear waves. The procedure is a fairly standard one (e.g. Mikhailovskii Reference Mikhailovskii1962; Stix Reference Stix1962). In §§ 4.2 and 4.3, the conductivity tensor is calculated explicitly for an equilibrium plasma described by a Maxwell–Boltzmann distribution in which the Larmor radii of all species are much smaller than the equilibrium gradient length scales. Some of the more technical steps in the derivation of the linear theory can be found in appendix A.
4.1 Perturbed distribution function and dispersion tensor
The calculation begins by decomposing the distribution function and electromagnetic fields into their equilibrium and fluctuating parts,
and substituting them into (3.1)–(3.4). Retaining terms up to first order in the perturbation amplitudes, the linearised version of the Vlasov equation (3.1) may be written as
where $\text{D}/\text{D}t$ is the time derivative taken along the unperturbed orbits given by (3.26). Equation (4.2) is solved by inverting the differential operator on the left-hand side and integrating along the phase-space trajectory $[\boldsymbol{r}^{\prime }(t^{\prime }),\boldsymbol{v}^{\prime }(t^{\prime })]$ terminating at $[\boldsymbol{r}(t),\boldsymbol{v}(t)]$ (see (3.26)):
In order to ease the notation, we henceforth drop the ‘0’ subscript on equilibrium quantities; because we are constructing a linear theory, this omission should cause no confusion.
To proceed, we assume that the equilibrium varies slowly on the scale of the fluctuations and adopt WKBJ plane-wave solutions:
and so (4.3) becomes
where $\unicode[STIX]{x1D70F}=t-t^{\prime }\geqslant 0$ is the new time variable and we have introduced the phase function
At this point, it is worth recalling that ${\mathcal{Z}}_{s}$ and ${\mathcal{E}}_{s}$ , and therefore the equilibrium distribution $F_{s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ , are constant along the orbit $[\boldsymbol{r}^{\prime }(t^{\prime }),\boldsymbol{v}^{\prime }(t^{\prime })]$ , i.e. ${\mathcal{E}}_{s}^{\prime }={\mathcal{E}}_{s}$ and ${\mathcal{Z}}_{s}^{\prime }={\mathcal{Z}}_{s}$ . We may use this fact alongside the chain rule to write
The phase-space derivatives of $F_{s}$ can then be pulled outside of the $\unicode[STIX]{x1D70F}$ -integral in (4.6) to yield, after integration by parts and some regrouping,
This equation describes the velocity-space response of the distribution function to the fluctuating electric (and, thereby, magnetic) field. The differences between this expression and that pertaining to a homogeneous plasma lie in three places. The first term in (4.9) represents the advection of the distribution function by the perturbed particle drifts. It is ultimately the cause of the instabilities we discuss in § 5. The second addition is the final term in parentheses, which becomes comparable to the usual $-T_{s}\unicode[STIX]{x2202}F_{s}/\unicode[STIX]{x2202}{\mathcal{E}}_{s}$ term when $k_{y}v_{ds}\sim \unicode[STIX]{x1D714}$ , i.e. the velocity associated with the guiding-centre drifts must be comparable to the phase speed in the $y$ direction of the waves. Finally, the phase function in (4.9) takes into account the particle drifts, and so the fluctuating electric (and thereby magnetic) fields will ultimately be sampled at spatial locations a distance more than a Larmor orbit away perpendicular to the mean magnetic field.
To complete the derivation of the general linear theory, we multiply (4.9) by $q_{s}\boldsymbol{v}$ , integrate over velocity space, and sum over species to find a linear relationship between the current density and the fluctuating electric field:
where the dispersion tensor
and $v_{A}\doteq B/\sqrt{4\unicode[STIX]{x03C0}\unicode[STIX]{x1D71A}}$ is the Alfvén speed. Setting the determinant of (4.11b ) to zero provides the sought-after dispersion relation.
4.2 Evaluation of the orbit and gyrophase integrals
To derive any practical use from (4.11), we must evaluate the integrals in (4.10a ). There are two main obstacles to doing so. First, we must choose a basis in which to express the conductivity and dispersion tensors. For the linear Vlasov theory of a homogeneous plasma, the basis typically used is defined by the magnetic field and the wavevector – the so-called Stix basis (Stix Reference Stix1962). For the problem investigated here, there is a preferred direction set by gravity and the equilibrium gradients of $F_{s}$ , and we therefore opt instead to use the $(x,y,z)$ coordinate system. While this complicates the appearance of simple wave motion in the dispersion tensor, it does greatly simplify the terms responsible for buoyancy and particle drifts, and affords a straightforward comparison with previous work on the magnetohydrodynamic stability of thermally stratified atmospheres. We thus decompose our wavevector as
for $\unicode[STIX]{x1D713}\doteq \tan ^{-1}(k_{z}/k_{y})$ , so that the phase function (4.7) becomes
The second obstacle is that the equilibrium distribution function $F_{s}$ that resides in the integrand of (4.10a ) is independent of gyrophase at fixed guiding centre, not at fixed position. As the conductivity tensor is determined by integrals over velocity $\boldsymbol{v}$ at fixed position $\boldsymbol{r}$ , rather than integrals over ${\mathcal{E}}_{s}$ at fixed ${\mathcal{Z}}_{s}$ , we must make both sides of (4.10a ) commensurate. To do so, we leverage the scale separation between the Larmor radii and the gradient length scales of the plasma by Taylor expanding $F_{s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ about the particle position $z$ and the new velocity-space variable
which is simply the particle energy ${\mathcal{E}}_{s}$ written out to zeroth order to $\unicode[STIX]{x1D70C}_{s}/H$ (cf. (3.22)). With $F_{s}$ then expressed as a function of our new phase-space variables $(z,\unicode[STIX]{x1D700}_{s})$ , the integration over gyrophase angle $\unicode[STIX]{x1D717}$ becomes straightforward.
Both of these tasks – writing (4.10) in the $(x,y,z)$ coordinate system to facilitate the orbit integral over $\unicode[STIX]{x1D70F}$ , and expanding the equilibrium distribution function about the phase-space coordinates $(z,\unicode[STIX]{x1D700}_{s})$ to ease the integration over gyrophase angle $\unicode[STIX]{x1D717}$ – are detailed in appendix A. Here, we spare the reader the details and simply state the result:
where $\text{d}^{3}\boldsymbol{v}=2\unicode[STIX]{x03C0}v_{\bot }\text{d}v_{\bot }\text{d}v_{\Vert }$ and $\ast$ denotes the complex conjugate. The velocities $\boldsymbol{u}_{n,s}$ and $\boldsymbol{w}_{n,s}$ are composed of combinations of $v_{\Vert }$ and $v_{\bot }$ with the $n$ th-order Bessel function $\text{J}_{n}(a_{s})$ and derivatives with respect to its argument $a_{s}\doteq k_{\bot }v_{\bot }/\unicode[STIX]{x1D6FA}_{s}$ , and are given by (A 4) and (A 5), respectively. Note that the partial $z$ -derivatives of $F_{s}$ in (4.15) are taken at fixed $\unicode[STIX]{x1D700}_{s}$ , not at fixed $v_{\Vert }$ and $v_{\bot }$ .Footnote 2
Up to this point, the only assumptions that have been made about the equilibrium distribution function are that it is isotropic in velocity space at fixed guiding-centre position, that it is independent of the $y$ coordinate and that its gradient length scales are much larger than the Larmor radii of all species. In order to carry out the remaining integrals in (4.15), we must choose a particular form for $F_{s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ . In the next section, we use (4.15) to compute the dispersion tensor for an isotropic Maxwell–Boltzmann distribution function.
4.3 Dispersion tensor and perturbed distribution for a Maxwell–Boltzmann equilibrium
A natural choice for our isotropic equilibrium distribution function is the quasi-Maxwell–Boltzmann distribution
where ${\mathcal{N}}_{s}$ is the number density of guiding centres. Besides being physically motivated by the systems of interest likely being near local thermodynamic equilibrium, (4.16) also has the convenient mathematical property that the integrals in (4.15) can be readily performed. To wit, the integrals over $v_{\Vert }$ may be written in terms of
which for $p=0$ is the standard plasma dispersion function (Fried & Conte Reference Fried and Conte1961), while those over $v_{\bot }$ can be calculated with the aid of the Watson (Reference Watson1966) relation
where $\unicode[STIX]{x1D6E4}_{n}(\unicode[STIX]{x1D706})\doteq I_{n}(\unicode[STIX]{x1D706})\exp (-\unicode[STIX]{x1D706})$ and $I_{n}$ is the $n$ th-order modified Bessel function. (Several particularly useful integrals derived from (4.17) and (4.18) are given in appendix B.)
As explained in § 4.2 and detailed in appendix A, we may use the scale separation between the Larmor radii and the gradient length scales of the plasma to write
in which
so that the zeroth-order equilibrium distribution function $F_{M,s}(z,\unicode[STIX]{x1D700}_{s})$ is equivalent to the standard (gyrotropic) Maxwell–Boltzmann distribution
where $v_{ths}\doteq (2T_{s}/m_{s})^{1/2}$ is the thermal speed of species $s$ . This is useful. Taking the first moment of (4.19) and using (3.21), (4.20) and (4.21), one can also show that the mean flow of the equilibrium plasma
In other words, the average equilibrium velocity is a superposition of the particle drifts, which are caused by the equilibrium electric field and gravity (see § 3.3), and the diamagnetic drift, which is driven by the equilibrium pressure gradient. The two balance in equilibrium, giving zero net flow. (In proving this result, it is important to remember that the partial $z$ -derivative of $F_{s}$ is taken at fixed $\unicode[STIX]{x1D700}_{s}$ .)
Our equilibrium distribution function now specified, we substitute (4.20) into (4.15), perform the integrals using (4.17) and (4.18) and plug the resulting conductivity tensor into (4.11b ) to find
The species- and Bessel-order-dependent tensors $\unicode[STIX]{x1D650}_{n,s}$ and $\unicode[STIX]{x1D652}_{n,s}$ , whose elements are given by (B 7) and (B 9), respectively, are functions of the $n$ th-order Landau pole $\unicode[STIX]{x1D701}_{n,s}\doteq (\unicode[STIX]{x1D714}+k_{y}v_{ds}-n\unicode[STIX]{x1D6FA}_{s})/|k_{\Vert }|v_{ths}$ and the dimensionless square of the perpendicular wavenumber $\unicode[STIX]{x1D6FC}_{s}\doteq (k_{\bot }\unicode[STIX]{x1D70C}_{s})^{2}/2$ . The perturbed distribution function may be calculated from (4.9) using (4.20). The general result is not particularly useful in the discussion that follows, and so we provide here only the leading-order (in $\unicode[STIX]{x1D70C}_{s}/H$ ) terms:
with $\unicode[STIX]{x2202}F_{s}/\unicode[STIX]{x2202}z$ being given by (4.20a ).
5 Results
The dispersion relation that results from taking the determinant of (4.24) is rather unwieldy, and it seems best not to write it down here in full. Fortunately, all of the interesting results we have found can be obtained by taking some relatively simple asymptotic limits, which are catalogued in table 1. To make contact with previous MHD calculations of the linear stability of a thermally stratified atmosphere, we first consider the long-wavelength (‘drift-kinetic’) limit, $k\unicode[STIX]{x1D70C}_{s}\ll 1$ (§ 5.1). Contact is made with the theory of kinetic MHD for a stratified plasma. For parallel propagating modes ( $k_{\bot }=0$ ), we obtain in § 5.2 the leading-order correction that incorporates FLR effects – the so-called gyroviscosity. Proceeding to yet smaller scales, in § 5.3 we allow for low-frequency fluctuations whose spatial scales perpendicular to the guide magnetic field are comparable to the Larmor radius of the particles. In this ‘gyrokinetic’ limit, the physics of drift waves becomes important, as does the damping of modes at perpendicular Larmor scales. In appendix C, we show that one may alternatively obtain the results of § 5.3 directly from the nonlinear gyrokinetic theory for a stratified atmosphere. We follow this with a discussion of our results as they relate to the ion- and electron-temperature-gradient drift-wave instabilities well known to the magnetic-confinement-fusion community (§ 5.4).
In the results that follow, we have taken the thermal-pressure scale height to be species independent, viz. $\text{d}\ln P_{s}/\text{d}z=\text{d}\ln P/\text{d}z$ . The same is true for the temperature scale height, $\text{d}\ln T_{s}/\text{d}z=\text{d}\ln T/\text{d}z$ . While we have generally allowed the equilibrium temperatures of the various species to differ (viz. $T_{i}\neq T_{e}$ ), in preparing the figures we have taken the plasma to be composed of electrons and a single species of singly charged ions ( $q_{i}=e$ , $m_{i}/m_{e}=1836$ ) with identical equilibrium temperatures ( $T_{i}/T_{e}=1$ ), unless otherwise explicitly stated. Without loss of generality, we take $k_{\Vert }>0$ . Because we have astrophysical applications in mind, we primarily focus on high- $\unicode[STIX]{x1D6FD}$ plasmas, the discussion in § 5.4 being a notable exception.
5.1 Drift-kinetic limit: $k\unicode[STIX]{x1D70C}_{s}\ll 1$
In the drift-kinetic limit, there is no practical difference between the particle position and the guiding-centre position. This simplifies things considerably. For one thing, the pressure tensor becomes diagonal in a coordinate system set by the magnetic-field direction (cf. (3.7)), so that
In addition, the higher-order contribution to the conductivity tensor, $\unicode[STIX]{x1D652}_{n,s}$ , may be completely dropped from the analysis, and $\unicode[STIX]{x1D650}_{n,s}$ takes on a relatively simple form in which contributions from $n=-1$ , $0$ and $1$ are all that is needed. Let us proceed.
We start by applying the low-frequency, long-wavelength ordering
with $\unicode[STIX]{x1D6FD}_{s}\doteq 8\unicode[STIX]{x03C0}n_{s}T_{s}/B_{0}^{2}\sim 1$ , to (4.24). (A subsidiary expansion in, e.g. large $\unicode[STIX]{x1D6FD}_{s}$ may be taken, so long as $\unicode[STIX]{x1D6FD}_{s}\nsim 1/\unicode[STIX]{x1D716}$ .) Note that this ordering precludes drift waves, since $k_{y}v_{ds}\sim k_{y}v_{ths}(\unicode[STIX]{x1D70C}_{s}/H)\ll \unicode[STIX]{x1D714}$ . Expanding the $Z_{p}(\unicode[STIX]{x1D701}_{n,s})$ functions about $\unicode[STIX]{x1D701}_{s}\doteq \unicode[STIX]{x1D714}/k_{\Vert }v_{ths}$ and using the lowest-order approximations $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{s})\simeq -\unicode[STIX]{x1D6E4}_{0}^{\prime }(\unicode[STIX]{x1D6FC}_{s})\simeq 1$ and $\unicode[STIX]{x1D6E4}_{1}(\unicode[STIX]{x1D6FC}_{s})\simeq \unicode[STIX]{x1D6FC}_{s}/2$ , the elements of the tensor $\unicode[STIX]{x1D650}_{n,s}$ (cf. (B 7)) simplify considerably. They are
Inserting (5.3) into (4.24) and setting its determinant to zero yields the general dispersion relation in the drift-kinetic limit, which we analyse for a variety of wavevector geometries in the next four subsections (§§ 5.1.1–5.1.4). Before doing so, it will serve useful to have readily available the long-wavelength limit of the perturbed distribution function (cf. (4.25)):
where
is the lowest-order mean perpendicular velocity – the (species-independent!) perturbed $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ flow. The appearance of the combination $(2\boldsymbol{v}_{\bot }\boldsymbol{\cdot }\unicode[STIX]{x1D6FF}\boldsymbol{u}_{\bot })(\unicode[STIX]{x2202}F_{s}/\unicode[STIX]{x2202}\unicode[STIX]{x1D700}_{s})$ in (5.4) suggests that, in the long-wavelength limit, it would be advantageous to redefine our energy variable to account for the kinetic energy of this flow,
such that the energies of all particles are measured in the frame of the $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ velocity. In other words, Alfvénic fluctuations do not change the form of the distribution function, but rather define the moving frame in which any changes to it are to be measured. Physically, this is because particles in a magnetised plasma adjust on a cyclotron time scale to take on the local $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ velocity; in a sense, $\unicode[STIX]{x1D70C}_{i}$ performs the role of the mean free path. This principle is what underlies Kulsrud’s formulation of kinetic MHD, in which the perpendicular particle velocities are measured relative to the $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ drift, the latter being governed by a set of MHD-like fluid equations rather than by a kinetic equation (Kulsrud Reference Kulsrud and Rosenbluth1964, Reference Kulsrud, Galeev and Sudan1983). In this frame, the perturbed distribution is gyrotropic at fixed position – a particularly useful result.
Interpreting the other terms in (5.4) is eased by shifting from $F_{s}(v_{\Vert },\unicode[STIX]{x1D700}_{s})$ to $f_{s}(v_{\Vert },w_{\bot })$ and by using Faraday’s law (4.5) to write the perpendicular components of $\unicode[STIX]{x1D6FF}\boldsymbol{E}$ in terms of $\unicode[STIX]{x1D6FF}B_{z}$ and $\unicode[STIX]{x1D6FF}B_{\Vert }$ :
Each of the terms on the right-hand side of (5.7) has a straightforward physical interpretation. The first term arises from particles sampling the equilibrium gradients as they stream along the perturbed magnetic field, or, equivalently, from the distribution function being advected in a Lagrangian sense upwards along a vertically perturbed magnetic-field line. The second term merely serves to enforce adiabatic invariance of $\unicode[STIX]{x1D707}_{s}$ , and would vanish if we were to have written our distribution function in terms of $\unicode[STIX]{x1D707}_{s}$ rather than $w_{\bot }$ . The remaining three terms – those on the second line of (5.7) – are due to the acceleration of particles along the magnetic field by the parallel electric force $q_{s}\unicode[STIX]{x1D6FF}E_{\Vert }$ , the divergence of the parallel pressure $-(1/n_{s})\unicode[STIX]{x1D735}\boldsymbol{\cdot }(\hat{\boldsymbol{b}}P_{\Vert s})\rightarrow -(\unicode[STIX]{x1D6FF}B_{z}/B)(1/n_{s})\text{d}P_{s}/\text{d}z$ , and the mirror force $-\unicode[STIX]{x1D707}_{s}\hat{\boldsymbol{b}}\boldsymbol{\cdot }\unicode[STIX]{x1D735}B\rightarrow -\unicode[STIX]{x1D707}_{s}\text{i}k_{\Vert }\unicode[STIX]{x1D6FF}B_{\Vert }$ . Each of these forces may lead to collisionless damping if a significant fraction of the particles in the distribution function is Landau resonant, $\unicode[STIX]{x1D714}\sim k_{\Vert }v_{\Vert }$ , thereby transferring the free energy stored in the lower-order moments of (5.7) (e.g. density, momentum, pressure) to the higher-order ‘kinetic’ moments (i.e. fine-scale structure in velocity space). Lower-order moments of particular interest in this section are the perturbed density of species $s$ ,
the perturbed perpendicular pressure of species $s$ ,
The final preparatory exercise before proceeding with the analysis of the drift-kinetic dispersion relation is to compute the parallel electric field found in (5.4). This is most easily accomplished by multiplying (5.8) by $q_{s}n_{s}$ , summing the result over species and using quasi-neutrality (3.2) to eliminate all but the Landau terms. The result is
This equation amounts to a statement of parallel pressure balance, taking into account the forcing out of large-pitch-angle particles (those with $\cos ^{-1}(v_{\Vert }/v)\sim \unicode[STIX]{x03C0}/2$ ) from regions of increased magnetic-field strength and the squeezing of particles upwards along vertically displaced magnetic-field lines by the pressure gradient.
The rapid establishment of quasi-neutrality by the parallel electric field (5.11) means, upon comparison with (5.5), that $\unicode[STIX]{x1D6FF}E_{\Vert }/\unicode[STIX]{x1D6FF}E_{\bot }\sim k\unicode[STIX]{x1D70C}_{s}\ll 1$ . Faraday’s law (4.5) then expresses the freezing of the magnetic flux into the mean perpendicular flow of the perturbed plasma,
a natural outcome at long wavelengths. Introducing the Lagrangian perpendicular displacement via $\unicode[STIX]{x1D6FF}\boldsymbol{u}_{\bot }\doteq \text{d}\unicode[STIX]{x1D743}_{\bot }/\text{d}t\rightarrow -\text{i}\unicode[STIX]{x1D714}\unicode[STIX]{x1D743}_{\bot }$ , (5.12) may be written
This particular form of the induction equation affords a better physical understanding of the second moments (5.9b ) and (5.10b ), as follows. Using (5.13) to replace $\unicode[STIX]{x1D6FF}B_{z}/B$ by $\text{i}k_{\Vert }\unicode[STIX]{x1D709}_{z}$ in (5.9b ) and (5.10b ), we find that the Lagrangian changes in the perpendicular and parallel temperatures as a fluid element of species $s$ is displaced are
5.1.1 Drift-kinetic limit: parallel propagation $(\boldsymbol{k}\Vert \hat{\boldsymbol{b}})$
The simplest wavevector geometry to analyse is $k_{\bot }=0$ , in which case (4.24) and (5.3) combine to give $\unicode[STIX]{x1D63F}^{(xx)}\sim O(1)$ , $\unicode[STIX]{x1D63F}^{(xy)}\sim O(\unicode[STIX]{x1D716})$ , $\unicode[STIX]{x1D63F}^{(yy)}\sim \unicode[STIX]{x1D63F}^{(zz)}\sim O(\unicode[STIX]{x1D716}^{2})$ , and $\unicode[STIX]{x1D63F}^{(xz)}=\unicode[STIX]{x1D63F}^{(zx)}=0$ . It is straightforward to show further that the contributions from $\unicode[STIX]{x1D650}_{n,s}^{(yz)}$ and $\unicode[STIX]{x1D650}_{n,s}^{(zy)}$ to the dispersion tensor vanish to leading order by quasi-neutrality (3.2), a fact that will become important again in § 5.2 when we consider gyroviscous contributions to the wave dispersion. The long-wavelength dispersion relation then reads
The first solution, $\unicode[STIX]{x1D63F}^{(zz)}=0$ , corresponds to Alfvén waves with frequencies $\unicode[STIX]{x1D714}=\pm k_{\Vert }v_{A}$ , which are polarized with $\unicode[STIX]{x1D6FF}\boldsymbol{B}$ along the $y$ -axis. They are unaffected by buoyancy, and are completely decoupled from the compressive fluctuations. The latter satisfy
where
The terms on the right-hand side of (5.16) correspond, respectively, to the effects of magnetic tension, buoyancy and collisionless damping by way of Landau-resonant particles undergoing acceleration by a combination of parallel electric fields and parallel gravity (both of which are proportional to the pressure gradient, and thus have been combined by introducing $\unicode[STIX]{x1D6F6}_{s}$ ). It is instructive to compare (5.16) with its more familiar Braginskii-MHD (i.e. collisional, magnetised) counterpart (cf. Balbus Reference Balbus2000; Kunz Reference Kunz2011), written in a rather auspicious guise:
where $\unicode[STIX]{x1D705}$ is the (electron-dominated) thermal diffusivity of a collisional plasma. At large enough wavenumbers such that their final terms (those proportional to $\text{d}\ln P_{s}/\text{d}z$ ) are negligible, both solutions return what now bears the moniker ‘magnetothermal instability’ (Balbus Reference Balbus2000, Reference Balbus2001):
which may be unstable in an atmosphere with $\text{d}\ln T_{s}/\text{d}z<0$ , provided the magnetic field is sufficiently subthermal (i.e. $\sqrt{\unicode[STIX]{x1D6FD}}\gg k_{\Vert }H\gg 1$ ).
That both the collisional and collisionless calculations return nearly the same result in the long-wavelength limit is encouraging – the MTI in its simplest form passes unscathed through to the kinetic case. But while, in both cases, the MTI is driven by the rapid transport of heat along perturbed magnetic-field lines, there is quite different physics at play in what governs that transport. For the collisional case, in which particles interact with one another on the collisional mean free path $\unicode[STIX]{x1D706}_{mfp}$ , the MTI grows maximally once $k_{\Vert }^{2}\unicode[STIX]{x1D706}_{mfp}H\gtrsim 1$ , that is, for wavelengths intermediate between the mean free path and the thermal-pressure scale height. The lower bound of $\unicode[STIX]{x1D706}_{mfp}$ is readily understood: fluid elements separated by a fluctuation wavelength must be able to equilibrate their temperatures via conduction, short-circuiting the usual adiabatic response and allowing, say, hot fluid elements to remain hot as they are convectively lifted upwards into cooler surroundings. The upper bound of $H$ results from a displaced fluid element needing to maintain pressure balance with its surroundings by radiating sound waves faster than the fluid element rises (i.e. $k_{\Vert }v_{thi}\gg \unicode[STIX]{x1D714}\sim v_{thi}/H$ ).Footnote 3 In the collisionless case, conduction is not set by particle–particle interactions, and sounds waves are collisionlessly damped. The application of the collisional MTI to a collisionless plasma is not obvious a priori.
In fact, the collisionless case is somewhat simpler than the collisional case, in that the same condition that allows displaced fluid elements to maintain pressure balance with their surroundings (that is, $\unicode[STIX]{x1D701}_{s}\ll 1$ ) also ensures that such displacements proceed approximately isothermally. Indeed, (5.9)–(5.14) with $k_{\bot }=0$ and $\unicode[STIX]{x1D701}_{s}\ll 1$ imply
Simply put, conduction and buoyancy in a collisionless plasma are inextricably linked. This is not true in a collisional, magnetised plasma, for which (cf. Kunz Reference Kunz2011)
In this case, pressure balance can be achieved regardless of whether conduction is rapid enough to ensure nearly isothermal displacements (just take $k_{\Vert }\unicode[STIX]{x1D706}_{mfp}\rightarrow 0$ ). By contrast, in a collisionless plasma, particles redistribute heat along field lines not by communicating thermodynamic information via collisions with other particles, but by free streaming along perturbed magnetic fields. This naturally occurs on roughly the same time scale on which a sound wave propagates and establishes pressure balance. If the perturbed field lines are not isothermal, it is for the same reason that a fluid element cannot maintain pressure equilibrium with its surroundings, that is, an appreciable number of particles are Landau resonant with the fluctuation.
Figure 1 displays the growth rates as calculated from (solid line) the collisionless drift-kinetic dispersion relation (5.16) and (dotted and dashed lines) the collisional Braginskii-MHD dispersion relation (e.g. equation (24) of Kunz Reference Kunz2011) with $k_{\bot }=0$ . All solutions asymptote to a value close to the maximum growth rate $\sqrt{-g\,\text{d}\ln T/\text{d}z}$ before being cutoff by magnetic tension at $k_{\Vert }H\approx \sqrt{\unicode[STIX]{x1D6FD}\,\text{d}\ln T/\text{d}\ln P}$ . Their approach to the maximum growth rate, however, is quite different. These trends may be retrieved by balancing the final (stabilizing) terms in (5.16) and (5.18) with the destabilizing temperature gradient, which gives $\unicode[STIX]{x1D714}\propto k_{\Vert }$ for the collisionless case and $\unicode[STIX]{x1D714}\propto k_{\Vert }^{2}$ for the collisional case. As explained above, this is related to the scale on which conduction can efficiently isothermalize displaced magnetic-field lines (cf. (5.21) and (5.21)).
5.1.2 Drift-kinetic limit: oblique propagation with $k_{\bot }=k_{z}$ $(\unicode[STIX]{x1D713}=\unicode[STIX]{x03C0}/2)$
Next, we relax our assumption that $k_{\bot }=0$ and allow for wavevectors with a vertical component: $\boldsymbol{k}=k_{\Vert }\hat{\boldsymbol{e}}_{x}+k_{z}\hat{\boldsymbol{e}}_{z}$ . The first thing to notice is that we no longer have $\unicode[STIX]{x1D6FF}B_{\Vert }=0$ , but rather $\unicode[STIX]{x1D6FF}B_{\Vert }/B=-\text{i}k_{z}\unicode[STIX]{x1D709}_{z}$ (see (5.13)), and so the Lagrangian change of the perpendicular temperature (5.14) as a fluid element is displaced upwards or downwards no longer identically vanishes. Instead, changes in temperature go hand-in-hand with changes in magnetic-field strength. Because the latter are collisionlessly damped (Barnes Reference Barnes1966), this makes the physics of the MTI much richer than in the collisional case, for which $\unicode[STIX]{x0394}T\simeq 0$ in the fast-conduction limit regardless of wavevector orientation. This is particularly true at large $k_{z}$ , as we now show.
For $k_{\bot }=k_{z}$ , we have $\unicode[STIX]{x1D63F}^{(xx)}\sim O(1)$ , $\unicode[STIX]{x1D63F}^{(xy)}\sim O(\unicode[STIX]{x1D716})$ and $\unicode[STIX]{x1D63F}^{(xz)}\sim \unicode[STIX]{x1D63F}^{(yy)}\sim \unicode[STIX]{x1D63F}^{(zz)}\sim O(\unicode[STIX]{x1D716}^{2})$ . (As in § 5.1.1, the contributions from $\unicode[STIX]{x1D650}_{n,s}^{(yz)}$ and $\unicode[STIX]{x1D650}_{n,s}^{(zy)}$ vanish by quasi-neutrality.) The leading-order dispersion relation is thus of the same form as (5.15), with the Alfvénic fluctuations $\unicode[STIX]{x1D714}=\pm k_{\Vert }v_{A}$ completely decoupled from the compressive ones. The dispersion relation for the latter may be manipulated into a form similar to the $k_{\bot }=0$ result (5.16),
the only substantive difference between the two being that $\unicode[STIX]{x1D6F6}_{s}$ (see (5.17)) is replaced by
The additional term in $\widetilde{\unicode[STIX]{x1D6F6}}_{s}$ that is proportional to $k_{z}^{2}H^{2}$ arises from the Barnes damping of magnetic-field-strength fluctuations, an effect we anticipate to stabilize the drift-kinetic MTI once that term becomes comparable to the destabilizing temperature gradient – that is, once $\unicode[STIX]{x1D701}_{i}k_{z}^{2}H^{2}\sim 1$ or $k_{z}H\sim \sqrt{k_{\Vert }H}$ . Just before that, where $\unicode[STIX]{x1D701}_{i}\ll (k_{z}/k_{\Vert })^{2}$ , the two terms in brackets in (5.22) must balance and the growth rate ought to decrease as $k_{z}^{-2}$ .
This is indeed what we find in figure 2(a), which shows the instability growth rate as a function of $k_{\Vert }$ and $k_{z}$ . As for the $k_{\bot }=0$ modes, the growth rate increases with $k_{\Vert }$ until it attains its maximum value ${\approx}\sqrt{-g\,\text{d}\ln T/\text{d}z}$ at $k_{\Vert }H\sim 100$ . Beyond $k_{\Vert }H\approx \sqrt{\unicode[STIX]{x1D6FD}\,\text{d}\ln T/\text{d}\ln P}$ , magnetic tension stabilizes the mode and it propagates as a pseudo-Alfvén wave undergoing weak Landau damping. As anticipated, increasing $k_{z}$ reduces the growth rate $\propto k_{z}^{-2}$ , strongly suppressing the drift-kinetic MTI once $k_{z}H\sim \sqrt{k_{\Vert }H}$ (demarcated by the white dashed line in the figure). Note that, in the collisional case, the limiting $k_{z}$ is set by the Braginskii viscosity and is instead given by $k_{z}H\sim \sqrt{H/\unicode[STIX]{x1D706}_{mfp}}$ , independent of parallel wavenumber (see § 4.2.1 and figure 5(b) of Kunz Reference Kunz2011). The difference is because, in the collisional case, wave damping is caused by the parallel diffusion of momentum (occurring at a rate ${\sim}k_{\Vert }^{2}\unicode[STIX]{x1D706}_{mfp}v_{ths}$ ) rather than by collisionless wave–particle interactions ( ${\sim}k_{\Vert }v_{ths}$ ).
5.1.3 Drift-kinetic limit: oblique propagation with $k_{\bot }=k_{y}~(\unicode[STIX]{x1D713}=0)$
The story gets a bit more complicated when the wavevector has a component perpendicular to both the mean magnetic field and the equilibrium gradients, $\boldsymbol{k}=k_{\Vert }\hat{\boldsymbol{e}}_{x}+k_{y}\hat{\boldsymbol{e}}_{y}$ . First, the contributions from $\unicode[STIX]{x1D650}_{n,s}^{(yz)}$ and $\unicode[STIX]{x1D650}_{n,s}^{(zy)}$ to the dispersion tensor no longer vanish. This is crucial, as these matrix elements are responsible for coupling Alfvén waves and slow modes, a feature heretofore lacking (and one that will reappear in § 5.2). Secondly, with $\unicode[STIX]{x1D63F}^{(xx)}\sim O(1)$ , $\unicode[STIX]{x1D63F}^{(xy)}\sim \unicode[STIX]{x1D63F}^{(xz)}\sim O(\unicode[STIX]{x1D716})$ and $\unicode[STIX]{x1D63F}^{(yy)}\sim \unicode[STIX]{x1D63F}^{(yz)}\sim \unicode[STIX]{x1D63F}^{(zz)}\sim O(\unicode[STIX]{x1D716}^{2})$ , all the elements of the dispersion tensor ultimately contribute to the dispersion relation, and the algebra gets considerably more thorny. Surprisingly, a fortuitous rearrangement of terms leads to a rather compact dispersion relation:
The first term on the left-hand side – the one in braces – should look familiar from § 5.1.1 (cf. (5.16)). It captures the effects of magnetic tension, buoyancy and collisionless damping on fluctuations polarized with $\unicode[STIX]{x1D6FF}\boldsymbol{B}$ in the $z$ direction. The second term on the left-hand side of (5.24) should also look familiar; it includes the contribution from Barnes damping of $\unicode[STIX]{x1D6FF}B_{\Vert }\neq 0$ fluctuations ( $\propto k_{y}^{2}(1+\unicode[STIX]{x1D6F6}_{s})$ ) that featured prominently in § 5.1.2. The right-hand side is new, and deserves some discussion.
Recall from §§ 5.1.1 and 5.1.2 that, for $k_{\bot }=0$ or $k_{\bot }=k_{z}$ , Alfvén waves polarized with $\unicode[STIX]{x1D6FF}\boldsymbol{B}$ in the $y$ direction are completely decoupled from all other fluctuations. This is because they change neither the density nor the temperature (see (5.8)–(5.10)), and therefore are not subject to buoyancy forces. This is clearly not the case when $k_{y}\neq 0$ . These fluctuations are now accompanied by compressions and rarefactions in the magnetic-field lines in order to preserve the divergence-free constraint, $\unicode[STIX]{x1D6FF}B_{\Vert }=-(k_{y}/k_{\Vert })\unicode[STIX]{x1D6FF}B_{y}$ . Accompanying these magnetic-field-strength fluctuations are changes in temperature and density (see (5.14)), which influence the behaviour of vertical displacements subject to gravity. As a result, all components of the magnetic field are coupled.
In particular, consider the limit $k_{y}^{2}H^{2}\gg k_{\Vert }H$ . The first bracket on the second line of (5.24) – the one associated with $y$ -polarized magnetic-field fluctuations – is then dominated by the $\propto k_{y}^{2}$ term. This corresponds to strong Barnes damping of $\unicode[STIX]{x1D6FF}B_{\Vert }$ fluctuations. In the previous section (§ 5.1.2), this damping was detrimental to the growth of the drift-kinetic MTI, since, there, $\unicode[STIX]{x1D6FF}B_{z}=-(k_{z}/k_{\Vert })\unicode[STIX]{x1D6FF}B_{\Vert }$ . However, in the present arrangement, the system has a mechanism to minimize the collisionless damping: simply counter the mirror force $-\unicode[STIX]{x1D707}_{s}\unicode[STIX]{x1D735}_{\Vert }\unicode[STIX]{x1D6FF}B_{\Vert }$ with the parallel gravitational force $m_{s}\boldsymbol{g}\boldsymbol{\cdot }\unicode[STIX]{x1D6FF}\hat{\boldsymbol{b}}$ . As a Landau-resonant particle accelerates downwards along a vertically perturbed magnetic-field line into a region of lower potential energy, arrange the magnetic field so that the particle is also entering a region of enhanced magnetic-field strength. The mathematical manifestation of this arrangement is the cancellation of $k_{y}^{2}$ from the Barnes-damping term in (5.24) with the $k_{y}^{2}$ on that equation’s right-hand side, the result being that
where
Comparing (5.25) with (5.16), we see that the collisionless damping (captured by the final term in those equations) is reduced from the $k_{\bot }=0$ case by a factor of $\simeq 2$ at large $k_{\Vert }H$ – that is, the growth rate is actually larger for $k_{y}\neq 0$ , independent of $k_{y}$ (so long as $k_{y}^{2}H^{2}\gg k_{\Vert }H\gg 1$ ).
This scenario is illustrated in figure 3, which shows the geometry of perturbed magnetic-field lines for an MTI-unstable mode with $k_{\Vert }H=k_{y}H\gg 1$ . The horizontal field fluctuations satisfy $\unicode[STIX]{x1D6FF}B_{\Vert }=-\unicode[STIX]{x1D6FF}B_{y}$ in order to preserve the solenoidality constraint (shown in the left panel), with vertical displacements of field lines anti-correlated with the strength of the perturbed field (indicated by the colouring of field lines on the right). As a result of this arrangement, the mirror force acts oppositely to the parallel gravitational force, reducing the parallel acceleration experienced by Landau-resonant particles. Perhaps the simplest way of understanding this advantageous field-line topology is to examine the energy of a $\unicode[STIX]{x1D707}$ -conserving particle in the frame of the perturbed $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ drift, $(1/2)m_{s}v_{\Vert }^{2}+\unicode[STIX]{x1D707}_{s}B+q_{s}\unicode[STIX]{x1D6F7}_{s}$ . The change in energy of this particle as it interacts with the wave is reduced by ensuring that increasing (decreasing) magnetic-field strengths occur in regions of decreasing (increasing) potential.
At the risk of belabouring this point, let us use (5.8) and (5.13) to write the evolution equation for the vertical displacement of a fluid element (neglecting magnetic tension):
Growth is maximized by reducing the amplitude of the final term, which is achieved by having magnetic-field-strength fluctuations $\unicode[STIX]{x1D6FF}B_{\Vert }$ anti-correlated with upward displacements $\unicode[STIX]{x1D709}_{z}>0$ into regions of larger potential ( $\text{d}\unicode[STIX]{x1D6F7}_{s}/\text{d}z>0$ ). This is precisely what is shown in figure 3.
5.1.4 Drift-kinetic limit: arbitrary wavevector orientation
Having now considered all simple wavevector geometries, we finally examine the general drift-kinetic dispersion relation, which, after some straightforward but tedious algebra, may be written in the following convenient form:
with $\widetilde{\unicode[STIX]{x1D6F6}}_{s}$ given by (5.23). Setting the first term in braces to zero returns the dispersion relation analysed in § 5.1.2; it includes the effects of collisionless damping due to the parallel electric field, the magnetic mirror force and the parallel gravitational force. The second term in brackets (the one on the second line of (5.28)) is identical to the second term in (5.24), and is related to fluctuations polarized with $\unicode[STIX]{x1D6FF}\boldsymbol{B}$ in the horizontal plane; it includes the contribution from Barnes damping of $\unicode[STIX]{x1D6FF}B_{\Vert }\neq 0$ fluctuations. The first term on the right-hand side, responsible for coupling these two branches, was discussed thoroughly in the previous section. The final term on the right-hand side is new, requiring both $k_{y}\neq 0$ and $k_{z}\neq 0$ .
Motivated by the results in the previous section, we take the limit $k_{y}^{2}H^{2}\gg k_{\Vert }H$ (while keeping $k_{z}\lesssim k_{y}$ ). Equation (5.28) becomes
which, aside from the tension term, is exactly the same as the $k_{z}=0$ solution (5.25)! The stabilization of the drift-kinetic MTI at moderate $k_{z}$ due to Barnes damping seen in § 5.1.2 has been removed by the same physics as was the focus of the previous section – the mirror force, by which Landau-resonant particles bleed energy from the magnetic-field-strength fluctuations, is offset by the divergence of the parallel pressure.
Such behaviour can be seen in figure 4, which displays the instability growth rate calculated from (5.28) as a function of $k_{y}H$ and $k_{z}H$ for $k_{\Vert }H=61.4$ (corresponding to the fastest-growing $k_{\bot }=0$ mode for $\unicode[STIX]{x1D6FD}_{i}=10^{6}$ ). At large $k_{y}$ , the growth rate is maximal and approximately constant with $k_{z}$ , at least until $k_{y}\sim k_{z}$ . At that point the ordering leading to (5.29) breaks down, and for $k_{y}^{2}/k_{z}^{2}\ll 1$ we have from (5.28)
to leading order in $k_{y}^{2}/k_{z}^{2}$ , i.e. the growth decreases as $\propto k_{z}^{-1}$ . The physical reason for this decrease is that, when $k_{z}\gg k_{y}$ , the mirror force can no longer be compensated adequately by the divergence of the parallel pressure, and Barnes damping wins out. The mathematical manifestation of this is an additional term in the evolution equation for the vertical displacement of a fluid element due to the perturbed perpendicular pressure (cf. 5.27); neglecting magnetic tension and pressure,
Minimizing the second term in the final equality of (5.31) by having upward displacements go hand-in-hand with weak magnetic-field strengths – an arrangement exploited by the $k_{y}^{2}H^{2}\gg k_{\Vert }H$ modes in § 5.1.3 – causes the final term in (5.31) to strongly damp the fluctuations when $k_{z}^{2}H^{2}\gg k_{\Vert }H$ .
5.2 Gyroviscous limit: $k_{\Vert }^{2}\unicode[STIX]{x1D70C}_{s}H\sim 1$
Having examined the long-wavelength limit, we now progress to smaller scales, namely, those intermediate between the atmospheric scale height and the species’ Larmor radii. To keep things simple at first, we focus only on fluctuations whose wavevectors are aligned with the mean magnetic field, $\boldsymbol{k}=k_{\Vert }\hat{\boldsymbol{e}}_{x}$ . A natural question to ask is what would happen to the growth rate of the MTI if the magnetic tension were insufficient to provide a parallel-wavenumber cutoff, thus allowing for growth at relatively small scales. Indeed, many numerical simulations of the (collisional) MTI presented in the published literature had adopted rather large initial plasma $\unicode[STIX]{x1D6FD}$ parameters $\unicode[STIX]{x1D6FD}\sim 10^{8}$ – $10^{12}$ (Parrish & Stone Reference Parrish and Stone2007; Parrish, Stone & Lemaster Reference Parrish, Stone and Lemaster2008; McCourt et al. Reference McCourt, Parrish, Sharma and Quataert2011; Parrish et al. Reference Parrish, McCourt, Quataert and Sharma2012), one of the motivations being to test whether the MTI in its nonlinear phase might act as an efficient dynamo capable of explaining the presence of the ${\sim}\unicode[STIX]{x1D707}\text{G}$ magnetic fields currently observed in nearby galaxy clusters (e.g. Carilli & Taylor Reference Carilli and Taylor2002). At such high $\unicode[STIX]{x1D6FD}$ , the ion-Larmor radius can become large enough for the characteristic scale at which FLR effects are important – roughly the geometric mean of the ion-Larmor radius and the scale height, viz.
– to be comparable to the size of some fluctuations, i.e. $k_{\Vert }^{2}\unicode[STIX]{x1D70C}_{i}H\sim 1$ . Because FLR effects introduce wave dispersion and thus are typically a stabilizing influence, one might inquire at what magnetic-field strength does such stabilization become more important than that provided by magnetic tension. The answer, as it turns out, is when $\unicode[STIX]{x1D6FD}_{i}\gtrsim H/\unicode[STIX]{x1D70C}_{i}$ or, for parameters characteristic of galaxy cluster outskirts,
Such a field is within the range of observationally allowed primordial magnetic-field strengths (e.g. Neronov & Vovk Reference Neronov and Vovk2010; Planck Collaboration et al. Reference Ade, Aghanim, Arnaud, Arroja, Ashdown, Aumont, Baccigalupi, Ballardini and Banday2015).
This short preamble concluded, we begin our reduction of the general dispersion relation (4.24) in this ‘gyroviscous’ limit by adopting the ordering
under which the components of the dispersion tensor satisfy $\unicode[STIX]{x1D63F}^{(xx)}\sim O(\unicode[STIX]{x1D716})$ , $\unicode[STIX]{x1D63F}^{(xy)}\sim \unicode[STIX]{x1D63F}^{(yy)}\sim \unicode[STIX]{x1D63F}^{(yz)}\sim \unicode[STIX]{x1D63F}^{(zz)}\sim O(\unicode[STIX]{x1D716}^{2})$ . The slight complication to obtaining these estimates is that the $Z_{0}(\unicode[STIX]{x1D701}_{n,s})$ function for $n=\pm 1$ must be expanded in its large arguments to third order in $k_{\Vert }\unicode[STIX]{x1D70C}_{s}$ , so that, e.g.
(Recall from § 5.1.1 that the first term in (5.35) ultimately vanishes by quasi-neutrality. A similar cancellation occurs for the leading-order term in $\unicode[STIX]{x1D650}^{(xy)}$ .) Taking these subtleties into account, the dispersion relation to leading order in $\unicode[STIX]{x1D716}$ is
Since $\unicode[STIX]{x1D63F}^{(xx)}$ is in general non-zero, the term in brackets must vanish, yielding
where we have used the fact that the ordering (5.34) precludes Landau resonances (i.e. $\unicode[STIX]{x1D701}_{s}\sim O(\sqrt{\unicode[STIX]{x1D716}})$ ). The first term in parentheses on the left-hand size of (5.37) corresponds to Alfvén waves ( $\unicode[STIX]{x1D714}=\pm k_{\Vert }v_{A}$ ) whose $\unicode[STIX]{x1D6FF}\boldsymbol{B}$ is polarized in the $y$ direction. The second term in parentheses is readily identifiable as the pseudo-Alfvén-wave branch, which, in the absence of stabilizing magnetic tension, is unstable to the MTI for $\text{d}\ln T/\text{d}z<0$ . Its magnetic-field fluctuations are polarized in the $z$ direction, thus enabling rapid field-aligned conduction for vertical displacements. As in §§ 5.1.3 and 5.1.4, these two branches are coupled by a non-zero right-hand side, but here the cause is ‘gyroviscosity’ (Kaufman Reference Kaufman1960; Braginskii Reference Braginskii1965; Ramos Reference Ramos2005) – an FLR-induced dissipationless cross-field transport of momentum that arises from Larmor-scale spatial variations in the guiding-centre $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ drifts. It is dispersive and stabilizing, an effect that can be seen clearly in figure 5(b), which shows the wave frequency versus parallel wavenumber for a variety of $\unicode[STIX]{x1D6FD}_{i}$ and thermal-pressure scale heights. For $\unicode[STIX]{x1D6FD}_{i}\gtrsim H/\unicode[STIX]{x1D70C}_{i}$ , the gyroviscosity is the dominant stabilizing effect (over the magnetic tension), providing a wavenumber cutoff
where $d_{i}\doteq (m_{i}c^{2}/4\unicode[STIX]{x03C0}e^{2}n_{i})^{1/2}=\unicode[STIX]{x1D70C}_{i}/\sqrt{\unicode[STIX]{x1D6FD}_{i}}$ is the ion skin depth. For $n_{i}=10^{-3}~\text{cm}^{-3}$ , $H=500~\text{kpc}$ and $\text{d}\ln T/\text{d}\ln P=1/3$ , (5.38) gives $(k_{\Vert }H)_{max}\approx 5\times 10^{4}\,(10^{12}/\unicode[STIX]{x1D6FD}_{i})^{1/4}$ .
The physics driving the gyroviscous stabilization of the MTI may be elucidated by examining the leading-order perturbed distribution function (4.25) in the frame of the perturbed $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ flow under the gyroviscous ordering (5.34):
The final term in this expression makes the perturbed distribution function $\unicode[STIX]{x1D6FF}f_{s}(v_{\Vert },w_{\bot })$ gyrophase dependent (a feature absent in the drift-kinetic limit), which endows the pressure tensor with off-diagonal elements:
Substituting this expression into the perturbed momentum equation and taking the limit $m_{e}/m_{i}\rightarrow 0$ , we find that the Lagrangian displacement perpendicular to the mean magnetic field $\unicode[STIX]{x1D743}_{\bot }$ evolves according to
The first term on the right-hand side – the buoyancy term – drives vertical displacements, which would grow exponentially (provided $\text{d}\ln T/\text{d}z<0$ ) were it not for the magnetic tension and the final (gyroviscous) term on the right-hand side. The latter rotates a vertical displacement into the horizontal plane (i.e. vertical momentum is transported in the cross-field $y$ direction). If $k_{\Vert }^{2}\unicode[STIX]{x1D70C}_{i}H\gtrsim 1$ , this rotation occurs on a faster time scale than buoyancy can grow the displacement and the MTI becomes stabilised.
A final note is perhaps in order regarding the gyroviscosity and the contributions to it that sometimes arise in the presence of temperature gradients, the so-called Mikhailovskii & Tsypin (Reference Mikhailovskii and Tsypin1971, Reference Mikhailovskii and Tsypin1984) terms (see also Catto & Simakov Reference Catto and Simakov2004). For an isotropic equilibrium distribution function, the gyroviscous contribution to the pressure tensor of species $s$ is
where the superscript $\text{T}$ denotes the matrix transpose and $\boldsymbol{Q}_{\bot s}=Q_{\bot s}\hat{\boldsymbol{b}}$ is the parallel flux of perpendicular heat (viz. $Q_{\bot s}=m_{s}\int \text{d}^{3}\boldsymbol{v}\,v_{\Vert }(v_{\bot }^{2}/2)f_{s}$ ). (We refer the reader to equation (6) and appendix A of Schekochihin et al. (Reference Schekochihin, Cowley, Rincon and Rosin2010) for a simple derivation of (5.42)). The $\unicode[STIX]{x1D735}\boldsymbol{u}_{s}$ terms, first obtained (in the collisional limit) by Braginskii (Reference Braginskii1965), simplify considerably for our problem. With $\boldsymbol{k}=k_{\Vert }\hat{\boldsymbol{b}}$ , they are $\text{i}k_{\Vert }(P_{s}/\unicode[STIX]{x1D6FA}_{s})(\hat{\boldsymbol{e}}_{x}\boldsymbol{\times }\unicode[STIX]{x1D6FF}\boldsymbol{u}_{\bot }\hat{\boldsymbol{e}}_{x}-\hat{\boldsymbol{e}}_{x}\hat{\boldsymbol{e}}_{x}\boldsymbol{\times }\unicode[STIX]{x1D6FF}\boldsymbol{u}_{\bot })$ , which, using flux freezing (5.12), returns precisely (5.40). The gyroviscous stress associated with these terms has a non-zero divergence, and therefore exerts a force on the plasma (the final term in (5.41)). The $\unicode[STIX]{x1D735}\boldsymbol{Q}_{\bot s}$ terms, on the other hand, vanish identically for $\boldsymbol{k}=k_{\Vert }\hat{\boldsymbol{b}}$ , and so heat-flux contributions to the gyroviscosity do not affect the modes investigated in this section. For arbitrary wavevector orientation, however, these terms become
where $\unicode[STIX]{x1D6FF}Q_{\bot s}=m_{s}\int \text{d}^{3}\boldsymbol{v}v_{\Vert }(v_{\bot }^{2}/2)\unicode[STIX]{x1D6FF}f_{s}$ . Its divergence is $k_{\Vert }\boldsymbol{k}_{\bot }(\unicode[STIX]{x1D6FF}Q_{\bot s}/\unicode[STIX]{x1D6FA}_{s})\boldsymbol{\times }\hat{\boldsymbol{e}}_{x}$ , and so the perturbed heat flux exerts a cross-field force on the plasma that is larger at smaller scales. One can show by directly computing
and taking the relevant limit $\unicode[STIX]{x1D701}_{s}\ll 1$ , that this force is only important at perpendicular scales of order the ion-Larmor scale. This is precisely our next destination.
5.3 Gyrokinetic limit: $k_{\Vert }\unicode[STIX]{x1D70C}_{s}\ll k_{\bot }\unicode[STIX]{x1D70C}_{s}\sim 1$
We have found that proceeding further to yet smaller parallel scales (e.g. $k_{\Vert }\unicode[STIX]{x1D70C}_{i}\sim 1$ ) is not a particularly fruitful venture, and so we now turn our attention to fluctuations with small (i.e. kinetic-scale) perpendicular wavelengths. To do so in an analytically tractable way, we adopt the gyrokinetic ordering (Rutherford & Frieman Reference Rutherford and Frieman1968; Taylor & Hastie Reference Taylor and Hastie1968; Antonsen & Lane Reference Antonsen and Lane1980; Catto, Tang & Baldwin Reference Catto, Tang and Baldwin1981; Frieman & Chen Reference Frieman and Chen1982):
which selects low-frequency fluctuations that are highly oblique with respect to the magnetic-field direction. Note that, in this ordering, the perpendicular length scale of the fluctuations $k_{\bot }^{-1}$ satisfies $k_{\bot }\unicode[STIX]{x1D70C}_{s}\sim k_{\Vert }H\sim 1$ , and so the fluctuations are permitted to have (perpendicular) scales of the order of the Larmor radius. As a result, $k_{y}v_{ds}\sim k_{\Vert }v_{ths}\sim \unicode[STIX]{x1D714}$ , i.e. drift waves (and the various instabilities typically associated with them) are allowed. In common with the long-wavelength limit taken in § 5.1, collisionless damping via the Landau resonance is captured (but cyclotron resonances are not).
The procedure for deriving the dispersion relation in the gyrokinetic limit is similar to that employed in the long-wavelength and gyroviscous limits: we apply the ordering (5.45) to the dispersion tensor (4.24), retain its leading-order terms, and set its determinant to zero. Unfortunately, the resulting equation is rather formidable, requiring myriad algebraic gymnastics to massage into a parsable form. Instead, we have found that a better route – one suggested by previous work on linear and nonlinear gyrokinetic theory – is to use (3.10a,b ) and (4.5) to shift our independent variables from the field components $(\unicode[STIX]{x1D6FF}E_{\Vert },\unicode[STIX]{x1D6FF}E_{y},\unicode[STIX]{x1D6FF}E_{z})$ to the potentials $(\unicode[STIX]{x1D711},A_{\Vert },\unicode[STIX]{x1D6FF}B_{\Vert })$ and then variously combine the three rows of the transformed dispersion tensor to obtain as our new rows what amounts to the quasi-neutrality constraint, the parallel component of Ampère’s law, and the perpendicular component of Ampère’s law (which, in gyrokinetics, is equivalent to a statement of perpendicular pressure balance). The determinant, and thus the dispersion relation, is of course identical following this profitable change of basis.
The technical details of these manipulations are not given here; instead, we provide in appendix C an alternate derivation of the linear gyrokinetic theory starting from a derivation of the nonlinear gyrokinetic theory of a thermally stratified atmosphere. Apart from some tedious bookkeeping and adroit consolidation of terms, the essential ingredients of the calculation are an expansion of the $Z_{p}(\unicode[STIX]{x1D701}_{n,s})$ functions about $\unicode[STIX]{x1D701}_{0,s}=(\unicode[STIX]{x1D714}+k_{y}v_{ds})/k_{\Vert }v_{ths}$ and several uses of the identity $\sum _{n=1}^{\infty }2\unicode[STIX]{x1D6E4}_{n}(\unicode[STIX]{x1D6FC}_{s})=1-\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{s})$ . After much effort, (4.11) becomes
where $\overline{\unicode[STIX]{x1D714}}_{i}\doteq \unicode[STIX]{x1D714}/k_{\Vert }v_{Ai}$ and we have employed the shorthand notation (cf. § 2.6 of Howes et al. Reference Howes, Cowley, Dorland, Hammett, Quataert and Schekochihin2006)
are the temperature-gradient and diamagnetic drift frequencies of species $s$ .Footnote 4
Setting the determinant of the $3\times 3$ matrix in (5.46) to zero gives the dispersion relation in the gyrokinetic limit:
The left-hand side of (5.49) contains two factors, the first corresponding to the Alfvén-wave branch and the second to the slow-wave branch. The right-hand side represents coupling of the two branches that occurs either at finite-Larmor radius (see § 2.6 of Howes et al. (Reference Howes, Cowley, Dorland, Hammett, Quataert and Schekochihin2006) for a thorough discussion of this effect in a homogeneous plasma) or in pressure-stratified plasmas when $k_{y}\neq 0$ . The latter has already been discussed as part of the long-wavelength limit in § 5.1.3. We now examine the short-wavelength regime. As for the previous sections, having the leading-order expression for the perturbed distribution function will serve useful:
where $\unicode[STIX]{x1D712}\doteq \unicode[STIX]{x1D711}-v_{\Vert }A_{\Vert }/c-\boldsymbol{v}_{\bot }\boldsymbol{\cdot }\boldsymbol{A}_{\bot }/c$ is the gyrokinetic potential (see (C 5)); its ring average $\langle \unicode[STIX]{x1D712}\rangle _{\boldsymbol{R}_{s}}$ is given by (C 29). Equation (5.50) may be obtained by applying the gyrokinetic ordering to (4.25) or, alternatively, by consulting appendix C.
In figure 6, we plot solutions of (5.49) for a hydrogenic plasma with $\unicode[STIX]{x1D6FD}_{i}=10^{4}$ , $\unicode[STIX]{x1D713}=0$ (i.e. $k_{z}=0$ ), and $\text{d}\ln T/\text{d}\ln P=1/3$ . At $k_{y}\unicode[STIX]{x1D70C}_{i}\ll 1$ , the growth rate remains roughly constant with $k_{y}\unicode[STIX]{x1D70C}_{i}$ , an effect discussed in § 5.1.3. As $k_{y}\unicode[STIX]{x1D70C}_{i}$ gets larger, the growth rate decreases and $\unicode[STIX]{x1D714}$ acquires a real part comparable to its imaginary part. The decrease in growth rate may be obtained from (5.49) by expanding in $\unicode[STIX]{x1D701}_{i},\unicode[STIX]{x1D6FC}_{i}\ll 1$ and $\unicode[STIX]{x1D701}_{e},\unicode[STIX]{x1D6FC}_{e}\rightarrow 0$ ; in this limit, ${\mathcal{B}}\simeq \unicode[STIX]{x1D6FC}_{i}$ , ${\mathcal{E}}\simeq -(3/2)\unicode[STIX]{x1D6FC}_{i}$ , and ${\mathcal{F}}\simeq (\unicode[STIX]{x1D714}_{Pi}\unicode[STIX]{x1D714}_{Ti}/\unicode[STIX]{x1D714}^{2})(1+T_{e}/T_{i})$ , and so (5.49) becomes
where we have employed the shorthand $\unicode[STIX]{x1D6EC}\doteq 2/\unicode[STIX]{x1D6FD}_{i}-{\mathcal{D}}+{\mathcal{C}}^{2}/{\mathcal{A}}$ to represent the slow-mode piece of (5.49). The first term in parentheses is readily identified as the Alfvénic MTI (cf. (2.3)). The second term in (5.51) captures the FLR stabilization, with the growth rate eventually falling off as ${\sim}(k_{\bot }\unicode[STIX]{x1D70C}_{i})^{-1}$ . As $k_{\bot }\unicode[STIX]{x1D70C}_{i}$ approaches unity, the Alfvénic MTI becomes more and more coupled to the collisionlessly damped slow-mode branch, and the dispersion relation becomes approximately $({\mathcal{D}}-2/\unicode[STIX]{x1D6FD}_{i}){\mathcal{F}}\simeq {\mathcal{E}}^{2}$ – a weakly Barnes-damped slow mode (i.e. $\unicode[STIX]{x1D714}\simeq -\text{i}k_{\Vert }v_{A}/\sqrt{\unicode[STIX]{x03C0}\unicode[STIX]{x1D6FD}_{i}}$ ) coupled by FLR effects ( ${\mathcal{E}}^{2}$ ) to a weakly unstable Alfvénic-MTI mode ( ${\mathcal{F}}$ ). This is a result of the ions gyro-averaging over Larmor-scale fluctuations and thus reducing the connection between the electromagnetic fluctuations and the ion distribution function. The destabilised Alfvén wave, for example, loses its character altogether, as the ions drift across the magnetic-field lines and the mode becomes dispersive.
On the other side of $k_{\bot }\unicode[STIX]{x1D70C}_{i}\sim 1$ , the growth rate begins to increase, eventually attaining values comparable to (and even exceeding!) those achieved at long wavelengths. This trend continues unabated until kinetic effects intervene at the electron-Larmor scale, where the wave is ultimately damped. To uncover the physics driving the sub-ion-Larmor-scale growth, we expand (5.49) in the region $k_{\bot }\unicode[STIX]{x1D70C}_{e}\ll 1\ll k_{\bot }\unicode[STIX]{x1D70C}_{i}$ , in which $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{i})$ and its derivatives asymptote to zero and $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{e})\simeq \unicode[STIX]{x1D6E4}_{1}(\unicode[STIX]{x1D6FC}_{e})\simeq 1$ , whence ${\mathcal{A}}\simeq 1+T_{i}/T_{e}$ , ${\mathcal{B}}\simeq 1-(T_{i}/T_{e})(\unicode[STIX]{x1D714}_{Te}/\unicode[STIX]{x1D714})$ , ${\mathcal{C}}\simeq 0$ , ${\mathcal{E}}\simeq -1$ and ${\mathcal{F}}\simeq (T_{i}/T_{e})(\unicode[STIX]{x1D714}_{Pe}/\unicode[STIX]{x1D714})[1+(\unicode[STIX]{x1D714}_{Te}/\unicode[STIX]{x1D714})]$ . As we are generally concerned with $\unicode[STIX]{x1D6FD}_{i}\gg 1$ , the coefficient ${\mathcal{D}}$ warrants some care. Usually, in such sub-ion-Larmor expansions, $\unicode[STIX]{x1D6FD}_{i}$ is taken to be order unity relative to the ion-to-electron mass ratio, and ${\mathcal{D}}$ can be safely taken to be $\simeq 0$ . Specifically, the ion contribution to ${\mathcal{D}}$ scales with $\unicode[STIX]{x1D6E4}_{0}^{\prime }(\unicode[STIX]{x1D6FC}_{i})\sim -1/\sqrt{8\unicode[STIX]{x03C0}\unicode[STIX]{x1D6FC}_{i}^{3}}$ , which is small, and the electron contribution scales with $(m_{e}/m_{i})^{1/2}$ , which is small as well. Indeed, this expansion is what, in a homogeneous plasma, leads to the usual dispersion relation describing kinetic-Alfvén waves (e.g. Kingsep, Chukbar & Yan’kov Reference Kingsep, Chukbar and Yan’kov1990):
But, with our interest in high- $\unicode[STIX]{x1D6FD}$ plasmas, $\unicode[STIX]{x1D6FD}_{i}$ can be ${\sim}m_{i}/m_{e}$ or perhaps even more, and the electron contribution to ${\mathcal{D}}$ may be ${\gtrsim}2/\unicode[STIX]{x1D6FD}_{i}$ (its main competition in (5.49)). For that reason, the mass-ratio expansion involved in taking $k_{\bot }\unicode[STIX]{x1D70C}_{e}\ll 1$ cannot be performed independently of taking $\unicode[STIX]{x1D6FD}_{i}\gg 1$ .
These complications taken into account, we retain both $1/\unicode[STIX]{x1D6FD}_{i}$ and ${\mathcal{D}}$ in (5.49) while adopting the asymptotic values of ${\mathcal{A}}$ , ${\mathcal{B}}$ , ${\mathcal{C}}$ , ${\mathcal{E}}$ and ${\mathcal{F}}$ given in the preceding paragraph. Equation (5.49) becomes
or
For $\unicode[STIX]{x1D6FC}_{i}\gg (m_{i}/m_{e})^{1/3}$ – a limit well satisfied in the regime of interest – ${\mathcal{D}}$ is dominated by its electron contribution $\simeq 2\text{i}\sqrt{\unicode[STIX]{x03C0}}(T_{e}/T_{i})(\unicode[STIX]{x1D714}/k_{\Vert }v_{the})[1-(\unicode[STIX]{x1D714}_{Te}/2\unicode[STIX]{x1D714})]$ . The first term in parentheses on the left-hand side of (5.54) then reduces to either $2/(1+T_{e}/T_{i})$ if $\unicode[STIX]{x1D6FD}_{i}\ll 1$ , $\unicode[STIX]{x1D6FD}_{i}$ if $1\ll \unicode[STIX]{x1D6FD}_{i}\ll (m_{i}/m_{e})^{1/2}$ , or $2\text{i}(T_{i}/T_{e})(k_{\Vert }v_{the}/\sqrt{\unicode[STIX]{x03C0}})/(2\unicode[STIX]{x1D714}-\unicode[STIX]{x1D714}_{Te})$ if $\unicode[STIX]{x1D6FD}_{i}\gg (m_{i}/m_{e})^{1/2}\gg 1$ . (Here, we have taken $k_{y}\unicode[STIX]{x1D70C}_{i}\sim k_{\Vert }H$ with respect to the subsidiary expansion in $\unicode[STIX]{x1D6FD}_{i}$ and mass ratio.) In all cases, a necessary (but not sufficient) condition for instability is
or, equivalently,
where $c_{s}\doteq (T_{e}/m_{i})^{1/2}$ is the ion-acoustic speed. With $\text{d}\ln P_{e}/\text{d}z<0$ , we require $\text{d}\ln T_{e}/\text{d}z<0$ (as in the drift-kinetic and gyroviscous limits). This is a kinetic-Alfvén drift wave, destabilised by rapid electron conduction along perturbed magnetic-field lines; we name this the electron magnetothermal instability, or eMTI, due to its reliance on the electron-temperature gradient and the relatively passive role played by the ion species. Indeed, at these scales, the ion response is nearly Boltzmann, $\unicode[STIX]{x1D6FF}f_{i}\approx -(e\unicode[STIX]{x1D711}/T_{i})f_{i}$ (cf. (5.50)).
We now obtain some approximate analytic solutions to (5.54). In the high- $\unicode[STIX]{x1D6FD}_{i}$ limit, equation (5.54) becomes
where $\unicode[STIX]{x1D70C}_{sound}\doteq c_{s}/\unicode[STIX]{x1D6FA}_{i}$ is the ion-sound radius. We may expand this solution in the parameter $(k_{y}\unicode[STIX]{x1D70C}_{e}/k_{\Vert }H)^{1/2}\sim (m_{e}/m_{i})^{1/4}\ll 1$ to find the leading-order expression
i.e. the mode is overstable, with a growth rate ${\sim}k_{y}\unicode[STIX]{x1D70C}_{sound}(c_{s}/H)\sqrt{k_{y}\unicode[STIX]{x1D70C}_{e}/k_{\Vert }H}$ .Footnote 5 Despite all the approximations, the growth rate agrees rather well with the numerical solution (see figure 6). The real part of (5.58) differs from that seen in the figure, which is not too surprising – $(m_{e}/m_{i})^{1/4}$ is not that small of a number! – but the accuracy of the solution improves as $m_{e}/m_{i}\rightarrow 0$ . This is shown in figure 7, which gives the solution to (5.49) for $m_{e}/m_{i}=1/183\,600$ . The blue dotted line, equation (5.58), is a clearly an excellent fit to both the real and imaginary parts of $\unicode[STIX]{x1D714}$ . On the other hand, the solution to (5.54) for $1\ll \unicode[STIX]{x1D6FD}_{i}\ll (m_{i}/m_{e})^{1/2}$ is approximately
In this regime, the kinetic-Alfvénic nature of the mode is readily apparent from the final term in the square root.
In all of these limits, the maximum growth rate of the eMTI is set by electron FLR effects (the mode is damped at $k_{y}\unicode[STIX]{x1D70C}_{e}\sim 1$ – see figure 6 b) and is generally larger than the growth rate for the (long-wavelength) kinetic MTI. This is shown explicitly in figure 8, which gives the respective maximum growth rates $\unicode[STIX]{x1D6FE}_{max}$ evaluated across the full range of $k_{\Vert }$ and $k_{y}$ ( $k_{z}=0$ ) for fixed temperature gradient and ion-to-electron temperature ratio. The difference is probably related to the relatively mundane role that the inertia-bearing ions play in the eMTI relative to the MTI (recall that the ion response is Boltzmann for the eMTI). It is the electrons which behave similarly in both cases: indeed, (5.50) with $(k_{y}\unicode[STIX]{x1D70C}_{e})^{2}\ll 1$ implies
which states that the Lagrangian change in a fluid element’s perpendicular electron temperature is nearly zero as the flux-frozen electrons carry the magnetic field upwards or downwards. While there appears to be an asymptotic eMTI maximum growth rate as $\unicode[STIX]{x1D6FD}_{i}\rightarrow \infty$ , which is roughly twice as large as the maximum MTI growth rate, the complicated dependence of $k_{\Vert }H$ and $k_{y}\unicode[STIX]{x1D70C}_{e}$ on $\unicode[STIX]{x1D6FD}_{i}$ at maximum eMTI growth precludes a straightforward analytic calculation for it.
We close this section by commenting on the transition between the drift-kinetic solution for $k_{y},k_{z}\neq 0$ explored in § 5.1.4 – see figure 4 in particular – and the gyrokinetic solution. Recall from that section that the growth rate at long wavelengths is maximal and roughly constant for $k_{y}/k_{z}\lesssim 1$ . In figure 9 we see that this trend continues to small wavelengths until $k_{y}\unicode[STIX]{x1D70C}_{i}\sim k_{z}\unicode[STIX]{x1D70C}_{i}\sim 1$ , at which point ion FLR effects sharply reduce the MTI growth rate. Thereafter, the eMTI emerges with an even larger growth rate, which is independent of $k_{z}$ until $k_{z}\unicode[STIX]{x1D70C}_{e}\sim 1$ . This independence is not particularly surprising: for $k_{\bot }^{2}\unicode[STIX]{x1D70C}_{i}^{2}\gg 1$ the gyrokinetic dispersion relation (5.49) only depends upon $k_{z}$ through $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{e})$ and its derivatives, the drift frequency depending only upon $k_{y}$ . As with the long-wavelength solution (5.29), the mirror force, by which Landau-resonant particles bleed energy from the magnetic-field-strength fluctuations, is offset by the divergence of the parallel pressure.
5.4 Relationship between the magnetothermal and the ion- and electron-temperature-gradient instabilities
At this point in the manuscript, those readers familiar with drift-wave instabilities in magnetic-confinement-fusion devices might be curious as to the relationship between the various incarnations of the magnetothermal instability elucidated above and the ion- and electron-temperature-gradient instabilities long known to plague fusion plasmas (the so-called ‘ITG’ and ‘ETG’). Indeed, both the MTI and ITG/ETG are convective instabilities driven by the free energy stored in the background temperature gradient. Their relationship is complicated somewhat by the usual restriction of ITG/ETG analyses to the electrostatic $\unicode[STIX]{x1D6FD}\rightarrow 0$ limit, in which the magnetic field is assumed rigid, and the geometrical complexity of tokamak magnetic fields, which introduces its own drift-wave instabilities apart from the ‘slab’ ITG/ETG instabilities (e.g. toroidal ITG/ETG, trapped electron modes).Footnote 6 Nevertheless, the transition between the ITG/ETG, on the one hand, and the MTI, on the other, as the plasma $\unicode[STIX]{x1D6FD}$ is increased, is an interesting question. In this section, we explore their connection.
Following conventional work on the slab ITG/ETG (Rudakov & Sagdeev Reference Rudakov and Sagdeev1961; Cowley et al. Reference Cowley, Kulsrud and Sudan1991), we begin by taking the electrostatic limit, in which the electric field $\boldsymbol{E}=-\unicode[STIX]{x1D735}\unicode[STIX]{x1D711}$ . The relevant dispersion relation is obtained by substituting the perturbed distribution function (4.3) with $\unicode[STIX]{x1D6FF}\boldsymbol{B}=0$ into the quasi-neutrality constraint (3.2), adopting plane-wave solutions, performing the resulting integrals (after conveniently setting $F_{s}=F_{M,s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ – see (4.16)) and demanding non-trivial solutions. The result is that
In the low-frequency limit, we need only retain the $n=0$ contributions to the summation above, and (5.61) reduces to ${\mathcal{A}}=0$ (cf. (5.47a )). We adopt this limit throughout the remainder of this section.
Explicit solutions of the dispersion relation ${\mathcal{A}}=0$ may be calculated in the analytically tractable limit $k_{\Vert }v_{the}\gg \unicode[STIX]{x1D714},\unicode[STIX]{x1D714}_{Pi}>k_{\Vert }v_{thi}$ , whence $Z_{0}(\unicode[STIX]{x1D701}_{0,i})$ and $Z_{2}(\unicode[STIX]{x1D701}_{0,i})$ can be expanded in their large argument. The result is
We now make some standard simplifying assumptions, which are known to be quantitatively imprecise but nevertheless afford a compact analytic solution whose qualitative physics may be readily deduced. Taking $\text{d}\ln n_{i}/\text{d}\ln T_{i}\ll (k_{\Vert }H/k_{y}\unicode[STIX]{x1D70C}_{i})^{2/3}\ll (T_{i}/T_{e})^{1/3}$ , we find that the entire right-hand side of (5.62) may be dropped. Expanding $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{i})\simeq -\unicode[STIX]{x1D6E4}_{0}^{\prime }(\unicode[STIX]{x1D6FC}_{i})\simeq 1$ , we find that
which is unstable for $\unicode[STIX]{x1D714}_{Ti}>0$ or, equivalently, $\boldsymbol{k}\boldsymbol{\cdot }(\hat{\boldsymbol{b}}\boldsymbol{\times }\unicode[STIX]{x1D735}T_{i})>0$ . This is the ITG instability, Doppler shifted by the equilibrium diamagnetic drift. The physics of the instability is explained in detail by Cowley et al. (Reference Cowley, Kulsrud and Sudan1991), to which we refer the reader. In brief, parallel flow along fields lines compresses the ions. An electrostatic potential is produced in order to enforce quasi-neutrality (the electron response is Boltzmann), which drives an $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ flow that advects cool ions into the compressed region. This lowers the pressure there, pulling in yet more ions along the field lines and reinforcing the original density enhancement. The result is a feedback loop that leads to instability. The key requirement here is for the mode frequency $\unicode[STIX]{x1D714}$ to be larger than the rate of pressure relaxation along field lines $k_{\Vert }c_{s}$ , so that the Lagrangian change in an ion fluid element’s pressure vanishes to leading order. Other requirements are that the density gradient not be too strong and $k_{y}\unicode[STIX]{x1D70C}_{i}$ not be too large.
These stabilizing effects can be seen in figure 10(a), which shows the growth rate of the electrostatic ITG in the $k_{\Vert }H$ – $k_{y}\unicode[STIX]{x1D70C}_{i}$ plane obtained from numerically solving (5.61) with $m_{e}/m_{i}=0$ (so that ETG is suppressed – see below). The dashed line shows the marginal stability threshold, which may be obtained analytically by setting $\text{Im}(\unicode[STIX]{x1D714})=0$ in the dispersion relation ${\mathcal{A}}=0$ , demanding that $\text{Im}({\mathcal{A}})=0$ (which determines the real frequency of the mode at marginal stability) and substituting the result back into ${\mathcal{A}}=0$ :
(Here we have taken $\text{Re}(\unicode[STIX]{x1D714})/k_{\Vert }v_{the}\ll 1$ .) Consider the numerator inside the square root on the right-hand side of (5.64), which must be positive for a sensible solution. This shows the standard result that, at long wavelengths (i.e. $\unicode[STIX]{x1D6FC}_{i}\rightarrow 0$ ), the density gradient can be no steeper than $\text{d}\ln n_{i}/\text{d}\ln T_{i}=1/2$ , or $\text{d}\ln T_{i}/\text{d}\ln P_{i}>2/3$ for instability. (In the fusion literature, this is expressed as $\unicode[STIX]{x1D702}_{i}\doteq \text{d}\ln T_{i}/\text{d}\ln n_{i}>2$ .) FLR effects change this slightly: using the asymptotic result $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{i})\approx 1/\sqrt{2\unicode[STIX]{x03C0}\unicode[STIX]{x1D6FC}_{i}}$ in (5.64) gives a critical density gradient $\text{d}\ln n_{i}/\text{d}\ln T_{i}\approx 1$ , or $\text{d}\ln T_{i}/\text{d}\ln P_{i}\gtrsim 1/2$ . (Note that our fiducial MTI-unstable atmosphere with $\text{d}\ln T_{i}/\text{d}\ln P_{i}=1/3$ is stable to ITG/ETG.)
Equilibrium gradients in the electron species also give rise to instability. A calculation analogous to the ITG may be performed at sub-ion-Larmor scales ( $\unicode[STIX]{x1D6FC}_{i}\gg 1$ ) by assuming $\unicode[STIX]{x1D714}/k_{\Vert }v_{the}>1$ and setting $\unicode[STIX]{x1D6E4}_{0}(\unicode[STIX]{x1D6FC}_{i})=\unicode[STIX]{x1D6E4}_{0}^{\prime }(\unicode[STIX]{x1D6FC}_{i})=0$ . The dispersion relation is identical to (5.62), but with the ion and electron subscripts reversed. Accordingly, the electron counterpart to (5.63) is
which is unstable for $\unicode[STIX]{x1D714}_{Te}>0$ or, equivalently, $\boldsymbol{k}\boldsymbol{\cdot }(\hat{\boldsymbol{b}}\boldsymbol{\times }\unicode[STIX]{x1D735}T_{e})<0$ . This is the ETG instability. The numerical solution to (5.62) for $m_{i}/m_{e}=1836$ is given in figure 10(b), which shows the ETG peaking at $k_{\Vert }H\sim 0.1$ and $k_{y}\unicode[STIX]{x1D70C}_{e}\sim 1$ .
The common characteristic of both the ITG and ETG instabilities is that the unstable mode drifts across stationary magnetic-field lines. This is allowed because, at the low values of $\unicode[STIX]{x1D6FD}$ required for the electrostatic approximation to hold, the electron skin depth is at least as large as the ion-Larmor radius, and so the ion (and electron) species can slip through the otherwise rigid field lines. By contrast, at large $\unicode[STIX]{x1D6FD}$ , the magnetic field is unfrozen only by electron FLR effects, and so a super-Larmor-scale displacement in the plasma carries the magnetic-field lines with it. It is this displacement of field lines, and the accompanying rapid streaming of particles along it, that both stabilises the plasma to ITG/ETG and destabilises it to the MTI/eMTI. While the physics is markedly different, the end result is the same: the free energy stored in the temperature gradient drives convective instability both at long wavelengths and, more vigorously, at Larmor scales.
6 Summary
Efforts to provide a proper accounting of kinetic effects in high- $\unicode[STIX]{x1D6FD}$ astrophysical plasmas are still in their infancy, particularly when it comes to the influence of plasma microphysics on the macroscale dynamics – convection being but one example. In this paper we have studied the linear stability properties of a collisionless, magnetised, thermally stratified atmosphere. This is an extension of previous work on collisional, magnetised, stratified atmospheres (Balbus Reference Balbus2000, Reference Balbus2001; Quataert Reference Quataert2008; Kunz Reference Kunz2011), which had shown them to be linearly unstable when the temperature increases in the direction of gravity. This instability – the magnetothermal instability (MTI) – is facilitated by rapid, field-aligned conduction obviating adiabatic fluid displacements. In this fast-conduction limit, slow-mode perturbations exhibit almost no (Lagrangian) change in their temperature; for $\text{d}\ln T/\text{d}z<0$ , it is thus energetically favourable for vertically displaced fluid elements to continue rising or falling in the gravitational potential. In certain situations, the accompanying field-aligned viscous momentum transport can drive Alfvénic fluctuations unstable as well.
This paper focuses instead on the case in which the heat and momentum transport that enables these instabilities is not due to collisional processes, but rather stems from the complex interplay between field-aligned particle streaming, collisionless wave damping, and finite-Larmor-radius effects. Because prior work on the MTI relied on the buoyant ion species being thermally equilibrated with the inertialess parallel-conducting electrons, the application of those results to a collisionless system is not by any means straightforward. To what extent they do carry over to the collisionless case is a particularly intriguing question in the context of the hot, diffuse plasma in the outskirts of galaxy clusters – the intracluster medium – in which the collisional mean free path can approach the thermal-pressure scale height. Radiatively inefficient black-hole accretion flows, such as that onto Sgr A $^{\ast }$ at the Galactic centre, are also a natural venue for application of our results, although the situation there is further complicated by differential rotation and the consequent kinetic magnetorotational instability (Quataert, Dorland & Hammett Reference Quataert, Dorland and Hammett2002; Sharma, Hammett & Quataert Reference Sharma, Hammett and Quataert2003; Heinemann & Quataert Reference Heinemann and Quataert2014; Islam Reference Islam2014; Quataert, Heinemann & Spitkovsky Reference Quataert, Heinemann and Spitkovsky2015).
Our results may be summarized as follows. First, at long wavelengths, the fluid MTI carries over relatively unscathed to the kinetic case, enjoying the same stability criterion ( $\text{d}\ln T/\text{d}z>0$ ) and maximum growth rate. The physics driving the requisite isothermalisation of perturbed field lines is, of course, quite different in each case, but that appears to be of little concern to buoyancy, which does its job regardless. That being said, it is rather important to note that conduction and buoyancy in a collisionless plasma are inextricably linked, a feature not present in a collisional fluid. In a collisionless plasma, approximately isothermal displacements are tandem with the ability of those displacements to maintain pressure balance with their surroundings. Both processes suffer from collisionless damping if an appreciable number of particles are Landau resonant with the fluctuation.
Second, because of adiabatic invariance, changes in temperature go hand-in-hand with changes in magnetic-field strength. Since the latter are collisionlessly damped, this makes the physics of the MTI much richer than in the collisional case, for which $\unicode[STIX]{x0394}T\simeq 0$ in the fast-conduction limit (provided $\boldsymbol{B}_{0}\bot \boldsymbol{g}$ ). One consequence is that fluctuations whose wavevectors have a component along gravity ( $k_{z}\neq 0$ ) are strongly damped once $k_{z}H\gtrsim \sqrt{k_{\Vert }H}$ . This damping can be mitigated by allowing $\boldsymbol{k}\boldsymbol{\cdot }(\boldsymbol{g}\boldsymbol{\times }\boldsymbol{B}_{0})\neq 0$ , which introduces an additional degree of freedom for the magnetic-field perturbations. This freedom allows the magnetic field to seek out favourable arrangements, in which field-strength fluctuations are anti-correlated with upward (downward) displacements into regions of greater (lesser) potential. The result is that the mirror force acts oppositely to the parallel gravitational force, and so the collisionless damping of field-strength fluctuations is reduced (see figure 3).
Third, for extremely sub-thermal magnetic fields characteristic of those found in primordial galaxy clusters and used in some published (fluid) simulations of the MTI, gyroviscosity – not magnetic tension – sets the (parallel) wavenumber cutoff. It does so by providing a dissipationless cross-field transport of momentum, which redirects buoyant motions unfavourably into the horizontal plane.
Fourth, and perhaps most importantly, we have shown that the plasma is overstable at sub-ion-Larmor scales ( $k_{\bot }\unicode[STIX]{x1D70C}_{i}\gg 1\gg k_{\bot }\unicode[STIX]{x1D70C}_{e}$ ) to a kinetic-Alfvén drift-wave instability – the electron MTI, or eMTI – which exhibits a growth rate even larger than does the long-wavelength MTI. It is driven by the electron-temperature gradient and relies on (approximately) Boltzmann ions and isothermal-along-perturbed-field-lines electrons. To give this last result some context, we briefly reviewed the slab ITG and ETG – two of the drift-wave instabilities long known to the magnetic-confinement-fusion community for driving the anomalous transport limiting tokamak performance. While our focus throughout the manuscript has been on applications to the high- $\unicode[STIX]{x1D6FD}$ intracluster medium, these drift-wave instabilities may also find applications in the hot and tenuous coronae of stars.
As a closing thought, it is tempting to speculate on the nonlinear evolution in a high- $\unicode[STIX]{x1D6FD}$ plasma of the instabilities investigated herein. One cannot go too far in extrapolating, of course, if only because the velocity-space anisotropy generated by adiabatic invariance as the instabilities grow will quickly drive secondary Larmor-scale instabilities (firehose, mirror), which have been shown to be efficient at isotropising the distribution function and thus providing an effective ‘viscosity’ (Kunz, Schekochihin & Stone (Reference Kunz, Schekochihin and Stone2014), Riquelme, Quataert & Verscharen (Reference Riquelme, Quataert and Verscharen2015), Sironi & Narayan (Reference Sironi and Narayan2015), Melville, Schekochihin & Kunz (Reference Melville, Schekochihin and Kunz2016); in the context of the $\unicode[STIX]{x1D6FD}\sim 1$ solar wind, e.g. Hellinger & Trávníček (Reference Hellinger and Trávníček2008, Reference Hellinger and Trávníček2015)). If the production of such Larmor-scale parasites also reduces the efficacy of electron conduction (as recently suggested by Komarov et al. (Reference Komarov, Churazov, Kunz and Schekochihin2016) and Riquelme, Quataert & Verscharen (Reference Riquelme, Quataert and Verscharen2016)), one might imagine a scenario in which the MTI, reliant upon the efficient field-aligned transport of heat, grows only as fast as allowed by the parasites it triggers. However, even without the complications introduced by MTI-driven firehose and mirror instabilities, the linear system itself presents an interesting puzzle. We have a plasma in which the same temperature gradient drives growth simultaneously at large scales (drift-kinetic MTI) and small scales (eMTI), the latter with growth rates larger than the former. This implies that the large-scale MTI will likely be dependent upon the saturation of the eMTI, particularly if the latter greatly impacts the nature and efficacy of electron heat transport. Such speculation awaits well-constructed kinetic simulations.
Acknowledgements
We thank I. Abel for sharing with us his expertise on the gyrokinetic theory of stratified plasmas; A. Schekochihin, G. Hammett, and S. Cowley for useful discussions on the ITG and ETG instabilities; and H. Latter and the anonymous referees for comments that lead to an improved presentation.
Appendix A. Details concerning the derivation of the conductivity tensor (4.15)
In this appendix, we flesh out some mathematical details from the derivation of the conductivity tensor that were mercifully expunged from § 4.2. The first task of that section was to perform the $\unicode[STIX]{x1D70F}$ integral in (4.10a ), a task that is greatly aided by using the Jacobi–Anger expansion
where $a_{s}\doteq k_{\bot }v_{\bot }/\unicode[STIX]{x1D6FA}_{s}$ is the argument of the $n$ th-order Bessel function $\text{J}_{n}$ , the prime denotes differentiation with respect to that argument, and $\ast$ denotes the complex conjugate. The integration is then straightforward, the result being that the conductivity tensor defined by (4.10) simplifies to
Note the denominator of the second term, which reveals that a fluctuation with frequency $\unicode[STIX]{x1D714}$ , Doppler shifted to the frame drifting at the $-\unicode[STIX]{x1D735}\unicode[STIX]{x1D6F7}_{s}\boldsymbol{\times }\hat{\boldsymbol{e}}_{x}$ velocity, may be Landau or cyclotron resonant; the former can lead to collisionless damping, a feature of critical importance for the behaviour of convective motions.
The next step is to carry out the integral over gyrophase angle in (A 2). Before doing so, recall that the equilibrium distribution function $F_{s}$ that resides in the integrand of (A 2) is independent of gyrophase at fixed guiding centre, not at fixed position. As explained in § 4.2, we circumvent this issue by Taylor expanding $F_{s}({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ about the particle position $z$ and the gyrophase-independent particle energy $\unicode[STIX]{x1D700}_{s}$ (see (4.14)). The practical result is that the partial derivatives of $F_{s}$ in (A 2) may be replaced by the expressions
The final result is given by (4.15).
Appendix B. Definitions and integrals related to the Maxwell–Boltzmann equilibrium (4.16)
To obtain (4.24), several $v_{\Vert }$ -, $v_{\bot }$ -, and Bessel-function-weighted integrals over the Maxwell–Boltzmann equilibrium distribution function (4.16) were performed using (4.17) and (4.18). In this appendix, we provide explicit expressions for these integrals, in addition to the tensors $\unicode[STIX]{x1D650}_{n,s}$ and $\unicode[STIX]{x1D652}_{n,s}$ for which they were needed.
First, we note that the Landau-type integral (4.17) with $p=0$ is just the standard plasma dispersion function $Z(\unicode[STIX]{x1D701})$ first tabulated by Fried & Conte (Reference Fried and Conte1961). The next six orders ( $p=1\cdots 6$ ) can be written in terms of $Z(\unicode[STIX]{x1D701})$ as follows (e.g. Goodrich, Durney & Grow Reference Goodrich, Durney and Grow1972):
The Bessel-function-weighted integrals over $v_{\bot }$ are more difficult to compute, as they require multiple differentiations of the Watson (Reference Watson1966) relation (4.18) with respect to its three parameters and some creative use of the recursion relations linking different orders $n$ of the modified Bessel function $\text{I}_{n}$ . With $a_{s}\doteq k_{\bot }v_{\bot }/\unicode[STIX]{x1D6FA}_{s}$ , $\unicode[STIX]{x1D6FC}_{s}\doteq (k_{\bot }\unicode[STIX]{x1D70C}_{s})^{2}/2$ , and $\unicode[STIX]{x1D6E4}_{n}(\unicode[STIX]{x1D6FC}_{s})\doteq \text{I}_{n}(\unicode[STIX]{x1D6FC}_{s})\exp (-\unicode[STIX]{x1D6FC}_{s})$ , some particularly useful integrals of this type are
are
where $\boldsymbol{w}_{n,s}$ and $\boldsymbol{u}_{n,s}^{\ast }$ are defined via (A 5) and (A 1), respectively. Once expanded with the integrals performed, (B 9) contributes a lot of terms. Fortunately, none of them enters into any of the limits we consider in the main text, and so we feel no pressing need to write out the components of $\unicode[STIX]{x1D652}_{n,s}$ explicitly. A precocious reader may compute them using the integrals (B 1)–(B 5), several pens and lots of paper.
Appendix C. Derivation of (5.46) from the gyrokinetic theory of a thermally stratified atmosphere
In this appendix, we show how one may instead obtain (5.46) by linearizing the nonlinear electromagnetic gyrokinetic theory of a thermally stratified atmosphere. We begin by providing a derivation of that theory. While similar derivations may be found elsewhere (e.g. Frieman & Chen Reference Frieman and Chen1982; Abel et al. Reference Abel, Plunk, Wang, Barnes, Cowley, Dorland and Schekochihin2013), our restriction to slab geometry removes many of the complications introduced by those authors’ applications to axisymmetric tokamaks and makes the theory rather easier grasp, which hopefully gives this appendix some pedagogical value. Our presentation also complements the derivation of ‘astrophysical gyrokinetics’ for a homogeneous plasma found in the appendix of Howes et al. (Reference Howes, Cowley, Dorland, Hammett, Quataert and Schekochihin2006), whose traditional approach, based on an order-by-order expansion of the Vlasov equation, we adopt. However, we make no attempt to formulate ‘global’ gyrokinetic equations that simultaneously describe the long and short scales in our model atmosphere, nor do we investigate the transport-time evolution of a gyrokinetic stratified atmosphere. Our aim here is simply to provide the reader with an alternative path to the linear gyrokinetic theory analysed in § 5.3.
C.1 Gyrokinetic ordering
The derivation of the nonlinear gyrokinetic theory begins by adopting the ordering
and expanding the distribution function in powers of $\unicode[STIX]{x1D716}$ :
where the subscript indicates the order in $\unicode[STIX]{x1D716}$ . The parameters $\unicode[STIX]{x1D6FD}_{s}$ and $\unicode[STIX]{x1D70F}$ are taken to be order unity in this expansion. We neglect any temporal evolution of the equilibrium quantities due to interactions with the fluctuating fields; these occur on the transport time scale, much longer than any time scale of interest in this paper. Note that (C 1) gives $k_{\bot }\unicode[STIX]{x1D70C}_{s}\sim 1$ , and thus perpendicular spatial structure in the fluctuations is allowed on Larmor scales. The difference between the ordering (C 1) and that used to reduce the linearized Vlasov equation in § 5.3 (viz., equation (5.45)) is that, here, while the amplitudes of the fluctuations are small, they are not infinitesimally small (relative to $\unicode[STIX]{x1D716}$ ). As a result, certain nonlinear effects (such as the advection of the perturbed distribution function by the fluctuating $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ drift) are retained. To avoid any confusion between linear and nonlinear quantities, we retain the subscript ‘0’ labelling the equilibrium fields.
Following the discussion in § 3.2, we take the equilibrium distribution function $F_{s}=F_{0s}+F_{1s}+F_{2s}+\cdots \,$ to be a function of the guiding-centre position $\boldsymbol{R}_{s}$ (3.23), the particle energy ${\mathcal{E}}_{s}$ (3.18), and the gyrophase angle $\unicode[STIX]{x1D717}$ . In the presence of electromagnetic fluctuations, these variables evolve according to
and $\text{d}\unicode[STIX]{x1D717}/\text{d}t=-\unicode[STIX]{x1D6FA}_{s}+O(\unicode[STIX]{x1D716}\unicode[STIX]{x1D6FA}_{s})$ , where we have used (3.10a,b ) to write the fluctuating electromagnetic fields in terms of the scalar potential $\unicode[STIX]{x1D711}$ , the vector potential $\boldsymbol{A}$ and the combination
Equation (C 5) is often referred to as the gyrokinetic potential, for reasons that will soon become clear. We also use the results of § 3.3 and decompose the velocity into its parallel, perpendicular and drift parts (see (3.21)):
where $v_{ds}=(c/B_{0})\text{d}\unicode[STIX]{x1D6F7}_{s}/\text{d}z$ is the equilibrium drift speed.
The formal expansion of (3.1)–(3.3) that results is worked out order by order. We begin with the Vlasov equation (3.1), ordered relative to $\unicode[STIX]{x1D714}f_{0s}$ .
C.2 Gyrokinetic equation
C.2.1 Minus-first order, $O(1/\unicode[STIX]{x1D716})$
The largest term in (3.1) corresponds physically to Larmor motion of the (zeroth-order) equilibrium particle distribution about the guiding centre:
Thus, $F_{0s}$ is independent of gyrophase (gyrotropic): $F_{0s}=F_{0s}(\boldsymbol{R}_{s},{\mathcal{E}}_{s})$ .Footnote 7
C.2.2 Zeroth order, $O(1)$
Proceeding to next order and decomposing the velocity into its parallel, perpendicular and drift parts (C 6), equation (3.1) becomes
Equation (C 8) contains both fluctuating and mean quantities. To isolate the latter, we average over space to obtain
further integrating over gyrophase and using (C 7), we find that $F_{0s}$ is independent of the $x$ coordinate (as promised in § 3.2). As a result, $F_{1s}$ , like $F_{0s}$ before it, is gyrotropic: $F_{1s}=F_{1s}(\boldsymbol{R}_{s},{\mathcal{E}}_{s})$ .
Returning to (C 8) and using our newfound knowledge about the equilibrium distribution function, the first-order perturbed distribution function is related to $F_{0s}$ via
Equation (C 10) admits a homogeneous solution and a particular solution. The homogeneous solution $h_{s}$ satisfies
where we have transformed the $\unicode[STIX]{x1D717}$ derivative taken at constant position $\boldsymbol{r}$ to one taken at constant guiding centre (see (3.23)). Thus, $h_{s}$ is independent of the gyrophase angle at constant guiding centre $\boldsymbol{R}_{s}$ (but not at constant position $\boldsymbol{r}$ ):
It represents the response of charged rings to the perturbed fields, and is thus referred to as the gyrokinetic response. The particular solution of (C 10) – the so-called adiabatic, or ‘Boltzmann’, response – is given by
It arises from the evolution of $F_{0s}$ under the influence of the perturbed electromagnetic fields. Indeed,
i.e. the Boltzmann response does not change the form of the equilibrium distribution function if the latter is written as a function of sufficiently precisely conserved particle invariants.
At this point in the derivation, we assume that $F_{0s}$ is independent of the $y$ coordinate, i.e. $F_{0s}=F_{0s}(Z_{s},{\mathcal{E}}_{s})$ , where $Z_{s}$ is the position of the guiding centre projected onto the $z$ axis. At this order, we may replace $Z_{s}$ with ${\mathcal{Z}}_{s}$ (see (3.22a )), and fold their ${\sim}O(\unicode[STIX]{x1D70C}_{s}/H)$ difference into $F_{1s}$ . Combining (C 12) and (C 14), the distribution function of species $s$ may thus be written
where we have absorbed the $O(\unicode[STIX]{x1D716}^{2}F_{0s})$ terms into $\unicode[STIX]{x1D6FF}f_{2s}$ . Next we derive an evolution equation for the gyrokinetic response $h_{s}$ .
C.2.3 First order, $O(\unicode[STIX]{x1D716})$
Collecting $O(\unicode[STIX]{x1D716})$ terms in (3.1) and decomposing $\unicode[STIX]{x1D6FF}f_{1s}$ into its adiabatic (C 13) and non-adiabatic (C 12) parts (cf. Antonsen & Lane Reference Antonsen and Lane1980; Catto et al. Reference Catto, Tang and Baldwin1981), we find that the gyrokinetic response evolves according to
where
is the velocity of the guiding centre required at this order (cf. C 3). Upon performing a ring average of (C 16) over $\unicode[STIX]{x1D717}$ at fixed $\boldsymbol{R}_{s}$ , defined, for any function $a(t,\boldsymbol{r},\boldsymbol{v})$ , as
we find that the entire right-hand side of (C 16) vanishes. This follows from the periodicity of $F_{2s}$ and $\unicode[STIX]{x1D6FF}f_{2s}$ in $\unicode[STIX]{x1D717}$ , the gyrotropy of $F_{0s}$ and $F_{1s}$ at fixed guiding centre, and from the fact that, for any arbitrary function $a(\boldsymbol{r})$ , the ring average $\langle \boldsymbol{v}_{\bot }\boldsymbol{\cdot }\unicode[STIX]{x1D735}_{\bot }a\rangle _{\boldsymbol{R}_{s}}=0$ (see equation (A21) of Howes et al. (Reference Howes, Cowley, Dorland, Hammett, Quataert and Schekochihin2006)). Thus, the ring-averaged (C 16) is
Using the decomposition $\unicode[STIX]{x1D6FF}\boldsymbol{B}=\unicode[STIX]{x1D735}A_{\Vert }\boldsymbol{\times }\hat{\boldsymbol{e}}_{x}+\unicode[STIX]{x1D6FF}B_{\Vert }\hat{\boldsymbol{e}}_{x}$ in (C 3) and retaining only first-order contributions, the ring-averaged guiding-centre velocity (C 17) is
The second term in (C 20a ) captures the (slow) $-\unicode[STIX]{x1D735}\unicode[STIX]{x1D6F7}_{s}\boldsymbol{\times }\boldsymbol{B}$ drift across the magnetic-field lines. The third term is due to the ring-averaged fluctuating $\boldsymbol{E}\boldsymbol{\times }\boldsymbol{B}$ drift. Finally, the fourth term represents the $\unicode[STIX]{x1D735}B$ drift; if we expand the ring average to lowest order in $\boldsymbol{v}\boldsymbol{\times }\hat{\boldsymbol{e}}_{x}/\unicode[STIX]{x1D6FA}_{s}$ , we obtain the familiar drift velocity
where $\unicode[STIX]{x1D735}B=\unicode[STIX]{x1D735}\unicode[STIX]{x1D6FF}B_{\Vert }$ is evaluated at the centre of the ring, $\boldsymbol{r}=\boldsymbol{R}_{s}$ .
Substituting (C 20) into (C 19), we finally obtain the gyrokinetic equation
where the Poisson bracket is defined by $\{U,V\}=\hat{\boldsymbol{e}}_{x}\boldsymbol{\cdot }\{(\unicode[STIX]{x2202}U/\unicode[STIX]{x2202}\boldsymbol{R}_{s})\boldsymbol{\times }(\unicode[STIX]{x2202}V/\unicode[STIX]{x2202}\boldsymbol{R}_{s})\}$ . In order to solve (C 23), we must relate the electromagnetic potentials to the gyrokinetic response. This is accomplished via the field equations.
C.3 Field equations
To obtain equations governing the electromagnetic potentials, we substitute (C 15) into the leading-order expansions of the quasi-neutrality constraint (3.2) and Ampère’s law (3.3). To first order, $O(\unicode[STIX]{x1D716})$ , (3.2) becomes
Because the electrostatic potential $\unicode[STIX]{x1D711}$ is a function of the spatial variable $\boldsymbol{r}$ , the velocity integral of $h_{s}$ must be performed at constant location $\boldsymbol{r}$ of the charges rather than at constant guiding centre $\boldsymbol{R}_{s}$ . This introduces a gyro-averaging operation dual to the ring average defined in (C 18):
Equation (C 24) may then be written as
Likewise, the parallel and perpendicular components of Ampère’s law become, respectively,
Together with the gyrokinetic equation (C 23), the field equations (C 26)–(C 28) provide a closed system that describes the nonlinear evolution of a gyrokinetic stratified atmosphere.
C.4 Fourier space
Our next task is to linearize (C 23) and substitute the resulting expression for the gyrokinetic response into the field equations. Before doing so, let us calculate the Fourier representation of the ring-averaging (C 18) and gyro-averaging (C 25) operations. Decomposing the gyrokinetic potential into plane waves, $\unicode[STIX]{x1D712}(t,\boldsymbol{r},\boldsymbol{v})=\sum _{\boldsymbol{k}}\unicode[STIX]{x1D712}_{\boldsymbol{k}}\exp (\text{i}\boldsymbol{k}\boldsymbol{\cdot }\boldsymbol{r})$ , the ring average of its Fourier coefficient is
where $a_{s}\doteq k_{\bot }v_{\bot }/\unicode[STIX]{x1D6FA}_{s}$ is the argument of the zeroth- ( $\text{J}_{0}$ ) and first-order ( $\text{J}_{1}$ ) Bessel functions. Similarly, we can Fourier decompose $h_{s}(t,\boldsymbol{R}_{s},v_{\Vert },v_{\bot })=\sum _{\boldsymbol{k}}h_{s\boldsymbol{k}}(t,v_{\Vert },v_{\bot })\exp (\text{i}\boldsymbol{k}\boldsymbol{\cdot }\boldsymbol{R}_{s})$ and perform its gyroaverages at constant $\boldsymbol{r}$ , viz.
Using these formulae, the field equations (C 26)–(C 28) expressed in Fourier space are
We now derive the linear gyrokinetic theory.
C.5 Linear theory
Returning to the gyrokinetic equation (C 23), we drop the nonlinear term $\propto \{\langle \unicode[STIX]{x1D712}\rangle _{\boldsymbol{R}_{s}},h_{s}\}$ and adopt plane-wave solutions to find
where we have safely omitted the ‘0’ subscript adorning the equilibrium quantities. (This is essentially equation (24) of Antonsen & Lane Reference Antonsen and Lane1980.) At this order in the gyrokinetic expansion, we may replace the phase-space variables $({\mathcal{Z}}_{s},{\mathcal{E}}_{s})$ by their leading-order counterparts $(z,\unicode[STIX]{x1D700}_{s})$ (see (4.14) and related discussion in § 4.2) and use (4.20), valid for a Maxwell–Boltzmann equilibrium distribution, to evaluate the phase-space derivatives of $F_{s}$ in (C 35). Substituting the resulting expression for $h_{s\boldsymbol{k}}$ into the field equations (C 32)–(C 34) and recalling the definitions of the temperature-gradient and diamagnetic drift frequencies,
we have
with $f_{M,s}=f_{M,s}(z,v_{\Vert },v_{\bot })$ given by (4.22). The final step in the derivation of the linear gyrokinetic theory is to replace the ring-averaged gyrokinetic potential in (C 37)–(C 39) by (C 29) and perform the resulting integrals using (4.17), (B 1a–d ), and (B 2a,b ), (B 3a,b ), (B 4a ). This leads identically to (5.46), as promised.