Scalar quantum kinetic theory for spin-1/2 particles: mean field theory

Starting from the Pauli Hamiltonian operator, we derive a scalar quantum kinetic equations for spin-1/2 systems. Here the regular Wigner two-state matrix is replaced by a scalar distribution function in extended phase space. Apart from being a formulation of principal interest, such scalar quantum kinetic equation makes the comparison to classical kinetic theory straightforward, and lends itself naturally to currently available numerical Vlasov and Boltzmann schemes. Moreover, while the quasi-distribution is a Wigner function in regular phase space, it is given by a Q-function in spin space. As such, nonlinear and dynamical quantum plasma problems are readily handled. Moreover, the issue of gauge invariance is treated. Applications (e.g. ultra-dense laser compressed targets and their diagnostics), possible extensions, and future improvements of the presented quantum statistical model are discussed.


Introduction
Quantum kinetic theory has a long history. In many respects, it all started with the seminal paper by Wigner in 1932 [1], see also Refs. [2,3], and the later developments of Moyal [4]. While the approach of Wigner has the advantage of being of interest for the interpretation of quantum mechanics [5], and also for the development of quantum optics (for an overview, see e.g. [6]), detailed calculations of material properties in condensed matter systems have relied to a large extent on either semiclassical techniques [7], in which the collisional operator in Boltzmann's equations involves quantum transition probabilities, or Green's function techniques [8,9], as well as diagrammatic techniques [10]. The theory of Baym and Kadanoff, as well as the works of Keldysh [11,12], has been successful in dealing with certain quantum transport phenomena. The theory contains memory effects (nonlocal terms, both in space and time), has a straightforward interpretation in terms of the different Green's functions, and the theory works well even on timescales shorter than the typical relaxation time of the system in question. However, the gap between classical plasma physics and quantum transport theory does not seem to have been bridged, probably due to reasons of formalism as well as a difference in application of the respective models. Moreover, while the Kadanoff-Baym equation gives a very good description of certain systems, it is perhaps not well-suited to some of the future applications of quantum kinetic theories, such as high intensity laser-plasma interactions [13], high energy density physics [14], and nonlinear collective quantum problems [15,16,17,18,19].
In particular, the field of quantum plasmas has recently attracted, a perhaps unexpected, interest in the field of laser plasmas [20,21,22], where high density ionized plasmas can be created in the laboratory. Moreover, the event of nano-devices and technology on sub-micron scales, such as quantum dots [23,24,25] and plasmonic components [26,27], has sparked the interest of many researchers of analyzing the dynamic and nonlinear properties of such systems. A recent result is that quantum effects in plasmas can be important in parameter regimes that for a long time have been considered purely classical [28].
The above discussion is mainly related to the statistical and dispersive behavior of unmagnetized quantum plasmas [29]. However, one intrinsic non-classical property of quantum systems is the spin. The magnetization that follows from the intrinsic spin, as well as that of orbital angular momentum, is of course the foundation for many important material properties [30]. Investigations of such condensed matter systems are often directed towards equilibrium properties, although the nonlinear dynamics of magnetization is sometimes interest and probed using the Landau-Lifshitz-Gilbert equation [31]. There are a variety of different physical systems where the spin can be of importance, such as metal alloys and semiconductors material for memory use [31], cold atom gases [32], and high density and high field astrophysical plasmas [14], to mention a few. Collective effects originating in the plasma particle species spin has therefore recently become an active field of research for fully ionized systems (see e.g. [19,33] and references therein), in particular in the nonlinear regime, where spin solitons [34] and ferromagnetic behavior in plasmas can be found [35]. Many of the studies presented in the literature have so far been of a theoretical nature, but it is not difficult to envision future applications to e.g. plasmonic devices [26] or femtosecond physics [36].
For the purpose of connecting classical plasma physics to the evolution of nonequilibrium quantum systems, utilization of quasi-distributions is of great value. First of all, the interpretation of the quasi-distribution function using ensemble averages of observables is in direct analogy with the classical case. It is even possible to directly construct a quasi-distribution, such as the Wigner function from measurements [37] (with the only information loss being the initial phase). Second, the quasi-distribution evolution follows from the quantum Liouville equation for the density operator, and gives a quantum analog of the Vlasov or Boltzmann equation. This may also render a quantum kinetic theory for the quasi-distribution function useful for adaption of classical numerical codes to the quantum regime. There are of course infinitely many ways to construct a quasi-distribution function, giving certain elementary requirements (see next section). However, a few quasi-distribution functions are more prominent in the literature than others. The best known quasi-distribution function is probably the Wigner distribution [1], but there are many others frequently used. In short, different definitions correspond to different operator ordering, so depending on the application different definitions are natural. For example, when considering optical coherence normally ordered operators occurs naturally and hence the Glauber-Sudarshan P-distribution [38,39] is a convenient choice. On the other hand, anti-normal ordered operators i.e. the Q-function or the more general Husimi function [40] are useful when dealing with quantum chaotic systems. For reviews of the subject see for example Refs. [5,41,42].
In this paper, we will construct a quasi-distribution function for a particle with spin-1/2 as a combination of a Wigner distribution for the position and momenta and the Q-function for the spin degree of freedom. Moreover, a quantum kinetic equation giving the evolution of this scalar distribution function, in the mean field or Hartree approximation, will be derived and applications to magnetized systems will be presented. A discussion of possible future applications and research directions will also be given.
The structure of the paper is as follows. In Sec. 2 we give a short overview of different quantum quasi-distributions. In Sec. 3 we consider the evolution equation for a density matrix for a spin-1/2 particle in an external electromagnetic field. In Section 4 we go on to derive a combined transformation for the phase space and spin variable. This transformation then renders an evolution equation for the system in extended phase space (x, p,ŝ) which is derived in Section 5. The extension to the mean field approximation is reviewed in Section 6 and in the following section we calculate the thermodynamic equilibrium density matrix for a set of N noninteracting particles. In Section 8 we consider the evolution equation in the long scale length limit and compare our results to previous semi-classical kinetic descriptions in the literature. In Section 9 we consider the linear solutions to the derived equations. Section 10 is devoted to a discussion of gauge properties and the fully gauge invariant evolution equation is presented. Finally we summarize the main results and discuss future development and applications in Section 11.

