A 14 Moments approximation

# Effect of bulk viscosity on Elliptic Flow near QCD phase transition

## Abstract

In this work, we examine the effect of bulk viscosity on elliptic flow taking into account the critical behavior of the EoS and transport coefficients near the QCD phase transition. We found that the dependence of is quantitatively changed by the presence of the QCD phase transition. Within reasonable values of the transport coefficients, decreases by a factor of at low ( GeV). However, for larger values of ( GeV), the interplay between the velocity of sound and transport coefficient near the QCD phase transition enhances . We further point out that Grad’s 14 moments approximation cannot be applied for the calculation of the one-particle distribution function at the freeze-out.

###### pacs:
47.10.-g,25.75.-q

## I Introduction

Among many novel discoveries in the studies of relativistic heavy-ion physics, one of the remarkable facts is the presence of strong collective flows, which are well described using the ideal hydrodynamic model (1). In this sense, hydrodynamics is expected to be an important tool to investigate the properties of the QCD matter. If such a scenario is established, we may further use the hydrodynamic model to infer the initial condition just after the collision from the final state observables. However, there are still several open questions in the hydrodynamic description of relativistic heavy ion collisions (2) and more extensive studies are necessary. The effect of dissipation is one of such problems.

The effect of shear viscosity in heavy ion collisions has been investigated by several authors (3); (4); (5); (6); (7); (9); (8), and it was found that the elliptic flow parameter can be significantly affected. On the other hand, the effect of bulk viscosity has not been analyzed in detail yet. Roughly speaking, while the shear viscosity acts as the resistance against the deformation of a fluid element, the bulk viscosity acts against the expansion or compression of the fluid. Thus, one may naturally expect that shear viscosity affects directly , which characterizes the spatial anisotropy of the dynamics. However, this does not mean that the bulk viscosity plays a secondary role. Because of the anisotropic expansion and compression of the produced matter, the bulk viscosity should alter the anisotropic flow as well.

Another important factor is the behavior of the bulk viscosity coefficient, , near the QCD phase transition temperature . Recently, lattice QCD simulations (10); (11) suggest that is strongly enhanced near , although the shear viscosity coefficient is suppressed (13). A similar enhancement will appear even below , as is shown in (14), where in the hadronic phase increases also quickly towards . Considering the fact that the fluid dynamics slows down near due to the small sound velocity , the critical behavior of discussed above may play a crucial role on .

As far as we know, the effect of the bulk viscosity on was calculated in (15) where the effect of the phase transition is also considered. As a matter of fact, a quantitative estimate of the effect of the phase transition in dissipative hydrodynamics is not straightforward because the above critical behavior of the bulk viscosity should be consistently constructed with the equation of state (EoS).

As a reliable EoS, the first choice would be the results of the lattice QCD calculation for the QGP phase and of the hadron resonance gas for the hadron phase. However, the connection of these two EoS’s is not trivial because no reliable information is available near . In particular, near is very sensitive to the way that the two phases are connected and it is important to check how these ambiguities influence .

The purpose of this paper is to investigate the effect of bulk viscosity on the elliptic flow within a 2+1 dimensional hydrodynamic modeling, incorporating the critical behaviors from and near . As was pointed out already (5); (8); (3), is affected also by the dissipative corrections on the one-particle distribution function at the freeze-out. However, we will not consider these corrections in this work because the 14 moments approximation, which is the method commonly employed to evaluate this correction, is inapplicable for the freeze-out values of bulk viscosity, as we will discuss later.

This paper is organized as follows. In Section II, we briefly review the equations of causal dissipative hydrodynamics in the hyperbolic coordinate system and present the SPH (smoothed particle hydrodynamics) parameterization for the numerical implementation. In Section III, we set up the parameters of our calculation, such as the initial condition, the EoS and the transport coefficients. We futher discuss the freeze-out process and show the intrinsic problem of the 14 moments approximation. The numerical results for the effects of bulk viscosity on are given in Section IV. Concluding remarks are given in Section V.

## Ii Relativistic Dissipative Hydrodynamics

### ii.1 Memory Function Method in a Hyperbolic Coordinate System

As is well known, the application of dissipative hydrodynamics in the relativistic regime requires some caution. Recently, it became clear that the instability of the hydrodynamic model can be induced by the violation of relativistic causality (16). Thus, the relativistic Navier-Stokes (Landau-Lifshitz) theory (17) is essentially unstable and not applicable for the complex systems generated in the relativistic heavy ion collisions. Several formulations of relativistic dissipative hydrodynamics have been proposed so far to solve the problem of acausality and instability (20).

