Relativistic kinetic theory for spin-1/2 particles: Conservation laws, thermodynamics, and linear waves

TL;DR

This paper develops a comprehensive relativistic kinetic model for spin-1/2 particles, detailing conservation laws, thermodynamics, and linear wave behavior, with implications for quantum plasma physics.

Contribution

It introduces a fully relativistic kinetic framework for spin-1/2 particles, including conservation laws, thermodynamic equilibrium, and linear dispersion relations, expanding the theoretical understanding of quantum relativistic systems.

Findings

Derived conservation laws for energy, momentum, and angular momentum.

Provided thermodynamic equilibrium distributions in various limits.

Calculated linear dispersion relations for coupled spin-momentum eigenfunctions.

Abstract

We study a recently derived fully relativistic kinetic model for spin-1/2 particles. Firstly, the full set of conservation laws for energy, momentum and angular momentum are given, together with an expression for the (non-symmetric) stress-energy tensor. Next, the thermodynamic equilibrium distribution is given in different limiting cases. Furthermore, we address the analytical complexity that arises when the spin- and momentum eigenfunctions are coupled in linear theory, by calculating the linear dispersion relation for such a case. Finally, we discuss the model and give some context by comparing with potentially relevant phenomena that are not included, such as radiation reaction and vacuum polarization.

Equations

$$f(\mathbf{x},\mathbf{p},\mathbf{s},t)$$

$$=\frac{1}{4\pi }\left({\delta }_{\alpha \beta }+\mathbf{s}\cdot {\mathbf{\sigma }}_{\alpha \beta }\right)W_{\beta \alpha }(\mathbf{x},\mathbf{p},t),$$

$$⟨\Phi ⟩:=\int d\Omega \Phi f$$

$$⟨\mathbf{x}|{\mathbf{\sigma }}_{\alpha \beta }{\rho }_{\beta \alpha }|\mathbf{x}⟩=3\int d\Omega \mathbf{s}f.$$

$$\int d\Omega {\nabla }_{s}f$$

$$\int d\Omega \mathbf{s}\times \mathbf{X}\cdot {\nabla }_{s}f$$

$$0$$

$$\mathbf{T}$$

$$\tilde{T}$$

$$\nabla \cdot \mathbf{E}$$

$$\nabla \cdot \mathbf{B}$$

$$\nabla \times \mathbf{E}$$

$$\nabla \times \mathbf{B}$$

$${\rho }_f$$

$${\mathbf{j}}_f$$

$$\mathbf{P}$$

$$\mathbf{M}$$

$$\mathbf{v}=\frac{\mathbf{p}}{ϵ}-{\mu }_B{\nabla }_p\mathbf{T}\cdot 3\mathbf{s}$$

$$W=\frac{1}{2}(E^2+B^2)+\int d\Omega \left(ϵ-3{\mu }_{B}m\frac{\mathbf{B}}{ϵ}\cdot \mathbf{s}\right)f$$

$$\begin{gathered} \mathbf{K}=\int d\Omega \left[ϵ+{\mu }_{B}m3\mathbf{s}\cdot \left(\frac{\mathbf{B}}{ϵ}-\frac{\mathbf{p}\times \mathbf{E}}{ϵ(ϵ+m)}\right)\right]\mathbf{v}f \\ +\mathbf{E}\times \mathbf{H}, \end{gathered}$$

$${\partial }_{t}W+\nabla \cdot \mathbf{K}=0.$$

$${\partial }_t⟨\mathbf{p}⟩={\partial }_t\int d\Omega \mathbf{p}f=\int d\Omega \mathbf{p}{\partial }_{t}f.$$

$$-\int d\Omega\,\mathbf{p}\partial_{t}f=\int d^{3}p\,\mathbf{p}\ \int d^{2}s\,\big{[}\mathbf{v}\cdot\nabla_{x}f\\ +q(\mathbf{E}+\mathbf{v}\times\mathbf{B})\cdot\nabla_{p}f+3\mu_{B}(\mathbf{s}\cdot\partial_{x_{i}}\mathbf{T})\partial_{p_{i}}f\big{]}.$$

$$\begin{gathered} q\int d\Omega p_i(E_j+{ϵ}_{jkl}{v}_k{B}_l){\partial }_{p_j}f \\ =-q\int d\Omega {\delta }_{ij}\left(E_j+{ϵ}_{jkl}{v}_k{B}_l\right)f \\ =-{\left({\rho }_{\text{f}}\mathbf{E}+{\mathbf{j}}_{\text{f}}\times \mathbf{B}\right)}_i \end{gathered}$$

$$\begin{gathered} 3{\mu }_B\int d\Omega p_i(s_k\cdot {\partial }_{x_j}{T}_k){\partial }_{p_j}f \\ =-3{\mu }_B\int d\Omega ({\delta }_{ij}{s}_k{\partial }_{x_j}{T}_k+p_i{s}_k{\partial }_{x_j}{\partial }_{p_j}{T}_k)f \\ =\int d\Omega (-3{\mu }_B{s}_k{\partial }_{x_i}{T}_k+p_i{\partial }_{x_j}{v}_j)f. \end{gathered}$$

$$3{\mu }_B\int d\Omega s_k({\partial }_{x_j}{T}_k)f=M_k{\partial }_{x_j}{B}_k+P_k{\partial }_{x_j}{B}_k$$

$$\begin{gathered} {\partial }_t⟨p_i⟩=-\int d\Omega \left[p_i{v}_j{\partial }_{x_j}+p_i({\partial }_{x_j}{v}_j)-3{\mu }_B{s}_j({\partial }_{x_i}{T}_j)\right]f+{\left({\rho }_{\text{f}}\mathbf{E}+{\mathbf{j}}_{\text{f}}\times \mathbf{B}\right)}_i \\ =-{\partial }_{x_j}\int d\Omega p_i{v}_{j}f+{\left({\rho }_{\text{f}}\mathbf{E}+{\mathbf{j}}_{\text{f}}\times \mathbf{B}\right)}_i+M_j{\partial }_{x_i}{B}_j+P_j{\partial }_{x_i}{E}_j \end{gathered}$$

