2.1. Reduced Vlasov–Maxwell equations

Using the reduced Hamiltonian (2.1), the reduced Vlasov–Maxwell equations are now expressed as follows. First, the reduced equations of motion are given in Hamiltonian canonical form as

(2.2)\begin{gather} \textrm{d}\bar{\boldsymbol{x}}/\textrm{d}t = \partial\bar{H}/\partial\bar{\boldsymbol{p}}= \bar{\boldsymbol{v}}+\partial\bar{\varPsi}/\partial\bar{\boldsymbol{p}}, \end{gather}

(2.3)\begin{gather}\textrm{d}\bar{\boldsymbol{p}}/\textrm{d}t ={-}\bar{\boldsymbol{\nabla}}\bar{H} ={-}e\bar{\boldsymbol{\nabla}}\varPhi + (e/c)\,\bar{\boldsymbol{\nabla}} \boldsymbol{A}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{v}}- \bar{\boldsymbol{\nabla}}\bar{\varPsi}. \end{gather}

If the reduced force equation (2.3) is written in terms of $\bar {\boldsymbol {v}}$, we find

(2.4)\begin{equation} m\,\frac{\textrm{d}\bar{\boldsymbol{v}}}{\textrm{d}t} = e\boldsymbol{E} + \frac{e}{c}\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\boldsymbol{\times}\boldsymbol{B} + \bar{\boldsymbol{\nabla}}\boldsymbol{E}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{\rm \pi}} + \bar{\boldsymbol{\nabla}}\boldsymbol{B}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{\mu}}, \end{equation}

where we used (2.2) on the right-hand side and the ponderomotive force

(2.5)\begin{equation}{-}\bar{\boldsymbol{\nabla}}\bar{\varPsi} ={-}\bar{\boldsymbol{\nabla}}\boldsymbol{A}\boldsymbol{\cdot}\left(-\frac{e}{c}\frac{\partial \bar{\varPsi}}{\partial \bar{\boldsymbol{p}}}\right) + \bar{\boldsymbol{\nabla}}\boldsymbol{E}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{\rm \pi}} + \bar{\boldsymbol{\nabla}}\boldsymbol{B}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{\mu}} \end{equation}

includes the reduced electric and magnetic dipole moments $(\bar {\boldsymbol {{\rm \pi} }}, \bar {\boldsymbol {\mu }}) \equiv (-\partial \bar {\varPsi }/\partial \boldsymbol {E}, -\partial \bar {\varPsi }/ \partial \boldsymbol {B})$ derived from the ponderomotive Hamiltonian. The reduced Vlasov equation is, therefore, expressed as

(2.6)\begin{equation} \frac{\partial \bar{f}}{\partial t} ={-}\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{\nabla}}\bar{f} - \frac{\textrm{d}\bar{\boldsymbol{p}}}{\textrm{d}t}\boldsymbol{\cdot}\frac{\partial \bar{f}}{\partial \bar{\boldsymbol{p}}}. \end{equation}

The reduced Maxwell equations, on the other hand, are expressed as

(2.7)\begin{gather} \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\boldsymbol{E} = 4{\rm \pi}\varrho \equiv 4{\rm \pi}\left(\bar{\varrho} - \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\bar{\mathbb{P}}\right), \end{gather}

(2.8)\begin{gather}\boldsymbol{\nabla}\boldsymbol{\times}\boldsymbol{B} - \frac{1}{c} \frac{\partial \boldsymbol{E}}{\partial t} =\frac{4{\rm \pi}}{c}\,\boldsymbol{J} \equiv \frac{4{\rm \pi}}{c} \left( \bar{\boldsymbol{J}} + \frac{\partial \bar{\mathbb{P}}}{\partial t} + c\boldsymbol{\nabla}\boldsymbol{\times}\bar{\mathbb{M}} \right), \end{gather}

with the source-free Maxwell equations

(2.9)\begin{equation} \left.\begin{aligned} \partial\boldsymbol{B}/\partial t + c\boldsymbol{\nabla}\boldsymbol{\times}\boldsymbol{E} = 0, \\ \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\boldsymbol{B} =0. \end{aligned} \right\} \end{equation}

In (2.7) and (2.8), the reduced charge and current densities $(\bar {\varrho }, \bar {\boldsymbol {J}})$ and the reduced polarization and magnetization $(\bar {\mathbb {P}},\bar {\mathbb {M}})$ are derived from the reduced Hamiltonian:

(2.10)\begin{equation} \int_{\bar{\boldsymbol{p}}}\,\bar{f}\,\partial\bar{H}(\bar{\boldsymbol{p}};\varPhi,\boldsymbol{A},\boldsymbol{E},\boldsymbol{B}) \equiv \partial\varPhi\bar{\varrho} - \partial\boldsymbol{A}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{J}}/c - \partial\boldsymbol{E}\,\boldsymbol{\cdot}\,\bar{\mathbb{P}} - \partial\boldsymbol{B}\,\boldsymbol{\cdot}\,\bar{\mathbb{M}}, \end{equation}

where the notation $\int _{\bar {\boldsymbol {p}}}$ indicates an integral over canonical-momentum space (as well as including a sum over particle species) and $\partial$ denotes either a space–time partial derivative $(\boldsymbol {\nabla }, \partial /\partial t)$ or an Eulerian variation $\delta$. Specifically, we find the definitions

(2.11)\begin{equation} \left(\bar{\rho}, \bar{\boldsymbol{J}}, \bar{\mathbb{P}} , \bar{\mathbb{M}} \right) \equiv \int_{\bar{\boldsymbol{p}}} \,\bar{f} \left( e, e\,\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}, - \frac{\partial \bar{\varPsi}}{\partial \boldsymbol{E}}, -\frac{\partial \bar{\varPsi}}{\partial \boldsymbol{B}}\right), \end{equation}

where contributions arise from reduced particles located at the field position (i.e. $\bar {\boldsymbol {x}} = \boldsymbol {x})$.

We also note that the reduced Maxwell equations (2.7) and (2.8) can be written in terms of the reduced Maxwell fields

(2.12)\begin{equation} \left. \begin{aligned} \bar{\mathbb{D}} = \boldsymbol{E} + 4{\rm \pi}\bar{\mathbb{P}}, \\ \bar{\mathbb{H}} = \boldsymbol{B} - 4{\rm \pi}\bar{\mathbb{M}}, \end{aligned} \right\} \end{equation}

as

(2.13)\begin{gather} \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\bar{\mathbb{D}} = 4{\rm \pi}\bar{\varrho}, \end{gather}

(2.14)\begin{gather}\boldsymbol{\nabla}\boldsymbol{\times}\bar{\mathbb{H}} - \frac{1}{c}\frac{\partial \bar{\mathbb{D}}}{\partial t} = \frac{4{\rm \pi}}{c} \bar{\boldsymbol{J}}, \end{gather}

which guarantees that the reduced charge conservation law

(2.15)\begin{equation} \frac{\partial \varrho}{\partial t} + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\boldsymbol{J} = \frac{\partial }{\partial t}\left(\bar{\varrho} - \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\bar{\mathbb{P}}\right) + \boldsymbol{\nabla}\boldsymbol{\cdot}\left( \bar{\boldsymbol{J}} + \frac{\partial \bar{\mathbb{P}}}{\partial t} + c\boldsymbol{\nabla}\boldsymbol{\times}\bar{\mathbb{M}} \right) = \frac{\partial \bar{\varrho}}{\partial t} + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\bar{\boldsymbol{J}} = 0 \end{equation}

follows directly from the charge conservation law.

2.2. Reduced energy–momentum conservation laws

Since the electromagnetic field $(\boldsymbol {E},\boldsymbol {B})$ is not split into equilibrium and perturbed components in the reduced Vlasov–Maxwell equations (2.6)–(2.8), the energy–momentum conservation laws derived for the reduced Vlasov equation (2.6) and the reduced Maxwell equations (2.7) and (2.8) (or (2.13) and (2.14)), with the source-free Maxwell equations (2.9), are direct consequences of the Noether theorem. Here, the reduced Noether equation (Brizard Reference Brizard2008) is first expressed as

(2.16)\begin{align} 0 &= \frac{\partial }{\partial t} \left[ \int_{\bar{\boldsymbol{p}}} \,\bar{f}\delta{\bar{\mathcal{S}}} - \delta\boldsymbol{A}\boldsymbol{\cdot}\frac{\bar{\mathbb{D}}}{4{\rm \pi} c} + \frac{\delta t}{8{\rm \pi}}\left( |\boldsymbol{E}|^{2} - |\boldsymbol{B}|^{2} \right) \right] \nonumber\\ &\quad +\boldsymbol{\nabla}\boldsymbol{\cdot}\left[ \int_{\bar{\boldsymbol{p}}}\, \bar{f}\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\delta \bar{\mathcal{S}} - \frac{1}{4{\rm \pi}}\left(\delta\varPhi\bar{\mathbb{D}} + \delta\boldsymbol{A}\boldsymbol{\times}\bar{\mathbb{H}}\right) + \frac{\delta\boldsymbol{x}}{8{\rm \pi}}\left( |\boldsymbol{E}|^{2} - |\boldsymbol{B}|^{2} \right) \right], \end{align}

where the Eulerian variations

(2.17)\begin{equation} \left. \begin{aligned} \delta\bar{\mathcal{S}} = \bar{\boldsymbol{p}}\,\boldsymbol{\cdot}\,\delta\boldsymbol{x} - \bar{H}\delta t, \\ \delta\varPhi = \boldsymbol{E}\,\boldsymbol{\cdot}\,\delta\boldsymbol{x} - c^{{-}1}\partial\delta\chi/\partial t, \\ \delta\boldsymbol{A} =\boldsymbol{E}c\delta t + \delta\boldsymbol{x}\boldsymbol{\times}\boldsymbol{B} + \boldsymbol{\nabla}\delta\chi, \end{aligned} \right\} \end{equation}

with the gauge term defined as $\delta \chi \equiv \varPhi c\delta t - \boldsymbol {A}\,\boldsymbol {\cdot }\,\delta \boldsymbol {x}$, are generated by the space–time virtual displacements $(\delta \boldsymbol {x}, \delta t)$. We note that the gauge term $\delta \chi$ appears naturally when the Euler variations $\delta \varPhi = -\delta t \partial \varPhi /\partial t - \delta \boldsymbol {x}\,\boldsymbol {\cdot }\,\boldsymbol {\nabla }\varPhi$ and $\delta \boldsymbol {A} = -\delta t\partial \boldsymbol {A}/\partial t - \delta \boldsymbol {x}\,\boldsymbol {\cdot }\,\boldsymbol {\nabla }\boldsymbol {A}$ are expressed in terms of the electric field $\boldsymbol {E} = -\boldsymbol {\nabla }\varPhi - c^{-1}\partial \boldsymbol {A}/\partial t$ and the magnetic field $\boldsymbol {B} = \boldsymbol {\nabla }\boldsymbol {\times }\boldsymbol {A}$.

We now remove the gauge-dependent terms by using the identity

(2.18)\begin{align} & -\frac{\partial }{\partial t}\left(\boldsymbol{\nabla}\delta\chi\boldsymbol{\cdot}\bar{\mathbb{D}}\right) + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\frac{\partial \delta\chi}{\partial t}\bar{\mathbb{D}} - c\boldsymbol{\nabla}\delta\chi\boldsymbol{\times}\bar{\mathbb{H}}\right)\nonumber\\ &\quad = \frac{\partial }{\partial t}\left(\delta\chi\boldsymbol{\nabla}\boldsymbol{\cdot}\bar{\mathbb{D}}\right) - \boldsymbol{\nabla}\boldsymbol{\cdot}\left[ \delta\chi\left( \frac{\partial \bar{\mathbb{D}}}{\partial t} - c\boldsymbol{\nabla}\boldsymbol{\times}\bar{\mathbb{H}}\right) \right], \end{align}

and, using the reduced Maxwell equations (2.13) and (2.14) and the gauge-independent term $\delta \bar{\mathcal {S}} + (e/c)\delta \chi = m\bar {\boldsymbol {v}}\,\boldsymbol {\cdot }\,\delta \boldsymbol {x} - \bar {K}\delta t$, the reduced Noether equation (2.16) yields the reduced energy–momentum conservation law

(2.19)\begin{equation} \frac{\partial }{\partial t}\left( \boldsymbol{\mathcal{P}}\boldsymbol{\cdot}\delta\boldsymbol{x} - \mathcal{E}\delta t \right) + \boldsymbol{\nabla}\boldsymbol{\cdot}\left( \boldsymbol{\mathsf{T}}\,\boldsymbol{\cdot}\,\delta\boldsymbol{x} - {\boldsymbol{S}} \delta t \right) = 0. \end{equation}

Here, the reduced energy–momentum densities