In this work, we will use the memory function method (18); (19). This theory is formulated in such a way that the magnitude of the bulk viscosity has naturally a lower bound and we can carry out stable numerical simulations even for ultra-relativistic initial conditions. Here, we present the basic equations of this formalism in the hyperbolic coordinate system.

For simplicity, we consider the case of vanishing baryon chemical potential where only the conservation of the energy and momentum is required. For a general metric , we have

 1√−g∂μ(√−gTμν)+ΓνλμTλμ=0, (1)

where is the Christoffel symbol and is the determinant of .

We use the Landau definition for the local rest frame and assume, as usual, that the thermodynamic relations are valid in this frame. In this work, we further ignore the shear viscosity. Then, the energy-momentum tensor is expressed as

 Tμν=(ε+p+Π)uμuν−(p+Π)gμν, (2)

where, , , and are the energy density, pressure, four velocity and bulk viscosity, respectively.

As was discussed in (19), one of the important factors to obtain a stable causal dissipative hydrodynamics is to consider the deformation of the fluid element in the hydrodynamic flow. We introduce a reference density defined by

 ˙σσ=θ=uμ;μ,

where is the expansion rate of a fluid element and denotes the covariante derivative. The above relation can also be expressed in the form of a continuity equation for ,

 1√−g∂μ(√−gσuμ)=0. (3)

The equation for the entropy production is written as,

 1√−g∂μ(√−gsuμ)=−ΠT1√−g∂μ(√−guμ), (4)

where is the entropy density. From this, we can define the thermodynamic force associated with the bulk viscosity as . In the memory function method, an irreversible current is induced by the thermodynamic force as

 J(τ)σ(τ)=−∫ττ0dτ′1τR(τ′)exp(−∫ττ′dτ′′τR(τ′′))F(τ′)σ(τ′)+J(τ0)σ(τ0)exp(−∫ττ0dτ′τR(τ′)). (5)

Thus, the equation for the bulk viscosity is obtained by setting and , and can be expressed in a differential form as follows

 τRuμ∂μ(Πσ)+Πσ=−ζσ1√−g∂μ(√−guμ), (6)

where is the relaxation time.

In hyperbolic coordinates Eqs. (1), (4) and (6) are explicitly given by

 γddτ((ε+p+Π)σui) = 1σ∂i(p+Π), (7) γddτ(sσ) = −ΠTσ(∂μuμ−γτ), (8) τRγddτΠ+Π = −(ζ+τRΠ)(∂μuμ−γτ). (9)

In this case ,where

 τ=√t2−z2,   η=12tanh(t+zt−z),   rT=(x,y), (10)

and . The three equations above are solved in the following calculations.

### ii.2 Smoothed Particles Hydrodynamics

To solve numerically the relativistic hydrodynamic equations, we use the Smoothed Particle Hydrodynamic (SPH) method (21). See (16); (19); (22) for the application to the heavy ion dynamics. For the sake of book keeping, here we just write down explicitly the SPH equations in hyperbolic coordinates leaving the detailed derivation to the references above.

In the SPH method, all the hydrodynamic variables are expressed in terms of a discrete set of Lagrangian coordinates together with a normalized kernel function for the smoothing procedure. The parameter represents the width of and serves as a cut-off parameter for short wavelength modes. For hyperbolic coordinates, this kernel function takes the form

 W(~r,h)=W(η,hη)W(rT,hT), (11)

with

 ∫W(η,hη)dη = 1, (12) ∫W(rT,hT)drT = 1, (13)

where the two width parameters, and , are introduced separately for and , respectively.

The basic procedure of the SPH scheme is to express a reference density and its current , where is the Lorentz factor, as

 τγσ→σ∗SPH(r,t)=NSPH∑α=1ναW(r−rα(t),h), (14) j→jSPH(r,t)=NSPH∑α=1ναdrα(t)dtW(r−rα(t),h), (15)

so that the continuity equation (3) is automatically satisfied (22). Eqs. (14) and (15) can be interpreted as if there are particles associated with each Lagrangian coordinate , with velocity , carrying a quantity for the reference density . These particles are called SPH particles. The following calculation is independent of the choice of and, hence, we set .

The dynamic variables in the SPH scheme are then,

 {rα,uα,(sσ)α,(Πσ)α; α=1,..,NSPH}, (16)