$$T_{ij}^{\text{e}}=⟨p_i{v}_j⟩$$

$$\begin{gathered} {\partial }_t(\mathbf{D}\times \mathbf{B}{)}_i+{\partial }_{x_j}{T}_{ij}^{\text{EM}} \\ =-{\left({\rho }_f\mathbf{E}+{\mathbf{j}}_f\times \mathbf{B}\right)}_i-M_k{\partial }_{x_i}{B}_k-P_k{\partial }_{x_i}{E}_k \end{gathered}$$

$$T_{ij}^{\text{EM}}=\frac{1}{2}\left(E^2+B^2-2\mathbf{M}\cdot \mathbf{B}\right){\delta }_{ij}-H_i{B}_j-E_i{D}_j.$$

$${\partial }_t{\left(⟨\mathbf{p}⟩+\mathbf{D}\times \mathbf{B}\right)}_i+{\partial }_{x_j}(T_{ij}^{\text{e}}+T_{ij}^{\text{EM}})=0.$$

$$:{\sigma }_i{\hat{v}}_j:={\sigma }_i\frac{{\hat{p}}_j}{m}-{\mu }_B{\partial }_{p_j}{\hat{T}}_i.$$

$$⟨\mathbf{x}|:{\sigma }_i{\hat{v}}_j:\hat{\rho }|\mathbf{x}⟩=\int d\Omega \left(3s_i\frac{p_j}{m}-{\mu }_B{\partial }_{p_j}{T}_i\right)f$$

$$\int d^{2}s{s}_i{s}_{j}f=\frac{1}{3}{\delta }_{ij}\int d^{2}sf$$

$$\int d^{2}s{s}_i{\partial }_{s_j}f=\frac{1}{4\pi }\int d^{2}s{s}_i\mathrm{Tr}[({\sigma }_j-s_j{s}_k{\sigma }_k)W]=0$$

Full text

I Introduction

Current and next-generation high-intensity laser facilities provide opportunities to study highly relativistic electron dynamics. Here the electron spin Hu and Keitel (1999); Walser and Keitel (2000); Walser et al. (2002); Wen et al. (2016) is of interest, as the electrons may spin-polarize Del Sorbo et al. (2017); Sorbo et al. (2018); Li et al. (2019), in just a few laser cycles. The ultra-strong magnetic fields present in astrophysical environments Harding and Lai (2006); Uzdensky and Rightley (2014); Mahajan et al. (2014) also enhance the significance of spin-related phenomena.

In less extreme settings, applications in, for example, spintronics Wolf et al. (2001), quantum wells Manfredi and Hervieux (2007), and plasmonics Atwater (2007) have also stimulated an interest in high density plasmas where quantum effects are significant Shukla and Eliasson (2011); Haas (2011).

Earlier spin-kinetic models have included effects such as the magnetic dipole force, magnetization currents, and spin precession Brodin et al. (2008); Zamanian et al. (2010); Asenjo et al. (2012); Hurst et al. (2014); Andreev (2015, 2017), but have been limited to particle velocities well below the speed of light in vacuum. However, in a previous paper Ekman et al. (2017), which we will refer to as Paper I below, a fully relativistic kinetic equation for spin-$\frac{1}{2}$ was presented, along with its couplings to Maxwell’s equations. This forms a model that can be used to describe self-consistent relativistic plasma dynamics, including spin effects.

Here, our purpose is to study some basic properties of the model, and also provide some tools needed for analytical studies. The organization of the paper is as follows: In Section II we give an overview of the model, including the assumptions made in the derivation. Next, in Section III, we derive conservation laws for energy, momentum, and angular momentum and give an expression for the stress energy tensor. It turns out that the stress-energy tensor is non-symmetric in our case, which is related to the presence of spin angular momentum. Section IV is devoted to a study of the thermodynamic background distribution for various cases and the associated background magnetization due to the electrons. In particular we present the expression for the relativistic and the non-relativistic case allowing for, in both cases, non-degenerate and degenerate electrons. Next, in Section V, we analyze linear theory in a homogeneous magnetized medium. It turns out due to relativistic effects the spin- and momentum eigenfunctions become coupled, and the standard solution procedure must be generalized. We give an example of how this can be done, and present a dispersion relation in the limiting case of wave propagation perpendicular to the magnetic field. Finally, in Section VI, the significance of our results and the applicability of the model are discussed.

II Overview of the model

The model comes from separating positive and negative energy solutions of the Dirac equation by means of a Foldy-Wouthuysen (F-W) transformation Foldy and Wouthuysen (1950); Silenko (2008). Since we are decoupling electrons and positrons, the physical condition of applicability is that pair production is negligible. Quantitatively, the fields should not be comparable to the critical fields, $E\ll E_{c}=m^{2}/q\hbar$ and similarly for $B$, and their typical scale lengths should be long compared to the Compton wavelength $\hbar/m$. We then take a gauge-invariant Wigner transformation Wigner (1932); Stratonovich (1956); Serimaa et al. (1986) (see Ref. Case (2008) for a pedagogical introduction to the Wigner formalism) to obtain an evolution equation for a $2\times 2$-matrix valued Wigner function $W_{\alpha\beta}(\mathbf{x}_{i},\mathbf{p}_{i},t)$. Here $\mathbf{x}_{i},\mathbf{p}_{i}$ denote that in a many-body system the Wigner function depends on all the positions and momenta. However, the BBGKY hierarchy applies, and neglecting collisions the evolution equation for the one-particle Wigner function is found from the one-particle Hamiltonian.

Next, we apply a spin transformation Zamanian et al. (2010),

[TABLE]

to obtain a scalar Wigner function (summation convention applied to the spin indices). In this formalism, densities in space are given by moments of the Wigner function over the momentum and spin variables $\mathbf{p}$ and $\mathbf{s}$. For example, the number density is $n=\int d^{3}pd^{2}s\,f$. Below, we will often use the notation $d\Omega=d^{3}pd^{2}s$ and

[TABLE]

where $\Phi$ is some function on phase space. The normalization of the spin transformation is such that the spin density is

[TABLE]

Also, it is easy to derive the identities

[TABLE]