(2.20)\begin{gather} \mathcal{E} = \int_{\bar{\boldsymbol{p}}} \bar{f}\bar{K} + \frac{1}{4{\rm \pi}}\boldsymbol{E}\,\boldsymbol{\cdot}\,\bar{\mathbb{D}} - \frac{1}{8{\rm \pi}} \left( |\boldsymbol{E}|^{2} - |\boldsymbol{B}|^{2} \right), \end{gather}

(2.21)\begin{gather}\boldsymbol{\mathcal{P}} = \int_{\bar{\boldsymbol{p}}}\, \bar{f}\,\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}m\bar{\boldsymbol{v}} + \frac{\bar{\mathbb{D}}\boldsymbol{\times}\boldsymbol{B}}{4{\rm \pi} c} \end{gather}

both include reduced polarization effects, with (2.21) displaying the Minkowski form ($\bar {\mathbb {D}}\boldsymbol {\times }\boldsymbol {B}/4{\rm \pi} c$) for the reduced electromagnetic momentum density. The reduced energy-density flux

(2.22)\begin{equation} {\boldsymbol{S}}=\int_{\bar{\boldsymbol{p}}} \,\bar{f}\,\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\bar{K} +\frac{c}{4{\rm \pi}}\boldsymbol{E}\boldsymbol{\times}\bar{\mathbb{H}}, \end{equation}

on the other hand, displays the Abraham form ($\boldsymbol {E}\boldsymbol {\times } c\bar {\mathbb {H}}/4{\rm \pi}$) for the reduced Poynting flux, while the reduced stress tensor

(2.23)\begin{equation} \boldsymbol{\mathsf{T}} = \int_{\bar{\boldsymbol{p}}}\, \bar{f}\,\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}m\bar{\boldsymbol{v}} + \frac{\boldsymbol{I}}{4{\rm \pi}} \left[ \frac{1}{2} \left(|\boldsymbol{E}|^{2} - |\boldsymbol{B}|^{2}|\right) + \boldsymbol{B}\,\boldsymbol{\cdot}\,\bar{\mathbb{H}} \right] - \frac{1}{4{\rm \pi}} \left(\boldsymbol{B}\bar{\mathbb{H}} + \bar{\mathbb{D}}\boldsymbol{E}\right) \end{equation}

is composed of the reduced Reynolds stress tensor, which includes the ponderomotive velocity $\partial \bar {\varPsi }/\partial \bar {\boldsymbol {p}} = \textrm {d}\bar {\boldsymbol {x}}/\textrm {d}t - \bar {\boldsymbol {v}}$, and the reduced Maxwell stress tensor, which includes polarization and magnetization corrections. We immediately see that the reduced stress tensor (2.23), which can be expressed as

(2.24)\begin{align} \boldsymbol{\mathsf{T}} & = \frac{\boldsymbol{I}}{8{\rm \pi}} \left( |\boldsymbol{E}|^{2} + |\boldsymbol{B}|^{2} \right) -\frac{1}{4{\rm \pi}}\left(\boldsymbol{E}\boldsymbol{E} + \boldsymbol{B}\boldsymbol{B}\right) + \int_{\bar{\boldsymbol{p}}} \bar{f} \left[ m\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\,\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t} + {\boldsymbol{I}} \left(\boldsymbol{B}\boldsymbol{\cdot}\frac{\partial \bar{\varPsi}}{\partial \boldsymbol{B}}\right) \right] \nonumber\\ &\quad -\int_{\bar{\boldsymbol{p}}} \,\bar{f} \left( m\,\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\frac{\partial \bar{\varPsi}}{\partial \bar{\boldsymbol{p}}} + \boldsymbol{B}\frac{\partial \bar{\varPsi}}{\partial \boldsymbol{B}} - \frac{\partial \bar{\varPsi}}{\partial \boldsymbol{E}}\boldsymbol{E} \right), \end{align}

is manifestly asymmetric as a result of ponderomotive, polarization and magnetization effects appearing in the last line of (2.24).

The apparent asymmetry of the reduced stress tensor (2.24) implies that the azimuthal angular momentum may not be conserved:

(2.25)\begin{equation} \frac{\partial \mathcal{P}_{\varphi}}{\partial t} +\boldsymbol{\nabla}\,\boldsymbol{\cdot}\,{\boldsymbol{T}}_{\varphi} = \boldsymbol{\mathsf{T}}^{\top}:\boldsymbol{\nabla}(\partial\boldsymbol{x}/\partial\varphi) \equiv \hat{\boldsymbol{\mathsf{z}}}\boldsymbol{\cdot}\left(\int_{\boldsymbol{p}} \,\bar{f}\bar{\boldsymbol{N}} \right), \end{equation}

unless the reduced torque $\bar {\boldsymbol {N}}$ vanishes identically. In (2.25), $\mathcal {P}_{\varphi } \equiv \boldsymbol{\mathcal {P}}\,\boldsymbol {\cdot }\,\partial \boldsymbol {x}/\partial \varphi$ is the azimuthal angular momentum density, ${\boldsymbol {T}}_{\varphi } \equiv \boldsymbol{\mathsf{T}}\,\boldsymbol {\cdot }\,\partial \boldsymbol {x}/\partial \varphi$ is the azimuthal angular momentum-density flux, $\boldsymbol{\mathsf{T}}^{\top }$ denotes the transpose of $\boldsymbol{\mathsf{T}}$ and, since the dyadic tensor $\boldsymbol {\nabla }(\partial \boldsymbol {x}/\partial \varphi )$ is antisymmetric, the reduced torque $\bar {\boldsymbol {N}}$ is expressed as

(2.26)\begin{equation} \bar{\boldsymbol{N}} \equiv \frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\boldsymbol{\times} m\bar{\boldsymbol{v}}- \left(\bar{\boldsymbol{\rm \pi}}\boldsymbol{\times}\boldsymbol{E} + \bar{\boldsymbol{\mu}}\boldsymbol{\times}\boldsymbol{B}\right), \end{equation}

which includes contributions from the electric and magnetic torques. The required symmetry of the reduced stress tensor (2.24) must, therefore, introduce constraints on the reduced polarization and magnetization, which force the reduced torque (2.26) to vanish identically.

2.3. Guiding-centre Vlasov–Maxwell equations

The apparent asymmetry of the guiding-centre stress tensor was independently noted by Pfirsch & Morrison (Reference Pfirsch and Morrison1985) and Similon (Reference Similon1985), using different variational formulations. It was recently shown by Brizard & Tronci (Reference Brizard and Tronci2016) and Sugama et al. (Reference Sugama, Matsuoka, Satake and Kanno2016), however, that the guiding-centre stress tensor is indeed explicitly symmetric. Here, we use the variational formulation of the guiding-centre Vlasov–Maxwell model of Brizard & Tronci (Reference Brizard and Tronci2016) to show that the guiding-centre torque, derived from the generic reduced torque (2.26), vanishes identically.

In the work of Brizard & Tronci (Reference Brizard and Tronci2016), which considers the simplest case $\boldsymbol {E} = 0$, the guiding-centre canonical momentum is defined as $\bar {\boldsymbol {p}} = (e/c)\boldsymbol {A} + \bar {p}_{\|}\,\hat {{\boldsymbol{\mathsf{b}}}}$, which implies that $\bar {\boldsymbol {v}} = (\bar {p}_{\|}/m)\,\hat {{\boldsymbol{\mathsf{b}}}}$ and the guiding-centre electric and magnetic dipole moments are

(2.27)\begin{equation} \left. \begin{aligned} \bar{\boldsymbol{\rm \pi}} = (e\hat{{\boldsymbol{\mathsf{b}}}}/\varOmega)\boldsymbol{\times} \textrm{d}\bar{\boldsymbol{x}}/\textrm{d}t = (e\hat{{\boldsymbol{\mathsf{b}}}}/\varOmega)\boldsymbol{\times}\partial\bar{\varPsi}/\partial\bar{\boldsymbol{p}},\\ \bar{\boldsymbol{\mu}} ={-}\bar{\mu}\hat{{\boldsymbol{\mathsf{b}}}} + \bar{\boldsymbol{\rm \pi}}\boldsymbol{\times}(\bar{p}_{\|}\hat{{\boldsymbol{\mathsf{b}}}}/mc). \end{aligned} \right\} \end{equation}

Hence, using (2.2) and (2.27), we easily verify that the guiding-centre torque (2.26) vanishes:

(2.28)\begin{equation} \bar{\boldsymbol{N}} = \frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\boldsymbol{\times} m\bar{\boldsymbol{v}} - \left( -\bar{\mu}\hat{{\boldsymbol{\mathsf{b}}}} + \bar{\boldsymbol{\rm \pi}}\boldsymbol{\times}\frac{\bar{p}_{\|}\hat{{\boldsymbol{\mathsf{b}}}}}{mc} \right)\boldsymbol{\times}\boldsymbol{B} = \left({-}m\bar{\boldsymbol{v}} + \bar{p}_{\|}\hat{{\boldsymbol{\mathsf{b}}}}\right) \boldsymbol{\times}\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t} \equiv 0, \end{equation}

and the guiding-centre stress tensor (2.23) is symmetric (Sugama et al. Reference Sugama, Matsuoka, Satake and Kanno2016; Brizard & Tronci Reference Brizard and Tronci2016):

(2.29)\begin{align} \boldsymbol{\mathsf{T}}_\textrm{gc} &= \frac{1}{4{\rm \pi}}\left(\frac{\boldsymbol{I}}{2}|\boldsymbol{B}|^{2} - \boldsymbol{B}\boldsymbol{B} \right) + \int_{\bar{\boldsymbol{p}}} \,\bar{f} \left[m\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t}\left(\frac{\textrm{d}\bar{\boldsymbol{x}}}{\textrm{d}t} -\frac{\partial \bar{\varPsi}}{\partial \bar{\boldsymbol{p}}}\right) + {\boldsymbol{I}}\left(\boldsymbol{B}\boldsymbol{\cdot}\frac{\partial \bar{\varPsi}}{\partial \boldsymbol{B}}\right) -\boldsymbol{B}\frac{\partial \bar{\varPsi}}{\partial \boldsymbol{B}} \right] \nonumber\\ &=\frac{1}{4{\rm \pi}}\left(\frac{\boldsymbol{I}}{2}|\boldsymbol{B}|^{2} - \boldsymbol{B}\boldsymbol{B} \right) + \boldsymbol{\mathsf{P}}_\textrm{CGL} + \int_{\bar{\boldsymbol{p}}}\, \bar{f} \left[ \bar{p}_{\|} \left( \frac{\partial \bar{\varPsi}}{\partial \bar{\boldsymbol{p}}}\,\hat{{\boldsymbol{\mathsf{b}}}} + \hat{{\boldsymbol{\mathsf{b}}}}\, \frac{\partial \bar{\varPsi}}{\partial \bar{\boldsymbol{p}}}\right) \right], \end{align}

where $\boldsymbol{\mathsf{P}}_\textrm {CGL} = \int _{\bar {\boldsymbol {p}}}\,\bar {f}[(\bar {p}_{\|}^{2}/m)\hat {{\boldsymbol{\mathsf{b}}}}\hat {{\boldsymbol{\mathsf{b}}}} + \bar {\mu } B({\boldsymbol {I}} - \hat {{\boldsymbol{\mathsf{b}}}}\hat {{\boldsymbol{\mathsf{b}}}})]$ is the symmetric Chew–Goldberger–Low pressure tensor and the ponderomotive velocity $\partial \bar {\varPsi }/\partial \bar {\boldsymbol {p}}$, which is assumed to be perpendicular to $\boldsymbol {B}$, represents the magnetic-drift velocity.

3.1. Gauge-free gyrocentre models

In the gyrocentre Hamiltonian model of Burby & Brizard (Reference Burby and Brizard2019), presented here in the drift-kinetic limit considered by Hirvijoki et al. (Reference Hirvijoki, Burby, Pfefferlé and Brizard2020), we find the definitions

(3.3)\begin{equation} (\varPhi_{1\textrm{gy}}, \boldsymbol{A}_{1\textrm{gy}}, \boldsymbol{\varPi}_\textrm{gy}) = (\varPhi_{1}, \boldsymbol{A}_{1}, 0), \end{equation}

where the perturbation fields are evaluated at the gyrocentre position ${\boldsymbol {X}}$ and the gyrocentre kinetic energy is

(3.4)\begin{align} K_\textrm{gy} &= \frac{p_{\|}^{2}}{2m} + \mu\left( B_{0} + \epsilon B_{1\|} + \frac{\epsilon^{2}}{2}\,\frac{|\boldsymbol{B}_{1}|^{2}}{B_{0}} \right) - \boldsymbol{\rm \pi}_\textrm{gc}\,\boldsymbol{\cdot}\,\epsilon\left( \boldsymbol{E}_{1} + (p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc)\boldsymbol{\times}\boldsymbol{B}_{1}\right) \nonumber\\ &\quad -\epsilon^{2}\frac{mc^{2}}{2B_{0}^{2}}\left| \boldsymbol{E}_{1} +(p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc)\boldsymbol{\times}\boldsymbol{B}_{1}\right|^{2}, \end{align}

