General thermodynamic equilibrium with axial chemical potential for the free Dirac field

We calculate the constitutive equations of the stress-energy tensor and the currents of the free massless Dirac field at thermodynamic equilibrium with acceleration and rotation and a conserved axial charge by using the density operator approach. We carry out an expansion in thermal vorticity to the second order with finite axial chemical potential $\mu_A$. The obtained coefficients of the expansion are expressed as correlators of angular momenta and boost operators with the currents. We confirm previous observations that the axial chemical potential induces non-vanishing components of the stress-energy tensor at first order in thermal vorticity due to breaking of parity invariance of the density operator with $\mu_A \ne 0$. The appearance of these components might play an important role in chiral hydrodynamics.


Introduction
The physics of chiral fluids has drawn remarkable attention lately. Much of the interest is driven by the quest of the Chiral Magnetic Effect (CME) in relativistic heavy ion collisions [1]. Besides, the discovery of a finite global polarization of the Λ hyperons at the Relativistic Heavy Ion Collider (RHIC) [2] in agreement with the predictions of local thermodynamic equilibrium [3] has demonstrated that the Quark-Gluon Plasma (QGP) has sufficient local vorticity and acceleration to induce polarization of produced particles. Because of the relation between spin and chirality, this phenomenon ought to have a close connection with the physics of CME.
Thermodynamic equilibrium coefficients include -amongst others -the well known proportionality coefficient between vector current and magnetic field in the CME, originally derived in ref. [17] assuming equilibrium. In fact, the systematic derivations of all constitutive relations at thermodynamic equilibrium, which is non-dissipative by construction, can be carried out with a functional approach [18,19], with finite temperature field thoery [20] or with the covariant density operator approach [21]. The latter was used to derive the constitutive relations of the stress-energy tensor in presence of acceleration and vorticity.
In this paper, we will take this approach to study the thermodynamic equilibrium of a chiral fluid endowed with an axial chemical potential (i.e. a non-vanishing mean axial charge) including acceleration and vorticity and we will present explicit calculation for the free massless Dirac field. We will see that the presence of a non-vanishing axial charge is responsible for additional parity-breaking terms in the constitutive equations of the stressenergy tensor and currents, which would not be otherwise present. These terms can be very important for the physics of CME and particularly for the the implementations of chiral magnetohydrodynamics. Indeed, general thermodynamic equilibrium relations are the fundamental building blocks of hydrodynamics once they are promoted to space-time dependent functions. Particularly, we will see that the axial chemical potential entails the appearance of new terms which are first-order in the gradient expansion.

Notation
In this paper, we use the natural units, with = c = k B = 1. The Minkowskian metric tensor is diag(1, −1, −1, −1); for the Levi-Civita symbol we use the convention ǫ 0123 = 1. Operators in Hilbert space will be denoted by a large upper hat, e.g. T (with the exception of Dirac field operator that is denoted by Ψ). The stress-energy tensor used to define Poincaré generators is always assumed to be symmetric with an associated vanishing spin tensor.

