A field–particle correlation analysis of a perpendicular magnetized collisionless shock

Abstract

Using the field–particle correlation technique, we examine the particle energization in a three-dimensional (one spatial dimension and two velocity dimensions; 1D-2V) continuum Vlasov–Maxwell simulation of a perpendicular magnetized collisionless shock. The combination of the field–particle correlation technique with the high-fidelity representation of the particle distribution function provided by a direct discretization of the Vlasov equation allows us to ascertain the details of the exchange of energy between the electromagnetic fields and the particles in phase space. We identify the velocity-space signatures of shock-drift acceleration of the ions and adiabatic heating of the electrons arising from the perpendicular collisionless shock by constructing a simplified model with the minimum ingredients necessary to produce the observed energization signatures in the self-consistent Vlasov–Maxwell simulation. We are thus able to completely characterize the energy transfer in the perpendicular collisionless shock considered here and provide predictions for the application of the field–particle correlation technique to spacecraft measurements of collisionless shocks.

1. Introduction

Shock waves, disturbances propagating faster than the largest local wave speed, are ubiquitous in space and astrophysical plasmas. From supernova remnants to the Earth's bow shock, we observe a variety of plasma environments in which these shock waves efficiently convert the bulk kinetic energy of their supersonic flows into other forms of energy, e.g. plasma heat, accelerated particles and electromagnetic radiation. Critically, a commonality amongst this variety of plasma environments is that they are weakly collisional, i.e. the collisional mean free path is much larger than the relevant plasma length scales, such as the gyroradius. Thus, this energy conversion must be mediated by collisionless interactions. Here, collisionless energy transfer refers to the myriad of mechanisms for transferring energy between the particles in the plasma and the electromagnetic fields, or vice versa, from wave–particle resonances to instabilities.

Diagnosing collisionless energy transfer is a grand challenge in plasma physics, as the plasma has many routes at its disposal for converting energy from one form to another. For collisionless shocks, a number of processes have been identified as potential energy transfer mechanisms, and the efficiency of each of these mechanisms is strongly dependent upon factors such as the shock geometry and the fast magnetosonic Mach number $M_{f} = U_{\text {shock}}/v_{f}$, where $U_{\text {shock}}$ is the shock velocity and $v_{f}$ is the fast magnetosonic wave velocity. In low-Mach-number shocks, dispersive radiation and wave–particle interactions provide an effective resistivity through the shock ramp and the requisite dissipation routes for the energy conversion of the shock (e.g. Kennel, Edmiston & Hada 1985; Balogh & Treumann 2013 and references therein).

As the shock velocity increases, particles can be reflected in the shock transition, further complicating the energy exchange (Schwartz, Thomsen & Gosling 1983; Burgess & Schwartz 1984; Scholer & Terasawa 1990; Guo & Giacalone 2013). Example energization mechanisms driven by particle reflection in high-Mach-number shocks include shock surfing acceleration (Sagdeev 1966; Sagdeev & Shapiro 1973; Lever, Quest & Shapiro 2001; Shapiro & Üçer 2003), shock-drift acceleration (Paschmann et al. 1982; Sckopke et al. 1983; Anagnostopoulos & Kaliabetsos 1994; Anagnostopoulos et al. 1998; Ball & Melrose 2001; Anagnostopoulos, Tenentes & Vassiliadis 2009; Park et al. 2013), diffusive shock acceleration (Fermi 1949, 1954; Blandford & Ostriker 1978; Ellison 1983; Blandford & Eichler 1987; Decker 1988; Malkov & Drury 2001; Caprioli, Amato & Blasi 2010), and the ‘fast Fermi’ mechanism (Leroy & Mangeney 1984; Wu 1984; Savoini, Lembége & Stienlet 2010). Additionally, the picture is not made simpler by the ways in which these different energization mechanisms may interact. The transition from particles gaining energy via shock-drift acceleration to particles gaining energy via diffusive shock acceleration (Caprioli, Pop & Spitkovsky 2014), the onset of upstream kinetic instabilities generated by the reflected particles (Schwartz et al. 1985, 1992; Schwartz 1995; Omidi, Eastwood & Sibeck 2010; Wilson et al. 2010, 2012; Turner et al. 2013; Wilson et al. 2013a, 2014a,b), the electromagnetic fluctuations resulting from these instabilities themselves contributing to the energetics of the shock via processes such as magnetic pumping (Lichko et al. 2017; Lichko & Egedal 2020) and the prospect of shock reformation due to the reflected particles, especially reflected ions (Leroy & Winske 1983; Kucharek & Scholer 1991; Giacalone et al. 1992), all complicate the energy exchange between the plasma and electromagnetic fields through the collisionless shock. Disentangling the competition between these processes remains challenging.

To ascertain the details of the energy transfer in collisionless shocks, we perform a first-principles, continuum kinetic simulation of a perpendicular collisionless shock and use the field–particle correlation technique (Klein & Howes 2016; Howes, Klein & Li 2017; Klein 2017; Klein, Howes & TenBarge 2017; Howes, McCubbin & Klein 2018; Chen, Klein & Howes 2019; Li et al. 2019; Horvath, Howes & McCubbin 2020; Klein et al. 2020) to characterize this energy exchange directly in phase space. We consider a reduced dimensionality and simplified geometry to isolate the available energization mechanisms available to the plasma, focusing on the energization mechanisms of shock-drift acceleration for the ions and adiabatic heating for the electrons. We emphasize that, while these processes have been studied previously using kinetic simulations and the particle-in-cell numerical method in higher dimensionality and greater generality (e.g. Park et al. 2013; Guo, Sironi & Narayan 2014a,b; Park, Caprioli & Spitkovsky 2015; Xu, Spitkovsky & Caprioli 2020), this study is the first direct diagnosis of the energy transfer in a collisionless shock in phase space and identification of the velocity-space signatures of shock-drift acceleration and adiabatic heating.

This Eulerian perspective—focusing on individual regions of phase space for determining the details of the energy exchange, in contrast to the more commonly used Lagrangian perspective of integrating particle trajectories to identify how individual particles are energized—is of high utility for interpreting spacecraft data. For example, using Magnetospheric Multiscale mission measurements of the electron distribution function in the Earth's turbulent magnetosheath and the field–particle correlation technique, Chen et al. (2019) found the velocity-space signature of electron Landau damping and determined that the observed turbulence was principally dissipating via electron Landau damping. In this regard, the work presented here is the beginning of a broader program of study to identify the velocity-space signatures of energization mechanisms in collisionless shocks and deploy the field–particle correlation technique for the analysis of energy exchange using in situ measurements of collisionless shocks.

To understand the observed velocity-space signatures and connect the resulting signatures to known mechanisms for plasma energization, we construct simplified analytical models for ions and electrons being energized by similar processes absent of the complications of a fully self-consistent shock, and compute the field–particle correlation on the particle distribution functions predicted by these idealized models. These simplified models allow us to proceed pedagogically and connect the two distinct pictures, the Eulerian point-of-view for identifying where in phase space the particles are being energized and the Lagrangian point-of-view for analysing how individual particles gain and lose energy. While significant intuition is gained from the Lagrangian perspective, this novel Eulerian perspective provided by the field–particle correlation technique has some advantages, chief among them is the ability to easily distinguish how different regions of phase space are being energized. We will show how the field–particle correlation technique allows us to easily separate the energy exchange occurring between the electromagnetic fields and the multi-component distribution functions (e.g. incoming beam versus reflected ions in the shock foot and ramp) which frequently characterize collisionless shock dynamics. We may thus distinguish the different effects the same electromagnetic fields are having on different parts of the distribution function, from how the cross-shock electric field decelerates the incoming bulk flow and accelerates reflected ions to how the motional electric field supporting the upstream $\boldsymbol {E} \times \boldsymbol {B}$ motion energizes both the reflected ions and bulk electrons.

The rest of the paper is organized as follows. In § 2 we provide details of the simulations performed, followed by a broad overview of the results of the simulations examining the shock structure in the electromagnetic fields and electron and ion distribution functions. We then present an overview of the central analysis tool of this paper, the field–particle correlation technique, in § 3. We apply the field–particle correlation technique to obtain the key results of the paper in §§ 4 and 5: the velocity-space signatures of (i) shock-drift acceleration of the ions and (ii) adiabatic heating of the electrons. We conclude in § 6 with a discussion of the implications of the results presented for spacecraft observations and future avenues of research applying field–particle correlations to a larger range of shock parameters.

2. Computational model and overview of results

To perform a self-consistent simulation of a perpendicular collisionless shock, we employ the continuum Vlasov–Maxwell solver in the Gkeyll framework (Juno et al. 2018; Hakim & Juno 2020). We emphasize that, unlike traditional particle based approaches such as the particle-in-cell method, Gkeyll directly discretizes the Vlasov–Maxwell system of equations on a phase-space grid to obtain a high-fidelity representation of the distribution function, free of the shot noise introduced by finite-sized particles. In other words, we solve the following system of equations with a grid-based method for every species $s$ in the plasma,

(2.1)\begin{equation} \left.\begin{array}{c} \displaystyle \dfrac{\partial f_s}{\partial t} + \boldsymbol{\nabla} \boldsymbol{\cdot} (\boldsymbol{v} f_s) + \boldsymbol{\nabla} _{\boldsymbol{v}} \boldsymbol{\cdot} \left( \dfrac{q_s}{m_s}[\boldsymbol{E} + \boldsymbol{v} \times \boldsymbol{B}] f_s \right) = 0, \\ \displaystyle \dfrac{\partial \boldsymbol{B}}{\partial t} + \boldsymbol{\nabla} \times\boldsymbol{E} = 0, \quad \epsilon_0\mu_0\dfrac{\partial \boldsymbol{E}}{\partial t} - \boldsymbol{\nabla} \times\boldsymbol{B} ={-}\mu_0\boldsymbol{J}, \\ \displaystyle \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{E} = \dfrac{\varrho_c}{\epsilon_0}, \quad \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{B} = 0, \end{array}\right\} \end{equation}

where $f_s = f_s(\boldsymbol {x}, \boldsymbol {v})$ is the particle distribution function for species $s$, $q_s$ and $m_s$ are the charge and mass of species $s$, respectively, $\boldsymbol {E} = \boldsymbol {E}(\boldsymbol {x})$ and $\boldsymbol {B} = \boldsymbol {B}(\boldsymbol {x})$ are the electric and magnetic fields,respectively, and the coupling between the electromagnetic fields and particles is given by velocity moments of the particle distribution function,

(2.2a,b)\begin{equation} \varrho_c = \sum_s q_s \int f_s \, \textrm{d} \boldsymbol{v}, \quad \boldsymbol{J} = \sum_s q_s \int \boldsymbol{v} f_s \, \textrm{d} \boldsymbol{v}, \end{equation}

i.e. the charge and current density. This approach has been previously leveraged in the study of electrostatic collisionless shocks (Pusztai et al. 2018; Sundström et al. 2019), allowing for a detailed study of the phase-space dynamics which results from the evolution of the shock. Comparisons to the particle-in-cell method for the study of kinetic instabilities have clearly demonstrated the advantages of using a continuum representation to eliminate discrete particle noise in the particle velocity distributions (Skoutnev et al. 2019; Juno et al. 2020).