where $\mathbf{X}$ is any vector independent of $\mathbf{s}$. These identities are used repeatedly below.

The evolution equation for the scalar Wigner function $f$ is

[TABLE]

where $\epsilon^{2}=\mathbf{p}^{2}+m^{2}$, $\mu_{B}=q\hbar/2m$ is the Bohr magneton and

[TABLE]

The system is closed with Maxwell’s equations, in units where $c=\epsilon_{0}=\mu_{0}=1$,

[TABLE]

where $\mathbf{P}$ and $\mathbf{M}$ are the polarization and magnetization densities, and $\rho_{f}$ and $\mathbf{j}_{f}$ are the free charge and current densities. These are given by

[TABLE]

and we refer to Paper I for the derivation. It follows from the evolution equation, Eq. 6, that the free charge is conserved, $\partial_{t}\rho_{f}+\nabla\cdot\mathbf{j}_{f}=0$, and we interpret

[TABLE]

as the function on phase space corresponding to the velocity – it is in fact the Weyl transform of the velocity operator $\hat{\mathbf{v}}=\frac{i}{\hbar}[\hat{H},\hat{\mathbf{x}}]$ given by the Heisenberg equation of motion. The spin-dependent term is related to the “hidden momentum” Shockley and James (1967); Shockley (1968); Coleman and Van Vleck (1968); Babson et al. (2009) of systems with magnetic moments. It is discussed further in Paper I, with additional references. Here we make some further observations:

Anomalous magnetic moment

The derivation of the model considered here is based on the Dirac theory where the gyromagnetic ratio, $g=2$, discarding the anomalous magnetic moment (AMM). However, making this assumption is not necessary for the Foldy-Wouthuysen transformation Silenko (2008), and we will allow for $g\neq 2$ below. Let us discuss the validity of this.

Including the AMM corresponds to adding a term $(g-2)F_{\mu\nu}\overline{\psi}\sigma^{\mu\nu}\psi$ to the Dirac equation and consequently new terms will appear in the Hamiltonian in the F-W representation, rather than just modifying $g$. One of those new terms has the form $(g-2)\bm{\sigma}\cdot\mathbf{B}$, i.e., unlike the Dirac magnetic moment, there is no $1/\varepsilon$ length contraction factor. Since $g-2\approx\frac{\alpha}{\pi}$ is small, the modification $g\mapsto g^{\prime}=2+\alpha/\pi$ should be satisfactory for a phenomenological description. Non-relativistically, such an analysis was carried out in Ref. Brodin et al. (2008) and later more generally in Ref. Lundin and Brodin (2010), finding new wave-particle resonances due to the mismatch between the electron cyclotron and the Larmor precession frequencies. We will continue in this spirit and consider $g-2\neq 0$ in Section V.

Radiation reaction

As it stands, the model does not include the effect of radiation reaction (RR, reviewed in, e.g., Refs. Di Piazza et al. (2012); Burton and Noble (2014)). In principle, a radiation reaction force could be added to Eq. 6 as a $\mathbf{F}_{\text{RR}}\cdot\nabla_{p}f$ term, possibly including the spin-dependence of the RR force Bhabha and Corben (1941). It is also now well known how to include RR in a particle-in-cell scheme Gonoskov et al. (2015). However, the ratio of the RR to the Lorentz force is $\eta=\alpha\gamma^{2}E/E_{\text{crit}}$ where $\gamma$ is the Lorentz factor and $\alpha$ the fine-structure constant. The strong RR regime $\eta\approx 1$ is only expected to be reached with next-generation $10\text{\,}\mathrm{P}\mathrm{W}$ laser facilities Di Piazza et al. (2012) (but see Refs. Cole et al. (2018); Poder et al. (2018) for a recent experiment using the Gemini laser), and by inspection there are clearly regimes where $\eta\ll 1$, but $\gamma$ is large enough that an $O(v^{2}/c^{2})$ treatment is inapplicable. Furthermore, even for strong fields, the spin-dependent non-RR forces may be comparable to the RR force or even dominate it Mahajan et al. (2014). However, the present model is not developed exclusively for strong fields and relativistic spin effects could be important, e.g., for relativistic temperatures (thermodynamic or Fermi).

III Conservation laws and stress-energy tensor

As a sanity check, the model should have conservation laws for energy, momentum, and angular momentum. In the previous paper, the conservation law for energy was given and discussed in relation to the Abraham-Minkowski dilemma; here we give the full set of conservation laws.

III.1 Energy and momentum

The total energy density is given by

[TABLE]

and with the energy flux vector

[TABLE]

where $\mathbf{H}=\mathbf{B}-\mathbf{M}$, we have a conservation law on divergence form

[TABLE]

To derive the conservation of momentum, we look at

[TABLE]

We substitute for $\partial_{t}f$ using the kinetic equation.

Because $\mathbf{p}$ is independent of $\mathbf{s}$, we can divide the integral as $\int d^{3}p\,\mathbf{p}\int d^{2}s\,\partial_{t}f$ and using the spin integral identities Eqs. 4 and 5 we arrive at

[TABLE]

The first term here is $\int d\Omega\,p_{i}v_{j}\nabla_{x_{j}}f$, however, the $x$-divergence cannot be taken outside the integral to write this as the divergence of a moment, because $\mathbf{v}$ depends on $\mathbf{x}$ through $\mathbf{T}$.

In the second term we integrate by parts. Since $\partial_{p_{j}}v_{k}$ is symmetric in its indices, we have

[TABLE]

which is the Lorentz force density. The spin-dependent part will be found below.

For the third term in Eq. 19, we have

[TABLE]

For the first term above, it is simple to establish that

[TABLE]

using the definitions of $\mathbf{M},\mathbf{P}$ in Eqs. 13 and 12. The second term in Eq. 21, containing $\partial_{x_{j}}v_{j}$, is what we need for the divergence of a moment with an $\mathbf{x}$-dependence. It is of course not a miracle that this happens, but a consequence of that the kinetic equation essentially has the form of a Poisson bracket between $f$ and the Hamiltonian.

To summarize,

[TABLE]

using, on the second line, Eq. 22. Here we can identify

[TABLE]