where they represent, respectively, position, velocity, entropy, and bulk viscosity associated with the -th SPH particle. The time evolution of these variables can be obtained from Eqs. (7), (8) and (9) by expressing them in the SPH representation. Writing the entropy density and the bulk viscosity as

 τγs → s∗SPH(r,t)=∑ανα(sσ)αW(r−rα(t),h), (17) Π → ΠSPH=∑ανα1γατ(Πσ)αW(r−rα(t),h), (18)

we obtain the equation for the bulk viscosity of the -th SPH particle

 τRαγαddτ(Πσ)α+(Πσ)α=−ζασα(∂μuμ−γτ)α, (19)

and the equation of motion Eq. (7) (22); (16)

 σ∗SPHddτ((ε+p+Π)ασαui α)=τNSPH∑β=1νβσ∗α⎛⎜ ⎜⎝pβ+Πβ(σ∗β)2+pα+Πα(σ∗α)2⎞⎟ ⎟⎠ ∂iW(|rα−rβ(t)|). (20)

The r.h.s. of Eq.(20) is nothing but the SPH representation of the gradients of pressure and bulk viscosity. We remark that the l.h.s. of Eq.(20) contains an implicit velocity dependence through the Lorentz factor in the thermodynamic variables. After some straightforward manipulations, we get

 Mmidumdτ=Fi, (21)

with

 Mji = γCδmi−Aγglmuiul, (22) Fi = Bui+∂i(p+Π), (23)

and

 A = ε+p−dwds(ΠT+s)−ζτR, (24) B = (γσ∗dσ∗dτ−γτ+12γulumdglmdτ)A+ΠτR, (25) C = ε+p+Π. (26)

In the above equations, we calculated the four divergence of the velocity as

 ∂μuμ=−γασ∗dσ∗dτ−gijuiγdujdτ−uiuj2γdgijdτ,

and, from the continuity equation,

 dσ∗dt=1σ∗ jSPH⋅∇σ∗SPH−∇⋅jSPH.

As mentioned in the introduction, we study the transverse dynamics near the central rapidity region, assuming the Bjorken scaling behavior in , that is, . Then the dynamics is restricted on the transverse plane and is simplified as

 (γCδij+AγuiTujT)dujTdτ=BuiT−∂i(p+Π), (27)

with

 B=(γσ∗dσ∗dτ−γτ)A+ΠτR. (28)

These are the SPH representations of our equations which will be solved.

## Iii Parameter Settings of the Model

### iii.1 Initial Conditions

One of the central issues of hydrodynamical simulations of relativistic heavy-ion collisions is to infer what are the real initial conditions attained in the instant of the collision. Here we just assume the initial conditions suggested by (23); (5) for the sake of comparison. The initial energy density is parameterized as

 ε0 = K⎧⎪ ⎪ ⎪⎨⎪ ⎪ ⎪⎩TA(x+b2,y)⎡⎢ ⎢ ⎢⎣1−⎛⎜ ⎜⎝1−σNNTB(x−b2,y)B⎞⎟ ⎟⎠B⎤⎥ ⎥ ⎥⎦ (29) +TB(x−b2,y)⎡⎢ ⎢ ⎢⎣1−⎛⎜ ⎜⎝1−σNNTA(x+b2,y)A⎞⎟ ⎟⎠A⎤⎥ ⎥ ⎥⎦⎫⎪ ⎪ ⎪⎬⎪ ⎪ ⎪⎭,

where for Au-Au collisions with the nucleon-nucleon cross section . The free parameter is fixed so that the maximum energy density is GeV for central collisions, as is done in (5). The term is the nuclear thickness function defined as

 TA(x,y)=∫dzρ01+exp[(r−RA)/λ], (30)

where , and , respectively. The impact parameter is chosen as .

### iii.2 Equation of State

We consider three different EoS’s; one for an ideal gas of massless quarks and the other two corresponding to different connections between the QGP and hadron gas phases, as shown in Fig. 1. The dash-dotted line corresponds to the ideal massless quark gas, which we refer to as EoS I. The dotted line (EoS II) and solid line (EoS III) are both constructed by connecting the EoS from lattice QCD calculations with three flavors for the QGP phase (11) and that of the hadron resonance gas for the hadron phase. However, as was mentioned in the introduction, it is not well established how the lattice EoS is connected to that of the hadron resonance gas EoS. In EoS II the two phases are connected smoothly for a wide interpolation domain of the temperature, MeV, whereas in EoS III the connetion is done within a narrower interpolation domain, MeV (see (12) for the connection method). Thus, compared with the EoS II, the EoS III is more faithful, both to the lattice QCD and the hadronic resonance gas values except in a narrow transition region near MeV. On the other hand, since the connection is rather abrupt, a kink-like structure appears at the QCD phase transition.