Thermodynamic equilibrium and chirality
In this section we will summarize the density operator approach and the inclusion of a conserved axial charge.
The covariant density operator formalism has been described in detail elsewhere, here we will briefly summarize it. The most general equilibrium density operator can be obtained by maximizing the total entropy S = −tr( ρ log ρ ) with the constraints of given energy, momentum and charge densities over some spacelike hypersurface Σ [21][22][23] and can be written as [23][24][25], see also [26,27], where T µν and j µ i are the symmetric stress-energy tensor and some conserved currents operators, β µ is the four-temperature vector such that T = 1/ β 2 is the proper comoving temperature and ζ i are the ratio of comoving chemical potentials and the temperature ζ i = µ i /T . The Eq. (2.1) describes a local thermal equilibrium; to obtain a global equilibrium distribution the operator (2.1) is required to be time independent, implying that the integrand is divergenceless. This is possible if [25] all ζ i are constants and the fourtemperature β is a Killing vector field: The general solution of Eq. (2.2) in Minkowsky space-time is where b is a constant four-vector and ̟ a constant antisymmetric tensor, dubbed thermal vorticity, that accounting for Eq. (2.2) can be expressed as an exterior derivative of the β field Recalling the definition of four-momentum operator P , conserved charges Q i The state characterized by this density operator is a generalization of global equilibrium, and since the parameters are globally defined we will refer to this condition as general global equilibrium.
As conserved currents, one usually have electric and baryon currents. However, for massless fermions there is one additional conserved current, the axial current j A , such that: which has an associated conserved axial charge where Σ is an arbitrary space-like hypersurface. Thus, for a system of massless fermions, one would deal with a general global thermodynamic equilibrium of this sort: In fact, the axial current is not conserved in presence of interactions because of the well known axial anomaly and the above continuity equation gets modified into (for QED): Nevertheless, in presence of an external field such that the right hand side of (2.9) vanishes, the axial charge is still constant and thermodynamic equilibrium is possible. If, on the other hand, we deal with dynamical gauge fields, the right hand side can be written as the divergence of the so-called Chern-Simons current: and its integral (the helicity, for the EM field) can be embodied in the definition of a gauge-invariant overall conserved axial charge including both fermion and gauge boson contribution [29]. In both the above cases, the density operator (2.8) is perfectly meaningful.
In this work, we will confine ourselves to the simplest case of free fermions, so that (2.6) is fulfilled. Thus, the obtained relations will not include contributions from possible external gauge fields, yet they will be very useful as a leading order expressions for a chiral fluid.
As a final remark, it is important to stress that β must be a time-like vector for the local mean values to be physical (see next Section). In general, for ̟ = 0 there is a 3D boundary in space-time separating a region where β is space-like; for an analysis of the effect of the boundary inside this region on rotating fermions see [28]. In the physical region, the four-temperature vector β defines a four-velocity u of the fluid -the β frame [21] or thermodynamic frame [18] - (2.10) As we will see, this frame is the most appropriate for thermodynamic equilibrium.

Mean value of local operators
Our goal is to calculate the mean values of physical quantities at thermodynamic equilibrium with the density operator ρ in Eq. (2.8): For this purpose, we will use the method already employed in refs. [26,30]. The first step is to shift the angular momentum operator at x by using the translational operator T(x) = exp[ ix · P ]: Thanks to this definition and Eq. (2.3) we can rewrite the density operator (2.8) as The second step is to make an expansion in thermal vorticity. Indeed, thermal vorticity (2.4) is adimensional in natural units and, under most circumstances (e.g. in relativistic heavy ion collisions) is ≪ 1. An expansion of (3.2) in ̟ will result in a good approximation of the mean values. At the leading order one has: where in the last step we have taken advantage of the translational invariance of the homogeneous operator ∝ exp[−β · P ]. The above expression shows that -as expected -for small deviations from homogeneous equilibrium, the mean value is well approximated by the homogeneous equilibrium value with four-temperature equal to its value in the same point x where the operator is evaluated. The corrections to the leading order term are obtained by using the expansion based on a path ordered method developed in ref. [23], in our case along the fluid velocity u (2.10). This method is also appropriate for Euclidean space-time. The resulting expression is given in terms of connected correlators: weighted with the homogeneous statistical operator At the second order in thermal vorticity we obtain [30] O where T τ (· · · ) is the T τ -ordered product defined as with p the permutation that orders τ by value: Thermal vorticity can be decomposed into two space-like fields α and w [26] much the same way as the electro-magnetic field tensor, by using the fluid four-velocity u: At global equilibrium, when β is a Killing vector, it can be shown that w and α are the vorticity and acceleration vectors divided by the local temperature T [26]: Furthermore, it is useful to define another four-vector: where ∆ µν = g µν − u µ u ν is the projector onto the hyperplane perpendicular to u. The four-vector γ is second-order in thermal vorticity and is, by construction, orthogonal to u, w and α. Similarly, the generators of the Lorentz group J can be decomposed as: The operators J and K are the generators of rotations and boosts with respect to the local frame defined by u. Making use of (3.5), we obtain whence the mean value expansion (3.4) becomes or, after defining the correlators the expansion in thermal vorticity of the mean value of a local operator O(x) can be written as: (3.9) Note that the above expansion applies to any local operator O, whether it is a scalar or a component of a tensor of any rank.

