Kubo formulas for thermodynamic transport coefficients

Uncharged relativistic fluids in 3+1 dimensions have three independent thermodynamic transport coefficients at second order in the derivative expansion. Fluids with a single global $U(1)$ current have nine, out of which seven are parity preserving. We derive the Kubo formulas for all nine thermodynamic transport coefficients in terms of equilibrium correlation functions of the energy-momentum tensor and the current. All parity-preserving coefficients can be expressed in terms of two-point functions in flat space without external sources, while the parity-violating coefficients require three-point functions. We use the Kubo formulas to compute the thermodynamic coefficients in several examples of free field theories.


Introduction
When a macroscopic system in equilibrium is subject to external fields, it reacts to the fields by adjusting its pressure, energy density, and other thermodynamic functions. If the system is subject to external electric and magnetic fields, the response is described by the electric and magnetic susceptibilities, which determine the electric permittivity and the magnetic permeability of matter. If the system is subject to external gravitational fields, the analogous gravitational susceptibilities determine the appropriate response of the free energy.
When the same system is perturbed out of equilibrium, the equilibrium susceptibilities can contribute to non-equilibrium phenomena. For example, the electric and magnetic susceptibilities determine the speed of light in matter, while the pressure determines the speed of sound in matter, through the equation of state.
Motivated by the applications of relativistic hydrodynamics to the hot sub-nuclear matter [1,2], we will be focusing on relativistic fluids in this paper. The most basic thermodynamic susceptibility is of course the pressure itself, which can be viewed as a response of the free energy to a diagonal metric perturbation thanks to the covariant definition of the energy-momentum tensor in relativistic systems.
For relativistic fluids in 3+1 dimensions in curved space, the transport coefficients that are thermodynamic in nature were first noticed in [3,4] at second order in the derivative expansion, though their thermodynamic nature was not fully appreciated at the time. For fluids in 2+1 dimensions, analogous thermodynamic transport coefficients already appear at first order in the derivative expansion [5]. They have been variously referred to in the literature as "thermodynamic response parameters", "thermodynamic transport coefficients", "thermodynamical hydrodynamic coefficients", "equilibrium hydrodynamic coefficients", or "non-dissipative transport coefficients". 1 There is a multitude of notations for these coefficients in the literature, and the translation between different conventions is not always straightforward.
While the thermodynamic transport coefficients were first noticed in the context of hydrodynamics [3][4][5], their connection with thermodynamics was not explored until [6,7]. These papers showed that the relevant coefficients in the constitutive relations follow from the equilibrium partition function, including the highly non-trivial constraints [5,8] demanded by the local positivity of entropy production. We will refer to the thermodynamic coefficients that appear in the constitutive relations as "thermodynamic transport coefficients", and to the coefficients in the equilibrium free energy as "thermodynamic susceptibilities". Thermodynamic transport coefficients are linear combinations of thermodynamic susceptibilities and their derivatives [6,7]. In the classification of non-dissipative transport coefficients in [9], thermodynamic transport coefficients correspond to class L.
We will be considering fluids with a conserved global U (1) charge, such as the baryon number. We will refer to the fluids that can be locally described as having a temperature and a chemical potential for the global U (1) charge as "charged fluids". The system can be coupled to the corresponding non-dynamical external U (1) gauge field, and to the nondynamical external metric. The thermodynamic susceptibilities then include the usual "electric" and "magnetic" susceptibilities, as well as the response of the free energy to the vorticity, to the Riemann curvature, the magneto-vortical response, etc. In 3+1 dimensions, there are nine such susceptibilities at two-derivative order [6]. These susceptibilities will appear in the constitutive relations and in equilibrium correlation functions for fluids in flat space without external U (1) fields.
As we will see later, the second-order thermodynamic transport coefficients in QCD at non-zero baryon number chemical potential (in flat space without external fields) are determined by five thermodynamic susceptibilities. For a parity-preserving conformal theory, the second-order thermodynamic transport coefficients at non-zero chemical potential (in flat space without external fields) are determined by three thermodynamic susceptibilities.
Kubo formulas for second-order thermodynamic transport coefficients were derived in [10][11][12] for uncharged fluids. Further, [12,13] evaluated these coefficients for noninteracting massless scalars, fermions, and gauge fields. Ref. [14] evaluated the seven parity-even thermodynamic susceptibilities of a charged fluid in a theory of free massless fermions, using a dimensionally reduced partition function in curved space. In addition to the above functional methods, [15] evaluated the thermodynamic second-order coefficients for charged fluids of free scalars and fermions using operator methods.
Our focus in this paper will be on the Kubo formulas for all nine second-order susceptibilities in 3+1 dimensions. We will write down the Kubo formulas for the susceptibilities, rather than for the thermodynamic transport coefficients. This is natural, as the thermodynamic transport coefficients are derived objects, while the susceptibilities are fundamental. We will find that the Kubo formulas for all seven parity-preserving susceptibilities can be written in terms of two-point correlation functions of the energy-momentum tensor and the U (1) current. In other words, using three-point functions to evaluate the thermodynamic transport coefficients as in [12,15] is not necessary, and two-point functions are sufficient. 2 Using two-point functions will hopefully allow for an easier evaluation of these transport coefficients on the lattice [12,17], and in holography [3,4,18]. Further, we write the free energy in a covariant form as in [7], which directly gives covariant expressions for the energy-momentum tensor and the U (1) current, generalizing the results of [7] for secondorder transport coefficients to charged fluids. (See [19] for the generalization of the results of [6] to charged fluids). We illustrate our Kubo formulas by evaluating the susceptibilities in a few examples of free field theories.
The paper is organized as follows. In section 2.1 we introduce the equilibrium generating functional for thermodynamic correlation functions, and define the energy-momentum tensor as well as the conserved current in terms of its variation under external sources. In section 2.2 we write the generating functional in terms of the nine independent susceptibilities appearing at second order in the derivative expansion. Section 2.3 briefly talks about the trace anomaly. Section 2.4 then provides the expressions for the energy-momentum tensor and the conserved current in terms of the second order susceptibilities for a charged fluid in the absence of external fields. We present the Kubo formulas in section 2.5, and express the susceptibilities in terms of equilibrium two-and three-point functions. In sections 3.1, 3.2 and 3.3 we evaluate the susceptibilities for free scalar, free Dirac fermion, and free gauge fields. We end with a discussion in section 4. Appendix A provides relations between our thermodynamic susceptibilities and the ones that have appeared previously in literature.