The temperature dependence of the pressure and are shown in Figs. 1 and 2, respectively. One can see that, although the behavior of the pressure is not drastically modified, is still very sensitive to the choice of connection; in EoS III has a narrow minimum reaching a much smaller value than that of EoS II. The nature of the transition of EoS III is almost that of the first order phase transition. If the transition were of first order, would vanish at the phase transition. Since the acceleration of the fluid is directly related to , this difference would affect significantly the hydrodynamic evolution of the system near the phase transition. Thus, we perform a comparative study of the behavior of calculated using these equations of state.

### iii.3 Transport Coefficients

The behavior of the transport coefficients and has not yet been established. There are several approaches to calculate transport coefficients. A well-known method is to use the Green-Kubo-Nakano (GKN) formula. Exactly speaking, the validity of the GKN formula is, however, not clear in the relativistic regime since it was developed for Newtonian fluids. Recently it has been shown that even in the presence of memory effects, the viscosity coefficients coincide with those of the GKN formula, at least, in the leading order approximation (24). Here, we use the result from lattice QCD (10), calculated with the GKN formula, for the QGP phase. For the hadron phase, we use the result obtained recently by (14). Then, the bulk viscosity coefficient shows a maximum at and starts to decrease almost exponentially as the system departs from this critical region. The coefficient vanishes in the high temperature limit in the QGP phase whereas in the low temperature limit (hadronic gas), it converges to a finite value. We fitted these results by analytic expressions and connected them by a parabola around MeV. The result is

 ζs=⎧⎪⎨⎪⎩A1x2+A2x−A3(0.995Tc≥T≥1.05Tc)λ1exp(−(x−1)/σ1)+λ2exp(−(x−1)/σ2)+0.001(T>1.05TC)λ3exp((x−1)/σ3)+λ4exp((x−1)/σ4)+0.03(T<0.995TC),

where . The fitted parameters are , , , , , , and . The Fig. 3 shows this fit and compares it with the results from (10); (14).

Inspired by the result for the relaxation time of the causal shear viscosity coefficient (24), we use the following parameterization for the relaxation time of the bulk viscosity,

 τR=ζϵ+pbΠ. (31)

In principle, the parameter is a function of thermodynamical quantities constrained by the causality condition (19). Here, for simplicity, we will always consider as a constant. Due to this assumption, the relaxation time also shows critical behavior as in the case of the bulk viscosity coefficient, as is shown in Fig. 4.

### iii.4 Freeze out and Kinetic Corrections

From extensive studies of the shear viscosity, it is known that there are two different origins for the dissipative corrections to the elliptic flow. One is the effect of viscosity on the collective evolution of the fluid itself. The other is the correction to the one-particle distribution function at the freeze-out hyper surface. It was shown that the latter correction is more important than the former for the shear viscosity (5); (8).

The correction to the one-particle distribution function is usually estimated using Grad’s 14 moments approximation. However, the applicability of this approach to the heavy-ion scenario is not so obvious. In the following, we will show that this procedure exhibits intrinsic problems at least for the bulk viscosity.

In the 14 moments approximation, the one-particle distribution function for bosons is expressed as

 f=f0+δf≈f0+f0(1+f0)(E0+D0(pμuμ)+B0(p2−4(pμuμ)2))Π, (32)

where refers to the Bose-Einstein distribution function and is the four momentum of a particle. In this expression, we neglected the contributions from the shear viscosity and the heat conduction. See Appendix A for details 1.

In Figs. 5 and 6, the coefficients of the viscous corrections , and are shown as a function of . For the sake of comparison, we plot also in Fig. 5, which is the coefficient corresponding to the viscous correction from the shear viscosity (see Eq. (66)). One can see that the bulk correction is one order of magnitude larger than the shear correction . By using these coefficients with MeV and MeV, we calculated the correction to the one-particle distribution function for different values of the bulk viscosity, which are shown in Fig. 7. In our simulations, the typical value of the bulk viscosity at the freeze-out hyper surface reaches fm, which means that is much bigger than itself. Then, the corrected one-particle distribution function can be even negative. The negative corrections come from the normalization conditions, Eqs. (57) and (58), which imply that

 ∫d3→K(2π)3δf = 0, (33) ∫d3→K(2π)3Eδf = 0. (34)