Historical note
Following the success of the classical theory of non-equilibrium statistical mechanics, it was natural to seek a similar theory for quantum systems in the late 20s and early 30s. However, while the classical Liouville equation generates trajectories in phase space as in a classical Hamilton-Jacobi theory, we in the quantum realm have to consider the Heisenberg uncertainty principle. This will not allow us to describe, as in classical systems, precise trajectories, but rather "smeared out" paths in what would be the corresponding phase space. Indeed, the attempts by de Broglie, Bohm, and others to give a close-to classical interpretation of the Schrödinger equation by Hamilton-Jacobi theory shows that, if one is inclined to stick to this interpretational scheme and extend this to statistical interpretations, one has to consider the wave function rather as an ensemble of (nonclassical [43,44]) trajectories (a similar conclusion can be drawn from path integral [45] as well as Ehrenfest techniques [46]), satisfying certain initial and boundary conditions. Thus, the introduction by Wigner of a quasi-distribution function (see below) was a natural step in the direction of relating measurements to classical transport theory. This is perhaps most obvious in the field of quantum optics, where phase space techniques since long has been widely used. Three main definitions of quasi-distributions can be found in this field, namely the Wigner function [1], the Husimi (or, equivalently, the Q-) function [40], and the Glauber-Sudarshan P-distribution [38,39]. Below we will give a short summary of some of the properties of the first two types of quasi-distribution functions (the P-distribution will not be used in the present work).

Basic requirements
Some basic requirements can be imposed on a quantum probability distribution function in phase space, in order for it to have a reasonable interpretation [44,47]. We denote the quantum (quasi-)distribution function by f (x, p) (for the moment, we drop the explicit time-dependence for notational convenience) for a given quantum stateρ of the system. Then the marginal distribution functions x|ρ|x and p|ρ|p should be related to f (x, p) according to and respectively. Moreover, we should require that the distribution function is positive definite, i.e.
However, it can be shown that the conditions (1)-(3) is not sufficient to uniquely determine a suitable quantum distribution function in phase space. In fact, Cohen [48] has shown that there are infinitely many function f (x, p) satisfying (1)-(3). A more complete list of properties that are desirable is found in [42], where expect for the three properties above, the additional properties that the distribution function is real, bilinear in the wave function and that the distribution functions for eigenstates of the Hamiltonian form a complete and orthogonal set. In fact, it can be shown that in general one cannot find a distribution function that satisfies all of (1)-(3) simultaneously, if one requires the distribution function to be bilinear in the wave function [75].
Though the above conditions are important when it comes to interpreting the distribution functions a perhaps more important condition is that it should be possible to calculate the expectation value of any operator. This condition is important since it means that all physically relevant information is included. To calculate the expectation value one first map the operator to the corresponding phase space functionÔ = O(x,p) → O(x, p), using the Weyl-correspondence, and then calculate the phase space average weighted by the distribution function The mapping from the operator space to phase space depends on which distribution function is used (see Ref. [42] for a details). Below we will collect the properties of two distribution functions of interest in our context, the Wigner distribution [1] and the Husimi function [40] (or Q-function [6,47]). These are also perhaps the most frequently encountered quantum probability distribution function in the literature (see Ref. [6] and [47] for further references and other prominent distributions used in quantum optics, such as the P-distribution of Glauber [38] and Sudarshan [39], and their interrelations). 3