Chiral thermodynamic coefficients
In this section, we identify and determine the thermodynamic coefficients up to second order in thermal vorticity for the stress energy tensor as well as vector and axial currents. We have seen in the previous section, that the mean values of local observables can be expanded in thermal vorticity (3.9) and that each term involves a coefficient expressed as connected correlator between the observables and the Lorentz group generators (3.8), caluculated with the homogeneous statistical operator ρ h (3.3) This operator features several symmetries: translations, rotations in the hyperplane perpendicular to β(x), and time-reversal T in the β(x) rest frame. However, as Q is odd under charge conjugation C and Q A is odd under parity P (see Table 1), ρ h is not invariant under C and P.
Similarly, rotational invariance implies that correlators in (3.8) will be non-vanishing only if the involved operators transform like the same irreducible vectors of the same representations of the rotation group. To study the impact of symmetries on the correlators (3.8) it is convenient to redefine: and to decompose both A and O in irreducible components under rotation: scalar, vector, symmetric traceless tensor, etc. The only non-vanishing correlators will be those matching corresponding components of the same kind and with the same transformation property under time-reversal.
In the homogeneous equilibrium density operator (4.1), the term µ A Q A is the only responsible for parity breaking, if the Hamiltionian is parity invariant. For this reason, all correlators (3.8) in which A and O transform with different sign under parity will be proportional to odd powers of µ A , because: These correlators will be henceforth dubbed as chiral correlators, and the corresponding coefficients chiral coefficients (see Table 1). Table 1. Transformation properties under discrete transformation for various operators in the rest frame. P is parity, T is time-reversal and C is charge conjugation transformation. T µν is the stress-energy tensor, j µ V is an (electric) vector current, j µ A is an axial current, φ is a scalar and φ A is a pseudoscalar, K µ and J µ are the Lorentz generators, Q and Q A are vector and axial charge respectively.
We refer to [26,30] for a detailed description of the method of calculating the coefficients; herein we briefly summarize it. The coefficients can be expressed as thermal expectation values in the rest frame where β µ = (1/T (x), 0) and they are function of T , ζ = µ/T and ζ A = µ A /T ; we denote the mean values in rest frame with a subscript T , that is: where H is the Hamiltonian. Consequently, the correlators (3.8) in the rest frame are written as: All the thermodynamic coefficients are linear combination of correlators (4.2) with an appropriate K, J and O components choice, and such expressions can be considered as "Kubo formulae" for the corresponding coefficients. From the Lorentz generators definition (2.5), it follows that the (4.2) are given by means of correlators of the symmetric stress-energy tensor and the operator O. In this work, we work out the (4.2) by using the imaginary time formalism. The shifted boost and angular momentum generators, here obtained by definition (3.7) in the rest frame, are given by Eq. (3.1) and, according to definition (2.5), arise from the spatial integration of the time evolved stress-energy tensor operator Thus, the correlators (4.2) are expressed as a linear combination of the following basic structure with suitable indices: where X n = (τ n , x n ).
We are now in a position to determine the expansion in thermal vorticity of conserved currents including the axial chemical potential µ A = ζ A T . But before that, it should be pointed out that any scalar or a pseudo-scalar operator, denoted as φ and φ A in Table 1, does not have chiral corrections up to second order. Indeed, the correlator between φ or φ A and J or K simply vanish at first-order because the latter are both vectors under rotation. Besides, it can be readily seen again from Table 1, that no second-order correlator can be formed with K and J that transform at the same time even under time reversal and odd under parity. Hence, there are no chiral-vortical corrections even at second order for scalar and pseudoscalar operators.