Here, a perpendicular shock refers to the orientation of the magnetic field with respect to the shock normal. Because the magnetic field is perpendicular to the shock normal, in one spatial dimension, we require only the two velocity components perpendicular to the magnetic field to fully describe the dynamics of the system, i.e. 1D-2V. The particular one-dimensional (1-D) geometry we choose is the one spatial coordinate along the shock normal in the $x$ direction, with the initial magnetic field in the $z$ direction, $\boldsymbol {B} (t = 0) = B_0 \hat {\boldsymbol {z}}$. For completeness, in this dimensionality and field geometry, the Vlasov–Maxwell system of equations is

(2.3)\begin{equation} \left. \begin{array}{c} \displaystyle \dfrac{\partial f_s}{\partial t} + v_x \dfrac{\partial f_s}{\partial x} + \dfrac{q_s}{m_s} \left( [E_x + v_y B_z] \dfrac{\partial f_s}{\partial v_x} + [E_y - v_x B_z] \dfrac{\partial f_s}{\partial v_y} \right ) = C[f_s], \\ \displaystyle \dfrac{\partial B_z}{\partial t} ={-}\dfrac{\partial E_y}{\partial x}, \\ \displaystyle \dfrac{\partial E_x}{\partial t} ={-}\dfrac{J_x}{\epsilon_0}, \quad \dfrac{\partial E_y}{\partial t} ={-}c^2 \displaystyle \dfrac{\partial B_z}{\partial x} - \dfrac{J_y}{\epsilon_0}, \end{array}\right\} \end{equation}

where we have added a collision operator $C[f_s]$ on the right-hand side of the Vlasov equation.¹

The electrons and ions are initialized with the same supersonic flow directed in the negative $x$ direction towards a reflecting wall, which leads to a shock wave that propagates in the positive $x$ direction in our simulation. Note that the particles reflect from the wall, but the ‘reflecting wall’ boundary condition for the electromagnetic fields is a conducting wall boundary condition in the traditional sense, with zero normal magnetic field and zero tangential electric field. This method of initialization is often called the ‘reflecting-wall’ set-up² and has been previously employed in numerous particle-in-cell studies of collisionless shocks (e.g. Papadopoulos, Wagner & Haber 1971; Spitkovsky 2005, 2007).

Detailed parameters are as follows: the reflecting wall for the particles and conducting wall for the electromagnetic fields are at $x = 0$, and plasma is injected with a copy boundary condition³ at $x = 25 d_i$, where $d_i = c/\omega _{pi}$ is the ion collisionless skin depth. Here, $c$ is the speed of light and $\omega _{pi} = \sqrt {e^2 n_0/\epsilon _0 m_i}$ is the ion plasma frequency. Note that the subscript $0$ denotes the upstream value, e.g. $n_0$ is the upstream density and $B_0$ is the upstream magnetic field magnitude. We use a reduced mass ratio between the ions and electrons, $m_i/m_e = 100$. The total plasma beta, $\beta = 2 \mu _0 n_0 (T_{e_0} + T_{i_0})/B_0^2 = 2$, with the ion beta, $\beta _i = 1.3$, and the electron beta, $\beta _e = 0.7$. Both the ions and electrons are non-relativistic, with $v_{te}/c = 1/16$, where $v_{ts} = \sqrt {2 T_{s_0}/m_s}$. With this choice of electron beta and $v_{te}/c$, the ratio of the electron plasma frequency, $\omega _{pe} = \sqrt {e^2 n_0/\epsilon _0 m_e}$, to the electron cyclotron frequency, $\varOmega _{ce} = -e B_0/m_e$, is $\omega _{pe}/\varOmega _{ce} \sim 13.4$. The in-flow velocity in the simulation frame to initialize the perpendicular, electromagnetic shock is $U_x = -3 v_A$, where $v_A = B_0/\sqrt {\mu _0 n_0 m_i}$ is the ion Alfvén speed. Note that the in-flow velocity is negative because the plasma initially flows in the negative $x$ direction. Because the plasma is initialized with a flow in a background magnetic field, we initialize the corresponding electric field to support this flow, $\boldsymbol {E} = -U_x \hat {\boldsymbol {x}} \times \boldsymbol {B} = U_x B_0 \hat {\boldsymbol {y}}$.

For the grid resolution in configuration space, we use $N_x = 1536$ grid cells, which correspond to ${\rm \Delta} x \sim d_e/6 \sim 3.7 \lambda _D$, where $d_e = c/\omega _{pe}$ and $\lambda _D = v_{te}/(\sqrt {2} \omega _{pe})$ are the electron inertial length and electron Debye length, respectively, and we employ piecewise quadratic Serendipity elements for the discontinuous Galerkin basis expansion (Arnold & Awanou 2011). In velocity space, the electron extents are $\pm 8 v_{te}$, and the ion extents are $\pm 16 v_{ti}$, with zero-flux boundary conditions at the velocity-space limits and $N_{v_x} = N_{v_y} = 64$ grid cells for both species, which correspond to ${\rm \Delta} v = v_{te}/4$ for the electrons and ${\rm \Delta} v = v_{ti}/2$ for the ions. The basis expansion in velocity space is also piecewise quadratic Serendipity elements. For further details about the algorithm and the choice of basis expansion, we refer the reader to Juno et al. (2018) and Hakim & Juno (2020).

We have run the simulation with a small amount of collisions to regularize velocity space. In this case, we choose an electron–electron collision frequency, $\nu _{ee} = 0.01 \varOmega _{ci}$, much less than the ion cyclotron frequency, $\varOmega _{ci} = e B_0/m_i$, with the ion–ion collision frequency correspondingly smaller based on the square root of the mass ratio, $\nu _{ii} = 0.001 \varOmega _{ci}$. Note that because the ions are hotter than the electrons, they should formally be even more collisionless than the electrons; however, these collisionalities are larger than typical solar wind collisionalities (Wilson et al. 2018) and are not chosen to be realistically small, but instead chosen to be just large enough to provide regularization of the velocity-space structure given finite velocity resolution. Details on the implementation of the collision operator and its conservation properties can be found in the paper by Hakim et al. (2019).

In figure 1, we show the electromagnetic fields and particle distribution functions for the electrons and ions in $x-v_x$ phase space after the perpendicular shock has formed and propagated through the simulation domain, $t_{end} = 11 \varOmega _{ci}^{-1}$. Although the downstream of the shock is fairly oscillatory as the energy injected into the plasma by the shock sloshes back and forth between the electromagnetic fields and particles, we can estimate the compression ratio of this low-Mach-number shock from the average downstream magnetic field to be roughly $r \sim 2.5$ (solid black line in figure 1(c)). With this estimate for the compression ratio, we calculate the shock velocity in the simulation frame to be $U_{\text {shock}} = U_x/(r-1) = 2 v_A$. Note that in this reflecting wall set-up, the simulation frame is equivalent to the frame in which the downstream plasma is at rest.

Figure 1. The $x$-electric field (a), $y$-electric field (b), $z$-magnetic field (c), ion distribution function integrated in $v'_y$ (d) and electron distribution function integrated in $v'_y$ (e) after the perpendicular shock has formed and propagated through the simulation domain, $t = 11 \varOmega _{cp}^{-1}$. Note that the distribution functions are plotted in the simulation frame $f_s(x, v_x', v_y')$ for each species $s$. We have marked an approximate transition from upstream of the shock to the shocked plasma (dashed–dotted line), and likewise an approximate transition from the shock to the downstream region (dashed line). To mark the oscillation of the electromagnetic fields and the sloshing of energy between the fields and particles in the downstream region, we have used a solid black line to mark the approximate compression of the magnetic field, along with $\boldsymbol {E} = 0$. We expect the $y$-electric field to roughly oscillate about zero in the frame of the simulation, as the ‘reflecting-wall’ set-up is performed in the frame of the downstream plasma, where the $\boldsymbol {E} \times \boldsymbol {B}$ velocity is zero.

Thus, combining the velocity of the incoming flow with the velocity of the shock in the downstream frame, this self-consistently produced perpendicular shock is a $M_A = 5$ shock, where $M_A$ is the Alfvén Mach number. Equivalently, using the definitions for the sound speed and magnetosonic speed given by

(2.4)\begin{gather} c_s = \sqrt{\frac{\gamma_i T_i + \gamma_e T_e}{m_i}}, \end{gather}

(2.5)\begin{gather}v_{f} = \sqrt{c_s^2 + v_A^2}, \end{gather}

where $\gamma _i = \gamma _e = 1 + 2/\textrm {VDIM} = 2$ because the simulation domain has two velocity dimensions, we find this shock has a fast magnetosonic Mach number, $M_{f} \approx 2.89$. With these plasma parameters and this magnetosonic Mach number, we note that this shock is supercritical $M_f> M_{f_{crit}} \simeq 2$ (Wilson 2016), similar to the Earth's bow shock, and thus bodes well for the ultimate goal of predicting velocity-space signatures of energization mechanisms in spacecraft observations of heliospheric shocks.

In this regard, we focus our attention now on the particle distribution functions and the phase-space structure generated through the shock. While the particle distribution functions in the $x-v_x$ phase space shown in figure 1 are illustrative of the dynamics through the shock, which shows a reflected population of ions in figure 1(d) and a clear compression of the electrons in figure 1(e), we can gain further insight into the dynamics of this shock by looking at the distribution function in $v_x-v_y$ at fixed points in configuration space through the shock. In figure 2, we plot the ion (a–f) and electron (g–i) distribution functions in velocity space at several points through the shock, from upstream through the ramp to downstream.

Figure 2. The ion (panels (a)–(f)) and electron (panels (g)–(i)) distribution functions in the simulation frame (downstream frame) $f_s(x, v_x', v_y')$, for each species $s$, plotted at various points through the shock at $t = 11 \varOmega _{cp}^{-1}$. As we move from upstream, $x = 24.5 d_i$, through the shock ramp, $x = 21.5 d_i$, we can identify the reflected ion population as well as a broadening of the electron distribution function.

As an example of the wealth of data contained in the distribution function, we draw special attention to the ion distribution function in the shock ramp. As we move from upstream into the shock, at the beginning of the ramp at $x=22.5d_i$, we begin to see a small population of reflected ions, forming a small ‘crescent’ distribution in the lower right quadrant of the $v_x-v_y$ space. Further up the ramp at $x=21.5d_i$, we observe that the incoming ion beam begins to be deflected by the fields in the shock transition, which generates a ‘boomerang’ distribution that smoothly connects the decelerated incoming ion beam with the reflected ion population. It is this reflected population, in agreement with previous studies of supercritical shocks (e.g. Ball & Melrose 2001; Balogh & Treumann 2013), which dominates the ion energization and provides a segue to our key result: diagnosing the velocity-space signatures of particle energization in this perpendicular electromagnetic shock. To obtain these velocity-space signatures, we now describe our tool of choice for our analysis of the high-quality distribution function data provided by the continuum kinetic simulation: the field–particle correlation technique.

3. The field–particle correlation technique

