Relativistic spin hydrodynamics with torsion and linear response theory for spin relaxation

Using the second law of local thermodynamics and the first-order Palatini formalism, we formulate relativistic spin hydrodynamics for quantum field theories with Dirac fermions, such as QED and QCD, in a torsionful curved background. We work in a regime where spin density, which is assumed to relax much slower than other non-hydrodynamic modes, is treated as an independent degree of freedom in an extended hydrodynamic description. Spin hydrodynamics in our approach contains only three non-hydrodynamic modes corresponding to a spin vector, whose relaxation time is controlled by a new transport coefficient: the rotational viscosity. We study linear response theory and observe an interesting mode mixing phenomenon between the transverse shear and the spin density modes. We propose several field-theoretical ways to compute the spin relaxation time and the rotational viscosity, via the Green-Kubo formula based on retarded correlation functions.


Introduction
Spin in the non-relativistic limit, where the spin-orbit coupling is suppressed, becomes an approximately conserved quantity (good quantum number). This is the case for the heavy quarks in QCD (e.g., charm quark); the resulting heavy-quark symmetry shows up as the degeneracy in the spectrum of particular hadrons [1][2][3][4][5]. Another familiar realization is a certain class of condensed matter systems (e.g., magnetic materials), which has motivated recent development of the hydrodynamic theory describing transport phenomena of spin, leading to a fruitful research area known as spintronics (see Ref. [6] for a review). On the other hand, spin in relativistic systems is not a conserved quantity even approximately -It is only a part of the total angular momentum, and is not conserved due to the spin-orbit coupling inherent in relativistic dynamics, whose corresponding transport theory has not been formulated until recently.
Spin angular momentum in the hydrodynamic regime of long time and distance scales is in general not a conserved quantity, due to exchanges with orbital angular momentum via spin-orbit coupling. However, the total angular momentum is conserved. As a result, an equilibrium state is characterized by a corresponding thermodynamic parameter [41] -the thermodynamic conjugate to the total angular momentum, which we refer to as the angular momentum "chemical" potential. In the absence of torsion, this chemical potential is equal to the thermal fluid vorticity [33,34], see Appendix A. In a globally rotating (torsionfree) equilibrium, the angular momentum chemical potential, i.e. the thermal vorticity, is constant in space. Since there is only one angular momentum chemical potential for the conserved total angular momentum, both spin and orbital angular momenta in equilibrium must be determined uniquely by that angular momentum chemical potential. That is, the spin polarization in equilibrium is fixed by the thermal vorticity, and their relation is one of the equilibrium thermodynamic properties.
On the other hand, since the system has a finite microscopic correlation length, the state of a local fluid element and its time evolution can only experience the local environment, and hence should be determined by the local thermal vorticity as well as other local thermodynamic parameters such as temperature and chemical potentials of conserved charges. This means that the fluid element in off-equilibrium relaxes to "local" equilibrium, with the thermal vorticity as one of the parameters characterizing local thermodynamic equilibrium. 1 In other words, spin polarization in the strict regime of hydrodynamics based on local equilibrium is not an independent variable, but is slaved to the conventional hydrodynamic variable of fluid velocity. Any deviation of spin polarization from its local equilibrium value is a non-hydrodynamic mode, which relaxes to zero with a finite relaxation rate that is determined by microscopic dynamics of the theory [19-21, 42, 43], similar to any other non-hydrodynamic modes. In a theory where the relaxation rate of spin polarization is much slower than other non-hydrodynamic modes, one can have two distinct hydrodynamic descriptions of spin polarization, depending on the time scale of interest compared to the relaxation time of spin polarization.
When the relevant time scale is much longer than the spin relaxation time, the spin is in local equilibrium with the thermal vorticity, and is not an independent degree of freedom. The system is described in a strict sense by hydrodynamics [44]. The true novelty of spin hydrodynamics in this regime is two-fold [13]: 1) the spin, or equivalently the fluid vorticity, affects the local thermodynamic laws used in hydrodynamics, as a second order gradient correction to the first law of thermodynamics, 2) the energy-momentum tensor has an anti-symmetric part which is proportional to the rate of change of the spin tensor, i.e., the fluid vorticity. One can formulate spin hydrodynamics in this regime based on the gradient expansion of thermal vorticity. In leading order it gives the ideal limit of spin hydrodynamics with no production of entropy [13], which is in accord with the existence of globally rotating equilibrium with a constant, but non-zero thermal vorticity. Interestingly, this ideal spin hydrodynamics has been shown to be equivalent, by a pseudo-gauge transformation [35,36,38,40], to conventional hydrodynamics based on the symmetric energy-momentum tensor with no spin degrees of freedom, with certain non-dissipative second order transport coefficients involving fluid vorticity [13]. The equivalence provides an important conceptual bridge between the two different formulations in the strict regime of hydrodynamics.
On the other hand, in the regime where the time scale of interest is comparable to or shorter than the spin relaxation time (but is still longer than relaxation time of other non-hydrodynamic modes), one can include spin polarization as an additional independent dynamical mode [11] in an extended hydrodynamic framework generally called Hydro+ [45]. The essential feature of this extended framework is that spin polarization as a non-hydrodynamic mode relaxes to its local equilibrium value, due to spin-orbit coupling inherent in the microscopic relativistic theory [11,21]. A new transport coefficient appears in this regime, dictated by the second law of thermodynamics [11,12,14,15]: the rotational viscosity, which determines the characteristic time scale for spin relaxation. In the present paper, we focus on this regime of spin hydrodynamics.
We address several important theoretical issues, founded upon the microscopic definition of the spin tensor in a background with torsion. First, the definition of spin current has been unclear in some of the previous works. Although the canonical spin current of Dirac fermions based on quantum field theory is totally anti-symmetric with respect to its three Lorentz indices, that was often not assumed in previous works. 2 Therefore, it is desirable to work with a microscopic definition of the spin-current operator from the underlying theory, possessing the appropriate anti-symmetric property of the spin-current. Second, the regime of applicability of relativistic spin hydrodynamics is rather subtle. A naive application of derivative expansion is problematic, and the proper expansion scheme needs to be clarified. Last but not least is the issue of how rapidly the spin density in a relativistic plasma relaxes to its equilibrium -one of the most crucial information relevant to the ongoing heavy-ion experiments. To answer this question, we need theoretical methods to evaluate the spin relaxation rate in the microscopic theory. One approach is the quantum kinetic theory with collisions [19-21, 42, 43]. In this work, we follow an alternative path that does not rely on weak-coupling approximation and leads to the Green-Kubo formula based on retarded correlation functions of spin observables, which is similar to those for other transport coefficients [46][47][48], but also with important differences due to non-conservation of spin density. 3 Our technique involves torsion which makes spin connection an independent source for the spin current and facilitates derivation of Green-Kubo formulae. 2 One can perform a pseudo-gauge transformation to bring the spin current into a desired form, which however changes the definition of the spin current. 3 Ref. [37] has recently derived a Green-Kubo formula based on the statistical operator approach, which has some overlap with our results.
Our approach is based on the combination of entropy-current analysis and linear response theory. Our starting point is the Ward-Takahashi identity for the local Lorentz symmetry in the underlying quantum field theory of Dirac fermions, such as QED and QCD. The first-order formalism in the background geometry with torsion naturally leads to a totally anti-symmetric spin current, which is not conserved due to spin-orbit coupling. As a result of the total anti-symmetric nature of the spin current, our relativistic spin hydrodynamics contains three non-hydrodynamic modes, corresponding to a spin vector that is the generator of spatial rotations.
We then perform the entropy-current analysis in a general background with torsion, which introduces a new transport coefficient, the rotational viscosity η s . The constitutive relations we find in a general torsionful background allow us to derive the Green-Kubo formula for the rotational viscosity, as well as other possible ways to compute it from the underlying quantum field theory.
Furthermore, the constitutive relations tell us that the local equilibrium value of the angular momentum chemical potential receives additional contributions coming from the background torsion. We study the linear response theory of our extended hydrodynamics with spin modes, and obtain the dispersion relations of dynamical modes. The relaxation rate of the spin modes is shown to be given by Γ s = 2η s /χ s , where χ s is the spin susceptibility (see (4.4) for the definition).
We consider the frequency scale ω Γ s , where these spin modes are treated as independent dynamical variables in our relativistic spin hydrodynamics. We also assume that ω Γ, where Γ is the relaxation rate of other non-hydrodynamic modes. This window of scale, where our extended hydrodynamics with additional spin modes is justified, is possible only when Γ s is parametrically smaller than Γ (see Figure 1). For example, in a weakly coupled quark-gluon plasma, when the mass M of fermions is much larger than the temperature T , Γ s is expected to be smaller than Γ by additional powers of T /M .
One of the interesting features we find in our linear response theory is the mixing between the transverse shear mode and the spin mode, that can be seen in Figure 1. This mode mixing arises because the transverse shear gradient ∂ z u x is a linear combination of shear tensor and fluid vorticity, the latter of which couples to spin modes.
The paper is organized as follows: In section 2, we review the definition of energymomentum tensor and spin current in quantum field theory of Dirac fermions, together with the associated Ward-Takahashi identity of local Lorentz symmetry. In section 3, we derive the constitutive relations of relativistic spin hydrodynamics in a torsionful curved background, where spin density is treated as an independent dynamical variable. In section 4, we develop the linear response theory of our spin hydrodynamics, which gives the dispersion relation of dynamical modes, as well as the derivation of the Green-Kubo formula for the spin relaxation rate. Section 5 is devoted to summary and discussion. In Appendices A and B, we present a discussion on the angular momentum in global equilibrium under rotation and the relation between the currents in the vierbein-spin connection formalism and the metric-affine connection formalism, respectively.
o v d f e 9 X 9 H s 1 / q f G q t r z C 3 o 0 c j 2 T f / n W V e f V w 8 O n 6 w 6 G y 2 t d t T R X y f 2 b z U M S c n 0 n n j L b P N N N q L X / 1 8 P w l O 5 8 Z r 0 / Q F T 1 z z k t q 0 D 0 n N a u v 2 n V a Z C 4 Q 5 k H J 3 8 f y / 4 y y c X r 7 n 0 + k R 1 k l 2 z F 8 p 5 x e r s g Z L a 5 T f 9 Z o 2 v 1 x C j Q c k / x / I b b M 6 k 5 L n U z N p s M r M U j a w L Y x j H F M 1 l H h k s Y x X 5 8 H b n u E R N G p H S 0 q K U b U i l l s g z j G 8 l r X w C m H + U V g = = < / l a t e x i t > ⇡ = wave number < l a t e x i t s h a 1 _ b a s e 6 4 = " K I 2 s + J 4 m m d e z + y u 6 2 U j f 9 B z h 5 M + p J E w / G P 8 O L / 4 A H / g E T 4 x E T L x 5 4 8 z o E K 6 X d a d p 9 8 + Z 9 v n 1 v 3 u y I w F F R z N j + W O b M u H X 2 3 M T 5 y Q t T F y 9 N z 8 x e L k V + K 7 R l 0 f Y d P 9 w S P J K O 8 m Q x V r E j t 4 J Q c l c 4 s i y a j / R 6 u S 3 D S P n e 0 7 g T y J r L G 5 6 q K 5 v H 6 C p V H 3 P X 5 d s z W Z Z j N O Z P G n l j Z M G M d X 8 t 4 7 8 4 + A 2 R n I d F 9 p 1 9 Z g f s G / v C f r K / p 2 o l p K F z 6 e B T 9 F g Z b E / v X d 3 8 M 5 J y 8 R n D 8 2 N q C C E w m u K e 3 a 6 U h 9 Y W Q x 3 u U 0 0 K a w z I o 6 u 1 e 3 z 7 1 Z u D z Q e F x e Q m + 8 h + Y Z 0 f 2 D 7 7 i p V 6 7 d / 2 p w 1 Z e E / q H j I v a V d c y t P D P i T U D Y X 6 d d T S s 6 M K Q r Q S X O u a f T + N 1 n F q I K t S s X I g K 0 e y H N q 4 E p p d 6 t L J q d G K Q 1 E N O j N Z 9 F Y p 6 s i j V Y f p C j p F w m Q l K J u q 8 W n l U X k J 7 N y / b D M F s d t H 7 K Y g O n 1 E J w W x 1 k e s G W J 4 d y R 2 Z F B 3 C i m 6 s 4 P 9 G c Q + H M n a 9 P 4 c 9 1 P P J b 0 R C f X T Z I 5 3 T / 7 / m + a k U V r K 5 e / m l j f u Z F d W z S 0 0 A d f g B t x C p X u w A k 9 g H Y q o / g J e w 1 t 4 N / 7 D m r L m r C u 9 0 M y Y Y e a g b 1 j X D w G L z y P n < / l a t e x i t > Figure 1. A schematic picture of frequency scales as well as sound, shear, and spin density modes. The spin hydrodynamics discussed in this work extends the regime of validity of hydrodynamics to ω Γ s , where Γ s is the relaxation rate for spin density. For ω Γ s , the spin is no longer an independent degree of freedom, and though spin hydrodynamics still applies, we also have pure hydrodynamics as an alternative and simpler description. The green and red lines show a mode mixing phenomenon between transverse shear and spin modes at the scale ω Γ s .