Stress-energy tensor
The expansion up to second order in thermal vorticity of the stress-energy tensor with axial chemical potential features additional terms with respect to the case ζ A = 0 [26,30]. Indeed, three new chiral coefficients appear at the first order in thermal vorticity, A, W 1 and W 2 : which can be obtained as: It is useful to express the previous coefficients in terms of the quantity, see Eq. (4.4), which will be computed in Section 6.1 for the massless Dirac field. By using the definition (2.5), we obtain The conservation equation for the mean value of the stress-energy tensor: implies that not all the coefficients in (4.5) are independent. Indeed, enforcing the continuity equation, the following constraints for second order coefficients are obtained [30]: where all derivatives are to be taken with fixed ζ = µ|β| and ζ A = µ A |β|. For the new first-order chiral coefficients we have: meaning that only W 1 and W 2 are really independent. In this work, however, we have computed all the chiral coefficients for a massless Dirac field using eqs. (4.6) and used the condition (4.8) as a consistency check.

Currents
Consider now the conserved current j µ V related to the charge Q (2.7) and its transformation properties in Table 1. The only non vanishing terms up to the second order in thermal vorticity expansion turn out to be: and the chiral coefficient W V reads: which is the conductivity of the so-called Chiral Vortical Effect (CVE). We can write W V by means of the general correlator: Taking the divergence of (4.9), one realizes that the mean vector current is conserved For a massless spinor field, using dimensional analysis, W V must be proportional to |β| −3 ensuring that the previous relation is always verified. Similarly, taking into account its transformation properties, the conserved axial current j µ A has the following expansion: (4.14) The so-called Axial Vortical Effect (AVE) shows up in the coefficient W A and has been discussed in detail in [30], see also [31] for complete non chiral corrections derived from a Wigner function ansatz. This coefficient does not vanish for µ A = 0 and it is thus independent of the anomaly (see discussion in Section 7). On the other hand, the proper chiral terms are given by the coefficients: and, in terms of the auxiliary correlators they can be written as: (4.16) Again, conservation equation of the mean axial current entails a condition for W A which is equivalent to that of the vector current and discussed in [30].

Discussion: β vs Landau frame in chiral hydrodynamics
As an interlude preceding the explicit calculation of the thermodynamic coefficients, it is worth discussing the question of the frame, which is relevant for the hydrodynamics of chiral fluids, that can be defined in general as fluids where the axial current is relevant. As has been mentioned, the four-velocity in our formulae has been defined starting from the four-temperature vector β which is a Killing vector in global equilibrium. This frame is called β or thermodynamic frame [18,21] and, as it can be seen from Eq. (4.5) it does not coincide with Landau frame, where u is by definition an eigenvector of the (symmetrized) stress-energy tensor. However, unlike in the non-chiral case, where the difference between Landau frame velocity and β frame velocity shows up only at second order in the gradients of u or β owing to the terms in u ν γ µ [26], see also [32] for an in-depth analysis; here the difference between the two velocities shows up at first order in the gradients because of the chiral terms in u ν w µ in (4.5) (see Appendix A for the transformation equations). This is especially important for chiral hydrodynamics because it implies that one cannot formulate it at first order without taking these equilibrium non-dissipative terms in u ν w µ into account. Even if one tried to eliminate them in the stress-energy tensor by going to the Landau frame, they would reappear in the constitutive equations of the axial current at first order in the gradients; this is just what happens, see Appendix A. For chiral-magneto-hydrodynamics, with dynamical electro-magnetic fields, the transition to Landau frame would be even more complicated because a change of the velocity definition implies a change of definition in the comoving electric field.
In conclusion, the β frame appears as the most appropriate to write down hydrodynamic equations, as advocated in [18,21].

