Coulomb collisions in strongly anisotropic plasmas I. Cyclotron cooling in electron–ion plasmas

Abstract

The behaviour of a collisional plasma that is optically thin to cyclotron radiation is considered, and the distribution functions accessible to it on the various time scales in the system are calculated. Particular attention is paid to the limit in which the collision time exceeds the radiation emission time, making the electron distribution function strongly anisotropic. Unusually for plasma physics, the collision operator can nevertheless be calculated analytically although the plasma is far from Maxwellian. The rate of radiation emission is calculated and found to be governed by the collision frequency multiplied by a factor that only depends logarithmically on plasma parameters.

1. Introduction

In the presence of a strong magnetic field, charged particles execute helical gyromotion around the magnetic field lines. The Lorentz force acts on particles in a direction perpendicular to both the magnetic field lines and the particle motion through them, causing the particles to accelerate, and in turn to release energy in the form of electromagnetic waves, ‘cyclotron radiation’. Plasmas with sufficiently large particle density, such as fusion plasmas, are optically thick to such emissions, meaning that any energy released through this process is simply reabsorbed back into the plasma. However, certain plasma systems are optically thin to cyclotron emission and this radiated energy can therefore be lost to the surroundings. These optically thin plasmas are the focus of this current work. It is the purpose of the present paper to show how these emissions lead to the manifestation of anisotropy in the plasma distribution function, and how the distribution evolves in these regimes.

In this work, paper (I) of the series, our contributions are: (a) to show that cyclotron emission results in strongly anisotropic distribution functions on the radiation time scale; (b) to calculate the evolution of the distribution function under collisional scattering which, in the absence of any radiation terms, acts to drive the plasma towards Maxwellian behaviour; (c) to show that this behaviour manifests itself under very general conditions; and (d) to present the applications and the limitations of this theory. In the companion work, paper II, we apply this theory to the first laboratory electron–positron plasma experiment. We begin in § 2 by briefly recapitulating the single-particle picture of an electron emitting cyclotron radiation, before using this description to explore the relaxation of the distribution function in a collisionless plasma. In § 3, collisions are introduced back into the picture and used to show that an initially isotropic distribution will relax to a distribution that is Maxwellian over the perpendicular velocity. In computing this distribution we introduce a free function $T_{\perp }$; in §§ 4 and 5, the evolution of this free function is explored. In § 6, we show that the conclusions drawn in the previous sections actually hold in quite general magnetic geometry; subsequently, we discuss some of the broad classes of applications in § 7.

2. The collisionless system in brief

Whilst we endeavour to provide a physical understanding of this theory in general magnetic geometries, let us first understand the basic premise before embarking on a more detailed calculation. To this end, we begin with a straight, constant magnetic field and, for the time being, neglect collisions in the plasma.

The fundamental aim of this section is thus to understand the following: (a) in the most simple case, cyclotron emission is responsible for an exponential decay of the perpendicular kinetic energy on the radiation time scale and (b) the inclusion of cyclotron emission, and the absence of collisions, results in a strongly anisotropic distribution function.

2.1. Cyclotron cooling

The emission of cyclotron radiation by particles in a plasma gives rise to a reaction force which must be included in the kinetic equation. The theory of relativistic plasmas, accounting for this radiation reaction by including the Abraham–Lorentz reaction force in the kinetic equation, has been developed by Andersson, Helander & Eriksson (2001), Hazeltine & Mahajan (2004) and subsequent authors. Here, we begin by specialising these results to a non-relativistic plasma.

The non-relativistic Abraham–Lorentz reaction force is given by

(2.1)\begin{equation} \boldsymbol{K} = \frac{{e}^{2}}{6{\rm \pi}\epsilon_{0}c^{3}}\dot{\boldsymbol{a}}, \end{equation}

where $\boldsymbol {a}$ is the acceleration vector and $e$ is the charge. Following the arguments of the aforementioned authors, it can be seen that, to leading order in a small-gyroradius expansion, the change in perpendicular energy, $w_{\perp }$, of a non-relativistic point charge as it accelerates in a magnetic field is given by Larmor's formula

(2.2a,b)\begin{equation} \left(\frac{{\mathrm{d}} w_{{\perp}}}{{\mathrm{d}} t}\right)_{\text{rad}} ={-} \frac{e^{2}a^{2}}{6{\rm \pi} \epsilon_{0}c^{3}}, \quad w_{{\perp}} = \frac{mv_{{\perp}}^{2}}{2}, \end{equation}

where $\boldsymbol {v}$ is the particle velocity vector and $m$ is the rest mass. Throughout this work, the subscripts $\perp$ and $\parallel$ indicate directions perpendicular and parallel to the magnetic field, respectively. The centripetal acceleration of a particle undergoing gyromotion is given by

(2.3)\begin{equation} a = \frac{eB}{m}v_{{\perp}}, \end{equation}

and hence it follows that

(2.4)\begin{equation} \left(\frac{{\mathrm{d}} w_{{\perp}}}{{\mathrm{d}} t}\right)_{\text{rad}} =- \frac{e^{4}B^{2}}{3{\rm \pi}\epsilon_{0}(mc)^{3}} w_{{\perp}}. \end{equation}

Simply put, the plasma will radiate its perpendicular kinetic energy on a time scale given by the radiation time:

(2.5)\begin{equation} \tau_{r} = \frac{3{\rm \pi}\epsilon_{0}(mc)^{3}}{e^{4}B^{2}}. \end{equation}

It is perhaps easiest to understand the influence of radiative emission when the change of energy is rewritten in terms of the perpendicular and parallel components of the velocity. We thus obtain

(2.6a,b)\begin{equation} \left(\frac{{\mathrm{d}} v_{{\perp}}}{{\mathrm{d}} t}\right)_{\text{rad}} ={-}\frac{v_{{\perp}}}{2\tau_{r}}, \quad \left(\frac{{\mathrm{d}} v_{{\parallel}}}{{\mathrm{d}} t}\right)_{\text{rad}} = 0. \end{equation}

Equations (2.6a,b) summarise the average effect of radiation reaction on the plasma. Fundamentally, it is these equations which will give rise to anisotropy in the distribution function.

It is pertinent to point out at this stage that collisions, whose effect is to isotropise the plasma, can mediate this cooling process due to the scattering of the velocity vector, as a result allowing the conversion of energy between perpendicular and parallel components. We will return to this discussion in § 3, following a discourse on the simpler collisionless case.

2.2. Collisionless kinetic equation

The kinetic equation governing the evolution of the distribution function $f(\boldsymbol {r},\boldsymbol {v},t)$ in the absence of collisions is given by

(2.7)\begin{equation} \frac{\partial f}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot} (\boldsymbol{v}f) + \boldsymbol{\nabla}_{v}\boldsymbol{\cdot} \left[\frac{q}{m}(\boldsymbol{E} + \boldsymbol{v}\times\boldsymbol{B})\,f + \left(\frac{{\mathrm{d}} \boldsymbol{v}}{{\mathrm{d}} t}\right)_{\text{rad}}f\right] = 0 , \end{equation}

where $\boldsymbol {\nabla }_{v}$ stands for the velocity derivative:

(2.8)\begin{equation} \boldsymbol{\nabla}_{v} = \hat{{\boldsymbol{x}}}\frac{\partial }{\partial v_{x}} + \hat{{\boldsymbol{y}}}\frac{\partial }{\partial v_{y}} + \hat{{\boldsymbol{z}}}\frac{\partial }{\partial v_{z}}. \end{equation}

Cyclotron emission is included in the kinetic equation through the inclusion of (2.6a,b) in the usual collisionless Vlasov equation.

We will immediately specialise to the case where the magnetic field is constant, $\boldsymbol {B} = B\hat {{\boldsymbol {b}}}$, and there is zero electric field, $\boldsymbol {E} = \boldsymbol {0}$. Note that, in this limit, we have assumed the plasma is homogeneous and there is no spatial dependence of the distribution function. We then write (2.7) in cylindrical $\boldsymbol {v}$-space coordinates $(v_{\perp },\alpha ,v_{\parallel })$ to obtain

(2.9)\begin{equation} \frac{\partial f}{\partial t} + \frac{1}{v_{{\perp}}}\frac{\partial }{\partial v_{{\perp}}}(v_{{\perp}}\dot{v}_{{\perp}}\,f) + \frac{\partial }{\partial \alpha} (\dot{\alpha}f) = 0, \end{equation}

