Kinetic theory for massive spin-1/2 particles from the Wigner-function formalism

We calculate the Wigner function for massive spin-1/2 particles in an inhomogeneous electromagnetic field to leading order in the Planck constant. Going beyond leading order in $\hbar$ we then derive a generalized Boltzmann equation in which the force exerted by an inhomogeneous electromagnetic field on the particle dipole moment arises naturally. Carefully taking the massless limit we find agreement with previous results. The case of global equilibrium with rotation is also studied. Finally, we outline the derivation of fluid-dynamical equations from the components of the Wigner function. The conservation of total angular momentum is promoted as an additional fluid-dynamical equation of motion. Our framework can be used to study polarization effects induced by vorticity and magnetic field in relativistic heavy-ion collisions.


I. INTRODUCTION
Relativistic heavy-ion collisions (HICs) create a new phase of hot and dense strong-interaction matter, the quark-gluon plasma (QGP) [see e.g. Ref. [1]]. The interaction rates between its constituents are sufficiently large that the matter rapidly reaches a state which can be described by fluid dynamics [2]. In non-central HICs the global angular momentum generates a non-vanishing vorticity of the QGP fluid. Furthermore, in such collisions a strong magnetic field is formed due to the electric current produced by the spectator protons constituting the colliding ions.
In the QGP, quarks can be considered as (nearly) massless fermions. The interplay between the chiral anomaly on the one hand and the magnetic field and the fluid vorticity on the other hand gives rise to novel transport phenomena called chiral effects. Two such phenomena are the chiral magnetic effect (CME) [3] and the chiral vortical effect (CVE) [4], where a charge current is induced along the direction of the magnetic field and the vorticity, respectively. Large-scale experimental efforts are currently under way to discover these phenomena in HICs [for a recent review, see Ref. [5]].
From the theoretical point of view, it is therefore mandatory to develop a theory which allows to study such transport phenomena in chiral fluids. One approach is chiral kinetic theory, which has been derived using various methods, e.g. the classical action [6][7][8][9][10][11], the Wigner function [12][13][14][15][16], and the world-line formalism [17][18][19]. In Refs. [12,13] it was shown that, using Wigner functions, one is able to recover the "side-jump" phenomenon first discussed in Refs. [10,11] in order to ensure total angular-momentum conservation in binary collisions. Furthermore, the inclusion of the chiral effects in fluid dynamics was studied in Refs. [4,20,21].
Another intriguing phenomenon occurring in the rotating QGP is that particles in the medium can be polarized in a way resembling the Einstein-de Haas [22] and Barnett effects [23]. Recently, the STAR Collabora-tion presented experimental evidence for the alignment of the spin of Λ hyperons with the global angular momentum in peripheral HICs [24]. This finding revealed, for the first time, the strong vortical structure of the QGP. Many theoretical works have explored spin-polarization mechanisms triggered by vorticity in HICs. In particular, the importance of the spin-orbit interaction [25][26][27] and the relation between spin polarization and thermal vorticity in local thermodynamical equilibrium have been studied [28][29][30][31]. A fluid-dynamical description, which includes the space-time evolution of the spin polarization, was proposed in Refs. [32][33][34][35]. However, this formulation is based on a specific choice for the energy-momentum and spin tensors. The physical implications of different sets of energy-momentum and spin tensors in fluid dynamics was investigated in Ref. [36].
Although there has been intense theoretical activity which has led to a deeper understanding of the transport properties of chiral matter, few studies have attempted to derive a covariant kinetic theory for massive particles using Wigner functions [37,38]. The aim of this paper is to fill this gap. We derive kinetic theory for massive spin-1/2 particles in an inhomogeneous electromagnetic field as a basis to study polarization effects in HICs. Our starting point is the covariant formulation of the Wigner function [39][40][41][42][43][44][45]. In order to solve the equations of motion for the Wigner function, we employ an expansion in the Planck constant and truncate at the lowest nontrivial order. This approximation is valid if the following two assumptions hold: (i) γ µ ∇ µ W ≪ mW , where W is the Wigner function, m is the particle mass, and ∇ µ represents gradients, and (ii) ≪ ∆R∆P , where ∆R is a spatial scale over which the electromagnetic field tensor varies significantly and ∆P a momentum scale over which the Wigner function varies significantly. Assumption (i) implies that the -expansion is effectively a gradient expansion. Assumption (ii) allows to truncate the power-series expansion of the Bessel functions entering the equations of motion of the Wigner function [41].
Under these assumptions, we first give an explicit derivation of the leading-order solution. Then, consid-ering the equation of motion for the Wigner function to first and second order in , we derive a generalized Boltzmann equation, where the external force acting on the particles is given by two contributions. The first one is the Lorentz force, which gives rise to the usual Vlasov term, and the second one is the Mathisson force [46], i.e., the force exerted on the particle's dipole moment in an inhomogeneous electromagnetic field. In our context, the dipole moment arises from the spin of the particle. We show how to take the massless limit, obtaining a result that agrees with previous works [12,13]. We also study the solution of the Boltzmann equation in the case of global equilibrium with rigid rotation. Finally, we derive fluid-dynamical equations of motion with spin degrees of freedom from the Wigner function using the canonical definitions of the energy-momentum and spin tensors. In accordance with previous works [32,36], the conservation of the total angular momentum is promoted as an additional fluid-dynamical equation, where the divergence of the spin tensor is related to the antisymmetric part of the energy-momentum tensor.
We use units c = k B = 1 throughout this paper. It is useful to explicitly keep Planck's constant , since it will be our power-counting parameter. The convention for the metric tensor is g µν = diag(+1, −1, −1, −1) and ǫ 0123 = −ǫ 0123 = +1 for the rank-four Levi-Civita tensor. We use the notation a µ b µ ≡ a · b for the scalar product of two four-vectors a µ , b µ and a · b for the corresponding scalar product of two spatial vectors a, b. A two-dimensional vector in spin space is denoted by a. The electromagnetic four-potential is A µ , where the electromagnetic charge is absorbed into its definition. We denote the dipolemoment tensor as Σ µν . This quantity corresponds to the spin tensor S µν of Refs. [10,11]. In this paper the term "spin tensor" is reserved for the rank-three Lorentz tensor S λ,µν .