Free massless Dirac field
We are now going to calculate the coefficients for the massless free Dirac field at finite temperature, vector and axial chemical potential. The basic quantity to compute is the grand-canonical partition function Z: For the non interacting massless Dirac field described by the Lagrangian the standard vector and axial currents read: with the corresponding conserved vector and axial charges (2.7) As it is known, the partition function can be written as a path integral of fields in Euclidean spacetime with antiperiodic boundary condition (ABC) Ψ(|β|, x) = −Ψ(0, x) and Euclidean gamma matricesγ µ , fulfilling the relation {γ µ ,γ ν } = 2δ µν . It is useful to define the following notation where p n = 2π(n + 1/2)/|β| with n ∈ Z are the Matsubara frequencies, and to introduce the right and left chiral chemical potentials and the corresponding projection operator where χ =R,L is assigned to be respectively +1, −1. Then, the Dirac propagator in the imaginary time reads [33]: where latin characters a, b, . . . denote spinorial indices, / P =γ µ P µ is the standard contraction between the (Euclidean) Dirac matricesγ µ and the (Euclidean) four-momenta P . The Euclidean canonical stress-energy tensor, see Eq. (4.3), equals where the i δ 0µ factor stems from Wick rotation. The Belinfante-symmetrized stress-energy tensor used to construct all Poincaré generators, is simply the symmetrization of the previous one: which can be expressed according to the point-splitting procedure as: where:

First order correlators
Here we outline the procedure to evaluate first and second order coefficients of a general operator O α...β (x) for free massless Dirac field theory. We can write a generic correlator (4.4) related to first order corrections in vorticity as where T µν is the symmetric stress-energy tensor coming from Lorentz generators (2.5) and the two operators connected correlator is In the same way as stress-energy tensor in (6.3), we use the point-splitting procedure to express the operator O as Thus, to compute (6.4), we first consider Now the connected correlator only concerns Dirac fields and can be computed using Wick's theorem With this result and exchanging the order of the two anti-commuting fieldsΨ(X 1 ) a Ψ(Y 2 ) d we are able to recreate a trace operation: Now, we express fermionic propagator in Fourier space as in (6.2) where, as previously, P + R/L = (p n + iµ R/L , p) and Q + R/L = (q m + iµ R/L , q). The derivatives, appearing in the point-splitted operators D and O, act on the propagators and are readily obtained: At this point, we take these derivatives, we send X 1 , X 2 → X and Y 1 , Y 2 → Y and we rename Q into −Q; note that Q + χ ′ = (q m + iµ χ ′ , q) goes to −Q − χ ′ = (iµ χ ′ − q m , −q). At the end, the correlator C(X, Y ) becomes where we defined The trace is to be carried out over spinorial indices by using the Euclideanγ matrices properties: tr (γ µγν ) = 4 δ µν , tr γ µ 1 . . .γ µ 2n+1 = 0, tr (γ kγλγµγν γ 5 ) = 4ǫ kλµν , tr (γ kγλγµγν ) = 4 δ kλ δ µν − 4 δ kµ δ λν + 4 δ kν δ λµ .
However, all the operators O that will be considered are analytic functions in the fourmomentum Q and P ; as a result, the function F χ will also be analytic in those variables. This enables us to analytically entail the sum over the Matsubara frequencies.
Before carrying out the thermal sum, we send Y to 0, as dictated by (6.4), and separate the momentum integral from the frequency sum in (6.1) As stated previously, F χ is an analytic function, hence the sum over frequencies can be carried out using the formula: where ω could be the modulus of p 2 or of q 2 and n F is the Fermi-Dirac distribution function: Introducing the notationP and similarly forQ(±), and taking advantage of Eq. (6.5), we have Note that we have added and subtracted the chemical potential in the exponential in such a way that the corresponding energy couples with the relating argument of the Fermi distribution function n F . This feature will be used later on. Finally, putting together and integrating by parts, we can write the coefficient (6.4) as: that is an integral over time and momenta of a known function. We can simplify this expression further with the following steps: • after the evaluation of F χ (P, Q) calculate S χ ; • then derive S χ (p, q, τ ) and replace q → −p; • integrate over solid angle Ω: d 3 p = p 2 dp dΩ first and over τ thereafter; • express all exponential factors with the Fermi distribution function using • express the derivative of n F (E) in terms of its powers • after the due simplifications, express powers of n F (E) in terms of its derivative: in this way a simple expression made only of n F (E) and its derivative is obtained; • finally, integrate by part n • sum up all the pieces with different χ.
This procedure always leads to expressions of the form where P N is a polynomial of N degree, η is +1 if the correlator C is even under charge transformation, otherwise it is −1, likewise η A = ± reflect parity transformation, and K is a numerical factor. The correlators in the form (6.7) can be easily integrated once we know the polynomial P N . The coefficients evaluated in [30] have either η = +, η A = + or η = −, η A = + and µ A = 0, that means µ R = µ L = µ and thus have the form So we can immediately obtain the dependence on chiral chemical potential for non chiral coefficients putting them in the form (6.7) and setting K =K/2, whereK is the numerical factor evaluated in [30]. The results obtained using this procedure are written in Section 6.3.