where we have used that $\dot {{v}}_{\parallel } = 0$. The final term in this equation can then be eliminated by averaging over the gyroangle $\alpha$.

It is fruitful to write the resulting equation in terms of the perpendicular energy of the particles, so as to obtain

(2.10)\begin{equation} \frac{\partial f}{\partial t}- \frac{w_{{\perp}}}{\tau_{r}}\frac{\partial f}{\partial w_{{\perp}}} = \frac{f}{\tau_{r}}. \end{equation}

This equation can then be easily solved via the method of characteristics; the general solution is given by

(2.11)\begin{equation} f(w_{{\perp}},w_{{\parallel}},t) = F_{0}(w_{{\perp}}\,\mathrm{e}^{t/\tau_{r}},w_{{\parallel}}) \mathrm{e}^{t/\tau_{r}}. \end{equation}

That is, the distribution function will be a function of the parallel and perpendicular kinetic energies, the latter of which will decay exponentially quickly on the radiation time scale.

As an illustrative example, if the initial distribution is Maxwellian, then

(2.12)\begin{equation} f(w_{{\perp}},w_{{\parallel}},t) = n \left(\frac{m}{2{\rm \pi} T}\right)^{3/2} \exp \left[ - \frac{1}{T} (w_{{\perp}}\,\mathrm{e}^{t/\tau_{r}} + w_{{\parallel}}) + \frac{t}{\tau_{r}}\right], \end{equation}

which corresponds to a bi-Maxwellian

(2.13)\begin{equation} f(w_{{\perp}},w_{{\parallel}},t) = n \left(\frac{m}{2{\rm \pi} T_{{\perp}}^{2/3}T_{{\parallel}}^{1/3}}\right)^{3/2} \exp \left(-\frac{w_{{\perp}}}{T_{{\perp}}} - \frac{w_{{\parallel}}}{T_{{\parallel}}}\right) \end{equation}

with $T_{\parallel } = T = \text {const.}$ and

(2.14)\begin{equation} T_{{\perp}}(t) = T\, \mathrm{e}^{{-}t/\tau_{r}}. \end{equation}

At this stage of our investigation, we have been able to accomplish the first fundamental aim of this paper; namely, we now understand the following. Firstly, in the most simple case, cyclotron emission is responsible for an exponential decay of the perpendicular energy on the radiation time scale whilst the parallel energy is kept constant. Secondly, the inclusion of radiative emission, and the absence of collisions, results in a strongly anisotropic distribution function on the radiation time scale.

2.3. Validity of the collisionless approach

A caveat here is that, of course, many plasmas are in fact collisional. It is therefore expected that eventually collisions will come into play and must be taken into account.

Following Braginskii (1965), the conventional definitions of the electron–ion and electron–electron collision times are given by

(2.15a,b)\begin{equation} \tau_{ei} = \frac{6\sqrt{2}{\rm \pi}^{3/2}\epsilon_{0}^{2}{m_{e}^{1/2}}T_{e}^{3/2}}{ Ze^{4} n_{e} \ln \varLambda}, \quad \tau_{ee} = Z\tau_{ei} \end{equation}

respectively, where $\ln \varLambda$ is the Coulomb logarithm. We will also assume bulk neutrality throughout this work so that $n_{e} = Zn_{i}$. From the above expressions, we can see that the collision time will decrease as the plasma cools.

We might therefore envision the following scenario. A plasma might begin in a regime where the radiation time $\tau _{r}$ is smaller than the initial collision time $\tau _{c} = \text {min}(\tau _{ee},\tau _{ei})$ and hence the collisionless theory will be applicable at first. However, the collision time itself will decrease as the plasma cools and hence this assumption will be violated after sufficient time has elapsed. As a result, collisional effects must be taken into consideration.

We thus imminently turn our attention to the collisional problem. However, before doing so, let us allow ourselves a small digression to say something about the nature of cyclotron radiation.

2.4. Cyclotron cooling and plasma waves

Before we embark on solving the collisional kinetic equation, let us first turn our attention to the nature of the electromagnetic waves in play. The cyclotron radiation that is being generated has to propagate through our optically thin plasma as some electromagnetic plasma wave (Krall & Trivelpiece 1986).

For the reasons outlined in the previous section, we declare an interest in optically thin plasmas where the radiation time is (at least initially) much shorter than the collision time. That is, we declare an interest in the limit

(2.16)\begin{equation} \frac{\tau_{r}}{\tau_{c}} = \frac{mn\ln \varLambda}{{\rm \pi} \epsilon_{0}B^{2}}\left(\frac{c}{v_{\text{th}}}\right)^{3} = \frac{\ln \varLambda}{\rm \pi} \left(\frac{\omega_{p}}{\varOmega_{c}}\right)^{2}\left(\frac{c}{v_{\text{th}}}\right)^{3} \ll 1, \end{equation}

where have introduced the two characteristic frequencies in the plasma, namely the plasma frequency and the cyclotron frequency:

(2.17a,b)\begin{equation} \omega_{p}^{2} = \frac{ne^{2}}{m\epsilon_{0}}, \quad \varOmega_{c} = \frac{eB}{m}. \end{equation}

Thus, it is easy to see from this ordering that electromagnetic waves with frequencies comparable to the cyclotron frequency are almost light waves:

(2.18)\begin{equation} \omega^{2} ({\sim} \varOmega_{c}^{2}) = \omega_{p}^{2} + k^{2}c^{2} \implies \omega \simeq kc. \end{equation}

That is, our cyclotron radiation is simply photon emission which can immediately escape the plasma because of the optical thinness.

We now turn our attention to the collisional problem.

3. Lowest-order collisional kinetic equation

Our aim now is to deduce the class of distribution functions to which a radiating collisional plasma can relax. When collisions are retained, the kinetic equation satisfied by the distribution function $f$ in straight field lines is

(3.1a,b)\begin{equation} \frac{\partial \,f}{\partial t} - \frac{1}{\tau_{r}} \frac{\partial}{\partial \mu}(\mu \,f) = {C}(\,f), \quad \mu = \frac{mv^{2}}{2B}, \end{equation}

where we have opted to write this equation in terms of the magnetic moment $\mu$ instead of the perpendicular velocity. The distribution function can now be changed both via radiation (the second term on the left-hand side) and also via collisions (the term on the right-hand side). The Landau collision operator for Coulomb interaction (Landau 1936) is ${C} = {C}_{ee} + {C}_{ei}$ with

(3.2)\begin{equation} C_{ab}(\,f) = \sigma_{ab} \boldsymbol{\nabla} \boldsymbol{\cdot} \int ff^{\prime} {\boldsymbol{\mathsf{U}}} \boldsymbol{\cdot} \left(\boldsymbol{\nabla}_{v} \ln \,f - \boldsymbol{\nabla}_{v}^{\prime}\ln \,f^{\prime}\right) {\mathrm{d}}^{3}\boldsymbol{v}^{\prime} \end{equation}

and

(3.3a,b)\begin{equation} \sigma_{ee} = \sigma = \frac{n_{e}e^{4}\ln \varLambda}{8{\rm \pi}\epsilon_{0}^{2}m_{e}^{2}}, \quad \sigma_{ei} = Z\sigma, \end{equation}

if $f$ is normalised so that

(3.4)\begin{equation} \int f \, {\mathrm{d}}^{3}\boldsymbol{v} = 1. \end{equation}

Here we have also introduced ${\boldsymbol{\mathsf{U}}}$, the second-rank tensor:

(3.5)\begin{equation} {\boldsymbol{\mathsf{U}}}(\boldsymbol{u}) = \frac{u^{2} {\boldsymbol{\mathsf{I}}} - \boldsymbol{uu}}{u^{3}}, \end{equation}

where $\boldsymbol {u} = \boldsymbol {v} - \boldsymbol {v}^{\prime }$ is the difference in velocity vectors between colliding particles and ${\boldsymbol{\mathsf{I}}}$ is the identity matrix.

The collision frequency is of order $\sigma /v_{\text {th}a}^{3}$ and we are interested in the limit ${\sigma \tau _{r}/v_{\text {th}a}^{3} \ll 1}$, where $v_{\text {th}a}$ is the thermal velocity of species $a$. That is, for the reasons outlined in § 2, we are declaring an interest in optically thin plasmas where the radiation time is (at least initially) much shorter than the collision time. The key idea here is that, in this limit, by the time collisions enter the picture our plasma will have already radiated a large fraction of its perpendicular energy through the cyclotron cooling process outlined earlier. In this limit $f$ will be strongly anisotropic, due to the cyclotron cooling process, and we write