Spin current in the first-order formalism
In this section, we review the definition of the energy-momentum tensor and the spin current from the viewpoint of the first-order (or Palatini) formalism for background spacetime [49]. We consider quantum field theory (QFT) in a torsionful (Einstein-Cartan) background geometry [50] and introduce currents and the Ward-Takahashi identities [51,52] associated with diffeomorphism, local Lorentz invariance and flavor symmetry.

QFT in a torsionful background
When the system enjoys a continuous global symmetry, it brings about a conserved current thanks to Noether's theorem [53]. Nevertheless, the naive definition of the Noether current often contains conceptual problems, an example of which is the gauge non-invariance of the canonical energy-momentum tensor for gauge theories. Although these problems may be resolved in an ad hoc way by using, e.g., the Belinfante improvement [54][55][56], it is not clear whether the improved one is the best definition of current operators (see, e.g., Refs. [57,58] and references therein for reviews on the definition of spin currents in gauge theories).
In order to avoid potential problems in the Noether current, one can instead utilize the definition of the gauge current in the following systematic way: we introduce the background gauge field that couples to the global symmetry of the system, and by taking the variation of the gauged action with respect to the background gauge field, one immediately gets a well-behaved current operator. A famous example is the symmetric and gauge-invariant energy-momentum tensor in general relativity, which is canonically defined by the variation of the matter action with respect to the metric.
In the same manner, we can define the spin current operator, although it is not a conserved current due to its possible exchanges with orbital angular momentum. The key concept here is the spin connection, which enables us to promote the Lorentz invariance to a local one. Let us consider QFTs in a background of torsionful curved spacetime, whose action reads as where we introduce a set of dynamical fields ϕ and background fields j, respectively. As an explicit example, we will consider QCD in this paper, but the whole discussion also applies to QED or any other theories including Dirac fermions. We consider the following QCD Lagrangian density is the field strength tensor for the gluon field a µ ≡ a a µ t a , with su(3) c generators t a , and the QCD coupling constant g. The dynamical fields are quarks and gluons: ϕ = {q, a µ }, and corresponding to global symmetries of QCD, a set of background fields is given by j = {eâ µ , ωâb µ , A µ }, where eâ µ (e μ a ) and ωâb µ = −ωbâ µ denote the (inverse) vierbein and spin connection with e ≡ det eâ µ , respectively, and A µ is a flavour gauge field, coupled to, e.g., the U(1) baryon symmetry (generalization to non-Abelian flavour symmetry is straightforward). We note that the anti-symmetry of ωâb µ implies metric compatibility, ∇ α g µν = 0. We also introduced the covariant derivative of the quark field as where γâ are the Dirac gamma matrices, and Σâb ≡ i[γâ, γb]/4 are the generators of the Lorentz group in the spinor representation. In this paper, we use Greek (hatted Latin) letters for coordinate (local Lorentz) indices. The Minkowski metric convention is chosen to be the mostly plus one: ηâb = diag(−1, +1, +1, +1).
It is worth emphasizing that the fundamental variables describing the background geometry are the vierbein and the spin connection (or contorsion which will be introduced shortly). We consider a torsionful background, where the vierbein and the spin connection are independent background fields. One can then define the space-time metric g µν and the affine (curved space-time) connection Γ ρ µν through where the second equation, or the so-called tetrad postulate, defines the affine connection from the given vierbein and spin connection. The anti-symmetric part of the tetrad postulate gives the first Cartan equation, which defines the torsion T ρ µν ≡ Γ ρ µν − Γ ρ νµ as follows: The introduction of the torsion is reasonable because it promotes the spin connection to an independent background field. This becomes manifest by solving the tetrad postulate for the spin connection as where we introduced the torsion-free spin connectionωâb µ (e) and the contorsion Kâb µ as follows [50]:ωâb where C µνρ is the Ricci rotation coefficient given by Note that the torsion-free partωâb µ (e) is fixed by the vierbein. Thus, introducing the spin connection as an independent field is equivalent to putting the system into a background with torsion. One observes that the quark couples to the background torsion through the spin connection.
It is worthwhile to emphasize that the gluon (or photon in the case of QED) is assumed to have no coupling to the torsion, in contrast to the quark field. This is a consequence of imposing SU(3) c gauge invariance of the action [50]. Note that we define the gluon field strength tensor with a simple partial derivative, not with the covariant derivative. At first glance, this choice does not seem compatible with diffeomorphism invariance, but it indeed is. To see this, notice that both covariant derivatives with or without torsion are consistent with diffeomorphism invariance, since their difference is proportional to the torsion tensor. However, one sees that the use of the covariant derivative with torsion for the field strength tensor spoils SU(3) c gauge invariance due to the last term appearing in Throughout the paper, we use ∇ µ as the covariant derivative including only the affine connection, and D µ as the one including both the affine and spin connections. It should also be pointed out that we have the symmetrized derivative in the fermion kinetic term in Eq. (2.2). This symmetrized form ensures that the Lagrangian is real-valued, and as a consequence, the resulting spin current is a Hermitian operator as it should be. These choices are crucial for the proper definition of the spin current, or equivalently, for our decomposition of the total angular momentum into spin and orbital angular momenta.