The Wigner function
The Wigner function for a quantum stateρ is defined as the Fourier transform of the two-point correlation function (i.e., density matrix). Thus, we accordingly have Through this definition of the Wigner function, we see that it satisfies the marginal distribution requirements (1) and (2). However, it does not satisfy the positivity criteria (3). The latter property then prevents a probability distribution interpretation. However, the negativity of the Wigner function is limited in the sense that the proper number density in physical space is n(x) = d 3 p f W (x, p) which is thus always positive. For a pure stateρ = |ψ ψ|, this definition gives One of the important properties of the Wigner function is that it cannot have too sharp peaks, expressed by a result of the noncommutativity between coordinate and momentum operators. The time evolution for the Wigner function in an external (analytic) potential V(x, t) is given by where the sin-function is defined in terms of its Taylor expansion in the case of analytic potentials, and we have used the indices x and p on the ∇ to denote its operation in phase space. To find the phase space function that corresponds to a given operator we must first express the operator in Weyl order [2], i.e. express in symmetric products ofx i and p i , i = 1, 2, 3, using the commutation relations and then substitutex i → x andp i → p. For example, calculating the average of the operatorx ip j , we havê where δ i j denotes the Kronecker delta function, and hence

The Husimi function
The Husimi function (see (11) below) is based on minimum uncertainty wave packets, and it does not satisfy (1) and (2) but is positive definite (thus satisfying (3)). As will be seen below, this allows probability distribution interpretation of the Husimi function; however, it gives a different greater uncertainty measure than expected through naive application of the Heisenberg uncertainty relation. These properties can be immediately understood from the following definition. For a given Wigner function, the Husimi function can be obtained through a Gaussian smoothing as where the parameter d sets the scale of the smoothing. While the Husimi function is positive definite, and produce the correct expectation values of observables, it satisfies an indeterminacy relation of the form for a quantum state (the latter being satisfied by the Wigner function). This results is due to the smoothing introduced in the definition of the Husimi function. The Husimi function does not give the probability for the particle to be at a certain phase space position, but rather the probability to find the particle in the minimum uncertainty state centered around the phase space point in question [49]. Introducing minimum uncertainty states |x 0 , p 0 which satisfies ∆x 2 ∆p 2 = 2 /4 one can write the Husimi function as However, as mentioned above, it can still be used to calculate any observable, but the operator ordering rule is more complicated than in the Wigner case and we will not consider this further here. The evolution of the Husimi equation can be found from (8) and the definition (11). It is fairly complicated (see [42]) and it is more convenient to compute the evolution of the Husimi distribution function by evaluating the Wigner function for all times through (8). This said, we note that although the evolution equation for the Husimi function is more complicated than the corresponding equation for the Wigner equation, it is sometimes the convenient choice. One such example is when considering chaotic system in which the phase space distribution function becomes very complicated. The Husimi function, being a Gaussian average, may then behave more regularly (see, e.g., Refs. [50,51,52]).