as the stress tensor for the electrons. Note that the relation between $\mathbf{p}$ and $\mathbf{v}$ also contains field variables, and there will be no clean separation of the total stress tensor into “field” and “particle” parts. This is what one would expect for an interacting theory.

We now need to find an appropriate Poynting vector and electromagnetic stress tensor to match this force density. From Maxwell’s equations, using various vector identities, one finds that

[TABLE]

where $\mathbf{D}=\mathbf{E}+\mathbf{P}$, and the tensor $T_{ij}^{\text{EM}}$ is defined by

[TABLE]

where the source term is precisely the negative of that in Eq. 23. We should note that while we have used notation to indicate that this is the “field” part of the stress tensor, it still contains particle variables through the magnetization and polarization.

Thus, we can express conservation of momentum as

[TABLE]

The stress tensor we have found is not symmetric, and because it is $\mathbf{E}\times\mathbf{H}$ that appears in the energy conservation law, the full stress-energy tensor is not symmetric under $0i\leftrightarrow i0$. This is usually considered a defect, for two reasons, discussed in more detail in Landau and Lifshitz vol. 2 Landau and Lifshitz (1987), §32–33 and §94.

The first reason is that to conserve angular momentum $\mathbf{r}\times\bm{\Pi}$, the stress tensor should be symmetric. However, this condition applies only if the angular momentum is entirely orbital angular momentum. If there is spin angular momentum, angular momentum can be conserved even with a stress tensor that is not symmetric, and we shall demonstrate explicitly that this is the case for our model, below.

Identifying the correct stress-energy tensor is a general problem for field theories Forger and Römer (2004), since to a tensor $T_{\mu\nu}$ with $\partial_{\mu}T^{\mu}{}_{\nu}=0$, we can add any 4-divergence $\partial_{\rho}\psi^{\rho}{}_{\mu\nu}$ with $\psi_{\rho\mu\nu}$ anti-symmetric in $\rho\mu$ to obtain a new tensor with vanishing 4-divergence. Hence, the integrals $\int d^{3}x\,T^{i0}$ – the actual conserved quantities – will be unchanged, but their densities in space are not uniquely defined without a further physical principle. As a familiar example, indeed, the tensor obtained by applying Noether’s theorem to the free electromagnetic field is neither symmetric nor even gauge invariant, but can be made so [][; Ch.~12.10.]Jackson. Still, merely requiring the tensor to be symmetric may not be enough to guarantee uniqueness Forger and Römer (2004).

The physical principle that ensures uniqueness is that the components of the stress-energy tensor are in principle observable: one needs simply a ruler precise enough to measure the curvature of spacetime. This is the second reason to prefer a symmetric stress-energy tensor: it is the source of gravitation, and the left-hand side of Einstein’s field equations is the symmetric Einstein tensor. This also provides a method to obtain the stress-energy tensor directly: it is the variation of the matter action with respect to the metric.

For manifestly Lorentz-covariant theories, the Belinfante-Rosenfeld construction Belinfante (1940); Rosenfeld (1940) gives an explicit form for $\psi_{\rho\mu\nu}$ to symmetrize the Noether tensor, but this construction relies on the structure of the Lorentz group. Because our model is in Hamiltonian form and written in terms of $\mathbf{E},\mathbf{B}$, we have broken manifest covariance and cannot use the Belinfante-Rosenfeld construction directly.

The Belinfante-Rosenfeld method, or the variation with respect to the metric, can be applied to the manifestly covariant Lagrangian in the Dirac representation, yielding a symmetric tensor $\overline{\psi}\hat{T}_{\mu\nu}^{\text{D}}\psi=\overline{\psi}(\gamma_{\mu}D_{\nu}+\gamma_{\nu}D_{\mu})\psi$, $D_{\mu}$ being the gauge covariant derivative. An alternative approach would therefore be to transform the operator $\hat{T}_{\mu\nu}^{\text{D}}$ to the Foldy-Wouthuysen representation. However, there are issues with gauge invariance to consider Goldman (1977), and this is outside the scope of the present paper.

III.2 Spin transport

To establish conservation of angular momentum, we consider $\partial_{t}\langle\mathbf{s}\rangle$. Here we must establish a few facts about the spin moments. Because products of Pauli matrices are reducible to Pauli matrices, there are no higher spin moments. That is, because $\sigma_{i}\sigma_{j}=\delta_{ij}+i\epsilon_{ijk}\sigma_{k}$, we must have $\int s_{i}s_{j}f\,d^{2}s\propto\delta_{ij}$ since moments of the Wigner function should correspond to symmetrically ordered operators.

In calculating the spin transport, we will want the moment corresponding to the operator $:\mathrel{\sigma_{i}\hat{v}_{j}}:$ where the colons indicate Weyl ordering. Using the Pauli matrices relation, we find that

[TABLE]

Because $f=\frac{1}{4\pi}\operatorname{Tr}[(1+\mathbf{s}\cdot\bm{\sigma})W]$, the moment corresponding to this operator is

[TABLE]

where the factor of $3$ is needed only in the first term. Also using the definition of $f$, we find

[TABLE]

and, importantly

[TABLE]

because all odd moments over the sphere vanish by reflection symmetry.

Thus, to find

[TABLE]

we again look at the terms in the kinetic equation. The Lorentz force term will not contribute anything: the electric force is independent of $\mathbf{p}$, so we get $\int d^{3}p\,\nabla_{p}f=0$; when integrating the magnetic force term by parts we get $\nabla_{p}\cdot(\mathbf{p}\times\ldots)=0$. For the first term, we get

[TABLE]

using the identities above.

The index structure in the spin torque term will be

[TABLE]

where the term $4$:th order in $s_{i}$ does not contribute because it is contracted with the anti-symmetric Levi-Civita. Thus the spin-torque term, unsurprisingly and reassuringly, gives the spin precession as in the Bargmann-Michel-Telegdi equation Bargmann et al. (1959):

[TABLE]

The term containing $\nabla_{x}\tilde{T}$ in the kinetic equation is treated similarly to the first. It gives an $x$-derivative on the function in Eq. 29 complementing that on $f$ in Eq. 31.