where $B_{1\|} \equiv \hat {{\boldsymbol{\mathsf{b}}}}_{0}\,\boldsymbol {\cdot }\,\boldsymbol {B}_{1}$ denotes the parallel component of the perturbed magnetic field $\boldsymbol {B}_{1}$. We note that the gauge-free model considered by Hirvijoki et al. (Reference Hirvijoki, Burby, Pfefferlé and Brizard2020) omits the guiding-centre electric dipole moment $\boldsymbol {{\rm \pi} }_\textrm {gc}$ in the gyrocentre kinetic energy (3.4) and, thus, the gyrocentre polarization and magnetization derived without this term are incomplete. We explicitly show in § 5, however, that this omission does not jeopardize the energy–momentum conservation laws.

Next, in the gyrocentre symplectic model of Brizard (Reference Brizard2020), we find

(3.5)\begin{equation} \left. \begin{aligned} (\varPhi_{1\textrm{gy}}, \boldsymbol{A}_{1\textrm{gy}}) = (\langle\varPhi_{1\textrm{gc}}\rangle, \langle\boldsymbol{A}_{1\textrm{gc}}\rangle), \\ \boldsymbol{\varPi}_\textrm{gy} = \epsilon\left(\langle\boldsymbol{E}_{1\textrm{gc}}\rangle + (p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc)\boldsymbol{\times}\langle\boldsymbol{B}_{1\textrm{gc}}\rangle\right)\boldsymbol{\times} e\hat{{\boldsymbol{\mathsf{b}}}}_{0}/\varOmega_{0}, \end{aligned} \right\} \end{equation}

where perturbation fields are evaluated at ${\boldsymbol {X}} + \boldsymbol {\rho }_{0}$, with $\langle \cdots \rangle$ denoting the standard gyroangle averaging (since the lowest-order guiding-centre gyroradius $\boldsymbol {\rho }_{0}$ depends on the gyrocentre gyroangle $\zeta$), and the gyrocentre kinetic energy is

(3.6)\begin{align} K_\textrm{gy} &= \frac{p_{\|}^{2}}{2m} + \mu\left( B_{0} + \epsilon\langle\langle B_{1\|\textrm{gc}}\rangle\rangle + \frac{\epsilon^{2}}{2}\,\frac{|\boldsymbol{B}_{1}|^{2}}{B_{0}} \right) \nonumber\\ &\quad + \epsilon^{2}\frac{mc^{2}}{2B_{0}^{2}}\left| \boldsymbol{E}_{1} + (p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc)\boldsymbol{\times}\boldsymbol{B}_{1}\right|^{2}. \end{align}

In (3.6), the finite-Larmor-radius effects are included only at first order in the perturbation expansion, with $\langle \langle \cdots \rangle \rangle$ denoting the gyrosurface averaging introduced by Porazik & Lin (Reference Porazik and Lin2011).

Previous symplectic gyrokinetic models considered either the parallel component $\langle A_{1\|\textrm {gc}}\rangle$ of the perturbed vector potential (Hahm et al. Reference Hahm, Lee and Brizard1988; Brizard Reference Brizard2017) or the inclusion of the perturbed $E\times B$ velocity (Wang & Hahm Reference Wang and Hahm2010a,Reference Wang and Hahmb; Leerink, Parra & Heikkinen Reference Leerink, Parra and Heikkinen2010), or both (Duthoit, Hahm & Wang Reference Duthoit, Hahm and Wang2014). In the present symplectic gyrokinetic model (3.5) and (3.6), the addition of the perturbed magnetic flutter momentum to the $E\times B$ momentum yields a covariant treatment of the electric dipole moment in the gyrocentre polarization and magnetization; see (3.23)–(3.26). In their guiding-centre treatment, Fan et al. (Reference Fan, Qin and Xiao2020) considered an extension of the Pfirsch & Morrison (Reference Pfirsch and Morrison1985) variational formulation by including higher-order guiding-centre corrections, where both electric and magnetic fields $(\boldsymbol {E}, \boldsymbol {B} = B\hat {{\boldsymbol{\mathsf{b}}}})$ are considered as variational fields.

3.2. Gyrocentre Euler–Lagrange equations

The gyrocentre Euler–Lagrange equations involving arbitrary variations in $({\boldsymbol {X}}, p_{\|}, J)$ are, respectively,

(3.7)\begin{gather} 0 = e\boldsymbol{E}_\textrm{gy}^{*} + \frac{e}{c}\dot{\boldsymbol{X}}\boldsymbol{\times}\boldsymbol{B}_\textrm{gy}^{*} - \dot{p}_{\|}\boldsymbol{\mathsf{b}}_\textrm{gy}^{*}, \end{gather}

(3.8)\begin{gather}0 = \dot{\boldsymbol{X}}\,\boldsymbol{\cdot}\,\boldsymbol{\mathsf{b}}_\textrm{gy}^{*} - \partial K_\textrm{gy}/\partial p_{\|}, \end{gather}

(3.9)\begin{gather}0 = \dot{\zeta} + \dot{\boldsymbol{X}}\,\boldsymbol{\cdot}\,\partial{\boldsymbol{P}}_\textrm{gy}/\partial J - \partial H_\textrm{gy}/\partial J, \end{gather}

where the effective gyrocentre electric field $\boldsymbol {E}_\textrm {gy}^{*}$ is defined as

(3.10)\begin{equation} e\boldsymbol{E}_\textrm{gy}^{*} \equiv{-}\boldsymbol{\nabla} H_\textrm{gy} - \frac{\partial {\boldsymbol{P}}_\textrm{gy}}{\partial t} = \epsilon e\boldsymbol{E}_{1\textrm{gy}} - \left( \frac{\partial \boldsymbol{\varPi}_\textrm{gy}}{\partial t} + \boldsymbol{\nabla} K_\textrm{gy} \right), \end{equation}

with $\boldsymbol {E}_{1\textrm {gy}} \equiv -\boldsymbol {\nabla }\varPhi _{1\textrm {gy}} - c^{-1}\partial \boldsymbol {A}_{1\textrm {gy}}/\partial t$. The effective gyrocentre magnetic field $\boldsymbol {B}_\textrm {gy}^{*}$ is defined as

(3.11)\begin{equation} \boldsymbol{B}_\textrm{gy}^{*} \equiv \boldsymbol{\nabla}\boldsymbol{\times}\left(\frac{c}{e}{\boldsymbol{P}}_\textrm{gy}\right) = \boldsymbol{B}_{0}^{*} + \epsilon\boldsymbol{B}_{1\textrm{gy}} + \boldsymbol{\nabla}\boldsymbol{\times}\left(\frac{c}{e}\boldsymbol{\varPi}_\textrm{gy}\right), \end{equation}

with $\boldsymbol {B}_{0}^{*} \equiv \boldsymbol {\nabla }\boldsymbol {\times }\boldsymbol {A}_{0}^{*}$ and $\boldsymbol {B}_{1\textrm {gy}} \equiv \boldsymbol {\nabla }\boldsymbol {\times }\boldsymbol {A}_{1\textrm {gy}}$, while

(3.12)\begin{equation} \boldsymbol{\mathsf{b}}^{*}_\textrm{gy} \equiv \partial{\boldsymbol{P}}_\textrm{gy}/\partial p_{\|} = \hat{{\boldsymbol{\mathsf{b}}}}_{0} + \partial\boldsymbol{\varPi}_\textrm{gy}/\partial p_{\|}. \end{equation}

We note that the effective gyrocentre electromagnetic fields satisfy the source-free Maxwell equations $\boldsymbol {\nabla }\,\boldsymbol {\cdot }\,\boldsymbol {B}_\textrm {gy}^{*} = 0$ and $\partial \boldsymbol {B}_\textrm {gy}^{*}/\partial t + c\boldsymbol {\nabla }\boldsymbol {\times }\boldsymbol {E}_\textrm {gy}^{*} = 0$.

The gyrocentre Euler–Lagrange equations (3.7) and (3.8) can also be written in Hamiltonian form as

(3.13)\begin{gather} \dot{\boldsymbol{X}} \equiv \left\{ {\boldsymbol{X}}, H_\textrm{gy}\right\}_\textrm{gy} = \boldsymbol{E}_\textrm{gy}^{*}\boldsymbol{\times}\frac{c\boldsymbol{\mathsf{b}}^{*}_\textrm{gy}}{B_{\|\textrm{gy}}^{**}} + \frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\frac{\boldsymbol{B}^{*}_\textrm{gy}}{B_{\|\textrm{gy}}^{**}}, \end{gather}

(3.14)\begin{gather}\dot{p}_{\|} \equiv \left\{ p_{\|}, H_\textrm{gy}\right\}_\textrm{gy} = e\boldsymbol{E}_\textrm{gy}^{*}\boldsymbol{\cdot}\frac{\boldsymbol{B}^{*}_\textrm{gy}}{B_{\|\textrm{gy}}^{**}}, \end{gather}

where $\{\,,\,\}_\textrm {gy}$ denotes the gyrocentre Poisson bracket and $B_{\|\textrm {gy}}^{**} \equiv \boldsymbol{\mathsf{b}}^{*}_\textrm {gy}\,\boldsymbol {\cdot }\,\boldsymbol {B}^{*}_\textrm {gy}$. We note that (3.13) and (3.14) satisfy the Euler–Lagrange identity:

(3.15)\begin{equation} \frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\dot{p}_{\|} = e\boldsymbol{E}_\textrm{gy}^{*}\,\boldsymbol{\cdot}\,\frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\frac{\boldsymbol{B}^{*}_\textrm{gy}}{B_{\|\textrm{gy}}^{**}} \equiv e\boldsymbol{E}_\textrm{gy}^{*}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}}, \end{equation}

which will be useful in our discussion of energy conservation. The gyrocentre equations (3.13) and (3.14) also satisfy the Liouville theorem:

(3.16)\begin{align} \frac{\partial B_{\|\textrm{gy}}^{**}}{\partial t} & = \frac{\partial \boldsymbol{\mathsf{b}}_\textrm{gy}^{*}}{\partial t}\,\boldsymbol{\cdot}\,\boldsymbol{B}_\textrm{gy}^{*} + \boldsymbol{\mathsf{b}}_\textrm{gy}^{*}\,\boldsymbol{\cdot}\,\frac{\partial \boldsymbol{B}_\textrm{gy}^{*}}{\partial t} = \frac{\partial }{\partial p_{\|}}\left(\frac{\partial {\boldsymbol{P}}_\textrm{gy}}{\partial t}\right)\boldsymbol{\cdot}\boldsymbol{B}_\textrm{gy}^{*} - \boldsymbol{\mathsf{b}}_\textrm{gy}^{*}\,\boldsymbol{\cdot}\,\boldsymbol{\nabla}\boldsymbol{\times}\left(c\boldsymbol{E}_\textrm{gy}^{*}\right) \nonumber\\ &={-}\frac{\partial }{\partial p_{\|}}\left[\boldsymbol{\nabla}\boldsymbol{\cdot}\left(H_\textrm{gy}\boldsymbol{B}^{*}_\textrm{gy}\right) + \dot{p}_{\|}\,B_{\|\textrm{gy}}^{**}\right] - \frac{\partial {\boldsymbol{P}}_\textrm{gy}}{\partial t}\,\boldsymbol{\cdot}\, \boldsymbol{\nabla}\boldsymbol{\times}\left(\frac{c\boldsymbol{\mathsf{b}}^{*}_\textrm{gy}}{e}\right) \nonumber\\ &\quad - \boldsymbol{\nabla}\boldsymbol{\cdot}\left[ B_{\|\textrm{gy}}^{**}\dot{\boldsymbol{X}} - \frac{\partial }{\partial p_{\|}}\left( H_\textrm{gy}\boldsymbol{B}^{*}_\textrm{gy}\right) \right] - \left( e\boldsymbol{E}_\textrm{gy}^{*} + \boldsymbol{\nabla} H_\textrm{gy}\right)\boldsymbol{\cdot} \boldsymbol{\nabla}\boldsymbol{\times}\left(\frac{c\boldsymbol{\mathsf{b}}^{*}_\textrm{gy}}{e}\right) \nonumber\\ & ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\left(B_{\|\textrm{gy}}^{**}\, \dot{\boldsymbol{X}}\right) - \frac{\partial }{\partial p_{\|}} \left( B_{\|\textrm{gy}}^{**}\,\dot{p}_{\|} \right), \end{align}

where we used (3.10).

3.3. Eulerian field variations of the gyrocentre Lagrangian

In the next section, we will need the Eulerian field variation of the gyrocentre Lagrangian (3.1) at a field point $\boldsymbol {x}$:

(3.17)\begin{equation} \delta L_\textrm{gy} = \left( \frac{e}{c}\epsilon\delta\boldsymbol{A}_{1\textrm{gy}} + \delta\boldsymbol{\varPi}_\textrm{gy}\right)\boldsymbol{\cdot}\dot{\boldsymbol{X}} - \left( e \epsilon\delta\varPhi_{1\textrm{gy}} + \delta K_\textrm{gy} \right), \end{equation}

where, in contrast to the works of Sugama et al. (Reference Sugama, Matsuoka, Nunami and Satake2021) and Fan et al. (Reference Fan, Qin and Xiao2020), the guiding-centre Lagrangian terms $(e/c)\boldsymbol {A}_{0}^{*}\,\boldsymbol {\cdot }\,\dot {\boldsymbol {X}} + J\dot {\zeta } - (p_{\|}^{2}/2m + \mu B_{0})$ are invariant in our gyrokinetic formalism. While Sugama et al. (Reference Sugama, Matsuoka, Nunami and Satake2021) considered the simplest guiding-centre representation (with $e\boldsymbol {A}_{0}^{*}/c = e\boldsymbol {A}_{0}/c + p_{\|}\hat {{\boldsymbol{\mathsf{b}}}}$), with field variations easily computed (e.g. $\delta \hat {{\boldsymbol{\mathsf{b}}}} = (\hat {{\boldsymbol{\mathsf{b}}}}\boldsymbol {\times }\delta \boldsymbol {B})\boldsymbol {\times }\hat {{\boldsymbol{\mathsf{b}}}}/B$), the higher-order guiding-centre model used by Fan et al. (Reference Fan, Qin and Xiao2020) requires complex expressions for the variations of the gyrogauge vector $\boldsymbol {R}_{0} = \boldsymbol {\nabla }\hat {\boldsymbol{\mathsf{1}}}\,\boldsymbol {\cdot }\,\hat {\boldsymbol{\mathsf{2}}}$, for example, in which the functional derivatives of all three unit vectors $(\hat {\boldsymbol{\mathsf{1}}},\hat {\boldsymbol{\mathsf{2}}}, \hat {{\boldsymbol{\mathsf{b}}}} = \hat {\boldsymbol{\mathsf{1}}}\boldsymbol {\times }\hat {\boldsymbol{\mathsf{2}}})$ need to be computed, although they are not explicitly calculated.

Here, the field variations are defined in terms of the generic functional derivatives

(3.18)\begin{equation} \left( \begin{matrix} \delta\varPsi_{1}({\boldsymbol{X}}) \\ \langle\delta\varPsi_{1}({\boldsymbol{X}} + \boldsymbol{\rho}_{0})\rangle \end{matrix} \right) \equiv \int_{\boldsymbol{x}} \delta\varPsi_{1}(\boldsymbol{x})\left(\begin{matrix} \delta^{3}({\boldsymbol{X}} - \boldsymbol{x}) \\ \left\langle\delta^{3}({\boldsymbol{X}} + \boldsymbol{\rho}_{0} - \boldsymbol{x})\right\rangle \end{matrix} \right), \end{equation}

where $\varPsi _{1}$ denotes an arbitrary component of the perturbed electromagnetic potentials or fields. We note that the second expression in (3.18) is valid only if the equilibrium (non-variational) magnetic field appears in the definition of the lowest-order gyroangle-dependent gyroradius $\boldsymbol {\rho }_{0}$. Hence, we find

(3.19)\begin{equation} \left(\frac{\delta \varPhi_{1\textrm{gy}}}{\delta \varPhi_{1}(\boldsymbol{x})}, \frac{\delta A^{i}_{1\textrm{gy}}}{\delta A^{j}_{1}(\boldsymbol{x})}, \frac{\delta E^{i}_{1\textrm{gy}}}{\delta E^{j}_{1}(\boldsymbol{x})}, \frac{\delta B^{i}_{1\textrm{gy}}}{\delta B^{j}_{1}(\boldsymbol{x})} \right) = \begin{cases} \left( \delta^{3}, \delta^{i}_{j}\delta^{3}, \delta^{i}_{j}\delta^{3}, \delta^{i}_{j}\delta^{3}\right) \\ \left( \langle\delta_\textrm{gc}^{3}\rangle, \delta^{i}_{j} \langle\delta_\textrm{gc}^{3}\rangle, \delta^{i}_{j}\langle\delta_\textrm{gc}^{3}\rangle, \delta^{i}_{j}\langle\delta_\textrm{gc}^{3}\rangle\right) \end{cases} \end{equation}

and

(3.20)\begin{equation} \left( \frac{\delta B_{1\|}({\boldsymbol{X}})}{\delta \boldsymbol{B}_{1}(\boldsymbol{x})}, \frac{\delta \langle\langle B_{1\|\textrm{gc}}\rangle\rangle}{\delta \boldsymbol{B}_{1}(\boldsymbol{x})}\right) = \left( \delta^{3}\hat{{\boldsymbol{\mathsf{b}}}}_{0}, \langle\langle\delta_\textrm{gc}^{3}\rangle\rangle\hat{{\boldsymbol{\mathsf{b}}}}_{0}\right), \end{equation}

with $\delta ^{3} \equiv \delta ^{3}({\boldsymbol {X}} - \boldsymbol {x})$ and $\delta _\textrm {gc}^{3} \equiv \delta ^{3}({\boldsymbol {X}} + \boldsymbol {\rho }_{0} - \boldsymbol {x})$ used in the gyrocentre models of Burby & Brizard (Reference Burby and Brizard2019) and Brizard (Reference Brizard2020), respectively, and $\delta ^{i}_{j}$ denotes the standard Kronecker delta.

In the gyrocentre model (3.4) of Burby & Brizard (Reference Burby and Brizard2019), we find

(3.21)\begin{gather} \epsilon^{{-}1} \frac{\delta K_\textrm{gy}}{\delta \boldsymbol{E}_{1}(\boldsymbol{x})} ={-}\delta^{3} \left[ \boldsymbol{\rm \pi}_\textrm{gc} + \epsilon\frac{mc^{2}}{B_{0}^{2}} \left(\boldsymbol{E}_{1} + \frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\boldsymbol{\times}\boldsymbol{B}_{1}\right) \right] \equiv{-}\delta^{3}\left(\boldsymbol{\rm \pi}_\textrm{gc} + \epsilon\boldsymbol{\rm \pi}_{2}\right), \end{gather}

(3.22)\begin{gather}\epsilon^{{-}1} \frac{\delta K_\textrm{gy}}{\delta \boldsymbol{B}_{1}(\boldsymbol{x})} = \delta^{3}\mu \left( \hat{{\boldsymbol{\mathsf{b}}}}_{0} + \epsilon\frac{\boldsymbol{B}_{1}}{B_{0}}\right) - \delta^{3}\left(\boldsymbol{\rm \pi}_\textrm{gc} + \epsilon\boldsymbol{\rm \pi}_{2}\right)\boldsymbol{\times}\frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}, \end{gather}

where the gyrocentre electric dipole moment $\boldsymbol {{\rm \pi} }_\textrm {gc} + \epsilon \boldsymbol {{\rm \pi} }_{2}$ includes the guiding-centre contribution $\boldsymbol {{\rm \pi} }_\textrm {gc}$ and its first-order gyrocentre correction $\boldsymbol {{\rm \pi} }_{2}$ (derived from the second-order gyrocentre Hamiltonian), while the intrinsic gyrocentre magnetic dipole moment $-\mu (\hat {{\boldsymbol{\mathsf{b}}}}_{0} + \epsilon \boldsymbol {B}_{1}/B_{0})$ is accompanied by the moving gyrocentre electric dipole moment contribution $(\boldsymbol {{\rm \pi} }_\textrm {gc} + \epsilon \boldsymbol {{\rm \pi} }_{2})\boldsymbol {\times } p_{\|}\hat {{\boldsymbol{\mathsf{b}}}}_{0}/mc$.

In the gyrocentre model (3.5) and (3.6) of Brizard (Reference Brizard2020), on the other hand, we find

(3.23)\begin{gather} \epsilon^{{-}1} \frac{\delta \boldsymbol{\varPi}_\textrm{gy}}{\delta \boldsymbol{E}_{1}(\boldsymbol{x})}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} = \langle\delta_\textrm{gc}^{3}\rangle\,\frac{e\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{\varOmega_{0}}\boldsymbol{\times}\dot{\boldsymbol{X}} \equiv \langle\delta_\textrm{gc}^{3}\rangle\boldsymbol{\rm \pi}_\textrm{gy}. \end{gather}

(3.24)\begin{gather}\epsilon^{{-}1} \frac{\delta \boldsymbol{\varPi}_\textrm{gy}}{\delta \boldsymbol{B}_{1}(\boldsymbol{x})}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} = \langle\delta_\textrm{gc}^{3}\rangle\left(\boldsymbol{\rm \pi}_\textrm{gy} \boldsymbol{\times}\frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\right), \end{gather}

and

(3.25)\begin{gather} \epsilon^{{-}1} \frac{\delta K_\textrm{gy}}{\delta \boldsymbol{E}_{1}(\boldsymbol{x})} = \epsilon\delta^{3}\frac{mc^{2}}{B_{0}^{2}} \left(\boldsymbol{E}_{1} + \frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\boldsymbol{\times}\boldsymbol{B}_{1}\right) \equiv \epsilon\delta^{3}\boldsymbol{\rm \pi}_{2}, \end{gather}

(3.26)\begin{gather}\epsilon^{{-}1} \frac{\delta K_\textrm{gy}}{\delta \boldsymbol{B}_{1}(\boldsymbol{x})} = \mu \left( \langle\langle\delta_\textrm{gc}^{3}\rangle\rangle\hat{{\boldsymbol{\mathsf{b}}}}_{0} + \epsilon\delta^{3}\boldsymbol{B}_{1}/B_{0}\right) + \epsilon\delta^{3}\left( \boldsymbol{\rm \pi}_{2}\boldsymbol{\times}\frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc} \right). \end{gather}

We note that the gyrocentre polarization and magnetization derived from (3.21) and (3.22) for the gyrokinetic model of Burby & Brizard (Reference Burby and Brizard2019) are explicitly truncated at first order in the perturbation amplitudes of the electric and magnetic fields $(\boldsymbol {E}_{1},\boldsymbol {B}_{1})$. Because the gyrocentre velocity (3.13) appears in the expressions (3.23) and (3.24) for the gyrokinetic model of Brizard (Reference Brizard2020), however, the corresponding gyrocentre polarization and magnetization contain contributions at higher orders in perturbation amplitude.

5.1. Gyrokinetic energy conservation law

Since the equilibrium magnetic field $\boldsymbol {B}_{0}$ is time-independent, the total energy associated with the gyrokinetic Vlasov–Maxwell equations (4.20)–(4.22) is conserved. We derive the energy conservation law from the gyrokinetic Noether equation by setting $\delta t \neq 0$ and $\delta \boldsymbol {x} = 0$ in (5.9), which yields the gyrokinetic energy conservation law:

(5.10)\begin{equation} \partial\mathcal{E}_\textrm{gy}/\partial t + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,{\boldsymbol{S}}_\textrm{gy} = 0, \end{equation}

where the gyrokinetic energy density is

(5.11)\begin{align} \mathcal{E}_\textrm{gy} & = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy} FK_\textrm{gy} + \frac{\epsilon}{4{\rm \pi}}\boldsymbol{E}_{1}\,\boldsymbol{\cdot}\,\mathbb{D}_\textrm{gy} - \frac{1}{8{\rm \pi}}\left( \epsilon^{2}|\boldsymbol{E}_{1}|^{2} - |\boldsymbol{B}|^{2} \right) \nonumber\\ & = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left[ \frac{p_{\|}^{2}}{2m} + \mu\left( B_{0} + \epsilon\langle\langle B_{1\|\textrm{gc}}\rangle\rangle + \frac{\epsilon^{2}}{2}\frac{|\boldsymbol{B}_{1}|^{2}}{B_{0}} \right) \right. \nonumber\\ &\quad \left.+\, \boldsymbol{E}_{1}(\boldsymbol{x})\boldsymbol{\cdot}\left(\frac{\delta \boldsymbol{\varPi}_\textrm{gy}}{\delta \boldsymbol{E}_{1}(\boldsymbol{x})} - \frac{\delta K_\textrm{gy}}{\delta \boldsymbol{E}_{1}(\boldsymbol{x})}\right)\right] + \frac{1}{8{\rm \pi}}\left( \epsilon^{2}|\boldsymbol{E}_{1}|^{2} + |\boldsymbol{B}|^{2} \right), \end{align}

while the gyrokinetic energy-density flux is