Quasi-distribution functions for spin
Similarly to the case for phase space it is possible to construct quasi-distribution functions for the spin degree of freedom. This has been done already in the 1950's by Stratonovic [53]. Later on the spin quasi-distribution functions were further developed and were applied to problems related to calculating correlation between spins [54,55,56]. The spin quasi-distribution function has also been discussed in connection with quantum scattering problems [57].
As in the case of the regular phase space variables x and p, there is no unique way to introduce a spin quasidistribution function. Scully and Wódkiewicz [58] give a very good review of the many different choices that can be considered. There are at least three different methods for defining spin distribution functions: delta distributions, distributions based on coherent states (Q and P) and Stratonovic distribution functions. However, the different outcomes of these choices overlap.
In this paper we will consider only the Q-function for spin which is defined as where |s is the state which has spin up in the direction of the unit vector s=s(θ, ϕ) often called a spin-coherent state [59,60]. Note that this is analogous to the definition of the Q-function in position/momentum space, the latter given by Eq. (14). As for the Husimi or Q-function in phase space, this distribution does not give the correct marginal distributions. This means that integrating over the ϕ angle does not leave the correct distribution function for the θ variable. However, it still contains all the information about the system and it can be used to calculate the expectation value of any observable, just like in the density matrix formalism. The mapping between spin operators and the corresponding spin-space functions will be considered in detail in Section 4. The main reason for choosing to work with this particular distribution function for the spin variable is that it is a function on the unit sphere and hence resemble the classical picture of a dipole moment, in fact, the evolution equation in the long scale length limit (see Eq. (63)) is almost identical to an equation derived previously from a semiclassical treatment of the spin [61]. The Q-function for the spin is also nonnegative which may be desirable in some cases.

The density matrix description
In order to derive our phase space model we here start from the density matrix description for a spin-1/2 particle. The basis states we will use are |x, α = |x ⊗ |α where |x is the state with position definitely at position x and |α is the state with spin up (α = 1) and spin down (α = 2) along the axis of quantization, which we here take to be in the z-direction. The density matrix in this basis is then where p i is the probability to have state ψ i and the greek letters denotes the spin indexes. Here ψ(x, 1) and ψ(x, 2) gives, respectively, the probability amplitude to have spin up and spin down. The Hamiltonian for a particle in an external electromagnetic field is given bŷ where m is the mass of the particle, q is the charge (for an electron q = −e < 0 where e is the elementary charge), µ is the magnetic moment of the particle, which for electrons is given by the (signed) Bohr magneton µ e = −e /(2m e ), A and V are the electromagnetic potentials, and B = ∇ x × A is the magnetic field. σ is the vector containing the three Pauli matrices as its components. With the axis of quantization in the z-direction they are given by We will use the notation denotes the component on row α and column β of σ (x) and similarly for the σ (y) and σ (z) matrices. The evolution equation for the density matrix can be derived from the Schrödinger equation for the wave function and its complex conjugate, giving the von Neumann equation Using the basis described above and the Hamiltonian (17) we obtain where we have used the Coulomb gauge ∇ x · A = 0. In general, the evolution equation of the diagonal terms ρ(x, α; y, α), α = 1, 2 are coupled via the off-diagonal terms. However, for static fields it is possible to obtain two decoupled equations for the diagonal elements, by orienting the axes so that the magnetic field is in the direction of the axis of quantization [62].