The Ward-Takahashi identities
The systems enjoying the Poincaré invariance and some flavour symmetries are equipped with conserved charges, given by the energy, momentum, total angular momentum, and flavour charges. Since we would like to discuss the dynamics of spin polarization, we need to identify the proper definition of the spin current, in order to decompose the total angular momentum into the orbital and the spin parts. The decomposition of the angular momentum tensor in quantum gauge theories is still a controversial subject (see, e.g., Refs. [57,58] and references therein for a review). 4 Based on the background vierbein and spin connection, we here introduce one way to define a gauge-invariant spin current, which eventually is shown to involve only the fermion sector. First of all, note that the introduction of background vierbein and spin connection promotes the Poincaré symmetry to local symmetries -diffeomorphism and local Lorentz invariance. As a consequence of these local symmetries, one can derive in the usual manner the covariant (non-)conservation laws as the Ward-Takahashi identities. The spatial components of the associated spin current, i.e. Σ ijk , turn out to be completely determined by the axial charge density J 0 5 . Also, it gives us a useful prescription to define the spin current through the variation of the QCD action (2.1)-(2.2) with respect to the spin connection. One can regard this procedure as a familiar way to define the gauge current. We thus define the "canonical" energy-momentum tensor Θ µâ and the spin current Σ µâb , as well as the flavour current J µ , by taking the variation with respect to eâ µ and ωâb µ , as well as A µ , respectively: (2.11) where subscripts ω and e indicate that they are fixed when we compute the functional derivative.
After a straightforward calculation, we find that the QCD action (2.1) gives us the following expression for our currents where {A, B} ≡ AB + BA is the anti-commutator. One finds that only the local (or on-site) quark field contributes to the spin current, so that our spin current is manifestly SU(3) c gauge-invariant. Moreover, thanks to the anti-commutator, Σ µâb are Hermitian gauge-invariant operators, as all physical observables should be. In other words, the spin current would not be a Hermitian operator, if the fermion kinetic term in the Lagrangian was not symmetrized. The absence of gluon contributions to the spin current results from the absence of the torsion in the gluon field strength tensor due to gauge invariance. In that sense, our requirement of SU(3) c gauge invariance in torsionful geometry fixes our choice of the decomposition of the angular momentum. There is another useful expression for the spin current, thanks to a gamma matrix identity {γ µ , Σ νρ } = ε µνρσ γ σ γ 5 with γ 5 ≡ −iγ0γ1γ2γ3, where the normalized totally antisymmetric tensor is defined as ε µνρσ ≡ µνρσ /e with 0123 = +1. This formula tells us that the spin current is totally anti-symmetric with respect to its three indices. One can then easily show that it is equal to the Hodge dual of the axial U(1) A current Jâ 5 as Therefore, the number of independent components of the spin current operator is only four, which are in one-to-one correspondence to the U(1) A current J µ 5 . This means that some components of the spin density, e.g. Σ 00i , are identically zero, and the other components such as Σ 0ij are equivalent to the spatial components of the axial current J i 5 . The remaining component of Σ ijk is completely determined by the axial charge density J 0 5 . Let us now derive the Ward-Takahashi identities that follow from local symmetries. The vital point is that the QCD action (2.1) remains invariant, i.e. δ χ S QCD = 0, under the following infinitesimal transformations: (2.14) We introduced a set of local infinitesimal parameters, χ ≡ {ξ µ , αâb, θ}, which generates the general coordinate and local Lorentz transformations -four-vector ξ µ for diffeomorphism, anti-symmetric tensor αâb = −αbâ for local Lorentz symmetry, and θ for the U(1) symmetry. On the other hand, a direct computation of the induced variation of the QCD action leads to another expression for δ χ S QCD with the arbitrary parameters χ. For δ χ S QCD = 0 to hold for arbitrary χ, their coefficients need to vanish, which gives us the Ward-Takahashi identities for the currents. Since the computation is a little complicated, we first consider the identities associated with the local Lorentz and U(1) symmetries, and then move on to that for diffeomorphism. The variation of the action under the local Lorentz transformation is given by where we used the equation of motion for dynamical variables (δS QCD /δq = δS QCD /δq = δS QCD /δa µ = 0). We here introduced a torsional contribution G µ ≡ T ν νµ , which results from the integration by parts (all surface terms are neglected) in the presence of the torsion as The local Lorentz invariance of the action tells us δ α S QCD = 0 holds for arbitrary positiondependent αâb. As a result, we obtain the following identity: Therefore, based on Eq. (2.17), the covariant divergence of our spin current is tied to the anti-symmetric part of our energy-momentum tensor. We regard this Ward-Takahashi identity as the equation of motion for the spin density. However, we also note that not all the components of Eq. (2.17) give equations of motion for the spin density. This is because our spin current is totally anti-symmetric, and thus, we only have three spin components as dynamical variables. As a result, the remaining three equations are shown to be constraints for a part of the energy-momentum tensor Θ µâ . We will use these constraints when we derive the constitutive relations with the help of the second law of local thermodynamics. A similar computation leads to the following identity attached to the U(1) symmetry: which gives a conservation law for the flavour current.
Being equipped with the above identities, we further evaluate the variation of the action induced by the general coordinate transformation as where we again used the equation of motion and performed the integration by parts. Here, we introduced the field strength tensor R a bµν attached to the spin connection, and F µν as (2.20) For δ ξ S QCD = 0 to hold for arbitrary ξ ν , we obtain the identity Therefore, our energy-momentum tensor has source terms coming from the spacetime torsion, curvature, and external U(1) gauge field. In summary, we have obtained the following Ward-Takahashi identities: which result from diffeomorphism, local Lorentz invariance and U(1) symmetry of the QCD action (2.1)-(2.2). When we are interested in QCD in the flat spacetime with vanishing U(1) gauge field, we turn off the background fields eâ µ , ωâb µ , and A µ , which reduces Eqs. (2.22) to These equations give the conservation laws for the canonical-like energy-momentum density, total angular momentum, and flavour charge (the first two were used in Ref. [11]). We note that the spin current Σ µâb is totally anti-symmetric with respect to its indices owing to the gamma matrix structure of Eq. (2.12) [recall Eq. (2.13)], which is not assumed in Ref. [11]. As a result of the anti-symmetric property, the number of independent degrees of freedom in the spin current is not six but three in the present work.
One remark here is that we defined the energy-momentum tensor Θ µâ by considering the variation with the fixed spin connection ωâb µ instead of the fixed contorsion Kâb µ . As a consequence, one finds that Θ µâ corresponds to the canonical-like energy-momentum tensor. If we introduce the energy-momentum tensor T μ a with the fixed contorsion Kâb µ , we have where we used the totally anti-symmetric property of the spin current to obtain the rightmost side. We will see that this energy-momentum tensor becomes useful when performing the linear-mode analysis in Sec. 4.

