Quantum Kinetic Theory of Spin Polarization of Massive Quarks in Perturbative QCD: Leading Log

TL;DR

This paper derives a quantum kinetic equation describing the spin polarization evolution of massive quarks in a QCD plasma, valid for arbitrary quark masses much larger than the Debye mass, at leading logarithmic order in perturbative QCD.

Contribution

It introduces a universal quantum kinetic equation for massive quark spin dynamics in perturbative QCD at leading log order, applicable for arbitrary large quark masses.

Findings

Collision terms are of order α_s^2 log(1/α_s).

Equation valid for quark masses much larger than the Debye mass.

Provides a framework for studying spin relaxation in QCD plasma.

Abstract

We present the quantum kinetic equation for spin polarization of massive quarks in leading log order of perturbative QCD, which describes time evolution of the spin density matrix in momentum space of a massive quark interacting with a background QCD plasma. We find that the time evolution operator of the spin density matrix, or the quantum kinetic collision terms, are universally of order $\alpha_s^2\log(1/\alpha_s)$ in terms of the QCD coupling constant $\alpha_s=g^2/(4\pi)$. Our quantum kinetic equation is valid for an arbitrary quark mass $m\gg m_D\sim gT$, where $m_D$ is the Debye mass, and can be used to study relaxation dynamics of spin polarization of massive quarks in perturbative QCD regime.

Tables

Table 1. Table 1: The coefficient functions j n L / T subscript superscript 𝑗 𝐿 𝑇 𝑛 j^{L/T}_{n} ( η p = 1 2 ​ ln ⁡ E p + p E p − p subscript 𝜂 𝑝 1 2 subscript 𝐸 𝑝 𝑝 subscript 𝐸 𝑝 𝑝 \eta_{p}=\frac{1}{2}\ln\frac{E_{p}+p}{E_{p}-p} is rapidity).

$j_0^L=\frac{p}{E_p}$	$j_0^T=\frac{{\eta }_p}{2}$
$j_{1/2}^L=-\frac{p}{2E_p}$	$j_{1/2}^T=-\frac{p{E}_p}{4m^2}$
$j_1^L=\frac{p^3}{3E_p^3}$	$j_1^T=\frac{{\eta }_p}{2}-\frac{p}{2E_p}$
$j_{3/2}^L=-\frac{2p^3}{E_p^3}$	$j_{3/2}^T=-\frac{p^3}{4m^2{E}_p}$
$j_2^L=\frac{p^5}{5E_p^5}$	$j_2^T=\frac{{\eta }_p}{2}-\frac{p}{2E_p}-\frac{p^3}{6E_p^3}$

Equations

$$\frac{\partial \hat{\rho }}{\partial t}=-\Gamma \cdot \hat{\rho }-\frac{i}{\hslash }[H_{\mathrm{eff}},\hat{\rho }],$$

$$H_{\mathrm{eff}}=-\frac{\hslash }{2}\mathbf{\sigma }\cdot (\mathbf{\omega }+Q\mathbf{B}),$$

$$\hat{\rho }({\mathbf{x}}_r,{\mathbf{p}}_r)=\int \frac{d^3{\mathbf{p}}_a}{(2\pi {)}^3}{e}^{i{\mathbf{p}}_a\cdot {\mathbf{x}}_r}\hat{\rho }({\mathbf{p}}_r,{\mathbf{p}}_a).$$

$$\hat{\rho }=\int \frac{d^3\mathbf{p}}{(2\pi {)}^3}\hat{\rho }(\mathbf{p}),$$

$$\hat{\rho }(\mathbf{p})=\sum _{s,s^{\prime }=\pm }|\mathbf{p},s⟩{\rho }_{s,s^{\prime }}(\mathbf{p})⟨\mathbf{p},s^{\prime }|,$$

$$\hat{\rho }(\mathbf{p})=\sum _{s,s^{\prime }=\pm }{\xi }_s(\mathbf{p}){\rho }_{s,s^{\prime }}(\mathbf{p}){\xi }_{s^{\prime }}^{†}(\mathbf{p}).$$

$$\mathbf{S}(\mathbf{p})=\mathrm{Tr}(\mathbf{S}\hat{\rho }(\mathbf{p}))=\frac{1}{2}\mathrm{Tr}(\mathbf{\sigma }\hat{\rho }(\mathbf{p})).$$

$$\hat{\rho }(\mathbf{p})=\frac{1}{2}f(\mathbf{p})+\mathbf{S}(\mathbf{p})\cdot \mathbf{\sigma }.$$

$$\mathbf{S}=\int \frac{d^3\mathbf{p}}{(2\pi {)}^3}\mathbf{S}(\mathbf{p}),n=\int \frac{d^3\mathbf{p}}{(2\pi {)}^3}f(\mathbf{p}).$$

$$\hat{\rho }(t)=⟨U_1(t,t_0)\hat{\rho }(t_0)U_2^{†}(t,t_0){⟩}_A,$$

$$G_{\mu \nu }^{(ij)}(q^0,\mathbf{q})=\int d^3\mathbf{x}dt{e}^{i(q^{0}t-\mathbf{q}\cdot \mathbf{x})}⟨A_{\mu }^{(i)}(\mathbf{x},t)A_{\nu }^{(j)}(0,0){⟩}_A,i,j=1,2(\mathrm{SK}\mathrm{contours})$$

$$G_{\mu \nu }^{(12)}(q^0,\mathbf{q})=n_B(q^0){\rho }_{\mu \nu }(q^0,\mathbf{q}),G_{\mu \nu }^{(21)}(q^0,\mathbf{q})=(n_B(q^0)+1){\rho }_{\mu \nu }(q^0,\mathbf{q}),$$

$$H_I=g\int d^3\mathbf{x}\bar{\psi }(\mathbf{x}){\gamma }^{\mu }{t}^a\psi (\mathbf{x})A_{\mu }^a(\mathbf{x}),$$

$$\gamma^{0}=\left(\begin{array}[]{cc}0&i{\bf 1}_{2\times 2}\\ i{\bf 1}_{2\times 2}&0\end{array}\right)\,,\quad\gamma^{i}=\left(\begin{array}[]{cc}0&i\sigma^{i}\\ -i\sigma^{i}&0\end{array}\right)\,,i=1,2,3\,,$$

$$|\bm{p},s\rangle\sim u(\bm{p},s)=\left(\begin{array}[]{c}\sqrt{E_{p}-sp}\,\xi_{s}(\bm{p})\\ \sqrt{E_{p}+sp}\,\xi_{s}(\bm{p})\end{array}\right)\,,$$

$$\psi (\mathbf{x})\sim \sum _{\mathbf{p},s}u(\mathbf{p},s)a_{\mathbf{p},s}{e}^{i\mathbf{p}\cdot \mathbf{x}}+h.c$$

$$\sum _{\mathbf{n}}\to V\int \frac{d^3\mathbf{p}}{(2\pi {)}^3}.$$

$$\psi (\mathbf{x})=\frac{1}{\sqrt{V}}\sum _{\mathbf{n},s}\frac{1}{\sqrt{2E_{p_{\mathbf{n}}}}}u({\mathbf{p}}_{\mathbf{n}},s)e^{i{\mathbf{p}}_{\mathbf{n}}\cdot \mathbf{x}}{a}_{{\mathbf{p}}_{\mathbf{n}},s},$$

$$A_{\mu }(\mathbf{x},t)=\frac{1}{V}\sum _{\mathbf{n}}{A}_{\mu }({\mathbf{p}}_{\mathbf{n}},t)e^{i{\mathbf{p}}_{\mathbf{n}}\cdot \mathbf{x}}+h.c.,$$

$$\{a_{{\mathbf{p}}_{\mathbf{n}},s},a_{{\mathbf{p}}_{{\mathbf{n}}^{\prime }},s^{\prime }}^{†}\}={\delta }_{\mathbf{n},{\mathbf{n}}^{\prime }}{\delta }_{s,s^{\prime }},$$

$$⟨A_{\mu }^{(i)}({\mathbf{p}}_{\mathbf{n}},t)A_{\nu }^{(j)}({\mathbf{p}}_{{\mathbf{n}}^{\prime }},t^{\prime }){⟩}_A=V{\delta }_{\mathbf{n},-{\mathbf{n}}^{\prime }}{G}_{\mu \nu }^{(ij)}({\mathbf{p}}_{\mathbf{n}},t),$$