The Wigner and Q transformation
The Wigner transformation for a spin-1/2 particle is given by where we have emphasized that for a particle with spin the Wigner transform must be taken for each spin matrix element of the density matrix separately. The Wigner transform of the spin density matrix has been calculated previously [57,62,63]. One approach is to consider the different components of the Wigner matrix W(x, p, α, β), for α = 1, 2 and β = 1, 2 and derive evolution equations for W(x, p, 1, 1) and W(x, p, 2, 2) which, as the for the density matrix, are in general coupled via the off-diagonal terms [62,64]. Another approach is to define a quasi-distribution function for the spin degree of freedom. This can, as have been discussed above, be done in a variety of different ways [55,56,58]. The way which is a direct generalization of the Wigner function is to consider two different spin components in two arbitrary directions s 1 and s 2 , corresponding to the two operators σ 1 and σ 2 , see Ref. [55]. Since the two operators in general do not commute, the values of s 1 and s 2 cannot be known simultaneously. This manifests itself in that the Wigner function W(s 1 , s 2 ), as for the corresponding case of position and momentum, can take on 6 negative values. Another possible choice of distribution function (corresponding to anti-normal operator ordering) is the Q-function. In the position/momentum space, this distribution function is the Gaussian averaged Wigner function and due to this it is positive definite. In optics, the Q-function can be measured directly [6]. The spin Q-function [58] gives the probability to measure the spin in a given direction and it is this we will here use to describe the spin degree of freedom.
To derive an evolution equation for the extended phase space distribution function f (r, p,ŝ), whereŝ is a unit vector (not an operator), we impose the following properties: should give the probability density to find the particle at position r with spin up in the direction ofŝ and, similarly, should give the probability to have momentum p and spin up in theŝ direction, a direct extension of the marginal distribution conditions (1) and (2). In order to derive the distribution function in the extended phase space we note that for a state ψ(x, α) we have the probabilities |ψ(x, 1)| 2 (|ψ(x, 2)| 2 ) to measure spin up (spin down) in the z-direction.
The corresponding density matrix is given by ρ(x, α; y, β) = ψ(x, α)ψ * (y, β). We can then write the probability to measure spin-up in the direction of the unit vectorŝ as where we have defined the (Hermitian) operatorP and where δ αβ denotes the Kronecker delta. As an example we consider the the probability to measure the spin in the directionŝ = −ẑ and we get get 2 α,β=1 as we expect (note that measuring spin-up in the −z direction is equivalent to measure spin-down in the z direction). The generalization to a statistical distribution of states is straightforward. Using the Wigner transform for the position and momentum and the spin transform discussed above we obtain the function which have the properties (21) and (22) stated above. The function f may also be written as where W is the 2 × 2 matrix with elements W(x, p, α, β). The normalization of the extended Wigner function is given by Hence we obtain a distribution function which is normalized over the allowed spin values if we redefine the operator in Eq. (24) asP In the Wigner formalism without spin, the density matrix is transformed into the Wigner function and the operators are transformed into phase space functions. For an operatorĝ = g(x,p) the corresponding phase space function is given by It can also be obtained by using Weyl ordering as described in subsection 2.3. With this function the expectation value of the operator is calculated as a phase space integral where f (x, p) andρ are related via a Wigner transform. In analogy with this, for a given operatorĥ acting on the spin degree of freedom, we define the corresponding spin-space function where h(α, β) denotes the (α, β) component of the operatorĥ. With this definition the expectation value of the operator is now calculated as an integral over the possible spin directions according to where we have used d 2ŝ s a s b = (4π/3)δ ab and σ i αβ σ i γδ = 2δ αδ δ βγ − δ αβ δ γδ . In doing the calculation above we have also used that the general form of a spin-operator isĥ = aI + b · σ where a and b may be dependent of position and momenta. Note that the definition of the spin space function, Eq. (32) implies that the spin operator σ is related to the spin unit vectorŝ according to σ → 3ŝ.
For operators depending on both the position and momentum and the spin degree of freedom the corresponding extended phase space function is obtained by doing both the transformations (30) and (32). The operator (29) can also be writtenP(ŝ) = |ŝ ŝ|, where |ŝ is the spin coherent state [59,60]. The definition (27) is then seen to coincide with the definition of the spin Q-function, see Eq. (15). The function f (x, p,ŝ) is hence a combination of a Wigner function in the phase space variables and the Q-function for the spin.

Equivalence with the density matrix formalism
The construction above contains the same information as the density matrix, and the distribution function can be used to calculate the expectation value of any observable. A more direct way to see the equivalence is to note that, for a given distribution function f (x, p,ŝ), it is possible to obtain the corresponding Wigner matrix as From this it is the possible to obtain the density matrix by taking the inverse Wigner transform, (see for example Ref. [76]).

Evolution equation
To derive the evolution equation for f (x, p,ŝ), the Wigner transform of Eq. (19) is calculated with the result (assuming that the fields and potentials are analytic functions) where functions of an operator is defined by its formal Taylor expansion and the left (right) arrow above the differential operators indicate that they act on the functions on the left (right). If the potentials have discontinuities the above equation can instead be written explicitly in the form of an integro-differential equation. Next we multiply by [δ βα + s · σ(β, α)]/2 and sum over α and β. The operators acting on the left hand side and the first four terms on the right hand side commute with the Pauli matrices and for these we obtain W(x, p, α, β) → f (x, p,ŝ). For the last two terms we use the property and also that σ * αβ = σ βα . After some straightforward calculations we get the evolution equation for the extended Wigner function , p,ŝ). (38) An advantage of writing the evolution equation in this form is that we may Taylor expand the trigonometric function to sufficient order in to obtain the semi-classical limit directly. Next we make a variable transformation in the evolution equation. The canonical momentum p is related to the velocity by v = (p − qA)/m. Changing variables from x, p and t to x, v, t we get where ∇ xi = ∂/∂x i and ∇ vi = ∂/∂v i . For the time derivative we get We can then write the full quantum-kinetic equation (38) as displaying the classical and semiclassical terms more explicitly on the left-hand side of the equation. We note that the terms on the right-hand side all are higher-order derivative corrections.