Derivation of spin hydrodynamics
In this section, we provide a derivation of hydrodynamics with spin polarization [11,12,14,15] (see also [28]). Based on the Ward-Takahashi identities and the second law of local thermodynamics [44], we specify a first-order constitutive relation of hydrodynamics with the totally anti-symmetric spin current. In addition to two (shear and bulk) viscosities and a flavour charge conductivity, we find one kinetic coefficient, called rotational viscosity, appearing in the anti-symmetric part of the energy-momentum tensor.

Setup
Let us first introduce the dynamical variables. We employ the Landau-Lifshitz frame with respect to the symmetric part of the energy-momentum tensor in defining the fluid fourvelocity u µ . Therefore, we introduce u µ and the associated energy density in the rest frame by where we denote the symmetric part of the energy-momentum tensor as Θ µâ (s) ≡ (Θ µâ + Θ μ a )/2. Owing to the normalization u µ u µ = −1, we can parameterize the four-velocity as We next introduce the spin density σâb = −σbâ and the flavour charge density n by σâb = −u µ Σ µâb and n = −u µ J µ . (3.2) As given in Eq. (2.12), we microscopically know the spin current for QFTs with Dirac fermions is totally anti-symmetric so that the spin density for such QFTs satisfies σâbuâ = 0. This means that we only have three independent components of spin density attached to the spatial rotation symmetry, while the other three associated with the boost symmetry are absent. Thus, our dynamical variables are three spin components in addition to four energy-momentum densities and one flavour charge density. As for the spin density, we also use the dual variable σâ defined by In short, we identify eight dynamical variables as , n, uâ, σâ (or σâb) satisfying ηâbuâub = −1, σâuâ = 0 = σâbuâ. We then introduce thermodynamic variables conjugate to the energy, spin, and flavour charge densities. For that purpose, we introduce a generalized entropy density s( , n, σâb), or equivalently s( , n, σâ), and introduce the inverse temperature β, flavour chemical potential µ, and spin chemical potential µâb (or its dual µâ) by In other words, we require the generalized first law of local thermodynamics in the rest frame including spin density: Note that the spin chemical potential satisfies µâuâ = 0 = µâbub. Our purpose is to express Θ µâ , Σ µâb , and J µ by using the eight dynamical variables listed in Eq. (3.4). We then rewrite these currents by performing the tensor decomposition as follows: where we defined a pressure p with a projection matrix ∆ μ a : ∆ μ a ≡ e μ a + u µ uâ satisfying ∆ μ a uâ = 0 = ∆ μ a u µ . (3.9) Note that the Landau-Lifshitz definition of the fluid velocity restricts δΘ µâ and δqâ to satisfy δΘ µâ uâ = 0 = u µ δΘ µâ and δqâuâ = 0, and the definition of n leads to u µ δJ µ = 0.
Equation (3.8) just gives one possible parameterization of Θ µâ , Σ µâb and J µ , and thus, δΘ µâ , δqâ, δσ, and δJ µ could contain nondissipative terms. We will, however, show that such a nondissipative term is absent in the constitutive relation for δΘ µâ , δqâ and δσ in the first-order spin hydrodynamics. It is worth emphasizing that the frequency scale of our main interest is in the spin hydrodynamic regime, at which the spin density shows its intrinsic dynamics as an independent dynamical variable. As we will show in Sec. 4.1, the spin density shows a relaxation dynamics with its characteristic relaxation rate Γ s . Thus, we are now focusing on the spin hydrodynamic regime ω = O(Γ s ), which allows us to find balanced terms in the resulting equation of motion of the spin density. We also note that we implicitly assume presence of a scale separation between spin and the other non-hydrodynamic modes characterized by the frequency scale Γ as Γ s Γ. This additional assumption is not what we can show in general setups, but what we expect to be true in a specific case such as the heavy-fermion mass limit.