$$H_I(t)=\frac{g}{V}\sum _{\mathbf{n},s}\sum _{{\mathbf{n}}^{\prime },s^{\prime }}\frac{1}{\sqrt{2E_{p_{\mathbf{n}}}}}\frac{1}{\sqrt{2E_{p_{{\mathbf{n}}^{\prime }}}}}\bar{u}({\mathbf{p}}_{\mathbf{n}},s){\gamma }^{\mu }u({\mathbf{p}}_{{\mathbf{n}}^{\prime }},s^{\prime })A_{\mu }({\mathbf{p}}_{\mathbf{n}}-{\mathbf{p}}_{{\mathbf{n}}^{\prime }},t)a_{{\mathbf{p}}_{\mathbf{n}},s}^{†}{a}_{{\mathbf{p}}_{{\mathbf{n}}^{\prime }},s^{\prime }},$$

$$|{\mathbf{p}}_{\mathbf{n}},s⟩\equiv a_{{\mathbf{p}}_{\mathbf{n}},s}^{†}|0⟩.$$

$$U_I(\Delta t,0)\approx 1-i{\int }_0^{\Delta t}dt{H}_I^{\mathrm{int}}(t)+(-i{)}^2{\int }_0^{\Delta t}dt{\int }_0^{t}d{t}^{\prime }{H}_I^{\mathrm{int}}(t)H_I^{\mathrm{int}}(t^{\prime })+\cdots ,$$

$$\hat{\rho }(\Delta t)$$

$$\hat{\rho }({\mathbf{p}}_{\mathbf{n}})=\sum _{s,s^{\prime }}{\xi }_s({\mathbf{p}}_{\mathbf{n}}){\rho }_{s,s^{\prime }}({\mathbf{p}}_{\mathbf{n}}){\xi }_{s^{\prime }}^{†}({\mathbf{p}}_{\mathbf{n}}),$$

$${\int }_0^{\Delta t}d{t}_1{\int }_0^{\Delta t}d{t}_2{G}^{(12)}(t_1-t_2)e^{i{q}^0(t_1-t_2)}={\int }_0^{\Delta t}d{t}_r{\int }_{-2t_r}^{2t_r}d{t}_a{G}^{(12)}(t_a)e^{i{q}^0{t}_a},$$

$$\frac{d}{dt}{\rho }_{s,s^{\prime }}(\mathbf{p},t)=g^2{C}_2(F)({\Gamma }_{\mathrm{cross}}+{\Gamma }_{\mathrm{self}\mathrm{energy}}),C_2(F)=\frac{N_c^2-1}{2N_c},$$

$${\Gamma }_{\mathrm{cross}}=\int \frac{d^3{\mathbf{p}}^{\prime }}{(2\pi {)}^3}\frac{1}{4E_p{E}_{p^{\prime }}}\sum _{s^{\prime \prime },s^{\prime \prime \prime }}[\bar{u}(\mathbf{p},s){\gamma }^{\mu }u({\mathbf{p}}^{\prime },s^{\prime \prime })]{\rho }_{s^{\prime \prime },s^{\prime \prime \prime }}({\mathbf{p}}^{\prime })[\bar{u}({\mathbf{p}}^{\prime },s^{\prime \prime \prime }){\gamma }^{\nu }u(\mathbf{p},s^{\prime })]G_{\mu \nu }^{(12)}(E_p-E_{p^{\prime }},\mathbf{p}-{\mathbf{p}}^{\prime })$$

$${\Gamma }_{\mathrm{self}\mathrm{energy}}=$$

Full text

**Quantum Kinetic Theory of Spin Polarization of Massive Quarks in Perturbative QCD: Leading Log **

Shiyong Li***e-mail: sli72@uic.edu and Ho-Ung Yee†††e-mail: hyee@uic.edu

Department of Physics, University of Illinois,

*Chicago, Illinois 60607 *

2019

We present the quantum kinetic equation for spin polarization of massive quarks in leading log order of perturbative QCD, which describes time evolution of the spin density matrix in momentum space of a massive quark interacting with a background QCD plasma. We find that the time evolution operator of the spin density matrix, or the quantum kinetic collision terms, are universally of order $\alpha_{s}^{2}\log(1/\alpha_{s})$ in terms of the QCD coupling constant $\alpha_{s}=g^{2}/(4\pi)$. Our quantum kinetic equation is valid for an arbitrary quark mass $m\gg m_{D}\sim gT$, where $m_{D}$ is the Debye mass, and can be used to study relaxation dynamics of spin polarization of massive quarks in perturbative QCD regime.

1 Introduction

The quark-gluon plasma created in off-central heavy-ion collisions is expected to have a sizable collective fluid rotation that is originated from the large initial angular momentum of the two colliding projectiles. A part of the orbital angular momentum that resides in this fluid motion, or vorticity, will then be transferred to spin angular momenta of quasi-particles by interactions [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]. In equilibrium, the resulting spin-dependent distribution function of quasi-particles can be shown to be equal to what one would have as if the energy of a particle of spin $\bm{S}$ was shifted by $\Delta E=-\bm{S}\cdot\bm{\omega}$ in thermal distribution[13, 14, 15], where $\bm{\omega}={1\over 2}\bm{\nabla}\times\bm{v}$ is the vorticity of fluid***An intuitive derivation of this fact, based on a detailed balance argument with total angular momentum conservation was given in Ref.[16].. Since $\bm{S}\sim{\cal O}(\hbar)$, the net spin angular momentum is quantum in nature, and is parametrically small (by $\hbar$) compared to macroscopic orbital angular momentum of the fluid. The recent experimental observation of spin polarization of $\Lambda$ baryons in off-central heavy-ion collisions [17] confirms the existence of this phenomenon that involves quantum spin of quasi-particles.

In the early stage of heavy-ion collision, the QCD plasma is presumably in its de-confined phase, where quarks and gluons are the basic degrees of freedom. The magnitudes of vorticity and magnetic field that polarize the quasi-particle spin are strongest at such an early stage. Some of the spin polarization of quarks and gluons in this phase can be transferred to that of hadrons after hadronization, which may affect the experimentally observed spin polarization of hadrons. Whether this effect survives hadronic phase dynamics depends on the relaxation dynamics of spin polarization in hadronic phase, as well as many other realistic conditions of heavy-ion collisions [18, 19]. As a first step, a reliable tracking of time-evolution of spin polarization of quarks and gluons within the de-confined phase itself would be a prerequisite in any quantitative theory prediction of spin polarization of observed hadrons. The aim of the present work is to address this problem, at least partly, in leading log order of perturbative QCD (pQCD). This also complements our previous work on the similar question in strongly coupled regime described by AdS/CFT correspondence [20].

In a time-dependent background such as heavy-ion collisions, the spin polarization of quasi-particles would naturally be driven off equilibrium. The time evolution of spin polarization would roughly be a competition between QCD dynamics that tries to relax the spin to equilibrium and the time-variation of backgrounds, such as vorticity and magnetic field, that drives the spin polarization of the system off equilibrium. If the time-variation of background is much slower than the characteristic relaxation time due to QCD interactions (which turns out to be $\tau_{R}\sim(\alpha_{s}^{2}\log(1/\alpha_{s})T)^{-1}$), the system would follow closely the instantaneous equilibrium at each time. On the opposite case, the system would deviate significantly from equilibrium, and the spin polarization should be determined by solving the dynamical equation for time-evolution of spin polarization.

The present work is a small step in formulating such a dynamical equation of spin polarization of quasi-particles in QGP phase in leading log order of pQCD, focusing only on the spin polarization of massive quark, that may be suitable for strange quark or more massive quark species. Specifically, we assume that the mass is of hard-scale, $m\gg m_{D}\sim gT$ ($g$ is the QCD coupling constant). We will see that this justifies a few simplifications we will detail below. Therefore, our results would not be applicable for light (u,d) quarks and gluons, the study of which we defer to a future work.

Schematically, the evolution of spin density matrix $\hat{\rho}$ of massive quark would take a form in linear order as (basically a “Lindblad equation”)

[TABLE]