Second order correlators
A second order correlator of the perturbative expansion (4.4) is the mean value of the local operator O α...β with two Lorentz generators. We will follow similar steps as those used for the first order correlators. We start by writing the second order generic correlator (4.4) in Euclidean space-time as where X = (τ 1 , x), Y = (τ 2 , y) and the connected correlator is given by Again, using the point splitting procedure and taking advantage of the Wick's theorem, we can split the mean value of (6.8) in two parts where Z = (τ 2 , z) and Following similar procedure of first order correlators, we express the two point function in Fourier space (6.2) and we can show that C 1 and C 2 are equal to where we denoted the momenta as P + χ = (p n + iµ χ , p), Q + χ = (q m + iµ χ , q), K + χ = (k l + iµ χ , k), and we defined the trace functions as: Once the form of O is explicitly given, these two trace functions are computed thanks to the well known gamma matrix traces. In general, we can say that as far as O is an analytic function of the derivative operator ∂ X , the two trace functions F 1 and F 2 are analytic functions on the four-momenta P , Q and K. Now, we send Z to 0 in Eq. (6.9), as required by (6.8), and we separate the momenta integrals from the frequencies sums 9 e i(p+k)·x e i(q−k)·y S χ (p, q, k, τ 1 , τ 2 ), (6.10) where we have defined Since τ 1 and τ 2 are always larger than zero in the integration (6.8), the sums over p n and q m are made using (6.5). For the k l sum, the result is still (6.5) when τ = τ 1 − τ 2 > 0, instead when τ = τ 1 − τ 2 < 0 the sum yields Reminding the definition (6.6) forP (±),Q(±),K(±) and for p ± χ , q ± χ and k ± χ , after all the sums, we find S χ (p, q, k, τ 1 , τ 2 ) = 1 8 p q k s 1 ,s 2 ,s 3 =± F 1 χ P (s 1 ),Q(s 2 ),K(s 3 ) S 1,χ , where θ is the Heaviside theta function. Now, we can take advantage of the formula to integrate over the coordinates x and y in the Eq. (6.8) and, using Eq. (6.10), we obtain: From this point, we can adopt the same algorithm listed in the previous section, with an additional time integration, to obtain simple expressions of the form (6.7).