Derivation of constitutive relations
To derive hydrodynamic equations with spin polarization, we need to express δΘ µâ , δqâ, δσ and δJ µ using our dynamical variables and background fields. For that purpose, we require a generalized second law of local thermodynamics, which require the existence of the entropy current ∃ s µ (β, µ, u µ , µâb) satisfying (∇ µ − G µ )s µ ≥ 0 for ∀ β, µ, u µ , and µâb. (3.10) Let us specify the condition following from the local second law (3.10). By parametrizing the entropy current as s µ = su µ + δs µ , we can express its divergence as The vital point here is that we can write the material derivative u µ ∇ µ of the entropy density s( , n, σâb) as 12) where we used the definition of the thermodynamic variables β, µ, and µâb. Noting that the left-hand side of this equation is a scalar, we put the covariant derivative D µ σâb to make the right-hand sides to be scalar too. We can rewrite the rightmost side of the equation (3.12) by using the Ward-Takahashi identities with the assumed constitutive relation (3.8). To do so, we take contractions of Eqs. (2.22) with βâ = βuâ, 1 2 βµâb and βµ, which results in where we used δΘ µâ uâ = 0 = δqâuâ to obtain the first equation. At this point, δqâ seems to enjoy the same status as δΘ µâ , δσ and δJ µ , whose constitutive relation should be determined independently of other parts. However, this is not the case: the constitutive relation for δqâ is determined by that for Σ µâb . To show this, recall that the Ward-Takahashi identity for the spin current contains three constraints rather than six dynamical equations due to the totally anti-symmetric property of our spin current. This becomes manifest by contracting Eq. (2.17) with ub and noticing that the resulting identity does not contain the time derivative in the right-hand side. As a consequence, we only need to derive the constitutive relation for δΘ µâ , δσ and δJ µ , the second of which will determine that for δqâ in accordance with the constraint (3.16). It is thus useful to rewrite (3.13) by substituting the constraint (3.16) as (3.17) The result so far is exact without assuming the derivative expansion. Let us then apply the derivative expansion to derive relativistic spin hydrodynamic equations, which capture leading dissipative dynamics of charge densities. For that purpose, we first need to specify our power counting scheme for the dynamical variables and background fields. In this paper, we employ the power counting scheme defined by O(∂ 0 ) = {β, u µ , µ, eâ µ } and O(∂ 1 ) = {µâb, σâb, ωâb µ }, (3.18) and derive the constitutive relations up to O(∂ 1 ). This allows us to keep the restricted number of terms appearing in the constitutive relations. To demonstrate this, using Eqs. (3.13)-(3.15), we rewrite the divergence of the entropy current as Here, we also introduced a collection of the higher-order terms ∆S = O(∂ 3 ), whose explicit form is given by We shall explain the reason why these terms are counted as O(∂ 3 ). First of all, our power counting scheme (3.18) immediately tells us the terms in the first and second lines are the third order in derivative. The first term in the first line, if considered, leads to a secondorder correction of the Gibbs-Duhem relation as + p = T s + µn + 1 2 µâbσâb (see also the following discussion). From the third line, one finds the form of δσ compatible with the second law as δσ = −T κ s ε µâbĉ uĉ D µ (βµâb) + Râb µd βd + D µ ub ∇âβ + u ν (D ν βâ − βdTdâ µ ) , (3.22) with a positive coefficient κ s > 0. Equation (3.22) represents a diffusive transport of the spin current and its response to the curvature tensor; namely, κ s is identified as a spin conductivity. Besides, this closes the constitutive relation for δqâ thanks to the constraint (3.16). Our power counting scheme tells us δσ = O(∂ 2 ) and δqâ = O(∂ 2 ), which proves that all terms in the last two lines of Eq. (3.21) are O(∂ 3 ). Note that the last term, if considered, leads to a second-order correction of the entropy current as δs µ = −βδq µ − βδσµ µ = O(∂ 2 ). As a result, we find δσ = 0 and δqâ = 0 in the first-order spin hydrodynamics.
Going back to Eq. (3.19), one finds the entropy current s µ and the constitutive relations for δΘ µâ and δJ µ compatible with the generalized second law of local thermodynamics (3.10) in the following way. Firstly, the first line of Eq. (3.19) vanishes by identifying where the first equation is the usual Gibbs-Duhem relation up to O(∂). We note that both equations have second-order derivative corrections present in ∆S as we discussed above.
The second line and the first term in the third line of Eq. (3.19) will always be positive by identifying the following constitutive relations: where one can decompose the viscous tensors η µ ν ab and (η s ) µ ν ab as η µ ν ab = 2η Here, all four kinetic coefficients η, ζ, η s , and κ n are assumed to take positive values. They are shear viscosity, bulk viscosity, rotational viscosity, and flavour conductivity. The positivity of these kinetic coefficients makes the entropy production rate a positive semidefinite quadratic form within the second-order derivative to satisfy the second law of local thermodynamics.
One can also rewrite the constitutive relation by decomposing the spin connection to the torsion-free and torsion parts as where the circled derivative contains the torsion-free spin connection. From this decomposition, one sees that the symmetric part of the energy-momentum tensor has no contribution from the torsion. This is because the contorsion Kâ bµ is anti-symmetric under the exchange of (µ,â) indices so that it disappears from the constitutive relation. Then, we obtain another expression for the constitutive relation as follows: (3.27) Equations (3.24) or (3.27) give the constitutive relations for the first-order spin hydrodynamics. In addition to two viscosities η and ζ, we have another viscosity known as the rotational viscosity η s . In summary, we find the constitutive relations for the first-order spin hydrodynamics in the torsionful curved background as follows: where the rank-four viscous tensors are decomposed as Eq. (3.25). Thus, the first-order spin hydrodynamics with one flavour current is equipped with three viscosities η, ζ, and η s and one conductivity κ n . Once we know their values (and also the equation of state), the Ward-Takahashi identities form a closed set of equations for the energy-momentum, spin, and flavor charge densities. This completes the phenomenological derivation of hydrodynamics with spin polarization. A few remarks are in order here. The present analysis uses the same strategy as Ref. [11], but we have some modifications that originate from the total anti-symmetry of the spin current. This difference eliminates three dynamical spin densities attached to the boost, and, as a result, one of the transport coefficients, the boost heat conductivity, found in Ref. [11], disappears. Also, the power counting scheme is a little different from that used in Ref. [15], which makes it difficult to give a direct comparison with their constitutive relation (see the discussion at the end of Sec. 4.1). 5 We here emphasize that our identification of the spin hydrodynamic regime allows a consistent truncation for the spin hydrodynamic equations. To see this, let us consider the Ward-Takahashi identity in the flat spacetime ∂ 0 Σ 0âb +∂ i Σ 0âb = −(Θâb −Θbâ). According to our power counting scheme (3.18), the first term on the left-hand side is counted as O(ω∂) while the right-hand side will be identified as O(Γ s ∂) (see Sec. 4.1). These are precisely balanced when we consider the spin hydrodynamic regime ω = O(Γ s ). Lastly, one could try to use the novel derivative expansion scheme introduced in Ref. [13], treating fluid vorticity and temperature gradients differently from other dissipative first order gradients such as shear and bulk tensors. This allows to rigorously extend the entropy-current analysis to second or higher order terms that involve vorticity and temperature gradients, producing new constraints for second or higher order transport coefficients [13]. The application of this scheme to the present analysis may, for example, justify (3.22) more rigorously.