(3.6)\begin{equation} \boldsymbol{v}_{{\perp}} = \epsilon\boldsymbol{x}_{{\perp}}, \quad \boldsymbol{u} = \boldsymbol{v}_{{\parallel}} - \boldsymbol{v}_{{\parallel}}^{\prime} + \epsilon (\boldsymbol{x}_{{\perp}} - \boldsymbol{x}_{{\perp}}^{\prime}), \quad \epsilon \ll 1. \end{equation}

Our strategy is then to expand the collision operators in terms of this small parameter.

It is important to point out here that $\epsilon$, effectively a measure of $v_{\perp }/v_{\parallel },$ might itself be expected to change in time, which could have ramifications later in our analysis. However, these concerns can be mitigated by noting that the window of time in which this theory is valid is only a few collision times. Thus we expect $\epsilon$ will not change appreciably and will certainly remain $\epsilon \ll 1$. This is explored further in § 5.3. In what follows, $\epsilon$ is treated as a constant.

We begin with the electron–electron collision operator. To lowest order in $\epsilon$ we obtain

(3.7)\begin{equation} C_{ee}(\,f) \simeq \frac{\sigma}{v_{{\perp}}}\frac{\partial}{\partial v_{{\perp}}}\boldsymbol{v}_{{\perp}}\boldsymbol{\cdot} \boldsymbol{\nabla}_{v} f \int \frac{f^{\prime}}{u} \, {\mathrm{d}}^{3} \boldsymbol{v}^{\prime}. \end{equation}

Great care must be taken when dealing with the integral in (3.7) which is divergent when $v_{\parallel } = v_{\parallel }^{\prime }.$ We may evaluate

(3.8)\begin{equation} \int \frac{f^{\prime}}{u} \, {\mathrm{d}}^{3}\boldsymbol{v} = \int \frac{f^{\prime} \, {\mathrm{d}}^{3}\boldsymbol{v}}{\sqrt{(v_{{\parallel}} - v_{{\parallel}}^{\prime})^{2} + \epsilon^{2}(\boldsymbol{x} - \boldsymbol{x}^{\prime})^{2}}} = 2 \int_{-\infty}^{\infty} f^{\prime}(v_{{\parallel}},v_{{\perp}}^{\prime}) | \ln \epsilon| 2{\rm \pi} v_{{\perp}}^{\prime}\,{\mathrm{d}}v_{{\perp}}^{\prime} \equiv g(v_{{\parallel}}) \end{equation}

since, by partial integration, for any suitably well-behaved $f$ we have

(3.9)\begin{equation} \int_{0}^{\infty} \frac{f(x)}{\sqrt{x^{2} + \epsilon^{2}}} \, {\mathrm{d}} x = \left[f(x)\ln (x + \sqrt{x^{2} + \epsilon^{2}})\right]_{0}^{\infty} - \int_{0}^{\infty} \frac{{\mathrm{d}} f}{{\mathrm{d}} x} \ln (x + \sqrt{x^{2} + \epsilon^{2}}) \, {\mathrm{d}} x, \end{equation}

where

(3.10)\begin{equation} \int_{0}^{\infty} \frac{{\mathrm{d}} f}{{\mathrm{d}} x} \ln (x + \sqrt{x^{2} + \epsilon^{2}}) \, {\mathrm{d}} x \end{equation}

remains finite as $\epsilon \rightarrow 0$. Thus, it follows that

(3.11)\begin{equation} C_{ee}(\,f) \simeq \frac{\sigma g(v_{{\parallel}})}{v_{{\perp}}} \frac{\partial }{\partial v_{{\perp}}} v_{{\perp}} \frac{\partial f}{\partial v_{{\perp}}}. \end{equation}

The electron–electron collision term is larger than the electron–ion collision term,

(3.12)\begin{align} C_{ei}(\,f) &= Z \sigma \frac{\partial }{\partial \boldsymbol{v}}\boldsymbol{\cdot} \int \frac{v^{2}{\boldsymbol{\mathsf{I}}} - \boldsymbol{vv}}{v^{3}}ff^{\prime}\left(\boldsymbol{\nabla}_{v} \ln \,f - \boldsymbol{\nabla}_{v}^{\prime} \ln \,f^{\prime} \right) {\mathrm{d}}^{3}\,\boldsymbol{v}^{\prime} \end{align}

(3.13)\begin{align} &= Z\sigma \boldsymbol{\nabla}_{v}\boldsymbol{\cdot} \left(\frac{v^{2}{\boldsymbol{\mathsf{I}}}-\boldsymbol{vv}}{v^{3}} \boldsymbol{\cdot} \boldsymbol{\nabla}_{v} f\right) \end{align}

(3.14)\begin{align} &= Z\frac{\sigma}{v_{{\perp}}}\frac{\partial }{\partial v_{{\perp}}}\left(\frac{\boldsymbol{v}_{{\perp}}\boldsymbol{\cdot} \boldsymbol{\nabla}_{v} f}{v_{{\parallel}}}\right), \end{align}

provided $Z = O(1).$ This can be clearly seen by noting that

(3.15)\begin{equation} \frac{C_{ei}}{C_{ee}} \sim Z \epsilon. \end{equation}

The lowest-order kinetic equation thus becomes

(3.16)\begin{equation} -\frac{1}{v_{{\perp}}}\frac{\partial }{\partial v_{{\perp}}}\left(\frac{v_{{\perp}}^{2}}{2}f\right) = \frac{\tau_{r}\sigma g(v_{{\parallel}})}{v_{{\perp}}} \frac{\partial }{\partial v_{{\perp}}} v_{{\perp}} \frac{\partial f}{\partial v_{{\perp}}}, \end{equation}

which can be integrated to give

(3.17a,b)\begin{equation} f(v_{{\parallel}},v_{{\perp}},t) = C(v_{{\parallel}},t) \exp \left(-\frac{mv_{{\perp}}^{2}}{2T_{{\perp}}(v_{{\parallel}},t)}\right), \quad T_{{\perp}} = 2\sigma \tau_{r} m g(v_{{\parallel}},t). \end{equation}

The integration constant $C$ is determined by the requirement that

(3.18)\begin{equation} N(v_{{\parallel}},t) = \int_{0}^{\infty} f 2{\rm \pi} v_{{\perp}} \, {\mathrm{d}} v_{{\perp}} = \frac{2{\rm \pi} T_{{\perp}}}{m}C = 4{\rm \pi} \sigma \tau_{r}g(v_{{\parallel}}) C \end{equation}

should equal $g(v_{\parallel })/2|\ln \epsilon |$, which gives $C = 1/(8{\rm \pi} \sigma \tau _{r}|\ln \epsilon |)$ and thus

(3.19)\begin{equation} f(v_{{\parallel}},v_{{\perp}},t) = \frac{1}{8{\rm \pi}\sigma\tau_{r}|\ln \epsilon|}\exp\left(-\frac{mv_{{\perp}}^{2}}{2T_{{\perp}}(v_{{\parallel}},t)}\right). \end{equation}

It is interesting to note that the distribution of particles over perpendicular velocities is Maxwellian, but not for the usual reason. Normally, this happens because of a balance between two terms in the collision operator, describing friction (drag) and energy diffusion, respectively. The former term slows particles down and the latter increases their average energy. The two terms balance exactly for a Maxwellian distribution. In (3.19), radiative energy loss has replaced the collisional friction, and the result is a Maxwellian in $v_{\perp }$ with a different (lower) perpendicular temperature than in the purely collisional case.

We note that the distribution function only depends on $v_{\parallel }$ through the function $T_{\perp }(v_{\parallel },t)$. We can determine this function by taking a moment of the kinetic equation.

4. Evolution of the perpendicular density

The moment $N$ defined in (3.18) can be seen to satisfy

(4.1)\begin{equation} \frac{\partial N}{\partial t} = \int_{0}^{\infty} C(\,f) 2{\rm \pi} v_{{\perp}} \, {\mathrm{d}} v_{{\perp}}. \end{equation}