Many-particle evolution equation
So far we have only considered one particle in an external electromagnetic field. To make a straightforward generalization to an N-body system we will consider the mean field approximation. In order to keep things simple we will neglect effects due to spin statistics (antisymmetry of the wavefunction). To a certain degree such effects can be incorporated by choosing appropriate initial conditions, see the next section. Introducing the many-particle density matrixρ 1...N it will satisfy the von Neumann equation The N-body HamiltonianĤ (N) in general includes interactions between the particleŝ whereĤ i = (p i − qA 0 (x i )) 2 /2m + qV 0 (x i ) is the Hamiltonian for particle i and contains the kinetic energy and the interaction with an external electromagnetic field (V 0 , A 0 ), andĤ i j is the interaction between particle i and j which we assume to be the full electromagnetic interaction between the particles. The interaction is hence obtained by solving Maxwell's equations. Following Ref. [66] we introduce the reduced density matrix in the thermodynamic limit (N, where V is the volume of the system and the trace includes summing over the spin degree of freedom. The normalization is given by 1 V s Tr 1,...sρ1...s = 1. (44) Note that this means that for the diagonal elements of the one-particle reduced density matrix ρ 1 (x, x) is proportional to the probability density to find any one of the N particles in position x independently of the positions of all the other particles. Expectation values of an s-body operator is given by The evolution equations for the reduced density matrix is given by the BBGKY-hierarchy [69] i where H (s) is obtained by changing N → s in Eq. (42). Considering the first order equation and introducing the two-particle correlation asρ 12 =ρ 1ρ2 +ĝ 12 we may write this as where theĤ MF = Tr 2Ĥ12ρ2 is the mean field which is found by solving Maxwell's equations self-consistently. The effects of particle-particle scattering is included in the correlation operatorĝ 12 . This, in turn, satisfies an equation that is coupled to the three particle correlations and so on. Here we will be mainly interested in the collective effects of the plasma and hence we will neglect the right hand side of Eq. (47), i.e. use the Hartree approximation. In order to include self-energy effects and ionization/recombination it is necessary to keep higher order correlations.
Comparing Eq. (47) with the corresponding equation for a single particle in an external electromagnetic field, Eq. (18), we note that they are formally the same. Thus in order to include systems of N-particles in the Hartree approximation we hence need to assume that the fields in the evolution equation are the self consistent fields and then keep in mind that the density matrix is now normalized according to Eq. (44). However, since we will not pursue the issue of the quantum BBGKY-hierarchy further, we may redefine the one-particle distribution function so that it has the normalization Trρ 1 = n 0 , so that for example x|ρ 1 |x = n(x) gives the mean density of particles at position x.
The mean field interactionĤ MF is obtained by coupling the equation to Maxwell's equations. The expression for the charge and current densities are then where we thus obtain a magnetization current contribution due to the spin (see Ref. [67]).