Linear response theory for spin hydrodynamics
In this section, we present the linear response theory of hydrodynamics with spin polarization. The purpose of this section is twofold: the first one is to specify the spectrum (or dispersion relation) of the modes described by spin hydrodynamics [59]. We will find that the rotational viscosity controls the (gapped) relaxation rate of the spin density. The second one is to find the field-theoretical expression for the rotational viscosity. Applying the linear response theory with the background (mechanical) perturbation [46][47][48]60], we provide the Green-Kubo formula for the rotational viscosity. Since the flavour currents (e.g., the baryon current) do not affect the dynamics of the spin density, we neglect them in the following discussion.

Linear-mode analysis and non-hydrodynamic spin mode
Let us consider a spinful relativistic fluid in flat space-time, where the Ward-Takahashi identities and the constitutive relations are simplified as It is not necessary to distinguish the curved and local Lorentz indices in flat space-time, so we remove hats over the indices in the following analysis. We perform the linear-mode analysis of the spin hydrodynamic equations (4.1) around the global equilibrium without background flow and spin density, that is, we consider small perturbations given by where we introduced a set of static/kinetic coefficients as The χ s is the spin susceptibility, and c s is the speed of sound. We work in the Fourier space with δÕ(k) ≡ d 4 x e iωt−ik·x δO(x), and we choose the momentum direction to be along the z-direction, k = (0, 0, k). Then, the linearized equations of motion become the block diagonal form as where we introduce the parallel and perpendicular matrices as Solving the characteristic equations, det A 3×3 = 0 and det A 4×4 ⊥ = 0, we find the dispersion relations for the modes in our spin hydrodynamics.
The longitudinal modes, described by A 3×3 , decompose into the decoupled sound and spin modes, whose dispersion relations read • One longitudinal spin mode : ω spin, (k) = −iΓ s .

(4.7)
It is clear that we have the relaxational longitudinal spin mode in addition to the gapless sound mode (see Figure 2(a)). In contrast, the transverse momentum (shear) and spin modes are coupled together, which gives rise to two pairs of the following dispersion relations (see Figure 2(b)): (4.8) Here, we call each branch as the shear and the transverse spin mode, because of their low-momentum behaviors • Two shear modes : ω shear (k) = −iγ ⊥ k 2 + O(k 4 ), • Two transverse spin modes : ω spin,⊥ (k) = −iΓ s − iγ s k 2 + O(k 4 ). (4.9) From the low-momentum behaviors in (4.7) and (4.9), we identify three non-hydrodynamic (gapped) spin modes, in addition to the hydrodynamic (gapless) sound and shear modes. All spin modes show the relaxation behavior with the characteristic time scale of Γ −1 s . This non-hydrodynamic behavior is originated from the fact that the spin angular momentum is not a conserved quantity. The relaxation of spin modes persists until the source term in the spin equation of motion vanishes, namely when the spin potential µ ab coincides with the fluid vorticity when it reaches the local equilibrium. In this strict hydrodynamic regime of ω Γ s , we thus find only the gapless sound and shear modes, as expected in the conventional relativistic hydrodynamics.
On the other hand, in the regime of our spin hydrodynamics (ω ∼ Γ s ), three components of spin density start to behave as independent dynamical modes. In contrast to the

ω(k)
Re ω shear (k) Im ω shear (k) Re ω spin,⟂ (k) Im ω spin,⟂ (k) Re ω sound (k) Im ω sound (k) Re ω spin,|| (k) Im ω spin,|| (k) decoupled longitudinal spin mode, the transverse spin modes feature an interesting coupled dynamics with the transverse shear mode. In fact, in high momentum region, the shear and transverse spin branches in the dispersion relations behave as It is notable that they interchange their momentum dependence between the low and the high momentum regimes: the shear branch approaches to a constant relaxation rate independent of momentum, while the spin branch shows a quadratic dependence on momentum (see Figure 2(b)). At the cross-over momentum scale, we see the level repulsion phenomenon between the two branches. We should also point out that the eigenvector of each branch is a linear combination of δπ x and δσ y (or δπ y and δσ x ), due to the off-diagonal mixing terms in the matrix A 4×4 ⊥ in Eq. (4.6). This mixing is a result of the fact that the transverse velocity gradient ∂ z δv x is a linear combination of shear tensor and vorticity, and the former couples to momentum fluctuation while the latter couples to the spin density. In particular, the strict hydrodynamic shear mode in low momentum is described by the energy-momentum tensor T μ a defined in Eq. (2.24), which is a linear combination of our energy-momentum tensor Θ µâ and the spin density. In fact, the flat space-time limit of Eq. (2.24) gives us and we find that the transverse fluctuation, e.g., δ T 0 x , is one of the eigenvectors in low momentum regime satisfying This shows that T 0â describes the gapless hydrodynamic modes in the strict hydrodynamic regime. However, we emphasize that we need Θ µâ to describe the dynamics of spin polarization, because its anti-symmetric part appears in the equation of motion for spin polarization.
The zero-momentum relaxation rate of the spin modes is controlled by Γ s = 2η s /χ s , and thus, the rotational viscosity η s is essentially the spin relaxation rate. Note that the finite-momentum part of the spin relaxation rate is different between the longitudinal and transverse components: the former does not depend on momentum, while the latter is given also by the rotational viscosity η s as γ s k 2 = ηs 2( 0 +p 0 ) k 2 . If we consider the nextleading-order correction of the spin current, the relaxation rate for the longitudinal spin density also shows a momentum dependence. While the finite-momentum dependence of the spin relaxation is an interesting subject, we focus on the zero-momentum behavior in the remainder of this paper. We also note that the parameter set used in Figure 2 (and also in Figure 1) is chosen to have a well-defined spin hydrodynamic regime. This is manifest in our choice of small Γ s , compared to the scale where hydrodynamic gradient expansion breaks down, i.e., ω c ∼ c 2 s /γ . This is a required assumption that ensures the existence of the spin hydrodynamic regime (or the spin hydro+ regime in Figure 1), in which only spin modes appear as an additional non-hydrodynamic mode in the time scale much larger than that of other non-hydrodynamic modes. This may be justified when, e.g., there exist sufficiently heavy fermions in the system.
Let us translate the above results into field-theoretic language. For that purpose, we introduce the retarded Green's function for arbitrary operatorsÂ andB as where Ô eq denotes the expectation value in homogeneous thermal equilibrium. Then, the dispersion relation (4.7) and (4.9) implies that the retarded Green's function for the spin density G σ i σ j R (ω, k = 0) has a pole at the location ω spin (k = 0) = −iΓ s . Besides, taking into account that the zero-frequency limit of the retarded Green's function gives the spin susceptibility lim we identify the low frequency-momentum behavior of the retarded Green's function for spin densities as where the ellipsis stands for higher-order corrections with respect to frequency and momentum, which is beyond our first-order spin hydrodynamics. As a result, we can directly extract the value of the spin relaxation rate Γ s , if we identify the location of the pole in the retarded spin-spin correlator. This give one direct way to evaluate the spin relaxation rate from the correlation function in quantum field theory. Before closing this section, we comment on our power counting scheme. Noting that the equation of motion for the spin density is given by Eq. (2.17), one might worry whether our constitutive relation is consistent with (naive) derivative expansion scheme. This is a natural concern, because the left-hand side of Eq. (2.17) is of second order in derivatives, while the right-hand side is of first order. However, the crucial point is that we need to distinguish spatial derivatives from time derivatives in the spin hydrodynamic regime we are working on, since our frequency scale is implicitly assumed to be ω Γ s , and therefore, time derivative should not be counted equally to spatial derivative. Then, there is no discrepancy in the number of derivatives between both sides of Eq. (2.17): they are both O(Γ s ∂), where ∂ denotes a spatial derivative. This point will also be important when we check the consistency between the linear-mode analysis and the Green-Kubo formula in the next section.