In solving this equation, we obtain the temperature $T_{\perp }$ through (3.18) and can thus ascertain the distribution function $f$ through use of (3.19). Let us first remark on the nature of (4.1). In deriving (4.1), we have integrated the kinetic equation (3.1a,b) over perpendicular velocities. In doing so, we have eliminated the radiation term and instead isolated the slow evolution of the distribution over parallel velocities. We can think of the evolution of an initial isotropic distribution function, in the collisional regime, as being partitioned into two disparate stages:

1. (i) Firstly, the distribution function will quickly radiate most of its perpendicular energy and relax to a form which is Maxwellian over perpendicular velocity as per (3.19).

2. (ii) Secondly, on a longer time scale, the distribution over parallel velocities will evolve according to (4.1).

It is this slow behaviour (ii) that we are now interested in probing. We will soon see that the precise nature of this evolution depends on the specific ordering satisfied by the plasma parameters.

In order to ascertain the time evolution of $N$, we thus require the following moments of the collision operators:

(4.2a,b)\begin{equation} I_{ee} \equiv \int C_{ee}(\,f) \, {\mathrm{d}}^{2} \boldsymbol{v}_{{\perp}}, \quad I_{ei} \equiv \int C_{ei}(\,f) \, {\mathrm{d}}^{2} \boldsymbol{v}_{{\perp}}. \end{equation}

It can be shown that the leading-order contributions to these two terms are given by

(4.3a,b)\begin{equation} I_{ee} ={-}\sigma |\ln \epsilon| \frac{\partial ^{2}N^{2}}{\partial v_{{\parallel}}^{2}}, \quad I_{ei} =Z\sigma \frac{\partial}{\partial v_{{\parallel}}}\left(\frac{N}{v_{{\parallel}}|v_{{\parallel}}|}\right). \end{equation}

The details of the calculation are relegated to appendix A. It then follows that (4.1) becomes

(4.4)\begin{equation} \frac{\partial N}{\partial t} ={-}\sigma|\ln\epsilon| \frac{\partial ^{2}N^{2}}{\partial v_{{\parallel}}^{2}} + 2Z\sigma\frac{\partial}{\partial v_{{\parallel}}}\left(\frac{N}{v_{{\parallel}}|v_{{\parallel}}|}\right) \end{equation}

or, more compactly, we can write

(4.5a–c)\begin{equation} \frac{\partial N}{\partial \tau} = \frac{\partial N}{\partial v_{{\parallel}}}\left(\frac{N}{v_{{\parallel}}|v_{{\parallel}}|} - \alpha \frac{\partial N^{2}}{\partial v_{{\parallel}}}\right), \quad \tau = \frac{3}{4}\sqrt{\rm \pi} v_{\text{th}}^{3} \frac{t}{\tau_{ei}}, \quad \alpha = |\ln \epsilon| \frac{1}{2Z}. \end{equation}

Once again, we have been able to successfully reduce our integro-differential evolution equation (4.1) into an equation involving only differential operators (4.5a–c). The two terms on the right-hand side of (4.5a–c) have very different character. The first one describes friction on the electrons from collisions with the ions, which causes the former to slow down in the parallel direction. As usual for Coulomb collisions, the collision frequency decreases with increasing speed, which here causes the singularity at $v_{\parallel } = 0.$ The second term describes electron–electron collisions, whose effect is to scatter the velocity in a peculiar way, causing ‘anti-diffusion’ in $v_{\parallel }$. This phenomenon can be understood as follows. Most electron–electron collisions occur between particles with similar velocities $v_{\parallel } \simeq v_{\parallel }^{\prime }$, and the local (in parallel velocity space) collision frequency is thus proportional to $N$, making this term quadratic in $N$ (particles with given $v_{\parallel }$ mostly collide with each other). The effect of such collisions is to convert parallel kinetic energy to perpendicular energy, and thus to increase $T_{{\perp }}$, which according to (3.18) leads to an increase in $N$. A local accumulation of particles somewhere in parallel velocity space will thus tend to grow at the expense of neighbouring regions, making the distribution function undergo ‘anti-diffusion’ in $v_{\parallel }$.

Despite its relatively simple form, (4.5a–c) is a nonlinear parabolic partial differential equation. It seems impossible to make analytical progress in full generality here, but nevertheless some insight into the system can be gleaned by appealing to various limiting forms.

5. Limiting forms of the perpendicular density

We note that, from (4.3a,b), it is easy to see that the moments of the collision operator (which gives rise to the different terms in (4.5a–c)) are ordered:

(5.1)\begin{equation} \frac{I_{ei}}{I_{ee}} \sim \frac{Z}{|\ln \epsilon|}. \end{equation}

This is to be compared and contrasted to the ordering of the collision operators themselves:

(5.2)\begin{equation} \frac{C_{ei}}{C_{ee}} \sim Z \epsilon. \end{equation}

Based on these different orderings, we can turn to solving equation (4.5a–c) in two distinct limits:

(5.3)\begin{gather} (A): \quad|\ln \epsilon| \ll Z \ll \frac{1}{\epsilon}, \end{gather}

(5.4)\begin{gather}(B): \quad Z \ll |\ln \epsilon| \ll \frac{1}{\epsilon}. \end{gather}

In both of these limits, we are justified in dropping the electron–ion collision operator at lowest order as discussed in § 3. This approximation gives rise to the rapid evolution of the distribution function towards (3.19) at lowest order.

5.1. Limit (A): neglecting the contribution from electron–electron collisions

In limit (A), we are justified in neglecting the contribution of $I_{ee}$ whilst retaining that of $I_{ei}$, noting that this is still consistent with the derivation of (3.19) providing the required ordering is satisfied. In this instance, (4.5a–c) becomes

(5.5)\begin{equation} \frac{\partial N}{\partial \tau} = \frac{\partial}{\partial v_{{\parallel}}}\left(\frac{N}{v_{{\parallel}}^{2}}\right). \end{equation}

The equation has been reduced to a quasilinear partial differential equation and so is amenable to solution by the method of characteristics. We thus obtain the general solution

(5.6)\begin{equation} N(v_{{\parallel}},\tau) = v_{{\parallel}}^{2}F\left(\frac{v_{{\parallel}}^{3}}{3} + \tau\right). \end{equation}

So, if we start from a Maxwellian distribution initially

(5.7)\begin{equation} N(v_{{\parallel}},0) = \left(\frac{m}{2{\rm \pi} T_{{\parallel}}}\right)^{1/2}\exp\left(-\frac{mv_{{\parallel}}^{2}}{2T_{{\parallel}}}\right), \end{equation}

then

(5.8)\begin{equation} N(v_{{\parallel}},\tau) = \left(\frac{m}{2{\rm \pi} T_{{\parallel}}}\right)^{1/2} \frac{v_{{\parallel}}^{2}}{(v_{{\parallel}}^{3} + 3\tau)^{2/3}} \exp\left(-\frac{m}{2T_{{\parallel}}}(v_{{\parallel}}^{3} + 3\tau)^{2/3}\right). \end{equation}

From this equation, we can also obtain an expression for the perpendicular temperature:

(5.9)\begin{equation} T_{{\perp}}(v_{{\parallel}},\tau) = 4\sigma \tau_{r}m|\ln \epsilon|\left(\frac{m}{2{\rm \pi} T_{{\parallel}}}\right)^{1/2} \frac{v_{{\parallel}}^{2}}{(v_{{\parallel}}^{3} + 3\tau)^{2/3}} \exp\left(-\frac{m}{2T_{{\parallel}}}(v_{{\parallel}}^{3} + 3\tau)^{2/3}\right). \end{equation}

It is important to recall that this solution was for $v_{\parallel } > 0$, which means that, despite appearances, both $T_{\perp }$ and $N$ are continuous. For $v_{\parallel } < 0$ we obtain the same solution as above with $v_{\parallel } \rightarrow -v_{\parallel }$.

5.1.1. Rate of energy loss

A quantity of direct experimental interest is the amount of power radiated by the plasma. This can be calculated by looking at the rate of change of energy loss:

(5.10)\begin{equation} P ={-}\frac{{\mathrm{d}} W}{{\mathrm{d}} t}, \end{equation}

where we have defined the total thermal energy stored in the plasma (per unit volume):

(5.11)\begin{equation} W = \int \tfrac{1}{2} mnv^{2} f \, {\mathrm{d}}^{3}\boldsymbol{v}. \end{equation}

