Bulk viscosity of strongly interacting matter in the relaxation time approximation
Alina Czajka, Sigtryggur Hauksson, Chun Shen, Sangyong Jeon, Charles Gale
aa r X i v : . [ nu c l - t h ] A p r Bulk viscosity of strongly interacting matter in the relaxation time approximation
Alina Czajka,
1, 2
Sigtryggur Hauksson, Chun Shen, Sangyong Jeon, and Charles Gale Department of Physics, McGill University, 3600 rue University, Montreal, Quebec H3A 2T8, Canada Institute of Physics, Jan Kochanowski University, Swietokrzyska 15 street, 25-406 Kielce, Poland Department of Physics, Brookhaven National Laboratory, Upton, New York 11973-500, USA (Dated: April 27, 2018)We show how thermal mean field effects can be incorporated consistently in the hydrodynamicalmodeling of heavy-ion collisions. The nonequilibrium correction to the distribution function resultingfrom a temperature-dependent mass is obtained in a procedure which automatically satisfies theLandau matching condition and is thermodynamically consistent. The physics of the bulk viscosityis studied here for Boltzmann and Bose-Einstein gases within the Chapman-Enskog and 14-momentapproaches in the relaxation time approximation. Constant and temperature-dependent masses areconsidered in turn. It is shown that, in the small mass limit, both methods lead to the same valueof the ratio of the bulk viscosity to its relaxation time. The inclusion of a temperature-dependentmass leads to the emergence of the β λ function in that ratio, and it is of the expected parametricform for the Boltzmann gas, while for the Bose-Einstein case it is affected by the infrared cutoff.This suggests that the relaxation time approximation may be too crude to obtain a reliable form of ζ/τ R for gases obeying Bose-Einstein statistics.
1. INTRODUCTION
The vibrant experimental programs pursued at theRelativistic Heavy Ion Collider (RHIC) and at the LargeHadron Collider (LHC) have ushered in a new era of ex-ploration of systems governed by the nuclear strong in-teraction. One of the remarkable features that emergedfrom investigating the physics of relativistic heavy-ioncollisions is the fact that the created systems could bemodeled theoretically by relativistic fluid dynamics [1, 2].This realization led to developments in the formulationof relativistic viscous hydrodynamics in which observableconsequences of the dissipative effects were isolated [3–13]. Currently, second-order viscous hydrodynamics pro-vides a description of the fluid behavior [11–14] whichremedies the main failure of the Navier-Stokes – or first-order – formulation: acausal signal propagation and nu-merical instabilities plaguing relativistic systems.While the hydrodynamic equations are universal andprovide a macroscopic picture of a relativistic fluid be-havior in terms of conservation laws, transport coeffi-cients are governed by the underlying microscopic theorywhich must be used for their extraction. Although thefirst applications of viscous hydrodynamics focused onthe shear viscosity, it has recently become clear that bulkviscosity also plays an important role in the evolution ofthe QGP system [15–17]. The calculation of bulk viscos-ity from first principles, however, remains a challengingproject. It is on this aspect that we concentrate in thispaper.The equations of the second-order hydrodynamics de-scribe very efficiently the expansion of the system pro-duced in heavy-ion collisions. This is a strong indicationthat the system must thermalize very rapidly, which inturn indicates that the system is strongly interacting atpresently achievable energies. Current estimates of thebulk viscosity of QCD are mainly based on the equationof state obtained from lattice QCD simulations [19, 20], or rely on empirical extractions based on simulations ofrelativistic nuclear collisions [15–18]. Application of lat-tice QCD findings [21–23] and hadron resonance gas re-sults [24, 25] made it possible to determine that the bulkviscosity is notably enhanced near the critical tempera-ture of the QCD phase transition while the shear viscosityis substantially decreased in this region [26, 27]. Further-more, the importance of bulk viscosity near the transitiontemperature region was shown to have a remarkable im-pact on the elliptic flow coefficient v [25, 28, 29] andother heavy-ion observables [15–18, 30, 31]. Recently,the behavior of bulk viscosity was also obtained from hy-drokinetic theory, which incorporates thermal noise [32].Despite the progress described above, there is still aneed to develop methods which provide a better insightin the effects of bulk viscosity at different energy scales.In particular, one may be interested in having a consis-tent analytical approach to bulk viscosity physics in theregime of very high temperatures. At this energy scalethe coupling constant is small and fundamental quantumfield theoretical tools can be used to study bulk viscositysystematically. Having a comprehensive fluid dynamicformulation of a weakly coupled gas may also providean essential benchmark for different approaches and phe-nomenological applications.In Refs. [33, 34] it was shown that quantum field theoryis equivalent, at least at leading order of perturbative ex-pansion, to kinetic theory. Later calculations then coulduse this efficient and intuitive kinetic theory frameworkto study transport phenomena; see [35–37]. It has alsoprovided a natural language to formulate fluid dynamicsconcepts. Within the kinetic approaches, the Chapman-Enskog and Grad’s 14-moment methods are commonlyemployed to study the nonequilibrium processes of afluid. They, however, rely on different treatments of thedistribution function. While the Chapman-Enskog the-ory deals directly with solving the Boltzmann equation[38], Grad’s approach is based on an expansion of thenonequilibrium function in terms of the powers of mo-menta [39]. To date, great progress has been made inextraction of different transport coefficients within dif-ferent theories. It seems, however, that the comprehen-sive analysis of transport processes in a system exhibit-ing conformal anomaly is not yet complete, especially incases involving a mean field interaction.A violation of conformal symmetry has a different im-pact on different transport coefficients. It does not affectshear viscosity much: its leading order behavior is dom-inated by the kinetic energy scale in weakly interactingsystems. On the other hand, the breaking of scale in-variance dominates the physics of bulk viscosity. Conse-quently, the behavior of bulk viscosity is largely deter-mined by the sources of conformal symmetry breaking:either the physical mass of plasma constituents or theCallan-Symanzyk β λ function, which fixes the coupling asa function of the energy scale [33]. The parametric formof bulk viscosity should then be dictated by the sources ofscale invariance breaking squared, as shown in Ref. [40]for QCD. The bulk viscosity of systems exhibiting a con-formal anomaly, due to the presence of a constant massonly, was later studied within the Chapman-Enskog ap-proach and the 14-moment approximation, mostly in therelaxation time approximation [41–43], and also withinother approaches [44]. Moreover, quasiparticle modelswere also examined for systems of various matter con-tent in Refs. [45–53].We observe, however, that there is still a need to re-visit a formulation of nonequilibrium fluid dynamics withthe mean field background. Such a formulation is essen-tial when one needs to include variable thermal massesconsistently in the equations of viscous hydrodynamics.Having the correct form of a nonequilibrium momentumdistribution is also critical while studying some aspects ofnuclear matter behavior phenomenologically, in particu-lar, when implementing the Cooper-Frye prescription inhydrodynamic simulations or examining electromagneticprobes in heavy-ion collisions [17, 54–56]. Furthermore,such a consistent approach allows for an exhaustive cal-culation of transport coefficients.The central part of this paper is devoted to derivationof the nonequilibium correction to the distribution func-tion where thermal effects are consistently included. Sub-sequently, it is shown how the correction influences thebulk viscosity behavior in the relaxation time approxima-tion. The analysis is done systematically and it comprisesdifferent cases, namely, formulation of equilibrium andnonequilibrium fluid dynamics and then computation ofthe ratio of bulk viscosity to relaxation time. A computa-tion is provided for gases of Boltzmann and Bose-Einsteinstatistics in both the Anderson-Witting model of theChapman-Enskog method and the 14-moment approxi-mation. The analysis performed in this paper is specificto single-component bosonic degrees of freedom. Con-sequently, when the explicit forms of the thermal massand the β λ function are needed, we will use those of thescalar λφ theory [33, 34]. The method developed here is not appropriate for a one-component system following aFermi-Dirac distribution function. Such a system wouldbe a system of noninteracting fermionic degrees of free-dom where the thermal mass and bulk viscosity cannotbe determined. To count fermions accurately one needsto consider a many-component system with the inclusionof bosons mediating the interaction. This is not donehere and is left for future work.The correction to the distribution function is foundby noticing that there is a twofold source of departurefrom equilibrium. First, there are hydrodynamic forcesthat generate a deviation in the distribution function δf ,that is, they change the functional form of the distri-bution function. The other source is related directly tointerparticle interactions, the effect of which is statisti-cally averaged and emerges as the mean field. Therefore,the correction is expressed by two terms; for the Bose-Einstein gas the correction is∆ f = δf − T dm dT f (1 + f ) E k R dKδf R dKE k f (1 + f ) . (1.1)For the description of quantities, see Table I. The ob-tained form of the correction allows one to formulate hy-drodynamic equations in a coherent way, where the Lan-dau matching condition and thermodynamic relations areguaranteed. Since the thermal mean field has a negligibleimpact on shear viscosity, we further concentrate on bulkviscosity dynamics, where the influence of the thermalbackground reveals itself through the Landau conditionand the speed of sound.We show that both the Chapman-Enskog and the 14-moment approaches lead to the same final expressions forthe ζ/τ R ratio in the small mass limit, where τ R is thebulk relaxation time. In general, temperature-dependentmass results in the emergence of the β λ function, whichdictates the very high temperature form of the ratio. Inthe Boltzmann case the ratio is ζ Boltz τ R ≈ T (cid:18) − c s (cid:19) (cid:18) π − m x πT (cid:19) , (1.2)where (1 / − c s ) is directly related to M c , the noncon-formality parameter; see Table I. This shows the ex-pected behavior of the source of scale invariance break-ing. One may observe that one factor of the scale invari-ance breaking parameter is introduced directly by theLandau matching, which comes from the small depar-ture from equilibrium. The other factor emerges as acorrection to the pressure given by purely equilibriumquantities, but not provided by the equation of state, asargued in [40]. For a system with Bose-Einstein statis-tics, the result is ζτ R ≈ T (cid:18) − c s (cid:19) π T m x − π − m m x !! . (1.3)The leading order term is not of the expected depen-dence because of the factor T /m x , which comes froman infrared cutoff. The same behavior is reflected if weneglect either the constant mass term or thermally af-fected quantities. Therefore, it rather indicates that therelaxation time approximation, which assumes that τ R isenergy independent, may not allow one to entirely cap-ture microscopic physics, in particular, of soft momentain quantum gases following a Bose-Einstein distributionfunction. A similar conclusion was reached in Ref. [40].The paper is organized as follows. In Sec. 2 the in-gredients of the effective kinetic theory are briefly sum-marized and the derivation of the noneqilibrium thermal correction is provided. Section 3 is devoted to the for-mulation of fluid dynamic basic equations with the meanfield background. In Sec. 4, the analysis of the ratioof bulk viscosity to relaxation time ratio is presentedin the Chapman-Enskog theory, within which we solvethe Anderson-Witting model. In Sec. 5 we use the 14-moment approximation to derive the evolution equationfor the bulk pressure and then to calculate the bulk vis-cosity over the relaxation time ratio and other transportcoefficients in the bulk channel in the relaxation time ap-proximation. Sec. 6 summarizes and concludes the work.Appendices contain some technical details. Description Equilibrium quantity Nonequilibrium quantityPhysical, zero-temperature mass of a particle m m Quasiparticle thermal mass m eq m th Quasiparticle mass m x = p m + m ˜ m x = p m + m Quasiparticle energy E k = p k + m x E k = √ k + ˜ m x Quasiparticle four-momentum k µ ≡ ( k , k ) = ( E k , k ) ˜ k µ ≡ (˜ k , k ) = ( E k , k )Lorentz invariant measure dK = d k / [(2 π ) E k ] d K = d k / [(2 π ) E k ]Distribution function (in the local rest frame) f = 1 / [ e βE k − β = 1 /T f = f + ∆ f Beta function for a coupling constant λ β λ = T dλ/dT = 3 λ / (16 π )Temperature dependence of the thermal mass T dm /dT = m + aT β λ , with a = 1 / M = ( − m + aT β λ ) / c : m eq ,c , f ,c = e − βE k , f c , m th ,c , a c = 1 / (8 π ), and M c .
2. NONEQUILIBRIUM DEVIATION FROM THEEQUILIBRIUM DISTRIBUTION FUNCTIONA. Boltzmann equation with the mean field effect
Kinetic theory provides an efficient classical descrip-tion of complex microscopic dynamics of an interactingmany-body system. It is a good alternative to quan-tum field theory to study transport phenomena in theweakly coupled limit dominated by quasi-particle dy-namics. By quasiparticles one means particles which,apart from zero temperature mass, gain additional ther-mal mass due to interactions with the medium: the effectof the mean field. They are characterized by a mean freepath which is much larger than the Compton wavelengthof the system’s constituents, and by a mean free time,which is much larger than the time between collisions[36]. The dynamics of quasiparticles is encoded in thephase-space distribution function which evolves accord-ing to the Boltzmann equation. We consider a system of uncharged thermally influ-enced particles of a single species for which the Boltz-mann equation reads(˜ k µ ∂ µ − E k ∇E k · ∇ k ) f = C [ f ] , (2.1)where C [ f ] is the collision term, f = f ( x, k ) is a distri-bution function of quasiparticles, and the second termof the left-hand side involves the force F = d k /dt = −∇E k . The quasiparticle four-momentum is defined as˜ k µ = (˜ k , k ), where ˜ k ≡ E k is the nonequilibrium energygiven by E k = p k + ˜ m x , (2.2)which is a time- and space-dependent variable since ˜ m x ≡ ˜ m ( x ) = m + m ( x ), where m is the physical mass and We use here such a notation that whenever x and k appear asarguments of a function, we mean x µ and ˜ k µ (or k µ in the caseof f ), respectively. m th ( x ) is the nonequilibrium thermal mass, which variesin time and space. Knowing the x dependence of theenergy, one may rewrite Eq. (2.1) as (cid:0) ˜ k µ ∂ µ − ∇ ˜ m x · ∇ k (cid:1) f = C [ f ] . (2.3)The central object of the kinetic theory is the phase-space density function f ( x, k ). What we assume aboutthe system is that its departure from the equilibriumstate is small, which, in turn, means that the processof system equilibration is controlled by a small deviationin the distribution function, which we denote as∆ f ( x, k ) = f ( x, k ) − f ( x, k ) , (2.4)where f ( x, k ) is the equilibrium Bose-Einstein distribu-tion function and, in a general frame, it has the form f ( x, k ) = 1exp[ u µ ( x ) k µ ( x ) β ( x )] − , (2.5)where β ≡ β ( x ) = 1 /T ( x ) with T ( x ) being the local tem-perature, and u µ ≡ u µ ( x ) is the fluid four-velocity. Thefour-velocity in the local rest frame is u µ = (1 , , , k µ = ( k , k ), where k component is the equilibrium x -dependent energy E k = p k + m x , (2.6)where the dependence of x enters through the mass m x ≡ m ( x ) = m + m ( x ) with m ( x ) being the equi-librium thermal mass, which is not the same as m ( x ),the nonequilibrium thermal mass. The Bose-Einsteindensity function in the fluid rest frame takes the form f ( x, k ) = 1exp (cid:0) E k ( x ) β ( x ) (cid:1) − . (2.7)Let us add that in the forthcoming parts we will be de-riving all equations for the Bose-Einstein gas, but theseequations may be analogously found for the classicalBoltzmann gas with the distribution function f ,c ( x, k ) = exp( − β ( x ) u µ ( x ) k µ ( x )) (2.8)and these will be briefly presented as well. Our aim isto reformulate the equations of the viscous hydrodynam-ics when the effect of fluctuating thermal mass is incor-porated. Therefore, we assume that thermal influenceon the process of the system equilibration is controlledby the nonequilibrium correction to the thermal mass,∆ m = m − m , which will be specified further. B. Form of ∆ f As stated earlier, in this work we study systems withdistribution functions that are perturbed from their equi-librium value. More specifically, the nonequilibriumphase space density can be written as f ( x, k ) = f th ( x, k ) + δf ( x, k ) (2.9) The first part, f th ( x, k ), still retains the local-equilibriumform of the distribution function, but the thermal masscontains the nonequilibrium corrections f th ( x, k ) ≡ f ( x, k ) | m + m ( x ) → m + m ( x )+∆ m ( x ) (2.10)= (cid:20) exp (cid:16)q k + m + m ( x ) + ∆ m ( x ) β ( x ) (cid:17) − (cid:21) − . The second part, δf ( x, k ), is a change in the functionalform of f ( x, k ) caused by hydrodynamic forces, or equiv-alently, nonvanishing gradients of energy and momentumdensities. The nonequilibrium correction ∆ f then hastwo parts, ∆ f ( x, k ) = f ( x, k ) − f ( x, k )= δf ( x, k ) + δf th ( x, k ) , (2.11)where, to the leading order in small change, δf th ( x, k ) = f th ( x, k ) − f ( x, k ) is δf th ( x, k ) = − f ( x, k ) (cid:0) f ( x, k ) (cid:1) ∆ m ( x )2 E k ( x ) β ( x ) , (2.12)which is obtained by expanding f th . Since ∆ m is thenonequilibrium deviation, it itself is going to be a func-tional of ∆ f . Hence, the equation∆ f = δf − βf (1 + f ) ∆ m E k (2.13)must be solved self-consistently for ∆ f . C. Form of ∆ m Recalling the basic foundations of effective kinetic the-ory, the analysis here relies heavily on findings within thescalar λφ theory, as provided in Refs. [33, 34], whichmakes the introduction of thermal corrections analyt-ically feasible. But the analysis presented here worksequally well whenever the equilibrium thermal mass hasthe form ∼ g n T , where g is the dimensionless couplingconstant and n is a positive integer. We intend to pro-vide an effective macroscopic framework to study weaklyinteracting systems, where the strength of interaction isdetermined by the coupling constant λ ≪
1. The cou-pling constant is scale (temperature) dependent and theanalysis performed here pertains only to the perturbativeregime. Within this approach the equilibrium thermalmass is found to be m = λ ( q )2 q , (2.14)where we have introduced the equilibrium scalar quantity q . The function q and its nonequilibrium counterpart q are defined through the corresponding distribution func-tions as q = Z dKf , (2.15) q = Z d K f. (2.16)For the definitions of the symbols, see Table I. Therefore,one can observe that Eq. (2.14) contains the couplingconstant λ ( q ), which is temperature dependent since q is temperature dependent.Throughout the analysis we always keep the assump-tion that all nonequilibrium quantities are slowly vary-ing functions of space points, which justifies that thenonequilibrium dynamics is governed by small deviationsof the quantities from their equilibrium values. There-fore, we further assume that the nonequilibrium thermalmass is a function of the scalar quantity q only. Thesame assumption is applied to the running coupling λ ( q ).Thus, the nonequilibrium thermal mass can be expandedas m ( q ) = m ( q + ∆ q ) = m ( q ) + ∆ m (2.17)with ∆ m = dm dq ∆ q. (2.18)The function q is uniquely defined by Eq. (2.16) andshould be obtained self-consistently from this equation.Hence to evaluate ∆ m , we need to find ∆ q which isitself a function of ∆ m . The deviation of the scalarquantity q can be written as∆ q = Z dKδf + ∂q ∂m ∆ m . (2.19)Equation (2.18) then takes the form∆ m = 11 − dm dq ∂q ∂m dm dq Z dKδf. (2.20)On the other hand both m and q are related by tem-perature, so that one can find dm dT = dm dq dq dT = dm dq (cid:18) β Z dKE k f (1 + f ) + dm dT ∂q ∂m (cid:19) . (2.21)Extracting further dm dq ∂q ∂m and inserting it toEq. (2.20) leads to∆ m = 2 T dm dT R dKδfβ R dKE k f (1 + f ) , (2.22)where we used dm /dT = 2 T dm /dT .Inserting Eq. (2.22) into Eq. (2.13), one gets∆ f = δf − T dm dT f (1 + f ) E k R dKδf R dKE k f (1 + f ) . (2.23)Analogously, the correction for the Boltzmann gas is∆ f c = δf c − T dm ,c dT f ,c E k R dKδf c R dKE k f ,c , (2.24) where the subscript c has been used to emphasize thatthe formula holds for the classical gas. Equations (2.23)and (2.24) are one of the main results of this paper.In previous analyses [45–47, 49–53], the second term inEq. (2.23) was missing or was incomplete. When apply-ing the Cooper-Frye formula in viscous hydrodynamics,it is ∆ f , not δf that should be used. D. Temperature dependence of the thermal mass
The thermal mass is a function of the scalar quantity q and is defined by Eq. (2.14). Its temperature dependenceis dictated by dm dT = λ ( q )2 dq dT + q dλ ( q ) dT . (2.25) q is one of the thermodynamic functions discussed in de-tail in Appendix B, and its leading order value is found tobe T /
12. Additionally, the second term in Eq. (2.25) en-codes the running of the coupling constant as a functionof the energy scale, which is the essence of the renormal-ization group β λ function, defined by β λ ≡ β ( λ ) = T dλ ( q ) dT . (2.26)It should be obtained using diagrammatic methods. Inthe case of scalar theory, β λ is positive and proportionalto λ . Collecting these contributions, one finds T dm dT = m + aT β λ . (2.27)where m = λT /
24 and a = 1 / f ,c instead of f . This gives q c = T / (2 π ) + O ( m x ), asgiven by Eq. (B.16), and it leads to T dm ,c dT = m ,c + a c T β λ , (2.28)where m ,c = λT / (4 π ) and a c = 1 / (8 π ).
3. EQUATIONS OF HYDRODYNAMICS WITHTHERMAL CORRECTIONSA. Local equilibrium hydrodynamics
First consider a system under strict local equilibrium.By that we mean that the functional form of the distribu-tion function is still f given in Eq. (2.5) or in Eq. (2.8),but the temperature as well as the thermal mass are x dependent. Such a system possesses a conserved stress-energy tensor of the form T µν = Z dKk µ k ν f − g µν U , (3.1)where the metric tensor we use is g µν =diag(1 , − , − , − U ≡ U ( x ) isthe mean-field contribution that guarantees the thermo-dynamic consistency of hydrodynamic equations and theconservation of energy and momentum, via the followingcondition: dU = q dm , (3.2)where q is the Lorentz scalar defined by Eq. (2.15).Since we study here a system with no conservedcharges, the Landau frame is a natural kinetic frameworkto define the four-velocity u µ via u µ T µν = ǫ u ν , (3.3)where the eigenvalue ǫ can be identified as the local en-ergy density. With this definition the energy-momentumtensor may be decomposed using two orthogonal projec-tions u µ u ν and ∆ µν = g µν − u µ u ν . The equilibriumenergy-momentum tensor becomes T µν = ǫ u µ u ν − P ∆ µν , (3.4)where P is the local thermodynamic pressure. The en-ergy density and the pressure are in turn given by ǫ = ¯ ǫ − U , (3.5) P = ¯ P + U , (3.6)where ¯ ǫ = (cid:10) ( u µ k µ ) (cid:11) , (3.7)¯ P = − (cid:10) ∆ µν k µ k ν (cid:11) (3.8)with the notation (cid:10) . . . (cid:11) = R dK . . . f . Let us pointout that the enthalpy is not changed by the mean field¯ ǫ + ¯ P = ǫ + P . One may also check that the definitionsof energy density (3.5) and pressure (3.6), together withthe condition (3.2), guarantee that the thermodynamicrelation T s = T dP dT = ǫ + P , (3.9)where s is the entropy density, is satisfied. B. Nonequilibrium hydrodynamics
The stress-energy tensor of fluid dynamics out of equi-librium takes the following form: T µν = Z d K ˜ k µ ˜ k ν f − g µν U, (3.10) which is formally the same as Eq. (3.1). The mean-fieldcorrection U must be now a function of q = R d K f only[34]. We emphasize that the formulation of the fluid hy-drodynamic framework with the thermal correction stillhas to conform with all assumptions that were made toprovide the effective kinetic theory, discussed in Sec. 2.In particular, such a description requires the system tobe sufficiently dilute and the quasiparticles’ mean freepaths to be much longer than the thermal width of itsconstituents, which is maintained when the strength ofinteraction is weak. Furthermore, to allow for validityof hydrodynamics, the system has to be characterizedby some macroscopic length scale at which macroscopicvariables, such as pressure and energy density, vary. Un-der these assumptions, a nonequilibrium hydrodynamicdescription applies to systems where departures of allquantities from their equilibrium values are character-ized by small corrections. Therefore, the nonequilibriumfunction U , in particular, may be expanded as U = U + ∆ U, (3.11)where ∆ U = dU dq ∆ q. (3.12)However, as discussed before and explicitly shown byEqs. (2.17) and (2.18), the thermal mass is also a func-tion of q only. Therefore, applying the relation (2.18) to(3.12), one finds ∆ U = q m . (3.13)As before, this is also the condition that U must sat-isfy to maintain the energy-momentum conservation law ∂ µ T µν = 0.The stress-energy tensor of the viscous hydrodynamics(3.10) may be next decomposed into the local equilibriumpart and the nonequilibrium deviation T µν = T µν + ∆ T µν , (3.14)where T µν is given by (3.4) and ∆ T µν carries all dynam-ical information needed in order to determine how thenonequilibrium system evolves into equilibrium. Notethat a separation of the viscous correction from the equi-librium part in Eq. (3.14) has been done not as a re-arrangement of Eq. (3.10) but rather as an expansionof the stress-energy tensor around its local equilibriumvalue. As shown in Appendix A, we have∆ T = Z dKE k ∆ f, (3.15)∆ T i = Z dKE k k i ∆ f, (3.16)∆ T ij = Z dKk i k j ∆ f − ∆ m Z dK k i k j E k f + δ ij ∆ m Z dKf , (3.17)where ∆ m and ∆ f are given by (2.22) and (2.23), re-spectively. Equations (3.15) and (3.16) shall dictate theform of the Landau matching condition, and Eq. (3.17)contains the definitions of the viscous corrections. C. Landau matching condition in the rest frame
The Landau matching is defined by the eigenvalueproblem u µ T µν = ǫu ν , (3.18)where ǫ is the energy density of the nonequilibrium stateincluding the thermal correction U . In the fluid restframe it comes down to two equations, corresponding tothe conditions on the energy density and the momentumdensity: T = ǫ, T i = 0 . (3.19)Under the Landau matching condition, the local equilib-rium is defined to have the same local energy and themomentum density∆ T = 0 , ∆ T i = 0 . (3.20)Using Eqs. (3.15) and (3.16) with the correction to thedistribution function ∆ f given by Eq. (2.23), we obtain∆ ǫ = Z dK (cid:20) E k − T dm dT (cid:21) δf, (3.21)0 = Z dK (cid:20) E k k i − T dm dT R dK ′ k ′ i f ( f + 1) R dK ′ E ′ k f ( f + 1) (cid:21) δf. (3.22)However, the second term in Eq. (3.22) vanishes becauseof rotational symmetry in equilibrium. Hence the Landaumatching conditions are Z dK (cid:20) E k − T dm dT (cid:21) δf = 0 , (3.23) Z dKE k k i δf = 0 . (3.24)The second condition indicates that δf cannot have avector component: it can only contain a spin 0 part anda spin 2 part. D. Shear-stress tensor and bulk pressure in thelocal rest frame
The shear tensor π ij and the bulk pressure Π are foundfrom Eq. (3.17) in the local rest frame, where Eqs. (2.22)and (2.23) are inserted. Then, as shown in Appendix A,one obtains ∆ T ij = Z dKk i k j δf. (3.25) We can reorganize (3.25) to separate the spin 0 part andthe spin 2 part as follows:∆ T ij = π ij + δ ij Π , (3.26)where π ij = Z dKk h i k j i δf, (3.27)Π = 13 Z dK k δf, (3.28)where k h i k j i = k i k j − k δ ij /
3. These coincide with thecommonly known forms of the shear-stress tensor andbulk pressure in the local rest frame.
E. General frame
In a general frame where the flow velocity u µ may bearbitrary, the energy-momentum tensor is T µν = Z dKk µ k ν f − g µν U + Z dK (cid:20) k µ k ν − u µ u ν T dm dT (cid:21) δf. (3.29)The Landau condition then becomes Z dK (cid:20) ( u µ k µ ) k ν − u ν T dm dT (cid:21) δf = 0 (3.30)and the viscous corrections are given by π µν = (cid:10) k h µ k ν i (cid:11) δ , Π = − (cid:10) ∆ µν k µ k ν (cid:11) δ , (3.31)where h . . . i δ ≡ R dK ( . . . ) δf . We have also used the no-tation A h µν i ≡ ∆ µναβ A αβ , where ∆ µναβ ≡ (∆ µα ∆ νβ +∆ µβ ∆ να − / µν ∆ αβ ) /
2. The definitions (3.31) have well-knownstructures, but the thermal mass that enters them is now x dependent and the Landau matching contains a correc-tion due to the temperature-dependent mass. These ar-guments are essential when one aims at examining trans-port properties of the medium.
4. NONEQUILIBRIUM CORRECTION IN THECHAPMAN-ENSKOG APPROACH
Chapman-Enskog theory provides a way to directlyfind the solution to the Boltzmann equation for near-equilibrium systems. Solving the full Boltzmann equa-tion, however, is formidable task. In this paper, we In Ref. [33], the energy-momentum tensor correction was writtendown incorrectly, but the mistake vanished with the impositionof the Landau matching condition, ensuring the validity of thesubsequent derivations. use the Anderson-Witting model [57] to find the explicitleading order solution. In this section, we focus on thebosonic quantum gas case. Treatment for the Boltzmanngas case is identical if one replaces f (1 + f ) with theBoltzmann factor f ,c . A. Solution of the Anderson-Witting equation inthe rest frame
With the medium-dependent thermal mass, theAnderson-Witting model is given by (cid:16) ˜ k µ ∂ µ − E k ∇E k · ∇ k (cid:17) f = − ( u · ˜ k ) τ R ∆ f, (4.1)where ˜ k µ = ( E k , k ). In the fluid cell rest frame u µ =(1 , , ,
0) and u · ˜ k = E k .To use the Chapman-Enskog method, we let f = f + f + f + · · · (4.2)where each f n contains only the n -th derivatives of thethermodynamic quantities and the flow velocity. Thefirst-order equation is obtained by identifying ∆ f = f in the right-hand side and using all other quantities intheir equilibrium forms (cid:18) k µ ∂ µ − ∂ i m ∂∂k i (cid:19) f ( x, k ) = − E k τ R ∆ f ( x, k ) , (4.3)where now k µ = ( E k , k ).Evaluating the left-hand side yields (cid:16) k µ ∂ µ − ∂ i m ∂∂k i (cid:17) f ( x, k ) = − βf ( x, k )(1 + f ( x, k )) × (cid:20)(cid:18) c s E k − T dm dT ! − k (cid:19) ( ∂ i u i ) − k h j k i i ∂ j u i (cid:21) , (4.4)where the equations of motion from the ideal hydrody-namics ∂ u i = ∂ i TT , (4.5) ∂ T = − T c s ∂ i u i (4.6)are used to remove time derivatives.The ∆ f in the right-hand side of the Anderson-Wittingmodel is just Eq. (2.23). Letting δf = f (1 + f ) φ , weget∆ f ( k ) = f ( k )(1 + f ( k )) × φ ( k ) − T E k dm dT R dKφ ( k ) f ( k )(1 + f ( k )) R dKE k f ( k )(1 + f ( k )) ! , (4.7)where the x dependence of all quantities is suppressedfor the sake of brevity. In previous derivations, the lastterm was missing [47, 52, 53]. Dividing φ into the shear and the bulk parts φ = φ s + φ b , and comparing Eqs. (4.4)and (4.7), the shear part of φ is trivially obtained as φ s ( k ) = − τ R T E k k h j k i i ∂ j u i , (4.8)since the angle integration over the spin-2 tensor k h j k i i vanishes. For the bulk part, letting φ b ( k ) = (cid:18) aE k + bE k (cid:19) ∂ i u i (4.9)and comparing Eqs. (4.7) and (4.4), we get a = τ R β (cid:18) c s − (cid:19) (4.10)and b = − M βτ R J , J , − T ( dm /dT ) J − , , (4.11)where we defined M = − m x − T dm dT ! . (4.12)With m ∝ λT , we have M = − (cid:0) m − aβ λ T (cid:1) , (4.13)where β λ is the coefficient function of the coupling con-stant renormalization group and a = O (1) depends onthe theory. The parameter M can be identified as the pa-rameter of nonconformality of the system (or the sourceof the conformal invariance violation). We have also in-troduced a notation for thermodynamic integrals, J n,q = a q Z dK ( u · k ) n − q ( − ∆ µν k µ k ν ) q f ( k )(1 + f ( k )) , (4.14)where a q = 1 / (2 q + 1)!!, which can be evaluated in thefluid cell rest frame. The bulk part of the leading or-der Chapman-Enskog solution of the Anderson-Wittingequation is then φ b ( k ) = τ R β ( ∂ i u i ) × (cid:18) ( c s − / E k − E k M J , J , − T ( dm /dT ) J − , (cid:19) . (4.15)To show that φ b ( k ) is in fact proportional to ( c s − / c s = dP /dTdǫ /dT = J , J , − ( T dm /dT ) J , , (4.16)where P and ǫ are the pressure and the energy densitygiven in Eqs. (3.6) and (3.5). Using the identities fromAppendix B 2, one can also show that13 − c s = − M J , J , − T ( dm /dT ) J , . (4.17)Hence finally φ b ( k ) = τ R β ( ∂ i u i )( c s − / × E k − E k J , − T ( dm /dT ) J , J , − T ( dm /dT ) J − , ! . (4.18)Equation (4.18) is another main result in this work. Thisequation slightly differs from the analogous one for theBoltzmann statistics shown in Ref. [17, 46, 49].In hydrodynamic simulations, it is practical to replacethe system expansion rate by the bulk viscous pressureusing the Navier-Stokes relation Π = − ζθ , which gives φ b ( k ) = β (cid:18) − Π ζ/τ R (cid:19) ( c s − / × E k − E k J , − T ( dm /dT ) J , J , − T ( dm /dT ) J − , ! . (4.19)Having given the solution of the Anderson-Wittingequation, one can also find ∆ f explicitly. InsertingEqs. (4.18) and (4.8) into (4.7) one finds∆ f ( k ) = f ( k )(1 + f ( k )) τ R β (cid:20) − ( ∂ j u i ) k h j k i i E k +( ∂ i u i )( c s − / (cid:18) E k − E k J , J , (cid:19)(cid:21) . (4.20)The phase space density correction ∆ f has a much sim-pler form than φ . However, for transport coefficient cal-culations, it is φ (equivalently δf ), rather than ∆ f , thatis needed. B. Energy conservation and Landau matching inthe Anderson-Witting case
By multiplying ˜ k ν = ( E k , k ) and integrating over d K ,the left-hand side of Anderson-Witting equation (4.1)turns into ∂ µ T µν , where the stress-energy tensor T µν isdefined in Eq. (3.10). Assuming that the mean-field con-tribution U satisfies ∂ µ U ( x ) = ∂ µ ˜ m x ( x )2 Z d K f ( x, k ) , (4.21)we get ∂ µ T µν = 0.Under the same condition, the right-hand side of theAnderson-Witting model within the Chapman-Enskogapproach must also vanish, − τ R Z dK E k k µ ∆ f = 0 , (4.22)to ensure energy-momentum conservation. This condi-tion for energy-momentum conservation is actually ex-actly the same as the Landau conditions we derived inSec. 3 C. Upon using ∆ f in Eq. (2.23) in the fluid restframe, these become0 = Z dK E k − T dm dT ! δf (4.23) and 0 = Z dKE k k i δf. (4.24)Eq. (4.24) is automatically satisfied by the δf = f (1 + f )( φ s + φ b ) obtained in the previous subsection since itdoes not contain a vector part. In the condition (4.23),the shear part φ s also vanishes because it contains a spin-2 tensor. Using Eqs. (4.18) and (4.14), it is easy to showthat the energy conservation and the Landau conditionare indeed fulfilled. This automatic fulfillment of theLandau condition for the quasiparticle case would nothave been possible if one missed the ∆ m correction in∆ f . C. The shear and the bulk viscosities in theAnderson-Witting model
The full leading order Chapman-Enskog solution to theAnderson-Witting model is given by Eq. (4.7) with φ s and φ b obtained above. The shear viscosity can be eval-uated by using Eq. (3.27) for π ij and Eq. (4.8) for φ s as π ij = 2 β τ R Z dK f (1 + f ) k E k σ ij , (4.25)where σ ij = − / ∂ i u j + ∂ j u i − / g ij ∂ k u k ). Identifying π ij = 2 ησ ij , we get ητ R = βJ , (4.26)and subsequently find the shear viscosity in the relax-ation time approximation, which was examined in fewpapers, see, for example, [41, 43, 49], and has the form ητ R = ǫ + P . (4.27)For the bulk viscosity, we start with Eq. (3.28)Π = Z dK k δf. (4.28)Using the Landau condition, Eq. (3.23), one getsΠ = M Z dKδf, (4.29)in which only the bulk part is relevant:Π = M Z dK f ( k )(1 + f ( k )) φ b ( k ) (4.30)with φ b ( k ) given by Eq. (4.18). Since Π = − ζ∂ i u i , onecan read off the ratio of bulk viscosity to the relaxationtime from Eq. (4.30) as ζτ R = βM (cid:18) J − , J , J , − T ( dm /dT ) J − , − J , J , J , − T ( dm /dT ) J , (cid:19) . (4.31)0The integrals present in Eq. (4.31) have been computedin Appendix B 2. Using them, one gets the value of theratio as ζτ R ≈ M π (cid:18) πT m x − (cid:18) − m m x (cid:19)(cid:19) . (4.32)For the application in relativistic viscous hydrodynamics,it is more useful to use the speed of sound. ApplyingEqs. (4.17) and (4.31), one can explicitly show that theratio is proportional to (1 / − c s ) , namely ζτ R ≈ T (cid:18) − c s (cid:19) π T m x − π − m m x !! . (4.33)Note the appearance of T /m x in the expression (4.33).This is in clear contrast to the Boltzmann statistics casewhich does not show such a behavior. The analysis forthe Boltzmann statistics case is identical to the analysisabove except that in place of J n,q we have I n,q = a q Z dK ( u · k ) n − q ( − ∆ µν k µ k ν ) q f ,c ( k ) , (4.34)where f ,c ( k ) = e − βk µ u µ . In this case, one gets ζ Boltz τ R ≈ T (cid:18) − c s (cid:19) (cid:18) π − m x πT (cid:19) . (4.35)The origin of this discrepancy is the fact that the Bose-Einstein factor behaves like f ( k ) ∼ T /E k in the infraredlimit, which makes the thermodynamic integral J − , inEq. (4.31) diverge in the m x → I − , doesnot. As a result, soft momenta govern the structure of ζ/τ R . However, since the calculation was performed inthe relaxation time approximation, which assumes that τ R is independent of energy, it may not capture the rightsoft physics. A similar behavior was seen in Ref. [40],where QCD bulk viscosity is studied. The authors claimthat the correct behavior of bulk viscosity is obtainedin the relaxation time approximation by neglecting theinfrared divergent term. But in principle there is no rea-son why this term should be ignored within the presentframework.Further, notice that starting from Eq. (4.4), the spin 0part (the bulk part) and the spin 2 part (the shear part)of the analysis are totally independent. Hence, it is pos-sible to generalize the leading order Anderson-Wittingequation as (cid:18) k µ ∂ µ − ∂ i m ∂∂k i (cid:19) f ( x, k )= − E k τ π ∆ f s ( x, k ) − E k τ Π ∆ f b ( x, k ) , (4.36)where ∆ f s and ∆ f b are the shear and bulk parts of ∆ f .In fact, when the dominant physical processes for theshear relaxation and the bulk relaxation are different, thisis the most natural form of the Anderson-Witting model. The analysis of this generalized Anderson-Witting modelfollows exactly the same route as for the single τ R , ex-cept that the shear viscosity and the bulk viscosity havedifferent relaxation times.As discussed in Refs. [33, 34], the dominant physicalprocesses for the shear relaxation and the bulk relaxationcan be indeed very different, and the bulk relaxation canbe dominated by the soft sector. Hence, the appearanceof T /m x is not entirely unnatural given that τ Π can havevery different m x dependence from τ π and the bulk relax-ation is dominated by the soft number-changing process. D. Comparison of ∆ f to previous works The phase space correction ∆ f in Eq. (4.20) ultimatelycomes from solving the first-order Chapman-Enskog ap-proximation. Hence, it should come as no surprise thatEq. (4.20) is consistent with similar results found in othersimilar works, provided that the right expression for thespeed of sound is used. For instance, in Ref. [49] onefinds that the bulk part of the phase space correction inthe Boltzmann case is derived to be∆ f R ( k ) = f ,c ( k ) φ R ( k ) (4.37)with φ R ( k ) = τ R β ( ∂ i u i ) (cid:18) ( c sR − / E k (4.38) − E k (cid:18) c sR m x T dm x dT − m x (cid:19)(cid:19) , where the speed of sound is c sR = (3 + zK ( z ) /K ( z )) − ,with z = m x /T and K n ( z ) being the modified Besselfunctions of the second kind. This φ R is different than φ b in Eq. (4.19) since φ R is a part of ∆ f while φ b is apart of δf . The phase space correction ∆ f R is, however,equivalent to the bulk part of ∆ f in Eq. (4.20) if one usesthe speed of sound expression (4.16) with J n,q → I n,q . Asmentioned above, this is as it should be since both aresolutions of the first-order Chapman-Enskog approxima-tion.The big difference between the previous treatmentsand ours is in computing the bulk viscosity. The bulkviscosity must be calculated using δf and not ∆ f as ex-plained in the previous section. If one uses ∆ f (or ∆ f R )instead of δf , the ratio ζ/τ R would be incorrectly calcu-lated.
5. TRANSPORT COEFFICIENTS IN THE14-MOMENT APPROXIMATION
When a system features a conformal anomaly, first-order transport coefficients reveal different sensitivity tothe source of the conformal symmetry violation, as ex-plicitly shown in the previous section. In particular,1shear viscosity is fully determined by the dominant en-ergy scale, which is the temperature T , and thus theshear viscosity over its relaxation time ratio behaves as T at leading order in the conformal symmetry break-ing, making the effects of scale anomaly negligible. Onthe other hand, bulk viscosity over the relaxation time isfully determined by the breaking of conformal symmetry.Such a difference makes it justified to omit the analysisof shear viscous effects and to evaluate first- and second-order transport coefficients related to bulk pressure, be-cause the additional term in Eq. (2.23) indeed concernsonly the scalar part. The analysis is performed at leadingorder in the conformal breaking parameter while includ-ing the thermal mass consistently.The bulk pressure is given by Eq. (4.29). Noting thatEq. (2.23) can be expressed as M ∆ f = M δf − T dm dT f (1 + f ) E k Π R dKE k f (1 + f ) , (5.1)one can rewrite Eq. (4.29) asΠ = ˜ M Z dK ∆ f, (5.2)where ˜ M = M J , J , − T ( dm /dT ) J − , . (5.3)To obtain the equation of motion for the bulk pressure,we first take the time derivative of Π,˙Π = ˙˜ M Z dK ∆ f + ˜ M (cid:20) Z dK ∆ ˙ f − ˙ m Z dK E k ∆ f (cid:21) , (5.4)where we adopted the notation ˙ A = u µ ∂ µ A for an arbi-trary quantity A , which reduces to the time derivative inthe rest frame of the fluid. From the Boltzmann equation (cid:16) ˜ k µ ∂ µ − E k ∇E k · ∇ k (cid:17) f = C [ f ] , (5.5)where C [ f ] is the collision integral, one finds u µ ∂ µ (∆ f ) = 1( u · ˜ k ) (cid:20) C [ f ] − ˜ k µ ∂ µ f − ˜ k µ ∇ µ ∆ f + 12 ∇ ˜ m x ∇ k f + 12 ∇ ˜ m x ∇ k ∆ f (cid:21) . (5.6)Inserting the expression (5.6) to Eq. (5.4) and keepingonly leading order terms, that is, terms which are evalu-ated with ˜ k → k , we have˙Π − C = − ˜ M (cid:20) − ˙ β (cid:16) J , − T ( dm /dT ) J − , (cid:17) + β θ (cid:16) J , − m x J − , (cid:17)(cid:21) + ˙˜ M ˜ M − θ ! Π − M (cid:16) ˙ m m x θ (cid:17) ρ − − M ρ µν − σ µν , (5.7) where θ ≡ ∇ µ u µ and σ µν = ∂ h µ u ν i is the Navier-Stokesshear tensor. In Eq. (5.7) we adopted the following no-tation for the collision term: C = ˜ M Z dK ( u · k ) − C [ f ] (5.8)and, for the irreducible moments, ρ n = h ( u α k α ) n i δ , ρ µνn = h ( u α k α ) n k h µ k ν i i δ . (5.9)Evaluating u ν ∂ µ T µν = 0 and implementing the for-mula (4.16) for the speed of sound squared, one obtains˙ β = Π θ − π µν σ µν J , − T ( dm /dT ) J , + c s βθ. (5.10)Next, calculating time derivatives ˙˜ M and ˙ m , Eq. (5.7)simplifies to˙Π − C = − β ˜ M (cid:20)(cid:18) − c s (cid:19) J , − T dm dT J − , ! + M J − , (cid:21) θ − (cid:16)
23 + 2 c s aT β λ M − A (cid:17) θ Π − π µν σ µν A + M ρ − θ − M ρ µν − σ µν , (5.11)where A = ˜ M J , − T ( dm /dT ) J − , J , − T ( dm /dT ) J , = c s −
13 (5.12)with the quantity ( c s − /
3) given by Eq. (4.17).To close Eq. (5.11) in terms of Π and π µν , one canapply the 14-moment approximation, which allows one toexpress the irreducible moments by Π and π µν as follows: ρ − = γ (0)2 Π , (5.13) ρ µν − = γ (2)2 π µν , (5.14)where the coefficients γ (0)2 and γ (2)2 are combinations ofdifferent thermal functions J n,q . Their particular formsare presented in Appendix C. Also, using the Anderson-Witting model for the collision term C [ f ] = − ( u · k ) ∆ fτ R , (5.15)where ∆ f is given by Eq. (2.23), the collision integralbecomes C = − Π τ R . (5.16)Applying the collision term in the relaxation time ap-proximation (5.16), the irreducible moments, Eqs. (5.13)and (5.14), and the relation for the speed of sound (4.17)to the evolution equation (5.11), one obtains˙Π + Π τ R = − ζθτ R − δ ΠΠ τ R θ Π + λ Π π τ R π µν σ µν , (5.17)2where ζτ R = βM (cid:20) J , J − , J , − T ( dm /dT ) J − , − J , J , J , − T ( dm /dT ) J , (cid:21) (5.18)is identical to the expression obtained in the Chapman-Enskog approach found in the previous section,Eq. (4.31). The remaining transport coefficients are δ ΠΠ τ R = 1 − c s + M γ (0)2 + 2 aT β λ M , (5.19) λ Π π τ R = 13 − c s − M γ (2)2 . (5.20)Converting M to the speed of sound and taking m → δ ΠΠ τ R ≈ T dm dT J − , J , ! + (cid:18) − c s (cid:19) + γ (0)2 J , J , − T dm dT ! (cid:18) − c s (cid:19) , (5.21) λ Π π τ R ≈ γ (2)2 J , J , − T dm dT !! (cid:18) − c s (cid:19) , (5.22)where γ (0)2 and γ (2)2 are calculated in Appendix C andare given by Eqs. (C.14) and (C.15), respectively. Wheninserted, one gets the leading orders of the coefficients, δ ΠΠ τ R ≈ π m eq T − π m T ! + (cid:18) − c s (cid:19) (cid:18) π Tm eq + 1 (cid:19) +0 . (cid:18) − c s (cid:19) T m , (5.23) λ Π π τ R ≈ . (cid:18) − c s (cid:19) , (5.24)where the numerical factors come from evaluating g (12 / ≈ .
97 and (1 + 12 g / ≈ .
05 with g and g given by Eqs. (C.16) and (C.17). As seen, the coef-ficient δ ΠΠ /τ R is affected by the soft physics even morestrongly than the bulk viscosity which is manifested bythe factors 1 /m eq and 1 /m .Repeating the same analysis for the Boltzmann gas,which leads simply to replacement of the thermody-namic functions J n,q → I n,q , one obtains the same valueof ζ Boltz /τ R as within the Chapman-Enskog approach,Eq. (4.35). The other two coefficients have the forms(5.21) and (5.22) with γ (0)2 and γ (2)2 given by Eqs. (C.19)and (C.20). The explicit expressions in the m → δ ΠΠ , Boltz τ R ≈ m , c T ! + 5 (cid:18) − c s (cid:19) − . (cid:18) − c s (cid:19) , (5.25) λ Π π, Boltz τ R ≈ . (cid:18) − c s (cid:19) , (5.26)where the numerical factors were found from 144 g c ≈− . g c ) ≈ . g c and g c writtenup below Eq. (C.20). One can also see from Eqs. (5.19)and (5.20) that when thermal quantities are neglectedand the constant mass is kept, we reproduce ζ Boltz /τ R , λ Π π, Boltz /τ R and the two first terms of δ ΠΠ , Boltz /τ R fromRef. [41].
6. SUMMARY AND CONCLUSIONS
In this paper we analyzed the influence of the meanfield on fluid dynamics in weakly interacting systems ofa single species, where all occurring masses are muchsmaller than the system’s temperature. Our main atten-tion was paid to proper determination of the form of thenonequilibrium correction to the distribution functionwhich depends on the mass varying as the temperaturevaries. The correction guarantees a consistent hydro-dynamic description which satisfies thermodynamic rela-tions and the conservation of energy and momentum andfurthermore gives an accurate fixing of the temperaturethrough Landau matching. The correction plays a centralrole in studying thermal dependence of bulk viscous dy-namics. Therefore, we further considered the Anderson-Witting model of the Chapman-Enskog approach andcomputed ζ/τ R of single-component Bose-Einstein andBoltzmann gases. We also derived the evolution equa-tion for the bulk pressure in the 14-moment approxima-tion and obtained relevant transport coefficients. Bothmethods provide the same result for ζ/τ R .The ratio ζ/τ R obtained for the Boltzmann statisticsbehaves as expected, that is, it is given by the noncon-formality parameter squared. When thermal effects areomitted, we reproduce the result from Refs. [41, 43]. Onthe other hand, for very high temperatures the ratio getsdominated by the β λ function. We also see that in spiteof breaking conformal invariance, bulk viscosity vanishesat some critical temperature where c s = 1 /
3. In the caseof the Bose-Einstein gas, we have shown that the leadingorder term of ζ/τ R is different than expected if we neglecteither the physical mass or thermal effects. The ratio inthis case is strongly redounded by the infrared physics,which introduces an additional energy-scale-dependentfactor T /m x . We suspect that the relaxation time ap-proximation used here does not include the entire mi-croscopic physics of a quantum gas, in particular, it isinsensitive to phenomena at the soft scale. Therefore, we3conclude that to compute the bulk viscosity over its re-laxation time for quantum gases of Bose-Einstein statis-tics, one needs to use more advanced methods and solvean integral equation. It can be done starting from eitherthe linearized Boltzmann equation or Kubo formulas, inwhich case note that the formula for the bulk relaxationtime was recently found [58]. ACKNOWLEDGMENTS
It is a pleasure to thank G. Denicol, J. Kapusta, andJ.-F. Paquet for useful discussions. This work is sup-ported in part by the Natural Sciences and EngineeringResearch Council of Canada, and by the US DOE, underContract No. DE-SC0012704. In addition, we gratefullyacknowledge support from the program Mobility Plus ofthe Polish Ministry of Science and Higher Education (A.C. ), from the Canada Council for the Arts through itsKillam Research Fellowship program (C. G.), and fromthe Goldhaber Distinguished Fellowship program fromBrookhaven Science Associates (C. S.).
Appendix A: Components of the energy-momentumtensor correction
The correction to the energy-momentum tensor is∆ T µν = Z dKk µ k ν δf + ∂T µν ∂m ∆ m (A.1)and its particular components are derived as follows. For µ = ν = 0, one gets∆ T = Z dKE k δf + ∆ m Z dKf − ∆ m β Z dKE k f (1 + f ) − ∆ U = Z dKE k ∆ f, (A.2)where the condition on ∆ U , given by (3.13), andEq. (2.13) have been used. Using Eq. (2.23) for ∆ f , onehas∆ T = Z dKE k ∆ f, = Z dKE k (cid:20) δf − T dm dT f (1 + f ) E k R dKδf R dKE k f (1 + f ) . (cid:21) = Z dK E k − T dm dT ! δf (A.3) Analogously, one gets the momentum density variation∆ T i = Z dKE k k i δf − ∆ m β Z dKk i f (1 + f )= Z dKE k k i ∆ f. (A.4)The stress tensor variation is∆ T ij = Z dKk i k j δf − ∆ m β Z dK k i k j E k f (1 + f ) − ∆ m Z dK k i k j E k f + δ ij ∆ U = Z dKk i k j ∆ f − ∆ m Z dK k i k j E k f + δ ij ∆ m Z dKf , (A.5)where Eq. (2.23) has been applied. Equations (A.3) -(A.5) correspond to Eqs. (3.15) - (3.17). Among all theseexpressions, ∆ T ij needs further simplifications to showhow one can obtain Eq. (3.25). The second and the thirdterms of the first line in Eq. (A.5) may be combined toget∆ T ij = Z dKk i k j δf − ∆ m Z dKk i k j ∂ E k (cid:18) f E k (cid:19) + δ ij ∆ U. (A.6)Next, using ∂ E k ( . . . ) = E k k ∂ k ( . . . ) and then integratingby parts leads to∆ T ij = Z dKk i k j δf − δ ij ∆ m Z dKf + δ ij ∆ U = Z dKk i k j δf, (A.7)where the condition (3.13) has been used. Appendix B: Details of the thermodynamic integrals1. Boltzmann statistics
Our strategy to evaluate the integrals with the Boltz-mann statistics is to use the integral representation ofthe modified Bessel functions of the second kind K n ( z ) = Z ∞ dθ cosh( nθ ) exp ( − z cosh θ ) , (B.1)We will also need the Bickley functions defined byKi r ( z ) = Z ∞ dθ exp ( − z cosh θ )(cosh θ ) r . (B.2)4We will need the following series in the small z limit K ( z ) ≈ z − z (cid:0) − γ E + ln 4 − z (cid:1) , (B.3) K ( z ) ≈ z −
12 + z (cid:0) − γ E + 2 ln 4 − z (cid:1) , (B.4) K ( z ) ≈ z − z + z , (B.5) K ( z ) ≈ z − z + 14 , (B.6) K ( z ) ≈ z − z + 1 z , (B.7)where γ E = 0 .
577 is the Euler constant and ln 4 = 1 . z were neglected. For theBickley function [59], we needKi ( z ) ≈ π − z (1 − γ E − ln( z/ | k | ≡ k = m x sinh θ , the thermodynamic func-tions I n,q , defined by (4.34) and evaluated in the fluidrest frame, can be expressed as I n,q ( T, z ) = a q T n +2 z n +2 π × Z ∞ dθ (cosh θ ) n − q (sinh θ ) q +2 exp ( − z cosh θ ) , (B.9)where z = m x T and a q = 1 / ((1 + 2 q )!!). Using cosh x =( e x + e − x ) / x = ( e x − e − x ) /
2, and the defini-tion (B.1), these integrals can be expressed in terms ofmodified Bessel functions of the second kind.Let us consider I , first. After the angle integral, wehave I , = 12 π Z ∞ dk k E k e − E k /T . (B.10)Using k = m x sinh θ , this becomes I , = m x π Z ∞ dθ sinh θ cosh θ e − z cosh θ , (B.11)where z = m x /T . By using cosh x = ( e x + e − x ) / x = ( e x − e − x ) / I , = − m x π (2 K ( z ) − K ( z ) − K ( z )) ≈ m x π (cid:18) z − z (cid:19) = 12 T π (cid:18) − z (cid:19) . (B.12)The other useful integrals are found in a similar way: I , ≈ T π (cid:18) − z (cid:19) , (B.13) I − , ≈ T π (cid:16) − zπ (cid:17) , (B.14) I , ≈ T π (cid:18) − z (cid:19) , (B.15) I , ≈ T π (cid:18) − z (cid:0) − γ + ln 4 − z (cid:1)(cid:19) , (B.16) where I , ≡ q is needed for the thermal mass evalua-tion. For ǫ + P , we have ǫ + P ≈ T π (cid:18) − z z (cid:19) . (B.17)
2. Bose-Einstein statistics
The thermodynamic integrals for the Bose-Einstein gasare defined by Eqs. (4.14): J n,q = a q Z dK ( u µ k µ ) n − q ( − ∆ µν k µ k ν ) q [ f (1 + f )] . (B.18)In the fluid rest frame and after the angle integrals, J n,q becomes J n,q = a q π Z ∞ dk k E k F n,q ( E k ) f ( E k )(1 + f ( E k ))(B.19)and F n,q ( E k ) = E n − qk k q = E n − qk ( E k − m x ) q . (B.20)Using ∂ k f ( E k ) = − kT E k f ( E k )(1+ f ( E k )) and integrat-ing by parts, we can rewrite the above as J n,q = a q T π Z ∞ dk f ( E k ) ∂ k ( kF n,q ( E k )) . (B.21)Changing the integration variable to E k , we further get J n,q = a q T π Z ∞ m x dE k G n,q ( E k ) f ( E k ) , (B.22)where G n,q ( E k ) = ( E k /k ) ∂ k ( kF n,q ( E k ))= ( E k − m x ) q E n − q − k (( n + 1) E k + m x (2 q − n )) p E k − m x (B.23)using k = p E k − m x .Our strategy to evaluate this integral is to separatethe high momentum contribution and the low momentumcontribution. We know how to evaluate Z ∞ m x dE k E lk f ( E k ) = T l +1 Z ∞ z dx x l ∞ X n =1 e − nx (B.24)where x = E k /T and z = m x /T in terms of the poly-logarithmic functions Li n ( z ). Hence, we first expandthe square root in m x /E k and identify the non-negativepower terms in E k . Denoting the collection of such termsas H n,q ( E k ), we then separate the integral as J n,q = a q T π Z ∞ m x dE k H n,q ( E k ) f ( E k )+ a q T π Z ∞ m x dE k ( G n,q ( E k ) − H n,q ( E k )) f ( E k )(B.25)5One can show that the reminder G n,q ( E k ) − H n,q ( E k ) = O (1 /E k ) for all n and q . Then expanding f in the small E k /T limit, f ( E k ) = TE k −
12 + E k T + O (cid:0) ( E k /T ) (cid:1) (B.26)we can keep the first three terms in the integrand tocalculate the soft contribution. This integral can usuallybe exactly evaluated.Let us consider J , . From Eq. (B.22), we have J , = T π Z ∞ m x dE k E k − m x E k p E k − m x f ( E k ) . (B.27)Expanding the square-root in powers of m x /E k , we get G , ( E k ) = 4 E k − m x E k + O (1 /E k ) . (B.28)We can then separate the hard and the soft parts J , = T π Z ∞ m x dE k (cid:0) E k − m x E k (cid:1) f ( E k )+ T π Z ∞ m dE k " E k − m x E k p E k − m x − E k + m x E k f ( E k ) . (B.29)Since the square bracket behaves like 1 /E k , we can use f ( E k ) ≈ T /E k − / E k / T to evaluate the secondintegral. It is J soft n,q ≈ T π (cid:18) z − z
16 + 7 z (cid:19) (B.30)with z = m x /T . The hard part is J hard3 , = T π Z ∞ m x dE k (cid:0) E k − m x E k (cid:1) f ( E k )= T π Z ∞ z dx (cid:0) x − z x (cid:1) e x − T π (cid:16)
24 Li ( e − z ) + 24 z Li ( e − z ) + 11 z Li ( e − z )+3 z Li ( e − z ) (cid:17) = T π (cid:18) π − π z − z z − z O (cid:0) z (cid:1)(cid:19) . (B.31)Adding the two yields J , = J hard3 , + J soft3 , ≈ T π (cid:18) π − π z z (cid:19) . (B.32)This formula works better than 1 part in 10 up to z = m x /T = 1.The usual way of evaluating Bose-Einstein integrals isto use modified Bessel functions of the second kind: J , = T π Z ∞ dk kE k E k − m x E k k ∞ X n =1 e − nE k /T = T m x π ∞ X n =1 (cid:18) K ( nz ) + 12 K ( nz ) (cid:19) . (B.33) Using the small x expansion of K n ( x ) and collecting onlythe terms converging under the infinite sum, we get J Bessel3 , ≈ T m x π ∞ X n =1 (cid:18) n z − n z (cid:19) = T π (cid:18) π − π z (cid:19) . (B.34)which gets only the first two terms.The useful integrals are then found using the formermethod: J , ≈ T (cid:18) − z π + 3 z π (cid:19) , (B.35) J , ≈ T π (cid:18) − z π (cid:19) , (B.36) J − , ≈ π Tz (cid:18) − zπ + z (cid:19) . (B.37)These formulas provide very good approximation up to z = 1. In evaluating J − , one would expect to use theBickley functions equivalently but this method does notwork because the sum cannot be easily evaluated, evenfor the leading behavior.For the enthalpy, one gets ǫ + P ≈ π T (cid:18) − z π (cid:19) (B.38)and for the thermal mass q ≈ T (cid:18) − zπ (cid:19) . (B.39) Appendix C: Irreducible moments
To express the irreducible moments of the distributionfunction one can apply the Grad’s 14-moment approxi-mation, where the correction to the distribution functionof the Bose-Einstein gas is a generalization of the Boltz-mann one, shown in [14, 41, 43], and takes the form δf = f (1 + f ) h E + B m x + D ( u · k ) − B ( u · k ) i Π+ f (1 + f ) B p α p β π αβ . (C.1)The coefficients E , B , D , and B are functions of m x , T , and u µ k µ and they read B = 12 J , , (C.2) D B = − J , J , − J , J , J , J , − J , J , ≡ − C , (C.3) E B = m x + 4 J , J , − J , J , J , J , − J , J , ≡ − C , (C.4) B = − C J , + 3 C J , + 3 J , + 5 J , , (C.5)6where terms related to the particle diffusion have beendropped. Therefore, the irreducible moments ρ − n and ρ µν − n can be expressed by Π and π µν as follows: ρ − n = γ (0) n Π , (C.6) ρ µν − n = γ (2) n π µν , (C.7)where the coefficients γ (0) n and γ (2) n are γ (0) n = ( E + B m x ) J − n, + D J − n, − B J − n, , (C.8) γ (2) n = J − n, J , . (C.9)Only γ (0)2 and γ (2)2 are needed here.Using the prescription shown in Appendix B for eval-uating relevant thermodynamics functions, one finds theleading order terms of coefficients E , B , and D , whichare E ≈ e z T , D ≈ d z T , B ≈ b z T , (C.10) where e = 48 π ( π − ζ (3) ζ (5))5(19 π ζ (3) − ζ (3) − π ζ (5)) ≈ . , (C.11) d = − π ( π ζ (3) − ζ (5))19 π ζ (3) − ζ (3) − π ζ (5) ≈ − . , (C.12) b = 3 π ( π − ζ (3))19 π ζ (3) − ζ (3) − π ζ (5) ≈ − . . (C.13)Therefore, the leading orders of γ (0)2 and γ (2)2 are γ (0)2 ≈ g z T , (C.14) γ (2)2 ≈ g T (C.15)where g = 32( π − ζ (3) ζ (5))5(19 π ζ (3) − ζ (3) − π ζ (5)) ≈ . , (C.16) g = ζ (3)20 ζ (5) ≈ .
06 (C.17)For the Boltzmann statistics, one needs to change inall integrals J n,q → I nq and M → M c . The leading orderresults for the classical gas are E ≈ π T , D ≈ − π T , B ≈ − π T , (C.18)which lead to γ (0)2 ≈ g c T , (C.19) γ (2)2 ≈ g c T , (C.20)where g c = − (5 + 12 γ E −
12 ln 2) ≈ − .
075 and g c =1 / [1] C. Gale, S. Jeon and B. Schenke, Int. J. Mod. Phys. A , 1340011 (2013)[2] U. Heinz and R. Snellings, Ann. Rev. Nucl. Part. Sci. ,123 (2013)[3] R. Baier, P. Romatschke, D. T. Son, A. O. Starinets andM. A. Stephanov, JHEP , 100 (2008)[4] B. Betz, D. Henkel and D. H. Rischke, J. Phys. G ,064029 (2009)[5] P. Romatschke, Int. J. Mod. Phys. E , 1 (2010)[6] P. Kovtun, J. Phys. A , 473001 (2012)[7] S. Jeon and U. Heinz, Int. J. Mod. Phys. E , 1530010(2015)[8] W. Florkowski, M. P. Heller and M. Spalinski, Rept.Prog. Phys. , 046001 (2018)[9] G. S. Denicol, J. Noronha, H. Niemi and D. H. Rischke, Phys. Rev. D , 074019 (2011)[10] T. Koide, E. Nakano and T. Kodama, Phys. Rev. Lett. , 052301 (2009)[11] G. S. Denicol, T. Koide and D. H. Rischke, Phys. Rev.Lett. , 162501 (2010)[12] G. S. Denicol, E. Molnar, H. Niemi and D. H. Rischke,Eur. Phys. J. A , 170 (2012)[13] G. S. Denicol, H. Niemi, E. Molnar and D. H. Rischke,Phys. Rev. D , 114047 (2012), G. S. Denicol, H. Niemi,E. Molnar and D. H. Rischke, Phys. Rev. D , 039902(E) (2015)[14] W. Israel and J. M. Stewart, Annals Phys. , 341(1979)[15] S. Ryu, J.-F. Paquet, C. Shen, G. S. Denicol, B. Schenke,S. Jeon and C. Gale, Phys. Rev. Lett. , 132301 (2015) [16] S. Ryu, J. F. Paquet, C. Shen, G. Denicol, B. Schenke,S. Jeon and C. Gale, Phys. Rev. C , 034910 (2018)[17] J. F. Paquet, C. Shen, G. S. Denicol, M. Luzum,B. Schenke, S. Jeon and C. Gale, Phys. Rev. C , 044906(2016)[18] J. F. Paquet, C. Shen, G. Denicol, S. Jeon and C. Gale,Nucl. Phys. A , 429 (2017).[19] J. N. Guenther, R. Bellwied, S. Bors´anyi, Z. Fodor,S. D. Katz, A. P´asztor, C. Ratti and K. K. Szab´o, Nucl.Phys. A , 720 (2017)[20] S. Borsanyi, EPJ Web Conf. , 01006 (2017)[21] P. Romatschke and D. T. Son, Phys. Rev. D , 065021(2009)[22] F. Karsch, D. Kharzeev and K. Tuchin, Phys. Lett. B , 217 (2008)[23] Y. Hama, R. P. G. Andrade, F. Grassi, O. Socolowski,Jr., T. Kodama, B. Tavares and S. S. Padula, Nucl. Phys.A , 169 (2006)[24] J. Noronha-Hostler, J. Noronha and C. Greiner, Phys.Rev. Lett. , 172302 (2009)[25] G. S. Denicol, T. Kodama, T. Koide and P. Mota, Phys.Rev. C , 064901 (2009)[26] A. Nakamura and S. Sakai, Phys. Rev. Lett. , 072305(2005)[27] L. P. Csernai, J. I. Kapusta and L. D. McLerran, Phys.Rev. Lett. , 152303 (2006)[28] G. S. Denicol, T. Kodama and T. Koide, J. Phys. G ,094040 (2010)[29] H. Song and U. W. Heinz, J. Phys. G , 064033 (2009)[30] P. Bozek, Phys. Rev. C , 054909 (2017)[31] A. Monnai, S. Mukherjee and Y. Yin, Phys. Rev. C ,034902 (2017)[32] Y. Akamatsu, A. Mazeliauskas and D. Teaney, Phys. Rev.C , 024902 (2018)[33] S. Jeon, Phys. Rev. D , 3591 (1995)[34] S. Jeon and L. G. Yaffe, Phys. Rev. D , 5799 (1996)[35] P. B. Arnold, G. D. Moore and L. G. Yaffe, JHEP ,001 (2000)[36] P. B. Arnold, G. D. Moore and L. G. Yaffe, JHEP ,030 (2003)[37] P. B. Arnold, G. D. Moore and L. G. Yaffe, JHEP ,051 (2003) [38] S. Chapman and T. G. Cowling, The Mathematical The-ory of Non-uniform Gases (Cambridge University Press,Cambridge, 1970), 3rd ed.[39] H. Grad, Comm. Pure Appl. Math , 331, (1949)[40] P. B. Arnold, C. Dogan and G. D. Moore, Phys. Rev. D , 085021 (2006)[41] G. S. Denicol, S. Jeon and C. Gale, Phys. Rev. C ,024912 (2014)[42] A. Jaiswal, R. Ryblewski and M. Strickland, Phys. Rev.C , 044908 (2014)[43] W. Florkowski, A. Jaiswal, E. Maksymiuk, R. Ryblewskiand M. Strickland, Phys. Rev. C , 054907 (2015)[44] X. G. Huang, T. Kodama, T. Koide and D. H. Rischke,Phys. Rev. C , 024906 (2011)[45] C. Sasaki and K. Redlich, Phys. Rev. C , 055207 (2009)[46] P. Chakraborty and J. I. Kapusta, Phys. Rev. C ,014906 (2011)[47] M. Bluhm, B. Kampfer and K. Redlich, Phys. Rev. C ,025201 (2011)[48] K. Dusling and T. Sch¨afer, Phys. Rev. C , 044909(2012)[49] P. Romatschke, Phys. Rev. D , 065012 (2012)[50] M. Albright and J. I. Kapusta, Phys. Rev. C , 014903(2016)[51] P. Chakraborty and J. I. Kapusta, Phys. Rev. C ,014907 (2017)[52] L. Tinti, A. Jaiswal and R. Ryblewski, Phys. Rev. D ,054007 (2017)[53] M. Alqahtani, M. Nopoush, R. Ryblewski and M. Strick-land, Phys. Rev. Lett. , 042301 (2017)[54] C. Shen, U. W. Heinz, J. F. Paquet, I. Kozlov andC. Gale, Phys. Rev. C , 024908 (2015)[55] C. Shen, J. F. Paquet, U. Heinz and C. Gale, Phys. Rev.C , 014908 (2015)[56] S. Hauksson, C. Shen, S. Jeon and C. Gale, Nucl. Part.Phys. Proc. , 169 (2017)[57] J. L. Anderson, H. R. Witting, Physica , 466-488(1974)[58] A. Czajka and S. Jeon, Phys. Rev. C , 064906 (2017)[59] J. M. Blair, C. A. Edwards and J. H. Johnson, Math.Comput.32