Obviously, the 14 moments approximation is consistent only when is much smaller than the equilibrium distribution function . That is, the magnitude of bulk viscosity should be, at least, smaller than fm at the freeze-out for pions. This value of bulk viscosity is very small and corresponds only to of the value of pressure of a hadron resonance gas at MeV. Usually, the procedure of freeze-out in hydrodynamics does not impose any constraint on the magnitude of the bulk viscosity and it is not clear that such small values of viscosity will be achieved. For example, in our simulations, the value of the ratio between the bulk viscosity and the pressure is around , which seems to be a small correction but is already too large for the 14 moments approximation. This situation is unchanged even for heavier particles ().

One may notice that Grad’s method was originally formulated for a single component fluid and, to discuss the appropriate correction to the one-particle distribution function, we have to generalize it to a multiple component fluid. Thus, one might argue that if we consider the multi-component fluid the limitation discussed above might be relaxed. In this case, each hadronic species can have its own bulk viscosity and pressure, and , where

 Π = ∑iΠ(i), (35) p = ∑ip(i). (36)

Here the index  represents the different hadronic species. In the context of the first order theory, the ratio should be proportional to the ratio of the coefficients , that is,

 Π(i)Π=ζ(i)ζ=s(i)s∼10−1. (37)

In the above expression,  is the entropy density and  is the entropy density contribution from the -th specie. Then one should use the ratio for each species, , to estimate the dissipative correction. From the discussion above, the magnitude of the bulk viscosity for pions, , should be limited to be smaller than fm. From this, we conclude that at MeV should be at maximum , instead of , for the moment method to be applicable. Even though the limit of applicability is relaxed, it is still not clear that such a constraint will always be satisfied in a realistic hydrodynamic scenario. Furthermore, it should be noted that the generalization of the 14 moments method to a multiple component fluid is not trivial, neither unique. There are several different solutions proposed (26); (27).

We conclude that the applicability of the 14 moments approximation is not clear in heavy-ion collisions. For this reason, we will not consider this correction to the single-particle distribution function in the following discussions.

## Iv Numerical results

We calculate the elliptic flow, , for pions as a function of the transverse momentum with a freeze out temperature of MeV by the Cooper-Frye method (22). Unless we say otherwise, we take . Results for constant values of the ratio / will be shown for the sake of comparison.

In Fig. 8, we show calculated with the ideal gas of massless quarks EoS (EoS I). Analogous to the case of the shear viscosity, decreases as increases, which is expected for smooth expanding cases, since the viscosity converts a part of the collective kinetic energy to heat. Even when we use the / with the critical behavior of the QCD phase transition, this surppression of is essentially the same, as shown by the dotted line. As a matter of fact, we found a direct correlation between and the entropy production. In Fig. 9, the entropy production is shown for different values of /, and one can notice that the larger the entropy production, the bigger the suppression of .

Next, we consider the effect of the QCD phase transition through the EoS. In Fig. 10, of viscous (dashed, dash-dotted and dotted lines) and ideal (solid line) fluids are shown for the EoS II. The dashed and dash-dotted lines correspond to constant / cases, and the dotted line represents the with the QCD phase transition. Although the general behavior is similar to the case of EoS I, the curvatures of for the corresponding cases are changed: is reduced at low but starts to increase at high . This effect at high is enhanced when exhibits the critical behavior.

To see this effect in detail, we show up to 3 GeV in Fig. 11. For GeV, starts to increase and eventually surpasses the value of the ideal case. No significant difference was observed for the EoS III compared to the case of EoS II, although this EoS displays a sharper transition from the QGP phase to the hadron phase. As discussed below, this interesting feature is the result of an interplay between the viscosity and the sound velocity.

First, note that the enhancement of at high is not observed for the EoS I case. Thus, this effect should be attributed to the presence of the QCD phase transition (or cross over). In Fig. 12, the temperature distribution calculated with EoS II at fm is shown. The solid and dashed lines represent the ideal and viscous results, respectively. One can easily see that, for the viscous case, the behavior of the fluid flow is not monotone in the region where has its minimum. This is because the flow of the internal matter, which has a higher temperature, tends to catch up the foregoing fluid elements, generating a non-monotonous velocity field in the radial direction. This behavior of the velocity field causes the bulk viscosity to become piled up, which heats the matter in this region. If this happens, the gradient of becomes dominant in the acceleration of the fluid compared to the pressure gradient (note that the acceleration is given by ) since the gradient of the pressure is proportional to . Such a mechanism works more efficiently in the direction where the initial acceleration is large, and in consequence, in the direction to increase the elliptic flow. Furthermore, this effect becomes effective only when significant collective flow is formed near the phase transition. Thus, the recovery of by this mechanism is expected for higher momentum particles. We consider that this is the reason for the counter-intuitive behavior of our result. Of course one should keep in mind that the applicability of the hydrodynamic approach becomes dubious for the description of such a high dynamics so that such a mechanism may not necessarily be attributed directly to experimentally observed values.