The energy stored can be seen to satisfy the equation

(5.12a,b)\begin{equation} \frac{{\mathrm{d}} W}{{\mathrm{d}} t} ={-}\frac{W_{{\perp}}}{\tau_{r}}, \quad W_{{\perp}} = \int \tfrac{1}{2}mnv_{{\perp}}^{2} f \, {\mathrm{d}}^{3}\boldsymbol{v}. \end{equation}

Thus, it suffices to calculate $W_{\perp }$ in order to find the power radiated by the plasma. This is given by evaluating

(5.13)\begin{align} W_{{\perp}} &= \int \frac{1}{2}mnv_{{\perp}}^{2} f \, {\mathrm{d}}^{3} \boldsymbol{v} \nonumber\\ &= \int \frac{mn}{8\sigma\tau_{r}|\ln\epsilon|} v_{{\perp}}^{3} \exp \left(-\frac{mv_{{\perp}}^{2}}{2T_{{\perp}}(v_{{\parallel}},t)}\right) {\mathrm{d}} v_{{\parallel}} \,{\mathrm{d}} v_{{\perp}} \nonumber\\ &= \int \frac{nT_{{\perp}}(v_{{\parallel}},t)^{2}}{4\sigma\tau_{r}m|\ln\epsilon|} \, {\mathrm{d}}v_{{\parallel}}, \end{align}

where we have used (3.19) to carry out the $v_{\perp }$ integral. We arrive at

(5.14)\begin{equation} \frac{W_{{\perp}}(y)}{W_{{\perp}}(0)} = \frac{2}{\sqrt{\rm \pi}}\int_{0}^{\infty} \frac{x^{4}}{(x^{3}+y)^{4/3}}\exp (-(x^{3} + y)^{2/3})\, {\mathrm{d}} x, \end{equation}

where we have normalised $v_{\parallel }$ and $\tau$ to the electron thermal velocity by writing

(5.15a–c)\begin{equation} v_{{\parallel}} = v_{\text{th}}x, \quad 3\tau = v_{\text{th}}^{3}y, \quad v_{\text{th}} = \sqrt{\frac{2T_{{\parallel}}}{m}}. \end{equation}

The integral (5.14) cannot be expressed in terms of standard mathematical functions. Instead, we must turn to numerically evaluating this function. Firstly, however, we can determine the behaviour of (5.14) in the long-time and short-time limits.

In the limit $y \gg 1$ we obtain (see appendix 1)

(5.16)\begin{equation} \frac{W_{{\perp}}(y)}{W_{{\perp}}(0)} \simeq \left(\frac{2}{3{\rm \pi}}\right)^{8/9} \varGamma \left(\frac{5}{3}\right) \left(\frac{t}{\tau_{ei}}\right)^{{-}7/9} \exp\left( -\left(\frac{9\sqrt{\rm \pi}}{4}\frac{t}{\tau_{ei}}\right)^{2/3}\right). \end{equation}

In the limit $y \ll 1$ we obtain (see appendix 2)

(5.17)\begin{equation} \frac{W_{{\perp}}(t)}{W_{{\perp}}(0)} \simeq 1 - \frac{\varGamma(5/3) \varGamma(2/3)}{\varGamma(7/3)} \left(\frac{128}{3{\rm \pi}}\right)^{1/6} \left(\frac{t}{\tau_{ei}}\right)^{1/3}. \end{equation}

In figure 1 we show the full numerical solution of (5.14) as well as the analytic solution for the long-term and short-term behaviour, given by (5.16) and (5.17), respectively. Some caution must be taken in interpreting this figure, and one must remember that the solution of (5.14) is only valid on time scales longer than the radiation time. Note that the quantity plotted on the vertical axis, $W_{\perp }(y)/W_{\perp }(0),$ is directly proportional to the emitted power as can be seen from (5.12a,b). The abscissa is simply a scaled time coordinate: $y = (9\sqrt {{\rm \pi} }/4)t/\tau _{ei}$.

Figure 1. The ratio of perpendicular energy at normalised time $y = (9\sqrt {{\rm \pi} }/4)t/\tau _{ei}$ to the initial perpendicular energy, a quantity that is proportional to the rate of energy loss in the plasma. The blue curve is the calculated rate of energy loss from (5.14). The red curve is the asymptotic solution for $y \gg 1$ given by (5.16). The purple curve is the asymptotic solution for $y \ll 1$ given by (5.17). This figure is for the perpendicular energy loss in limit (A).

Although at a first glance it may seem as though the energy loss rate would be faster than exponential, this is of course not the case when the solution is restricted to the region where the orderings are valid.

5.2. Limit (B): neglecting the contribution from electron–ion collisions

In limit (B), when $n_{i} Z^{2} / n_{e} \ll | \ln \epsilon |$, we are justified in neglecting the contribution from $I_{ei}$ and obtain a relatively simple equation:

(5.18)\begin{equation} \frac{\partial N}{\partial \tau} ={-}\alpha \frac{\partial ^{2} N^{2}}{\partial v_{{\parallel}}^{2}}. \end{equation}

This is a nonlinear parabolic partial differential equation which essentially describes backwards diffusion in $v_{\parallel }$ space. Indeed, this is simply the backwards heat equation in a medium where the diffusion coefficient varies linearly with temperature. As such, this problem is ill-posed.

A simple solution exists (given suitable initial data) but the general solution will not depend continuously on the auxiliary data as it describes a reverse diffusion process meaning that arbitrarily small perturbations in the initial conditions will be amplified and can lead to even infinitesimal-time singularities. This means that it is possible to prepare arbitrarily small perturbations to the initial conditions (likely in any given Sobolev norm), which are sufficiently wiggly in higher derivatives to ensure the solution blows up before any given time $t > 0$.

In order to find this singular solution, we consider the problem

(5.19a,b)\begin{equation} \frac{\partial N}{\partial \tau} - \beta \frac{\partial}{\partial v_{{\parallel}}}\left(N\frac{\partial N}{\partial v_{{\parallel}}}\right) = 0, \quad \beta ={-} 2\alpha \end{equation}

with

(5.20)\begin{equation} \int_{-\infty}^{\infty} N(v_{{\parallel}}, \tau) \, {\mathrm{d}} v_{{\parallel}} = 1. \end{equation}

Equation (5.19a,b) admits a similarity solution. To see this, we introduce a dilation transformation:

(5.21a–c)\begin{equation} z = \delta^{a}v_{{\parallel}}, \quad s = \delta^{b}\tau, \quad v(s,z) = \delta^{c}N (\delta^{{-}a}z,\delta^{{-}b}s). \end{equation}

Our problem then becomes

(5.22)\begin{equation} \delta^{b-a} \frac{\partial v}{\partial s} - \beta \delta^{2a} \frac{\partial ^{2}}{\partial z^{2}}\left(\frac{1}{2}\delta^{{-}2c}v^{2}\right) = 0. \end{equation}

From this, we see that (5.19a,b) is invariant under the dilation transformation provided $b -c = 2(a-c)$, i.e. $c = 2a - b$.

We therefore pose the self-similar ansatz

(5.23a,b)\begin{equation} N(v_{{\parallel}},\tau) = \tau^{(2a-b)/c}y(\xi), \quad \xi = \frac{v_{{\parallel}}}{\tau^{a/b}}. \end{equation}

The condition (5.20) then gives

(5.24)\begin{equation} \int_{-\infty}^{\infty} N(v_{{\parallel}},\tau) \, {\mathrm{d}} v_{{\parallel}} = 1 \implies \tau^{3(a/b)-1} \int_{-\infty}^{\infty} y(\xi) \, {\mathrm{d}} \xi = 1. \end{equation}

This condition must hold for all $\tau > 0$ and hence we have $3a=b$. So we can write

(5.25a,b)\begin{equation} N(v_{{\parallel}},\tau) = \tau^{{-}1/3}y(\xi), \quad \xi = \frac{v_{{\parallel}}}{\tau^{1/3}}. \end{equation}

It is of course no coincidence to have arrived at this particular self-similar ansatz. We could have also posited this solution simply by noting that the quantity $\xi = v_{\parallel }^{3}/\tau$ is dimensionless, which we could have deduced from (5.6).

The derivatives then transform as

(5.26)\begin{gather} \frac{\partial N}{\partial \tau} ={-}\frac{1}{3}\tau^{{-}4/3} (y + \xi y^{\prime}), \end{gather}