Results
First, we consider both the canonical and symmetric stress energy tensor. In the Euclidean space-time they are given by and they give different results for A, W 1 , W 2 , G 1 , G 2 , i.e. the coefficients related to the non symmetric part of the expansion The non chiral coefficients can be ensued from [30] as described at the end of Section 6.1 and are All non-chiral coefficients reported in this section are also listed in Table 2. Instead, the chiral coefficients are given by Eq.s (4.6), while accounting for (4.7) can be calculated using the methods described above; they are found to be We can readily check that the coefficients (6.11) satisfy the relation (4.8). All chiral coefficients reported in this section are also listed in Table 3. The Euclidean version of the Dirac field electric vector current j µ V =Ψγ µ Ψ is The non chiral coefficients of the decomposition Table 2. The non-chiral thermodynamic coefficients up to second order in thermal vorticity of the stress-energy tensor, electric and axial current for a free massless Dirac field (see Eqs. (4.5),(4.9),(4.14) for definitions). Here we use T = 1/|β|, µ = ζ T and µ A = ζ A T . Table 3. The chiral thermodynamic coefficients up to second order in thermal vorticity of the stressenergy tensor, electric and axial current for a free massless Dirac field (see Eqs. (4.5),(4.9),(4.14) for definitions). Here we use T = 1/|β|, µ = ζ T and µ A = ζ A T .
are obtained using the Eq.s (4.10) while for the CVE conductivity W V (4.11) evaluating the correlators in Eq. (4.12) we obtain the well-known result The same is done for the axial current which has similar expansion on thermal vorticity The only non chiral coefficient is the axial vortical effect conductivity W A , for which we have [30] The chiral coefficients are again evaluated accounting for the relations (4.15) and (4.16) As expected for dimensional analysis, the coefficients W V (6.12) and W A (6.13) respectively fulfill the relations (4.13) and (4.17).