In conclusion, after multiplying by $\hbar/2$, giving the spin angular momentum, we find

[TABLE]

III.3 Angular momentum

Let $\bm{\Pi}=\mathbf{D}\times\mathbf{B}+\langle\mathbf{p}\rangle$ be the total linear momentum, and define the orbital angular momentum $\mathbf{L}=\mathbf{r}\times\bm{\Pi}$. Then

[TABLE]

and we see that the orbital angular momentum would be conserved, were $T_{kj}$ symmetric.

The source of orbital angular momentum is $-\epsilon_{ijk}T_{jk}$. For our stress tensor, this is

[TABLE]

using that $\mathbf{H}=\mathbf{B}-\mathbf{M},\mathbf{D}=\mathbf{E}+\mathbf{P}$, so only the magnetization/polarization contributes to the cross product.

Now,

[TABLE]

so

[TABLE]

using the Jacobi identity for the cross product.

We then recognize that

[TABLE]

But this source term is exactly the negative of the spin torque in (33), so the total angular momentum, orbital plus spin, is conserved. Hence, the reasons to prefer a symmetric stress-energy tensor, discussed above, are not relevant in our case.

IV The thermodynamic background Wigner function

The aim is to find the thermodynamic background Wigner function for electrons in a constant magnetic field $B_{0}\mathbf{\hat{z}}$, and also to compute the associated background magnetization. We will divide the treatment into the non-relativistic and the relativistic regime. Before we look into specific cases, we apply spherical coordinates in spin-space (where $\theta_{s}$ is the angle with the $z$-axis) and introduce the division of the equilibrium Wigner function $f_{0}$ into its spin-up and spin-down parts (see e.g. Ref. Zamanian et al. (2010)) according to

[TABLE]

and

[TABLE]

where $n_{0\pm}$ is the number density of the spin-up and spin-down populations respectively, and we write $n_{0}=n_{0+}+n_{0-}$. The factors $(1\pm\cos\theta_{s})$ correspond to the projection operators $1\pm\sigma_{z}$.

We note that non-relativistically, the background magnetization $\mathbf{M}_{0}=\mu_{B}\int\mathbf{s}f_{0}(\mathbf{p},\mathbf{s})\,d\Omega=M_{0}\mathbf{\hat{z}}$ is given by $M_{0}=\mu_{B}n_{0+}-\mu_{B}n_{0-}=R\mu_{B}n_{0}$ where the thermodynamic factor $R$ is defined by

[TABLE]

i.e., it is a functional of the distribution. Below, the value of $R$, i.e., the degree of spin-polarization, will be presented for a few specific cases and the relativistic generalization for the background magnetization will be given.

IV.1 The non-relativistic regime

For the non-relativistic case the characteristic kinetic energy $E_{k}$ should be much smaller than the electron rest mass energy. $E_{k}$ is the thermal energy $k_{B}T$ or the (non-relativistic) Fermi energy $E_{F}=(3\pi^{2}n_{0})^{2/3}\hbar^{2}/2m$, whichever is larger. Moreover, we assume the Zeeman energy $\mu_{B}B_{0}$ to fulfill $\mu_{B}B_{0}\ll E_{k}$. If this is fulfilled Landau quantization is not significant, which means that the energy states are continuous to a good approximation. The general non-relativistic expression for the Wigner function has been computed in Ref. Zamanian et al. (2010) (see Eq. (59) in that work). In the limit of $\mu_{B}B_{0}\ll E_{k}$, this expression reduces to

[TABLE]

where the normalization constant $C$ can be chosen such as to fulfill $\int f_{0}(\mathbf{p},\mathbf{s})\,d\Omega=n_{0}$. The chemical potential $\mu_{c}$ coincides with the Fermi energy for $T\ll T_{F}$, and for the opposite ordering it suffices to know that $\mu_{c}$ is large and negative. For $T\gg T_{F}$ we can use $\exp(-\mu_{c}/k_{B}T)\gg 1$ and get

[TABLE]

where the normalization condition gives us

[TABLE]

For the non-degenerate case given by (45), it is easily confirmed that the $R$-factor is given by the textbook result

[TABLE]

Next we consider the fully degenerate case where $T=0$. As a result, the number densities of spin-up and down particles are given by the volumes of the respective Fermi spheres, limited by the Fermi-momentum $p_{F\pm}=\sqrt{2m(E_{F}\mp\mu_{B}B_{0})}$ for spin-up and down states respectively. Thus we have

[TABLE]

The approximation made in (48) is a Taylor expansion to first order in $\mu_{B}B_{0}/k_{B}T_{F}$, in line with the general condition $\mu_{B}B_{0}\ll E_{k}$, which is needed to avoid the complications related to strong Landau quantization.

IV.2 The relativistic regime

Since $\mu_{B}B_{0}\ll E_{k}$ still holds relativistically, Landau quantization is not an issue. As a result, most of the previous section can be copied by simply replacing the velocity dependence with a momentum-dependence, the non-relativistic kinetic energy with $(\gamma-1)m$, replacing $\mu_{B}B_{0}$ with $\mu_{B}B_{0}/\gamma$ and now using the relativistic Fermi energy $E_{F}=\left((m)^{2}+c^{4}p_{F}^{2}\right)^{1/2}-m$ where $p_{F}^{2}=(3\pi^{2}n_{0})^{2/3}\hbar^{2}$. Here we may neglect the corrections of the relativistic factor due to the magnetic dipole energy Eq. 14 (there is no electric field in equilibrium, so the $\mathbf{E}\times\mathbf{s}$ term vanishes) and use $\gamma=(1+p^{2}/m^{2}c^{2})^{1/2}$. In the nondegenerate case (large negative chemical potential) we immediately get

[TABLE]

The corresponding thermodynamic factor thus becomes

[TABLE]

Since $\mu_{B}B_{0}$ is a correction term, we may Taylor expand the exponentials, in which case we get

[TABLE]

which agrees with the non-relativistic expression in the limit of $(\gamma-1)\ll 1$.

For a fully degenerate relativistic system we again use $f_{0}=F_{0+}(1+\cos\theta_{s})+F_{0-}(1-\cos\theta_{s})$ but with