Equilibrium generating functional
We consider a macroscopic system that has degrees of freedom which couple to the external metric g µν and to an external Abelian gauge field A µ . The matter in equilibrium is described by the generating functional for equilibrium (zero-frequency) correlation functions W [g, A]. Such generating functionals have been discussed starting from [6,7], with applications to relativistic hydrodynamics. We follow the presentation of [7,20,21]. Equilibrium is characterized by a timelike Killing vector which we denote by V . The coordinates in which V µ = (1, 0) correspond to the matter at rest. 3 The matter velocity, temperature, and the chemical potential are defined as where β 0 is a constant which sets the normalization of temperature, and Λ V is a gauge function which ensures that the chemical potential is gauge invariant. Denoting the Lie derivative by £, the conditions for being in equilibrium are For systems with a finite correlation length, the equilibrium generating functional is extensive in the thermodynamic limit, and can be written as where the density F[g, A] is a local function of the external sources g µν and A µ . We will assume that the microscopic theory has no chiral anomalies, and so W [g, A] is gauge-and diffeomorphism-invariant. We define the energy-momentum tensor T µν and the current J µ in the standard fashion, The diffeomorphism-and gauge-invariance of W [g, A] imply, respectively, where the U (1) gauge field strength is F µν = ∂ µ A ν − ∂ ν A µ . When the external sources g and A vary in space on length scales much longer than the correlation length, the density F[g, A] can be written as a derivative expansion of the external sources. The problem of finding the generating functional then boils down to finding the gauge-and diffeomorphisminvariants made out of the metric, the gauge field, and the quantities in (2.1), up to a given order in derivatives. Before we start writing down the invariants that appear in (2.3), it is worth emphasizing the identities which follow from the fact that the system is in equilibrium. The equilibrium conditions (2.2) imply that the fluid velocity, temperature, and the chemical potential defined by (2.1) are not arbitrary, but rather must obey Here the acceleration is a µ ≡ u λ ∇ λ u µ , the projector ∆ µν ≡ g µν + u µ u ν projects onto the space orthogonal to u µ , the shear tensor is σ µν ≡ ∆ µα ∆ νβ (∇ α u β + ∇ β u α − 2 3 ∆ αβ ∇·u), and the electric field is E µ ≡ F µν u ν . The first equation in (2.6) says that T and µ are time-independent in equilibrium. The second equation in (2.6) says that the gravitational potential induces a temperature gradient. This is a consequence of Tolman's law [22] (equilibrium temperature is proportional to 1/ √ −g 00 in the appropriate coordinates). The third equation in (2.6) says that the electric field induces a charge gradient. This is a formal way to express the phenomenon of electric screening. Alternatively, if (2.6c) were not true, there would be entropy production due to the electrical conductivity. The last equation in (2.6) says that the expansion and shear must vanish in equilibrium. If it were not so, there would be entropy production due to the bulk and shear viscosities. We will find it convenient in what follows to use the electromagnetic polarization tensor. As the density F[g, A] is local and gauge-invariant, one can formally consider it to be a function of A µ and the field strength F µν . We then have which defines the current J µ f and the anti-symmetric polarization tensor M µν . Of course, the exact way how one chooses to consider W [A] as a function of A µ and F µν is ambiguous. This ambiguity is the ambiguity of separating the charge/current into the components corresponding to "bound charges" and "free charges". While J µ f and M µν are ambiguous, the total current J µ defined by (2.4) is not, and is given by The first term can be called the current of free charges, and the second term the current of bound charges. A convenient choice of fixing the ambiguity in the definition of J µ f is to use (2.6c) to trade the derivatives of the chemical potential in the density F[g, A] for the electric field. This gives J µ f = ρu µ , where ρ ≡ ∂F/∂µ defines the density of free charges. Then to all orders in the derivative expansion. Note that M µν = −M νµ , and the bound current does not contribute to the conservation equation (2.5b). See [21] for more details about the electric and magnetic contributions to M µν .