(5.27)\begin{gather}\frac{\partial ^{2}}{\partial v_{{\parallel}}^{2}}\left(\frac{1}{2}N^{2}\right) = \frac{1}{2}\tau^{{-}2/3} (2yy^{\prime\prime} + 2 y^{\prime2} ) \tau^{{-}2/3}, \end{gather}

where the prime notation now denotes differentiation with respect to $\xi$.

Equation (5.19a,b) becomes

(5.28)\begin{equation} \beta y y^{\prime\prime} + \beta y^{\prime 2} + \tfrac{1}{3}y + \tfrac{1}{3}\xi y^{\prime} = 0, \end{equation}

or equivalently

(5.29)\begin{equation} 3\beta (yy^{\prime})^{\prime} + (\xi y)^{\prime} = 0. \end{equation}

Direct integration gives

(5.30)\begin{equation} 3\beta yy^{\prime} + \xi y = k_{0}. \end{equation}

The integration constant $k_{0}$ is chosen such that the flux of particles $N ({{\mathrm {d}} N}/{{\mathrm {d}} v_{\parallel }})$ is continuous at the edge of the support where $N$ goes to zero. In our case, this is the same constant that one gets from imposing zero flux at the symmetry point because the flux has to vanish either when $N=0$ (at the edge) or when ${{\mathrm {d}} N}/{{\mathrm {d}} v_{\parallel }} = 0$ (at a symmetry point).

This equation can then be integrated once again:

(5.31)\begin{equation} y(\xi) ={-}\frac{1}{6\beta}\xi^{2} + k_{1}. \end{equation}

In order to satisfy the requirement that $N \rightarrow 0$ as $| v_{\parallel } | \rightarrow \infty ,$ we write the solution as

(5.32)\begin{equation} N(v_{{\parallel}},\tau) = \begin{cases} \dfrac{\tau^{{-}1/3}}{6\beta}\left(A^{2} - \left(\dfrac{v_{{\parallel}}}{\tau^{1/3}}\right)^{2}\right) & |v_{{\parallel}}| < A\tau^{1/3} \\ 0 & |v_{{\parallel}}| > A\tau^{1/3}, \end{cases} \end{equation}

where $A = \sqrt [3]{9\beta /2}$ so as to satisfy (5.20).

We can thus obtain a weak, in the sense that the first derivative is discontinuous, solution of (4.5a–c) in the limit where electron–ion collisions are neglected. It is imperative to mention that this solution holds only when $\tau < 0.$ In essence, what we have done in this section is to solve an ill-posed backwards diffusion equation (5.18) by solving instead a well-posed forward diffusion equation (5.19a,b) and running the solution backwards in time. The rationale in adopting this approach is to motivate the correct boundary conditions to apply to the equation.

As this problem is ill-posed, our self-similar solution requires very specific initial conditions. An arbitrarily small perturbation (in any given Sobolev norm) can cause the solution to fail and develop an infinitesimal-time singularity. Of course the self-similar solution is still of value, and the formation of singularities is somewhat fictitious since our assumption that $\boldsymbol {\nabla }_{{\parallel }} \ll \boldsymbol {\nabla }_{{\perp }}$ will certainly break down before any singularity arises. The solution obtained is of the correct general form (i.e. it attains a maximum for $v_{\parallel } = 0,$ and satisfies $N \rightarrow 0$ as $|v_{\parallel }|\rightarrow 0$) as to possibly be approached from more general initial conditions. Such solutions are not uncommon in plasma physics and the generation of sharp structures by resistive or ambipolar diffusion has also been investigated (Low 1973).

It is also fitting to remark here that including higher powers of the expansion parameter $\epsilon$ in our expansion of the collision operator (i.e. in the calculation of $I_{ee}$) would usually bring higher-order derivatives of $N$ into (5.18) which would regularise the partial differential equation. This is a possible avenue of further exploration with the caveat that the underlying equations would then likely be rendered too difficult to solve with analytical techniques.

As before, one can calculate the perpendicular temperature:

(5.33)\begin{equation} T_{{\perp}}(v_{{\parallel}},\tau) = \begin{cases} 2\sigma\tau_{r}m|\ln\epsilon|\dfrac{\tau^{{-}1/3}}{3\beta}\left(A^{2} - \left(\dfrac{v_{{\parallel}}}{\tau^{1/3}}\right)^{2}\right) & |v_{{\parallel}}| < A\tau^{1/3} \\ 0 & |v_{{\parallel}}| > A\tau^{1/3}. \end{cases} \end{equation}

5.2.1. Rate of energy loss

Following the previous subsection, one can calculate the rate at which energy is radiated from the plasma by first calculating $W_{\perp }$ and then appealing to equation (5.12a,b). This is given by

(5.34)\begin{equation} W_{{\perp}} = \int \frac{nT_{{\perp}}(v_{{\parallel}},t)^{2}}{4\sigma\tau_{r}m|\ln \epsilon|} \, {\mathrm{d}}v_{{\parallel}} = 4 mn \sigma \tau_{r} |\ln \epsilon| \left(\frac{3}{4\tau A^{3}}\right)^{2} \int_{-\infty}^{\infty} (A^{2} \tau^{2/3} - v_{{\parallel}}^{2})^{2}\, {\mathrm{d}}v_{{\parallel}}, \end{equation}

which can be evaluated exactly to give

(5.35)\begin{equation} W_{{\perp}} = 8 m n\sigma \tau_{r}|\ln \epsilon| \left(\frac{3}{4\tau A^{3}}\right)^{2} \int_{0}^{A\tau^{1/3}}(A^{2}\tau^{2/3} - v_{{\parallel}}^{2}) \, {\mathrm{d}}v_{{\parallel}} = \frac{12}{5} \frac{mn \sigma \tau_{r} |\ln \epsilon|}{A\tau^{1/3}}. \end{equation}

This leads to an energy loss rate of

(5.36)\begin{equation} -\frac{{\mathrm{d}} W}{{\mathrm{d}} t} = \frac{12}{5} \left(\frac{1}{2Z}\right)^{1/3} \frac{mn\sigma^{2/3}|\ln \epsilon|}{At^{1/3}}, \end{equation}

which can be written in terms of the Braginskii electron–electron collision time $\tau _{ee}$ as

(5.37)\begin{equation} -\frac{{\mathrm{d}} W}{{\mathrm{d}} t} ={-} \frac{6}{5}{\rm \pi}^{1/3} n_{e}T_{{\parallel} e} |\ln \epsilon|^{2/3} \left(\frac{1}{t\tau_{ee}^{2}}\right)^{1/3}. \end{equation}

At first glance, this equation might appear somewhat peculiar. We must first of course recall that our solution was valid only when time is run backwards, and hence the term $t$ appearing on the right-hand side is negative; as a result, the total energy does indeed decrease as expected.

5.3. Validity of the collisional approach; a self-accelerating process

We have found the solution to (4.5a–c) in two different limits based on the plasma parameters, thus elucidating the plasma distribution function in the collisional regime (3.19). Throughout this work, a number of orderings have been invoked and we thus turn our attention to the validity of our solution; specifically, when it is valid to treat $\epsilon$ as constant. Thus far we have focused on the evolution of the plasma on time scales $t$ with

(5.38)\begin{equation} 0 < \tau_{r} < t \sim \tau_{c}. \end{equation}

In this regime, we treat $\epsilon$ as constant seeing as though any variation is very small. However, as per our earlier remarks in § 2.3, we must take care to recall that $\tau _c$ is not constant and will itself decrease as the plasma cools via cyclotron emission. Our entire process (radiative cooling and collisional scattering) is therefore self-accelerating. Collisional scattering converts parallel kinetic energy to perpendicular kinetic energy. This perpendicular energy is then radiated, lowering the temperature of the plasma and in turn increasing the collisional scattering rate. As a result, the plasma will only remain in this regime for a few collision times, before entering a regime where

(5.39)\begin{equation} \tau_{c} < \tau_{r} < t. \end{equation}

In this final regime, any remaining parallel kinetic energy will be converted to perpendicular kinetic energy and then radiated. The total energy of the plasma will thus decay exponentially quickly on the radiation time scale.

In essence, we can partition the evolution of the plasma energy into three distinct regimes:

1. (i) Initially, $0 < \tau _{r} < t < \tau _{c}.$ In this regime the behaviour is collisionless as described in § 2. The perpendicular energy decays exponentially quickly on the radiation time scale.

2. (ii) After some time, $\tau _{r} < t \sim \tau _{c}$ and collisional scattering becomes important. This regime is what we have studied in §§ 3–5. This is a self-accelerating process and the plasma will only remain in this regime for a few collision times.

3. (iii) Eventually, $\tau _{c} < \tau _{r} < t.$ In this regime, the distribution function is isotropic (Maxwellian), and any remaining parallel kinetic energy will be converted to perpendicular kinetic energy via collisional scattering and then radiated.

We have now accomplished our second aim; namely, we have been able to calculate the evolution of the distribution function under collisional scattering. Unusually for plasma physics, this was possible even though the distribution function is far from Maxwellian. Thus far, we have made the restrictive assumption of straight-field-line geometry. We will now ask ourselves whether we can relax this assumption.

6. General magnetic geometry

Lurking within this theory was an assumption of straight-field-line geometry which considerably simplified the derivations involved. It is clear that this assumption will generally not hold in any plasma of physical interest, so let us address this point.

We begin in general magnetic geometry, although we will later find it necessary to discuss trapped and circulating particles. The kinetic equation becomes

(6.1)\begin{equation} \frac{\partial f}{\partial t} + v_{{\parallel}}\frac{\partial f}{\partial l} - \frac{\partial }{\partial \mu}\left(\frac{\mu}{\tau_{r}}f\right) - \frac{\mu}{m}\boldsymbol{\nabla}_{{\parallel}}B \frac{\partial f}{\partial v_{{\parallel}}} = {C}_{ee} + {C}_{ei}, \end{equation}

where $l$ parametrises the length of a magnetic field line and we retain the assumption that the parallel electric field is negligible. Here we can exploit the existence of two time scales by expanding in the small parameter

(6.2)\begin{equation} \epsilon_{1} \equiv \frac{\tau_{b}}{\tau_{r}} \ll 1, \end{equation}

where $\tau _{b} \sim L/v$ denotes the bounce time, i.e. the typical time it takes to travel the macroscopic distance $L$ along the field. Thus, writing $f = f_{0} + f_{1} + \cdots$, one obtains at leading order

(6.3)\begin{equation} v_{{\parallel}} \frac{\partial f_{0}}{\partial l} - \frac{\mu}{m}\boldsymbol{\nabla}_{{\parallel}}B \frac{\partial f_{0}}{\partial v_{{\parallel}}} = 0. \end{equation}

One can solve this equation by making the change of variables $(l,\mu ,v_{\parallel }) \rightarrow (l,\mu ,w)$, where we have introduced the total particle energy

(6.4)\begin{equation} w = \frac{mv_{{\parallel}}^{2}}{2} + \mu B, \end{equation}

so that the kinetic equation becomes at leading order

(6.5)\begin{equation} v_{{\parallel}}\left(\frac{\partial f_{0}}{\partial l} + \mu \boldsymbol{\nabla}_{{\parallel}}B \frac{\partial f_{0}}{\partial w}\right)- \frac{\mu}{m}\boldsymbol{\nabla}_{{\parallel}}B \left(mv_{{\parallel}}\frac{\partial f_{0}}{\partial w}\right) = v_{{\parallel}}\left(\frac{\partial f_{0}}{\partial l}\right)_{w} = 0, \end{equation}

which has the solution for the lowest-order distribution function

(6.6)\begin{equation} f_{0} = f_{0}(\mu,w). \end{equation}

At the next order, the kinetic equation becomes

(6.7)\begin{equation} \frac{\partial f_{0}}{\partial t} - \frac{1}{\tau_{r}} \frac{\partial }{\partial \mu}\left(\mu f\right) - \mu\frac{B}{\tau_{r}}\frac{\partial f}{\partial w} + v_{{\parallel}}\left(\frac{\partial f_{1}}{\partial l}\right)_{w} = C_{ee}(\,f_{0}) + C_{ei}(\,f_{0}). \end{equation}

6.1. The bounce-averaged kinetic equation

In order to remove the dependence on $f_{1}$, we define the bounce average of a function $Q(\mu ,w,l)$ by

(6.8)\begin{equation} \bar{Q}(\mu,w) = \int Q(\mu,w,l) \left.\frac{{\mathrm{d}} l}{v_{{\parallel}}} \right/\int \frac{{\mathrm{d}} l}{v_{{\parallel}}}, \end{equation}

with

(6.9)\begin{equation} v_{{\parallel}} = \sqrt{\frac{2}{m}(w-\mu B)}. \end{equation}

The integral in (6.8) is taken between consecutive bounce points, defined by $v_{\parallel } = 0$, for trapped particles. For circulating particles, the integral is taken once around a field line if the field line is closed. If the field line is not closed, but instead traces out a magnetic surface, as in a stellarator or tokamak, then the bounce average for circulating particles is given by

(6.10)\begin{equation} \bar{Q}(\mu,w) =\lim_{L \rightarrow \infty} \int_{{-}L}^{L} Q(\mu,w,l) \left.\frac{{\mathrm{d}} l}{v_{{\parallel}}} \right/\int_{{-}L}^{L} \frac{{\mathrm{d}} l}{v_{{\parallel}}}. \end{equation}

Taking the bounce average of (6.7) then yields

(6.11)\begin{equation} \frac{\partial f_{0}}{\partial t} - \overline{\frac{1}{\tau_{r}}\frac{\partial }{\partial \mu}(\mu f_{0})} - \mu \overline{\left(\frac{B}{\tau_{r}}\right)}\frac{\partial f_{0}}{\partial w} = \overline{C_{ee}(\,f_{0})} + \overline{C_{ei}(\,f_{0})}. \end{equation}

Now, here we must be careful when bounce-averaging as $\tau _{r} \propto B^{-2}$. Let us write

(6.12)\begin{equation} \tau_{r} = \tau_{0}\left(\frac{B}{B_{0}}\right)^{{-}2}, \end{equation}

where $B_{0}$ and $\tau _{0}$ appearing on the right-hand side of this equation are constants. We thus arrive at the bounce-averaged kinetic equation

(6.13)\begin{equation} \frac{\partial f_{0}}{\partial t} - \frac{1}{\tau_{0}} \overline{\frac{B^{2}}{B_{0}^{2}}}\frac{\partial }{\partial \mu}\left(\mu f_{0}\right) - \frac{\mu}{\tau_{0}} \overline{\frac{B^{3}}{B_{0}^{2}}} \frac{\partial f_{0}}{\partial w} = \overline{C_{ee}(\,f_{0})} + \overline{C_{ei}(\,f_{0})}. \end{equation}

Since $\overline {B^{2}}$ and $\overline {B^{3}}$ are, in general, complicated functions of $\mu /w,$ one cannot hope to solve this equation in full generality. However, we can deduce the distribution function on a significantly long time scale $\tau$ which is larger than the radiation time and comparable to the collision time $\tau _{r} \ll \tau \lesssim \tau _{c}$. In this instance, the perpendicular energy will have been radiated away to leading order and we have

(6.14)\begin{equation} \lambda B \ll 1, \end{equation}

where we have introduced $\lambda = \mu /w$. Trapped particles will thus be absent from the population. In this regime, one can neglect the third term in (6.13) and thus obtain

(6.15)\begin{equation} \frac{\partial f_{0}}{\partial t} - \frac{1}{\tau_{0}}\overline{\frac{B^{2}}{B_{0}^{2}}}\frac{\partial }{\partial \mu}\left(\mu f_{0}\right) = \overline{C_{ee}(\,f_{0})} + \overline{C_{ei}(\,f_{0})}. \end{equation}

Moreover, we can also evaluate

(6.16)\begin{equation} \overline{\frac{B^{2}}{B_{0}^{2}}} = \frac{1}{B_{0}^{2}}\oint \left.\frac{B^{2}(l)}{\sqrt{1-\lambda B}} \,{\mathrm{d}} l \right/ \oint\frac{{\mathrm{d}} l}{\sqrt{1-\lambda B}} = \left.\frac{1}{B_{0}^{2}}\oint B^{2} \,{\mathrm{d}} l \right/\oint {\mathrm{d}} l \equiv 1, \end{equation}

where the last equivalence follows from the choice to simply define $B_{0}^{2}$ to be the average of $B^{2}$ over a field line.