with a linear relaxation operator $\Gamma$, and the effective one-particle Hamiltonian $H_{\rm eff}$ in a 2-dimensional spin space, that may include vorticity, $\bm{\omega}$, and magnetic field, $\bm{B}$, in somewhat phenomenological way as

[TABLE]

where $Q$ is the electromagnetic charge of quark. In general, the $\Gamma$ should also depend on vorticity and magnetic field in such a way that the equilibrium spin density matrix is given by $\hat{\rho}_{eq}=e^{-H_{\rm eff}/T}\approx{\bf 1}-H_{\rm eff}/T$, at least in linear order in $\bm{\omega}$ and $\bm{B}$. The magnetic field, for example, should modify the wave function and energy spectrum of quark states in the computation of $\Gamma$ that we describe in the following sections. We would expect to have an expansion of $\Gamma$ in small $\bm{\omega}$ and $\bm{B}$ as $\Gamma=\Gamma_{0}+\Gamma_{1}+\cdots$ where $\Gamma_{1}$ is linear in $\bm{\omega}$ or $\bm{B}$. In this work, we present our result for the leading relaxation operator $\Gamma_{0}$, that corresponds to the case of vanishing $\bm{\omega}$ and $\bm{B}$, and the computation of $\Gamma_{1}$ is planned in the future.

Although $\Gamma_{0}$ is not sufficient to describe the spin polarization in a time-varying vorticity and magnetic field, it can still be used to compute spin-related correlation functions in a linear response theory. At least, $\Gamma_{0}$ can describe how a spin density matrix, initially polarized, relaxes to the unpolarized one (the identity operator in spin space), when vorticity and magnetic field cease to exist.

In general, the density matrix $\hat{\rho}$ is defined in the phase space $(\bm{x},\bm{p})$ in addition to spin space. A convenient way to think about it is in the language of Schwinger-Keldysh contour. The position and momentum operators in forward and backward time contours (labeled as 1 and 2 respectively) satisfy the commutation relations $[\bm{x}^{i}_{1},\bm{p}^{j}_{1}]=i\hbar\delta^{ij}$ and $[\bm{x}^{i}_{2},\bm{p}^{j}_{2}]=-i\hbar\delta^{ij}$. In terms of “ra” variables where $r={1\over 2}(1+2)$ and $a=1-2$, the only non-vanishing commutators are $[\bm{x}^{i}_{r/a},\bm{p}^{j}_{a/r}]=i\hbar\delta^{ij}$, especially $\bm{x}_{r}$ and $\bm{p}_{r}$ commute with each other. This allows us to introduce a wave function $\hat{\rho}(\bm{x}_{r},\bm{p}_{r})$ which is our density matrix in phase space. Since $x_{r}$ and $p_{a}$ are conjugate variables, the wave function (or density matrix) in momentum space $\hat{\rho}(\bm{p}_{r},\bm{p}_{a})=\hat{\rho}(\bm{p}_{1},\bm{p}_{2})$ is related to the density matrix in phase space $\hat{\rho}(\bm{x}_{r},\bm{p}_{r})$ by a Fourier (or Wigner) transform

[TABLE]

We will assume that the density matrix in phase space $\hat{\rho}(\bm{x},\bm{p})$ is a slowly varying function on space $\bm{x}$, compared to a microscopic scale of QCD interactions, usually set by the mean free path $l_{mfp}\sim(\alpha_{s}^{2}T)^{-1}$. This means that we can consider $\bm{x}_{r}$ as constant in the computation of relaxation operator $\Gamma$ in (1.1). This is translated to a smallness of $\bm{p}_{a}\sim\partial_{\bm{x}}\ll l_{mfp}^{-1}$ by (1.3), that is, the density matrix in momentum space $(\bm{p}_{1},\bm{p}_{2})$ is nearly diagonal in momentum variables. If we neglected spin degrees of freedom, these diagonal elements would correspond to a usual distribution function $f(\bm{p})$ in momentum space. In our computation of $\Gamma$, we therefore work with diagonal elements in the density matrix in momentum space, defined by $\hat{\rho}(\bm{p}_{1},\bm{p}_{2})\equiv(2\pi)^{3}\delta(\bm{p}_{1}-\bm{p}_{2})\hat{\rho}(\bm{p}_{1})$. This is justified as long as we don’t care about the advective terms in $\bm{x}$ in quantum kinetic equation, but focus only on local “collision terms” of $\Gamma$ in (1.1) in the spatial homogeneous limit. Note that we still keep a full spin matrix of $\hat{\rho}(\bm{p})$ in the spin space.

The reason why we need to keep full quantum correlation of spin degrees of freedom in the density matrix $\hat{\rho}(\bm{p})$ is that the two spin states are degenerate in energy and the quantum correlation time is arbitrarily large $\tau_{q}\sim{\hbar\over\Delta E}\to\infty$. Even in the presence of background vorticity or magnetic field, the energy shift is $\Delta E\sim\bm{S}\cdot\bm{\omega}$ or $\Delta E\sim Q\bm{S}\cdot\bm{B}$ ($Q$ is charge) which is ${\cal O}(\hbar)$ since $\bm{S}\sim{\cal O}(\hbar)$. The quantum correlation time for spin is $\tau_{q}\sim{\hbar\over\Delta E}\sim{\cal O}(\hbar^{0})$, and it is in the classical time scale, that is usually described by a classical Boltzmann kinetic theory. Therefore, quantum correlation of spin should be considered even in a regime of kinetic theory, and hence “quantum kinetic theory”.

A more fundamental treatment of (1.1), that should require a significantly larger effort in the future, would involve a complete analysis of spin density matrix in the full phase space $(\bm{x},\bm{p})$. The free streaming, collision-less quantum kinetic equation for this case was recently studied in Ref [21, 22, 23, 24]. Our study in this sense can be viewed as providing the collision term in leading log of pQCD. However, our result should be improved in this case, including spatial gradient effects in the collision terms (we are restricting to the homogeneous limit, as described in the previous paragraph). This is because the vorticity is a spatial gradient of background fluids. Another way to understand this is that the orbital angular momentum of background fluids can in general be transferred to spin angular momentum of the massive quarks we are looking at. In this more fundamental picture, since total angular momentum has to be conserved, a loss of spin angular momentum of the massive quark (described by $\Gamma$) must be compensated by a gain of angular momentum in the background fluid. This gain term will be shared among all quasi-particles of the background, dominantly light quarks and gluons. Compared to the massive quark we consider, these other degrees of freedom is much larger, and the spin gain is diluted and its back reaction to the equation for $\hat{\rho}$ will be suppressed compared to the loss term. More importantly, the gain in angular momentum of background will be shared between orbital and spin angular momenta. Since spin is smaller than orbital by $\hbar$, the most of gain will go to the orbital angular momentum with a change of vorticity, $\Delta\bm{\omega}$. Since the spin is $\bm{S}\sim{\cal O}(\hbar)$, the change $\Delta\bm{\omega}\sim{\cal O}(\hbar)$, and its effect to the dynamics of $\bm{S}$ via $\Delta E=-\bm{S}\cdot\bm{\omega}\sim{\cal O}(\hbar^{2})$ is higher order in small $\hbar$. Based on this consideration, we neglect possible “gain terms” in our quantum kinetic equation (1.1). In essence, we treat the background as a spin reservoir that can absorb any change of spin angular momentum in $\hat{\rho}$ of “dilute” massive quarks, without any back reaction of the absorbed spin to the evolution of $\hat{\rho}$ itself. This is justified as long as we track the spin polarization of dilute massive quarks only, without caring about those of light quarks and gluons.

The characteristic relaxation rate of spin polarization of massive quarks will be shown to be of order $\Gamma_{0}\sim\alpha_{s}^{2}\log(1/\alpha_{s})T$ where $\alpha_{s}=g^{2}/(4\pi)$†††It can be shown that this is in fact true for light quarks and gluons as well, as it is universal for all soft $t$-channel processes [27].. This is of the same microscopic relaxation rate that governs other transport coefficients, such as shear viscosity or charge conductivities. See Refs.[25, 26] for similar observations, but in terms of usual scattering rate picture, i.e. considering only diagonal elements of the density matrix. These contributions to $\Gamma_{0}$ arise from soft $t$-channel gluon exchange of momentum $q$ in the scatterings with background hard thermal particles, where the log comes from a range $m_{D}\sim gT\ll q\ll T$. If the quark was light, there would also exist soft $t$-channel quark exchange contribution of the same leading log order, making conversion of a quark to a gluon [27, 28]. See Figure 1.