[TABLE]

Thus the $R$-factor is determined by

[TABLE]

where $p_{F\pm}$ is determined from

[TABLE]

where we have used that the relativistic factor near the Fermi surface fulfills $\gamma_{F}\simeq 1+E_{F}/m$.

Simplifying Eq. 53 using $\mu_{B}B_{0}\ll m$ and $\mu_{B}B_{0}\ll E_{F}$ results in the final expression

[TABLE]

which coincides with the non-relativistic expression in the limit $E_{F}\ll m$.

It should be noted that the expression $M_{0}=R\mu_{B}n_{0}$ requires that each particle contributes with $\pm\mu_{B}$ to the magnetic moment. Since this is only true in the rest frame, we must compensate. If the magnetic dipole energy is large this would be complicated, but for $\mu_{B}B_{0}\ll m$ a simplified calculation suffices. The result is that the magnetization is reduced in proportion to the gamma-factor, i.e.,

[TABLE]

where $f_{0p}$ is the reduced momentum distribution function that would apply in the absence of spin-polarization (i.e. a relativistic Fermi-Dirac or a Synge-Jüttner distribution for the two cases considered above). The above expression gives the correct magnetization to first order in an expansion $\mu_{B}B_{0}/m$. It should be stressed that the expression for $M_{0}$ becomes significantly more complicated if higher order contributions are required.

Physically we can note that relativistic effects decrease the total background magnetization in two ways. Firstly, for relativistic particles the energy difference between the spin-up and spin-down states (in the laboratory frame) is smaller, reducing the difference in number density between the spin-up and down states (see e.g. Eq. 49). Secondly, the contribution to the magnetic dipole moment in the laboratory frame is smaller for a moving particle, reducing the magnetization a second time (see e.g. Eq. 56).

V Dispersion relation for linear waves

The enhanced complexity of Eq. 6, compared to non-relativistic theories, introduces some technical obstacles already in linearized theory. Our purpose here is to address these difficulties, and present a general calculation method.

After linearization Eq. 6 can be written as

[TABLE]

Here $\theta_{s}$ is the angle between the spin and $\mathbf{B}_{0}$. Some of the difficulties in Eq. 57, as compared to the classical Vlasov equation follow from $\gamma=\epsilon/m$ being a function of momentum. However, since it is well-known from classical relativistic theory theory how to deal with this, we will focus on some of the other subtleties, and use the lowest order approximations, $\epsilon/m\approx 1$, $\nabla_{p}(1/\epsilon)\approx-2\mathbf{p}/m^{3}$.

From now on we consider a restricted geometry, namely transverse waves with

[TABLE]

where $\mathbf{k}$ is the wave vector. The perturbed quantities follow the plane-wave ansatz according to $f_{1}=\tilde{f}_{1}e^{i(\mathbf{k}\cdot\mathbf{x}-\omega t)}$. We express the momentum of the particles in cylindrical coordinates $(p_{\perp},\varphi_{p},p_{z})$ and the spin in spherical coordinates, $\mathbf{s}=\sin\theta_{s}\cos\varphi_{s}\mathbf{e}_{x}+\sin\theta_{s}\sin\varphi_{s}\mathbf{e}_{y}+\cos\theta_{s}\mathbf{e}_{z}$.

Furthermore we make an expansion of $\tilde{f}_{1}$ in eigenfunctions of the operators of the right-hand side of Eq. 57 Lundin and Brodin (2010)

[TABLE]

where

[TABLE]

where $J_{\tau}$ is the Bessel function of first kind and $\omega_{ce}=qB_{0}/m$ is the electron cyclotron frequency. Here, Greek summation indices take integer values from $-\infty$ to $+\infty$ and we will suppress the argument of the Bessel functions, as it is always $k_{\perp}p_{\perp}/m\omega_{ce}$. Only $\alpha=\pm 1$ will contribute in Eq. 58. From now on, we will additionally use Latin summation indices that take only the values $\pm 1$ to distinguish between the two types of sums.

Using this eigenfunction expansion together with

[TABLE]

where $\omega_{cg}=(g/2)\omega_{ce}$ in Eq. 57, we get

[TABLE]

where RHS is the right-hand side of Eq. 57, simplified according to the above assumptions.

The differential equation Eq. 61 is relatively hard to solve analytically. However the troublesome operator $\cos\theta_{s}-\sin\theta_{s}\partial_{\theta_{s}}$ appears with a factor $\mu_{B}B_{0}/m$, which is of order of unity for magnetars, but is much smaller for other known environments Harding and Lai (2006). Hence we can solve this equation for $f_{1}$ by using perturbation theory. Expanding $f_{1}=f_{10}+f_{11}$, where $f_{11}$ is first order in $\mu_{B}B_{0}/m$, $f_{10}$ is the solution to Eq. 61, keeping only zeroth order terms in $\mu_{B}B_{0}/m$ in the right-hand side:

[TABLE]

The first order term $f_{11}$ is then obtained by taking $f_{1}=f_{10}$ in all first-order terms in Eq. 61:

[TABLE]

In these expansions,

[TABLE]

Now with $f_{1}$ expressed in terms of $f_{0}$, we can split $f_{0}$ into its spin parts, using the notation $f_{0\pm}=\frac{1}{4\pi}F_{0\pm}(1\pm\cos{\theta_{s}})$ as in the previous section, Eq. 41. This will give additional sums over $\pm$, analogous to sums over the Latin indices already used. In determining the dispersion relation, we calculate the total current density $\mathbf{J}=\mathbf{j}_{f}+\frac{\partial\mathbf{P}}{\partial t}+\nabla\times\mathbf{M}$ (see Eqs. 11, 12 and 13). Using Maxwell’s fourth equation, Eq. 9d, the dispersion relation is given by

[TABLE]

where $\chi^{f}_{\beta}$ is the contribution from the free current density

[TABLE]

$\chi^{M}_{\beta}$ is the contribution from the magnetization current density

[TABLE]

and $\chi^{P}_{\beta}$ is the contribution from the polarization current density

[TABLE]