So far, all the calculations presented were implemented fixing the relaxation time parameter . This is because the dependence for the EoS II is negligibly small. However, the dependence can be observed for an EoS with a sharper phase transition such as EoS III. In Fig. 13, we show the elliptic flow calculated with the EoS III for different relaxation times, and . One can clearly see a dependence on the parameter at high GeV although the curves of are independent of the values of for small .

## V Concluding remarks

In this work, we examined the effect of bulk viscosity on taking into account the critical behavior of the EoS and transport coefficients near the QCD phase transition. We found that the dependence of is quantitatively changed by the presence of the QCD phase transition. Within reasonable values of the transport coefficients, decreases by a factor of at low ( GeV). However, for larger values of ( GeV), the interplay between and near the QCD phase transition generates more collective flow, with even larger values than those of the ideal fluid. These effects can depend on the value of the transport coefficient when the phase transition is more abrupt. Of course one should keep in mind that the applicability of the hydrodynamic approach is not guaranteed for the high dynamics, but this is an interesting feature and deserves further investigation.

All the calculations of here were done without dissipative corrections to the one-particle distribution function. As discussed in this paper, the naive application of the usual 14 moment calculation for the bulk viscosity leads to unacceptable values for the corrections to the distribution function. It should be noted that the 14 moments approximation is not an unique choice to solve Grad’s method and there are other approaches to calculate the non-equilibrium distribution function (25). It may be possible that the problem presented here can be solved using alternative approaches. Further investigation on this aspect is required urgently.

Although not explored in the present paper, one interesting aspect observed in these studies should be mentioned. For much larger values of , say , and smaller (), we found some peculiar behavior of at large , exhibiting large fluctuations for small changes in the transport coefficients. We identified that such fluctuations are related to the emergence of instabilities in the fluid evolution and that they are amplified when a sharp phase transition (such as a first order phase transition) is present. We also found that these instabilities cannot be controlled even after introducing the additional viscosity discussed in (28). Although the shear viscosity was not considered here, we believe that these instabilities are an intrinsic phenomenon that appears for large bulk viscosity dynamics, such as the Burgers instabilities (16); (29); (30). If this is true, they should reflect in the observed . For instance, the event-by-event fluctuations in would become very large at high . Since this effect is closely related to the nature of the phase transition, it would be interesting to verify this experimentally.

We acknowledge illuminating discussions from T. Hirano and A. Monnai and thank J. Noronha for constructive comments on the manuscript. This work has been supported by CNPq, FAPERJ, CAPES, PRONEX and the Helmholtz International Center for FAIR within the framework of the LOEWE program (Landesoffensive zur Entwicklung Wissenschaftlich- Okonomischer Exzellenz) launched by the State of Hessen.

## Appendix A 14 Moments approximation

In this appendix, we discuss the Grad’s method with the 14 moments approximation (31). In this approximation, the non-equilibrium one-particle distribution function is assumed as

 f = (exp(−y)+a)−1, (38) y = (ϵ+α0)+Kμ(ϵμ−β0uμ)+KμKνϵμν, (39)

where the parameters , and represent the deviation from equilibrium up to second order in momentum and is for bosons, for fermions and for a Boltzmann gas. The index 0 is applied for the quantities in equilibrium. Without loss of generality is assumed to be traceless,

 ϵμμ=0. (40)

For small momentum, the non-equilibrium one-particle distribution function can be expanded around the equilibrium state

 f≈f0+f0(1−af0)(y−y0)+O(2), (41)

where

 f0 = (exp(α0(x,t)−β0(x,t)Kμuμ(x,t))+a)−1, (42) y−y0 = ϵ+Kμϵμ+KμKνϵμν≪1. (43)

For a rarefied gas, the conserved number current and the energy-momentum tensor are expressed as

 Nμ =∫d3→K(2π)3EKμf(K,T,μ), (44) Tμν =∫d3→K(2π)3EKμKνf(K,T,μ). (45)