(5.12)\begin{equation} {\boldsymbol{S}}_\textrm{gy} = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\dot{\boldsymbol{X}}K_\textrm{gy} + \frac{c}{4{\rm \pi}}\epsilon\boldsymbol{E}_{1}\boldsymbol{\times}\mathbb{H}_\textrm{gy}, \end{equation}

where the polarization and magnetization $(\mathbb {P}_\textrm {gy}, \mathbb {M}_\textrm {gy})$ are defined in (4.13) and (4.14), with $\mathbb {H}_\textrm {gy}$ defined in (4.19). In addition, we note that the gyrokinetic polarization and magnetization $(\mathbb {P}_\textrm {gy}, \mathbb {M}_\textrm {gy})$ include the full gyrocentre velocity $\dot {\boldsymbol {X}}$ defined in (3.13), which is expressed in terms of the effective electric and magnetic fields (3.10)–(3.11). We also note that, as shown by Burby et al. (Reference Burby, Brizard, Morrison and Qin2015), the gyrokinetic Vlasov–Maxwell Hamiltonian functional is naturally derived from the gyrokinetic energy density (5.11).

The explicit proof of gyrokinetic energy conservation, which applies to both gauge-free gyrokinetic models (Burby & Brizard Reference Burby and Brizard2019; Brizard Reference Brizard2020) considered here, proceeds as follows. First, we begin with the partial time derivative of the gyrokinetic energy density (5.11):

(5.13)\begin{align} \frac{\partial \mathcal{E}_\textrm{gy}}{\partial t} &= \int_{\boldsymbol{P}} \left[ \frac{\partial (\mathcal{ J}_\textrm{gy}F)}{\partial t}K_\textrm{gy} + \mathcal{J}_\textrm{gy}F\left( \frac{\partial \boldsymbol{E}_{1}}{\partial t}\boldsymbol{\cdot}\frac{\partial K_\textrm{gy}}{\partial \boldsymbol{E}_{1}} + \frac{\partial \boldsymbol{B}_{1}}{\partial t}\boldsymbol{\cdot}\frac{\partial K_\textrm{gy}}{\partial \boldsymbol{B}_{1}} \right) \right] \nonumber\\ &\quad + \frac{\epsilon\boldsymbol{E}_{1}}{4{\rm \pi}}\boldsymbol{\cdot}\frac{\partial \mathbb{D}_\textrm{gy}}{\partial t} + \epsilon\frac{\partial \boldsymbol{E}_{1}}{\partial t}\,\boldsymbol{\cdot}\,\mathbb{P}_\textrm{gy} + \frac{\boldsymbol{B}}{4{\rm \pi}}\,\boldsymbol{\cdot}\,\epsilon\frac{\partial \boldsymbol{B}_{1}}{\partial t}, \end{align}

where we expanded the term $\partial K_\textrm {gy}(\boldsymbol {E}_{1},\boldsymbol {B}_{1})/\partial t$ and used the definition (4.19) for $\mathbb {D}_\textrm {gy}$. Using the phase-space divergence form (4.20) of the gyrokinetic Vlasov equation, the first term on the right-hand side can be expressed as

(5.14)\begin{equation} \int_{\boldsymbol{P}} \frac{\partial (\mathcal{J}_\textrm{gy}F)}{\partial t}K_\textrm{gy} ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\left(\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}FK_\textrm{gy}\dot{\boldsymbol{X}}\right) +\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left( \frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\dot{p}_{\|} + \dot{\boldsymbol{X}}\,\boldsymbol{\cdot}\,\boldsymbol{\nabla} K_\textrm{gy}\right), \end{equation}

while, using the definitions (4.13) and (4.14) of the gyrokinetic polarization and magnetization, the gyrokinetic kinetic terms in (5.13) can be expressed as

(5.15)\begin{equation} \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\frac{\partial K_\textrm{gy}}{\partial t} ={-}\epsilon\left(\frac{\partial \boldsymbol{E}_{1}}{\partial t}\boldsymbol{\cdot}\mathbb{P}_\textrm{gy} + \frac{\partial \boldsymbol{B}_{1}}{\partial t}\boldsymbol{\cdot}\mathbb{M}_\textrm{gy} \right) + \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left( \frac{\partial \boldsymbol{\varPi}_\textrm{gy}}{\partial t}\boldsymbol{\cdot}\dot{\boldsymbol{X}}\right). \end{equation}

By combining these expressions, (5.13) becomes

(5.16)\begin{align} \frac{\partial \mathcal{E}_\textrm{gy}}{\partial t} & ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\left(\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}FK_\textrm{gy}\,\dot{\boldsymbol{X}}\right) + \frac{\mathbb{H}_\textrm{gy}}{4{\rm \pi}}\,\boldsymbol{\cdot}\,\epsilon\frac{\partial \boldsymbol{B}_{1}}{\partial t} + \frac{\epsilon\boldsymbol{E}_{1}}{4{\rm \pi}}\boldsymbol{\cdot}\frac{\partial \mathbb{D}_\textrm{gy}}{\partial t} \nonumber\\ &\quad + \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\left[\frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\dot{p}_{\|} +\dot{\boldsymbol{X}}\boldsymbol{\cdot}\left( \boldsymbol{\nabla} K_\textrm{gy} + \frac{\partial \boldsymbol{\varPi}_\textrm{gy}}{\partial t}\right) \right], \end{align}

where we introduced the definition (4.19) for $\mathbb {H}_\textrm {gy}$. Next, we use Faraday's law (4.25) to write

(5.17)\begin{align} \frac{\mathbb{H}_\textrm{gy}}{4{\rm \pi}}\,\boldsymbol{\cdot}\,\epsilon\frac{\partial \boldsymbol{B}_{1}}{\partial t} ={-}\frac{c\mathbb{H}_\textrm{gy}}{4{\rm \pi}}\,\boldsymbol{\cdot}\,\boldsymbol{\nabla}\boldsymbol{\times}\epsilon\boldsymbol{E}_{1} ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\left( \frac{c}{4{\rm \pi}}\epsilon\boldsymbol{E}_{1}\boldsymbol{\times}\mathbb{H}_\textrm{gy}\right) - \frac{\epsilon\boldsymbol{E}_{1}}{4{\rm \pi}}\,\boldsymbol{\cdot}\, c\boldsymbol{\nabla}\boldsymbol{\times}\mathbb{H}_\textrm{gy}, \end{align}

so that (5.16) becomes

(5.18)\begin{align} \frac{\partial \mathcal{E}_\textrm{gy}}{\partial t} + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,{\boldsymbol{S}}_\textrm{gy} &={-}\frac{\epsilon\boldsymbol{E}_{1}}{4{\rm \pi}}\boldsymbol{\cdot}\left( c\boldsymbol{\nabla} \boldsymbol{\times}\mathbb{H}_\textrm{gy} -\frac{\partial \mathbb{D}_\textrm{gy}}{\partial t} \right) \nonumber\\ &\quad + \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\left[ \frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\dot{p}_{\|} + \dot{\boldsymbol{X}}\boldsymbol{\cdot}\left( \boldsymbol{\nabla} K_\textrm{gy} + \frac{\partial \boldsymbol{\varPi}_\textrm{gy}}{\partial t}\right) \right], \end{align}

where we reconstructed the gyrokinetic energy-density flux (5.12) on the left-hand side of (5.18). Lastly, we use the identity derived from (3.10):

(5.19)\begin{equation} \boldsymbol{\nabla} K_\textrm{gy} + \frac{\partial \boldsymbol{\varPi}_\textrm{gy}}{\partial t} = e\left( \epsilon\boldsymbol{E}_{1\textrm{gy}} - \boldsymbol{E}_\textrm{gy}^{*}\right), \end{equation}

and we use the macroscopic gyrokinetic Maxwell equation (4.22), with

(5.20)\begin{equation}{-}\frac{\epsilon\boldsymbol{E}_{1}}{4{\rm \pi}}\boldsymbol{\cdot}\left( c\boldsymbol{\nabla}\boldsymbol{\times}\mathbb{H}_\textrm{gy} - \frac{\partial \mathbb{D}_\textrm{gy}}{\partial t} \right) ={-}\epsilon\boldsymbol{E}_{1}\,\boldsymbol{\cdot}\,\boldsymbol{J}_\textrm{gy} ={-}\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\left(\epsilon e\boldsymbol{E}_{1\textrm{gy}}\boldsymbol{\cdot}\dot{\boldsymbol{X}} \right), \end{equation}

to obtain

(5.21)\begin{equation} \frac{\partial \mathcal{E}_\textrm{gy}}{\partial t} + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,{\boldsymbol{S}}_\textrm{gy} = \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\left( \frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\dot{p}_{\|} - \dot{\boldsymbol{X}}\,\boldsymbol{\cdot}\, e\boldsymbol{E}_\textrm{gy}^{*} \right). \end{equation}

Using the Euler–Lagrange identity (3.15), the right-hand side of (5.21) is shown to vanish identically and we readily recover the exact gyrokinetic energy conservation law.

5.2. Gyrokinetic Noether momentum equation

Because the equilibrium magnetic field $\boldsymbol {B}_{0}$ considered in standard gyrokinetic Vlasov–Maxwell theory is spatially non-uniform (i.e. it serves to magnetically confine charged particles in accordance with the guiding-centre approximation), a general gyrokinetic Vlasov–Maxwell momentum conservation law does not exist. Indeed, according to the Noether theorem, the gyrokinetic Vlasov–Maxwell momentum is conserved only in directions corresponding to symmetries of the equilibrium magnetic field. Before we derive the gyrokinetic angular momentum conservation law associated with an axisymmetric equilibrium magnetic field, we wish to show that the gyrokinetic Noether momentum equation, from which our exact angular momentum conservation law will be derived, is consistent with the gyrokinetic Vlasov–Maxwell equations (4.20)–(4.22).

We begin with the gyrokinetic Noether momentum equation derived by setting $\delta t = 0$ and $\delta \boldsymbol {x} \neq 0$ in (5.9):

(5.22)\begin{equation} \frac{\partial \boldsymbol{\mathcal{P}}^{*}_\textrm{gy}}{\partial t} + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\boldsymbol{\mathsf{T}}^{*}_\textrm{gy} = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left( \frac{e}{c}\boldsymbol{\nabla}\boldsymbol{A}_{0}^{*}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} + \nabla^{\prime}\boldsymbol{\varPi}_\textrm{gy}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} - \nabla^{\prime}K_\textrm{gy}\right) - \boldsymbol{\nabla}\boldsymbol{B}_{0}\boldsymbol{\cdot}\frac{\boldsymbol{B}}{4{\rm \pi}}, \end{equation}

where the gyrokinetic canonical momentum density is defined as

(5.23)\begin{equation} \boldsymbol{\mathcal{P}}^{*}_\textrm{gy} = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left( \frac{e}{c}\boldsymbol{A}_{0}^{*} + \boldsymbol{\varPi}_\textrm{gy}\right) + \frac{\mathbb{D}_\textrm{gy}}{4{\rm \pi} c}\boldsymbol{\times}\epsilon\boldsymbol{B}_{1} \end{equation}

and the gyrokinetic canonical stress tensor is defined as

(5.24)\begin{align} \boldsymbol{\mathsf{T}}^{*}_\textrm{gy} & = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F \dot{\boldsymbol{X}}\left( \frac{e}{c}\boldsymbol{A}_{0}^{*} + \boldsymbol{\varPi}_\textrm{gy}\right) - \frac{\epsilon}{4{\rm \pi}} \left( \mathbb{D}_\textrm{gy}\boldsymbol{E}_{1} + \boldsymbol{B}_{1}\mathbb{H}_\textrm{gy}\right) \nonumber\\ &\quad +{\boldsymbol{I}} \left[ \frac{1}{8{\rm \pi}} \left( \epsilon^{2}|\boldsymbol{E}_{1}|^{2} -|\boldsymbol{B}|^{2} \right) + \frac{\epsilon}{4{\rm \pi}}\boldsymbol{B}_{1}\,\boldsymbol{\cdot}\,\mathbb{H}_\textrm{gy}\right], \end{align}

where ${\boldsymbol {I}}$ denotes the identity matrix. We note that, while the gyrokinetic stress tensor (5.24) is manifestly not symmetric, the conservation of the gyrokinetic angular momentum will follow exactly from (5.22).

We would now like to show that (5.22) is an exact consequence of the gyrokinetic Vlasov–Maxwell equations (4.20)–(4.22). We begin with the partial time derivatives of the first two terms in the gyrokinetic canonical momentum density (5.23):