Green-Kubo formula in the first-order formalism
As shown in the previous section, the rotational viscosity η s gives the relaxation rate of spin polarization. Our primary question is how to evaluate it from an underlying QFT such as QCD. We formulate in this section the linear-response theory based on the first-order formalism, which gives the Green-Kubo formula to evaluate the spin relaxation rate from the field-theory Green's functions.
First of all, we note that the spin density is a non-hydrodynamic mode with a characteristic relaxation rate Γ s , and therefore, it disappears in the strict hydrodynamic regime ω Γ s in our setup, where the local equilibrium state does not support vorticity. Thus, when writing the equations of spin hydrodynamics, we implicitly assume that we are working on the frequency scale given by ω Γ s (recall Figure 1). In other words, depending on the frequency scale of interest to us, the constitutive relation of the spin hydrodynamics is further reduced to be the usual hydrodynamics as where Γ represents the relaxation scale of other non-hydrodynamic modes. Here, we specified the constitutive relation in the true hydrodynamic regime ω Γ s by equating the source term of the spin equation of motion with zero. As a result, the small frequency limit which we will use below is not the strict ω → 0 limit, but should be regarded as a constrained limit Γ s ω Γ. This is in sharp contrast to the usual ω → 0 limit used in the Green-Kubo formulae for other transport coefficients for conserved charges.
We then introduce a mechanical perturbation [60], and identify the rotational viscosity from the linear response of the system to the background perturbation. The basic strategy is to compute the perturbation of currents in first order of the background torsion, and to match the result with the constitutive relation in our spin hydrodynamics, (3.24) or (3.27). From the first equation of (3.27), one sees that the bulk and shear viscosities are captured by the torsion-free spin connection, which results in the usual Green-Kubo formula for them. Thus, we focus on the rotational viscosity in the following discussion.
The vital point for the subsequent discussion is to identify Eqs. (3.24) or (3.27) as a non-equilibrium expectation value of the quantum operator δΘ µâ in a torsionful curved background. We here use the vierbein and contorsion as independent sources, while it is also possible to use the vierbein and spin connection, instead. Similarly to the previous section, we assume that the system is initially in homogeneous thermal equilibrium with β = const., uâ = (1, 0) t , and µâ µ = 0. On the time scale of ω Γ s , the spin chemical potential µâ µ remains zero even after introducing background perturbations, as we explained in the previous paragraph. This simplifies the constitutive relation for the anti-symmetric part of the energy-momentum tensor as This equation serves as a basis to obtain the Green-Kubo formula for the rotational viscosity η s . Below, we present two ways to compute the rotational viscosity, whose equivalence will be shown by the Ward-Takahashi identity. Firstly, we turn off the background vierbein by putting eâ µ = δâ µ , while keeping the non-vanishing contorsion, so that one has only the second term in Eq. (4.17). On the other hand, we can also compute the same δΘ µâ e.K in a different manner based on QFT. To describe this in the underlying quantum theory, we start from the initial density operator at past infinity t = −∞ given by the equilibrium Gibbs distribution whereĤ denotes the Hamiltonian of, e.g., QCD in the absence of background fields. Then, noting that Eq. (4.17) indicates that the contorsion Kb 0ν plays a role of a relevant background field, we perturb the system by adiabatically turning on the small contorsion and compute the induced response of the anti-symmetric part of Θ µâ . Since the contorsion couples to the spin current in the microscopic field theory, one can describe the situation by considering the time-evolution generated by a time-dependent Hamiltonian We shall then evaluate the first-order perturbation induced byĤ ext (t). A usual computation easily gives the energy-momentum tensor Θ µâ (x) K as where we use the definition of the retarded Green's function (4.13) in the second line. The first term with the subscript "eq" on the right-hand side represents the equilibrium expectation value in the absence of the background contorsion. Thus, after applying the Markov approximation, which essentially means a small frequency limit, we find the linear response result for δΘ µâ (x)| (a) = Θ µâ (x)| (a) K − Θ µâ (x)| (a) eq as Note that we use the constrained limit, lim Γs ω Γ , that was discussed previously, since in the strict ω → 0 limit, both sides are expected to vanish due to relaxation dynamics of spin density. By comparing this result with Eq. (4.17), we find the rotational viscosity given by This gives one Green-Kubo formula for the rotational viscosity.
As we mentioned previously, there is another way to compute the rotational viscosity. This is because Eq. (4.17) also involves the vierbein through the torsion-free spin connection ωb ν0 (e). Thus, considering the linear response induced by the anti-symmetric part of the vierbein instead of the contorsion, we can also compute the rotational viscosity. However, we need to be careful at this point because the operator coupled to the vierbein is notΘ µâ butT μ a defined in Eq. (2.24) when we treat the vierbein and contorsion as independent sources. This second way to compute the rotational viscosity suggests that we can use the torsion-free background to compute the rotational viscosity, but the identification of Θ µâ will be a little complicated (see also Appendix B).
We can, however, derive the same result while avoiding such a direct computation with the vierbein perturbation. To do so, we simply need to use the Ward-Takahashi identity for the spin current. In flat space-time, we find the following identity in the Fourier space:  This looks more familiar because it takes a similar form to the Green-Kubo formula for the shear and bulk viscosities. We can show that the constrained limit lim Γs ω Γ is indeed necessary in these Green-Kubo formulas. To see this, using the Ward-Takahashi identity (4.23), we can replace Θ xŷ in the first Green-Kubo formula (4.22) with the spin density, which leads to another expression for the rotational viscosity: (ω, k)). Indeed, the constrained limit of Eq. (4.15) leads to (ω, k) does not give the rotational viscosity, because it simply vanishes.
While the above Green-Kubo formula utilizes the constrained limit, the strict ω → 0 limit can also provide yet another formula for the rotational viscosity. Taking the imaginary part of Eq. (4.15), it is straightforward to obtain We emphasize that lim ω→0 in this equation is the usual zero-frequency limit. Interestingly, the rotational viscosity η s appears in the denominator of the equation, in contrast to the other Green-Kubo formula that we derive in this section. This can be understood from the requirement ω Γ s in this limit, which probes infrared singularities that should be regularized by a finite Γ s . The Eq. (4.27) may provide another practical way to evaluate the rotational viscosity in quantum field theory.