Derivative expansion
We next specify the derivative counting. We choose the counting scheme in which the metric is g µν ∼ O(1), so that the Riemann tensor is O(∂ 2 ). Similarly, the temperature is T ∼ O(1). If the matter in question has degrees of freedom that carry "electric" charges (as would be in a conductor), the chemical potential is also µ ∼ O(1). The equilibrium condition (2.6c) then requires that the electric field is E µ ∼ O(∂). In an insulator, on the other hand, µ is not a relevant thermodynamic variable, and one can take E µ ∼ O(1). We will be considering conducting matter without macroscopic O(1) magnetic fields, and will take A µ ∼ O(1). At zeroth order in derivatives we then have only two invariants, T and µ. Thus the generating functional is where p(T, µ) is some function of T and µ (which is in fact the pressure), and the dots denote the terms of order O(∂) and higher. The functional form of p(T, µ) is to be determined from the microscopic theory. The energy-momentum tensor and the current which follow from the definitions (2.4) are where ǫ ≡ −p + T ∂p/∂T + µ∂p/∂µ is the energy density, n ≡ ∂p/∂µ is the charge density, and again the dots denote the terms of order O(∂) and higher. The conservation laws (2.5) are satisfied identically, simply because the above T µν and J µ were obtained from a gaugeand diffeomorphism-invariant generating functional. Let us specialize to 3+1 dimensions for definiteness. Then at order O(∂), there are no invariants that could appear in the generating functional. At order O(∂ 2 ), we write the generating functional as where the dots denote the terms of order O(∂ 3 ) and higher. The coefficients f n (T, µ) are the second-order thermodynamic susceptibilities (sometimes called thermodynamic transport coefficients) which need to be determined from the microscopic theory, just like the pressure. Finally, s n are the two-derivative invariants made out of the metric, the gauge field, and the quantities in (2.1), such as ∇ 2 T , a µ a µ , R, F µν F µν etc. The invariants must be such that they do not vanish in equilibrium. We will find it convenient to write the invariants in terms of the magnetic field B µ ≡ 1 2 ǫ µναβ u ν F αβ and the vorticity vector These are also true out of equilibrium. The vorticity tensor ω µν ≡ 1 2 ∆ µα ∆ νβ (∇ α u β −∇ β u α ) is related to the vorticity vector by ω µν = − 1 2 ǫ µνρσ u ρ Ω σ , so that ω µν ω µν = 1 2 Ω 2 . Not all invariants are independent: for example, (2.10b) shows that the Ω·a term in the generating functional may be absorbed into the E·Ω term after an integration by parts W n/a n/a 2 4 3 4 n/a 4 n/a and a redefinition of the f n coefficients. 4 Similarly, (2.10a) shows that the E·Ω term in the generating functional may be absorbed into the B·a and B·E terms after an integration by parts and a redefinition of the f n coefficients. Further, in equilibrium we have where R µν is the Ricci tensor. As a result, the u µ R µν u ν term in the generating functional may be absorbed into the Ω 2 , a 2 , and E·a terms after an integration by parts and a redefinition of the f n coefficients. The independent second-order invariants in the generating functional were classified in [6], using a dimensionally reduced formulation. There are seven independent invariants in a parity-preserving theory, and nine independent invariants in a parity-violating theory. If the matter degrees of freedom do not couple to the gauge field, there are only three independent invariants. We choose the independent invariants as listed in Table 1. This fixes the definition of the susceptibility coefficients f n (T, µ) in the generating functional (2.9). The table also lists the conformal weights w of the invariants under the Weyl rescaling of the metric g µν →g µν = e −2ϕ g µν . The quantity Φ has conformal weight w if under the Weyl rescaling Φ →Φ = e wϕ Φ. The zeroth-order invariants T and µ have w = 1. The acceleration transforms inhomogeneously,ã µ = a µ − ∂ µ ϕ. In a conformal theory, the generating functional W [g, A] must be invariant under the Weyl rescaling. While the invariants R, a 2 , and E·a do not have well-defined weights, in d + 1 spacetime dimensions the combination is Weyl-invariant up to a boundary term, for f (T, µ) = T d−1 F (µ/T ). Thus in a conformal theory in 3+1 dimensions, f 1 = T 2 F (µ/T ), f 2 = 6f 1 , f 7 = −6∂f 1 /∂µ, and f 9 = 0. In par-ticular, thermal equilibrium of a neutral conformal fluid in 3+1 dimensions is characterized by two independent second-order susceptibility coefficients f 1 and f 3 .