Discussion and summary
Many of the coefficients presented here have been obtained with different methods. The expectation values of the stress-energy tensor and the currents for free massless Dirac field evaluated here are in agreement with the exact results in the case of pure rotation and vanishing chemical potentials of [20]. The CVE and AVE have been discussed at length over the past decade (see [1,34,35]), and several calculations of W V (6.12) and W A (6.13) were presented in literature.
As it was mentioned in Section 4.2, the AVE does not vanish for ζ A = µ A /T = 0, as it is apparent from eq. (6.13). This means that the AVE does not need chiral imbalance to show up and persists for a perfectly chirally balanced system. We interpret this feature -at least in the form in eq. (6.13) -as an evidence that AVE is not related to the anomalous axial current divergence. As it was found by A. Vilenkin [36][37][38] and lately highlighted in ref. [30], an axial current proportional to vorticity arises because it is allowed by the symmetry of the density operator (this was also recognized in [39]) with rotation and, indeed, the very same expression (6.13) is found for free massless fermions, without coupling to either external or dynamical gauge fields which make up the anomalous divergence of the axial current. We understand that this statement goes against the conventional interpretation in literature, and yet some other authors have lately cast doubts about the anomalous origin of the AVE [40,41]. We hope that our viewpoint will contribute to clarify the point.
As far as the stress-energy tensor is concerned, the consequence of parity breaking induced by ζ A , is an energy flux q µ along the vorticity w µ , see Eq. (4.5): To our knowledge, this kind of coefficient was obtained first by [36] as an energy flux resulting from neutrinos emitted by a rotating black hole; indeed, our result of W Can 2 for a free Dirac field (6.11) perfectly agrees with the one in [36](there is a factor 1 2 of a difference because Vilenkin considered only left-handed neutrinos). Lately, the same result for this coefficient was obtained by [7,8] using holographic techniques, by [11][12][13] in chiral kinetic theory and in [8,9] evaluating Kubo formulae in finite temperature field theory.
It is also worth noting that thermodynamic equilibria with vorticity imply different mean values for the canonical and symmetric stress-energy tensor [42,43]. This is seen here for the coefficients A, W 1 , W 2 , G 1 and G 2 . Particularly, the coefficient A vanishes if the stress-energy tensor is symmetric but not for the canonical.
The second order coefficients N A α , N A w and G A , see eqs (4.14) and (4.15) are indeed newly obtained and they appear as quantum corrections of the axial current in presence of acceleration and rotation. They are the axial counterpart of the second-order equilibrium corrections of the vector current discussed in [30]. N A α , N A w are corrections along the fluid velocity and hence modify the axial charge density, while G A yields a chiral flow along the four-vector γ µ defined in Eq. (3.6).
To summarize, we have studied global thermodynamic equilibrium with acceleration and vorticity with axial charge chemical potential inducing a macroscopic parity breaking. We carried out an expansion of the stress-energy tensor at the second order in the thermal vorticity tensor, including acceleration and vorticity and obtained the constitutive equations at equilibrium of the stress-energy tensor as well as of the vector and axial currents. These equations may be phenomenologically relevant for the physics of the Quark Gluon Plasma and especially for its hydrodynamic modelling in presence of chirality imbalance. eigenvector of the stress energy tensor: where ρ L is both the eigenvalue and the energy density in the Landau frame. Being a four-vector, the Landau velocity u L can be expressed in terms of the tetrad {u, α, w, γ} defined in Eq.s (3.5) and (3.6): where |w| = √ −w 2 , |γ| = −γ 2 , |α| = √ −α 2 and a, b, c and d are four unknown constants such that u µ L u Lµ = 1, i.e.
Furthermore, since for a vanishing thermal vorticity we expect that the two thermodynamic frames coincide, for ̟ = 0 we must have that a = 1 and ρ L = ρ. This means that the leading term of a and ρ L is zeroth order in thermal vorticity. Now, taking the stress-energy tensor mean value at the second order expansion on thermal vorticity and contracting it with the landau velocity (A.1), we obtain where for simplicity we defined ρ eff ≡ ρ − α 2 U α − w 2 U w and p eff ≡ p − α 2 D α − w 2 D w . Equating this expression to ρ L u ν L , with u L given by (A.1), we obtain the five equations in five unknown variables (a, b, c, d, ρ L ) that diagonalize the stress energy tensor (at second order in thermal vorticity): To solve it, we can simplify the equations furtherer. Indeed, we can write the first Eq. of (A.2) as b a = ρ eff − ρ L W 1 In the same way, we can write the third Eq. of (A.2) isolating b/a: hence, equating the previous two, we obtain We can also note that in the second equation of (A.2) the term in d is higher order and can be omitted, as a consequence, after the replacement of Eq.s (A.6) and (A.7), the second equation of (A.2) becomes a third grade equation for ρ L : W 1 W 2 (ρ L +p eff )|w| 2 +(ρ L +p eff )(ρ L −ρ eff )(ρ L +p eff +W |w| 2 )+W 1 G 1 A|w| 2 |α| 2 +A(ρ eff −ρ L )|α| 2 = 0.
We find that the only solution of this Eq. reproducing ρ L = ρ for ̟ = 0 is: Replacing this solution in (A.6), we obtain b a = W 2 ρ + p |w| + O(̟ 3 ) and replacing it in (A.7) gives From these two, we learn that b is actually first order on thermal vorticity and c is second order, thus from (A.5) it follows that d is in fact third order and can be set to zero. Consequently, the last equation to solve is a 2 − b 2 − c 2 = 1, that gives In conclusion, we found that the relation between the Landau fluid velocity u L and the β-frame fluid velocity u is or, reverting it, u is given by This transformation at first order in vorticity and the following relation between CVE conductivity in Landau frame with W 2 and W V were also pointed out in [8]. At this point, taking advantage of (A.8), we can express the stress-energy tensor mean value at second order in thermal vorticity (4.5) in the Landau frame: T µν = Aǫ µνκλ α κ u Lλ − AW 2 ρ + p ǫ µνκλ α κ w λ + (W 1 − W 2 )u ν L w µ − p eff ∆ µν L + ρ eff − W 2 2 ρ + p |w| 2 u µ L u ν L + Aα µ α ν + W (ρ + p) + W 1 W 2 ρ + p w µ w ν + AW 2 ρ + p u µ L γ ν + (G 2 − G 1 )(ρ + p) + AW 2 ρ + p u ν L γ µ .
Furthermore, using (A.8), we can also write the vectorial current (4.9) in the Landau frame Instead, replacing (A.8) in (4.14), we obtain the axial current in the Landau frame