The spin polarization of light quarks can be transferred to that of gluons and vice versa by these conversion processes, and a complete picture for light quarks would have to involve the spin density matrices of both light quarks and gluons. For the process of Figure 1 to happen, the exchanged quark should be the same species of the incoming quark. Our assumption of a hard-scale mass $m\gg gT$ for the massive quark implies $q\gtrsim m\gg gT$, which makes a soft $q$ quark exchange impossible for the massive quark. This justifies the absence of quark-gluon conversion process of Figure 1 for massive quark at leading log in our study.

2 Time evolution of spin density matrix in Schwinger-Keldysh formalism

We consider the Hilbert space of one-quark state interacting with background QCD plasma degrees of freedom in a finite temperature $T$. This reduced description is justified as long as the occupation number of quark per unit quantum state (given by $({\rm Number\,of\,quarks})\times(2\pi\hbar)^{3}/(d^{3}\bm{x}d^{3}\bm{p})$ in phase space, usually called the distribution function $f(\bm{x},\bm{p})$) is much less than unity (or “dilute” Boltzmann limit), so that quantum statistics of Pauli blocking is negligible. We assume that our massive quark species satisfies this condition, either by $m\gtrsim T$ due to thermal Boltzmann suppression, or at least in early stages in heavy-ion collisions when the massive quarks are scarce. A convenient basis of states for our purpose is $\{|\bm{p},\pm\rangle\}$ of a momentum $\bm{p}$ and helicity $h=\pm 1/2$ (meaning that the spin state is an eigenstate of the spin angular momentum along $\hat{\bm{p}}=\bm{p}/p$ ($p\equiv|\bm{p}|$) with the eigenvalue $\pm\hbar/2$, that is, $(\hat{\bm{p}}\cdot\bm{\sigma})|\bm{p},\pm\rangle=\pm|\bm{p},\pm\rangle$). As explained in the introduction, we consider a density matrix that is (approximately) diagonal in the momentum variable, that is sufficient for describing the local quantum collision term $\Gamma_{0}$ in the quantum kinetic equation. We therefore have the density matrix per unit volume as (we set $\hbar=1$ from here without much confusion)

[TABLE]

where $\hat{\rho}(\bm{p})$ is a $2\times 2$ spin density matrix at a fixed momentum $\bm{p}$. More explicitly, we have

[TABLE]

in bra-ket notation, with a set of four functions in momentum space, $\rho_{s,s^{\prime}}(\bm{p})$.

It is important to recall a phase ambiguity of the basis states $|\bm{p},s\rangle\to e^{i\phi(\bm{p},s)}|\bm{p},s\rangle$ with an arbitrary choice of $\phi(\bm{p},s)$, which is reflected to the compensating phase ambiguity of $\rho_{s,s^{\prime}}(\bm{p})\to e^{-i(\phi(\bm{p},s)-\phi(\bm{p},s^{\prime}))}\rho_{s,s^{\prime}}(\bm{p})$, such that the density matrix $\hat{\rho}$ is unambiguous. We will be careful about this ambiguity in our computation, such that our final quantum kinetic equation in terms of physical spin polarization is well defined free of this phase ambiguity.

One way of fixing the phase ambiguity is to work universally in the basis of z-component of spin operator. Then, the helicity $s/2$ state has an explicit 2-component spinor representation, $\xi_{s}(\bm{p})$, satisfying $(\hat{\bm{p}}\cdot\bm{\sigma})\xi_{s}(\bm{p})=s\xi_{s}(\bm{p})$, with normalization $\xi_{s}^{\dagger}\xi_{s}=1$. The density matrix in this basis is then an explicit $2\times 2$ matrix, given by

[TABLE]

We emphasize again that $\xi_{s}(\bm{p})$ and hence $\rho_{s,s^{\prime}}(\bm{p})$ are each phase ambiguous, but the net density matrix $\hat{\rho}(\bm{p})$ is free of ambiguity. The spin operator in this basis is $\bm{S}={1\over 2}\bm{\sigma}$ (with $\hbar=1$), and the spin polarization density in momentum space from $\hat{\rho}(\bm{p})$ is then given by

[TABLE]

Recalling that ${\rm Tr}(\hat{\rho}(\bm{p}))$ is the usual number distribution $f(\bm{p})$ that appears in the conventional semi-classical Boltzmann equation, we can express the density matrix in this basis as

[TABLE]

Note that $f(\bm{p})$ and $\bm{S}(\bm{p})$ are physical quantities and are independent of our basis choice: we will present our quantum kinetic equation for $f(\bm{p},t)$ and $\bm{S}(\bm{p},t)$ ($t$ is time). The total spin polarization and the number of quarks per unit volume is then given by

[TABLE]

The density matrix $\hat{\rho}(t)$ and its time evolution is most naturally described in the Schwinger-Keldysh formalism. The path integral in the time-forward contour (labeled as contour 1) gives the unitary time evolution of the ket part of the density matrix, and that in the time-backward contour (contour 2) gives the complex conjugate evolution of the bra part of the density matrix,

[TABLE]

where $U_{1,2}(t,t_{0})={\cal P}e^{-i\int_{t_{0}}^{t}dt^{\prime}\,H_{1,2}(t^{\prime})}$ are the unitary time evolutions in the contours 1 and 2 respectively, and $\langle\cdots\rangle_{A}$ means the thermal path integral average of background degrees of freedom of the QCD plasma. Since our system is an open system and is interacting with background degrees of freedom (more precisely, through the soft-scale color gauge field $A_{\mu}(x)$ in leading log order), the $H_{1}(t)$ depends on the operator of the background degrees of freedom in contour 1 that couples our system to the background, and it is time-dependent in general due to time-dependence of that operator (i.e. the color gauge field in contour 1, $A^{(1)}_{\mu}(\bm{x},t)$). The same is true for $H_{2}(t)$ and $A^{(2)}_{\mu}(\bm{x},t)$. The average $\langle\cdots\rangle_{A}$ in the above then involves the thermal correlation functions of $A^{(1)}$ and $A^{(2)}$ in the Schwinger-Keldysh contours (the two-point functions in our leading order computation). In the frequency space, these correlation functions satisfy the KMS relations. More explicitly, defining

[TABLE]

what we will need later are the relations,

[TABLE]

in terms of the gluon spectral density $\rho_{\mu\nu}\equiv i(G^{R}_{\mu\nu}-(G^{R}_{\nu\mu})^{*})=-2\,{\rm Im}[G^{R}_{\mu\nu}]$, where $G^{R}$ is the retarded two-point function and the last equality holds only for symmetric case that is true in parity (P)-even background that we assume‡‡‡See Refs.[29, 30] for an introduction to a possible anti-symmetric part, that is called “P-odd spectral density”, in the presence of background axial charge., and $n_{B}(q^{0})=1/(e^{\beta q^{0}}-1)$ is the Bose-Einstein distribution. In our leading log computation, these correlation functions include the well-known 1-loop Hard-Thermal-Loop (HTL) self-energy, the imaginary part of which gives the non-vanishing spectral density in soft t-channel space-like momenta, that represents the scatterings with background thermal particles by cutting the 1-loop, while the real part regulates the infrared divergence in these t-channel scatterings by (real-time) screening effects due to background thermal particles (see the section 3 for a more detailed review on this).

The Hamiltonian in our one-quark picture is a sum of the free kinetic energy, $H_{0}$, and the QCD interaction with background gluon fields, $H_{I}$. The interaction Hamiltonian arises from the field theory Hamiltonian

[TABLE]

where $\psi(\bm{x})$ is the quark field operator, and $A^{a}_{\mu}(\bm{x})$ is the gluon field with color index $a$ ($t^{a}$ are the color generators). We choose our convention as

[TABLE]

and $\bar{\psi}\equiv-\psi^{\dagger}\gamma^{0}$. In this convention, the quark spinor of momentum $\bm{p}$ and helicity $h=\pm 1/2$, that shares the same phase ambiguity as the state $|\bm{p},\pm\rangle$, is explicitly given by

[TABLE]