Note that only in $\chi^{f}_{\beta}$ do we have contributions from $f_{11}$, this is because both $B_{\rho}$ and $C_{\rho}$ vanish in $f_{11}$ when using $f_{0\pm}=\frac{1}{4\pi}F_{0\pm}(1\pm\cos{\theta_{s}})$.

Eq. 67 contains both classical and quantum modes. Since we are interested in the effects of the spin on the dispersion relation, we consider a regime where the spin effects are comparable to the classical ones. Considering

[TABLE]

in this regime the denominators $\omega-\omega_{ce}\beta$, $\,\omega-\omega_{ce}\beta\pm\omega_{cg}$ and $\omega-\omega_{ce}\beta\pm\Delta\omega_{ce}$ are minimized for $\beta=0,\pm 1$, and [math] respectively. Keeping only these values of $\beta$ in the summation over $\beta$ and focusing on the short Larmor radius regime, i.e., $kp_{\perp}/m\omega_{ce}\ll 1$, we expand the Bessel functions to second order. We refer to Appendix A for more details.

Furthermore, to compute the integrals in Eq. 67 explicitly, we must specify the background distribution function. Considering a Fermi-Dirac distribution for low temperatures $T\ll T_{F}$, where $T_{F}$ is the Fermi-temperature, the dispersion relation in equation Eq. 67 is

[TABLE]

where $\delta=g/2-1$. For $k\rightarrow 0$, $\omega\rightarrow\Delta\omega_{ce}$ as expected. The equation above shows that $\omega$ diverges from $\Delta\omega_{ce}$ for larger wave numbers. This divergence is proportional to the square of the Compton wavelength, but since the solution is in the short Larmor radius regime, $k$ should satisfy $k\ll m\omega_{ce}/p_{\perp}$.

VI Discussion

The present quantum kinetic relativistic model for electrons, first derived in Ref. Ekman et al. (2017), is based on the Dirac equation, where a phenomenological adjustment of the spin $g$-factor has been made to account for the QED-contribution. It generalizes the relativistic Vlasov equation to include the spin dynamics, while still allowing for a fully relativistic motion. In comparison with the models found in Ref. Zamanian et al. (2010) it is extended to contain spin-orbit interaction, a contribution to the polarization current associated with the spin, and Thomas precession. As compared to the model of Ref. Asenjo et al. (2012) it allows for fully relativistic motion (i.e., relativistic gamma factors not close to unity). Here we have continued studying the model, including the full set of conservation laws for energy, momentum and angular momentum, as well as an expression for the stress-energy tensor. Moreover, the thermodynamic background Wigner function and the background magnetization are given for the degenerate and non-degenerate cases, in both the relativistic and non-relativistic cases.

A complication in the present theory when it comes to practical calculations is the non-trivial relation between the momentum variable and the velocity, involving the spin state. To illustrate how this can be handled analytically, treating the spin-dependence perturbatively, we have calculated the linear dispersion relation for the case of perpendicular propagation across an external magnetic field. Finally a simple limiting case with $\omega\approx\Delta\omega_{ce}$ has been presented.