Trace anomaly
In a conformal theory, quantum effects give rise to the conformal anomaly [23]. For a conformal field theory subject to external gravitational and electromagnetic fields, the trace of the energy-momentum tensor is Here a and c are dimensionless coefficients that depend on the degrees of freedom of the theory. For example, for a free theory of N S real scalars, N F Dirac fermions, and N V vector fields, one has The coefficient b 0 is the coefficient of the leading-order beta function for the electromagnetic coupling e used to couple the theory to the external gauge field: , where M is the renormalization scale 5 . For example, for a free theory of n s complex scalars with charges q s,k (k = 1, . . . , n s ) and n f Dirac fermions with charges q f,i (i = 1, . . . , n f ), one has which gives the standard one-loop QED beta-function.
The gravitational contributions to the trace anomaly in (2.13) are fourth order in derivatives. As we are only interested in the generating functional (and hence the energymomentum tensor) up to second order in derivatives, we will ignore these contributions (see [25] for a general discussion in an arbitrary number of dimensions). The electromagnetic contribution in (2.13) is, however, second order in derivatives, and must emerge from the generating functional. This can be accounted for if the coefficients f 4 and f 6 are not themselves Weyl-invariant. In fact, the trace anomaly places constraints on the form of these coefficients.
Indeed, consider the "4,6" part of the generating functional, For the corresponding energy-momentum tensor T µν 4,6 one finds where f ′ n ≡ T f n,T + µf n,µ , and the comma denotes the derivative with respect to the argument that follows. For the trace anomaly, this has to match −b 0 5 The external electromagnetic field Aµ can be introduced by minimally coupling the fundamental matter fields of the theory to Aµ (without factors of e) and adding the kinetic term − 1 4e 2 F 2 µν to the action of the theory. The electromagnetic field becomes non-dynamical as e → 0. The last term in (2.13) describes the violation of scale invariance due to the renormalization-group running of e, caused by the charged matter fields of the theory. We have found the discussion in Section 2 of [24] helpful.
This is solved by where the integration constant M can be interpreted as the renormalization scale. This explicitly shows that f 4 and f 6 shift under the Weyl rescaling, due to the ln T terms. The total "electromagnetic" part of the generating functional is then where e 2 (M ) is the renormalized coupling, and F 4 (µ/T ), F 6 (µ/T ) are the renormalized susceptibilities. The renormalization-group equation for 1/e 2 (M ) ensures that M d dM W EM = 0, i.e. the generating functional does not depend upon the renormalization scale, as the case should be.

