1 Introduction
An interesting and fundamental question in plasma physics is how thermal equilibrium is determined in an essentially collisionless plasma. In such a system, there is no reason a priori to assume that the relative temperatures of the component species will be equal: on the contrary, there are numerous collisionless heating mechanisms that have been identified – shocks, magnetic reconnection, cyclotron resonance (Cranmer, Field & Kohl Reference Cranmer, Field and Kohl1999) and various forms of turbulent heating (Quataert & Gruzinov Reference Quataert and Gruzinov1999; Howes et al. Reference Howes, Cowley, Dorland, Hammett, Quataert and Schekochihin2008; Chandran et al. Reference Chandran, Li, Rogers, Quataert and Germaschewski2010; Barnes, Parra & Dorland Reference Barnes, Parra and Dorland2012), to name a few – and each of these mechanisms may drive temperature separation instead of equilibration. Indeed, evidence from observations of collisionless space and astrophysical plasmas, e.g. the solar wind and accretion flows onto compact objects, suggests that more massive charged particles may be preferentially heated; cf. Schmidt et al. (Reference Schmidt, Rosenbauer, Shelly and Geiss1980), Collier et al. (Reference Collier, Hamilton, Gloeckler, Bochsler and Sheldon1996), Kohl et al. (Reference Kohl, Noci, Antonucci, Tondello, Huber, Gardner, Nicolosi, Strachan, Fineschi and Raymond1997), Quataert & Gruzinov (Reference Quataert and Gruzinov1999).
A thorough treatment of this problem would require exploration of a parameter space spanning a wide range of plasma
$\unicode[STIX]{x1D6FD}$
(ratio of plasma to magnetic pressure), electron–ion temperature ratio, energy injection mechanisms, etc. while including all potentially relevant physics present from the macroscopic to microscopic space–time scales. As such a study is infeasible, we choose to focus on a single heating mechanism – turbulent Joule heating – in a simplified system: a two-component plasma immersed in a straight, homogenous magnetic field with a cross-field density gradient. This system supports electron drift waves and is known to be susceptible to the so-called ‘universal’ instability (Galeev, Oraevsky & Sagdeev Reference Galeev, Oraevsky and Sagdeev1963; Krall & Rosenbluth Reference Krall and Rosenbluth1965; Helander & Plunk Reference Helander and Plunk2015; Landreman, Antonsen & Dorland Reference Landreman, Antonsen and Dorland2015). The universal instability serves as an energy injection mechanism for turbulence at the ion Larmor scale that gives rise to plasma fluctuations that drive cross-field transport and Joule heating (or cooling) of the component plasma species. For sufficiently infrequent collisions, this turbulent heating and transport determines the thermal equilibrium of the plasma.
To guide our study, we note that the universal instability is derived analytically in the limit of disparate ion and electron thermal speeds. As noted in Helander (Reference Helander2014), the universal instability is eliminated for an equal temperature pair plasma, or, more generally, when the ion and electron thermal speeds are equal. In the following sections we analytically and numerically calculate the dependence of the universal instability and resultant turbulent heating on the electron–ion temperature ratio for two cases: pair plasmas and conventional mass-ratio plasmas. The former case is relevant for proposed experiments in which a small population of positrons is to be added to a pure electron plasma of a different temperature (Pedersen & Boozer Reference Pedersen and Boozer2002); the latter may provide insight into how astrophysical plasmas can achieve disparate electron and ion temperatures. We introduce the gyrokinetic model used for the analysis in § 2. In § 3 we calculate linear growth rates and use them to obtain quasilinear heating and cross-field flux estimates. We provide linear, quasilinear and nonlinear simulation results in § 4 to support the analytical calculations before concluding.
2 Model system
We consider a collisionless plasma immersed in a straight, homogeneous magnetic field
$\boldsymbol{B}=B\hat{\boldsymbol{z}}$
and with a fixed density gradient perpendicular to the field in the
$x$
-direction. The fixed density gradient could be, e.g. the result of gravitational equilibrium in an astrophysical plasma or of an external particle source in a laboratory plasma. We restrict our attention to electrostatic fluctuations whose frequency is small compared with the Larmor frequency and adopt the gyrokinetic ordering:
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D716}\doteq {\displaystyle \frac{\unicode[STIX]{x1D714}}{\unicode[STIX]{x1D6FA}_{s}}}\sim {\displaystyle \frac{k_{\Vert }}{k_{\bot }}}\sim {\displaystyle \frac{\unicode[STIX]{x1D70C}_{s}}{L_{n}}}\sim {\displaystyle \frac{\unicode[STIX]{x1D6FF}f_{s}}{f_{s}}}\sim {\displaystyle \frac{q_{s}\unicode[STIX]{x1D711}}{T_{s}}}\ll 1, & \displaystyle\end{eqnarray}$$
where
$\unicode[STIX]{x1D711}$
is the electrostatic potential fluctuation,
$f_{s}$
is the distribution function for species
$s$
,
$\unicode[STIX]{x1D6FF}f$
is the fluctuating component of
$f$
,
$\unicode[STIX]{x1D714}$
is the characteristic frequency of the fluctuations,
$k_{\Vert }$
and
$k_{\bot }$
are the associated wavenumbers along and across the mean field,
$\unicode[STIX]{x1D6FA}_{s}$
is the Larmor frequency for species
$s$
,
$\unicode[STIX]{x1D70C}_{s}$
its thermal Larmor radius,
$T_{s}$
its temperature,
$q_{s}$
its charge and
$L_{n}$
is the mean density gradient scale length.
Applying these orderings to the Fokker–Planck equation and averaging over the rapid gyration of particles about the mean magnetic field results in the gyrokinetic equation,
$$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\unicode[STIX]{x2202}h_{s}}{\unicode[STIX]{x2202}t}}+v_{\Vert }{\displaystyle \frac{\unicode[STIX]{x2202}h_{s}}{\unicode[STIX]{x2202}z}}+{\displaystyle \frac{c}{B}}\{\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}},h_{s}\}=-q_{s}{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}}}{\unicode[STIX]{x2202}t}}{\displaystyle \frac{\unicode[STIX]{x2202}F_{0,s}}{\unicode[STIX]{x2202}E_{s}}}+{\displaystyle \frac{c}{B}}{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}}}{\unicode[STIX]{x2202}y}}{\displaystyle \frac{\unicode[STIX]{x2202}F_{0,s}}{\unicode[STIX]{x2202}x}}+C[h_{s}], & \displaystyle\end{eqnarray}$$
where
$h_{s}=\unicode[STIX]{x1D6FF}f_{s}-q_{s}\unicode[STIX]{x1D711}(\unicode[STIX]{x2202}F_{0,s}/\unicode[STIX]{x2202}E_{s})$
is the non-Boltzmann piece of
$\unicode[STIX]{x1D6FF}f_{s}$
,
$E_{s}=m_{s}v^{2}/2$
is the kinetic energy of species
$s$
,
$m_{s}$
is species mass,
$t$
is time,
$v_{\Vert }$
is the parallel component of the particle velocity,
$c$
is the speed of light,
$F_{0,s}$
is the mean component of
$f_{s}$
,
$\{\}$
indicates a Poisson bracket,
$\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}}$
is the average of
$\unicode[STIX]{x1D711}$
over Larmor angle at fixed guiding centre position
$\boldsymbol{R}_{s}=\boldsymbol{r}+\boldsymbol{v}\times \hat{\boldsymbol{z}}/\unicode[STIX]{x1D6FA}_{s}$
,
$\boldsymbol{r}$
is particle position and
$C[h_{s}]$
represents the effect of Coulomb collisions on species
$s$
.
The gyrokinetic system is closed by coupling to Poisson’s equation:
$$\begin{eqnarray}\displaystyle & \displaystyle 4\unicode[STIX]{x03C0}\mathop{\sum }_{s}q_{s}\int \text{d}^{3}v\unicode[STIX]{x1D6FF}f_{s}=-\unicode[STIX]{x1D6FB}_{\bot }^{2}\unicode[STIX]{x1D711}, & \displaystyle\end{eqnarray}$$
with
$\unicode[STIX]{x1D6FB}_{\bot }^{2}=\unicode[STIX]{x2202}_{x}^{2}+\unicode[STIX]{x2202}_{y}^{2}$
. If the Debye length
$\unicode[STIX]{x1D706}_{s}\doteq (T_{s}/4\unicode[STIX]{x03C0}n_{s}q_{s}^{2})^{1/2}$
is much smaller than the electron Larmor radius, the right-hand side of equation (2.3) can be neglected. In this limit Poisson’s equation reduces to the quasineutrality constraint that the total charge density of the plasma is zero. The ratio of Debye length to Larmor radius can be written as
$\unicode[STIX]{x1D706}_{s}^{2}/\unicode[STIX]{x1D70C}_{s}^{2}=(1/2\unicode[STIX]{x1D6FD}_{s})(v_{\text{th},s}/c)^{2}$
, with
$\unicode[STIX]{x1D6FD}_{s}=8\unicode[STIX]{x03C0}p_{s}/B^{2}$
. Thus Debye length effects are only significant in Poisson’s equation if the plasma is relativistic and/or if
$\unicode[STIX]{x1D6FD}_{s}$
is sufficiently small (
${\lesssim}v_{\text{th},s}^{2}/c^{2}$
). The latter condition is consistent with the non-relativistic, electrostatic approximation we employ.
Note that even when
$\unicode[STIX]{x1D706}_{s}^{2}/\unicode[STIX]{x1D70C}_{s}^{2}\sim 1$
, the displacement current appearing in Ampere’s law is negligible in the gyrokinetic ordering compared to the plasma current. Only when taking the divergence of Ampere’s law, i.e. when considering charge conservation, must the displacement current be retained for
$\unicode[STIX]{x1D706}_{s}^{2}/\unicode[STIX]{x1D70C}_{s}^{2}\sim 1$
.
2.1 Turbulent heating and heat transport
By taking fluid moments of the Fokker–Planck equation and closing the system with the gyrokinetic ordering, one obtains an equation for the slow evolution of the mean temperature
$T_{s}$
(Sugama et al.
Reference Sugama, Okamoto, Horton and Wakatani1996; Hinton & Waltz Reference Hinton and Waltz2006; Howes et al.
Reference Howes, Cowley, Dorland, Hammett, Quataert and Schekochihin2006; Barnes et al.
Reference Barnes, Abel, Dorland, Goerler, Hammett and Jenko2010a
; Abel et al.
Reference Abel, Plunk, Wang, Barnes, Cowley, Dorland and Schekochihin2013):
$$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{3}{2}}n_{s}{\displaystyle \frac{\text{d}T_{s}}{\text{d}t}}+{\displaystyle \frac{\unicode[STIX]{x2202}Q_{s}}{\unicode[STIX]{x2202}x}}=H_{s}, & \displaystyle\end{eqnarray}$$
where the cross-field turbulent heat flux
$Q_{s}$
and turbulent heating
$H_{s}$
are given by
$$\begin{eqnarray}\displaystyle & \displaystyle Q_{s}={\displaystyle \frac{1}{V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}v\left({\displaystyle \frac{m_{s}v^{2}}{2}}-{\displaystyle \frac{3}{2}}T_{s}\right)\overline{(h_{s}v_{E})} & \displaystyle\end{eqnarray}$$
and
$$\begin{eqnarray}\displaystyle & \displaystyle H_{s}={\displaystyle \frac{1}{V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}vq_{s}\overline{\left(h_{s}{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}}}{\unicode[STIX]{x2202}t}}\right)}. & \displaystyle\end{eqnarray}$$
Here
$V$
is the volume of the region over which the spatial integration is performed, the subscript
$\boldsymbol{R}$
on the velocity integration indicates that it is carried out at fixed guiding centre,
$v_{E}=-(c/B)\unicode[STIX]{x2202}\unicode[STIX]{x1D711}/\unicode[STIX]{x2202}y$
is the
$x$
-component of the
$E\times B$
drift velocity and the overline indicates an average over time scales long compared to the fluctuation time
$1/\unicode[STIX]{x1D714}$
but short compared to the equilibrium time scale. Note that we have neglected collisional temperature equilibration as we are considering systems with collisional mean free path much longer than any other scales of interest.
An alternative expression for
$H_{s}$
is obtained by integrating equation (2.6) by parts in time (Candy Reference Candy2013):
$$\begin{eqnarray}\displaystyle & \displaystyle H_{s}={\displaystyle \frac{1}{2V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}vq_{s}\overline{\left(h_{s}{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}}}{\unicode[STIX]{x2202}t}}-{\displaystyle \frac{\unicode[STIX]{x2202}h_{s}}{\unicode[STIX]{x2202}t}}\langle \unicode[STIX]{x1D711}\rangle _{\boldsymbol{R}_{s}}\right)}. & \displaystyle\end{eqnarray}$$
Substitution of Poisson’s equation (2.3) in (2.7) immediately indicates that the net (species-summed) turbulent heating is zero in the absence of an external energy injection mechanism; i.e. for a two-component plasma,
$\overline{H_{i}}=-\overline{H_{e}}$
. If trace minority ions are present, one can obtain mass- and charge-dependent scalings for their turbulent heating rates relative to the main ions (Barnes et al.
Reference Barnes, Parra and Dorland2012). Expressing equations (2.5) and (2.7) in terms of Fourier modes with
$h_{s}=\sum _{\boldsymbol{k}}\tilde{h}_{\boldsymbol{k},s}(v_{\Vert },v_{\bot },t)\exp (\text{i}\boldsymbol{k}\boldsymbol{\cdot }\boldsymbol{R}_{s})$
gives
$$\begin{eqnarray}\displaystyle & \displaystyle Q_{s}=\text{i}{\displaystyle \frac{c}{B}}\mathop{\sum }_{\boldsymbol{k}}k_{y}\int \text{d}^{3}v\text{J}_{0}(\unicode[STIX]{x1D6FC}_{k_{\bot },s}){\displaystyle \frac{m_{s}v^{2}}{2}}\overline{(\tilde{h}_{\boldsymbol{k},s}\tilde{\unicode[STIX]{x1D711}}_{\boldsymbol{k}}^{\ast })}\doteq \mathop{\sum }_{\boldsymbol{ k}}\tilde{Q}_{\boldsymbol{k},s} & \displaystyle\end{eqnarray}$$
and
$$\begin{eqnarray}\displaystyle & \displaystyle H_{s}={\displaystyle \frac{q_{s}}{2}}\mathop{\sum }_{\boldsymbol{k}}\int \text{d}^{3}v\text{J}_{0}(\unicode[STIX]{x1D6FC}_{k_{\bot },s})\overline{\left({\displaystyle \frac{\unicode[STIX]{x2202}\tilde{\unicode[STIX]{x1D711}}_{\boldsymbol{k}}^{\ast }}{\unicode[STIX]{x2202}t}}\tilde{h}_{\boldsymbol{k},s}-\tilde{\unicode[STIX]{x1D711}}_{\boldsymbol{k}}^{\ast }{\displaystyle \frac{\unicode[STIX]{x2202}\tilde{h}_{\boldsymbol{k},s}}{\unicode[STIX]{x2202}t}}\right)}\doteq \mathop{\sum }_{\boldsymbol{k}}\tilde{H}_{\boldsymbol{k},s}, & \displaystyle\end{eqnarray}$$
where
$\unicode[STIX]{x1D6FC}_{k_{\bot },s}\doteq k_{\bot }v_{\bot }/\unicode[STIX]{x1D6FA}_{s}$
and
$\text{J}_{0}$
is a Bessel function of the first kind.
2.2 Connection to particle transport
The effect of cross-field particle transport is encapsulated in the continuity equation:
$$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\unicode[STIX]{x2202}n_{s}}{\unicode[STIX]{x2202}t}}+{\displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6E4}_{s}}{\unicode[STIX]{x2202}x}}=0, & \displaystyle\end{eqnarray}$$
where
$n_{s}$
is particle density and
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E4}_{s}={\displaystyle \frac{1}{V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}v\overline{(h_{s}v_{E})} & \displaystyle\end{eqnarray}$$
is the cross-field particle flux. Expanding
$\unicode[STIX]{x1D711}$
and
$h$
in terms of Fourier modes gives
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E4}_{s}=\text{i}{\displaystyle \frac{c}{B}}\mathop{\sum }_{\boldsymbol{k}}k_{y}\int \text{d}^{3}v\text{J}_{0}(\unicode[STIX]{x1D6FC}_{k_{\bot },s})\overline{(\tilde{h}_{\boldsymbol{k},s}\tilde{\unicode[STIX]{x1D711}}_{\boldsymbol{k}}^{\ast })}\doteq \mathop{\sum }_{\boldsymbol{k}}\tilde{\unicode[STIX]{x1D6E4}}_{\boldsymbol{k},s}. & \displaystyle\end{eqnarray}$$
The particle flux
$\unicode[STIX]{x1D6E4}$
can be related to the turbulent heating of (2.6) by multiplying the gyrokinetic equation (2.2) by
$\unicode[STIX]{x1D711}$
and averaging over the phase space:
$$\begin{eqnarray}\displaystyle & \displaystyle H_{s}=-{\displaystyle \frac{1}{V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}v\left({\displaystyle \frac{T_{s}h_{s}}{F_{0,s}}}C[h_{s}]\right)+T_{s}{\displaystyle \frac{\unicode[STIX]{x2202}\ln n_{s}}{\unicode[STIX]{x2202}x}}\unicode[STIX]{x1D6E4}_{s}. & \displaystyle\end{eqnarray}$$
Note that the effect of collisions is retained in (2.13), but collisional temperature equilibration is neglected in the energy equation (2.4). This is because the turbulent fluctuations undergo filamentation in phase space (Schekochihin et al.
Reference Schekochihin, Cowley, Dorland, Hammett, Howes, Plunk, Quataert and Tatsuno2008, Reference Schekochihin, Cowley, Dorland, Hammett, Howes, Quataert and Tatsuno2009; Tatsuno et al.
Reference Tatsuno, Dorland, Schekochihin, Plunk, Barnes, Cowley and Howes2009; Barnes, Dorland & Tatsuno Reference Barnes, Dorland and Tatsuno2010b
; Hatch et al.
Reference Hatch, Jenko, Bratanov and Navarro2014) and thus tiny collisional deflections in velocity have immediate impact on the turbulent fluctuations – even when the amount of energy that is exchanged in the course of the deflection is negligible. At the small phase space scales where collisional dissipation occurs, the dissipation is effectively diffusive (Abel et al.
Reference Abel, Barnes, Cowley, Dorland, Hammett and Schekochihin2008; Barnes et al.
Reference Barnes, Dorland and Tatsuno2010b
, Reference Barnes, Parra and Dorland2012). This gives rise to heating that is positive definite for each species. Furthermore, Poisson’s equation (2.3) applied to our two-component plasma dictates that
$\unicode[STIX]{x1D6E4}_{i}=\unicode[STIX]{x1D6E4}_{e}\doteq \unicode[STIX]{x1D6E4}$
and
$n_{i}=n_{e}\doteq n$
. Substituting these expressions in (2.13) and enforcing
$H_{i}=-H_{e}$
results in the constraint that
$\unicode[STIX]{x1D6E4}(\unicode[STIX]{x2202}\ln n/\unicode[STIX]{x2202}x)\leqslant 0$
and thus
$|H_{s}|\leqslant T_{s}|\unicode[STIX]{x1D6E4}(\unicode[STIX]{x2202}\ln n/\unicode[STIX]{x2202}x)|$
. Comparing this inequality with (2.4) and (2.10), one finds that the time scale associated with turbulent particle transport is always at least as fast as that associated with turbulent heating. One would therefore expect that when turbulent heating drives temperatures apart (as we find below), the effect would be bounded by the rate of transport down the density gradient. For systems supported by an external plasma source, e.g. a magnetically confined plasma, this may be set by the rate at which plasma is introduced.
3 Linear analysis
If we consider small amplitude perturbations, we may neglect the quadratic nonlinearity and carry out a linear analysis of the gyrokinetic equation. Upon assuming solutions of the form
$h_{s}=\sum _{\boldsymbol{k}}{\hat{h}}_{\boldsymbol{k},s}(v_{\Vert },v_{\bot })\exp (\text{i}\boldsymbol{k}\boldsymbol{\cdot }\boldsymbol{R}_{s}-\text{i}\unicode[STIX]{x1D714}t)$
, one obtains
$$\begin{eqnarray}\displaystyle & \displaystyle {\hat{h}}_{\boldsymbol{k},s}={\displaystyle \frac{q_{s}\hat{\unicode[STIX]{x1D711}}_{\boldsymbol{k}}}{T_{s}}}\text{J}_{0}(\unicode[STIX]{x1D6FC}_{k_{\bot },s})\left({\displaystyle \frac{\unicode[STIX]{x1D714}+\unicode[STIX]{x1D714}_{\ast ,s}}{\unicode[STIX]{x1D714}-k_{\Vert }v_{\Vert }}}\right)F_{M,s}, & \displaystyle\end{eqnarray}$$
where
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D714}_{\ast ,s}={\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{s}v_{\text{th},s}}{2L_{n,s}}}{\displaystyle \frac{q_{s}}{|q_{s}|}}, & \displaystyle\end{eqnarray}$$
$v_{\text{th},s}^{2}=2T_{s}/m_{s}$
,
$m_{s}$
is species mass,
$\unicode[STIX]{x1D70C}_{s}=v_{\text{th},s}/|\unicode[STIX]{x1D6FA}_{s}|$
and
$1/L_{n,s}=-\unicode[STIX]{x2202}\ln n_{s}/\unicode[STIX]{x2202}x$
.
From Poisson’s equation (2.3),
$$\begin{eqnarray}\displaystyle & \displaystyle \int \text{d}^{3}v(h_{i}-h_{e})={\displaystyle \frac{en}{T_{e}}}(1+\unicode[STIX]{x1D70F}-\unicode[STIX]{x1D706}_{e}^{2}\unicode[STIX]{x1D6FB}_{\bot }^{2})\unicode[STIX]{x1D711}, & \displaystyle\end{eqnarray}$$
where
$\unicode[STIX]{x1D70F}\doteq T_{e}/T_{i}$
, and we have restricted our attention to a single ion species with proton charge
$e$
. Substitution of (3.1) into (3.3) results in the dispersion relation
$$\begin{eqnarray}\displaystyle \displaystyle \unicode[STIX]{x1D716}(\unicode[STIX]{x1D714},\boldsymbol{k}) & = & \displaystyle 1+\unicode[STIX]{x1D70F}+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2}+\left(\unicode[STIX]{x1D701}_{e}-{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{e}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{e})Z(\unicode[STIX]{x1D701}_{e})\nonumber\\ \displaystyle & & \displaystyle +\,\unicode[STIX]{x1D70F}\left(\unicode[STIX]{x1D701}_{i}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})Z(\unicode[STIX]{x1D701}_{i})=0,\end{eqnarray}$$
where
$\unicode[STIX]{x1D701}_{s}\doteq \unicode[STIX]{x1D714}/k_{\Vert }v_{\text{th},s}$
,
$$\begin{eqnarray}\displaystyle & \displaystyle Z(x)\doteq \text{i}\sqrt{\unicode[STIX]{x03C0}}\text{e}^{-x^{2}}~\text{erfc}(-\text{i}x) & \displaystyle\end{eqnarray}$$
is the plasma dispersion function,
$\text{erfc}$
is the complementary error function and
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6F6}(x)\doteq \exp \left(-{\displaystyle \frac{x^{2}}{2}}\right)\text{I}_{0}\left({\displaystyle \frac{x^{2}}{2}}\right), & \displaystyle\end{eqnarray}$$
with
$\text{I}_{0}$
a modified Bessel of the first kind. Note that we have used quasineutrality to set
$L_{n,i}=L_{n,e}\doteq L_{n}$
.
3.1 Quasilinear energy exchange
Using equation (3.1) for
${\hat{h}}_{\boldsymbol{k},s}$
in (2.9) we obtain a quasilinear approximation for the energy exchange:
$$\begin{eqnarray}\displaystyle \displaystyle {\hat{H}}_{\boldsymbol{k},s} & \doteq & \displaystyle {\displaystyle \frac{\tilde{H}_{\boldsymbol{k},s}}{nT_{i}}}{\displaystyle \frac{|L_{n}|}{v_{\text{th},i}}}{\displaystyle \frac{T_{i}^{2}}{e^{2}|\tilde{\unicode[STIX]{x1D711}}_{\boldsymbol{k}}|^{2}}}\nonumber\\ \displaystyle & = & \displaystyle \unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{s})k_{\Vert }|L_{n}|{\displaystyle \frac{T_{i}}{T_{s}}}\text{Re}(\unicode[STIX]{x1D701}_{i})\text{Im}\left(\left(\unicode[STIX]{x1D701}_{s}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{s}}{2k_{\Vert }L_{n}}}{\displaystyle \frac{q_{s}}{|q_{s}|}}\right)Z(\unicode[STIX]{x1D701}_{s})\right).\end{eqnarray}$$
It is straightforward to verify that summing this expression over species and using the dispersion relation of (3.4) leads to zero net heating.
This quasilinear estimate provides information about linear phase relationships that indicate the sign of the heating contribution as a function of wavelength. It predicts neither the spectrum nor the saturated heating amplitude in steady state. As we are primarily interested in the sign of the heating, it is enough to make some assumptions about the turbulent heating spectrum; in particular, we assume that the steady-state heating has the same sign as the quasilinear estimate for the mode(s) with the largest linear growth rate. In the following subsections we obtain numerical and approximate analytical solutions for the mode frequencies and associated quasilinear heating in various limits.
3.2 Comparable thermal speeds
The dispersion relation (3.4) has been analysed thoroughly for both pair plasmas and conventional mass-ratio plasmas with equal temperatures; cf. Helander & Connor (Reference Helander and Connor2016), Mischenko et al. (Reference Mischenko, Zocco, Helander and Könies2018). Here, we begin by showing that there is no instability, and thus no turbulent heating, in the general case of
$v_{\text{th},e}=v_{\text{th},i}$
. To find the condition for marginal stability, we seek solutions for which
$\unicode[STIX]{x1D6FE}\doteq \text{Im}(\unicode[STIX]{x1D714})=0$
. In this case, the plasma dispersion function simplifies to
$Z(x)=\sqrt{\unicode[STIX]{x03C0}}\exp (-x^{2})(\text{i}-\text{erfi}(x))$
, with
$\text{erfi}$
the imaginary error function and
$x=\text{Re}(\unicode[STIX]{x1D701})$
. The constraint
$\text{Im}(\unicode[STIX]{x1D716})=0$
then gives
$$\begin{eqnarray}\displaystyle & \displaystyle \left(x_{e}-{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{e}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{e})\exp (-x_{e}^{2})=-\unicode[STIX]{x1D70F}\left(x_{i}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})\exp (-x_{i}^{2}). & \displaystyle\end{eqnarray}$$
Substituting this expression into the constraint
$\text{Re}(\unicode[STIX]{x1D716})=0$
gives
$$\begin{eqnarray}\displaystyle & \displaystyle 0=1+\unicode[STIX]{x1D70F}+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2}+\unicode[STIX]{x1D70F}\left(x_{i}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})\unicode[STIX]{x03C0}^{1/2}\exp (-x_{i}^{2})(\text{erfi}(x_{e})-\text{erfi}(x_{i})).\quad & \displaystyle\end{eqnarray}$$
When
$v_{\text{th},e}=v_{\text{th},i}$
, then
$x_{e}=x_{i}$
, and equation (3.9) has no solution. Such a plasma is therefore either always stable or always unstable, independent of wavenumber and density gradient. As there is no instability for zero density gradient, the plasma must therefore be always stable for
$v_{\text{th},e}=v_{\text{th},i}$
.
Next, we seek marginally stable solutions when
$v_{\text{th},e}=v_{\text{th},i}(1+\unicode[STIX]{x1D6FF})$
, with
$|\unicode[STIX]{x1D6FF}|\ll 1$
. In this limit, the constraint equation (3.9) becomes
$$\begin{eqnarray}\displaystyle & \displaystyle 0\approx 1+{\displaystyle \frac{k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2}}{2}}-x_{i}\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})\left(x_{i}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6FF}, & \displaystyle\end{eqnarray}$$
with solutions given by
$$\begin{eqnarray}\displaystyle & \displaystyle x_{i}=-{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{4k_{\Vert }L_{n}}}\pm \sqrt{\left({\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{4k_{\Vert }L_{n}}}\right)^{2}+{\displaystyle \frac{1+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2}/2}{\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})\unicode[STIX]{x1D6FF}}}}. & \displaystyle\end{eqnarray}$$
In order for the solutions to be real (as we assumed when we considered marginal stability), the term inside the square root must be positive definite. This constraint is satisfied when
$\unicode[STIX]{x1D6FF}>0$
or when
$(k_{y}\unicode[STIX]{x1D70C}_{i})^{2}\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})|\unicode[STIX]{x1D6FF}|>16(k_{\Vert }L_{n})^{2}(1+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2}/2)$
. The latter constraint can always be satisfied for sufficiently long parallel wavelengths. Consequently, an unbounded system can be unstable for all finite values of
$\unicode[STIX]{x1D6FF}$
; i.e. for all plasmas with
$v_{\text{th},i}\neq v_{\text{th},e}$
.
We now proceed to obtain the sign of the turbulent heating driven by instabilities with
$\unicode[STIX]{x1D6FF}<0$
and
$\unicode[STIX]{x1D6FF}>0$
, respectively. Since
$|\unicode[STIX]{x1D6FF}|\ll 1$
, we can approximate
$\unicode[STIX]{x1D701}_{s}\approx x_{s}$
in the quasilinear heating expression (3.7) to obtain:
$$\begin{eqnarray}\displaystyle & \displaystyle {\hat{H}}_{\boldsymbol{k},i}=k_{\Vert }|L_{n}|x_{i}\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})\left(x_{i}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{2k_{\Vert }L_{n}}}\right){\displaystyle \frac{\sqrt{\unicode[STIX]{x03C0}}}{2}}\exp (-x_{i}^{2}). & \displaystyle\end{eqnarray}$$
When
$\unicode[STIX]{x1D6FF}<0$
,
$x_{i}$
given by (3.11) satisfies
$-k_{y}\unicode[STIX]{x1D70C}_{i}/2k_{\Vert }L_{n}<x_{i}<0$
. Substituting this range of
$x_{i}$
values in (3.12) results in the constraint
${\hat{H}}_{\boldsymbol{k},i}<0$
. When
$\unicode[STIX]{x1D6FF}>0$
,
$x_{i}$
given by (3.7) satisfies
$x_{i}>0$
or
$x_{i}<-k_{y}\unicode[STIX]{x1D70C}_{i}/2k_{\Vert }L_{n}$
. Substituting this range of
$x_{i}$
values in (3.12) results in the constraint
${\hat{H}}_{\boldsymbol{k},i}>0$
. Combining these two constraints gives the general expression
$\text{sgn}({\hat{H}}_{\boldsymbol{k},i})=\text{sgn}(\unicode[STIX]{x1D6FF})$
. Thus we see that for small deviations from stability, the quasilinear turbulent heating acts to equalize the ion and electron thermal speeds and thus stabilize the mode.
For the case of pair plasmas (
$m_{i}=m_{e}$
) our quasilinear analysis indicates that turbulent heating driven by the electron drift wave acts to equalize the positron and electron temperatures. We can also use a symmetry of the gyrokinetic–Poisson system of equations to see how pair plasma heating depends on the temperature ratio. In particular, we consider how the equations are modified under an interchange of the electron and positron temperatures. First we note that the average over Larmor angle is unaffected by the interchange, as the Larmor radius of each particle is independent of temperature. Denoting the solutions when
$\unicode[STIX]{x1D70F}=\unicode[STIX]{x1D70F}_{0}$
as (
$h_{s}$
,
$\unicode[STIX]{x1D711}$
) and the solutions when
$\unicode[STIX]{x1D70F}=1/\unicode[STIX]{x1D70F}_{0}$
as (
$h_{s}^{\leftrightarrow }$
,
$\unicode[STIX]{x1D711}^{\leftrightarrow }$
), we have
$$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{\unicode[STIX]{x2202}h_{s}^{\leftrightarrow }}{\unicode[STIX]{x2202}t}}+v_{\Vert }{\displaystyle \frac{\unicode[STIX]{x2202}h_{s}^{\leftrightarrow }}{\unicode[STIX]{x2202}z}}+{\displaystyle \frac{c}{B}}\{\langle \unicode[STIX]{x1D711}^{\leftrightarrow }\rangle _{\boldsymbol{ R}},h_{s}^{\leftrightarrow }\}=-{\displaystyle \frac{q_{s^{\prime }}}{T_{s^{\prime }}}}{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}^{\leftrightarrow }\rangle _{\boldsymbol{ R}}}{\unicode[STIX]{x2202}t}}F_{0,s^{\prime }}+{\displaystyle \frac{c}{B}}{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}^{\leftrightarrow }\rangle _{\boldsymbol{ R}}}{\unicode[STIX]{x2202}y}}{\displaystyle \frac{\unicode[STIX]{x2202}F_{0,s^{\prime }}}{\unicode[STIX]{x2202}x}}+C[h_{s}^{\leftrightarrow }] & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$
and
$$\begin{eqnarray}\displaystyle & \displaystyle \int \text{d}^{3}v(h_{i}^{\leftrightarrow }-h_{e}^{\leftrightarrow })=\left(\sqrt{\unicode[STIX]{x1D70F}_{0}}+\sqrt{{\displaystyle \frac{1}{\unicode[STIX]{x1D70F}_{0}}}}-{\displaystyle \frac{T_{0}}{4\unicode[STIX]{x03C0}e^{2}}}\unicode[STIX]{x1D6FB}^{2}\right){\displaystyle \frac{e\unicode[STIX]{x1D711}^{\leftrightarrow }}{T_{0}}}, & \displaystyle\end{eqnarray}$$
with
$T_{0}\doteq (T_{e}T_{i})^{1/2}$
,
$s^{\prime }=i$
when
$s=e$
, and
$s^{\prime }=e$
when
$s=i$
.
There are two sets of solutions to these equations that can be obtained via symmetry considerations:
$h_{s}^{\leftrightarrow }(x,y,z,v_{\Vert },v_{\bot },t)=-h_{s^{\prime }}(-x,-y,z,v_{\Vert },v_{\bot },t)$
with associated potential
$\unicode[STIX]{x1D711}^{\leftrightarrow }(x,y,z,t)=\unicode[STIX]{x1D711}(-x,-y,z,t)$
, and
$h_{s}^{\leftrightarrow }(x,y,z,v_{\Vert },v_{\bot },t)=h_{s^{\prime }}(x,-y,z,v_{\Vert },v_{\bot },t)$
with associated potential
$\unicode[STIX]{x1D711}^{\leftrightarrow }(x,y,z,t)=-\unicode[STIX]{x1D711}(x,-y,z,t)$
. They can be verified by direct substitution into equations (3.13) and (3.14), noting that the interchange
$T_{i}\leftrightarrow T_{e}$
results in
$F_{0,i}\leftrightarrow F_{0,e}$
. As the interchanged solutions differ from the original solutions merely by sign and spatial orientation, the associated linear growth rates remain unchanged. For pair plasmas growth rates are thus symmetric under interchange of
$T_{i}$
and
$T_{e}$
, as seen in figure 2.
Substituting the interchanged solutions into the heating expression (2.6) gives
$$\begin{eqnarray}\displaystyle \displaystyle H_{s}^{\leftrightarrow } & = & \displaystyle {\displaystyle \frac{1}{V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}vq_{s}\overline{\left(h_{s}^{\leftrightarrow }{\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}^{\leftrightarrow }\rangle _{\boldsymbol{ R}_{s}}}{\unicode[STIX]{x2202}t}}\right)}\nonumber\\ \displaystyle & = & \displaystyle {\displaystyle \frac{1}{V}}\int \text{d}^{3}R\int _{\boldsymbol{ R}}\text{d}^{3}vq_{s^{\prime }}\overline{\left(h_{s^{\prime }}(\pm x,-y,z,v_{\Vert },v_{\bot },t){\displaystyle \frac{\unicode[STIX]{x2202}\langle \unicode[STIX]{x1D711}(\pm x,-y,z,t)\rangle _{\boldsymbol{R}_{s^{\prime }}}}{\unicode[STIX]{x2202}t}}\right)}\nonumber\\ \displaystyle & = & \displaystyle H_{s^{\prime }}=-H_{s},\end{eqnarray}$$
where the last equality follows from the fact that the species-summed heating is zero. Thus the heating of each species of a pair plasma is anti-symmetric under interchange of
$T_{i}$
and
$T_{e}$
, and there can be no turbulent heating when
$T_{i}=T_{e}$
.
3.3 Disparate thermal speeds
For plasmas with
$m_{i}\gg m_{e}$
, the ions and electrons will have disparate thermal speeds
$v_{\text{th},i}\ll v_{\text{th},e}$
for temperature ratios satisfying
$m_{e}/m_{i}\ll \unicode[STIX]{x1D70F}\ll m_{i}/m_{e}$
. In this case, we can look for solutions to the dispersion relation that satisfy
$|\unicode[STIX]{x1D701}_{e}|\ll 1\ll |\unicode[STIX]{x1D701}_{i}|$
. With this restriction, the plasma dispersion functions appearing in (3.4) can be greatly simplified. In particular, we use
$$\begin{eqnarray}\displaystyle & \displaystyle Z(\unicode[STIX]{x1D701}_{e})\approx \text{i}\sqrt{\unicode[STIX]{x03C0}} & \displaystyle\end{eqnarray}$$
and
$$\begin{eqnarray}\displaystyle & \displaystyle Z(\unicode[STIX]{x1D701}_{i})\approx -1/\unicode[STIX]{x1D701}_{i}, & \displaystyle\end{eqnarray}$$
giving the following approximate dispersion relation:
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D701}_{i}(1+\unicode[STIX]{x1D70F}+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2})+\unicode[STIX]{x1D701}_{i}\left(\unicode[STIX]{x1D701}_{e}-{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{e}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{e})\text{i}\sqrt{\unicode[STIX]{x03C0}}-\unicode[STIX]{x1D70F}\left(\unicode[STIX]{x1D701}_{i}+{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{i}}{2k_{\Vert }L_{n}}}\right)\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})=0. & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$
We next consider wavelengths much shorter than the ion Larmor radius but much longer than the electron Larmor radius; i.e.
$k_{\bot }\unicode[STIX]{x1D70C}_{i}\gg 1\gg k_{\bot }\unicode[STIX]{x1D70C}_{e}$
. In this limit we can approximate
$\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{e})\approx 1$
and
$\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{i})\approx 1/(\sqrt{\unicode[STIX]{x03C0}}k_{\bot }\unicode[STIX]{x1D70C}_{i})$
in (3.18) to obtain
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D701}_{i}(1+\unicode[STIX]{x1D70F}+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2})+\unicode[STIX]{x1D701}_{i}\left(\unicode[STIX]{x1D701}_{e}-{\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{e}}{2k_{\Vert }L_{n}}}\right)\text{i}\sqrt{\unicode[STIX]{x03C0}}-{\displaystyle \frac{\unicode[STIX]{x1D70F}}{\sqrt{\unicode[STIX]{x03C0}}}}\left({\displaystyle \frac{\unicode[STIX]{x1D701}_{i}}{k_{\bot }\unicode[STIX]{x1D70C}_{i}}}+{\displaystyle \frac{k_{y}}{k_{\bot }}}{\displaystyle \frac{1}{2k_{\Vert }L_{n}}}\right)=0.\quad & \displaystyle\end{eqnarray}$$
Seeking solutions for which
$\unicode[STIX]{x1D701}_{i}\sim k_{\bot }\unicode[STIX]{x1D70C}_{e}/k_{\Vert }L_{n}$
allows us to neglect the
$\unicode[STIX]{x1D701}_{e}$
and
$\unicode[STIX]{x1D701}_{i}/k_{\bot }\unicode[STIX]{x1D70C}_{i}$
terms. The resultant solution for
$\unicode[STIX]{x1D714}$
is
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D714}={\displaystyle \frac{\unicode[STIX]{x1D70F}k_{\Vert }v_{\text{th},i}}{\sqrt{\unicode[STIX]{x03C0}}}}{\displaystyle \frac{2k_{\Vert }L_{n}(1+\unicode[STIX]{x1D70F}+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2})\text{sgn}(k_{y})+\text{i}\sqrt{\unicode[STIX]{x03C0}}k_{\bot }\unicode[STIX]{x1D70C}_{e}}{(1+\unicode[STIX]{x1D70F}+k_{\bot }^{2}\unicode[STIX]{x1D706}_{e}^{2})^{2}(2k_{\Vert }L_{n})^{2}+\unicode[STIX]{x03C0}k_{\bot }^{2}\unicode[STIX]{x1D70C}_{e}^{2}}}. & \displaystyle\end{eqnarray}$$
In figure 1 we show an example comparison of the exact and analytical expressions for
$\unicode[STIX]{x1D714}$
for a plasma with
$k_{\bot }\unicode[STIX]{x1D706}_{e}=0$
and
$m_{i}/m_{e}=337\,824$
(corresponding to the mass of tungsten).
Figure 1. Comparison of exact (solid lines) and approximate analytical (dashed lines) growth rates for the case
$\unicode[STIX]{x1D70F}=1$
and
$m_{i}/m_{e}=337\,824$
(tungsten ions).
Equation (3.20) indicates that for
$k_{\Vert }>0$
there is an instability with peak growth rate
$\unicode[STIX]{x1D6FE}$
at wavelengths
$\boldsymbol{k}_{m}$
satisfying
$(\unicode[STIX]{x2202}\unicode[STIX]{x1D6FE}/\unicode[STIX]{x2202}k_{\Vert })|_{\boldsymbol{k}_{m}}=(\unicode[STIX]{x2202}\unicode[STIX]{x1D6FE}/\unicode[STIX]{x2202}k_{\bot })|_{\boldsymbol{k}_{m}}=0$
. These constraints give
$$\begin{eqnarray}\displaystyle & \displaystyle {\displaystyle \frac{2}{\sqrt{\unicode[STIX]{x03C0}}}}k_{\Vert ,m}|L_{n}|(1+\unicode[STIX]{x1D70F}+k_{\bot ,m}^{2}\unicode[STIX]{x1D706}_{e}^{2})=k_{\bot ,m}\unicode[STIX]{x1D70C}_{e}, & \displaystyle\end{eqnarray}$$
with the largest growth rate for
$k_{\bot ,m}\unicode[STIX]{x1D706}_{e}=0$
. Using these results in (3.20) gives
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D714}_{m}={\displaystyle \frac{\unicode[STIX]{x1D70F}}{4\sqrt{\unicode[STIX]{x03C0}}(1+\unicode[STIX]{x1D70F})}}{\displaystyle \frac{v_{\text{th},i}}{|L_{n}|}}(\text{i}+\text{sgn}(k_{y}L_{n})), & \displaystyle\end{eqnarray}$$
with
$\unicode[STIX]{x1D714}_{m}\doteq \unicode[STIX]{x1D714}(\boldsymbol{k}_{m})$
the complex frequency evaluated at the wavevector
$\boldsymbol{k}_{m}$
that maximizes the linear growth rate.
Plugging equation (3.22) for
$\unicode[STIX]{x1D714}_{m}$
into (3.7) for turbulent heating and using the appropriate approximations for
$\unicode[STIX]{x1D6F6}(k_{\bot }\unicode[STIX]{x1D70C}_{s})$
and
$Z(\unicode[STIX]{x1D701}_{s})$
gives
$$\begin{eqnarray}\displaystyle \displaystyle {\hat{H}}_{\boldsymbol{k},i} & = & \displaystyle -{\hat{H}}_{\boldsymbol{k},e}\approx {\displaystyle \frac{k_{y}\unicode[STIX]{x1D70C}_{e}}{k_{\Vert }v_{\text{th},i}}}\text{Re}(\unicode[STIX]{x1D714}){\displaystyle \frac{|L_{n}|}{L_{n}}}{\displaystyle \frac{\sqrt{\unicode[STIX]{x03C0}}}{4}}\nonumber\\ \displaystyle & = & \displaystyle {\displaystyle \frac{\unicode[STIX]{x1D70F}}{8\sqrt{\unicode[STIX]{x03C0}}}}>0.\end{eqnarray}$$
So ions are heated and electrons are cooled; i.e. the instability acts to equalize
$v_{\text{th},e}$
and
$v_{\text{th},i}$
and thus stabilize the mode.
Figure 2. Normalized linear growth rates maximized over
$k_{\Vert }$
and
$k_{y}$
(a) and corresponding turbulent heating (b) from GS2 simulations as a function of electron–ion temperature ratio for an electron–positron plasma. The quasilinear turbulent ion heating
${\hat{H}}_{i}$
(black circles and line) given by (3.7) is weighted by the linear growth rate
$\unicode[STIX]{x1D6FE}|L_{n}|/v_{0}$
to qualitatively estimate saturated fluctuation amplitudes. The turbulent ion heating
$(H_{i}|L_{n}|/nT_{0}v_{0})(|L_{n}|/\unicode[STIX]{x1D70C}_{0})^{2}$
obtained from nonlinear simulations is also given (red squares). Here
$\unicode[STIX]{x1D70C}_{0}=v_{0}/\unicode[STIX]{x1D6FA}_{i}$
, with
$v_{0}$
the geometric mean of the ion and electron thermal speeds.
4 Simulation results
We now provide numerical data to verify our analytical predictions. All simulations were conducted using the local, Eulerian gyrokinetic code
$\text{GS2}$
(Kotschenreuther, Rewoldt & Tang Reference Kotschenreuther, Rewoldt and Tang1995; Dorland et al.
Reference Dorland, Jenko, Kotschenreuther and Rogers2000) with kinetic electrons and a single ion species immersed in a straight, uniform magnetic field and with a cross-field density gradient. The fluctuations are constrained to be purely electrostatic, and we take
$k_{\bot }\unicode[STIX]{x1D706}_{e}=0$
. Each simulation uses 32 points along the magnetic field direction (
$z$
), and the velocity space is sampled on a polar grid (Barnes et al.
Reference Barnes, Dorland and Tatsuno2010b
), with 12 points in speed
$v$
and 16 points in pitch angle
$v_{\Vert }/v$
.
First we consider the case of an electron–positron plasma. The linear growth rates, maximized over
$k_{\Vert }$
and
$k_{\bot }$
, are plotted against
$T_{e}/T_{i}$
in figure 2. They are normalized by
$|L_{n}|/v_{0}$
, with
$v_{0}\doteq \sqrt{v_{\text{th},i}v_{\text{th},e}}$
, in order to make manifest the expected symmetry of the growth rates with respect to the interchange of
$T_{i}$
and
$T_{e}$
. As predicted, the plasma is stable only when
$T_{e}=T_{i}$
, and the growth rates are symmetric about
$T_{e}/T_{i}=1$
.
Also shown in figure 2 are the quasilinear heating estimates defined in (3.7) and the heating determined from nonlinear simulations as a function of
$T_{e}/T_{i}$
. The nonlinear simulations employed 47 modes (after de-aliasing) in both directions perpendicular to the magnetic field (
$x$
and
$y$
). The box sizes in both
$x$
and
$y$
were
$20\unicode[STIX]{x03C0}\unicode[STIX]{x1D70C}_{i}$
. The
$z$
the box size was
$40\unicode[STIX]{x03C0}L_{n}$
for the
$\unicode[STIX]{x1D70F}=0.5$
and
$\unicode[STIX]{x1D70F}=2$
cases and
$60\unicode[STIX]{x03C0}L_{n}$
for the
$\unicode[STIX]{x1D70F}=10$
case. A small amount of hyperviscosity, as well as a gyrokinetic collision operator (Abel et al.
Reference Abel, Barnes, Cowley, Dorland, Hammett and Schekochihin2008; Barnes et al.
Reference Barnes, Abel, Dorland, Ernst, Hammett, Ricci, Rogers, Schekochihin and Tatsuno2009) that includes the Landau test-particle operator and a model field-particle operator that conserves number, momentum and energy were used in all nonlinear simulations. The collisionality used for both species was
$\unicode[STIX]{x1D708}_{ii}|L_{n}|/v_{\text{th},i}=0.0005$
, with
$\unicode[STIX]{x1D708}_{ii}=\sqrt{2}\unicode[STIX]{x03C0}ne^{4}\ln \unicode[STIX]{x1D6EC}/(T_{i}^{3/2}m_{i}^{1/2})$
. All simulations were run until a steady state was reached, and time averages of the turbulent energy exchange were taken during the steady-state interval. An example time trace, along with the statistical average used, is shown in figure 3. As we showed analytically, the heating (both quasilinear and nonlinear) is antisymmetric about
$T_{e}/T_{i}=1$
, with the heating acting to equilibrate the electron and positron temperatures and thus shut off the linear instability.
Figure 3. Typical time traces (here for
$T_{e}/T_{i}=0.5$
) of the normalized turbulent heating from a nonlinear GS2 simulation, along with their time-averaged values (black dashed lines). Here
$\unicode[STIX]{x1D70C}_{0}=v_{0}/\unicode[STIX]{x1D6FA}_{i}$
, with
$v_{0}$
the geometric mean of the ion and electron thermal speeds.
Figure 4. Normalized linear growth rates maximized over
$k_{\Vert }$
and
$k_{y}$
(a) and corresponding turbulent heating (b) from GS2 simulations as a function of electron–ion temperature ratio for an electron–proton plasma. The quasilinear turbulent ion heating
${\hat{H}}_{i}$
(black circles and line) given by (3.7) is weighted by the linear growth rate
$\unicode[STIX]{x1D6FE}|L_{n}|/v_{\text{th},i}$
to qualitatively estimate saturated fluctuation amplitudes. The turbulent ion heating
$(H_{i}|L_{n}|/nT_{i}v_{\text{th},i})(|L_{n}|/\unicode[STIX]{x1D70C}_{i})^{2}$
obtained from nonlinear simulations is also given (red squares).
We next show the growth rates and turbulent ion heating as a function of electron–ion temperature for an electron–proton plasma in figure 4. Again we see that the plasma is unstable for
$v_{\text{th},e}\neq v_{\text{th},i}$
. Furthermore, the ion turbulent heating is positive definite for
$v_{\text{th},e}>v_{\text{th},i}$
, in agreement with our approximate analytic result.
5 Conclusions
The analytical and numerical results shown in this paper indicate that turbulent heating driven by the electron drift wave instability present in an inhomogeneous, magnetized plasma acts to stabilize the mode. Stabilization occurs when the ion and electron thermal speeds are equal. For a conventional mass-ratio plasma with
$T_{i}\sim T_{e}$
, this leads to the ions being heated and the electrons cooled; for a pair plasma, the turbulent heating acts to equalize the ion and electron temperatures.
While turbulent heating acts to stabilize the mode, it is not the only stabilization mechanism in the system. We showed that the instability drives particle transport that flattens the driving density gradient on a time scale at least as fast as the heating influences the electron–ion temperature ratio. Consequently, for the turbulent heating to fully set the thermal equilibrium there must be additional physics that fixes the density gradient; e.g. an external density source or a lowest order equilibrium set by gravitational forces. In this case, thermal equilibrium would correspond to
$T_{i}/T_{e}=\sqrt{m_{i}/m_{e}}$
.
Acknowledgements
The authors would like to thank F. I. Parra, A. A. Schekochihin and A. Zocco for useful discussions. M.B. was supported in part by STFC grant ST/N000919/1. The authors also acknowledge the use of ARCHER through the Plasma HEC Consortium EPRSC grant number EP/L000237/1 under project e281-gs2 and the use of the EUROfusion High Performance Computer (Marconi-Fusion) under project MULTEI.