Before we go on with the technical aspects, let us first discuss some specific systems where the current model is of particular interest, i.e. plasmas where the particle motion is relativistic and the electrons are spin-polarized. As described in Section IV, the thermodynamic background distribution has a degree of spin polarization of the order $\mu_{B}B_{0}/k_{B}T$ for a non-degenerate system and of the order $\mu_{B}B_{0}/k_{B}T_{F}$ for degenerate systems. For astrophysical objects such as pulsars ($B\sim${10}^{8}\text{,}\mathrm{T} at the surface) and magnetars ($B\sim${10}^{10}{10}^{11}\text{\,}\mathrm{T} at the surface) Uzdensky and Rightley (2014); Harding and Lai (2006), we typically have non-relativistic temperatures $T\sim${10}^{6}\text{,}\mathrm{K} in the atmosphere Page *et al.* ([1996](#bib.bib47)). This corresponds to a large electron spin-polarization $0.1<\mu_{B}B_{0}/k_{B}T<1$ even at a considerable distance from the neutron star surface. While the electron temperature may not be of relativistic magnitude, the particle interactions with the pulsar fields nevertheless induce relativistic motion in the atmosphere Mahajan *et al.* ([2013](#bib.bib48)); Da Costa and Kahn ([1991](#bib.bib49)). For electrons at the pulsar or magnetar surface, on the other hand, the degeneracy makes the particles strongly relativistic with a Fermi temperature of the order $T_{F}\sim${10}^{12}\text{\,}\mathrm{K} Lattimer and Prakash (2004), which, depending on the magnetic field, corresponds to a degree of spin polarization $0.0001<\mu_{B}B_{0}/k_{B}T_{F}<0.1$.

Thus for electrons belonging to the neutron star surface, we note that spin polarization is pronounced for magnetars but not as much for pulsars. Turning to laboratory applications, laser-plasma interactions with dense targets are of particular interest Shi et al. (2018); Lindman (2010); Wagner et al. (2004); Del Sorbo et al. (2017); Sorbo et al. (2018); Li et al. (2019); Marklund and Shukla (2006); Di Piazza et al. (2012); Stefan et al. (2011); Brodin et al. (2010). While such systems are not initially spin-polarized, theory aided by particle-in-cell simulations have predicted Shi et al. (2018) that quasi-static magnetic field with a field strength of the order $B_{0}\sim${10}^{5}\text{,}\mathrm{T}$$ can be formed. It should be stressed that these estimates are fully consistent with experiments, see e.g. Refs Lindman (2010); Wagner et al. (2004). Before the final stage of laser-plasma compression, where signficant electron heating may occur, this corresponds to a spin polarization of the order $0.01<\mu_{B}B_{0}/k_{B}T_{F}<0.1$, at the same time as the electrons are relativistic due to the quiver motion in the laser field Marklund and Shukla (2006); Di Piazza et al. (2012).

Here we have mainly been concerned with systems that, besides being relativistic, have a thermodynamic background distribution that is spin polarized. However, it is worth noting that interaction with a dynamical electromagnetic fields by itself can lead to spin polarization, see e.g. Refs Del Sorbo et al. (2017); Sorbo et al. (2018); Li et al. (2019); Stefan et al. (2011); Brodin et al. (2010). Moreover, we stress that the current theory deviates from simpler models to some extent already for modest gamma factors, see the discussion of hidden momentum in Paper I.

The current model is of particular interest for strong fields. However, classical or QED corrected radiation reaction Nikishov and Ritus (1964a, b); Harvey and Marklund (2012) is not covered and neither is strong field vacuum polarization due to QED included. Thus it is of interest to establish the regime of validity for the model, i.e. to find the conditions when the spin contributions are the most important extensions of the classical physics. First it should be noted that the non-relativistic spin-contributions are not necessarily dependent on strong fields, in particular the ratio of the magnetic dipole force over the Lorentz force scales as $\hbar k^{2}/m\omega$ for a plane wave field. Thus, in case we have a short-wavelength plasma perturbation with a low-frequency, the magnetic dipole force should be included independent of the field strength. While effects of this nature have been covered already in previous non-relativistic works Zamanian et al. (2010); Hurst et al. (2014) or weakly relativistic theories Asenjo et al. (2012); Stefan and Brodin (2013), the theory presented here is needed in case the thermodynamic temperature or the Fermi temperature is relativistic. In the absence of relativistic temperatures, however, the strong field effects of the present theory, whose relative importance is proportional to $\mu_{B}B/m$, need to be compared with the effects of radiation reaction and vacuum polarization.

The contribution from classical radiation reaction relative to the Lorentz force scales as the dimensionless parameter $R=2a_{0}^{2}\gamma\alpha\omega/3m$ where $a_{0}=qE_{0}/(\omega m)$ is the normalized vector potential, $E_{0}$ being the peak field strength and $\omega$ the laser frequency; and $\alpha$ is the fine-structure constant. Thus the relative importance of radiation reaction in relation to relativistic spin effects, in fields of strength $E$ and $B$, is given by

[TABLE]

For definiteness, if the relativistic particle velocities are induced by strong laser fields, such that $E/B\sim 1$ and $\gamma\sim a_{0}$, the condition for the spin effects to dominate over radiation reaction becomes

[TABLE]

Thus for moderately relativistic fields, the relativistic contribution to the spin dynamics is more important than radiation reaction, whereas for strong enough laser fields, the ordering is the opposite. However, for other types of field-configurations (i.e. not due to lasers), in particular in the vicinity of strongly magnetized objects (e.g. pulsars), we may have $E/B\ll 1$. As seen from Eq. 73, in this case the spin-relativistic effects may dominate over radiation reaction even for $\gamma\gg 10$.

Next, we note that our treatment assumes $\mu_{B}B/m<1$, which corresponds to a field magnitude below the critical field strength, i.e. the value for which pair production is exponentially suppressed. In case pair production does not take place, vacuum polarization scales as $(\alpha/90\pi)E^{2}/E_{\mathrm{cr}}^{2}$ Dittrich and Gies (2000); Marklund and Shukla (2006). Thus, provided the plasma is not dilute 111 Particle contributions enter a dispersion relation proportional to a factor $\omega_{p}^{2}/\omega^{2}$, and hence are suppressed for low densities or high frequencies.

, vacuum polarization is smaller than the spin-relativistic effects by a factor of the order

[TABLE]

where we have put $E=B$ to simplify the expression, which is justified since electric and magnetic fields induce effects of the same order.

Due to the complexity of the kinetic model, it is rather difficult to analyse for problems beyond linear waves in homogeneous media. Still, due to the richness of the physics included, the present theory can be useful in many different contexts. Firstly, Eqs. 6 to 13 can serve as a starting point for deriving simpler models, using e.g. a moment expansion Haas et al. (2010); [][; especiallyPaper~IV.]stefan2014models of the kinetic evolution equation. Secondly, the present theory can be used to determine the region of validity of many simpler models. Thirdly, potentially there are new phenomena to be found from Eqs. 6 to 13 even for the special case of linear homogeneous plasmas, as we have only covered a few of those possibilities here. Finally, the mathematical structure of the present theory can help shedding light over long-standing problems such as the Abraham-Minkowski controversy 222Discussed in Paper I Ekman et al. (2017), see also references therein..

In this paper we have written the model with the spin as an independent variable, which is convenient for analytical calculations. One can instead not take the spin transform and formulate the model in terms of the four independent components of the Hermitian matrix $W_{\alpha\beta}$, analogous to Ref. Hurst et al. (2014). In this case, the computational cost of simulating the model would be that of simulating four models for spin-less particles, plus some extra cost due to additional terms, which is not prohibitive compared to simulating the Vlasov equation.

Acknowledgements.

All authors would like to acknowledge financial support by the Swedish Research Council, grant number 2016-03806. G. Brodin, H. Al-Naseri and J. Zamanian would also like to acknowledge support by the Knut and Alice Wallenberg Foundation through the PLIONA project.

Appendix A The Harmonic Limit

In obtaining the harmonic dispersion relation in Eq. 72 from Eq. 67, several approximations have been done. Working with $\omega\approx\Delta\omega_{ce}$, gives us possibilities for minimizing the values of the denominators in Eq. 67 for certain values of $\beta$. Considering only the leading terms where the denominators are minimized according to

[TABLE]

where we set $\beta=0$, $\pm 1$ and [math] respectively. Note that when using $\beta=0$ in $\chi_{f}$ the contribution from $f_{11}$ in Eq. 67 becomes zero.

Considering the short Larmor radius regime where the argument of the Bessel functions is small, i.e., $kp_{\perp}/m\omega_{ce}\ll 1$, we can make a second order Taylor expansion of the Bessel functions. Furthermore, we consider the case where the background distribution is Fermi-Dirac with the low temperature $T\ll T_{F}$. Hence we have, $F_{0\pm}4\pi p_{F\pm}^{3}/3\approx n_{0}$. Eq. 67 is now

[TABLE]

Since we are looking for $\omega\approx\Delta\omega_{ce}$, we consider the $\omega^{2}-\Delta\omega_{ce}^{2}$-root and approximate $\omega^{2}$ in the right hand side of Eq. 76 to $\Delta\omega_{ce}^{2}$. We have finally the dispersion relation in Eq. 72.