Thermodynamic equilibrium density matrix
As an example we calculate the extended phase space distribution function for a system of N non-interacting particles in a constant magnetic field which are in thermodynamic equilibrium at temperature T . Assuming that the magnetic field is B = B 0ẑ and using the Landau gauge A = (−yB 0 , 0, 0) we can obtain the eigenstates of the Hamiltonian, Eq. (17), where φ n is the n'th harmonic oscillator wave function given by where H n are the Hermite polynomials [68]. Furthermore, we have introduced the spinors χ a (α) which satisfies σ z χ a (α) = aχ a (α) with a = ±1. The energy levels corresponding to Eq. (51) are given by Note that the energy is independent of the momentum in the x-direction so that the energy levels are degenerate. The thermal equilibrium density matrix at temperature T is given bŷ where k B is Boltzmann's constant and the partition function is given by Considering the one-particle density matrix, it can be written ρ(x, α; y, β) = p x ,n,p z ,a p p x ,n,p z ,a ψ p x ,n,p z ,a (x, α)ψ * p x ,n,p z ,a (y, β)χ a (α)χ † a (β), where the probability for the state with quantum numbers (p x , n, p y , a) is given by where µ c is the chemical potential. The Wigner transform of the density matrix for an harmonic oscillator has been calculated by Ref. [47]. Using their result we can calculate the Wigner transform to be where L n denotes the Laguerre polynomials [68]. Calculating the spin-transform, Eq. (27), of this and also changing variables to v = (p − qA)/m we finally obtain f (x, v,ŝ) = n,a n 0 (−1) n 2π(2π ) 3 1 + a cos θ s e (E n,pz ,a −µ c )/T + 1 exp where we have also multiplied by n 0 to obtain the chosen normalization (see the previous section). Note that the argument appearing in the exponential and the Laguerre polynomials is just the kinetic energy of the motion perpendicular to the magnetic field. However, as opposed to the classical case, v x and v y are non-commuting quantities and cannot be determined simultaneously. To verify that Eq. (59) indeed is a solution to the Wigner equation we note that for stationary solutions in the given choice of magnetic field the equation can be written and we see that in fact any spatially homogenous function solves this. The expression (59) contains Landau-quantization, spin splitting of energy states and Fermi-Dirac statistics. For cases where the chemical potential µ c is large, and the difference between nearby Landau levels is smaller than the thermal energy, the velocity distribution approaches the classical Maxwellian. An important quantum mechanical result that remains in this limit is that the probability distribution of the spin up and down populations scales as 1 + cos θ s and 1 − cos θ s , respectively. Thus, in the above regime the distribution can be approximated by where F ± are Maxwellian distributions. The ratio F + /F − in thermodynamic equilibrium is For small chemical potential, when Fermi-Dirac statistics applies to F ± , we can still have the form (61), but the ratio F + /F − is velocity dependent. Actually, even in the absence of thermodynamical equilibrium Eq. (61) is the most general time independent, homogenous expression for the distribution function in a constant magnetic field.

Long scale length limit
To obtain the long scale limit we Taylor expand the trigonometric operators to order , which applies if the characteristic scale lenghts are longer than the thermal de Broglie length. Thus we henceforth neglect higher order terms in , such as 2 4 ∇ 2 x V(x) · ∇ 2 p f . The evolution equation then becomes Making a variable change from x, p and t to x, v, t according to Eqs. (39) the second term in Eq. (39a) above will combine with other terms in Eq. (62) to produce the magnetic field term in the Lorentz force. The last term in the time derivative Eq. (39c) will combine with the gradient of the scalar potential ∇V to produce the electric field E = −∇V + ∂ t A. The evolution equation then takes the form Note that the last term contains derivatives both with respect to the velocity v and the spinŝ. The equation above has already been studied in [61] with the last term missing due to semi-classical approximations. It is there shown to give rise to new oscillation modes due to the anomalous magnetic moment of the electron. A similar equation to has also been studied in [72] where it is investigated whether spin may be of importance in magnetic confined fusion experiments. The difference between Eq. (63) and the semi-classical case (with the last term missing) is due to the fact that the quantum mechanical probability distribution is always spread out, as follows from Eq. (29). In order to demonstrate this we consider the distribution function for a single particle, that at a time t has a given spin state, pointing in the directionê, whereê is a unit vector. As follows from Eq. (29), the corresponding distribution function, which is smeared out in spin space, can be written f (x, v,ŝ) = F(x, v)(1 +ê ·ŝ)/4π. If we average over all spin directions, the last term combine with the magnetic dipole term according to 1 4π where we stress that 2/3 of the contribution comes from the latter term. The full evolution equation for this reduced distribution function can be written where B = |B| and we have allowed the spin directionê to be slowly variying, following the variations of the magnetic field direction, i.e.ê =B(x, t). The equation above is usefull when the spin state of each particle is conserved for a sufficiently long time.
For a semi-classical treatment of the spin we would expect that the probability to measure the spin in the direction s given that the spin is inê-direction is given by f cl (ŝ) = δ(ŝ −r). However, as can be seen from the above equation, the classical limit of a particle with spin in ther-direction is not a particle with definite magnetic moment in thê r-direction but a statistical distribution of spins in all directions (exceptŝ = −r which has zero probability). For a magnetized electron plasma, the magnetization is given by ∇ × M = ∇ × µ σ where the expectation value is taken with respect to the spin degree of freedom. In the quantum model developed here the spin is given by The factor 3 will account for the fact that probability to find the spin in a certain direction is smeared out over the whole unit sphere. In a classical treatment of the spin variable the corresponding integral contains no factor 3, but instead the distribution function is a delta function of the spin and hence the same result can be obtained. The latter is the model used in Ref. [61].