The energy-momentum tensor and the current
Let us now write down the energy-momentum tensor that follows from the generating functional (2.9). This was done in [6,7] for neutral matter, and in [19] for charged matter, in a dimensionally reduced formulation. Here we will write the energy-momentum tensor in the covariant form, decomposing T µν with respect to the fluid velocity u µ as T µν = Eu µ u ν + P∆ µν + Q µ u ν + Q ν u µ + T µν . (
As an example, consider matter that has a global U (1) charge (so that one can introduce the corresponding chemical potential), but which is not subject to any external electric and magnetic fields coupled to that U (1) current. An example is QCD at finite (or zero) baryon number chemical potential. A straightforward (and tedious) calculation gives the coefficients of the energy-momentum tensor (2.15) in terms of the three susceptibilities f n (T, µ) as where again f ′ n ≡ T f n,T + µf n,µ , f ′′ n ≡ T 2 f n,T,T + 2µT f n,T,µ + µ 2 f n,µ,µ , and the comma subscript denotes the partial derivative with respect to the argument that follows. The leading-order energy density is ǫ = −p + T p ,T + µp ,µ , as before. Equations (2.15) and (2.18) give the energy-momentum tensor of a relativistic fluid in hydrostatic equilibrium, up to O(∂ 2 ) terms beyond the perfect fluid approximation. This generalizes the result of [7] to non-zero µ.
The polarization vectors which determine the equilibrium current (2.17) are As an example, consider parity-invariant matter that has a global U (1) charge (so that one can introduce the corresponding chemical potential), but which is not subject to any external electric and magnetic fields coupled to that U (1) current. Again, QCD at finite (or zero) baryon number chemical potential would be an example. The charge density and the spatial current in (2.16) are then where n ≡ ∂p/∂µ is the zeroth-order charge density. We see that the thermodynamics of QCD with a baryon number chemical potential is specified by the pressure p(T, µ) at zeroth order in derivatives, as well as by the five susceptibilities f i (T, µ), with i = 1, 2, 3, 5, 7 at second order in derivatives. For a conformal theory, f 2 = 6f 1 , f 7 = −6∂f 1 /∂µ, f 9 = 0, and hence one only needs three susceptibilities f 1 , f 3 , and f 5 to specify the two-derivative thermodynamics of conformal matter not subject to external electromagnetic fields.

Kubo formulas
The above equilibrium expressions for T µν [g] and J µ [g, A] allow for a straightforward computation of equilibrium (zero frequency) correlation functions of the corresponding operators. In order to write down the correlation functions for matter at rest in flat space, we choose V α = (1, 0), and take the external sources as g µν = η µν + δg µν (k)e ik·x , A λ = µ 0 δ 0 λ + δA λ (k)e ik·x . The equilibrium two-point functions G AB of two operators A and B are then defined by varying the corresponding equilibrium one-point functions with respect to the source, The locality of the derivative expansion implies that the two-point functions are at most quadratic in k, with the coefficients of the O(k 2 ) terms determined by the susceptibilities f n . The energy-momentum tensor (2.15), (2.18) implies the following Kubo formulas in terms of the above zero-frequency correlation functions: There are of course other ways to write the Kubo formulas for f 1,2,3 which follow from the rotation invariance of the two-point functions. 6 Similarly, the equilibrium current (2.16), (2.17), (2.19) gives the following Kubo formulas Thus all seven parity-preserving thermodynamic susceptibilities admit Kubo formulas in terms of equilibrium two-point functions in flat space without external fields. The paritybreaking susceptibilities f 8 and f 9 do not appear in the above linearized analysis, but can 6 For example, (2f1−T f1,T −µf1,µ) = 1 4 lim k→0 As expected, in a CFT with T µ µ = 0 one recovers the constraints f1 = T 2 F (µ/T ), f2 = 6f1.  Figure 1. The one-loop diagram contributing to the T µν (k)T ρσ (−k) correlation function for fermions. The two vertices denote the two energy-momentum tensor insertions, the lines are the fundamental field propagators, and p denotes the loop momentum which needs to be integrated out. The hat inT µν signifies that it is an operator made out of the fundamental fields. Similar diagrams arise for free scalar as well as for free gauge fields.
be expressed in terms of equilibrium three-point functions in flat space without external fields.
In order to find the three-point functions, we expand the equilibrium T µν and J µ to quadratic order in small fluctuations h αβ (x), A λ (x). Note that we don't need to solve for the conservation of T µν in our setup -the conservation laws (2.5) are satisfied identically in equilibrium, as a consequence of the diffeomorphism invariance of W [g, A], and V being a Killing vector. For example, let us take h tt (z) and h tx (y) as the only non-vanishing perturbations. We then find As another example, let us take h tt (z) and A x (y) as the only non-vanishing perturbations. We then find Taking the variation of the one-point functions (2.29), (2.30) with respect to the sources h αβ , A α , one finds that the susceptibilities f 9 and f ′ 8 + f 9,µ are given in terms of the second derivatives of the appropriately defined three-point functions G J t T tx T tt (p, k) and G J t J x T tt (p, k), respectively.

Free fields
Let us now use the above Kubo formulas to evaluate the thermodynamic susceptibilities for non-interacting quantum fields in 3+1 dimensions. The energy-momentum tensor and the current are quadratic in the fields, hence the two-point functions can be evaluated from the diagram schematically shown in Figure 1. The diagram can be evaluated by the standard methods of equilibrium thermal field theory in flat space [26,27]. The integration over the intermediate momenta in the loop will give rise to ultraviolet divergences which can be regulated by introducing a high-momentum cutoff scale Λ. We will assume that the cutoff dependence is removed by the standard zero-temperature renormalization, and will only report the temperature-dependent (and cutoff-independent) contributions to the thermodynamic susceptibilities f n . [As an example, the zero-temperature contribution to f 1 (which we will not indicate explicitly) gives rise to the renormalization of Newton's constant by the quantum fluctuations of the matter fields. Similarly, the zero-temperature contributions to f 4 and f 6 give rise to electric charge renormalization.]

Scalars
We start with the massless real scalar field. The action is [28] where R is the Ricci scalar as before, and the dimensionless parameter ξ specifies the coupling to curvature. The scalar field is minimally coupled for ξ = 0, and conformally coupled for ξ = 1/6. The energy-momentum tensor of the theory is found by varying the action with respect to the metric, where G µν = R µν − 1 2 g µν R is the Einstein tensor. We use the hat to distinguish the microscopicT µν (which depends on the fundamental fields) from the macroscopic T µν defined by (2.4) (which only depends on the temperature, chemical potential, and the fixed external sources). The explicit metric dependence in theT µν operator will give rise to contact terms in two-point functions upon taking the metric variation. Schematically, We will be computing two-point correlation functions ofT µν in flat space, by taking the metric as g µν = η µν + δg µν . In this case φ 2 does not vary in space, and ∂ φ 2 vanishes. The term (∂φ) 2 contributes to the internal energy of the scalar field in equilibrium, and gives a constant momentum-independent contribution to the two-point function. Thus for the purpose of computing the thermodynamic susceptibilities (which appear as O(k 2 ) contributions to the two-point function), the only relevant contact term arises from φ 2 ∂ 2 g. The variation can be written as with the coefficients A, B, C that are local functions of the field φ. Expanding the Einstein tensor, we find C µν,αβ,ρσ = 1 2 ξφ 2 P µν,αβ,ρσ , with where the parentheses denote symmetrization (with the 1/2). Note that P µν,αβ,ρσ = P αβ,µν,ρσ . The "variational" two-point function G T µν T αβ defined by (2.21) is then related to the standard "operator" two-point function T µνT αβ by 7 The terms in the right-hand side of (3.1) may be evaluated diagrammatically by the standard rules of equilibrium thermal field theory in flat space in the Matsubara formalism.
The contact term contributes a simple "bubble" diagram with The susceptibilities f 1 , f 2 , f 3 may be computed from the Kubo formulas (2.22), (2.23), (2.24), using the diagrams schematically shown in Fig. 1. Performing the integral over the intermediate momenta, we find for the temperature-dependent contributions The rest of the susceptibilities f n vanish for the real scalar field. For a complex scalar field at µ = 0, the above f 1,2,3 get multiplied by a factor of 2. Minimally coupling the complex scalar field to the external gauge field A µ gives according to the general result (2.14). The rest of the susceptibilities f n vanish at µ = 0 by charge conjugation and parity.

Dirac fermions
We now consider a massless Dirac fermion field at µ = 0. The action is given by [28] Here γ µ are the spacetime dependent Dirac γ-matrices, γ µ (x) = e µ a (x) γ a , where e µ a (x) is the vierbein field, and γ a , a = 0, 1, 2, 3 are the usual position independent flat space γ-matrices. The Clifford algebras satisfied by the γ and γ matrices are The covariant derivative acting on the fermion field is given by where σ ab ≡ 1 4 [γ a , γ b ], and ω ab µ is the spin connection, ω ab µ = 1 2 e aν ∂ µ e b ν − ∂ ν e b µ − 1 2 e bν ∂ µ e a ν − ∂ ν e a µ + 1 2 e aν e bρ ∂ ρ e c ν − ∂ ν e c ρ e cµ . The energy-momentum tensor iŝ There are no terms in the energy-momentum tensor with two derivatives of the metric, hence there are no contact terms analogous to the ones we had for the scalar field. Hence we have where the right-hand side may be evaluated with the flat-space energy-momentum tensor (replacing the covariant derivatives in (3.4) with partial derivatives) by the standard rules of equilibrium thermal field theory in the Matsubara formalism. The Euclidean propagator is D(iω n , k) = / k[−(iω n ) 2 +k 2 ] −1 , with / k = γ 0 E ω n +γ·k, with γ 0 E = iγ 0 and ω n = (2n+1)πT , with n integer. We get Minimally coupling the Dirac current,Ĵ µ = −Ψγ µ Ψ, to the external gauge field A µ gives according to the general result (2.14). The rest of the susceptibilities f n vanish at µ = 0 by charge conjugation and parity.

Gauge fields
We now give results for the thermodynamic susceptibilities of a free U (1) gauge field. The action for the theory is given by with the energy-momentum tensor Once again there are no contact term contributions to the two-point function, and the relation (3.5) between the variational and operator definitions of the correlation function is valid. Proceeding in the same way as in [13] and evaluating the one-loop diagram similar to figure 1 by using the Euclidean propagator for the gauge field in the Feynman gauge, D µν (iω n , k) = δ µν [−(iω n ) 2 + k 2 ] −1 with ω n = 2πnT , one finds for the thermodynamic susceptibilities with the other susceptibilities vanishing.

Discussion
The emphasis of this note was on the Kubo formulas for thermodynamic susceptibilities that appear at two-derivative order in the constitutive relations of the energy-momentum tensor and of the global U (1) current (such as the baryon number current in QCD). Our work is close in spirit to [7,12]. Our main result is that all parity-even thermodynamic susceptibilities can be computed in terms of equilibrium two-point functions, while the earlier literature gave most of the susceptibilities in terms of equilibrium three-point functions. Explicitly, the Kubo formulas for parity-even susceptibilities are given by the equations (2.22) -(2.28), and the Kubo formulas for parity-odd susceptibilities are given by (2.29), (2.30). The Kubo formulas are applicable both at zero and non-zero chemical potential. The two-point functions are the zero-frequency equilibrium correlation functions in flat space, and can in principle be evaluated by the Euclidean methods, such as using lattice gauge theory. In fact, Ref. [17] has already performed a lattice evaluation of f 1 (or rather κ ≡ −2f 1 ) in the SU (3) Yang-Mills theory. We hope that the Kubo formulas derived in this paper will be useful for explicit calculations of the thermodynamic susceptibilities on the lattice as well as by holographic methods [3,4,18] in strongly interacting quantum field theories.
The equilibrium constitutive relations for the energy-momentum tensor and the current are written down in Eqs. (2.15), (2.18), and (2.16), (2.20), for a fluid not subject to external electromagnetic fields. 8 The constitutive relations are written in the "thermodynamic frame" [7], which means that the fluid velocity in equilibrium is aligned with the timelike Killing vector, according to the definition (2.1). In principle, one can redefine the thermodynamic variables and write down the constitutive relations in the "Landau-Lifshitz frame", which corresponds to a redefinition of T , µ, and u α , so that the expressions for E, N , Q µ written in terms of the new variables are made to look like E = ǫ, N = n, Q µ = 0.
While doing so is fine at order O(∂ 2 ), transforming to the Landau-Lifshitz frame (or any other frame) will also make the second-order susceptibilities f n (T, µ) appear at O(∂ 3 ) and higher in the constitutive relations, confusing their true two-derivative nature. 9 This makes any frame other than the thermodynamic frame ill-suited for a systematic understanding of thermodynamic contributions to the constitutive relations of T µν and J µ . Of course, the expectation values of T µν [g, A], J µ [g, A] and the corresponding correlation functions are physical objects, and do not depend on one's choice of "frame". 8 See [29] for fluids subject to a magnetic field (as would be relevant for magneto-hydrodynamics) and [21] for fluids subject to both electric and magnetic fields (as would be relevant for polarized fluids). 9 Similarly, in 2+1 dimensions, writing the constitutive relations in the Landau-Lifshitz frame will make the O(∂) thermodynamic susceptibilities also appear at O(∂ 2 ) and higher in the constitutive relations. The same comment applies to chiral anomalies: the thermodynamic frame is the frame in which the anomalous contributions to the constitutive relations only appear at one-derivative order [6,20]. Writing down the constitutive relations in any other frame (such as the Landau-Lifshitz frame) will make anomalous terms appear at O(∂), O(∂ 2 ), and higher in the constitutive relations, confusing their true one-derivative nature.
where as before f ′ n ≡ T f n,T +µf n,µ , and the comma subscript denotes the partial derivative with respect to the argument that follows. Note that the issue of "frame" transformations does not arise here, and the constitutive relations following from the partition function of [14] are in the same thermodynamic frame as ours.