Substituting Eq. (41) into them, we have

 Nμ =Iμ0+ϵJμ0+Jμν0ϵν+Jμνλ0ϵνλ, (46) Tμν =Iμν0+ϵJμν0+Jμνλ0ϵλ+Jμνλρ0ϵλρ, (47)

where,

 Iα1…αn0 ≡∫d3→K(2π)3EKα1…Kαnf0, Jα1…αn0 ≡∫d3→K(2π)3EKα1…Kαnf0(1−af0). (48)

The 0-th order terms are identified as the equilibrium currents,

 Iμ0 =Nμ0, (49) Iμν0 =Tμν0. (50)

All the remaining terms come from the dissipative corrections. The moments (48) can be expanded as

 Iμ0 =I10uμ, Iμν0 =I20uμuν−I21Δμν, Jμ0 =J10uμ, Jμν0 =J20uμuν−J21Δμν, Jμνλ0 =J30uμuνuλ−J31uμΔνλ−J31uλΔνμ−J31uνΔλμ Jμνλρ0 =J40uμuνuλuρ−J41uμuνΔλρ−J41uμuλΔνρ−J41uμuρΔλν−J41uνuλΔμρ −J41uνuρΔμλ−J41uλuρΔμν+J42ΔμνΔλρ+J42ΔμλΔνρ+J42ΔμρΔνλ, (51)

where

 Inq =1(2q+1)!!∫d3→K(2π)3EEn−2q(k2)qf0, (52) Jnq =1(2q+1)!!∫d3→K(2π)3EEn−2q(k2)qf0(1−af0). (53)

In hydrodynamics, the dissipative corrections in the currents are parameterized in terms of the bulk viscosity, the shear viscosity and the heat conductivity. The variables , and can be determined in terms of these hydrodynamic variables with the following set of constraints

 qν=ΔμνNμ, (54) πμν=Pμναβ(Tαβ−Tαβ0), (55) Π=−13Δμν(Tμν−Tμν0), (56) uμ(Nμ−Nμ0)=0, (57) uν(Tμν−Tμν0)=0. (58)

The first three relations define the irreversible currents while the last two restrictions come from the Landau choice for the local equilibrium reference frame. If we were using the Eckart frame, for example, we would have to use, instead of (57) and (58),

 (Nμ−Nμ0)=0, (59) uμuν(Tμν−Tμν0)=0. (60)

The conditions (54), (55) and (56) imply,

 qμ =−J21Δμνϵν−2J31Δμνuλϵνλ, (61) πμν =2J42Pμνλρϵλρ, (62) Π =ϵJ21+J31uλϵλ+J41uλuρϵλρ−53J42Δλρϵλρ. (63)

A general solution can be constructed using the following ansatz,

 ϵ =E0Π, (64) ϵλ =D0Πuλ+D1qλ, (65) ϵλρ =B0(Δλρ−3uλuρ)Π+B1(uλqρ+uρqλ)+B2πλρ. (66)

In order to know the form of the one-particle distribution function (41), we have to determine the form of the coefficients , , , , and . This can be done by substituting these general solutions into the equations discussed above,

 (J21+J31)D1+(J31+J41)B1 = −1, (67) B2 = 12J42, (68) J21E0+J31D0+(3J41+5J42)B0 = 1, (69) J10E0+J20D0+(J30+J31)3B0 = 0, (70) J31D1+J41B1 = 0, (71) J20E0+J30D0+(J40+J41)3B0 = 0. (72)

The first three equations, (67), (68) and (69), come directly from (61), (62) and (63), respectively. The last three equations, (70), (71), and (72), are consequences of Eqs. (57) and (58). The solution of this set of equations is

 B0 = 1−3C1J21−3C2J31−3J41−5J42, (73) B1 = J31J41J21−J31J31, (74) B2 = 12J41, (75) D03B0 = −4J31J20−J41J10J30J10−J20J20≡−C2, (76) D1 = −J41J41J21−J31J31, (77) E03B0 = m2+4J31J30−J41J20J30J10−J20J20≡−C1, (78)

where is the mass of the particle.

In this paper, we are interested only in the bulk viscosity. Then, the form of the distribution function is

 f=f0+f0(1−af0)(E0+D0(pμuμ)+B0(p2−4(pμuμ)2))Π, (79)

We remark that the 14 moments approximation was formulated for a single component gas. The generalization of this method for multiple components systems has not yet been done.

### Footnotes

1. The correction to the one-particle distribution function due to the bulk viscosity is also discussed in (8). Their expression, however, does not include the scalar and linear terms in momentum to the correction of the distribution function, that is, they neglect terms like and . However, this form of distribution function breaks the Landau conditions and satisfies only one of the normalization conditions, Eqs. (33) and (34).