(5.25)\begin{align} \frac{\partial }{\partial t}\left( \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\frac{e}{c} \boldsymbol{A}_{0}^{*}\right) &={-}\boldsymbol{\nabla}\boldsymbol{\cdot} \left(\int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\dot{\boldsymbol{X}}\frac{e}{c}\boldsymbol{A}_{0}^{*}\right)\nonumber\\ &\quad + \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy} F \left( \frac{e}{c}\boldsymbol{\nabla}\boldsymbol{A}_{0}^{*}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} + \dot{p}_{\|}\, \hat{{\boldsymbol{\mathsf{b}}}}_{0} - \frac{e}{c}\dot{\boldsymbol{X}}\boldsymbol{\times}\boldsymbol{B}_{0}^{*}\right), \end{align}

(5.26)\begin{align} \frac{\partial }{\partial t}\left( \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\boldsymbol{\varPi}_\textrm{gy}\right) &={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\left(\int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\dot{\boldsymbol{X}}\boldsymbol{\varPi}_\textrm{gy} \right) + \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left( \boldsymbol{\nabla}\boldsymbol{\varPi}_\textrm{gy}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} -\boldsymbol{\nabla} K_\textrm{gy} \right) \nonumber\\ &\quad + \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left[ \dot{p}_{\|}\frac{\partial \boldsymbol{\varPi}_\textrm{gy}}{\partial p_{\|}} + e \left(\epsilon\boldsymbol{E}_{1\textrm{gy}} - \boldsymbol{E}_\textrm{gy}^{*}\right) - \dot{\boldsymbol{X}}\boldsymbol{\times}\boldsymbol{\nabla}\boldsymbol{\times}\boldsymbol{\varPi}_\textrm{gy} \right], \end{align}

where we used the phase-space divergence form (4.20) of the gyrokinetic Vlasov equation, followed by integrations by parts, and used (5.19) to write $\partial \boldsymbol {\varPi }_\textrm {gy}/\partial t$. By combining these two expressions, we obtain

(5.27)\begin{align} \frac{\partial }{\partial t}\left[ \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left(\frac{e}{c}\boldsymbol{A}_{0}^{*} + \boldsymbol{\varPi}_\textrm{gy}\right)\right] & ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\left[ \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\dot{\boldsymbol{X}}\left(\frac{e}{c}\boldsymbol{A}_{0}^{*} + \boldsymbol{\varPi}_\textrm{gy}\right)\right]\nonumber\\ &\quad + \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left[ \left(\frac{e}{c}\boldsymbol{\nabla}\boldsymbol{A}_{0}^{*} + \boldsymbol{\nabla}\boldsymbol{\varPi}_\textrm{gy}\right)\boldsymbol{\cdot}\dot{\boldsymbol{X}} - \boldsymbol{\nabla} K_\textrm{gy} \right] \nonumber\\ &\quad +\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left[ \epsilon \left( e\boldsymbol{E}_{1\textrm{gy}} + \frac{e}{c}\dot{\boldsymbol{X}}\boldsymbol{\times}\boldsymbol{B}_{1\textrm{gy}} \right) \right] \nonumber\\ &\quad + \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left[ \dot{p}_{\|}\boldsymbol{\mathsf{b}}_\textrm{gy}^{*} - \left( e\boldsymbol{E}_\textrm{gy}^{*} + \frac{e}{c}\dot{\boldsymbol{X}}\boldsymbol{\times}\boldsymbol{B}_\textrm{gy}^{*} \right) \right], \end{align}

where the last line vanishes as a result of the gyrocentre Euler–Lagrange equation (3.7). Next, we take the partial time derivative of the third term in (5.23):

(5.28)\begin{align} \frac{\partial }{\partial t}\left( \frac{\mathbb{D}_\textrm{gy}}{4{\rm \pi} c}\boldsymbol{\times}\epsilon\boldsymbol{B}_{1} \right) & = \frac{1}{4{\rm \pi}}\left( \frac{1}{c}\frac{\partial \mathbb{D}_\textrm{gy}}{\partial t}\boldsymbol{\times}\epsilon\boldsymbol{B}_{1} + \mathbb{D}_\textrm{gy}\boldsymbol{\times}\frac{\epsilon}{c}\frac{\partial \boldsymbol{B}_{1}}{\partial t}\right) \nonumber\\ &= \left(\boldsymbol{\nabla}\boldsymbol{\times}\mathbb{H}_\textrm{gy}\right)\boldsymbol{\times}\frac{\epsilon\boldsymbol{B}_{1}}{4{\rm \pi}} - \frac{\mathbb{D}_\textrm{gy}}{4{\rm \pi}}\boldsymbol{\times}\left(\boldsymbol{\nabla}\boldsymbol{\times}\epsilon \boldsymbol{E}_{1} \right) - \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left( \frac{e}{c}\dot{\boldsymbol{X}}\boldsymbol{\times} \epsilon\boldsymbol{B}_{1\textrm{gy}}\right) \nonumber\\ &=\boldsymbol{\nabla}\boldsymbol{\cdot}\left[ \frac{\epsilon}{4{\rm \pi}} \left(\boldsymbol{B}_{1}\mathbb{H}_\textrm{gy} + \mathbb{D}_\textrm{gy}\boldsymbol{E}_{1}\right) - \frac{\boldsymbol{I}}{4{\rm \pi}} \left( \epsilon\boldsymbol{B}_{1}\,\boldsymbol{\cdot}\,\mathbb{H}_\textrm{gy} + \frac{\epsilon^{2}}{2}|\boldsymbol{E}_{1}|^{2} - \frac{1}{2} |\boldsymbol{B}|^{2} \right) \right] \nonumber\\ &\quad -\epsilon\left(\boldsymbol{\nabla}\boldsymbol{E}_{1}\,\boldsymbol{\cdot}\,\mathbb{P}_\textrm{gy} +\boldsymbol{\nabla}\boldsymbol{B}_{1}\,\boldsymbol{\cdot}\,\mathbb{M}_\textrm{gy}\right) -\boldsymbol{\nabla}\boldsymbol{B}_{0}\boldsymbol{\cdot}\frac{\boldsymbol{B}}{4{\rm \pi}} \nonumber\\ &\quad - \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \epsilon\left( e\boldsymbol{E}_{1\textrm{gy}} + \frac{e}{c}\dot{\boldsymbol{X}}\boldsymbol{\times} \boldsymbol{B}_{1\textrm{gy}}\right). \end{align}

When we combine (5.27) and (5.28), we obtain

(5.29)\begin{align} \frac{\partial \boldsymbol{\mathcal{P}}^{*}_\textrm{gy}}{\partial t} + \boldsymbol{\nabla}\,\boldsymbol{\cdot}\,\boldsymbol{\mathsf{T}}^{*}_\textrm{gy} &= \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left[ \left(\frac{e}{c}\boldsymbol{\nabla}\boldsymbol{A}_{0}^{*} + \boldsymbol{\nabla}\boldsymbol{\varPi}_\textrm{gy}\right)\boldsymbol{\cdot}\dot{\boldsymbol{X}} - \boldsymbol{\nabla} K_\textrm{gy} \right] \nonumber\\ &\quad -\epsilon\left(\boldsymbol{\nabla}\boldsymbol{E}_{1}\,\boldsymbol{\cdot}\,\mathbb{P}_\textrm{gy} +\boldsymbol{\nabla}\boldsymbol{B}_{1}\,\boldsymbol{\cdot}\,\mathbb{M}_\textrm{gy}\right) -\boldsymbol{\nabla}\boldsymbol{B}_{0}\boldsymbol{\cdot}\frac{\boldsymbol{B}}{4{\rm \pi}}, \end{align}

where

(5.30)\begin{align} -\epsilon\left(\boldsymbol{\nabla}\boldsymbol{E}_{1}\,\boldsymbol{\cdot}\,\mathbb{P}_\textrm{gy} +\boldsymbol{\nabla}\boldsymbol{B}_{1}\,\boldsymbol{\cdot}\,\mathbb{M}_\textrm{gy}\right) &={-}\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left(\boldsymbol{\nabla}\boldsymbol{E}_{1}\boldsymbol{\cdot}\frac{\delta \boldsymbol{\varPi}_\textrm{gy}}{\delta \boldsymbol{E}_{1}} + \boldsymbol{\nabla}\boldsymbol{B}_{1}\boldsymbol{\cdot} \frac{\delta \boldsymbol{\varPi}_\textrm{gy}}{\delta \boldsymbol{B}_{1}} \right)\boldsymbol{\cdot}\dot{\boldsymbol{X}} \nonumber\\ &\quad +\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left(\boldsymbol{\nabla}\boldsymbol{E}_{1} \boldsymbol{\cdot}\frac{\delta K_\textrm{gy}}{\delta \boldsymbol{E}_{1}} + \boldsymbol{\nabla}\boldsymbol{B}_{1}\boldsymbol{\cdot}\frac{\delta K_\textrm{gy}}{\delta \boldsymbol{B}_{1}} \right) \nonumber\\ &\equiv{-}\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F \left[ \left(\boldsymbol{\nabla}\boldsymbol{\varPi}_\textrm{gy} - \nabla^{\prime}\boldsymbol{\varPi}_\textrm{gy}\right)\boldsymbol{\cdot}\dot{\boldsymbol{X}} - \left(\boldsymbol{\nabla} K_\textrm{gy} - \nabla^{\prime}K_\textrm{gy}\right)\right]. \end{align}

By inserting these terms in (5.29), we recover the gyrokinetic Noether momentum equation (5.22).

We note that, while the gyrokinetic Noether momentum equation (5.22) is not a gyrokinetic conservation law, it can be used directly to obtain a gyrokinetic momentum transport equation (e.g. in the parallel direction) by taking its projection in the desired direction. For example, the gyrokinetic canonical parallel-momentum transport equation is expressed as

(5.31)\begin{align} \frac{\partial \mathcal{P}^{*}_{\textrm{gy}\|}}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\boldsymbol{\mathsf{T}}^{*}_\textrm{gy}\boldsymbol{\cdot}\hat{{\boldsymbol{\mathsf{b}}}}_{0}\right) &= \boldsymbol{\mathsf{T}}^{*\top}_\textrm{gy}:\boldsymbol{\nabla}\hat{{\boldsymbol{\mathsf{b}}}}_{0} - \hat{{\boldsymbol{\mathsf{b}}}}_{0}\boldsymbol{\cdot}\boldsymbol{\nabla}\boldsymbol{B}_{0}\boldsymbol{\cdot}\boldsymbol{B}/4{\rm \pi}\nonumber\\ &\quad + \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy} F\,\hat{{\boldsymbol{\mathsf{b}}}}_{0}\boldsymbol{\cdot}\left( \frac{e}{c}\boldsymbol{\nabla}\boldsymbol{A}_{0}^{*}\boldsymbol{\cdot}\dot{\boldsymbol{X}} + \nabla^{\prime}\boldsymbol{\varPi}_\textrm{gy}\boldsymbol{\cdot}\dot{\boldsymbol{X}} - \nabla^{\prime}K_\textrm{gy}\right), \end{align}

where $\mathcal {P}^{*}_{\textrm {gy}\|} \equiv \boldsymbol{\mathcal {P}}^{*}_\textrm {gy}\,\boldsymbol {\cdot }\,\hat {{\boldsymbol{\mathsf{b}}}}_{0}$ and $\boldsymbol{\mathsf{T}}_{\textrm {gy}}^{*\top }$ denotes the transpose of the gyrokinetic stress tensor (5.24). The gyrokinetic canonical parallel-momentum transport equation (5.31) can be transformed into a simpler form as the $p_{\|}$-moment of the gyrokinetic Vlasov equation:

(5.32)\begin{equation} \frac{\partial }{\partial t}\left(\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\,p_{\|} \right) +\boldsymbol{\nabla}\boldsymbol{\cdot}\left(\int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\dot{\boldsymbol{X}}p_{\|} \right) = \int_{\boldsymbol{P}}\mathcal{J}_\textrm{gy}F\dot{p}_{\|}, \end{equation}

where the gyrocentre parallel force $\dot {p}_{\|}$ is defined by (3.14). See Brizard & Tronko (Reference Brizard and Tronko2011) for the explicit transformation from (5.31) to (5.32) for the case of the gyrokinetic Vlasov–Poisson equations. We note that the parallel contraction of the gyrokinetic stress tensor $\boldsymbol{\mathsf{T}}^{*}_\textrm {gy}\,\boldsymbol {\cdot }\,\hat {{\boldsymbol{\mathsf{b}}}}_{0}$ on the left-hand side of (5.31) contains the gyrokinetic Maxwell stress tensor term $-\mathbb {D}_\textrm {gy}\epsilon E_{1\|}/4{\rm \pi}$, which plays a central role in the electrostatic gyrokinetic Vlasov–Poisson model of McDevitt et al. (Reference McDevitt, Diamond, Gürcan and Hahm2009) in discussing toroidal rotation driven by the gyrocentre polarization $\mathbb {P}_\textrm {gy}$. In particular, McDevitt et al. (Reference McDevitt, Diamond, Gürcan and Hahm2009) show how this polarization contribution can be retrieved from a perturbation expansion (up to fourth order) of the right-hand side of (5.32) through a $\delta F$-decomposition of the gyrocentre Vlasov distribution. Our gyrokinetic canonical parallel-momentum transport equation (5.31), in contrast, explicitly exhibits the complete gyrocentre polarization and magnetization effects in a full-$F$ gyrokinetic Vlasov–Maxwell theory.