Summary and Outlook
In this paper, we addressed several theoretical issues in relativistic hydrodynamics with spin polarization. First, we clarified the definition of the energy-momentum tensor and spin current for quantum field theories with Dirac fermions (such as QED and QCD) based on the first-order formalism in a background geometry with torsion. The resulting spin current is totally anti-symmetric with respect to its three indices, and we performed the consistent entropy current analysis to derive the first-order constitutive relations in the torsionful background. The total anti-symmetry of the spin current leads to three dynamical degrees of freedom arising from spin polarization, and we find one kinetic coefficient known as the rotational viscosity, which controls the relaxation rate for the spin density. Second, we presented a number of ways to evaluate the spin relaxation rate, or equivalently the rotational viscosity, from the retarded correlation functions of spin-related operators; the several Green-Kubo formulae for the rotational viscosity are derived with the help of mechanical perturbations described by the time-dependent Hamiltonian. We also point out a subtlety of the low-frequency limit in the obtained Green-Kubo formula, which originates from the non-conservation of spin angular momentum that is inherent in relativistic spin hydrodynamics.
Although we presented a rigorous framework for the first-order relativistic spin hydrodynamics, we must also point out the underlying assumptions for its validity, which have not been clearly explained in the previous literature. As was already pointed out in the old paper by Martin, Parodi, and Pershan [61], we cannot model-independently assume that the relaxational spin mode is in the near hydrodynamic regime. In other words, we cannot always find the regime for the spin hydrodynamics we discuss, which is defined by Γ s ω Γ, while we use this assumption to derive the Green-Kubo formula for the rotational viscosity. This is in sharp contrast to the strict hydrodynamic modes whose gapless property is protected by conservation laws, so that the presence of the hydrodynamic regime is always guaranteed in the k → 0 limit. There is a certain similarity between spin hydrodynamics and Müller-Israel-Stewart theory of second-order relativistic hydrodynamics [62][63][64][65] in that both introduce additional degrees of freedom. Neither of these two theories could be a universal effective description unless the additional degrees of freedom are parametrically slow. However, when an approximate spin symmetry emerges by considering, e.g., a large mass limit for fermions (e.g., charm quarks in QCD, and protons and/or electrons in QED), Γ heavy s for heavy-fermion spin is expected to be parametrically smaller than Γ for other non-hydrodynamic modes, which guarantees the existence of the spin hydrodynamic regime.
As an outlook let us point out a few research directions along the lines of the present paper. The most important one is to evaluate the rotational viscosity (or the spin relaxation rate) for QCD or gauge theories. The Green-Kubo formula for the rotational viscosity derived in this paper enables us to compute the spin relaxation rate in the underlying quantum field theory. For that purpose, similarly to the computation of the shear viscosity, we may rely on thermal perturbation theory for the weakly-coupled regime of QCD [66,67], or the AdS/CFT correspondence for the strongly-coupled limit of QCD-like theories [68,69]. The results would have direct implications for the theoretical study of spin polarization in the QGP created in relativistic heavy-ion collision experiments.
Our analysis in this paper employs the entropy-current analysis to derive the constitutive relations. Another interesting direction is to explore the non-equilibrium statistical operator method [70][71][72][73][74][75][76][77][78][79][80][81][82][83]. This will provide a complementary starting point to discuss spin transport with the help of thermal perturbation theory, which describes the relaxation process based on density operators (see, e.g., Ref. [37] for a recent attempt to derive the Green-Kubo formula with the statistical operator method). Besides, it will also be interesting to investigate possible anomalous transport phenomena in the spin transport, with the help of the path-integral formulation of the local thermodynamic functional [84]. In fact, non-dissipative anomalous currents such as the chiral magnetic and vortical effects [85][86][87][88][89][90][91] are captured by the local thermodynamic functional, or equivalently called the hydrostatic partition function [92][93][94][95][96][97][98][99][100][101]. Guided by the result obtained in Ref. [84], an extra amount of spin density, which corresponds to the deviation of the spin potential from its equilibrium value of thermal vorticity, would lead to the presence of thermal torsion emerging in the imaginary-time formalism. We speculate that this geometric description has a promising potential to provide a solid basis for systematic investigation of possible anomalous transport phenomena associated with the spin current. We leave all these interesting questions to future work.

A Global equilibrium under rotation
A global equilibrium state is characterized by a time-like Killing vector β µ (see, e.g., Refs. [102,103]). In the flat spacetime, it is a solution to the Killing equation and has the form β µ = β c (Ω µ ν x ν + b µ ) with a constant anti-symmetric tensor Ω µν = −Ω νµ , a constant time-like vector b µ , and an overall constant factor β c , which correspond, respectively, to angular velocity of rotation, the linear velocity of the center of rotation and the inverse temperature at the center. For instance, the fluid rigidly rotating around the z-axis is described by where r 2 = x 2 +y 2 . It can be checked explicitly that Eq.(A.2) is an exact solution of the ideal relativistic hydrodynamic equations. This is a consequence of the fact that conservation of angular momentum implies existence of equilibrium (nondissipative) rotating solution.
Note that thermal equilibrium in the fluid requires that the local temperature T = 1/β varies in such a way that T /γ is constant, where γ = 1/ √ 1 − Ω 2 r 2 is the relativistic time dilation factor. The solution given by Eq.(A.2) exists only for the angular velocity Ω smaller than 1/R, where R is the maximum radial size of the rotating system, i.e., r < R. This is obvious directly from Eq. (A.2), since β becomes complex otherwise, and physically, since the velocity ΩR of the outer edge of the system cannot become superluminal lest it violates causality.
We show in the following that it is physically impossible to spin the fluid fast enough so that it would violate causality, since no finite amount of angular momentum could make the fluid rotate rigidly with angular velocity Ω exceeding 1/R. The reason for that is that the temperature of the outer edge of the system approaches infinity when Ω → 1/R and, as a result, the edge carries an infinite amount of angular momentum, as we shall now show. The orbital angular momentum along the z-axis is given by J z = d 3 r x Θ 0y − y Θ 0x = d 3 r wu 0 (xu y − yu x ) = Ωβ 2 c d 3 r wβ −2 (x 2 + y 2 ) , (A.4) where w = + p is the enthalpy. For simplicity let us consider a rotating cylinder of radius R and length L and perform the integration in Eq.(A.4) in cylindrical coordinates. Since the local temperature T (r) = 1/β(r) is a monotonous function of the radial coordinate r = x 2 + y 2 , it is convenient to change the integration variables from r to T : Since thermodynamic stability requires w to be a non-decreasing function of T , the integral diverges when the edge temperature T (R) → ∞, i.e., when Ω → 1/R. Therefore, an infinite amount of angular momentum is needed to reach the causality upper bound for the angular velocity of rigid rotation Ω → 1/R.

B Currents in metric-affine connection formulation
The vierbein-spin connection formulation used in the main text is a useful tool to define the spin current and the canonical energy-momentum tensor. Since our energy-momentum tensor has a possible anti-symmetric component, one may be interested in how it appears in the formulation relying on the metric. We here provide a way to define the currents used in this paper with a help of the first-order metric-connection formulation (see also, e.g., Appendix G of Ref. [95]).
Let us now use the metric g µν and affine connection Γ µ νρ , instead of eâ µ and ωâ µ b , as independent backgrounds, which also allows a possible contribution from the non-vanishing torsion. While QCD enjoys diffeomorphism invariance as δ ξ S QCD = 0 even if we use these backgrounds, its variation with respect to the metric and affine connection leads to a different set of currents. We then introduce the following currents as an alternative to Eq. (2.11): , s µν ρ (x) ≡

(B.4)
This relation allows us to investigate dynamics of Θ µν and Σ µνρ , using t µν and s λµν that are defined by the variation with respect to the metric and affine connection.