where $p=|\bm{p}|$, $E_{p}=\sqrt{p^{2}+m^{2}}$ and $(\hat{\bm{p}}\cdot\bm{\sigma})\xi_{s}(\bm{p})=s\xi_{s}(\bm{p})$ ($s=\pm 1$). This explicit expression will be used in our computation of spin-dependent transition amplitudes. A quick way to see why the (arbitrary) phase of $|\bm{p},s\rangle$ is identical to that of $u(\bm{p},s)$ is to note the field operator $\psi(\bm{x})$ expanded as

[TABLE]

where $a_{\bm{p},s}$ is the annihilation operator of one-quark state. Noting that $|\bm{p},s\rangle\sim a^{\dagger}_{\bm{p},s}|0\rangle$, the phase ambiguity of $|\bm{p},s\rangle$ (or equivalently, $a^{\dagger}_{\bm{p},s}$) is precisely identical to the phase ambiguity of $u(\bm{p},s)$, such that the $\psi(\bm{x})$ operator and its conjugate entering the interaction Hamiltonian $H_{I}$ are unambiguous. Ultimately, this phase ambiguity becomes that of the 2-component spinor $\xi_{s}(\bm{p})$ in (2.15). As we use $H_{I}$ in our computation with a consistent use of $\xi_{s}(\bm{p})$ in both $H_{I}$ and the definition of density matrix (2.6), our result for $f(\bm{p})$ and $\bm{S}(\bm{p})$ is free of this ambiguity.

Since we need to keep the normalization of one-quark state correctly when discussing the density matrix, it is most convenient to work in a finite volume $V$ with discrete spectrum of states and then take an infinite volume limit. The momentum space becomes discrete $\bm{p}_{\bm{n}}$ with integer-valued label vector $\bm{n}$, and the infinite volume limit is

[TABLE]

The fields are expanded as (including only the quark sector, neglecting anti-quarks)

[TABLE]

and

[TABLE]

where $a_{\bm{p}_{\bm{n}},s}$ is the annihilation operator of one-quark state $|\bm{p}_{\bm{n}},s\rangle$, with unit normalization, that is

[TABLE]

and $A_{\mu}(\bm{p}_{\bm{n}},t)$ is defined such that it has the two-point correlation functions as

[TABLE]

with the usual infinite volume correlation function $G^{(ij)}_{\mu\nu}(\bm{p},t)$, so that $\langle A^{(i)}(\bm{x},t)A^{(j)}(\bm{x}^{\prime},t^{\prime})\rangle_{A}$ has the correct infinite volume limit (that is, independent of the volume $V$ as a local correlation function). Then, the one-quark Hamiltonian from $H_{I}$ becomes

[TABLE]

where the color structure is omitted for notational simplicity. The normalized one-quark states are created by

[TABLE]

To obtain the time evolution equation of the density matrix from (2.10) in perturbation theory of $H_{I}$, we work in the interaction picture of $H_{0}$ and need to expand $U_{1,2}$ to quadratic order in $H_{I}$, since one-point functions of gluon fields vanish, $\langle A\rangle_{A}=0$, and the first non-vanishing correlation functions are the two-point functions. Recall the interaction picture: $U(\Delta t,0)=U_{0}(\Delta t)U_{I}(\Delta t,0)$ where $\Delta t$ is the time step we are considering, $U_{0}(\Delta t)=e^{-iH_{0}\Delta t}$ is the free evolution, and we set the initial time as $t_{0}=0$ without loss of generality. Then, we have the interaction picture evolution as

[TABLE]

with $H^{\rm int}_{I}(t)=U_{0}(t)^{\dagger}H_{I}(t)U_{0}(t)$. Using these in (2.10), we have

[TABLE]

where the second term in the first line comes from “cross” combination obtained by one $H_{I}$ from $U_{1}$ and one from $U_{2}$, while the last two lines are “self-energy” contributions coming from quadratic expansions in $H_{I}$ in each $U_{1}$ and $U_{2}$. We need the both types of contributions in order to make sure the probability conservation of the density matrix, that is, ${\rm Tr}(\hat{\rho}(\Delta t))={\rm Tr}(\hat{\rho}(0))$ to this order, that can be checked easily using the definition of correlation functions in Schwinger-Keldysh contours. See Figure 2 for diagrammatic representation of these contributions. The $H_{I}^{{\rm int}(1)}$ and $H_{I}^{{\rm int}(2)}$ are the Hamiltonians (2.22) with $A^{(1)}$ and $A^{(2)}$ fields, respectively. As they are linear in $A$ fields, their correlation functions in (2.25) are proportional to the two-point correlation functions of background gluon fields, defined in (2.11).

Since our density matrix in (2.4) is diagonal in momentum space, and hence in energy spectrum of $H_{0}$, it commutes with $H_{0}$ and the first term in (2.25) is simply $\hat{\rho}(0)$. Writing $\hat{\rho}(\bm{p}_{\bm{n}})$ in the helicity basis (2.5) in terms of $\rho_{s,s^{\prime}}(\bm{p}_{\bm{n}},t)$ (now with time-dependence), and using the explicit form of $H_{I}$ in (2.22), the above evolution equation (2.25) can easily be translated to those of $\rho_{s,s^{\prime}}(\bm{p}_{\bm{n}},t)$. From the identification $|\bm{p}_{\bm{n}},s\rangle=a^{\dagger}_{\bm{p}_{\bm{n}},s}|0\rangle$, it can be seen that the phase ambiguity of $\rho_{s,s^{\prime}}(\bm{p}_{\bm{n}})$ via that of $|\bm{p}_{\bm{n}},s\rangle$ will cancel in the expression of the physical density matrix (2.6),

[TABLE]

due to the fact that $u(\bm{p}_{\bm{n}},s)$ appearing in $H_{I}$ is proportional to $\xi_{s}(\bm{p}_{\bm{n}})$ and it shares the same phase with the $a^{\dagger}_{\bm{p}_{\bm{n}},s}$ (so that $H_{I}$ is unambiguous) and hence $|\bm{p}_{\bm{n}},s\rangle$.

When the time step $\Delta t$ of evolution is much larger than the correlation time of the gluon two-point functions, the cross and self energy terms in (2.25) will be linear in $\Delta t$, and the resulting evolution equation of density matrix will be of first order in time. In our leading log order, the dominant contribution to these terms come from the HTL contribution to the gluon two-point functions, with soft frequency-momenta in a range $gT\ll q\ll T$. This gives an estimate for the correlation time $\tau_{c}\lesssim(gT)^{-1}$ for leading log. As long as $\tau_{c}\ll\Delta t\ll 1/\delta E=\infty$, where $\delta E$ is the energy difference between the two quantum spin states, this is a valid description of the time evolution.

More explicitly, for the cross term as an example, the two time integrals with the gluon two-point function can be written schematically as

[TABLE]

with a change of variable $t_{r}=(t_{1}+t_{2})/2$ and $t_{a}=t_{1}-t_{2}$, and $q^{0}$ is some combination of energies of states (see the following expressions in this section). As $G^{(12)}(t_{a})$ decays fast beyond $t_{a}>\tau_{c}$, we can extend the range of $t_{a}$ integral to $[-\infty,+\infty]$ for most of $t_{r}$ values in $[0,\Delta t]$ when $\Delta t\gg\tau_{c}$: this gives the leading term linear in $\Delta t$. Then the above becomes $G^{(12)}(q^{0})\Delta t$ where $G^{(12)}(q^{0})$ is the Fourier transform of $G^{(12)}(t)$. A similar manipulation can be done for self energy terms to get the leading linear term in $\Delta t$.

After some algebra with these ingredients, and taking an infinite volume limit, we obtain a well-defined evolution equation for $\rho_{s,s^{\prime}}(\bm{p},t)$. Writing

[TABLE]

the cross contribution is given by

[TABLE]

and the self-energy contribution is a sum

[TABLE]

In the Appendix, we prove the following important simplification in the self energy term, due to the rotational symmetry and P-even nature of thermal QCD background: the $\bm{p}^{\prime}$ integral appearing in the self energy term

[TABLE]

is proportional to $\delta_{s,s^{\prime\prime\prime}}$, that is, it is non-zero only when the initial and final spins are the same, $s=s^{\prime\prime\prime}$, and moreover, the value doesn’t depend on $s=\pm$. Physically, what it means is that the self energy can’t flip the longitudinally polarized spin (i.e. helicity, $s/2$) due to rotational invariance of the background, and the self energy can’t depend on the sign of helicity $s/2$ either, since $s$ flips under parity transformation. With this, the two terms in the self energy (2.30) nicely combine to give