5.3. Gyrokinetic angular momentum conservation law

Assuming now that the equilibrium magnetic field $\boldsymbol {B}_{0}$ is axisymmetric (i.e. $\partial \boldsymbol {B}_{0}/\partial \varphi = \hat {\boldsymbol{\mathsf{z}}}\boldsymbol {\times }\boldsymbol {B}_{0}$), we derive the gyrokinetic canonical angular momentum conservation law by taking the scalar product of (5.22) with $\partial \boldsymbol {x}/\partial \varphi$ (i.e. $\delta \boldsymbol {x} = \delta \varphi \partial \boldsymbol {x}/\partial \varphi$), where the toroidal angle $\varphi$ is associated with rotations about the $z$ axis. Hence, the toroidal canonical angular momentum density $\mathcal {P}_{\textrm {gy}\varphi }^{*} \equiv \boldsymbol{\mathcal {P}}_{\textrm {gy}}^{*}\,\boldsymbol {\cdot }\,\partial \boldsymbol {x}/\partial \varphi$ satisfies the Noether canonical angular momentum equation

(5.33)\begin{align} \frac{\partial \mathcal{P}_{\textrm{gy}\varphi}^{*}}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\boldsymbol{\mathsf{T}}_{\textrm{gy}}^{*}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) &= \boldsymbol{\mathsf{T}}_{\textrm{gy}}^{*\top}\boldsymbol{:}\boldsymbol{\nabla}\left(\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) - \frac{\partial \boldsymbol{B}_{0}}{\partial \varphi}\boldsymbol{\cdot}\frac{\boldsymbol{B}}{4{\rm \pi}} \nonumber\\ &\quad + \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left( \frac{e}{c}\,\frac{\partial \boldsymbol{A}_{0}^{*}}{\partial \varphi}\,\boldsymbol{\cdot}\,\dot{\boldsymbol{X}} + \frac{\partial^{\prime}\boldsymbol{\varPi}_\textrm{gy}}{\partial\varphi}\boldsymbol{\cdot}\dot{\boldsymbol{X}} - \frac{\partial^{\prime}K_\textrm{gy}}{\partial\varphi}\right). \end{align}

Under the assumption that the equilibrium magnetic field is axisymmetric, we have $\partial B_{0}/\partial \varphi \equiv 0$ and we will use the identity $\partial \hat {{\boldsymbol{\mathsf{b}}}}_{0}/\partial \varphi \equiv \hat {\boldsymbol{\mathsf{z}}}\boldsymbol {\times }\hat {{\boldsymbol{\mathsf{b}}}}_{0}$, so that $\boldsymbol {B}\,\boldsymbol {\cdot }\,\partial \boldsymbol {B}_{0}/\partial \varphi = \epsilon \boldsymbol {B}_{1}\boldsymbol {\cdot }(\hat {\boldsymbol{\mathsf{z}}}\boldsymbol {\times }\boldsymbol {B}_{0})$.

Instead of merely assuming that the right-hand side of (5.33) is zero, we now systematically show how the various terms do cancel each other out to yield an exact conservation law. Before we begin, however, we note that the first term vanishes identically if the gyrokinetic stress tensor (5.24) is symmetric (i.e. $\boldsymbol{\mathsf{T}}_{\textrm {gy}}^{*\top } = \boldsymbol{\mathsf{T}}_{\textrm {gy}}^{*}$), which is expected (and required) when there is no separation between dynamical fields and equilibrium fields, e.g. in guiding-centre Vlasov–Maxwell theory (Brizard & Tronci Reference Brizard and Tronci2016). In the present case, however, the asymmetry of the gyrokinetic stress tensor (5.24) is necessary in order to cancel the additional terms on the right-hand side of (5.33).

We now proceed with the proof that the right-hand side of (5.33) is zero for the gauge-free model of Burby & Brizard (Reference Burby and Brizard2019), where $\boldsymbol {\varPi }_\textrm {gy} \equiv 0$, and present the results for the gauge-free model of Brizard (Reference Brizard2020). First, we note that since the dyadic tensor $\boldsymbol {\nabla }(\partial \boldsymbol {x}/\partial \varphi ) = \hat {R}\,\hat {\varphi } - \hat {\varphi }\,\hat {R}$ is antisymmetric (where $R \equiv |\partial \boldsymbol {x}/\partial \varphi |$), only the antisymmetric part of $\boldsymbol{\mathsf{T}}_{\textrm {gy}}^{*\top }$ contributes in the first term of (5.33):

(5.34)\begin{align} \boldsymbol{\mathsf{T}}_{\textrm{gy}}^{*\top}\boldsymbol{:}\boldsymbol{\nabla}\left(\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) & = \hat{\boldsymbol{\mathsf{z}}}\boldsymbol{\cdot}\left[ \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left(\dot{\boldsymbol{X}}\boldsymbol{\times}\frac{e}{c}\boldsymbol{A}_{0}^{*}\right) -\frac{\epsilon}{4{\rm \pi}} \left(\mathbb{D}_\textrm{gy}\boldsymbol{\times} \boldsymbol{E}_{1} +\boldsymbol{B}_{1}\boldsymbol{\times}\mathbb{H}_\textrm{gy}\right) \right] \nonumber\\ & = \hat{\boldsymbol{\mathsf{z}}}\boldsymbol{\cdot}\left[ \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left(\dot{\boldsymbol{X}}\boldsymbol{\times}\frac{e}{c}\boldsymbol{A}_{0}^{*}\right) +\epsilon\boldsymbol{E}_{1}\boldsymbol{\times}\mathbb{P}_\textrm{gy} + \epsilon\boldsymbol{B}_{1}\boldsymbol{\times}\mathbb{M}_\textrm{gy}\right] \nonumber\\ &\quad -\frac{\hat{\boldsymbol{\mathsf{z}}}}{4{\rm \pi}}\boldsymbol{\cdot}(\epsilon\boldsymbol{B}_{1}\boldsymbol{\times}\boldsymbol{B}_{0}), \end{align}

where we used the dyadic identities ${\boldsymbol {I}}:\boldsymbol {\nabla }(\partial \boldsymbol {x}/\partial \varphi ) = \boldsymbol {\nabla }\boldsymbol {\cdot }(\partial \boldsymbol {x}/\partial \varphi ) = 0$ and ${\boldsymbol {VW}}: \boldsymbol {\nabla }(\partial \boldsymbol {x}/\partial \varphi ) \equiv \hat {\boldsymbol{\mathsf{z}}}\boldsymbol {\cdot }({W}\boldsymbol {\times }{\boldsymbol {V}})$, which holds for an arbitrary pair of vectors $({\boldsymbol {V}},{\boldsymbol {W}})$. Next, the last term is

(5.35)\begin{equation} \frac{\partial^{\prime}K_\textrm{gy}}{\partial\varphi} = \hat{\boldsymbol{\mathsf{z}}}\boldsymbol{\cdot} \left[\epsilon\mu\hat{{\boldsymbol{\mathsf{b}}}}_{0}\boldsymbol{\times}\langle\langle \boldsymbol{B}_{1\textrm{gc}}\rangle\rangle + \frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\boldsymbol{\times}\left(\epsilon\boldsymbol{B}_{1}\boldsymbol{\times}\boldsymbol{\rm \pi}_\textrm{gy} \right) + \boldsymbol{\rm \pi}_\textrm{gy}\boldsymbol{\times}\epsilon\left( \boldsymbol{E}_{1} + \frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\boldsymbol{\times}\boldsymbol{B}_{1}\right) \right], \end{equation}

where $\boldsymbol {{\rm \pi} }_\textrm {gy} \equiv \boldsymbol {{\rm \pi} }_\textrm {gc} + \epsilon \boldsymbol {{\rm \pi} }_{2}$. Lastly, we write $\partial \boldsymbol {A}_{0}^{*}/\partial \varphi = \hat {\boldsymbol{\mathsf{z}}}\boldsymbol {\times }\boldsymbol {A}_{0}^{*}$ and, after some cancellations, (5.33) becomes

(5.36)\begin{equation} \frac{\partial \mathcal{P}_{\textrm{gy}\varphi}^{*}}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\boldsymbol{\mathsf{T}}_{\textrm{gy}}^{*}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) = \hat{\boldsymbol{\mathsf{z}}}\,\boldsymbol{\cdot}\,\epsilon\left(\boldsymbol{E}_{1}\boldsymbol{\times}\mathbb{P}_\textrm{gy} +\boldsymbol{B}_{1}\boldsymbol{\times}\mathbb{M}_\textrm{gy}\right) - \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\frac{\partial^{\prime}K_\textrm{gy}}{\partial\varphi}, \end{equation}

where

(5.37)\begin{align} \hat{\boldsymbol{\mathsf{z}}}\,\boldsymbol{\cdot}\,\epsilon\left(\boldsymbol{E}_{1}\boldsymbol{\times}\mathbb{P}_\textrm{gy} +\boldsymbol{B}_{1}\boldsymbol{\times}\mathbb{M}_\textrm{gy}\right) & = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\hat{\boldsymbol{\mathsf{z}}}\,\boldsymbol{\cdot}\,\epsilon\left( \boldsymbol{E}_{1}\boldsymbol{\times}\boldsymbol{\rm \pi}_\textrm{gy} - \mu \langle\langle\boldsymbol{B}_{1\textrm{gc}}\rangle\rangle\boldsymbol{\times}\hat{{\boldsymbol{\mathsf{b}}}}_{0} \right) \nonumber\\ &\quad + \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\hat{\boldsymbol{\mathsf{z}}}\boldsymbol{\cdot}\left[ \boldsymbol{B}_{1}\boldsymbol{\times}\left(\boldsymbol{\rm \pi}_\textrm{gy}\boldsymbol{\times}\frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\right) \right]. \end{align}

Upon further cancellations, (5.36) becomes

(5.38)\begin{equation} \frac{\partial \mathcal{P}_{\textrm{gy}\varphi}^{*}}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\boldsymbol{\mathsf{T}}_{\textrm{gy}}^{*}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\hat{\boldsymbol{\mathsf{z}}}\,\boldsymbol{\cdot}\,{\boldsymbol{N}}, \end{equation}

where the gyrocentre torque

(5.39)\begin{equation} {\boldsymbol{N}} \equiv \epsilon\frac{p_{\|}}{mc}\left[ \boldsymbol{B}_{1}\boldsymbol{\times}\left(\boldsymbol{\rm \pi}_\textrm{gy}\boldsymbol{\times}\hat{{\boldsymbol{\mathsf{b}}}}_{0}\right) + \boldsymbol{\rm \pi}_\textrm{gy}\boldsymbol{\times}\left(\hat{{\boldsymbol{\mathsf{b}}}}_{0}\boldsymbol{\times} \boldsymbol{B}_{1}\right) + \hat{{\boldsymbol{\mathsf{b}}}}_{0}\boldsymbol{\times}\left( \boldsymbol{B}_{1}\boldsymbol{\times}\boldsymbol{\rm \pi}_\textrm{gy}\right)\right] \equiv 0 \end{equation}

vanishes according to the Jacobi identity

(5.40)\begin{equation} {\boldsymbol{U}}\boldsymbol{\times}({\boldsymbol{V}}\boldsymbol{\times}{\boldsymbol{W}}) + { \boldsymbol{V}}\boldsymbol{\times}({\boldsymbol{W}}\boldsymbol{\times}{\boldsymbol{U}}) + {\boldsymbol{W}}\boldsymbol{\times}({\boldsymbol{U}}\boldsymbol{\times}{\boldsymbol{V}}) \equiv 0 \end{equation}

for the double vector product of any three arbitrary vector fields $({\boldsymbol {U}}, {\boldsymbol {V}}, {\boldsymbol {W}})$. For the gauge-free model of Brizard (Reference Brizard2020), the gyrocentre torque

(5.41)\begin{equation} {\boldsymbol{N}} \equiv \sum_{i = 1}^{3}\left[{\boldsymbol{U}}_{i}\boldsymbol{\times}({\boldsymbol{V}}_{i}\boldsymbol{\times}{ \boldsymbol{W}}_{i}) +{\boldsymbol{V}}_{i}\boldsymbol{\times}({\boldsymbol{W}}_{i}\boldsymbol{\times}{\boldsymbol{U}}_{i}) + {\boldsymbol{W}}_{i}\boldsymbol{\times}({\boldsymbol{U}}_{i}\boldsymbol{\times}{\boldsymbol{V}}_{i}) \right] \equiv 0 \end{equation}

also vanishes as a result of the Jacobi vector identity (5.40), where