Spin induced damping of Alfvén waves
As an example of the usefulness of Eq. (63) we will consider shear Alfvén like waves in the linear limit. First we divide the variables as f = f 0 + f 1 and B = B 0 + B 1 , with in which case the linearized electron equation can be written

Generalized L and R waves
As a further example in the linear regime we linearize the full evolution equation, Eq. (40). To simplify the algebra we look at waves propagating parallel to a background magnetic field. Following the steps from the last subsection it is straightforward to derive the dispersion relation Figure 1: The dependence of the normalized growth rate γ = γ cl − γ sp /(γ cl + γ sp ) on the normalized wavenumber k = k z v thi /ω ci , for n 0 = 10 24 Thus, to lowest order in we haveẼ so that (cf. (40)) The gauge invariant Wigner function has a modified Weyl correspondence which is well suited for calculating fluid moments. In order to obtain the phase space O(x, v) function which corresponds to an operator O(x,v), all products of the operatorsx andv ≡ [p − qA(x)]/m are first ordered in a symmetric form using the commutation relationx and p and then the substitutionx → x andv → v is taken (details can be found in [65]).

Summary and discussion
In the present paper we have derived an evolution equation, Eq. (40), for a quasi distribution function of electrons, based on a Wigner transformation of the density matrix, together with a spin operator contracting the 2 × 2 Wignermatrix to a scalar function f (x, p, s). The free current and the magnetization can be directly computed from the quasi distribution function, and hence Eq. (38) (or the gauge invariant alternative, Eq. (83)) together with Maxwell's equations with the sources (49), and (50), form a closed set. The present theory has the advantage to include the full quantum dynamics in a single equation, and provide an immediate path between the classical and quantum descriptions. For macroscopic scale lengths longer than the characteristic de Broglie wavelength, the kinetic equation greatly simplifies. In particular, the semi-classical kinetic theory put forward in Ref. [61] is recovered, but with some small but significant deviations as shown in (63). The difference between the semi-classical theory and our result follows from the smeared out probability distribution of the spin.
In order to illustrate the theory, examples of linear wave propagation solving the full quantum theory, Eq. (40), as well as the long wavelength limit, Eq. (63) are given. An interesting result is that the wave damping can be much affected by the non-classical terms even in a supposedly classical temperature and density regime, although the real part of the wave frequency then is always well approximated by the classical Vlasov theory. The reason is that the spin terms give raise to new types of wave particle resonances.
Proper initial conditions for the quasi distribution function can be found by computing the Wigner transformation of the density matrix in the thermodynamical ground state. The result for the simple but important special case of a magnetic field is given, see Eq. (59), in which case the energy levels are Landau quantized and split due to the two spin states [74].
The present quantum theory can be used for a broad range of parameters, extending the applicability to regimes of high densities, strong magnetic fields, low temperatures and short scalelengths, that are not covered by the classical Vlasov equation. However, there is still much room for improvements. In particular, when the strong coupling parameter Γ is increased, collisional effects become important [71]. A schematic view of the different plasma regimes is given in Fig. 2. Furthermore, the present theory does not account for relativistic effects that is crucial in e.g. laser plasma interaction. Removal of these restrictions, as well as a more complete evaluation of the present theory, constitutes interesting projects for future research. Figure 2: Various plasma regimes in the temperature-density parameter space are illustrated. The dotted line is given by the strong coupling parameter Γ = E p /k B T = 1, where E p = e 2 n 1/3 0 /4πε 0 is the potential energy due to the nearest neighbor. For larger densities this parameter is repaced by Γ F = E p /k B (T + T F ), since the average kinetic energy of the particles is given by the Fermi energy rather than the thermal energy. The curve Γ F = 1 is illustrated by the dashed curve, and the strong coupling region, which is shaded, occurs below this line. In this region our model is not directly applicable, since collisions has not been taken into account. For comparison we have also drawn the lines ω p /k B T (the dotted grey line) and the line T F /T (the dotted-dashed grey line) which measures the importance of wave function dispersion and the Fermi pressure, respectively. As a rough estimate, the quantum regime is below either of these lines. Note, however, that spin effects can sometimes be important even above these lines [28]