Hence, we arrive at

(6.17)\begin{equation} \frac{\partial f_{0}}{\partial t} - \frac{1}{\tau_{0}}\frac{\partial }{\partial \mu}(\mu f_{0}) = \overline{C_{ee}(\,f_{0})} + \overline{C_{ei}(\,f_{0})}. \end{equation}

6.2. The bounce-averaged collision operator

We have already shown that at leading order, the contribution from electron–ion collisions is formally smaller than that from electron–electron collisions provided $n_{i}Z^{2}/n_{e} = O(1).$ Let us thus restrict our attention to the electron–electron collision operator.

We know that at leading order

(6.18)\begin{equation} C_{ee}(\,f_{0}) \simeq \frac{\sigma g(v_{{\parallel}})}{v_{{\perp}}}\frac{\partial }{\partial v_{{\perp}}}v_{{\perp}} \frac{\partial f_{0}}{\partial v_{{\perp}}}. \end{equation}

We must be cautious as $v_{\perp }$ varies along the orbit and should not be used as a coordinate in the bounce-averaged equation. When the collision operator is instead written using the magnetic moment $\mu$ as a coordinate, we obtain

(6.19)\begin{equation} \overline{C_{ee}(\,f_{0})} \simeq 2\sigma m\overline{\frac{g(v_{{\parallel}})}{B}} \frac{\partial }{\partial \mu}\left(\mu\frac{\partial f_{0}}{\partial \mu}\right) \end{equation}

and thus

(6.20)\begin{equation} \frac{\partial f_{0}}{\partial t} - \frac{1}{\tau_{0}}\frac{\partial }{\partial \mu}(\mu f_{0}) = 2 \sigma m \overline{\frac{g(v_{{\parallel}})}{B}} \frac{\partial}{\partial \mu}\left(\mu \frac{\partial f_{0}}{\partial \mu}\right). \end{equation}

It is of course no accident that this equation looks remarkably similar to the equation obtained in the straight-field-line limit. One can make the analogy exact and state that the leading-order dynamics of both systems are governed by the differential equation provided we associate the $B^{2}$ appearing in Larmor's formula, and the function $g(v_{\parallel })$ arising from the electron–electron collision operator, with their averages over a magnetic field line.

That is, on time scales $\tau$ with

(6.21)\begin{equation} \tau_{b} \ll \tau_{r} \ll \tau \lesssim \tau_{c} \end{equation}

the distribution function will become strongly anisotropic in general magnetic geometry.

7. Applicability of this theory

The theory developed in the preceding sections applies to any optically thin plasma where the collision time exceeds the radiation emission time, which is always true if the density is sufficiently small. Two broad classes of applications can be envisaged.

7.1. Cyclotron sources in the laboratory

Non-neutral plasmas, specifically plasmas consisting of charged particles with a single sign of charge can be confined in Penning–Malmberg traps (Dubin & O'Neil 1999; Danielson et al. 2015). Such a trap consists of a vacuum region inside an electrode structure consisting of a stack of hollow, metal cylinders. A uniform axial magnetic field is then applied to inhibit particle motion in the radial direction. Voltages must also be imposed on the end electrodes to prevent particle loss in the magnetic field direction.

It is frequently useful to compress plasmas radially; for instance, to increase the plasma density. This is usually accomplished by applying a torque on the plasma using rotating electric fields, the so-called ‘rotating wall technique’ (Anderegg, Hollmann & Driscoll 1998). Very long confinement times (of the order of hour or days) can be achieved using these techniques, making their use highly desirable. Particle cooling is often necessary to maintain good confinement by mitigating the heating caused by the torque using the rotating wall method. In the case of electrons or positrons, if the magnetic field is sufficiently strong, the particles will cool by cyclotron radiation (O'Neil 1985). An example of where cyclotron cooling is employed successfully to this end is in the production of anti-hydrogen, where this process is used cool pure electron plasma to sub-eV temperatures (Amoretti et al. 2002).

Another major application of the theory developed in this work will be for the first laboratory experiment, currently under development, to create and confine the first terrestrial electron–positron plasmas in the laboratory. This is done by first accumulating positrons from a powerful source and then injecting these into a pure electron plasma confined by the dipolar magnetic field of a current-carrying circular coil, so that a stationary, quasi-neutral electron–positron plasma is formed (Sunn Pedersen et al. 2012). This system will satisfy the necessary conditions of being optically thin and with a radiation time which should be initially shorter than the collision time. However, there is an additional complication in that the plasma will be so strongly magnetised that Coulomb collisions are no longer effectively described by the Landau operator. The theory of scattering in strongly magnetised anisotropic pair plasmas is developed in the companion paper (II).

7.2. Synchrotron sources in astrophysics

This work was built on the fact that a charged particle moving in a magnetic field radiates energy. At non-relativistic energies, the focus of this paper, this process is called cyclotron cooling. At relativistic velocities it is known as synchrotron radiation.

Synchrotron sources are ubiquitous and the emission of relativistic and ultra-relativistic electrons gyrating in a magnetic field is a process that dominates much of high-energy astrophysics. Indeed, it is known that synchrotron radiation is responsible for the non-thermal optical and X-ray emission observed in the Crab Nebula (Rees & Gunn 1974) and that pulsars are strong synchrotron sources (Sturrock 1971).

Of course, the physics of these systems should also be treated carefully using a model taking relativistic effects into account in the modelling of the radiation term and the collision terms. This is an area of active research and a potential axis along which this current work could be further developed. There also likely exist analogous systems in astrophysics where the electrons are subrelativistic but nevertheless are still strongly radiating and weakly collisional in the sense discussed here. Such a system would of course still be required to satisfy

(7.1)\begin{equation} \frac{\tau_{r}}{\tau_{c}} = \left(\frac{\omega_{p}}{\varOmega_{c}}\right)^{2}\left(\frac{c}{v_{\text{th}}}\right)^{3} \ll 1. \end{equation}

The easiest way to satisfy this is of course to allow relativistic plasmas, but we might also envision some exotic astrophysical plasma which is sufficiently nebulous (very low density) but also sufficiently strongly magnetised (very strong magnetic fields) as to allow this condition to be met in a non-relativistic plasma. In such systems, the theory in this paper will be directly relevant.

However, these is another caveat here which must be carefully considered. We have presented here a mechanism through which plasmas can become strongly anisotropic in velocity space. In actuality, there is a plethora of instabilities which could act on astrophysical plasmas and restore isotropy on time scales shorter than the collision time. It seems to be that the question of (in)stability could depend sensitively on plasma parameters and that considerations of isotropy-restoring instabilities do not necessarily preclude the types of systems studied in this work.

8. Conclusions

In this paper, we have developed a theory for collisional scattering in strongly anisotropic plasmas. Such plasmas arise due to the emission of radiation when charged particles move in magnetic fields, which leads to rapid depletion of the perpendicular kinetic energy. We have derived equations which describe the evolution of the electron distribution function in such plasmas. Unusually for plasma physics, the collision operator could be calculated analytically, albeit only to logarithmic accuracy, although the distribution function is far from Maxwellian.

It was found that in such strongly anisotropic populations, the evolution of the lowest-order distribution function is dominated by electron–electron collisions unless $n_{i}Z^{2}/n_{e} \gg 1$. Such collisions lead to a distribution that is Maxwellian in $v_{\perp }$ for any value of $v_{\parallel }$.

The distribution over $v_{\parallel }$ can be ascertained from an equation governing the density (in $v_{\parallel }$ space) $N,$ the integral of the distribution function over perpendicular velocities. We found that that this quantity satisfies an ill-posed, nonlinear, parabolic partial differential equation, reminiscent of the backwards diffusion equation.

This equation can be solved in two limits: firstly, in the case where the contribution from electron–ion collisions can be neglected, in which case a similarity solution is found; and secondly, in the case where the contribution from electron–electron collisions can be neglected and the equation can be solved by integration along the characteristic curves. In this latter limit, the equation is well-posed. In both cases, the remaining energy is radiated on the time scale of the ordinary collision frequency divided by $|\ln \epsilon |$.

Several candidates for areas of application of this theory were presented. These included both astrophysical applications and experimental applications. The latter of these applications will be developed in a companion paper (II).

Editor William Dorland thanks the referees for their advice in evaluating this article.

Declaration of interests

The authors report no conflict of interest.