### References

1. See for example, P. Huovinen and P.V. Ruuskanen, Ann. Rev. Nucl. Part. Sci. 56, 163 (2006); Y. Hama, T. Kodama and O. Socolowski, Braz. J. Phys. 35, 24 (2005), Jean-Yves Ollitrault, Euro. J. Phys. 29, 275 (2008) and references therein.
2. T. Kodama, T. Koide, G. Denicol and Ph. Mota, Int. J. Mod. Phys. E16, 763-776 (2007).
3. P. Romatschke and U. Romatschke, Phys. Rev. Lett. 99, 172301 (2007); R. Baier, P. Romatschke and U. Wiedemann, Nucl. Phys. A782, 313 (2007); R. Baier and P. Romatschke, Eur. Phys. J. C51, 677-687 (2007).
4. U. Heinz, H. Song and A. Chaudhuri, Phys. Rev. C73, 034904 (2006).
5. H. Song and U. Heinz, Phys. Rev. C77, 064901 (2008).
6. H. Song and U. Heinz, Phys. Rev. C78, 024902 (2008).
7. A. Chaudhuri, J. Phys. G35, 104015 (2008).
8. K. Dusling and D. Teaney, Phys. Rev. C 77, 034905 (2008).
9. D. Molnar and P. Huovinen, J. Phys. G35, 104125 (2008).
10. F. Karsch, D. Kharzeev and K. Tuchin, Phys. Lett. B663, 217-221 (2008).
11. Y. Aoki et al (WHOT-QCD Collaboration), arXiv:0810.4742 (2008).
12. Hama et al., Nucl.Phys.A774,169 (2006)
13. A. Nakamura and S. Sakai, Phys. Rev. Lett. 94, 072305 (2005).
14. J. Noronha-Hostler, J. Noronha and C. Greiner, arXiv:0811.1571 (2008).
15. H. Song and U. Heinz, arXiv:0812.4274 (2008).
16. G. Denicol, T. Kodama, T. Koide and Philipe Mota, J. Phys. G36 035103 (2009).
17. L. D. Landau and E. M. Lifshitz, Fluid Mechanics (Pergamon, New York, 1959).
18. T. Koide, G. Denicol, Philipe Mota and T. Kodama, Phys. Rev. C75, 034909 (2007).
19. G. S. Denicol, T. Kodama, T. Koide and Ph. Mota, J. Phys. G35, 115102 (2008).
20. I. Müller, Living Rev. Rel. 2, 1 (1999), D. Jou, J. Casas-Vázquez and G. Lebon, Rep. Prog. Phys. 51, 1105 (1988); Rep. Prog. Phys. 62, 1035 (1999), I. Müller, Z. Phys. 198, 329 (1967), W. Israel, Ann. Phys. (N.Y.) 100, 310 (1976), J. M. Stewart, Proc. Roy. Soc. A357, 59 (1977), W. Israel and J. M. Stewart, Ann. Phys. (N.Y.) 118, 341 (1979), H. Grad, Commun. Pure Appl. Math. 2, 331 (1949).
21. L.B. Lucy ApJ. 82, 1013 (1977), J.J. Monaghan, Annu. Rev. Astron. Astrophys. 30, 543 (1992).
22. C. E. Aguiar, T. Kodama, T. Osada e Y. Hama, J. Phys. G27, 75 (2001).
23. T. Hirano, K. Tsuda, Phys. Rev. C66, 054905 (2002)
24. T. Koide, E. Nakano and T. Kodama, arXiv:0901.3707 [hep-th].
25. M. A. York and G. D. Moore, arXiv:0811.0729 [hep-ph].
26. M. Prakash et al, Phys. Rep. 227, 321 (1993).
27. A. Monnai and T. Hirano, arXiv:0903.4436 [nucl-th].
28. G. S. Denicol, T. Kodama, T. Koide and Philipe Mota, Phys. Rev. C78, 034901 (2008).
29. G. Torrieri and I. Mishustin, 2008 Phys. Rev. C78, 021901(R).
30. J.M. Burgers, Adv. Appl. Mech. 1, 171 (1948).
31. W. Israel and J. M. Stewart, Ann. Phys. (N.Y.) 118, 341 (1979); A. Muronga, Phys. Rev. C76, 014910 (2007); B. Betz, D. Henkel and D. H. Rischke, arXiv:0812.1440 [nucl-th].