[TABLE]

in terms of the Fourier transform of the correlation function $G^{(21)}_{\mu\nu}$, where we use an identity for Schwinger-Keldysh two-point functions,

[TABLE]

to combine the two $t_{a}$ integrals.

The appearance of $G^{(12)}$ in the cross term and $G^{(21)}$ in the self energy is a reflection of the generic feature that ensures the thermal detailed balance with the KMS relation (2.12). In fact, using the identity

[TABLE]

one can easily check from (2.29) and (2.32) that the equilibrium thermal density matrix

[TABLE]

with any fugacity constant $z=e^{\mu/T}$ makes the sum of $\Gamma_{\rm cross}$ and $\Gamma_{\rm self\,energy}$ terms vanishes. This equilibrium is also equivalently described by

[TABLE]

The rest of the paper presents key elements of our computation of the integrals in (2.29) and (2.32) in leading log order. Readers who are interested in only the final results can go straight to section 4.

3 Leading log integrals with arbitrary quark mass

We first consider the evaluation of the cross term (2.29). Since the physical density matrix that is free of phase ambiguity that we discussed is $\hat{\rho}(\bm{p})=\sum_{s,s^{\prime}}\xi_{s}(\bm{p})\rho(\bm{p})_{s,s^{\prime}}\xi_{s^{\prime}}^{\dagger}(\bm{p})$, we consider this object. Using also that

[TABLE]

we can express (2.29) in terms of the unambiguous $\hat{\rho}(\bm{p})$. Since we will focus only on the soft $\bm{q}=\bm{p}-\bm{p}^{\prime}\sim gT$ regime that produces the leading log, we change the integration variables from $\bm{p}^{\prime}$ to $\bm{q}$,

[TABLE]

What needs to be computed in the resulting integrand is the following spinor summation, contracted with the gluon two-point function,

[TABLE]

where $\bm{p}^{\prime}=\bm{p}-\bm{q}$.

We work in the Coulomb gauge, where the gluon two-point functions are written in terms of the longitudinal and transverse spectral densities, $\rho_{L}$ and $\rho_{T}$ respectively, (not to be confused with the density matrix)

[TABLE]

with

[TABLE]

where $q\equiv|\bm{q}|$, and the transverse projection operator has only spatial components as

[TABLE]

From the explicit expression of spinor (2.15), we have

[TABLE]

where $\sigma^{\mu}=({\bf 1},\bm{\sigma})$ and $\bar{\sigma}^{\mu}=({\bf 1},-\bm{\sigma})$. Noting that

[TABLE]

is the projection operator to the helicity $s/2$ state in spin space, the spinor sum in (3.39) becomes a summation of various terms of the type,

[TABLE]

The computation of these spinor sum is challenging, but is doable with some efforts utilizing several properties of the projection operators, such as

[TABLE]

for any operator $A$. Note that these expressions are free of phase ambiguity of $\xi_{s}(\bm{p})$, due to the reasons as explained before.

Since $G^{(12)}_{\mu\nu}(q^{0},\bm{q})=n_{B}(q^{0})(\delta_{\mu 0}\delta_{\nu 0}\rho_{L}(q^{0},q)+\Pi^{T}_{\mu\nu}(\bm{q})\rho_{T}(q^{0},q))$, we present our results of computation of (3.39) in terms of longitudinal (that involves $\rho_{L}$) and transverse ($\rho_{T}$) gluon parts, respectively. The longitudinal part of (3.39) is given as, omitting the obvious factor of $n_{B}(q^{0})\rho_{L}(q^{0},q)$ ($q^{0}\equiv E_{p}-E_{p^{\prime}}$), and writing $\hat{\rho}(\bm{p})$ in terms of the physical quantities (see (2.8)), $\hat{\rho}(\bm{p})={1\over 2}f(\bm{p})+\bm{S}(\bm{p})\cdot\bm{\sigma}$,

[TABLE]

The transverse part of (3.39), that is proportional to $n_{B}(q^{0})\rho_{T}(q^{0},q)$, is obtained after a long computation as

[TABLE]

We have checked the validity of the above results at least in two special limits: 1) $\bm{p}^{\prime}=\bm{p}$ limit (treating $\hat{\bm{q}}$ arbitrary), and 2) massless ($m=0$) limit. In both limits, the spinor sum (3.39) reduces to $s^{\prime\prime}=s$ and $s^{\prime\prime\prime}=s^{\prime}$ cases only, due to the fact that $\bar{u}(\bm{p},s)\gamma^{\mu}u(\bm{p}^{\prime},s^{\prime\prime})$ vanishes in these limits unless $s=s^{\prime\prime}$, which can be easily checked from the explicit spinor expression (2.15). Using this fact, one can compute (3.39) in these limits directly, and then can compare with the above results in the same limits. The longitudinal part is easy to compare, but the comparison of the transverse part needs some non-trivial identities. In the limit 1), one needs the following “dyad” identity ($\otimes$ is a dyad product of row and column vectors)

[TABLE]

for any two unit vectors $\hat{\bm{p}}$ and $\hat{\bm{q}}$, where $\bm{I}$ is the $3\times 3$ identity dyad (matrix). In the case of the limit 2), one needs a more non-trivial identity that we checked by Mathematica,

[TABLE]

for any three unit vectors $\hat{\bm{p}}$, $\hat{\bm{p}}^{\prime}$ and $\hat{\bm{q}}$. These agreements give us confidence on the validity of the above spinor sum results.

The computation of spin sum in the self energy (2.32) is simpler. First note that the self energy term has a simple structure

[TABLE]

that is, it is a constant ($\gamma$) times of the identity operator in both the spin and momentum space. The “damping rate” $\gamma$ is (recall $\bm{p}^{\prime}=\bm{p}-\bm{q}$)

[TABLE]

Recalling that this expression doesn’t depend on $s$ (due to parity invariance as proved in Appendix), it turns out to be easier to compute the spin sum by expressing it as

[TABLE]

removing any reference to $s$. For the longitudinal gluon contribution ($\mu=\nu=0$), the spin sum becomes (recall $G^{(21)}_{\mu\nu}(q^{0},\bm{q})=(n_{B}(q^{0})+1)(\delta_{\mu 0}\delta_{\nu 0}\rho_{L}(q^{0},q)+\Pi^{T}_{\mu\nu}(\bm{q})\rho_{T}(q^{0},q))$)

[TABLE]

and for the transverse gluon contribution, we obtain

[TABLE]

After (3.39) and the spin sum in (3.53) is computed, what remains is to compute the $\bm{q}$ integrations in (2.29) and (3.53) to leading log order, with the gluon spectral densities $\rho_{L/T}(q^{0},q)$ given by the well-known HTL contributions. For completeness, they are given by

[TABLE]

where $m_{D}\sim gT$ is the Debye mass, and

[TABLE]

where

[TABLE]

As in the computations of shear viscosity and conductivities in massless limit [27], as well as in the computation of diffusion constant of heavy quark limit [31], we find that the leading log contribution comes from the same soft $q$ regime, for arbitrary quark mass, where the log arises from the range $gT\ll q\ll T$. Physically, this contribution represents the t-channel scatterings with background thermal particles of hard scale ($p\sim T$) with a soft gluon exchange of momentum $(q^{0},\bm{q})$. We emphasize that these HTL contributions include only the thermal background gluons and light quarks, which means that we don’t include the scatterings with the other massive quarks present in the background plasma. This is the same “diluteness” assumption we explained at the beginning of section 2. The Debye mass in this case is given by

[TABLE]

where $N_{F}=2$ is the number of light flavors.

We follow the known steps of computing $\bm{q}$ integration in leading log order [27, 32, 33, 28, 29]. A first step in this $\bm{q}$ integration is to make a change of variable from the azimuthal angle $\cos\theta_{pq}$ between $\bm{p}$ and $\bm{q}$ to the energy transfer $q^{0}=E_{p}-E_{p^{\prime}}$ (recall $\bm{p}^{\prime}\equiv\bm{p}-\bm{q}$), where they are related by

[TABLE]