(5.42)\begin{equation} \left. \begin{aligned} ({\boldsymbol{U}}_{1},{\boldsymbol{V}}_{1},{\boldsymbol{W}}_{1}) = \left(\epsilon\langle\boldsymbol{B}_{1\textrm{gc}}\rangle, \boldsymbol{\rm \pi}_\textrm{gy}, p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc\right), \\ ({\boldsymbol{U}}_{2},{\boldsymbol{V}}_{2},{\boldsymbol{W}}_{2}) = \left(\epsilon\boldsymbol{B}_{1}, -\epsilon\boldsymbol{\rm \pi}_{2}, p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc\right), \\ ({\boldsymbol{U}}_{3},{\boldsymbol{V}}_{3},{\boldsymbol{W}}_{3}) = \left(\dot{\boldsymbol{X}}, \epsilon\langle\boldsymbol{E}_{1\textrm{gc}}\rangle + p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}/mc\boldsymbol{\times} \epsilon\langle\boldsymbol{B}_{1\textrm{gc}}\rangle,e\hat{{\boldsymbol{\mathsf{b}}}}_{0}/\varOmega_{0}\right). \end{aligned} \right\} \end{equation}

5.4. Gyrokinetic angular momentum conservation in axisymmetric tokamak plasmas

Hence, we have explicitly proved that the gyrokinetic canonical angular momentum conservation law

(5.43)\begin{equation} \frac{\partial \mathcal{P}_{\textrm{gy}\varphi}^{*}}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\boldsymbol{\mathsf{T}}_{\textrm{gy}}^{*}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) = 0 \end{equation}

follows exactly from the gyrokinetic Vlasiov–Maxwell equations. We now evaluate this equation in axisymmetric tokamak geometry, in which the tokamak magnetic field is $\boldsymbol {B}_{0} = B_{0\varphi }(\psi )\boldsymbol {\nabla }\varphi + \boldsymbol {\nabla }\varphi \boldsymbol {\times }\boldsymbol {\nabla }\psi$, where $\psi$ denotes the magnetic poloidal flux and the toroidal component $B_{0\varphi }(\psi )$ is a flux function. In (5.43), the total toroidal angular momentum density

(5.44)\begin{equation} \mathcal{P}^{*}_{\textrm{gy}\varphi} = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\left( P_{\textrm{gc}\varphi}^{*} + \boldsymbol{\varPi}_\textrm{gy}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} \right) + \frac{\mathbb{D}_\textrm{gy}}{4{\rm \pi} c}\boldsymbol{\times}\epsilon\boldsymbol{B}_{1}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} \end{equation}

is the sum of three groups of terms.

The first group in (5.44) is defined as the gyrocentre moment of the guiding-centre toroidal angular momentum:

(5.45)\begin{equation} P_{\textrm{gc}\varphi}^{*} \equiv \frac{e}{c}\boldsymbol{A}_{0}^{*}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} ={-}\frac{e}{c}\psi + p_{\|}\,b_{0\varphi} - J \left[ 2b_{0z} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\frac{1}{2 B_{0}}\boldsymbol{\nabla}\psi\right)\right], \end{equation}

which contains higher-order guiding-centre corrections (Tronko & Brizard Reference Tronko and Brizard2015). In a careful numerical analysis of the exact particle orbits of energetic ions in a tokamak magnetic field, Belova, Gorlenkov & Cheng (Reference Belova, Gorlenkov and Cheng2003) have shown that the higher-order guiding-centre corrections to the lowest-order guiding-centre toroidal angular momentum $P_{\textrm {gc}\varphi }^{*} = -(e/c)\psi + p_{\|}b_{0\varphi } + \cdots$ play a crucial role in the guiding-centre toroidal angular momentum law (i.e. in the absence of electromagnetic field perturbations). We note that it is a common practice to extract the dominant guiding-centre contribution from $-(e/c)\psi$ using the identity

(5.46)\begin{equation}{-}\frac{\partial }{\partial t}\left( \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\frac{e}{c}\psi \right) - \boldsymbol{\nabla}\boldsymbol{\cdot}\left( \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}\dot{\boldsymbol{X}}F\frac{e}{c}\psi \right) ={-} \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\frac{e}{c}\dot{\psi} \equiv{-}\frac{1}{c}J_\textrm{gy}^{\psi}, \end{equation}

where the radial velocity $\dot {\psi } \equiv \dot {\boldsymbol {X}}\,\boldsymbol {\cdot }\,\boldsymbol {\nabla }\psi$ is expressed in terms of the gyrocentre velocity $\dot {\boldsymbol {X}}$:

(5.47)\begin{equation} \dot{\psi} = \boldsymbol{\nabla}\psi\boldsymbol{\cdot}\left( \boldsymbol{E}_\textrm{gy}^{*}\boldsymbol{\times}\frac{c\boldsymbol{\mathsf{b}}_\textrm{gy}^{*}}{B_{\textrm{gy}\|}^{**}} + \frac{\partial K_\textrm{gy}}{\partial p_{\|}}\,\frac{\boldsymbol{B}_\textrm{gy}^{*}}{B_{\textrm{gy}\|}^{**}} \right). \end{equation}

Hence, we may now define $P_{\textrm {gc}\varphi } \equiv P_{\textrm {gc}\varphi }^{*} + (e/c)\psi$, and thus (5.43) becomes

(5.48)\begin{equation} \frac{\partial \mathcal{P}_{\textrm{gy}\varphi}}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\boldsymbol{\mathsf{T}}_{\textrm{gy}}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}\right) = \frac{1}{c}J_\textrm{gy}^{\psi}, \end{equation}

where the toroidal angular momentum density (5.44) is now defined with $P_{\textrm {gc}\varphi }$.

The second group in (5.44), which appears because of the symplectic momentum perturbation $\boldsymbol {\varPi }_\textrm {gy}$, contains the toroidal components of the perturbed $E\times B$ velocity and magnetic flutter momentum:

(5.49)\begin{align} \boldsymbol{\varPi}_\textrm{gy}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} &= \epsilon\left(\langle\boldsymbol{E}_{1\textrm{gc}}\rangle\boldsymbol{\times}\frac{e\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{\varOmega_{0}} + \frac{p_{\|}}{B_{0}}\langle\boldsymbol{B}_{1\bot\textrm{gc}}\rangle\right) \boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} \nonumber\\ &= \frac{e}{B_{0}\varOmega_{0}} \left( \langle\boldsymbol{E}_{1\textrm{gc}}\rangle + \frac{p_{\|}\hat{{\boldsymbol{\mathsf{b}}}}_{0}}{mc}\boldsymbol{\times} \langle\boldsymbol{B}_{1\textrm{gc}}\rangle\right) \boldsymbol{\cdot}\boldsymbol{\nabla}\psi, \end{align}

which can be expressed in terms of the radial component of the perturbed gyrocentre force, where we used the tokamak identity

(5.50)\begin{equation} \boldsymbol{B}_{0}\boldsymbol{\times}\partial\boldsymbol{x}/\partial\varphi =\boldsymbol{\nabla}\psi. \end{equation}

The third group in (5.44) contains the toroidal component of the Minkowski electromagnetic momentum (Abiteboul et al. Reference Abiteboul, Garbet, Grandgirard, Allfrey, Ghendrih, Latu, Sarazin and Strugarek2011):

(5.51)\begin{equation} \frac{\mathbb{D}_\textrm{gy}}{4{\rm \pi} c}\boldsymbol{\times}\epsilon\boldsymbol{B}_{1}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} = \frac{1}{4{\rm \pi} c}\left[\left(\epsilon\boldsymbol{E}_{1} + 4{\rm \pi}\mathbb{P}_\textrm{gy}\right) \boldsymbol{\times} \epsilon\boldsymbol{B}_{1}\right]\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}. \end{equation}

The partial time derivative of this term can be directly obtained from the toroidal component of (5.28). We note that, in the electrostatic limit (i.e. in the absence of magnetic field perturbations), we recover the flux-averaged gyrokinetic toroidal angular momentum density previously derived (without guiding-centre corrections, i.e. $P_{\textrm {gc}\varphi } = p_{\|}b_{0\varphi }$) (Hahm et al. Reference Hahm, Diamond, Gürcan and Rewoldt2007; Scott & Smirnov Reference Scott and Smirnov2010; Abiteboul et al. Reference Abiteboul, Garbet, Grandgirard, Allfrey, Ghendrih, Latu, Sarazin and Strugarek2011; Brizard & Tronko Reference Brizard and Tronko2011).

Finally, we proceed with a flux-surface average (Brizard & Tronko Reference Brizard and Tronko2011):

(5.52)\begin{equation} [\kern-1pt[\cdots]\kern-1pt] \equiv \frac{1}{\mathcal{V}}\oint(\cdots)\mathcal{J}_{\psi}\,\textrm{d}\vartheta\,\textrm{d}\varphi, \end{equation}

where $\mathcal {V}(\psi ) \equiv \oint \mathcal {J}_{\psi }\,\textrm {d}\vartheta \,\textrm {d}\varphi$ is the surface integral of the magnetic coordinate Jacobian $\mathcal {J}_{\psi } \equiv (\boldsymbol {\nabla }\psi \boldsymbol {\times }\boldsymbol {\nabla }\theta \,\boldsymbol {\cdot }\,\boldsymbol {\nabla }\varphi )^{-1} = 1/B_{0}^{\theta }$. The flux-surface average (5.52) satisfies the property

(5.53)\begin{equation} [\kern-1pt[\boldsymbol{\nabla}\boldsymbol{\cdot}{\boldsymbol{C}}]\kern-1pt] \equiv \frac{1}{\mathcal{V}}\,\frac{\partial }{\partial \psi}\left( \mathcal{V}[\kern-1pt[ \boldsymbol{C}\,\boldsymbol{\cdot}\, \boldsymbol{\nabla}\psi]\kern-1pt] \right)\end{equation}

for any vector field $\boldsymbol {C}$. In a time-independent axisymmetric tokamak geometry, we note that $\partial /\partial t$ also commutes with magnetic surface averaging. The magnetic surface-averaged gyrokinetic canonical angular momentum conservation law (5.48) becomes

(5.54)\begin{equation} \frac{\partial [\kern-1pt[\mathcal{P}_{\textrm{gy}\varphi}]\kern-1pt]}{\partial t} + \frac{1}{\mathcal{V}}\,\frac{\partial }{\partial \psi}\left( \mathcal{V}\left[\left[ T^{\psi}_{\textrm{gy}\varphi}\right]\right]\right) = \frac{1}{c}\,[\kern-1pt[ J_\textrm{gy}^{\psi}]\kern-1pt], \end{equation}

where $T^{\psi }_{\textrm {gy}\varphi } \equiv \boldsymbol {\nabla }\psi \,\boldsymbol {\cdot }\,\boldsymbol{\mathsf{T}}_\textrm {gy}^{*}\,\boldsymbol {\cdot }\,\partial \boldsymbol {x}/\partial \varphi$ is defined as

(5.55)\begin{equation} T^{\psi}_{\textrm{gy}\varphi} = \int_{\boldsymbol{P}} \mathcal{J}_\textrm{gy}F\dot{\psi}\left( P_{\textrm{gc}\varphi} + \boldsymbol{\varPi}_\textrm{gy}\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} \right) - \frac{\epsilon}{4{\rm \pi}} \boldsymbol{\nabla}\psi\boldsymbol{\cdot} \left( \mathbb{D}_\textrm{gy}\boldsymbol{E}_{1} +\boldsymbol{B}_{1}\mathbb{H}_\textrm{gy}\right) \boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}, \end{equation}

where we have used $\boldsymbol {\nabla }\psi \boldsymbol {\cdot }\partial \boldsymbol {x}/\partial \varphi = 0$ and $\dot {\psi }$ is given in (5.47). We note that, using the tokamak identity (5.50), the third term in (5.55), which contains the polarization term derived by McDevitt et al. (Reference McDevitt, Diamond, Gürcan and Hahm2009) in the parallel limit, can be expressed as

(5.56)\begin{equation} \boldsymbol{\nabla}\psi\boldsymbol{\cdot}\left( \frac{\epsilon}{4{\rm \pi}}\mathbb{D}_\textrm{gy}\boldsymbol{E}_{1} \right)\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi} = \frac{\partial \boldsymbol{x}}{\partial \varphi}\boldsymbol{\cdot}\left[ \frac{\epsilon}{4{\rm \pi}}\left(\mathbb{D}_\textrm{gy}\boldsymbol{\times}\boldsymbol{B}_{0}\right)\boldsymbol{E}_{1} \right]\boldsymbol{\cdot}\frac{\partial \boldsymbol{x}}{\partial \varphi}, \end{equation}

and similarly for the fourth term. Similar terms have appeared in the toroidal angular momentum transport analysis of Parra & Catto (Reference Parra and Catto2010b).