From combining the Vlasov equation and Maxwell's equations, we can obtain a conservation equation for the total energy of the kinetic plasma (e.g. Klein et al. 2017),

(3.1)\begin{equation} W = \int \textrm{d} \boldsymbol{x} \left(\frac{\epsilon_0}{2} |\boldsymbol{E}|^2 + \frac{1}{2\mu_0}|\boldsymbol{B}|^2\right) + \sum_s \int \textrm{d} \boldsymbol{x} \int \textrm{d} \boldsymbol{v} \frac{1}{2}m_s v^2 f_s. \end{equation}

The first integral represents the energy of electromagnetic fields in the plasma and the second accounts for the combined microscopic kinetic energy⁴ of all plasma species $s$. In the absence of particle collisions, the net microscopic kinetic energy of a given plasma species may only be changed through collisionless interactions between the particles of that species and the electromagnetic fields.

To explore the energy transfer between fields and particles, we define the phase-space energy density for a particle species $s$ by $w_s(\boldsymbol {x},\boldsymbol {v},t) \equiv m_s v^2 f_s(\boldsymbol {x},\boldsymbol {v},t)/2$ in the non-relativistic limit. Multiplying the Vlasov equation by $m_s v^2/2$, we obtain an expression for the rate of change of this phase-space energy density,

(3.2)\begin{equation} \frac{\partial w_s(\boldsymbol{x},\boldsymbol{v},t)}{\partial t} ={-} \boldsymbol{v}\boldsymbol{\cdot} \boldsymbol{\nabla} w_s - q_s\frac{v^2}{2} \boldsymbol{E} \boldsymbol{\cdot} \frac{\partial f_s}{\partial \boldsymbol{v}} - q_s\frac{v^2}{2} (\boldsymbol{v} \times \boldsymbol{B}) \boldsymbol{\cdot} \frac{\partial f_s}{\partial \boldsymbol{v}}. \end{equation}

This equation describes the mechanisms that govern how the energy density in the 3D-3V phase space $(\boldsymbol {x},\boldsymbol {v})$ evolves, where each term has a clear physical interpretation.

The first term on the right-hand side of (3.2) describes how $w_s(\boldsymbol {x},\boldsymbol {v},t)$ changes due to particle advection from other spatial regions, which gives rise in fluid theory to the energy change through pressure forces and heat fluxes.⁵ Because this term describes the advection of particle kinetic energy as particles move from one spatial position to another, when integrated over the full plasma volume, this term yields zero net change of the total kinetic energy of particle species $s$, $\mathcal {W}_s= \int \textrm {d} \boldsymbol {x} \int \textrm {d} \boldsymbol {v} \frac {1}{2}m_s v^2 f_s$. The third term on the right-hand side of (3.2) describes the magnetic forces on the particles. Although this term can move kinetic energy from one location in velocity space to another, when integrated over all velocity space, this term does zero net work on the particles, as expected for the magnetic force.

The second term on the right-hand side of (3.2) describes the work done on the plasma species $s$ by the electric field. When (3.2) is integrated over all velocity space and all physical space to obtain the rate of change of the total kinetic energy $\mathcal {W}_s$ of a particle species $s$, the first and third terms have zero net contribution (Klein & Howes 2016; Howes et al. 2017), which yields

(3.3)\begin{equation} \frac{\partial \mathcal{W}_s}{\partial t} ={-} \int \textrm{d} \boldsymbol{x} \int \textrm{d} \boldsymbol{v} q_s \frac{v^2}{2} \frac{\partial \,f_{s}}{\partial \boldsymbol{v}} \boldsymbol{\cdot} \boldsymbol{E} = \int \textrm{d} \boldsymbol{x}\left( \int \textrm{d} \boldsymbol{v} q_s \boldsymbol{v} \,f_{s} \right) \boldsymbol{\cdot} \boldsymbol{E}= \int \textrm{d} \boldsymbol{x} \boldsymbol{j}_s \boldsymbol{\cdot} \boldsymbol{E}, \end{equation}

This expression makes clear that the change in species energy $\mathcal {W}_s$ arises from work done on that species by the electric field, $\boldsymbol {j}_s \boldsymbol {\cdot } \boldsymbol {E}$.

In our exploration of particle energization at collisionless shocks, we choose to focus on the second term in (3.2) to investigate the energization of the particles by the electric field. The form of that term demonstrates that the rate of particle energization can be computed at a single-point in physical space $\boldsymbol {x}_0$ by measuring the electric field at that position $\boldsymbol {E}(\boldsymbol {x}_0)$ and the particle velocity distribution at the same position $f_s(\boldsymbol {x}_0,\boldsymbol {v})$. This fundamental fact underlies the field–particle correlation (FPC) technique (Klein & Howes 2016; Howes et al. 2017; Klein et al. 2017), where the unnormalized correlation (essentially a time average) of the product of the electric field $\boldsymbol {E}(\boldsymbol {x}_0)$ and a term that depends on the particle velocity distribution $f_s(\boldsymbol {x}_0,\boldsymbol {v})$ over some correlation interval $\tau$ is computed by