The $q^{0}$ has a maximum $q^{0}_{\rm max}$ (minimum $q^{0}_{\rm min}$) when $\theta_{pq}=0$ ($\pi$), and

[TABLE]

Note that $q^{0}_{\rm min}$ is different from $-q^{0}_{\rm max}$ by a term proportional to $q^{2}$, that is present only in the massive case. We will see that this $q^{2}$ correction to the $q^{0}$ integration range, that is absent in massless case, gives rise to the same leading log contribution to the final result, so it is important to keep it to this order. From this, we can convert $\bm{q}$ integration in (2.29) and (3.53) into an integration of two variables $(q^{0},q)$,

[TABLE]

where (spinor sum) is the spinor part that we computed above, and the integration of polar angle around $\hat{\bm{p}}$ axis gives $(2\pi)$, after we make the spinor part to be independent of polar angle, exploiting rotational symmetry, since $\rho_{L/T}(q^{0},q)$ depends on $\bm{q}$ only via $q=|\bm{q}|$. Specifically, we can replace

[TABLE]

One then computes $(q^{0},q)$ integration with the gluon spectral densities $\rho_{L/T}(q^{0},q)$ to the desired order that produces the leading log in the final result. Since the leading log comes from soft $(q^{0},q)\sim gT$ regime, one expands the (spinor sum) part in power series of soft $(q^{0},q)\ll(p,E_{p},m)\sim T$: it is sufficient to keep only up to linear order in $(q^{0},q)$, as higher powers give higher order terms in $g$. The total sum of the cross and self energy terms, (2.29) and (3.51) respectively, takes a form of (3.62) after this expansion, where the (spinor sum) part has a structure of

[TABLE]

with two functions $C^{0,1}(q^{0}/q)$ on $q^{0}/q$. There are two important features in this result, that makes the leading log contribution possible:

1. In principle, since the (spinor sum) contains $n_{B}(q^{0})\sim T/q^{0}$ for $q^{0}\sim gT\ll T$ (see (2.12)), the expansion could start from $(1/q^{0})C^{-1}(q^{0}/q)$, instead of $C^{0}(q^{0}/q)$. In fact, both the cross and the self energy terms start from this order, but their sum cancels to this order. If this cancellation was absent, the final result of spin evolution rate would have been dominated by this order, which gives $g^{2}\log(1/g)$, instead of $g^{4}\log(1/g)$ that we find. This contribution would come from the ultra-soft range $g^{2}T\ll q\ll gT$, which represents “small angle” scatterings, contrary to our range $gT\ll q\ll T$ for $g^{4}\log(1/g)$ that represents “large angle” scatterings. This cancellation is important also in computations of shear viscosity and charge conductivities (but not in “color” conductivity [34, 35]), and has been shown to be related to conservation Ward identities of energy-momentum and charge currents [33]. The same cancellation we observe in our spin density matrix suggests it may be related to angular momentum conservation. This cancellation also adds confidence that our computation of spin sum is correct. 2) In principle, the functions $C^{0,1}(q^{0}/q)$ could be any function on $q^{0}/q$, but they turn out to be even functions on $q^{0}/q$, which is crucial to have the final $g^{4}\log(1/g)$ rate. This feature is important due to the fact that the spectral densities, $\rho_{L/T}$, are odd functions on $q^{0}$. If $q^{0}_{\rm min}$ was precisely equal to $-q^{0}_{\rm max}$ (as in the massless case), the $q^{0}$ integral of $C^{0}(q^{0}/q)$ would have vanished, and the first non-vanishing result would come from the next order $q^{0}C^{1}(q^{0}/q)$ term, which gives $g^{4}\log(1/g)$ rate. Due to the $q^{2}$ correction to the $q^{0}$ integration range, $q^{0}_{\rm max/min}$, in our massive case (see (3.61)), the $C^{0}(q^{0}/q)$ integral does contribute, but since this correction is one higher order than the leading range, the result is of the same order as the one from $q^{0}C^{1}(q^{0}/q)$, that is, the same $g^{4}\log(1/g)$ rate.

Let us define $q^{0}$ integrals of spectral densities, that are needed in our computation described above,

[TABLE]

In the massless case, only integer $n$ survives due to $q_{\rm min}=-q_{\rm max}$ in that case, and they can be computed by a sum-rule technique [32, 33, 29], utilizing analytic property of the spectral densities. In our massive case, and only for the range of our interests, $gT\ll q\ll T$, the spectral densities can be simplified to produce the results for $J_{n}^{L/T}$ to our desired order

[TABLE]

where $m_{D}\sim gT$ is the Debye mass. One can easily check that in the massless limit, these spectral densities produce the same results for $J_{n}^{L/T}$ from the sum-rule technique. These expressions are obtained from the full expressions, (3.56) and (3.57), by using the hierarchy $m_{D}\sim gT\ll(q^{0},q)\ll T$. It is convenient to write $J^{L/T}_{n}$ as

[TABLE]

in terms of the dimensionless coefficient functions $j^{L/T}_{n}$ on $(p,E_{p},m)$, after extracting the dependence on $q$ explicitly as above. By explicit evaluations, we find them as in Table 1.

Note that the half-integer $n$ cases are needed for the contributions from $C^{0}(q^{0}/q)$ as we explained above, which exist only in the massive case. After doing $q^{0}$ integration using these formula, one finally performs $q$ integration in (3.62) to get the leading log result, where $q$ ranges in $m_{D}\ll q\ll T$: these boundaries come from the fact that our expression for the integrand is valid only in this range. The log arises from

[TABLE]

4 Quantum kinetic equation for spin polarization of massive quarks

After a lengthy, but straightforward computation that we describe in the previous section, we present our final result for the time-evolution of the spin density matrix in momentum space, $\hat{\rho}(\bm{p})={1\over 2}f(\bm{p})+\bm{S}(\bm{p})\cdot\bm{\sigma}$, in leading log order of $g^{4}\log(1/g)$. We write these evolution equations as

[TABLE]

where $\Gamma_{f}$ and $\bm{\Gamma}_{S}$ are diffusion-like differential operators in momentum space, that contain up to second order derivatives in $\bm{p}$. The $\Gamma_{f}$ is given by ($\bm{\nabla}_{p}\equiv\partial/\partial\bm{p}$ and $(\bm{p}\cdot\bm{\nabla}_{p})^{2}f\equiv\bm{p}\cdot\bm{\nabla}_{p}\left(\bm{p}\cdot\bm{\nabla}_{p}f\right)$)

[TABLE]

This result passes a very non-trivial test of the expected detailed balance: one can check that $\Gamma_{f}=0$ when $f(\bm{p})=f^{\rm eq}(\bm{p})=ze^{-E_{p}/T}$ for any constant $z$. It should be emphasized that this check is satisfied irrespective of the values of $j^{L/T}_{n}$, because the detailed balance is independent of details of the spectral densities that determine $j^{L/T}_{n}$. This gives us confidence that our computation is correct. We also note that our result for $\Gamma_{f}$ that provides the local collision term in leading log order, together with free streaming advection term in Boltzmann equation, can be used to compute several conventional transport coefficients, such as shear viscosity and electric conductivity, arising from dilute massive quarks.

For the spin polarization part, we obtain ($i=1,2,3$ denotes a spatial index for vector)

[TABLE]

Note that $f(\bm{p})$ and $\bm{S}(\bm{p})$ do not mix with each other in these equations.

There is a highly non-trivial test of the above result in the massless limit. Note that our computation doesn’t include the quark-gluon conversion processes that becomes of the same order in the massless limit, and also our values of $J^{L/T}_{n}$ do not have the correct massless limit, so the massless limit of the above result should not be taken as the true result for the massless case. What we are testing is the “consistency” of the above equations with the “chirality conservation” in the massless limit, and this test is a kinematical one, and should hold true for each scattering processes included, independent of details of spectral densities, that is, the values of $J^{L/T}_{n}$. In the massless limit, the negative helicity state ($s=-1$, left-handed) and positive helicity state ($s=+1$, right-handed) are decoupled, and do not mix by gauge interactions. The spin density matrix should then take the following decoupled form

[TABLE]