II. EQUATIONS FOR THE WIGNER FUNCTION FOR MASSIVE FERMIONS
The Wigner function is defined as the Fourier transform of the two-point correlation function [40], Here, x 1 and x 2 are the space-time coordinates of two different points, with y µ ≡ x µ 1 −x µ 2 and x µ ≡ (x µ 1 +x µ 2 )/2. The gauge link is defined as In this paper, A µ will be treated as an external, classical field (otherwise, the gauge link would need to be path-ordered). The particular choice of path for the integration between x 1 and x 2 ensures that p µ is the kinetic momentum. Note that the factors 2π in the denominator in Eq. (1) belong to the phase-space volume and do not participate in the -counting employed throughout this paper. Starting from the Dirac equation and its adjoint, where D µ ≡ ∂ xµ + i A µ is the covariant derivative, one can derive the kinetic equation for the Wigner function as [40] ( Here one has defined the operator with the generalized space-time derivative and momentum operators where ∆ ≡ 2 ∂ p · ∂ x and F µν = ∂ µ x A ν − ∂ ν x A µ is the electromagnetic field-strength tensor. We should emphasize that in Eq. (4) the space-time derivative ∂ x contained in ∆ only acts on F µν , but not on the Wigner function. The functions j 0 (x) = sin x/x and j 1 (x) = (sin x − x cos x)/x 2 are spherical Bessel functions. If we assume that the particles only interact with the classical electromagnetic field but not among themselves (which, in the language of kinetic theory, is the limit of the collisionless Boltzmann-Vlasov equation), Eq. (4) is exact and contains the full dynamics of the Wigner function.
In order to derive a kinetic equation for massive spin-1/2 particles, it is advantageous to decompose the Wigner function in terms of a basis formed by the 16 independent generators of the Clifford algebra {1, γ 5 , γ µ , γ 5 γ µ , σ µν }, with γ 5 ≡ iγ 0 γ 1 γ 2 γ 3 and σ µν ≡ The coefficients F , P, V µ , A µ , and S µν are real functions of the phase-space coordinates x, p and correspond to the scalar, pseudo-scalar, vector, axial-vector, and tensor components of the Wigner function. Some of them have an obvious physical meaning [47]. For example, V µ is the fermion four-current and A µ is related to the spin density. Using the trace properties of the Dirac matrices, the coefficients in Eq. (8) are given by Replacing W in Eq. (4) by the decomposition (8), we find the following complex-valued equations: where Decomposing these equations into their real and imaginary parts, we obtain a set of coupled equations which determine the coefficients in the decomposition (8) of the Wigner function. The real parts read 2 and the imaginary parts are In the next sections, we will explicitly solve Eqs. (11) - (20) to zeroth order in , and then derive kinetic equations which the general solution has to fulfill up to first order in .

III. ZEROTH-ORDER SOLUTION
To zeroth order in , the operator K µ = p µ and Eq. (4) reduces to The solution is given by [37,48] where are the contributions from positive and negative energies, respectively. Here, f + rs (x, q) ≡ a † (q, r)a(q, s) and f − sr (x, q) ≡ b † (q,r)b(q, s) are the distribution functions for fermions and anti-fermions, respectively, which are in general matrices in spin space. The spin indices label spin states parallel, r, s = +, or anti-parallel, r, s = −, to the quantization direction s in the rest frame of the particle, respectively.
This spin quantization direction can in principle be chosen arbitrarily. However, the most convenient choice is to quantize the spin with respect to the polarization direction [37,41]. In other words, we choose a spin basis in which the new distribution functionsf ± rs are diagonal, i.e.,f In App. A we demonstrate that such a choice is always possible, at the expense of introducing space-time dependent spinors, cf. Eq. (A7). We will also use the diagonal basis in the calculation of the contributions of higher order in in the following sections. As shown in App. A, the spin quantization direction n (0)µ is given by Here, n ± is the spin quantization direction in the rest frame of the particle/anti-particle [cf. Eq. (A6)] and E p = p 2 + m 2 . The spin quantization direction n ± transforms as an axial vector under Lorentz boosts and parity transformations. We show in App. A that n ± depends in general on p and x, thus n ±µ is defined locally. The vector n (0)µ is aligned with the polarization direction and agrees with the classical spin vector, i.e., as we will see later, it obeys the classical equation for spin precession in an electromagnetic field, the so-called Bargmann-Michel-Telegdi (BMT) equation [49]. Moreover, n (0)µ fulfills p · n (0) = 0 (which can be seen using Eqs. (A9) and (A10) and applying the Dirac equation for the u-and v-spinors as well as the identitȳ u(p, r)γ 5 u(p, s) =v(−p, r)γ 5 v(−p, s) = 0).
Equations (22) -(24) represent the solution obtained in Ref. [48] for vanishing electromagnetic fields. However, this is also the solution for non-vanishing electromagnetic fields, since the form of Eq. (4) remains the same. The momentum variable p µ is then the kinetic (and not the canonical) momentum.
Closer inspection of Eq. (4) reveals that Eq. (22) with Eqs. (23), (24) is also a solution to Eq. (4) at arbitrary order in , if γ·∇W (0) = 0 and γ µ F µν ∂ pν W (0) = 0 (because then the -dependence of the operator K µ vanishes). In the absence of electromagnetic fields, one at least needs to require that γ · ∂ x W (0) = 0. In the full solution, i.e., the solution to all orders in , the momentum variable q is no longer equal to the kinetic momentum p. This is obviously not the case for Eqs. (23), (24), since they are proportional to ∼ δ 4 (p ∓ q), see also the discussion in Ref. [48]. Now we easily obtain the coefficients of the decomposition (8) using Eqs. (9) and (25). We find with and where e = ±, f (0)e s are the distribution functions in the diagonal basis, and the dipole-moment tensor is defined as for the proof, see App. A.
From their definitions (6), (7), we observe that the operators ∇ µ and Π µ can be expanded in a series of powers in 2 . In order to derive the semi-classical limit, we may truncate these series at order 0 and 2 , respectively, where We also expand the functions F , P, V µ , A µ , S µν into power series in , e.g., Inserting these expansions into Eqs. (11) - (20) and then comparing order by order in one can get a set of equations which we will analyze up to second order in in the remainder of this section.

A. Zeroth order in
We first analyze the on-shell conditions (33) for the scalar and tensor components to leading order in and show that the direct calculation of the Wigner function to this order presented in Sec. III is consistent with these conditions. To order O( 0 ), Eq. (33) reads where we have used p ν S (0) νµ = 0, which is the constraint equation (18) to zeroth order in . The general solution of the above equations reads where µν A (0) are up to now arbitrary functions which do not have singularities at p 2 = m 2 . We also demand that they go to zero sufficiently fast for large momenta (in order to neglect boundary terms when performing an integration by parts). Comparing to the previous section, we can identify V (0) with the spinsymmetric combination (29) and A (0) with spin-antisymmetric combination (30) of the zeroth-order distribution function, as well as Σ (0) µν with the dipole-moment tensor, which satisfies p µ Σ (0) µν = 0 in order to fulfill Eq. (18). In order to be consistent with Eq. (31), we demand With the help of Eq. (32) we can now write down the remaining components of the Wigner function to leading order in , It is straightforward to check that our solutions (37), (38) satisfy Eqs. (16) - (20). All zeroth-order solutions are on mass-shell and agree with the results from the direct calculation of the Wigner function in Sec. III.

B. First order in
The starting point for our analysis of the contributions of next-to-leading order in is again the on-shell equation (33). The O( ) part reads where we used p µ S (0) νµ = 0 and the relation which follows from Eq. (18) to first order in . Here the leading-order functions S µν and F (0) have been obtained in the previous subsection. The solutions to Eq. (39) can in general be written as Here,Σ µν is, up to a factor m, the on-shell part of the first-order dipole moment. The functions V (1) andΣ (1) µν will be determined from the kinetic equations that we will derive below. The function V (1) can be identified as the O( ) correction to the spin-symmetric combination of the distribution function. Using Eq. (40), we derive a constraint forΣ Expanding all quantities in Eq. (32) into power series in , to O( ) we obtain Inserting the zeroth-and first-order solutions from Eqs. (37), (38), and (41), we can derive the first-order pseudoscalar, vector, and axial-vector functions, is the first-order on-shell correction to n (0) µ A (0) . To first order in , the constraints (16), (20) read They lead to the kinetic equations of the particle distributions and the dipole moment to zeroth order in ; for details see App. C, As we have shown in the previous subsection, the zeroth-order kinetic equations are derived from the firstorder constraint equations. In order to obtain the firstorder kinetic equations, we focus on the second-order parts of Eqs. (16), (20), with the operator Π After some calculation (cf. App. C), one derives the following kinetic equations, Multiplying the second equation (49) by − 1 2m ǫ αβµν p β and using Eq. (45), we obtain a kinetic equation forn whereF µν ≡ 1 2 ǫ µναβ F αβ is the dual field-strength tensor.

V. BOLTZMANN EQUATION FOR SPIN-1/2 PARTICLES
In this section, we will derive the Boltzmann equation for spin-1/2 particles. Noting that we combine Eqs. (37) and (41) and use Eqs. (29) and (30) to obtain to order O( ) . Thus, the on-shell condition is to first order in modified to Since the energies of particles and anti-particles differ in sign, the distribution functions have support on different mass shells. Thus, they have to fulfill Eqs. (47) and (49) separately [41]. Then, using Eqs. (47) and (49), we find to order O( ) The last term can be shown to vanish. To this end, we choose the solution of the third equation (47) where the term multiplying the δ function is zero (this choice will be also made in the remainder of this paper, whenever employing Eqs. (47) and (49) in the calculations). Then, the last term in Eq.
Equation (55) is the main result of the present paper. In the following, we discuss the massless and the classical limit, as well as some consequences for global equilibrium and fluid dynamics.

VI. MASSLESS LIMIT
In this section, we explain how to obtain the massless limit of the currents V µ and A µ , cf. Eqs. (38) and (44). The crucial step is to replace the dipole-moment tensor (31) for m = 0 by the corresponding one for m = 0. Obviously, this cannot be achieved simply by taking the limit m → 0 in Eq. (31).
For massive particles, the dipole-moment tensor as well as the particle's position are uniquely defined in the rest frame. The Pauli-Lubanski operator is defined as [50] whereP µ ≡ i D µ is the (kinetic) momentum operator.
On the classical level, Σ µν is the intrinsic angularmomentum tensor about the center of mass. In a relativistic theory, the center of mass of a particle is framedependent. In order to have a frame-independent definition of Σ µν , one requires p µ Σ µν = 0 as a gauge condition. This requirement identifies the dipole-moment tensor (57) as the intrinsic angular-momentum tensor about the center of mass in the rest frame of the particle [51].
For massless particles there is no rest frame, thus both the position (in the classical case the center of momentum) and the dipole-moment tensor can at first be defined in an arbitrary frame, which makes them framedependent. For massless particles, the polarization vector n µ is always parallel to the momentum p µ . Thus, the requirement p µ Σ µν = 0 can no longer be used as a gauge condition, since Eq. (57) automatically satisfies this constraint. [In the massless limit, one also needs to change the normalization of the spinors toūu = 2|p| [48].] If we choose the dipole-moment tensor to be defined in a frame characterized by a time-like four-vector u µ , we can choose the gauge condition u µ Σ µν = 0 [11]. Consequently, the frame vector u µ must assume the role of p µ in Eq. (57). Moreover, since n µ and p µ are parallel for massless particles, the momentum p µ can assume the role of n µ in Eq. (57). Finally, in order to obtain the massless case we need to replace the normalization factor 1/m in Eq. (57). The energy of a massive particle in its rest frame is p 0 rf = p 2 . If the particle is on the mass-shell, this is equivalent to p 0 rf = m. The energy of a massless particle in the rest frame of u µ , however, is p 0 u = p · u. Thus, it is natural to replace the normalization 1/m in Eq. (57) by 1/(p · u). We emphasize that this replacement can only be done in the presence of a δ-function which sets the rest-frame energy equal to the mass m. The explicit expression for the dipole-moment tensor in the massless case is then given by which agrees with the definition of the "spin tensor" in Ref. [11]. This tensor corresponds classically to the intrinsic angular momentum about the center of momentum as seen from the frame where u µ = (1, 0, 0, 0) and will have the quantum-mechanical properties of an angular-momentum operator in that frame.
With this knowledge, we can make the transition between the Wigner functions of massive and massless particles. For zero fermion mass, Eqs. (15) and (20) decouple. By defining right-and left-handed currents J µ χ ≡ 1 2 (V µ + χA µ ), χ = ± for right-/left-handed particles, we have to order These equations have been solved in Refs. [12,14,15], with the result where f χ is the distribution function for right-/lefthanded fermions and u µ is the four-velocity of an arbitrary frame. We remark that in the massive case, s describes "spin up" or "spin down", which corresponds to positive or negative helicity in the massless limit (mn µ → p µ ). On the other hand, the currents above are defined for given chirality χ. Since helicity and chirality are identical for massless particles, but opposite for massless anti-particles, the relation between chirality χ and spin/helicity s is χ = e s with e = ± representing particles/anti-particles.
To obtain the massless limit of our solutions, we replace the massive dipole-moment tensor by the massless one, Σ (0)µν → Σ µν u . In order to obtain the vector current for the massless case from Eq. (44), we need to consider the term ∼ ∇ (0)ν Σ (0) µν . We first pull the constant factor 1/m out of the derivative and then replace δ(p 2 −m 2 )/m = δ(p 2 −m 2 )/ p 2 → δ(p 2 )/(p·u). Finally, replacing p µ /m → u µ , mn µ → p µ in this term we find In Ref. [12] the frame-vector u µ is assumed to be independent of space-time coordinates. In order to compare to the solution found in that reference, we adopt the same assumption. Evaluating the derivatives, contracting the ǫ-tensors, and using p 2 δ ′ (p 2 ) = −δ(p 2 ), we find from Eqs. (38) and (61) where V = V (0) + V (1) . Note that V (1) depends on the frame vector u µ such that the whole expression (62) is frame independent [11,14,16]. To obtain the axialvector current in the massless case from Eqs. (38) and (44), we note that the general solution of Eq. (42) reads where the first and second terms depend on arbitrary time-like unit vectors u µ and v µ , respectively. Here, one makes use of the first equation (47) to see that the constraint (42) is fulfilled. Inserting Eq. (63) into Eq. (44), and replacing the zeroth order dipole-moment tensor Σ (0) µν by Σ u,µν , we find where A ≡ A (0) + A (1) , with A (1) dependent on u µ . Note that, in order for the above axial current to be frame-independent, the function A (1) cannot depend on v µ . Adding and subtracting Eqs. (62) and (64), we recover the result in Eq. (60).

VII. CLASSICAL LIMIT
In this section, we show that Eq. (55) gives rise to the first and second Mathisson-Papapetrou-Dixon (MPD) equations [46,52] as well as to the BMT equation [49], which were derived for classical, extended, spinning particles with non-vanishing dipole moment. Comparing Eq. (55) to the generic form of the collisionless relativistic Boltzmann-Vlasov equation where f s is the distribution function and F µ s is the external force, we find that in our case i.e., the external force is given as the sum of the Lorentz force and the Mathisson force [46]. In Refs. [46,52], the kinetic equation for particles with classical dipole moment m µν was derived. Our results agree with the latter, setting with Bohr's magneton µ B ≡ e /(2m), where e is the electric charge, and the gyromagnetic ratio g = 2, as expected for Dirac particles with spin 1/2. The classical dipole moment does not depend on the kinetic momentum, so that after dividing by m Eq. (65) becomeṡ where the four-velocity of the particles isẋ µ ≡ dx µ /dτ = p µ /m, where τ is the world-line parameter, anḋ With Eq. (66), this is the first MPD equation [46,52].
The evolution of the dipole-moment tensor is given by the third equation (47), which can be rewritten as where we usedΣ with F µ s given by Eq. (66) to zeroth order. Equation (70) is identical to the second MPD equation [46,52]. Using Eq. (31), we obtain Inserting Eq. (70) and contracting with ǫ ρσµν yields Contracting with p ρ and using Eqs. (66), (69) to zeroth order in , we conclude that This is the BMT equation for classical spin precession in an electromagnetic field [49].

VIII. GLOBAL EQUILIBRIUM
Equation (55) determines the single-particle distribution function f ± s in a general non-equilibrium state. A special solution is obtained in global equilibrium, which we will consider in this section.
A necessary condition for equilibrium is vanishing entropy production. In this case, the distribution function must have the form f eq s = (e gs + 1) −1 , with g s being a linear combination of the collisional invariants, namely, charge, kinetic momentum p µ , and total angular momentum which is the sum of orbital angular momentum L µν = x [µ p ν] and spin angular momentum, which to first order is given by the dipole-moment tensor s 2 Σ (0)µν . (Also the canonical momentum π µ is conserved in a collision and could be used instead of the kinetic momentum. Here, we will at first use the kinetic momentum, since it is independent of space-time coordinates, as well as gaugeindependent.) Thus, Here, b µ (x), a s (x), and Ω µν (x) are Lagrangian multipliers, which can depend on x. Since J µν s is anti-symmetric, the symmetric part of Ω µν can be dropped without loss of generality.
Let us consider the case of global equilibrium with rigid rotation. Using Eq. (76), Eq. (77) can be written as where In global equilibrium, the Boltzmann equation (55) needs to be fulfilled. From the part of Eq. (55) proportional to the derivative of f eq s we obtain where we used Eq. (47). This equation is fulfilled, if which makes the terms in the first and second line of Eq.
(80) vanish. The terms in the third line of Eq. (79) can be shown to vanish if b µ is constant, since then Ω µν is equal to the thermal vorticity, i.e., For the proof, one also employs the relation which can be proven with the help of the homogeneous Maxwell equations and Eq. (80). These equilibrium conditions agree with those found in the classical case [52]. It is amusing to note that, without electromagnetic fields, the tensor Ω µν does not need to be equal to the thermal vorticity.
We introduce the Lie derivative of A µ along the direction of β λ as Choosing a gauge in which L β A µ = 0, we can rewrite Eq. (80) as Defining the function g s becomes Here, π µ ≡ p µ +A µ is the canonical momentum, U µ is the fluid velocity, β ≡ 1/T is the inverse temperature, and µ s is the chemical potential for particles with spin s (for anti-particles, we need to reverse the sign of the chemical potential). This form of g s , and thus the distribution function f eq s , agrees in the massless and field-free limit to the one suggested in Ref. [11]. Moreover, recalling the definition (A12), (A13) of the dipole-moment tensor one can prove that the distribution function agrees with the one proposed in Ref. [29] to first order in if µ + = µ − and if the electromagnetic field vanishes. However, f eq s does not agree with the distribution function used in Ref. [13]. In that work, the authors use the kinetic momentum (which they call canonical).
The part of Eq. (55) which is proportional to f ± s vanishes if µ + = µ − to zeroth order in . In the presence of a spin imbalance, ∆µ ≡ µ + − µ − = 0, it only vanishes if The reason that global equilibrium with spin imbalance can in general not be realized for massive particles is that in this case the axial-vector current is only conserved if the pseudo-scalar function P = 0, see also Eq. (12).
To zeroth order in , the distribution function is given by with We define the dual thermal vorticity tensor asω µν ≡ 1 2 ǫ µναβ ω αβ . Now we calculate the vector current by inserting the distribution function (75) into the equation µ , cf. Eqs. (38), (44). With where we used L β A α = 0, Taylor-expanding and noting that where we used p·n = 0 and δ(p 2 −m 2 )p 2 = δ(p 2 −m 2 )m 2 , we find The current given by Eq. (93) contains contributions which are not parallel to p µ . To first order in , particles are not transported parallel to their momenta. The term containingF µν in Eqs. (93) is caused by off-shell effects and describes the vector current induced by electromagnetic fields, which yields the analogue of the CME in the case of non-zero mass. On the other hand, the term containingω µν describes the current induced by vorticity and thus gives the analogue of the CVE. We furthermore calculate the axial-vector current. In order to do so, it is convenient to decompose the tensor Σ (1) µν introduced in Eq. (41) in the following way, The tensor Ξ µν is anti-symmetric and satisfies p µ Ξ µν = 0. On the other hand, χ µν represents the dipole moment induced by the gradients of the distribution function since, according to Eq. (42), it satisfies Inserting V (0) into Eq. (95) and using Eq. (90) we can derive the following constraint for χ µν , where we adopted the short-hand notation The most general solution of Eq. (96) can be written as where κ 1,2 are arbitrary coefficients which satisfy κ 1 + κ 2 = 1, and v µ is an arbitrary vector such that v · p = 0. Other possible terms which vanish when being contracted with the momentum are absorbed into Ξ µν . The decomposition into χ µν and Ξ µν is not unique, but allows for the transformations with C being an arbitrary function of x and p. For any value of κ 2 , we can apply the above transformation with C ≡ κ 2 to Eq. (98), which yields Thus, we can set κ 2 = 0 without loss of generality. In other words, it is always possible to isolate the contribution proportional to ω µν in the decomposition forΣ (1) µν . This decomposition will assume a physical meaning when looking at the kinetic equation for S µν .
Inserting Eq. (94) into Eq. (49), we obtain Noting that p·∇ (0) χ µν = 0 and using Eqs. (82) and (100), we find that the χ µν -dependent part vanishes and which is the second MPD equation for Ξ µν . This part of the dipole-moment corresponds, together with the zeroth-order dipole moment, to the classical spin precession in electromagnetic fields. We now derive from Eq. (44) the full axial-vector part of the Wigner function up to first order in , i.e., By looking at the different terms in Eq. (103), we identify three contributions to the axial-vector current in the massive case. The first term in the second line and the term in the last line describe the spin precession in the presence of an electromagnetic field according to the BMT equation. We remark that the function Ξ µν is not specified and has to be determined through Eq. (102). The second term in the second line gives rise to the axial current in the direction of the vorticity, which is the analogue of the axial chiral vortical effect (ACVE). Finally, the term in the third line describes the axial current along the magnetic field, which is the analogue of the chiral separation effect (CSE). In the case where µ + = µ − , i.e, f − , the only non-zero contributions are the analogue of the ACVE and CSE (even though not specified the last term contains the difference between the distribution functions of spin up and spin down). This result is in agreement with previous works, where only a vanishing spin imbalance was considered [31,37].

IX. FLUID-DYNAMICAL EQUATIONS
In this section, we present the equations of motion of the fluid-dynamical variables, i.e., of the net particlenumber current and the energy-momentum tensor. We also give an equation for the spin tensor, which supplements these equations in the case of spin-1/2 particles.
The net particle-number current is defined as Inserting the zeroth-and first-order solutions (38), (44) into Eq. (104) we obtain where dP ≡ d 4 p δ(p 2 − m 2 ). Equation (16) represents the conservation law for the vector component of the Wigner function. Integrating this equation over kinetic 4-momentum, we immediately obtain the conservation law for the net particle-number current, where we assumed that F µν is independent of p ν and V µ vanishes sufficiently rapidly for large momenta, which ensures the vanishing of a boundary term.
The Lagrangian operator for a Dirac spinor in an electromagnetic field is [41] From the Lagrangian we can derive the canonical energymomentum tensor as follows, where we have separated the total energy-momentum tensor into three parts: the gauge-invariant matter part T µν mat , the part containing the interaction between gauge potential and matter current, T µν int , and the electromagnetic part T µν em , Note that none of these are in general symmetric under µ ↔ ν. The total energy-momentum tensor is conserved ∂ xν T µν = 0, which can be checked using the Dirac and Maxwell equations. However, the matter part is not conserved, This equation can be derived by acting ∂ xµ on V µ in the definition of T µν mat , cf. first equation (109), then using Eq. (16), and finally integrating by parts. Inserting Eqs. (38), (44) into the energy-momentum tensor, we get The total canonical angular momentum tensor is calculated as follows, :ψ{γ λ , σ µν }ψ : The first two terms, x µ T λν − x ν T λµ , can be interpreted as the orbital angular-momentum tensor. The remaining terms constitute the spin angular-momentum tensor, which can be further separated into a matter and a field part. The spin tensor of matter can be defined as [29] S λ,µν mat (x) ≡ 1 4 :ψ{γ λ , σ µν }ψ : With the help of Eq. (20) we find, to any order in , where we assumed that boundary terms vanish. Thus, the spin of matter is not conserved separately. To zeroth order in , T µν mat is symmetric according to Eq. (111), thus both sides of Eq. (114) vanish. To first order in , both sides are non-zero. Inserting the zeroth-order Wigner function into Eq. (113) we obtain The above expressions for the energy-momentum and spin tensor emerge directly from Noether's theorem and thus correspond to the canonical ones. However, one can obtain different sets of tensors by applying pseudogauge transformations that keep the conservation laws for energy-momentum and spin. It has been shown that using different sets of tensors related through this pseudogauge freedom is not equivalent and leads to different measurable quantities [36]. Note that with Eq. (114), we can also prove that ∂ xν T µν mat = F µα J α , the form of the equation of motion for the matter energy-momentum tensor given in Refs. [53,54].

X. CONCLUSIONS
In this paper we have derived kinetic theory for massive spin-1/2 particles in an inhomogeneous electromagnetic field starting from the covariant formulation of the Wigner function. Carrying out an expansion in and truncating it at first order, we found a general solution of the equations of motion. We showed how to consistently take the massless limit and demonstrated agreement with previous works, which describe the CME and CVE. One of the crucial results of our work is the derivation of the collisionless Boltzmann equation for particles that carry a dipole moment due to their spin. We also recovered well-known results in the classical limit. The external force acting on the particles is the sum of the Lorentz force and the Mathisson force, i.e., the first MPD equation. The time evolution of the dipole moment follows the second MPD equation, and the spin polarization precesses according to the BMT equation. Moreover, as an example, we studied the case of a rigidly rotating fluid in global equilibrium. In particular, we found the conditions that the Lagrange multipliers related to the conservation of charge, energy, momentum, and angular momentum have to satisfy in order for the distribution function to be a solution of the Boltzmann equation. Finally, fluiddynamical equations of motion are provided, in which the spin tensor is included among the evolved densities.
A straightforward extension of this work could be the inclusion of a collision term into our generalized Boltzmann equation and the derivation of the equations of motion for dissipative relativistic magneto-hydrodynamics for spin-1/2 particles. This could be achieved using the method of moments, following Refs. [53,54], where this has already been done for spin-0 particles. Another potential extension would be the derivation of a transport equation starting from the equal-time Wigner-function formalism [42].
In the rest frame, the standard spinors u and v are given as [50] u(0, +) = √ 2m (1, 0, 0, 0) T , Note that u(0, +) corresponds to a particle with spin parallel to the z-direction, while v(0, +) corresponds to an anti-particle with spin anti-parallel to the z-direction. We diagonalize the distribution functions f (0)e in the rest frame, with D e being 2 × 2 matrices in spin space, where d e ± are the eigenvectors of b e · σ corresponding to the eigenvalues ±, respectively, where n e ≡ b e / √ b e · b e is the unit vector along the direction of b e . Note that the distribution functions f (0)e in general depend on the space-time coordinates x µ , thus the transformation matrices D e as well as n e are defined locally. We then define the following spinors, which can be derived by rotating the standard ones, The spinorsũ(x, 0, ±) andṽ(x, 0, ±) now correspond to particles/anti-particles with spin parallel/anti-parallel to n ± . Using Eqs. where n ±µ is given by Eq. (27). We rewrite the axialvector current as where the vector n (0) µ (x, p) and the distribution function A (0) (x, p) are determined by Eqs. (26) and (30), respectively.
zeroth-and first-order solutions we obtain In order to derive the kinetic equation for the firstorder dipole-moment tensor, we first need V (2) µ , which is calculated by expanding Eq. (32) into a series in and identifying the 2 term, Inserting this, as well as A µν from Eqs. (37) and (41) into Eq. (C8) and using the above commutators, one obtains the kinetic equation forΣ (1) µν .