(3.4)\begin{equation} C_{\boldsymbol{E}}(\boldsymbol{x}_0,\boldsymbol{v},t,\tau) = \frac{1}{\tau} \int_{t-\tau/2}^{t+\tau/2} \left[{-}q_s\frac{v^2}{2} \frac{\partial f_s(\boldsymbol{x}_0,\boldsymbol{v},t')}{\partial \boldsymbol{v}} \right] \boldsymbol{\cdot} \boldsymbol{E}(\boldsymbol{x}_0,t')\, \textrm{d} t'. \end{equation}

The resulting correlation $C_{\boldsymbol {E}}(\boldsymbol {x}_0,\boldsymbol {v},t,\tau )$ directly measures the rate of change of phase-space energy density at position $\boldsymbol {x}_0$ as a function of 3V particle velocity space $\boldsymbol {v}$, which produces a velocity-space signature that is characteristic of the mechanism of energization and can be used to identify a particular, locally-occurring energization process, e.g. Landau damping (Howes et al. 2017; Klein et al. 2017) and cyclotron damping (Klein et al. 2020). We note that as part of this identification, further analysis may be required to ascertain certain details; for example, if one obtains velocity-space signatures corresponding to the presence of Landau damping, the resonant velocity that the velocity-space signature is concentrated around is necessary to determine what wave modes are Landau damping in the plasma, as one can find similar structure whether a Langmuir wave (Howes et al. 2017) or kinetic Alfvén wave (Klein et al. 2017; Horvath et al. 2020) is undergoing Landau damping. However, even with this caveat, a key advantage of the FPC method to diagnose particle energization is that it requires only measurements at a single spatial point $\boldsymbol {x}_0$ to determine the energization by the electric field. An appropriately instrumented single spacecraft mission can provide the requisite full 3V particle velocity distribution $f_s(\boldsymbol {x}_0,\boldsymbol {v},t)$ and electric field $\boldsymbol {E}(\boldsymbol {x}_0,t)$ at the spacecraft position $\boldsymbol {x}_0$ as a function of time. Thus, the velocity-space signatures determined here using kinetic numerical simulations, our key results, may be directly sought using spacecraft observations.

In the case of particle energization as a consequence of the dissipation of weakly collisional plasma turbulence, the rate of particle energization represented by the second term in (3.2) generally includes two distinct contributions: (i) an often large-amplitude oscillatory component that leads to zero net energization, which is associated with undamped wave motion, and (ii) a typically smaller amplitude secular component that corresponds to the net collisionless transfer of energy from the fields to the particles (Klein & Howes 2016; Howes et al. 2017). By an appropriate choice of the correlation interval $\tau$, the oscillatory energy transfer is largely eliminated by time averaging, which exposes the secular energy transfer associated with the collisionless damping of the turbulent fluctuations. For the perpendicular collisionless shock in this study, the shock is quasi-stationary in the shock–rest frame of reference, with smooth electromagnetic fields through the shock as seen in figure 1, and thus we need not time average the correlation, but instead take the instantaneous correlation (the limit $\tau \rightarrow 0$). We will thus suppress the dependence of the correlation on $\tau$ henceforth. We note that the FPC with $\tau =0$ is simply the instantaneous rate of change of the phase-space energy density, $\partial w_s/\partial t$, due to work done on the particles by the electric field. If kinetic instabilities were to arise upstream or within the shock transition region, or if the shock itself were to become non-stationary, then it is likely that taking a correlation interval $\tau$ longer than either the unstable wave period or the shock reformation time would be necessary to recover a meaningful velocity-space signature of the net particle energization.

In addition, we adopt two final modifications of the FPC analysis that are well suited for the study of collisionless shocks: (i) we separate the contributions to the rate of energization by the different components of the electric field, $E_x$ and $E_y$; and (ii) we replace $v^2$ in (3.4) by the component associated with the electric field, e.g. using $v_x^2$ for the correlation using $E_x$. We refer the reader to Appendix A for a discussion of the validity and usefulness of this transformation. Therefore, the form of the FPCs implemented here for a position $x=x_0$ in our 1D-2V Gkeyll simulation is given by

(3.5)\begin{gather} C_{E_x} (x_0, v_x, v_y, t) ={-}q_s \frac{v_x^2}{2} E_x(x_0,t) \frac{\partial f_s(x_0, v_x, v_y, t)}{\partial v_x}, \end{gather}

(3.6)\begin{gather}C_{E_y} (x_0, v_x, v_y, t) ={-}q_s \frac{v_y^2}{2} E_y(x_0,t) \frac{\partial f_s(x_0, v_x, v_y, t)}{\partial v_y}. \end{gather}

An issue which cannot be overemphasized in performing the FPC analysis of a collisionless shock is making a judicious choice of the frame of reference in which to calculate (3.5) and (3.6) (Goodrich & Scudder 1984). We choose to evaluate the correlations in the frame of reference in which the shock is at rest (the shock–rest frame, unprimed variables), as opposed to the frame of reference of the simulation, in which the plasma is at rest downstream of the shock (the simulation frame or downstream frame, primed variables).⁶ For clarity, the shock velocity in the simulation frame is given by ${\boldsymbol {U}}_{\textrm {shock}} = U_{\textrm {shock}} {\hat {\boldsymbol {x}}} = 2 v_A {\hat {\boldsymbol {x}}}$. It is critical not only that the velocity coordinates are transformed to the shock–rest frame, $\boldsymbol {v}=\boldsymbol {v}'-\boldsymbol {U}_{\text {shock}}$, but also that the electromagnetic fields are appropriately Galilean transformed to the shock–rest frame,

(3.7)\begin{equation} \boldsymbol{E} = \boldsymbol{E}' + \boldsymbol{U}_{\text{shock}} \times \boldsymbol{B}', \end{equation}

and $\boldsymbol {B}=\boldsymbol {B}'$.

We note that our discussion and application of the FPC to the collisionless shock is principally concerned with how the plasma is energized via the electromagnetic fields and focuses on the phase-space dynamics governed by the electric field term in the Vlasov equation. As mentioned previously, plasmas additionally convert bulk kinetic to thermal energy, and vice versa, via other terms in (3.2) such as the $\boldsymbol {v} \boldsymbol {\cdot } \boldsymbol {\nabla }$ term, which gives rise to pressure forces and heat fluxes. To explore how these other physical mechanisms impact the flow of energy through 3D-3V phase space, one can perform complementary correlations with these other terms. Correlating with the magnetic term in the Lorentz force allows the determination of how the magnetic field leads to changes in $w_s(\boldsymbol {x},\boldsymbol {v},t)$ as a function of velocity $\boldsymbol {v}$—e.g. energy can be moved between different degrees of freedom by the magnetic field, even though the net energy change (integrated over velocity space) must always be zero. Similarly, if spatial gradients of $f_s(\boldsymbol {x},\boldsymbol {v},t)$ are available, the velocity-space signatures of the work done on the particles by the pressure tensor can be determined. Of course, computing the total rate of change of the phase-space energy density $w_s(\boldsymbol {x},\boldsymbol {v},t)$ at a particular point in configuration and velocity space requires all the terms of (3.2).

Our focus here, however, is on the term in the Vlasov equation which produces net energization of a plasma species $s$, the electric field term in (3.2). In fact, as shown in Appendix B, these additional terms, such as the $\boldsymbol {v} \times \boldsymbol {B}$ term, can have a cancellation effect on the evolution of the phase-space energy density, so that the net energization due to, for example, an $\boldsymbol {E} \times \boldsymbol {B}$ drift is identically zero, as it should be. As such, we are well justified in formulating the FPC to focus only on the net energization and avoid obfuscating the signatures of energization with the additional motion of phase-space energy density arising from these other terms in the Vlasov equation. For a further discussion of energization versus energy conversion, we refer the reader to Appendix C.

4. Field–particle correlation analysis: ions

4.1. Velocity-space signature of ion energization

In figure 3(a), the ion distribution function includes both a component from the incoming beam of ions upstream, as well as the aforementioned ‘crescent’ population of reflected ions. For the $E_x$ contribution to the FPC, $C_{E_x}(v_x,v_y)$, at position $x = 22.9 d_i$, figure 3(b) shows that the incoming ion beam is being acted upon strongly by the cross-shock electric field $E_x$, but that $E_x$ has little effect on the reflected ion population at this position. On the other hand, for the $E_y$ contribution to the FPC, $C_{E_y}(v_x,v_y)$, at position $x = 22.9 d_i$, figure 3c) shows that the reflected ions principally interact with this component of the electric field, i.e. the motional electric field which supports the incoming $\boldsymbol {E} \times \boldsymbol {B}$ flow.

Figure 3. The ion distribution function $f_i(v_x,v_y)$ (a), and the $C_{E_x}$ (b) and $C_{E_y}$ (c) components of the FPC, (3.5) and (3.6), computed at $x = 22.9 d_i$ from the self-consistent Gkeyll simulation. Note that the FPC is computed in the shock–rest frame. While the bulk incoming ions are slowed down by the cross-shock electric field, $E_x$, we see the distribution of reflected ions gain energy due to the motional electric field, $E_y$, which supports the incoming supersonic $\boldsymbol {E} \times \boldsymbol {B}$ flow.

To understand this visual representation of the rate of ion energization over velocity space, recall that the FPC determines the rate of change of the phase-space energy density of a particular plasma species, $w_s(\boldsymbol {x},\boldsymbol {v},t) = m_s |\boldsymbol {v}|^2 f_s(\boldsymbol {x},\boldsymbol {v},t)/2$, due to the electric field. The phase-space energy density of the ions, $w_i$, can only change if the number of ions in that volume of phase space changes. Therefore, the nested blue and red crescents in figure 3(c) indicate that ions are accelerated by $E_y$ from the blue region to the red region. Conservation of particle number requires that the number of ions lost from the blue region is the same as the number gained in the red region, but because the red region is at higher velocity $v_y$, the net effect, obtained by integrating $C_{E_y}$ over velocity space $(v_x,v_y)$, is an increase in the ion phase-space energy density $w_i$. We also note that the observed $C_{E_y}$ signature is a larger amplitude than the observed $C_{E_x}$, such that $E_y$ dominates the energy exchange at this particular point in space. Furthermore, the FPC method computes the rate of change of energy density, so the rate of energization per ion in the low-density population of reflected ions is much higher in amplitude than the loss of energy per ion by the much more dense incoming beam.

As a first attempt to understand this signature, consider that the gradient length scale of the collisionless shock in our simulation is $L_{\text {shock}} \sim \rho _i$, where $\rho _i = v_{t_i}/\varOmega _{ci}$ is the ion Larmor, or gyro-, radius. Therefore, ions encountering this gradient in the magnetic field will not necessarily have closed orbits and smoothly transition downstream. Depending on an ion's gyrophase when it encounters this magnetic field gradient, the ion's new Larmor orbit may cause the ion to move back upstream, where the magnetic field magnitude is smaller. The increased Larmor radius of this reflected ion in the upstream region then allows the ion to gain energy along the motional electric field supporting the incoming $\boldsymbol {E} \times \boldsymbol {B}$ motion. This energization of the reflected ion population via $E_y$ is consistent with the well-known energization mechanism, shock-drift acceleration (Paschmann et al. 1982; Sckopke et al. 1983; Anagnostopoulos & Kaliabetsos 1994; Anagnostopoulos et al. 1998; Ball & Melrose 2001; Anagnostopoulos et al. 2009; Park et al. 2013). However, to understand why shock-drift acceleration would produce the particular velocity-space signature observed in figure 3(c), we turn to a simplified analytic model to connect the well-known Lagrangian picture for shock-drift acceleration with the new Eulerian perspective granted by the FPC.

4.2. Shock-drift acceleration in an idealized perpendicular shock

We consider now a simplified reduction of the electromagnetic fields observed in our self-consistent simulation to a step function in the magnetic field,

(4.1)\begin{equation} B_z(x) = \left\{ \begin{array}{@{}cc} B_{u} & x \ge 0\\ B_{d} & x < 0 \end{array} \right., \end{equation}

with amplitude jump $B_{d}/B_{u}=4$. We will also continue to work exclusively in the shock–rest frame, where to a good approximation the motional electric field, $E_y$, is a constant through the entire shock. The value of the constant $E_y$, as well as the ion and electron plasma betas, are chosen so that the shock velocity is similar to the self-consistent simulation, $M_A = 4.9$ and $M_f = 3.0$. This reduced model corresponds to the limit $L_{\text {shock}}/\rho _i \ll 1$ and allows us to decompose the ion motion more easily between upstream and downstream gyro- and $\boldsymbol {E} \times \boldsymbol {B}$ motion. To mimic the geometry of the self-consistent simulation, we take $E_y < 0$ and $B_z > 0$ so that the inflow $\boldsymbol {E} \times \boldsymbol {B}$ is in the negative $x$ direction.

In figure 4, we plot (a) the trajectory of an ion in the $(x,y)$ plane and (b) its corresponding trajectory in $(v_x,v_y)$ velocity space in the shock–rest frame, where the colours indicate the corresponding segments of the trajectory. In the upstream region at $x>0$ (black), the black circle centred about the upstream $\boldsymbol {E} \times \boldsymbol {B}$ velocity (black star) corresponds to the Larmor orbit of the ion about the upstream inflow velocity in the $(v_x,v_y)$ plane. Upon first crossing the magnetic discontinuity to $x<0$, the ion changes to a Larmor gyration in the $(v_x,v_y)$ plane (blue) about the downstream $\boldsymbol {E} \times \boldsymbol {B}$ velocity (green star). In the larger amplitude downstream perpendicular magnetic field, the radius of the Larmor motion in the $(x,y)$ plane in the shock–rest frame is reduced (blue).

Figure 4. (a) Real space trajectory of an ion as it traverses the shock front and (b) velocity-space trajectory.

Depending on the ion's gyrophase when the ion crosses the magnetic discontinuity, the ion passes back upstream to $x>0$, and once again undergoes a Larmor orbit in the $(v_x,v_y)$ plane (red) about the upstream $\boldsymbol {E} \times \boldsymbol {B}$ velocity (black star). In this segment of the trajectory (red), the ion gains perpendicular energy in the shock–rest frame, graphically represented by the distance in velocity space of the ion from the origin of the $(v_x,v_y)$ plane. Finally, the ion will eventually cross back into the downstream region to $x<0$ (green), and resumes its Larmor orbit in the $(v_x,v_y)$ plane (green) about the downstream $\boldsymbol {E} \times \boldsymbol {B}$ velocity (green star). Without any additional crossings of the magnetic discontinuity, the ion will simply $\boldsymbol {E} \times \boldsymbol {B}$ drift downstream.

In the segment of the trajectory where the ion can gain energy, it is the motional electric field, $E_y$, that is doing positive work on the ion, exactly like in our self-consistent simulation. We note that this ion's dynamics—the reflection due to the magnetic gradient and energy gain from its traversal upstream and alignment with the motional electric field—is the well-known single-particle picture of shock-drift acceleration. In fact, this picture in velocity space, where a single ion gains energy via this reflection by a magnetic gradient, has been previously noted (Gedalin 1996a). We wish now to connect this Lagrangian perspective on how a single ion gains energy from this reflection off a magnetic gradient to the Eulerian point-of-view we have from the FPC.

4.3. Velocity-space signature of shock-drift acceleration

To connect the single-particle picture of shock-drift acceleration with how a distribution of ions is energized, we employ the Vlasov-mapping technique (Scudder et al. 1986; Kletzing 1994; Hull et al. 1998; Hull & Scudder 2000; Hull et al. 2001; Mitchell & Schwartz 2013, 2014), described in Appendix D, to determine the velocity distribution function in our simplified model for the electromagnetic fields through the shock. We show in figure 5 the reconstructed ion distribution $f_i(v_x,v_y)$ (a) at $x = 0.4 d_i$ and the corresponding FPC $C_{E_y}(v_x,v_y)$ (b) computed from the motional electric field, $E_y$, and gradients of this reconstructed distribution function. In addition, we repeat figures 3(a) and (c), for reference in comparing the distribution function and generated velocity-space signature between the simplified model and self-consistent simulation.

Figure 5. Comparison of the reconstructed ion distribution function (a) and $C_{E_y}$ component of the FPC (b) computed from this reconstruction to the self-consistently produced ion distribution function (c) and $C_{E_y}$ component of the FPC (d) from the Gkeyll simulation. Using the Vlasov-mapping technique, we can connect the single-particle orbits (overplotted white (a) and black (b) lines) to the distribution function dynamics. Here, $C_{E_y}$ integrated over velocity space is net positive, which means the observed velocity-space signature corresponds to an energization process. We identify this particular velocity-space signature as the signature of shock-drift acceleration, energization of the reflected ions via the motional electric field in the upstream, via the connection between where in velocity space a single ion is energized and the specific region of velocity space where the strongest energy exchange is occurring.

In the reconstructed distribution function from the idealized model, we identify, in addition to the incoming upstream population centred at the upstream $\boldsymbol {E} \times \boldsymbol {B}$ velocity, a component of reflected particles that have returned upstream, exactly like in the self-consistent simulation. Overplotted on the ion distribution function and computed FPC from the Vlasov-mapping technique is the trajectory in $(v_x,v_y)$ for the ion analysed in figure 4, which shows that this reflected population and velocity-space signature are coincident with the red segment of the trajectory in figure 4. Integrating this field–particle correlation over velocity space simply yields the net rate of work done by $E_y$, $\int C_{E_y}(v_x,v_y) dv_x dv_y = j_yE_y$, and we find the integration to be positive. We thus identify the whole population of reflected ions as experiencing net energization, with the velocity-space signature of this energization process, shock-drift acceleration, given by figure 5(b,d).

We have now connected the Lagrangian picture of shock-drift acceleration with the Eulerian picture provided by the FPC technique, and we conclude this section noting that while shock-drift acceleration has been studied extensive theoretically and numerically (e.g. Gedalin 1996a,b, 1997; Gedalin, Newbury & Russell 2000; Park et al. 2013; Guo et al. 2014a,b; Park et al. 2015; Gedalin et al. 2018; Xu et al. 2020), the velocity-space signature of shock-drift acceleration provides a new perspective on the energization of the ions in phase space via this process. In both cases, we understand that a portion of the distribution of ions is reflected via the magnetic field gradient and return upstream, where they can gain energy via the motional electric field. Although the single-particle trajectory in phase space guides our understanding of where we expect the ions to be gaining energy, using the FPC technique enables us to see clearly the exact region of phase space in which ions are being energized via shock-drift acceleration.

This confirmation of the velocity-space signature of shock-drift acceleration in a self-consistent simulation is a vitally important step for the comparison to measured velocity-space signatures of energy exchange using in situ spacecraft measurements; however, it is also interesting that the velocity-space signature of shock-drift acceleration is unchanged between the idealized model and a self-consistent simulation given the additional physics of the self-consistent simulation: the finite shock width and cross-shock electric field. We explore the reasons for the excellent agreement despite these two key differences between the simulation and the idealized model in Appendix E, where we find the cross-shock electric field assists in reflecting ions, which allows the ions to traverse further back upstream and gain additional energy via shock-drift acceleration. Thus, while the combination of the finite shock width and cross-shock electric field quantitatively changes the population of ions that are reflected, the qualitative signature of energization in velocity space via shock-drift acceleration remains unchanged.

5. Field–particle correlation analysis: electrons

5.1. Velocity-space signature of electron energization

We now examine the energization of the electrons by the simulated perpendicular collisionless shock. Similar to figure 3 for the ions, figure 6(a) shows the electron distribution function $f_e(v_x,v_y)$, figure 6(b) shows the $C_{E_x}$ and figure 6(c) shows the $C_{E_y}$ components of the FPC, (3.5) and (3.6). As shown in figure 1(e), the electron distribution broadens through the entire shock ramp, so we plot in figure 6 the results of the FPC analysis at $x_B = 21.8 d_i$, where the cross-shock electric field peaks. In figure 6(a), the centre of the distribution in the shock–rest frame is displaced away from the origin to $v_x<0$ and $v_y<0$ due to the particle drifts in the varying electric and magnetic fields through the shock ramp. Because the thermal width of the electron velocity distribution is much larger than the net drift of the distribution in the shock-rest frame, computing the $C_{E_x}$ and $C_{E_y}$ correlations using (3.5) and (3.6) leads to the qualitative ‘two-lobed’ velocity-space signatures observed in figure 6(b,c). The small drifts—i.e. $|v_x/v_{te}| \ll 1$—lead to a slight asymmetry of the two-lobed structure, so we overplot contours of constant $C_{E_{x,y}}$ to make these slight asymmetries more visually apparent. Although the gain (red) and loss (blue) of electron energy largely cancels out upon integration over velocity space $(v_x,v_y)$, asymmetries in the two lobes lead to a non-zero net energization.

Figure 6. The electron distribution function $f_e(v_x,v_y)$ (a), and the $C_{E_x}$ (b) and $C_{E_y}$ (c) components of the FPC, (3.5) and (3.6), computed at $x_B = 21.8 d_i$ from the self-consistent Gkeyll simulation. Note that the FPC is computed in the shock–rest frame. The contours on the FPC plots, which are the same in both $C_{E_x}$ and $C_{E_y}$, make clear that $E_y$ leads to a net loss of electron energy, whereas $E_x$ yields a net increase of electron energy.

In figure 6(b), the asymmetry of the velocity-space signature leads to a net energization of the electrons by the cross-shock component of the electric field $E_x$. Additionally, in contrast to the shock-drift acceleration of the ions by the motional electric field $E_y$ seen in figure 3(c), we see in figure 6(c) that the electrons experience a net loss of energy due to the $E_y$ component of the electric field. As with the ion analysis, we now turn to an idealized model for the electron dynamics through the shock layer to understand these velocity-space signatures for the electron energization.

5.2. Adiabatic heating in an idealized perpendicular shock

Unlike the ions for which $L_{\text {shock}} \lesssim \rho _i$, the electron gyroradius is much smaller than the gradient length scale of the shock, $L_{\text {shock}} \gg \rho _e$, and we thus expect the electrons to stay well magnetized through the shock. Therefore, we adopt a simplified model for the electron dynamics through the shock by taking a shock ramp with a linearly increasing magnetic field over a length $L= 2 d_i$ and mass ratio $m_i/m_e=1836$, which satisfy the limit $L_{\text {shock}} \gg \rho _e$. The other parameters are the same as the simple shock model used for the ion analysis in § 4.2: a magnetic field increase of $B_{d}/B_{u}=4$, a constant and uniform motional electric field $E_y<0$ in the shock–rest frame, with the same ion and electron plasma betas such that the shock velocity is comparable to the self-consistent simulation $M_A = 4.9$ and $M_f = 3.0$. Note that this idealized model for the electrons has no cross-shock component of the electric field, $E_x=0$; the implications of this choice are discussed at length in the upcoming subsections.

In this model, the increase in the magnetic field magnitude through the ramp leads to a steady decrease in the $\boldsymbol {E} \times \boldsymbol {B}$ velocity as the plasma flows through the shock transition. In addition, the gradient of the magnetic field magnitude in the shock-normal direction induces a $\boldsymbol {\nabla } B$ drift in the $+y$ direction. In figure 7(a), we plot the profile of the perpendicular magnetic field $B_z(x)$ (blue) and the motional electric field $E_y(x)$ (red) along the shock normal direction, and figure 7(b) shows the trajectory of an electron in the $(x,y)$ plane as it flows through the shock ramp over $0 \le x/d_i \le 2$. The trajectory plot shows clearly the $\boldsymbol {\nabla } B$ drift in the $+y$ direction. As the electrons flow through the shock ramp over $0 \le x/d_i \le 2$ and undergo a $\boldsymbol {\nabla } B$ drift in the $+y$ direction, figure 7(c) shows that the net energization for a distribution of electrons $j_y E_y$ is positive.

Figure 7. (a) Profiles along the shock normal direction of the perpendicular magnetic field $B_z$ (blue) and the motional electric field $E_y$ (red), (b) trajectory in the $(x,y)$ plane of an electron as it traverses the finite-width ramp in the magnetic field and (c) the rate of work done by the electric field on the distribution of particles $j_y E_y$.

In the region where the perpendicular magnetic field changes magnitude, $0 \le x/d_i \le 2$, the $\boldsymbol {\nabla } B$ drift is anti-aligned with the motional electric field, thus leading to net energization of electrons. In fact, the rate of energization of the electrons by the $\boldsymbol {\nabla } B$ drift in the motional electric field is precisely that needed to conserve the first adiabatic invariant of the electron, i.e. the magnetic moment $\mu = m_e v_\perp ^2/2 B_z$. This energization via conservation of the electron's adiabatic invariant is thus often referred to as adiabatic heating.

We can show the relationship between the energization via the $\boldsymbol {\nabla } B$ drift and the conservation of the electron's magnetic moment by considering the change in the perpendicular kinetic energy of the electrons,

(5.1)\begin{equation} \frac{\textrm{d} m_e v_\perp^2/2}{\textrm{d} t}=q_e u_{\boldsymbol{\nabla} B} E_y,\end{equation}

where the magnitude of the $\boldsymbol {\nabla } B$ drift in the $+y$ direction is given by

(5.2)\begin{equation} u_{\boldsymbol{\nabla} B}= \frac{m_e v_\perp^2}{2 q_e B_z} \left(\frac{1}{B_z} \frac{\partial B_z}{\partial x} \right).\end{equation}

For the static fields in this model, the total time derivative is dominated by the $\boldsymbol {E} \times \boldsymbol {B}$ velocity, $\textrm {d} / \textrm {d} t = \partial /\partial t + u_x \partial /\partial x = u_{E \times B} \partial /\partial x$. Substituting $u_{E \times B}=E_y/B_z$, we can manipulate (5.1) to obtain

(5.3)\begin{equation} \frac{\partial }{\partial x} \frac{m_e v_\perp^2}{2B_z}= \frac{\partial \mu }{\partial x}=0,\end{equation}

proving that the electron's magnetic moment $\mu$ is conserved. As before, we now wish to connect this single-particle, Lagrangian picture of adiabatic heating with the Eulerian picture provided by the FPC technique.

5.3. Cross-shock electric field impact on velocity-space signature of adiabatic heating

We again employ the Vlasov-mapping technique, as in § 4.3, to reconstruct the electron distribution function through the idealized shock model and the resulting FPC velocity-space signatures. Figure 8(a) presents the resulting electron distribution function $f_e(v_x,v_y)$ and figure 8(b) shows the $C_{E_y}$ correlation at position $x_A = 1.8 d_i$ in the model, where the velocity-integrated energy transfer rate $j_y E_y$ is positive, as seen in figure 7(c). Note that the small shift (relative to $v_{te}$) of the distribution to $v_y>0$, due to the $\boldsymbol {\nabla } B$ drift, leads to an asymmetry in the velocity-space signature because of the $v_y^2$ weighting in (3.6) for $C_{E_y}$.⁷ Although a large part of the energy transfer rate represented by this two-lobed velocity-space signature cancels out upon integration over $v_y$, the slight asymmetry leads to a net positive energization of the electrons, which yields $j_y E_y>0$, as plotted in figure 7(c).

Figure 8. Comparison of the reconstructed electron distribution function (a) and $C_{E_y}$ component of the FPC (b) computed from this reconstruction to the self-consistently produced electron distribution function (c) and $C_{E_y}$ component of the FPC (d) from the Gkeyll simulation. Unlike with the ion comparison presented in figure 5, the model $C_{E_y}$ displays the opposite asymmetry from the $C_{E_y}$ computed from the self-consistent simulation. Even though the signatures are qualitatively similar, this first computation of $C_{E_y}$ suggests that $E_y$ is responsible for a net loss of energy through the shock.

We compare the predicted velocity-space signature of electron adiabatic heating from the idealized shock model in figure 8(b) to the Gkeyll simulation results, where we plot in figure 8(c) the $f_e(v_x,v_y)$ and figure 8(d) the $C_{E_y}$ from the simulation at position $x_B=21.8 d_i$. While the correlation $C_{E_y}$ from the simulation has the same qualitative, two-lobed structure indicative of drift energization of the electrons, the net drift of the electron velocity distribution has $v_y<0$, and therefore the resulting asymmetry in $C_{E_y}$ has the opposite overall sign, which leads to a net loss of energy for the electrons due to $E_y$. How can we reconcile these apparently contradictory results for the electron energization by the motional electric field $E_y$, particularly given that we have shown in (5.3) that the $\boldsymbol {\nabla } B$ drift in the $+y$ direction leads to adiabatic heating of the electrons?

We can check if the adiabatic invariant of a distribution of electrons,

(5.4)\begin{equation} \mu_e= \frac{1}{\displaystyle \int f_e \, \textrm{d} \boldsymbol{v}} \left\{ \int \frac{m_e [ (v_x - u_x)^2 + (v_y - u_y)^2 ]}{2 B_z} f_e \, \textrm{d} \boldsymbol{v} \right\},\end{equation}

is constant through the shock, and whether our intuition about how the electron temperature should increase through a magnetic field gradient has merit. Note that in the 1D-2V geometry of the simulated perpendicular shock where both velocity coordinates are perpendicular to the magnetic field, the adiabatic invariant of the distribution of electrons can be simplified to

(5.5)\begin{equation} \mu_e = \frac{T_{{\perp},e}}{B_z}, \end{equation}

where $T_{\perp ,e}$ is the temperature of the electrons perpendicular to the magnetic field. Indeed, as shown in figure 9, the adiabatic invariant for the distribution of electrons, $\mu _e$, is well conserved through the shock transition.

Figure 9. The electron adiabatic invariant, $\mu _e = T_\perp /B_z$ (blue solid), the electron temperature (red solid) and the magnetic field (black dashed) normalized to their value upstream and plotted through the shock. The electron temperature rises commensurate with the compression of the magnetic field such that the electron $\mu _e$ is well conserved through the shock.

This result suggests that the cross-shock electric field is complicating the electron dynamics and energy exchange through the shock. Here we can leverage the Vlasov-mapping model to explore the effect of the cross-shock electric field on the energetics of the electrons. Figure 10(b) presents a comparison of the electron trajectories, figure 10(c) shows the rate of energization of the electrons by the electric field $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ and figure 10(d) presents the cumulative electron energization integrated from upstream $\int _{x_{\textrm {up}}}^x \,\textrm {d} \boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ for two models: (i) a ‘full model’ (solid) which integrates electron trajectories in the self-consistently produced electromagnetic fields in figure 10(a); and (ii) a ‘zero $E_x$ model’ (dashed), in which we artificially set the cross-shock electric field to zero.

Figure 10. (a) Electromagnetic fields approximated from the self-consistent Gkeyll simulation. (b) Example electron trajectories for full model (black) and zero $E_x$ model (blue), which show qualitatively different drifts in the $y$ direction. (c) Rate of work done by the components of the electric field, $j_yE_y$ (red) and $j_xE_x$ (blue) for the full model (solid) and zero $E_x$ model (dashed), along with total $\boldsymbol {j} \boldsymbol {\cdot } \boldsymbol {E}$ (black). Note that the total energization (black, solid and dashed) is the same for both cases. (d) Cumulative work done integrated from upstream $\int _x \boldsymbol {j} \boldsymbol {\cdot } \boldsymbol {E}$. The electrons experience adiabatic heating in both cases, although the detailed mechanisms of energization involve qualitatively different drifts.

The trajectories in figure 10(b) show a clear qualitative difference between the zero $E_x$ model and the full model: in the zero $E_x$ model (blue), the transverse drift through the shock ramp is relatively weak and in the $+y$ direction, whereas the full model (black) yields a much larger amplitude drift in the opposite direction. This qualitative difference in the transverse drift direction explains the stark differences in $C_{E_y}$ between the idealized model in figure 8(b) and the simulation in figure 8(d).

Looking at the rate of electron energization by the electric field in figure 10(c), we indeed see that $j_{ye} E_y$ (red dashed) is positive for the zero $E_x$ model (and is the only means of energization of the electrons because $E_x=0$), but $j_{ye} E_y$ is negative for the full model (red solid). However, when the energization by $j_{xe} E_x$ (blue solid) is combined with $j_{ye} E_y$ (red solid), the net rate of energization $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ of the two models is exactly the same (black solid and red dashed overlap). Therefore, although the dynamics of the electrons differ qualitatively in the presence or absence of the cross-shock electric field, their energization is the same in either case. To explain this puzzling finding, we exploit the limit $L_{\text {shock}} \gg \rho _e$ to execute a guiding-centre drift analysis of the electron energization.

5.4. Guiding-centre drift analysis of electron energization

In the idealized model presented in § 5.2, there are only two drifts: an $\boldsymbol {E} \times \boldsymbol {B}$ drift in the $x$ direction, due to a constant $E_y$ through the magnetic field ramp, and a $\boldsymbol {\nabla } B$ drift in the $y$ direction, due to the linearly increasing magnetic field. The introduction of a cross-shock electric field, along with the transition from a single-particle picture to a distribution of particles, adds several new drifts to the full list of potentially dynamically important drifts. We now not only have an $\boldsymbol {E} \times \boldsymbol {B}$ drift that has a component in the $x$ direction, due to the motional electric field $E_y$ supporting the incoming supersonic flow, but also a new $y$ component, due to the cross-shock electric field $E_x$. For a distribution of electrons the $\boldsymbol {\nabla } B$ drift in the $y$ direction is modified from its single particle form,

(5.6)\begin{equation} u_{\boldsymbol{\nabla} B}^{SP} = \frac{m_e v_{{\perp},e}^2}{2 q_e B^2_z} \frac{\partial B_z}{\partial x} \hat{\boldsymbol{y}}, \end{equation}

to

(5.7)\begin{equation} u_{\boldsymbol{\nabla} B} = \frac{1}{q_e n_e} \frac{p_{{\perp},e}}{B^2_z} \frac{\partial B_z}{\partial x} \hat{\boldsymbol{y}}; \end{equation}

where $p_{\perp ,e}$ is the electron perpendicular pressure,

(5.8)\begin{equation} p_{{\perp},e} = \frac{1}{2} m_e \int (\boldsymbol{v}_\perp{-} \boldsymbol{u}_{e \perp})^2 f_e \, \textrm{d} \boldsymbol{v}. \end{equation}

We now must also consider a polarization drift in the $x$ direction, and we note that the total time derivative as the electrons flow through the shock is dominated by the convective contribution $\textrm {d} / \textrm {d} t = \partial /\partial t + \boldsymbol {U} \boldsymbol {\cdot } \boldsymbol {\nabla } \simeq U_x \partial /\partial x$, which gives

(5.9)\begin{equation} u_{\textrm{d} E/\textrm{d} t} = \frac{1}{\varOmega_{cs} B_z} \frac{\textrm{d} \boldsymbol{E}}{\textrm{d} t} \sim \frac{1}{\varOmega_{cs} B_z} U_x \frac{\partial E_x}{\partial x} \hat{\boldsymbol{x}}. \end{equation}

Finally, we have the magnetization drift, $\boldsymbol {\nabla } \times \boldsymbol {M}$, a bulk drift in the $y$ direction due to the increasing density through the shock ramp,

(5.10)\begin{equation} u_{\boldsymbol{\nabla} \times \boldsymbol{M}} = \frac{\boldsymbol{\nabla} \times \boldsymbol{M}}{q_e n_e} =\frac{1}{q_e n_e} \boldsymbol{\nabla} \times \left(- \frac{p_{{\perp}, e}}{B_z} \right) = \frac{1}{q_e n_e} \frac{\partial }{\partial x} \left(\frac{p_{{\perp},e}}{B_z} \right)\hat{\boldsymbol{y}}, \end{equation}

where the magnetization vector $\boldsymbol {M}$ (Hazeltine & Waelbroeck 1998) is given by

(5.11)\begin{equation} \boldsymbol{M} ={-}p_\perp \frac{\boldsymbol{B}}{|\boldsymbol{B}|^2}. \end{equation}

Although there can, in principle, be a polarization drift in the $y$ direction due to the variation in $E_y$ through the shock, in the shock rest frame $E_y$ changes very little, so this drift is negligible. We note that the $\boldsymbol {\nabla } B$ drift and magnetization drift can be combined to form the diamagnetic drift,

(5.12)\begin{align} u_{\mbox{diamag}} & = u_{\boldsymbol{\nabla} B} + u_{\boldsymbol{\nabla} \times \boldsymbol{M}} = \frac{1}{q_e n_e} \left[ \frac{p_{e \perp}}{B^2_z} \frac{\partial B_z}{\partial x} + \frac{\partial }{\partial x} \left(\frac{p_{e \perp}}{B_z} \right) \right], \nonumber\\ & = \frac{1}{q_e n_e} \frac{1}{B_z} \frac{\partial p_{e \perp}}{\partial x} ={-}\frac{1}{q_e n_e} \frac{\boldsymbol{\nabla} p_{e \perp} \times \boldsymbol{B}}{|\boldsymbol{B}|^2}. \end{align}

This generalization from the single-particle picture to a distribution of particles is somewhat subtle and often dubbed Spitzer's paradox from the early work on plasma equilibria in a fusion context (Spitzer 1952, 1962; Qin et al. 2000). While the concept of pressure is ill-defined for a single particle, magnetic field gradients must inevitably be balanced by pressure gradients in equilibrium because the $\boldsymbol {\nabla } B$ drift depends on the particles’ velocity, and thus different parts of the distribution of particles will experience different $\boldsymbol {\nabla } B$ drifts. As a consequence, the pressure of the plasma must change through the magnetic field gradient, presuming of course that the plasma is magnetized. To understand this generalization from the single-particle picture to a distribution of particles, a detailed derivation of these drifts from the first moment of the Vlasov equation can be found in Appendix F.

We plot these drifts in the shock-rest frame in figure 11(a) as a function of the distance through the shock, which shows that there is a clear ordering of the magnitude of these drifts: both components of the $\boldsymbol {E} \times \boldsymbol {B}$ drift are dominant, the $\boldsymbol {\nabla } B$ and magnetization $\boldsymbol {\nabla } \times \boldsymbol {M}$ drifts provide equal contributions at a smaller amplitude, and the polarization drift is yet smaller. To demonstrate that the sum of these drifts indeed completely describes the bulk flow of the plasma in the shock–rest frame, we show in figure 11(b) that the first moment of the electron distribution agrees well with the sum of the drifts for both $x$ and $y$ components.

Figure 11. (a) A comparison of the major drifts through the shock evaluated in the shock–rest frame of reference: (i) $\boldsymbol {E} \times \boldsymbol {B}$ drift in $x$ (black) and $y$ (green dashed), (ii) the $\boldsymbol {\nabla } B$ drift in $y$ (blue), (iii) the magnetization drift in $y$ (red dashed), and (iv) the polarization drift in $x$ (magenta dashed–dotted). We check that these drifts sum to the total first moment computed from the electron distribution function (b) as well as determine how each of these drifts contributes to the overall energy exchange, $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ (c) and compare the $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ computed from these drifts to the total $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ computed from moments of the electron distribution function. Note that in comparing how each of these drifts contributes to $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$, we plot the polarization drift multiplied by $E_x$ (green dashed–dotted), the $\boldsymbol {\nabla } B$ and magnetization drifts multiplied by $E_y$ in the shock-rest frame (blue and red dashed, respectively), and the total $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ arising from both components of the $\boldsymbol {E} \times \boldsymbol {B}$ flow (black), as we expect the total energization due to the $\boldsymbol {E} \times \boldsymbol {B}$ flow to be zero.

The rate of energization of the electrons in the drift approximation is equal to the rate of work done by the electric field on the drifting electrons, given by

(5.13)\begin{equation} \boldsymbol{j}_{d,e} \boldsymbol{\cdot} \boldsymbol{E} = q_e n_e \boldsymbol{U}_d \boldsymbol{\cdot} \boldsymbol{E}, \end{equation}

where $\boldsymbol {U}_d$ is the total drift motion, which can be decomposed into the contributions by the individual drifts identified above in, e.g. (5.7)–(5.11). This analysis follows in the vein of studies of electron energization in reconnection using the guiding centre approximation (Dahlin, Drake & Swisdak 2014, 2015, 2016, 2017; Drake et al. 2019). However, we emphasize again that we have generalized from a single-particle picture to a distribution function picture as part of our goal to understand the energization of the electrons in phase space, and thus our equivalent guiding centre energization equation has additional terms, such as the energization of the plasma via the magnetization drift (5.10), which do not appear in the guiding centre energy equation of a single particle, e.g. (1) in Dahlin et al. (2014) and Dahlin et al. (2016).⁸ We reiterate that the perspective provided by the Lagrangian point-of-view, wherein one considers the energization of individual particles, has significant merit, but that the Eulerian perspective has its own advantages, and to obtain the Eulerian point-of-view we must generalize from single particles to distributions of particles.

In figure 11(c), we plot the rate of electron energization by each of these four drifts in the shock–rest frame. The first key takeaway from this figure is that, although the two components of the $\boldsymbol {E} \times \boldsymbol {B}$ drift dominate the total drift motion, there is no net work done by the $\boldsymbol {E} \times \boldsymbol {B}$ drift as we expect because $(\boldsymbol {E} \times \boldsymbol {B}) \boldsymbol {\cdot } \boldsymbol {E}=0$. So, while the energization arising from one component of the $\boldsymbol {E} \times \boldsymbol {B}$ drift may appear large, summed over all components the energization must be zero, as shown in figure 11(c). In this regard, separately plotting the contributions to the rate of electron energization $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ by the different components of the $\boldsymbol {E} \times \boldsymbol {B}$ drift can be somewhat misleading, as both components of the $\boldsymbol {E} \times \boldsymbol {B}$ flow are much larger than the other drifts.

That the full $\boldsymbol {E} \times \boldsymbol {B}$ drift leads to zero net energization of the particles also explains the puzzling finding in our single-particle modelling of the electron energization shown in figure 10(c). Although the $j_{e_x} E_x$ and $j_{e_y} E_y$ from the full model were larger than and significantly different from the zero $E_x$ model, when summed they yielded the same net rate of energization of the electrons as the zero $E_x$ model. This cancellation is exactly the result of the two components of ‘energization’ from the $\boldsymbol {E} \times \boldsymbol {B}$ flow cancelling, as we know they must. The remaining net energization is then solely from the other drifts and their alignment with the motional electric field $E_y$.

In figure 11(d), we check that the total rate of energization of the electrons by the electric field in the shock–rest frame, $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$, agrees with the sum of the energization by the $\boldsymbol {\nabla } B$, magnetization and polarization drifts, finding good agreement. Note that for this comparison, we are computing the electron current in the shock–rest frame from the electron distribution function, i.e.

(5.14)\begin{equation} \boldsymbol{j}_e = \int q_e (\boldsymbol{v}' - U_{\text{shock}} \hat{\boldsymbol{x}}) f_e \, \textrm{d}\boldsymbol{v}'. \end{equation}

After the comparison between the full model and zero-$E_x$ model in figure 10(c,d) revealed that the total $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ was roughly equivalent between the two models, when the zero-$E_x$ model only had energy exchange, due to the alignment of the $\boldsymbol {\nabla } B$ drift with the motional electric field, we might have anticipated that the only energy gain arose from this same adiabatic heating process from the idealized model in § 5.2. Importantly though, we see from the drift analysis in figure 11(c) that the model in § 5.2 must be generalized to the case when a distribution of particles is drifting.

While the energy gain by a single electron $\boldsymbol {\nabla } B$ drifting in the model fields is exactly the energy gain required for that single electron's magnetic moment $\mu$ to be conserved, in the self-consistent simulation it is not only the $\boldsymbol {\nabla } B$ drift that ensures the electron distribution's adiabatic invariant $\mu _e$ is well conserved in the shock. We also have an equal contribution to the energization from the magnetization drift. Together, as shown by (5.12), the $\boldsymbol {\nabla } B$ and magnetization drifts are equivalent to the diamagnetic drift, $u_{\mbox {diamag}} = u_{\boldsymbol {\nabla } B} + u_{\boldsymbol {\nabla } \times \boldsymbol {M}}$, so another perspective on the electron energization via adiabatic heating is the adiabatic invariant of a distribution of electrons, $\mu _e$ in (5.4), is conserved due to the alignment of the diamagnetic drift and the motional electric field, $E_y$.

Two important questions remain: (i) what is responsible for the slight disagreement between the energy gain arising from the $\boldsymbol {\nabla } B$ and magnetization drifts and $\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ from velocity moments of the electron distribution function in figure 11(d)?; and (ii) even if the physics of adiabatic heating is the same with the electrons gaining energy via the alignment of drifts with the motional electric field, is the velocity-space signature of adiabatic heating the same as that predicted by the Vlasov-mapping model in figure 10(b)? To answer the first question, we have performed a more realistic mass ratio simulation, $m_i/m_e = 400$, in Appendix G where we find better agreement between the energy gain arising from the $\boldsymbol {\nabla } B$ and magnetization drifts and $\,\boldsymbol {j}_e \boldsymbol {\cdot } \boldsymbol {E}$ from velocity moments of the electron distribution function. Thus, the small disagreement between these methods of measuring the electron's energy gain in the $m_i/m_e = 100$ simulation is simply a result of the fact that the electron gyroradius is not asymptotically smaller than the shock's extent. To answer the second question, we seek a means of eliminating the $\boldsymbol {E} \times \boldsymbol {B}$ component of the energy exchange in the FPC.

5.5. Velocity-space signature of adiabatic electron heating

If the large $\boldsymbol {E} \times \boldsymbol {B}$ flows are polluting the analysis of the overall exchange of the energy, when fundamentally the electron heating principally arises from the alignment of the $\boldsymbol {\nabla } B$ and magnetization drifts with the motional electric field, we return to the velocity-space signature plotted in figures 6(c) and 8(d) to determine how we might remove the contribution from the large $\boldsymbol {E} \times \boldsymbol {B}$ flows to the FPC signal. Guided by the fact that the transverse electric field component $E_y$ governs the adiabatic heating of the electrons, we seek to eliminate the large contribution to the rate of energization associated with the $y$-component of the $\boldsymbol {E} \times \boldsymbol {B}$ drift (which is ultimately cancelled by energization associated with its $x$-component, as we show in Appendix B). Therefore, we transform to a frame of reference moving in the transverse direction at the same velocity as the $y$-component of the $\boldsymbol {E} \times \boldsymbol {B}$ drift, ${\boldsymbol {U}}_{td}= -E_x/B_z {\hat {\boldsymbol {y}}}$. We define the transverse drift frame of reference: (i) in the $x$ direction, or shock-normal direction, the shock is at rest; (ii) in the $y$ direction, or transverse direction, the frame moves at a velocity equal to the $y$-component of the local $\boldsymbol {E} \times \boldsymbol {B}$ drift in the shock–rest frame.

Critically, the electric field transforms from the shock–rest frame ($sf$, unprimed) to the transverse drift frame ($td$, double primed), as

(5.15)\begin{equation} \boldsymbol{E}''= \boldsymbol{E} + \boldsymbol{U}_{td} \times \boldsymbol{B} = \boldsymbol{E} + \left(-\frac{E_x}{B_z}{\hat{\boldsymbol{y}}}\right) \times B_z {\hat{\boldsymbol{z}}} = \boldsymbol{E} -E_x {\hat{\boldsymbol{x}}} = E_y {\hat{\boldsymbol{y}}}, \end{equation}

where we have assumed $E_z=0$ in this 1D-2V perpendicular shock. Therefore, in the transverse drift frame, the cross-shock electric field is zero, $E_x''=0$, and the motional electric field is unchanged from the shock–rest frame, $E_y''=E_y$. The $y$-component of the velocity coordinate in the transverse drift frame is

(5.16)\begin{equation} v_y''= v_y - U_{td,y} = v_y + \frac{E_x}{B_z}. \end{equation}

For completeness, the $x$ velocity coordinate and magnetic field are unchanged from the shock–rest frame, $v_x''=v_x$ and $\boldsymbol {B}''=\boldsymbol {B}$.

Although the transverse drift frame changes with position through the shock as the cross-shock electric field changes along the normal direction, one can determine this frame of reference at any position from a local, single-point measurement of the electric field. In contrast, other drifts generally depend on gradients, and therefore cannot be uniquely specified using only single-point measurements. The benefit of the transverse drift frame of reference is not only that the rate of energization associated with the total $\boldsymbol {E} \times \boldsymbol {B}$ drift is equal to zero, which is true in any frame of reference, but also that the rates of energization associated with each component of the $\boldsymbol {E} \times \boldsymbol {B}$ drift are separately zero.⁹

Because $E_x''=0$ in the transverse drift frame, we need only examine the $y$ correlation, $C_{E_y''}(v_x'',v_y'')$, to explore the energization of the electrons. We plot in figure 12(a) the correlation $C_{E_y''}(v_x'',v_y'')$ in the transverse drift frame at position $x_B=21.8 d_i$, and compare it with the corresponding correlation $C_{E_y}(v_x,v_y)$ in the shock–rest frame shown in figure 12(b) (a repeat of figure 6c).

Figure 12. A comparison of the FPC from $E_y$, where $v_x$ and $v_y$ are shifted to the shock-rest frame and local $\boldsymbol {E} \times \boldsymbol {B}$ frame respectively, i.e. the transverse drift frame, $C_{E''_y}(v''_x, v''_y)$ (a) to our previous computation of the FPC from $E_y$ using only a frame transformation in $v_x$ to the shock–rest frame, $C_{E_y}(v_x,v_y)$ (b). Note that panel (b) is a repeat of figure 6(c) and panel (d) is a repeat of figure 8. While the previous correlation, $C_{E_y}$, suggested the electrons were losing energy in this degree of freedom, the newly transformed correlation, $C_{E''_y}$, demonstrates that the electrons are in fact gaining energy once the $\boldsymbol {E} \times \boldsymbol {B}$ motion in this degree of freedom is subtracted. This energization mechanism, caused by an alignment of the $\boldsymbol {\nabla } B$ and magnetization drifts with the motional electric field, $E_y$, is the same mechanism responsible for energizing the electrons in the idealized model, and the velocity-space signature for this adiabatic heating now exactly matches the results of the idealized model presented in §§ 5.2 and 5.3.

We note two things immediately from figure 12. First, transformation to the transverse drift frame has switched the sign of the net rate of electron energization (from integrating the correlation over velocity space) relative to the case in the shock–rest frame, made clear by counting the overplotted equally spaced contours in the blue-red, two-lobe structure of each case. Second, the correlation $C_{E_y''}(v_x'',v_y'')$ in the transverse drift frame is strikingly similar to the correlation found in figure 8(b) computed from the reconstructed distribution function from the idealized model. Both points are perhaps unsurprising, as this transformation has removed the energy exchange in the $v_y$ degree of freedom due to the $\boldsymbol {E} \times \boldsymbol {B}$ flow, and we are left with a similar energization mechanism found in the idealized model: alignment of the $\boldsymbol {\nabla } B$ and magnetization drifts with the motional electric field, $E_y$. The alignment of these two combined drifts (constituting the diamagnetic drift) with the motional electric field is exactly what is required for the electron distribution's adiabatic invariant $\mu _e$ to be well-conserved and for the electrons to gain energy through the increasing magnetic field of the perpendicular shock. Figure 12(a) contains the second key result of this paper, the velocity-space signature of adiabatic electron heating.

It is worth emphasizing that transformation to the transverse drift frame not only revealed the same velocity-space signature for adiabatic heating as in the idealized model analysed in §§ 5.2 and 5.3, but also allowed the extraction of an energization signature which was buried in the large background energy exchange from the $\boldsymbol {E} \times \boldsymbol {B}$ flow. Ultimately, the adiabatic heating via alignment of the $\boldsymbol {\nabla } B$ and magnetization drifts with the motional electric field was masked by the large oscillation of energy between the $v_x$ and $v_y$ degrees of freedom due to the two components of the $\boldsymbol {E} \times \boldsymbol {B}$ motion. The need to carefully consider the frame of reference in which the energization analysis is performed, especially due to the impact the frame of reference choice has on the cross-shock electric field, has been pointed out in previous studies (Goodrich & Scudder 1984).

One must therefore exercise extreme caution in attributing energization to a particular drift when separately considering the work done by the different components of the electric field. In figure 10(c), for example, it would be easy to attribute erroneously the energization of the electrons to the cross-shock electric field $E_x$. Rather, the energization arising from the $x$-component of the $\boldsymbol {E} \times \boldsymbol {B}$ drift and $E_x$ is exactly cancelled by a loss of energy arising from the $y$-component of the $\boldsymbol {E} \times \boldsymbol {B}$ drift and $E_y$. The net result is a much smaller positive rate of energization due to the motional electric field $E_y$ and the remaining drifts in the transverse direction.

That the addition of the magnetization drift did not qualitatively alter the velocity-space signature of adiabatic electron heating can be understood by considering the general qualitative features of drift energization. In the drift approximation for a warm electron distribution relevant to most heliospheric plasmas, it is often true that the magnitude of the drift (not including the dominant $\boldsymbol {E} \times \boldsymbol {B}$ drift, which cannot energize particles) is much less than the thermal velocity of the electrons, $U_{d,e} \ll v_{te}$. In this case, the centre of the drifting velocity distribution is offset from the origin of velocity space, but this offset will be much smaller than the thermal width of the distribution. When the correlation is computed by taking the velocity derivatives of $f_e$ and multiplying by the appropriate component of the electric field weighted by $|\boldsymbol {v}|^2$, the result will generally produce a two-lobed velocity-space signature, qualitatively similar to that shown in figure 12(a).

A careful consideration of the drifts is essential to interpret properly the mechanism responsible for the adiabatic heating. However, much of the interest in collisionless shocks is focused not on adiabatic heating and instead on identifying mechanisms of non-adiabatic heating. In fact, because the observed temperature increase in the electrons is entirely in the perpendicular temperature to conserve the electron distribution's adiabatic invariant, even if the electron response is initially adiabatic, this anisotropy will itself be a source of instabilities, such as the whistler anisotropy instability (Gary 1993; Gary, Liu & Winske 2011; Wilson et al. 2013b, 2020), which will further complicate the energy exchange.

Likewise, the transverse drift frame relies on the electrons being magnetized so that their $\boldsymbol {E} \times \boldsymbol {B}$ motion is well-defined. While the transformation is sensible here because $L_{\text {shock}} \gg \rho _e$, there are many observations of collisionless shocks where the shock ramp is not asymptotically larger than the electron gyroradius (Hobara et al. 2010) and thus would warrant caution in the application of the transverse drift frame to reveal any masked energization signatures such as the velocity-space signature of adiabatic heating in figure 12(a). In addition, the transverse-drift frame is a non-inertial frame because the $\boldsymbol {E} \times \boldsymbol {B}$ velocity is changing through the shock. While the transformation to the instantaneous transverse-drift frame at a single point in configuration space in a simulation is a perfectly reasonable frame transformation, care will be required in applying this technique to the analysis of spacecraft data, which inevitably average over a small volume of configuration space.

Nevertheless, the velocity-space signatures of the mechanisms that govern non-adiabatic heating are likely to be entirely distinct from the simple two-lobed appearance of adiabatic heating, and so will be easy to distinguish using a field–particle correlation analysis. In addition, it is useful to first characterize the ‘background’ signature of adiabatic heating, as represented by the typical velocity-space signature shown in figure 12(a). In cases where the adiabatic response produces a temperature anisotropy that drives instabilities, the energetic response, and thus the electron's velocity-space signature through the shock, is likely to be characterized by a combination of both adiabatic and non-adiabatic signatures.

6. Summary and future outlook

This paper presents the first attempt to characterize the energy exchange in a collisionless shock using the field–particle correlation technique (Klein & Howes 2016; Howes et al. 2017). We have examined a self-consistent perpendicular collisionless shock using the continuum Vlasov–Maxwell solver in the Gkeyll simulation framework and our results can be summarized as follows:

1. (i) we have identified the velocity-space signatures of shock-drift acceleration of ions in figure 5(d) and adiabatic heating of electrons in figure 12(a);

2. (ii) using simplified models of single-particle motion through idealized models of the electromagnetic fields through the shock transition, we identified the conditions under which we expect to observe these velocity-space signatures for these energization processes.

   1. (a) We determined that the velocity-space signature of shock-drift acceleration can be seen clearly in a reflected ion population and is robust to the presence of the finite shock width and cross-shock electric field that arise in the self-consistent simulation.

   2. (b) For the electrons, we determined it was critical to eliminate the energization associated with the separate components of the $\boldsymbol {E} \times \boldsymbol {B}$ drift by transforming to the transverse drift frame of reference, as the large $\boldsymbol {E} \times \boldsymbol {B}$ drifts in both $x$ and $y$ due to the incoming supersonic flow and cross-shock electric field significantly obscured the effect of the energization of the electrons via other drifts;

3. (iii) finally, we observed a general strength of this method of analysis using the field–particle correlation technique, which goes beyond identification of where $\boldsymbol {j} \boldsymbol {\cdot } \boldsymbol {E}$ is positive and which components are positive, and furthermore captured the subtleties of a generalization from a single-particle picture to a distribution of particles.

   1. (a) In the case of the ions, a cursory glance at the $x$ component of $\boldsymbol {j} \boldsymbol {\cdot } \boldsymbol {E}$, which is overall negative owing to the slowing down of the bulk distribution, obscures the role the cross-shock electric field plays in increasing the reflected fraction of ions—see Appendix E for further details.

   2. (b) In the case of electrons, just using the components of $\boldsymbol {j} \boldsymbol {\cdot } \boldsymbol {E}$ would completely miss the actual source of energization, as one might expect that the positive $x$ component of $\boldsymbol {j} \boldsymbol {\cdot } \boldsymbol {E}$ corresponds to net energization via the cross-shock electric field. In fact, it is not the cross-shock component of the electric field $E_x$ which leads to the observed energization of the electrons through the shock, but rather the motional electric field $E_y$. Additionally, while in the single-particle picture the energization is simply the alignment of said motional electric field with the $\boldsymbol {\nabla } B$ drift, the generalization to a distribution of electrons introduces an additional drift, the magnetization drift, which when combined with the $\boldsymbol {\nabla } B$ drift, forms the diamagnetic drift. When aligned with the motional electric field, these two drifts provide the necessary energization for the adiabatic invariant of the electron distribution to be well conserved through the shock.

The work presented here is only the beginning of a program of studies to determine, in general, how we may be able to leverage the full information contained in the particle velocity distribution function to ascertain the details of the energy exchange in a collisionless shock. While historically the Lagrangian perspective of examining how individual particles gain and lose energy has led to enormous improvements in our understanding of the dynamics and energetics of collisionless shocks, this complementary Eulerian approach, directly analysing the energy exchange in phase space using the field–particle correlation technique, has significant value for interpreting both simulation and spacecraft data. Especially when advances in spacecraft instrumentation provide ever higher resolution and higher cadence particle velocity distribution measurements of collisionless shocks (Chen et al. 2018; Goodrich et al. 2018), the time is ripe for fully exploiting the information contained in phase space to provide a deeper understanding of the mechanisms of particle energization at a collisionless shock.

Further studies of higher dimensional, higher magnetosonic Mach number, and more general geometry collisionless shocks are of the utmost importance. As reviewed in the introduction, there is a large variety of processes not considered in this study that have been studied previously as potential energization mechanisms, from shock surfing acceleration to diffusive shock acceleration. As with the body of work using the field–particle correlation technique for analysing dissipation via resonant processes, we will require a systematic study of all of these processes if we have any hopes of distinguishing their velocity-space signatures. We may expect certain energization processes, such as diffusive shock acceleration and the fast Fermi process, which also rely on particle reflection, may have qualitatively similar velocity-space signatures to shock-drift acceleration but still contain all the requisite information to identify the particular process locally occurring, just as we can use information, such as the velocity around which a resonant wave process is identified, to characterize the particular waves which are resonantly energizing the plasma.

We expect, as we move to higher dimensionality, higher magnetosonic Mach number and more general shock geometry, our analysis will also be further complicated by upstream kinetic instabilities such as those observed in the Earth's bow shock. The velocity-space signatures of these instabilities are of equal importance to characterize and study using the field–particle correlation technique. While the instantaneous field–particle correlation technique employed in this study was well-suited to the impulsive ion energization via shock-drift acceleration and the steady electron energization via adiabatic heating, we may require finite-time correlations to characterize the energy exchange within the upstream fluctuations and turbulence of the shocks in exact analogy with previous field–particle correlation studies (Klein et al. 2017, 2020; Horvath et al. 2020). Nevertheless, an exciting frontier awaits in applying the field–particle correlation technique to the distribution functions produced by more realistic collisionless shock simulations and classifying the observed velocity-space signatures of particle energization.