as a sum of positive and negative helicity chiral quark contributions, where ${\cal P}_{\pm}$ are nothing but the spin projection operators to the two decoupled helicity states. The $f_{\pm}(\bm{p})$ are the number distribution functions of chiral quarks of helicity $s/2=\pm 1/2$ in momentum space. In a parity-even background that we are considering, $f_{+}$ and $f_{-}$ should satisfy the same evolution equation. Writing the above density matrix as

[TABLE]

we see the correspondence to our variables $f(\bm{p})$ and $\bm{S}(\bm{p})$ as,

[TABLE]

Since $f_{+}$ and $f_{-}$ satisfy the same evolution equation, the two functions $f(\bm{p})$ and $f_{s}(\bm{p})$ should satisfy the same equation as well. This means that our above result, when we take the massless limit while keeping $J^{L/T}_{n}$ arbitrary, should pass the following non-trivial tests: 1) the evolution equation for $\bm{S}(\bm{p})$ must admit a consistent Ansatz, $\bm{S}(\bm{p})=f_{s}(\bm{p})\hat{\bm{p}}$, and 2) the resulting evolution equation for $f_{s}(\bm{p})$ must be the same as the one for $f(\bm{p})$ in the massless limit. Both tests require non-trivial cancellations between various terms in (4.71), and it is amusing to check that the tests are satisfied by our results (LABEL:finalF) and (4.71): the $f({\bm{p}})$ and $f_{s}(\bm{p})$ satisfy the same evolution equation in the massless limit with

[TABLE]

This also means that the massless limit allows a broader set of equilibria as

[TABLE]

with arbitrary chiral fugacity constants $z_{\pm}=e^{\mu_{\pm}/T}$. With these remarkable checks satisfied, we become confident that the results (LABEL:finalF) and (4.71) are correct.

Using the explicit values of $j^{L/T}_{n}$ given in Table 1, we have the following expression for the $\Gamma_{f}$,

[TABLE]

It can be checked again that the detailed balance condition is satisfied with $f^{\rm eq}(\bm{p})=ze^{-E_{p}/T}$. For $\bm{\Gamma}^{i}_{S}$, we have

[TABLE]

5 Discussion

Our work is a small step toward a more complete picture of quantum kinetic theory of spin dynamics in perturbative QCD plasma. It is important to extend our work, going beyond the spatial homogeneous limit. This would introduce a spin density matrix $\hat{\rho}(\bm{x},\bm{p})$ that depends on both position and momentum, or equivalently $\hat{\rho}(\bm{p}_{1},\bm{p}_{2})$ which is non-diagonal in momentum space, as explained in the introduction. One could expect that the resulting quantum kinetic equation would look like

[TABLE]

where $\Gamma$, the quantum kinetic collision term, is what we compute in this work. As discussed in the introduction, this picture has to be improved by including 1) gradient corrections to the collision term $\Gamma$ (that is, corrections that involve $\partial_{\bm{x}}\hat{\rho}$), in order to allow spin-orbital angular momentum exchange, and 2) the effects of background electromagnetic fields in both free streaming and collision terms. The 2) for the free streaming part was recently studied in Ref.[21, 22, 23, 24] for massive quarks, extending the recent development of chiral kinetic theory for massless chiral quarks [36, 37, 38, 39, 40, 41, 42], for which a Berry’s curvature in momentum space due to spin projection plays a critical role. The 1) is also expected to be intimately related to the “side-jump” phenomenon in chiral kinetic theory [43, 44, 45]. We hope to make further progress on these important goals in a near future.

The spin polarization in local equilibrium of conventional hydrodynamics description is fixed by hydrodynamic variables, such as temperature and vorticity, and it is not an independent hydrodynamic variable. However, since the total angular momentum including spin has to be conserved, one may think of formulating a hydrodynamics description of possible interplay between spin and orbital angular momenta. There has been recent development in this “spin hydrodynamics” in relativistic regime [46, 47, 48]. Since spin is of order $\hbar$, one can think of this as a ${\cal O}(\hbar)$ quantum correction to the conventional classical hydrodynamics description. Some of the transport coefficients in this spin hydrodynamics [48] should in principle be determined by the quantum kinetic theory that we aim to construct.

Acknowledgment

We thank Yukinao Akamatsu, Gokce Basar, Misha Stephanov and Derek Teaney for helpful discussions. This work is supported by the U.S. Department of Energy, Office of Science, Office of Nuclear Physics, with the grant No. DE-SC0018209 and within the framework of the Beam Energy Scan Theory (BEST) Topical Collaboration.

Appendix: A simplification by rotational symmetry

In this appendix, we show that the following integral is non-zero only when $s=s^{\prime}$, and the value doesn’t depend on $s$ either,

[TABLE]

when the gluon two-point function, $G_{\mu\nu}$, has a decomposition to longitudinal and transverse parts in Coulomb gauge,

[TABLE]

where the only non-zero elements of $\Pi^{T}_{\mu\nu}(\bm{q})$ is

[TABLE]

Let us first show this for the longitudinal case. From the explicit expressions for $\gamma^{\mu}$ in (2.14), and the spinor in (2.15) that we reproduce here

[TABLE]

we have (recall $\bar{u}=-u^{\dagger}\gamma^{0}$)

[TABLE]

where we use

[TABLE]

Now, consider the $\bm{p}^{\prime}$ integral in (5.80), and introduce spherical coordinates $(p^{\prime},\theta^{\prime},\phi^{\prime})$, where $\theta^{\prime}$ is the azimuthal angle between $\bm{p}^{\prime}$ and $\bm{p}$, and $\phi^{\prime}$ is the polar angle around the perpendicular plane to $\bm{p}$. It is easy to see, due to rotational symmetry of other parts in the integrand around $\phi^{\prime}$, that only $\phi^{\prime}$ dependence appears in the spinor sum (5.89), as transverse components of $\bm{p}^{\prime}$ with respect to $\bm{p}$. Since these transverse part of $\bm{p}^{\prime}$ will integrate to zero after $\phi^{\prime}$ integration, it is clear that $\bm{p}^{\prime}\cdot\bm{\sigma}$ in (5.89) will be replaced by something proportional to $\bm{p}\cdot\bm{\sigma}$ after $\phi^{\prime}$ integration. Therefore, (5.80) becomes after $\bm{p}^{\prime}$ integration

[TABLE]

with some constants $A,B,C$. As $u(\bm{p},s)\propto\xi_{s}(\bm{p})$ and $(\bm{p}\cdot\bm{\sigma})\xi_{s}(\bm{p})=sp\xi_{s}(\bm{p})$, and $\xi_{s}^{\dagger}(\bm{p})\xi_{s^{\prime}}(\bm{p})=\delta_{s,s^{\prime}}$, we conclude that (5.80), that is (5.91), is non-zero only when $s=s^{\prime}$. To show that the value doesn’t depend on $s=s^{\prime}$, we use the explicit spinor (5.83) to evaluate (5.91) to obtain

[TABLE]

which is indeed independent of the choice of $s$.

The proof in the transverse case is more complicated, but the idea is the same.

[TABLE]

Upon $\bm{p}^{\prime}$ integration, we again have $\bm{p}^{\prime}\cdot\bm{\sigma}$ become proportional to $\bm{p}\cdot\bm{\sigma}$, and recalling that $\bm{q}=\bm{p}-\bm{p}^{\prime}$, we have

[TABLE]

Due to symmetry of $\phi^{\prime}$ and $\bm{q}=\bm{p}-\bm{p}^{\prime}$, the $\hat{\bm{q}}^{i}\hat{\bm{q}}^{k}$ part will become, after $\phi^{\prime}$ integration, a linear combination of $\delta^{ik}$ and $\bm{p}^{i}\bm{p}^{k}$. Obviously, these two structures are the only possible rank-2 structures with $\bm{p}$, since the only available vector after $\bm{p}^{\prime}$ integration is $\bm{p}$. The $\bm{p}^{i}\bm{p}^{k}$ piece doesn’t contribute to the above due to $\epsilon^{jkl}\bm{p}^{l}$, while the $\delta^{ik}$ piece results in

[TABLE]

that is, the same $\bm{p}\cdot\bm{\sigma}$ structure. Therefore, (LABEL:sumapp3) becomes, after $\bm{p}^{\prime}$ integration in (5.80),

[TABLE]

with constants $A^{\prime},B^{\prime},C^{\prime}$, that is the same structure we obtain in the longitudinal polarization case (see (5.91)), and hence the same conclusions follow.