Analytical formulas, general properties and calculation of transport coefficients in the hadron gas: shear and bulk viscosities.

Analytical formulas, general properties and calculation of transport coefficients in the hadron gas: shear and bulk viscosities.

Oleg Moroz Bogolyubov Institute for Theoretical Physics, Kiev, Ukraine
July 12, 2019
Abstract

Elaborated calculations of the shear and the bulk viscosities in the hadron gas, using the ultrarelativistic quantum molecular dynamics (UrQMD) model cross sections, are made. These cross sections are analyzed and improved. A special treatment of the resonances is implemented additionally. All this allows for better hydrodynamical description of the experimental data. The previously considered approximation of one constant cross section for all hadrons is justified. It’s found that the bulk viscosity of the hadron gas is much larger than the bulk viscosity of the pion gas while the shear viscosity is found to be less sensitive to the hadronic mass spectrum. The maximum of the bulk viscosity of the hadron gas is expected to be approximately in the temperature range with zero chemical potentials. This range covers the critical temperature values found from lattice calculations. We comment on some important aspects of calculations of the bulk viscosity, which were not taken into account or were not analyzed well previously. Doing this, a generalized Chapman-Enskog procedure, taking into account deviations from the chemical equilibrium, is outlined. Some general properties, features, the physical meaning of the bulk viscosity and some other comments on the deviations from the chemical equilibrium supplement this discussion. Analytical closed-form expressions for the transport coefficients and some related quantities within a quite large class of cross sections can be obtained. Some examples are explicitly considered. Comparisons with some previous calculations of the viscosities in the hadron gas and the pion gas are done.

bulk viscosity, shear viscosity, hadron gas, pion gas, kinetics, transport theory, chemical equilibrium, chemical freeze-out, crossover
pacs:
25.75.-q, 24.10.Pa, 47.45.Ab, 51.20.+d

I Introduction

The bulk and the shear viscosity coefficients are transport coefficients which enter in the hydrodynamic equations, and thus are important for studying of nonequilibrium evolution of any thermodynamic system.

There are two more additional reasons to study the shear viscosity. The first one is the experimentally observed minimum of the ratio of the shear viscosity to the entropy density near the liquid-gas phase transition for different substances, which may help in studying of the quantum chromodynamics phase diagram and finding of the location of the critical point Csernai:2006zz (); lacey ()111Fireballs, created in heavy ion collisions, have finite sizes and finite times of existence of their thermalized part. This puts important restrictions on detection of the critical fluctuations of thermodynamic functions Stephanov:1999zu (). Because of this it’s also important to consider nonequilibrium dissipative corrections and nonequilibrium phenomenons like critical slow down/speed up.. Such a minimum was observed in theoretical results in several models, see e. g. chakkap (); Dobado:2009ek (). For a counterexample see Chen:2010vf () and references therein. The second reason is the calculation of the in strongly interacting systems, preferably real ones, to compare physical inputs which provide small values of the . The conjectured lowest bound222In Danielewicz:1984ww () the bound coming from the Heisenberg uncertainty principle was obtained for the . However, it was obtained using a formula, which is justified in rarified gases with short-range interactions. It’s well known already from the nonrelativistic kinetic theory that dense gases get corrections over the particle number densities (see e. g. landau10 (), Section 18), corresponding to more than binary collisions, and in very dense gasses this bound can be quite inaccurate. In liquids and other substances the mechanism of appearance of the shear viscosity may be different (see Schafer:2009dj () for a review). In particular, the shear viscosity of water can be very well described by a phenomenological formula with an exponential dependence on the inverse temperature, see e. g. Sengers (). adscftbound () was violated with different counterexamples. For some reasonable ones see Buchel:2008vz (); Sinha:2009ev (). Also see the recent review Cremonini:2011iq (). The bulk viscosity, being very sensitive to violation of the equation of state and being connected with fluctuations through the fluctuation-dissipation theorem Callen:1951vq (), can have a maximum near a phase transition Kharzeev:2007wb (); Karsch:2007jc (); chakkap (); Dobado:2012zf (). In chakkap () and Dobado:2012zf () sharp maxima were observed in the bulk viscosity and the ratio in the linear -model for the vacuum mass . Decreasing the vacuum mass the maximum eventually disappears. Any maximum of the was not observed in the large-N limit of the linear -model in the Dobado:2012zf (). Any maximum of the was not observed in the large-N limit of the -dimensional Gross-Neveu model FernandezFraile:2010gu () (see also Sec. II for comments).

Whether one uses the Kubo333The Kubo formulas are distinguished from the Greet-Kubo formulas e. g. in Muronga:2003tb (); Kadanoff (). formula or the Boltzmann equation one faces nearly the same integral equation for the transport coefficients jeon (); jeonyaffe (); Arnold:2002zm (). The preferable way to solve it is the variational (or Ritz) method. Due to its complexity the relaxation time approximation is used often in the framework of the Boltzmann equation. Though this approximation is inaccurate, does not allow to control precision of approximation and can potentially lead to large deviations. The main difficulty in the variational method is in calculation of collision integrals. To calculate any transport coefficient in the lowest order approximation in a mixture with a very large number of components (like in the hadron gas) one would need to calculate roughly 12-dimensional integrals if only the elastic collisions are considered. Fortunately, it’s possible to simplify these integrals considerably and perform these calculations in a reasonable time.

This paper contains calculations of the shear and the bulk viscosity coefficients for the hadron gas using the (corrected, see Sec. III) UrQMD cross sections. The calculations are done in the framework of the Boltzmann equation with the classical Maxwell-Boltzmann statistics, without medium effects and with the ideal gas equation of state. The Maxwell-Boltzmann statistics approximation allows one to obtain some relatively simple analytical closed-form expressions. Originally the calculations in the same approximations for the hadron gas but with one constant cross section for all hadrons were done in Moroz:2011vn (). The deviations in the worst cases are relatively small. In that paper some analytical formulas of the viscosities for 1-, 2- (explicitly) and N-component (up to solution of the matrix equation) gases with constant cross sections were obtained. Analogical formulas can be written down for quite a large class of non-constant cross sections, in particular, for the ones which appear in the chiral perturbation theory. The final expressions may become somewhat more cumbersome; anyway this is better than numerical integration at least in the speed of the computation. Explicit formulas for the viscosities with the elastic pion-pion isospin averaged cross section and somewhat more general one are obtained in the present paper. The results of the Moroz:2011vn () are partially reproduced in the present paper, improving the text and adding more detailed explanations. The presented calculations can be considered as quite precise ones at low temperatures where the elastic collisions dominate and the equation of state is close to the ideal gas equation of state. At higher temperatures the calculations with the total cross section are expected to give the qualitative description.

For comparison the calculations of the viscosities are performed for the pion gas (throughout the paper the chemical potentials are equal to zero if else is not stated). The results are relatively close to the results in prakash (); davesne (). There the calculations are made in the same approximations except for the davesne (), where the Bose-Einstein statistics is used instead of the Maxwell-Boltzmann one. The discrepancies from the used classical statistics are not large at zero chemical potential and become larger as the chemical potential grows (see Sec. IV for the errors and comments). The comparison is made with the results of the prakash (), see fig. 1. The discrepancies up to a factor of for the bulk viscosity and up to a factor of for the shear viscosity come most probably from somewhat different elastic plus the quasielastic (through the intermediate -resonance) cross section of the prakash ()444The author could not reproduce this plotted total cross section by its formula. In fact it was approximately 2.6 times larger. But the plotted total cross section is quite close to the isospin averaged (corrected) UrQMD total cross section. A notable deviation is only at , when the UrQMD cross section becomes times smaller. At the UrQMD cross section is a little larger instead. (the averaging over the scattering angle is expected to give small errors; see also comments below the formula (128)). The minima of the shear viscosities near are attributed to the peaks from the -resonances in the cross sections. It’s not noticeable in the figure for the dash-dotted line. Nonzero values of the bulk viscosities and theirs maxima are solely due to the masses of the pions. The paper prakash () implements also the isospin averaged current algebra elastic cross sections. These cross sections can be reproduced in the lowest order in the chiral perturbation theory Scherer:2002tk (). They obviously have quite large deviations from the experimental data at high enough energies and wrong asymptotic dependence, which can be seen from the comparison of them with the isospin averaged elastic plus the quasielastic experimental cross sections in the prakash (). The elastic cross sections are rather close to the constant Bass:1998ca (); Bleicher:1999xi ().

In several papers the bulk viscosity was calculated for the pion gas, using the chiral perturbation theory (or the unitarized chiral perturbation theory) and some other approaches, with quite large discrepancies between the quantitative results. In fernnicola () the calculations were done by the Kubo formula in a rough approximation. There the number-changing processes were neglected too, and the non-vanishing value of the bulk viscosity is obtained due to a trace anomaly and the pions’ masses. At small temperatures, where the effects of the trace anomaly are small, the magnitude of the bulk viscosity is large in compare to the results of this paper and the prakash (); davesne (). For example, at () it’s larger approximately in 39 (8) times than the bulk viscosity in this paper. The maximal values differ in approximately 24 times. In lumoore () the calculations are done in the framework of the Boltzmann equation and have a divergent dependence of the for because of remained weak number-changing processes (at the bulk viscosity is nearly 57 times larger than the bulk viscosity calculated in this paper). This dependence should change at low enough temperatures, or higher ones for the pion gas, see Sec. II for explanations. Joining the results of the calculations at low and high temperatures, the function may turn out to be not continuous at the middle temperatures (which is not a physical effect, see Sec. II), and the smooth function is to be obtained through some interpolation. In dobado () the bulk viscosity was calculated in the framework of the Boltzmann equation with the ideal gas equation of state and only the elastic collisions taken into account. The Inverse Amplitude Method was used to get the scattering amplitudes of the pions. The quantitative results are close to the results in this paper (discrepancies up to a factor of ). In chenwang () the calculations are done in the framework of the Boltzmann equation for the massless pions. There the bulk viscosity increases rapidly so that the ratio increases with the temperature.

Calculations of the shear viscosity in the hadron gas with a large number of components were done in Gorenstein:2007mw (), using some approximate phenomenological formula, and in toneev (), using the relaxation time approximation. These results are in good agreement with the calculations of this paper. Hence, as long as the ratio calculated in the toneev () for the free massive pion gas is times larger (in the temperature range with the deviations growing as the temperature decreases) than the one calculated in this paper, one can suspect that the difference comes from the bulk viscosity because of the used relaxation time approximation555In the relaxation time approximation the bulk viscosity source term is treated somewhat differently: the becomes proportional to the integral of the squared source term (times some functions of momentum) and not to the square of the integrated source term (times some functions of momentum). Note that in the fernnicola () the used formula has this relaxation time approximation form. Also there the source term is the one of a system with the inelastic processes. These facts could help to understand the enlarged values of the bulk viscosity. Not small quantitative discrepancies can be noticed between the calculations of the chakkap () and the Dobado:2012zf (). and likely not conserved particle numbers at low temperatures, provided that the SHMC model’s cross sections, used in the toneev (), don’t have large deviations from the UrQMD cross sections or the experimental data, which seems to be the case. Also note that the results in the toneev () for the free particles and the SHMC model don’t differ very much. These facts may explain why the of the hadron gas in the toneev () is times larger (in the temperature range ) than the calculated in this paper. At the low temperature and the vanishing chemical potentials it is 11.3 times more (at the same temperature the factor is 8.2 for the case of the pion gas). In nhngr () the calculation of the bulk viscosity is done for the hadron gas (with an excluded-volume equation of state) with the masses less than using some special formula, obtained though some ansatz Kharzeev:2007wb (). Its quantitative accuracy has not been clarified. The ratio in the nhngr () deviates from the of this paper up to a factor of 1.8 in the temperature range and is different on at .

Also the shear viscosity has been calculated using the Kubo formula (or the Green-Kubo formula) in a gas of mesons and their resonances Muronga:2003tb (). There the UrQMD simulations are performed to calculate the energy-momentum tensor, used in the calculations by the Kubo formula. The in the Muronga:2003tb () is times smaller then the for the hadron gas in this paper. At it is times smaller. In Muroya:2004pu () similar calculations, using the URASiMA event generator, are done for the shear viscosity with close results.

The structure of the paper is the following. A misleading viewpoint on the bulk viscosity, connected with the inelastic processes, is commented on in Sec. II together with some properties, features and physical meaning of the bulk viscosity. In that section some questions concerning the deviations from the chemical equilibrium are addressed too. Sec. III contains some comments on the constant cross sections, which are used in approximating calculations, and some other general comments on cross sections. Also it contains a description of the UrQMD cross sections, which are used in the main calculations, together with their analysis, corrections and the consequences of the corrections for the freeze-out temperatures. The applicability of the used through the paper approximations is discussed in Sec. IV. The system of the Boltzmann equations, its solution and formal expressions of the transport coefficients can be found in Sec. V. The numerical calculations for the hadron gas are presented in Sec. VI. In Sec. VII.1 analytical results for the single-component gas are presented. In particular, an analytical expression for the first order single-component shear viscosity coefficient with constant cross section, found before in anderson (), is corrected while the bulk viscosity coefficient remains the same. The nonequilibrium distribution function in the same approximation is written down. Also the viscosities with some non-constant cross sections are written down. Some analytical results for the binary mixture with constant cross sections are considered in Sec. VII.2. Integrals of source terms needed for the calculation of the transport coefficients can be found in Appendix A. The general entropy density formula can be found in Appendix B. It is used in the numerical calculations for the hadron gas. Transformations of collision brackets, being the 12-dimensional integrals which enter in the viscosities, and some analytical formulas for them can be found in Appendix C. The closed-form expressions for collision rates, mean free paths and mean free times are included in Appendix D.

Ii Some features and properties of the bulk viscosity

First, it should be reminded that the transport coefficients are defined as coefficients next to their gradients in the formal expansion of the energy-momentum tensor and the charge density flows over the gradients of the thermodynamic functions or the flow velocity (see e. g. landau6 (), Section 136). The Kubo formulas are not definitions of the transport coefficients, as one might think. They may introduce some assumptions. In particular, the Kubo formulas in the form as in the jeon () have zero frequency and zero momentum limits, which neglect finite size and finite time effects. Zero momentum limit implies the thermodynamical limit. This limit is needed to avoid possible nonphysical contributions from inappropriate choice of a current and an ensemble Kadanoff (). The Kubo formulas in the form as in Muroya:2004pu (); Kubo () suppose thermal equilibrium in the initial moment of time . So that any infinite space-time scale cannot be connected with the transport coefficients by their definitions.

The Boltzmann equations will be used in what follows. In the case of the elastic collisions they can be derived from the Liouville equation666The Boltzmann equations can also be derived for the case of the inelastic collisions from some physical considerations, see groot () (Chap. I, Sec. 2). in the approximations ( is the effective radius of two particle interactions between the particles of the species and ) that is for rarified gases with short-range interactions777This is the case of interest. Coulomb interactions can be neglected in heavy ion collisions at all the considered energies in this paper.. Also the linear integral equations for the viscosities and other transport coefficients, derivable from the Boltzmann equation, can be obtained (with some corrections) from the perturbative calculations for quantum field theories at finite temperature (including the inelastic processes) using the Kubo formulas jeon (); Gagnon:2006hi (); Gagnon:2007qt (), which justifies application of the Boltzmann equation when the inelastic processes are present.

The bulk viscosity can reveal itself only when there is a nonzero divergence of the flow velocity. This nonequilibrium perturbation should not be confused with another possible independent perturbation (as was done in several papers, some of which are mentioned below; the roots of the misleading viewpoint, perhaps, can be found in landau6 (), Section 81). Namely, this is the homogeneous perturbation. It can be both the chemical and the kinetic one888The inclusion of this kinetic perturbation is similar to the inclusion of the chemical one so that it is omitted for simplicity below. Usually this perturbation should fade first because also the inelastic processes are responsible for the relaxation of the momentum spectra. However, see comments for -dimensional systems below.. Then it can be generalized and made dependent on the coordinate. It just should not be proportional to any gradient. Then the constraints of the local conservation laws should be imposed on these perturbations. The perturbations for the particle numbers should be such that don’t violate conservation of all charges. Considering the case of homogeneous chemical perturbation in a gas with fixed volume, one concludes that the temperature should change with time, being some energy per particle. So that energy conservation should be obtained varying the temperature. Adding an infinitesimal correction to the temperature one gets a perturbation of the form . Such perturbations don’t contribute to all collision integrals. To describe purely chemical perturbations they have to be chosen in the form of the momentum-independent terms (except for the terms), otherwise there will be contributions from the elastic collision integrals. Such perturbations can be considered as chemical potentials-like ones (being small, one can expand the distribution functions over them and get these momentum-independent terms) with the arguments for maximization of the entropy. To find the evolution of these terms they should be separated. Let’s write this in some formulas. Multicomponent gas distribution functions with the leading perturbations can be represented in the form (detailed definitions can be found in Sec. V.1)

 fk=f(0)k(1+~φk)(1+φk)≈f(0)k(1+~φk+φk),|~φk|≪1,|φk|≪1, (1)

where are the perturbation due to the gradients and are the chemical perturbations999Note that if the -th species have conserved particle numbers, then the nonequilibrium chemical potential is nonphysical or redefining the usual (thermodynamic) chemical potential.. Following steps of Sec. V.1, one can get the following linearized equations from the Boltzmann equations:

 pμk(UμD+∇μ)f(0)k+f(0)kpμk(UμD+∇μ)~φk≈−f(0)kLk[φk]−f(0)kLinelk[~φk], (2)

where and is the sum of the linearized collision integrals (divided on the ) of all the processes and of the inelastic processes correspondingly. The 2-nd order gradients and the squared 1-st order gradients are neglected in the l. h. s. because they are of the next order101010The question of validity of this expansion over the gradients (which coincides with the usual order counting in the formal expansion over the gradients in the hydrodynamics) for some profiles is not discussed in this paper. and should be cancelled in the next iteration by the next corrections to the distribution functions. Also the smallness of the is used. If the spatial covariant gradients (at the initial moment of time) are of the same order as the gradients of the thermodynamic functions or the flow velocity, then the terms in the l. h. s. should be retained111111It’s a reasonable assumption in the case when the hydrodynamical description is applicable. For example, the chemical perturbations can be a result of a fast previous expansion (faster than the chemical equilibration). Then the inhomogeneities of the chemical perturbations should correlate with the inhomogeneities of the thermodynamic functions, the flow velocity or it’s divergence.. The covariant temporal derivatives are needed to describe the temporal evolution of the . Then the equations (2) can be split onto the separate equations for the and the

 pμkUμD~φk+pμk∇μ~φk≈−Linelk[~φk], (3)
 pμkUμDf(0)k+pμk∇μf(0)k≈−f(0)kLk[φk]. (4)

The equations (2) can be split within the framework of the perturbation theory over the gradients. Let’s consider also the condition in the (1). Then neglecting the in the (1) and repeating the steps of Sec. V.1, one can get the following linearized equations:

 pμkUμD~φk≈−Linelk[~φk]. (5)

The equations (5) are precise in the homogeneous case (the approximation is only from the linearization). The 1-st order gradients and the are neglected. Then using the (5) and the (2), one can get

 (pμkUμD+pμk∇μ)f(0)k+f(0)kpμk∇μ~φk≈−f(0)kLk[φk]. (6)

Solving the system of equations (5) in the local rest frame, one gets the leading exponential fading dependencies on time121212If the expansion rate is much larger than the collision rates of the inelastic processes (e. g. because of a substantial decrease of the temperature), then the chemical perturbations should enlarge instead. If the r. h. s. of the (3) is smaller than the second term of the l. h. s., then one can consider another approximation, when the -th species particle numbers are conserved. Then the chemical perturbation becomes an addition to the thermodynamic chemical potential. (in a covariant form this should be an explicit space-time dependence). Such dependencies were obtained in some previous studies, see e. g. Matsui:1985eu (); Song:1996ik (). The equations (6) are different from the ones obtained from the common Chapman-Enskog procedure (see e. g. groot (), Chap. V) because of the terms. The contributions from the small chemical perturbations can be neglected in the considered order in the transport coefficients because they are multiplied on the 1-st order gradients. The terms can be cancelled, introducing terms proportional to the into the terms. If the spatial distributions of the are such that are of the 2-nd or a higher order, then the can be neglected. This assumption or approximation is used in the calculations of this paper. In the linear response theory one can also introduce independent small chemical perturbations with the same conclusions for the 1-st order transport coefficients and find evolution of the perturbations with time.

Note that the deviation from the chemical equilibrium itself is not necessarily a source of the bulk viscosity, as is stated in Paech:2006st (). If the bulk viscosity is not equal to zero only because of the particles’ masses and they are tended to zero, the bulk viscosity source term and the bulk viscosity tend to zero even if there are inelastic processes (see the end of Sec. V.1). In the Paech:2006st () the independent chemical perturbations and the perturbations due to the gradients were just connected through the perturbations of particle numbers, and the bulk viscosity became proportional to the chemical relaxation time. Formally infinite chemical relaxation time doesn’t imply any divergencies in the chemical perturbations , but rather approximation of conserved particle numbers. Note that the dependence on the strength of the inelastic processes is different for the chemical perturbations and the perturbations due to the gradients. Increasing the strength of the inelastic processes the chemical relaxation time decreases. And the gradients’ relaxation time increases, because the transport coefficients, at least in rarified gases with short-range interactions, roughly speaking, are inversely proportional to the integrated cross sections (in an ideal liquid the gradients’ relaxation time is infinite). What happens with the bulk viscosity if the inelastic processes become weaker is discussed below.

Making the inelastic processes weaker in compare to the elastic ones, the bulk viscosity eventually gets a formal dominant contribution from them because of the approximate zero mode(s) jeon (), connected with possible conservation of particle number(s)131313If the particles involved into the inelastic processes are massive, then the formal dominant contribution is the exponential one over the temperature and grows as the temperature decreases. If the particles are massless or approximately massless, as in high-temperature QCD Arnold:2006fz (), then a more complicated situation can occur, and one may need to compare some differences of processes’ rates (and not just equilibrium collision rates), arising in the collision matrix ( in assignments of the Arnold:2006fz ()). Under the same pair of used test-functions (indexed by , in the ), and for the same pair of particle species, smaller differences of processes’ rates can be neglected. Comparing among different pairs of test-functions the smallest nonzero contributions dominate, or rather as can be obtained directly from the inverted collision matrix. . As long as it’s clear that the bulk viscosity is not responsible for the chemical equilibration, it’s also clear that there may be the approximation of conserved particle numbers if the momentum spectrum, as well as the gradients, can relax by means of only the elastic collisions (which is usually the case) and the elastic processes make a dominant contribution to the collision rates. The question is only at what concrete temperature does this approximation sets in. Let’s make an illustrative example of what nonphysical contributions one can get from formally remained weak inelastic processes. Consider infinitely weak inelastic processes and the perturbation of the flow velocity such that the energy-momentum tensor gets a sizable contribution from the bulk viscosity term, not large in compare to the pressure (cf. (20), (27)) to remain the perturbation theory applicable. Then it’s obvious that this contribution is not physical because it is created by the practically absent processes and the infinitesimal perturbation of the flow velocity. Instead, this system is practically described by the equilibrium thermodynamic functions. This also answers positively the question whether the thermodynamic chemical potential can be introduced for approximately conserved particle number in principle. As far as the author knows, the first correct comment (albeit somewhat inaccurate) on this issue can be found in the jeonyaffe (). However, note that in fact there is no divergent mean free paths, corresponding to the inelastic processes (IMFP) in this case. They are cut by the mean free paths, corresponding to the elastic processes (and the overall collision rate have the dominant contribution from the elastic collisions). So that it may be not necessary for the chemical relaxation time to be much larger than any relevant time scale (like the gradients relaxation time or the time of existence of the thermal part of the system) to switch off the inelastic processes. That’s why a criterion based on comparison of collision rates of elastic and inelastic processes can be considered to switch off the inelastic processes. Such a comparison is done in the UrQMD studies of the hadron gas in Bleicher:2002dm () (see Sec. IV for farther discussions). According to Goity (), the chemical relaxation time of the processes in the pion gas is much larger than the thermal relaxation time. And e. g. at the chemical relaxation time is equal to , which is larger than the typical lifetime of the thermal part of the expanding fireball (see e. g. Bleicher:2002dm ()). So that it’s the inelastic processes which should be neglected in the pion gas at or even higher temperatures, which wasn’t done in the lumoore (). To show importance of the gradients relaxation time, let’s consider the following possible case. Let’s consider the only perturbation - propagating sound wave, perturbed in a point. It’s possible for the IMFP to be much larger than the gradients relaxation size (on which the wave can be considered as damped) and be much smaller than the system’s size at the same time. Then, the bulk viscosity cannot be defined by the IMFP in this case, because it enters in the sound attenuation constant. Thus, the gradients relaxation size and time are cutting parameters. Note that they exists even in infinite systems considered during infinite time interval.

The bulk viscosity source terms increases substantially if particle numbers are not conserved (cf. (120), (121); in mixtures these particle numbers should also be not small). This reflects additional fluctuations from not conserved particle numbers. Though the inelastic processes have to be effective enough to consider the approximation of not conserved particle numbers. Perhaps, the point at which the bulk viscosities in the different approximations cross can provide a criterion for switching on/off the inelastic processes. If this is not so, then one would have to make some interpolation in the intermediate region141414Perhaps, the bulk viscosity calculated without constant test-functions (except for zero modes of the inelastic collision integrals, used to conserve charges) can provide a good interpolation.. Note that e. g. in the calculations by the Kubo formulas through the direct calculations of the energy-momentum tensor as in the Muronga:2003tb () it’s not needed to use the approximation of conserved or not conserved particle numbers (which defines the number of independent thermodynamic chemical potentials, through which the chemical potentials of all particles are expressed, cf. (18)). There the energy-momentum tensor should be a smooth function of time and the thermodynamic functions as long as the inelastic processes fade smoothly. Then the bulk viscosity should be a smooth function of the temperature and particles’ chemical potentials regardless of the number of the independent chemical potentials.

In the Arnold:2006fz () a bottleneck for the relaxation to equilibrium characterized by the bulk viscosity due to the weakest processes’ rates is assumed. Instead, there are rather dominant contributions from some test-functions151515Not a bottleneck from some perturbations, because one actually doesn’t have a choice in the form of the momentum dependence of the perturbations corresponding to the transport coefficients. The kinetic perturbation can be of different forms of the momentum dependence. (as is commented in the footnote 13), which should not be specially treated though, except for the ones which are the approximate zero modes making a dominant contribution. A similar dominance161616Another similar dominance can exist from particle species interacting weakly with all particles. is present also in other transport coefficients, in particular, when there is only one type of processes. Although in QCD at high enough temperatures the equilibrium elastic collisions rate is parametrically the largest one171717The estimate can be easily inferred from Arnold:2000dr ()., , because of cancellations the momentum transfer takes place with the rate , which is parametrically smaller than the particle number change rate from the effective ”” processes. This provides an example when the equilibrium collision rates may differ substantially from the relevant collision rates. The ”” processes provide small chemical relaxation time in compare to the thermal relaxation time, which justifies the approximation of not conserved particle numbers and the enhancement of the bulk viscosity from the source terms at least at small enough , whereas the contributions to the bulk viscosity from the collision integrals of the ”” processes are suppressed at small enough (the inelastic processes are not suppressed, but they are of the order ). To avoid misunderstanding it may be mentioned that taking the total collision rate of the ”” processes as formally infinite by taking the corresponding matrix elements as formally infinite ones, one gets zero bulk viscosity and zero mean free paths as long as both the gluons and quarks take part in these processes (see also footnote 24).

In the case of a -dimensional single-component gas the elastic collisions cannot result in the relaxation of the momentum spectra and, hence, cannot stimulate the system to evolute towards equilibrium181818There are forward scatterings and momentum interchange. As long as the particles are not distinguishable the momentum interchange from the elastic collisions is equivalent to the forward scatterings or absence of the elastic collisions at all.. As a result, the exponentially divergent bulk viscosity was obtained in the paper FernandezFraile:2010gu (). Considering again the example about the infinitely small perturbation of the flow velocity and assuming also a finite size of the system, it’s again obvious that the weak inelastic processes may make nonphysical contributions (in this case the mean free path is formally cut by the system’s size). If this is the case, then the hydrodynamical description becomes inapplicable, and might use simulations of particles’ collisions or the Boltzmann equations in the approximation without collisions (on a time scale much smaller than the chemical relaxation time). If the -dimensional description is only an approximate one (that is with small angle elastic scatterings in higher dimensions), the relaxation of the momentum spectrum by the elastic collisions should be considered. And if a -dimensional gas has at least two components with different masses, then a nontrivial momentum exchange in the elastic collisions is possible. This results in the possibility of the relaxation of the momentum spectra by only the elastic collisions Cubero ().

Let’s summarize this section with formulation of the physical meaning of the bulk viscosity. The bulk viscosity reflects deviation of the value of the pressure from its local equilibrium value (as can be seen from the (27)), appearing when the system expands/compresses, because of the delay in the equilibration. The bulk viscosity is not responsible for the restoration of the chemical or the kinetic equilibria - it’s responsible for the relaxation of the divergence of the flow velocity. If there are inelastic processes, then the particle numbers also get nonequilibrium contributions (cf. (14), (37), (68)) such that the charge is conserved locally (cf. (74))191919One should keep in mind that while studying the chemical perturbations through the thermodynamic functions first the contributions from the transport coefficients’ terms should be subtracted.. Though these contributions together with the contribution to the pressure may become nonphysical because of the approximate zero modes (if such ones appear in the calculations). The magnitude of the bulk viscosity changes from theory to theory. Under some quite general assumptions a nonzero value of the bulk viscosity can be connected with violation of the scale invariance due to a nonzero value of the energy-momentum tensor Coleman:1970je (); Callan:1970ze (). Of course, the beta function can contribute to the energy-momentum tensor and the bulk viscosity too jeon ().

Iii The hard core interaction model and the UrQMD cross sections

In a non-relativistic classical theory of particle interactions there is a widespread model, used in approximate calculations, called the hard core repulsion model or the model of hard spheres with some radius . For its applications to the high-energy nuclear collisions see Gorenstein:2007mw () and references therein. The differential scattering cross section for this model can be inferred from the problem of scattering of point particle on the spherical potential if and if landau1 (). In this model the differential cross section is equal to . To apply this result to the gas of hard spheres with the radius one can notice that the scattering of any two spheres can be considered as the scattering of the point particle on the sphere of the radius , so that one should take . The total cross section is obtained after integration over the angles of the which results in the . For collisions of hard spheres of different radiuses one should take or replace the on the :

 σtot,kl=π(rk+rl)2. (7)

The relativistic generalization of this model is the constant (not dependent on the scattering energy and angle) differential cross sections model.

The hard spheres model is classical, and connection of its cross sections to cross sections, calculated in any quantum theory, is needed. For particles, having a spin, the differential cross sections averaged over the initial spin states and summed over the final ones will be used.202020It’s assumed that particle numbers of the same species but with different spin states are equal. If this were not so then in approximation, in which the spin interactions are neglected and probabilities to have certain spin states are equal, the numbers of the particles with different spin states would be approximately equal in the mean free time. With equal particle numbers their distribution functions are equal too. This allows one to use the summed over the final states cross sections in the Boltzmann equations. If colliding particles are identical and their differential cross section is integrated over the momentums (or the spatial angle to get the total cross section) then it should be multiplied on the factor to cancel double counting of the momentum states. These factors are exactly the factors next to the collision integrals in the Boltzmann equations (30). The differential cross sections times these factors will be called the classical differential cross sections.

The UrQMD cross sections are used in the numerical calculations of Sec. VI212121Very high energy dependence of any used UrQMD cross section is not important because of the exponential suppression . The used cross sections were cut on the and were continued by a corresponding constant continuously at higher energies. At small enough momentums there is another somewhat weaker suppression. The momentum space density of each particle provides suppression. This may (partially) suppress some deviations from the experimental data of some UrQMD cross sections (like for the pair) at . To estimate at what temperatures some discrepancies in cross sections can appear one can equate the to the sum of the averaged one-particle energies (23) of the two colliding particles.. These cross sections are described in Bass:1998ca (); Bleicher:1999xi (). More details can be found in the UrQMD program codes. Below there is some description mainly of what is different or new.

The UrQMD cross sections are averaged over the initial spin states and summed over the final ones. As long as the UrQMD cross sections are total ones (integrated over the scattering angle), the factors are already absorbed into them (in what follows only such cross sections will be considered in this section tacitly). Dividing them on the , one gets the classical differential cross sections, averaged over the scattering angle.

The UrQMD codes (version 1.3) were modified to get accurately tabulated (with a step of ) cross sections. Resonances’ masses and widths (they are tuned in their uncertainty regions to describe the experimental data better), used in the UrQMD codes, have somewhat different values than the ones in the Bass:1998ca (). Influence of variation of these parameters was studied in Gerhard:2012fj (). The UrQMD codes implement somewhat different averaging of the c. m. momentums over the resonances’ masses222222Averaged powers of the momentums are used, not powers of the averaged momentums. than in the papers Bass:1998ca (); Bleicher:1999xi (). It was found that using the resonance dominating cross sections from the papers Bass:1998ca (); Bleicher:1999xi () some of these cross sections could have a large rise at small c. m. momentums if constant widths are used in the calculations of the averaged c. m. momentums in the energy dependent widths. So that one should be aware of this fact232323It may be mentioned that one should be also aware of possible differences in storing of the floating point numbers in different programming languages or while using different compilers.. The UrQMD codes have a low energy cut-off at (and a similar one over the c. m. momentum if triggered) for the resonance dominating cross sections, and no large low energy rise was found there.

An important ingredient of the UrQMD model is the Additive Quark Model (AQM), which is used for unknown cross sections. Universality of hadrons, based on jet quenching arguments, is used to support this model. This model describes the experimentally known cross sections well at sufficiently high energies. Application of this model is better than elimination of the corresponding hadrons, which is the same as equating their all cross sections to zero and, hence, exclusion of their contributions from the thermodynamic functions (infinite mean paths, no thermalization).

At this point an interruption should be made to consider some important questions related to different types of the UrQMD cross sections. These different types are used due to several reasons and are the following: the elastic cross section(s) (ECS(s)), the elastic plus the quasielastic cross section(s) (EQCS(s)), the total cross section(s) (TCS(s)) and the previous two types with enhanced in some way resonances’ cross sections (index ”2” is appended in the abbreviations).

Of course, the system of the Boltzmann equations would have a solution with any of these cross sections. Usage of the ECSs is completely self-consistent as long as only the elastic collision integrals are used in the calculations of the viscosities. However, there are reasons to consider also the EQCSs. Exactly this type of cross sections, being averaged over the isospin, is implemented in prakash (). The quasielastic cross sections can be used as rightful contributions to the ECSs in the approximation that the 4-momentum of the intermediate resonance does not change (the effects of the exclusion of the resonances as independent particles are considered in Sec. VI). The mean free paths of the intermediate resonances without contributions of the decays, being not equal to zero, also introduce some errors, which are neglected. The EQCSs conserve particle numbers, which is consistent with the only elastic collision integrals, implemented in the calculations. There are also some additional arguments for the usage of these cross sections. From the phenomenological considerations one can take into account shortening of the mean free paths (or enlarging of the collision rates) due to creation of the resonances. In other words, there would have to be contributions from the inelastic collision integrals next to the elastic collision integrals, and they are taken into account approximately by the contributions from the quasielastic cross sections.

Resonances are not just intermediate particles, and they can collide with other particles. They make not negligibly small contribution to the thermodynamic functions and the viscosities. So that they are also included in the calculations as independent particles with their parameters and corresponding collision integrals. They would have to have shortening of their mean free paths from their decays and contributions from the inelastic collision integrals too. These contributions may be taken into account from the following collision rate considerations. A resonance’s decay rate can be approximately replaced with just its total width. Then, given a resonance, one would have to redistribute its width (that is not changing the whole collision rate containing the contribution of the decay rate) in such a way that the cross section of the collision of this resonance with a resonance of the same species gets an addition242424This enhancement leads to the shortening of the mean free paths of the resonances of only this species, as needed. In the formal limit of this infinitely large enhancement other collision integrals can be neglected and the Boltzmann equation for this species decouples. Then from the solution of the Boltzmann equation for a single-component gas (see Sec. VII.1) one concludes that the nonequilibrium perturbation to the distribution function of this species vanishes in this limit. Note that infinitely strong interactions also with particles of all other particle species would result in zero transport coefficients.. Using an approximate expression for the collision rates (in the nonrelativistic approximation, applicable in this case) from Appendix D, one easily finds the addition (where is the width) to the . Such cross sections seem to be the most physically preferable ones because they take into account more realistic mean free paths than in the previous case while not violating the conservation of the particle numbers too.

The TCSs are used to take into account even larger shortening252525Not the largest one. The effect of the enhancements of the resonances’ TCSs is of for the bulk viscosity and of for the shear viscosity so that TCS2s are additionally considered. of the mean free paths than in the case of the EQCSs. However, such cross sections introduce some inconsistency, implying that the conservation of the particle numbers is violated. As long as there may be contributions from some partial cross sections to the UrQMD ECSs or the EQCSs which were not taken into account (see below), the TCSs can be used as the upper bounds for the ECSs and the EQCSs. However, it’s expected that these bounds are excessively high. If so, the TCSs (rather TCS2s) can be considered not only as the approximation taking into account real mean free paths but also as some measure of deviation from the approximation of only the elastic and the quasielastic collisions with the following arguments. If the TCSs were approximately equal to the ECSs or the EQCSs, or the numbers of particles with large inelastic cross sections were small, then one could expect small errors due to the negligibility of the inelastic collisions.

Continuing the discussion of the details of the UrQMD cross sections, it should be mentioned that the UrQMD TCSs are the most reliable ones. The sum of the partial cross sections is not always equal to the TCSs by their construction. If this is the case, then some partial cross sections are rescaled depending on their reliability262626This information, including some other information about the cross sections, is stored in the array SigmaLn of the file blockres.f.

The magnitudes of the partial cross sections, implemented in the UrQMD codes, are used to determine, what a partial cross section to choose in a given collision, using a random number generator. Among these partial cross sections there are the ECSs. Exactly these ECSs are used in the present calculations. However, if a partial cross section with a string excitation is chosen in a given collision, there is a probability to end up with the elastic collision if the is too small. These contributions to the ECSs are not calculated and are not added to the ECSs. Also the string excitations can, possibly, end up with creation of a resonance. Contributions to the EQCSs from the string excitations are taken into account partially (see below).

The ECSs, if not known from the experiment, are taken in the form of some extrapolations, discussed below, or the AQM is used. The normalization on the corresponding TCSs can change the ECSs notably. The meson meson (MM) ECSs are equal to . The meson baryon (MB) ECSs are equal to the AQM rescaled experimental cross sections. But after the normalization they become equal to zero in the resonances dominated energy range (approximately below ). The anti-baryon baryon () ECSs are equal to the AQM rescaled experimental cross sections. Other ECSs are equal to the AQM ECSs.

Before discussing the quasielastic cross sections first let’s write for convenience the resonance dominated cross sections formula for a reaction . Correcting a typo and rewriting it in a somewhat different form than in Bass:1998ca (); Bleicher:1999xi (), one gets

 σijtot(√s)=∑RgRgigjπp2cmΓ2R,totbR→ij(MR−√s)2+Γ2R,tot/4,bR→ij≡|⟨ji,mi,jj,mj||JR,MR⟩|2ΓR→ijΓR,tot, (8)

where is the partial energy-dependent width of the decay of the resonance into particles of types and without specification of their isospin projection, is the total energy-dependent width of the decay of the resonance , is the spin degeneracy factor, is the energy-dependent branching ratio. The squared Clebsch-Gordan coefficients allow to specify the branching ratio for the pair of the particles with concrete isospin projections. The squared Clebsch-Gordan coefficients should be normalized in such a way that they give unity after summation over all isospin projections in a given multiplet. This formula represents contributions from all possible resonances through which the reaction can take place. Now it’s easy to write down the cross sections for the quasielastic scatterings:

 σijquasi(√s)=∑RgRgigjπp2cmΓ2R,totb2R→ij(MR−√s)2+Γ2R,tot/4. (9)

One more multiplier takes into account the fact that a resonance decays only into the pair and represents the probability of this decay.

The TCS is not described by the formula (8) completely, and a partial cross section, attributed to the s-channel strings excitations, is added in the UrQMD model to fit the TCS to the experimental data. In the UrQMD model this s-channel strings cross section is added also to other strange meson nonstrange baryon TCSs when annihilation is possible due to the quark content. From comparison with the experimental data for the ECS Beringer () (actually it’s believed to be the EQCS because smaller peaks from the resonances are reproduced there) it was found that the half of the s-channel strings cross section is enough to describe well this experimental cross section. Then the half of the s-channel strings cross section is added to other strange meson nonstrange baryon EQCSs when annihilation is possible. These contributions from the strings excitations are the most low energetic ones. They are the only contributions from the strings excitations which are added. The next in the energy scale possible contributions to the EQCSs may be in the cross sections. In other pairs the string excitations appear approximately from .

There is an important omission, found in the UrQMD codes (present also in the last version 3.3). The function fcgk returns incorrect (two times smaller) values of the squared Clebsch-Gordan coefficients for the resonances dominated cross sections in some cases. The first case is for the pairs of unflavored mesons from the same multiplet with the isospin . For example, the function fcgk returns for the only possible isospin decomposition of the to the pair, because the states and are counted as different ones. As a result, the peak from the -resonance becomes two times smaller than e. g. in the -resonance isospin decomposition. The second less important case is for the pairs of unflavored mesons with the isospin and anti-nucleons. The third even less important case is for the pair and it’s charge conjugate.

Let’s make some comments on the errors what the above-mentioned omissions cause in some quantities at zero chemical potentials, which in turn demonstrate sensitivity to different changes in the cross sections. The errors in the viscosities with the ECSs are less than . The errors in the shear viscosity with the EQCSs (the TCSs) reach () at . Outside the temperature range the errors reach (). The errors in the bulk viscosity with the EQCSs (the TCSs) reach () at . Outside the temperature range the errors reach (). The errors in the total number of collisions per unit time per unit volume (using the TCSs and including the decay rates) reach (at ). Outside the temperature range the errors reach . In view of the errors for the total number of collisions the kinetic freeze-out temperatures found in the UrQMD studies Bleicher:2002dm () should decrease, becoming closer to the experimentally extracted ones (see Heinz:2007in () and references therein). The chemical freeze-out temperature may change in a less extent. This is because both the inelastic and the quasielastic processes’ cross sections (like of the quasielastic collision of the pair and of the reaction ) increase, so that the temperature at which the inelastic processes cease to be dominant may almost not change.

It’s observed that some of the UrQMD detailed balance cross sections (e. g. for the pair) are not symmetric under the particle interchange. This is because the function W3j, calculating the Wigner symbols, doesn’t return zero in some cases. Namely, the selection rule for the sum is not included. In principle, such omission could result in negative values of the essentially non-negative viscosities but, as long as only small fraction of cross sections is affected, this omission has caused only negligibly small errors in the viscosities. But e. g. the error in the TCS is approximately .

Also some fixes of the UrQMD cross sections are made. It’s found that the UrQMD ECS has large deviations from the experimental data Beringer () in the range (the UrQMD cross section reaches in the region ). To fit this cross section to the experimental data it is replaced with the AQM ECS in the range and is interpolated smoothly with the sine function in the range with the cross section being equal to at . This replacement is also applied to other MB ECSs, when annihilation is not possible due to the quark content.

The next fix is for the BB ECSs. It’s found that the UrQMD ECS has quite large deviations form the experimental data Beringer () too. To fit this cross section to the experimental data it is replaced with the AQM ECS in the range and interpolated smoothly with the sine function in the range with the cross section being equal to the AQM TCS at . This replacement is also applied to other BB ECSs.

Some other found lacks result in negligible errors in the viscosities. However, errors in the corresponding mean free paths and possible other quantities may be not negligible ones. Two of such lacks can be mentioned. The first one is the following. The and cross sections are fitted to the experimental data. And their charge conjugates are calculated using general formulas and so cause deviations up to for . The second lack is the following. In some not large energy regions with the resonance dominated cross sections are equal to zero for some small numbers of pairs because there is no resonances which could be created by this pair. These regions are replaced by a constant continuously.

Let’s also comment on the deviations from the fixes described in the last four preceding paragraphs. The altogether deviations in the viscosities and the total number of collisions with the TCSs are less than . The altogether deviations in the viscosities with the ECS or the EQCS are in the range . The largest contribution is from the MB cross sections’ fixes. At the deviations are less than (the temperatures above are not studied).

Iv Conditions of applicability

Before proceeding forth first the applicability of the Boltzmann equation and of the calculations of the transport coefficients should be clarified.

Although the Boltzmann equations are valid for any perturbations of the distribution functions they should be slowly varying functions of the space-time coordinates to justify that they can be considered as functions of macroscopic quantities like the temperature, the chemical potentials or the flow velocity or, in other words, that one can apply thermodynamics locally. Then one can make the expansion over the independent gradients of the thermodynamic functions and the flow velocity (the Chapman-Enskog method), which vanish in equilibrium. Smallness of these perturbations of the distribution functions in compare to their leading parts ensures the validity of this expansion and that the gradients are small272727The magnitudes of thermodynamic quantities can also be restricted by this condition or, conversely, not restricted even if transport coefficients diverge. See also Sec. VII.1 of this paper. The smallness of the shear and the bulk viscosity gradients can also be checked by the condition of smallness of the (27) in compare to the (20). Of course, the next corrections should be small too.. Because these perturbations are inversely proportional to coupling constants one can say that they are proportional to some product of particles’ mean free paths and the gradients. So that, in other words, the mean free paths should be much smaller than characteristic lengths, on which the macroscopic quantities change considerably282828 It’s clear that the mean free paths should be smaller than the system’s size too. .

As is discussed in Sec. II, the inelastic processes may need addition treatment in the calculations of the bulk viscosity. There is a need to specify reasonable conditions when the inelastic processes can be neglected. One could use the following reliable criterion, which takes into account both the particle number densities and the intensity of the interactions:

 ∫t2t1dt∫V(t)0d3x∑n∈ all% channels˜Rinelk′,n<1, (10)

where is the system’s volume, is the number of reactions of particles of the -th species292929Primed indexes run over the particle species without regard to their spin states. This assignment is clarified more in Sec. V.1. over all channels per unit time per unit volume (analog of (299)), is chosen to satisfy the inequality, and is equal to the moment of time at which the divergence of the flow velocity is relaxed (if this time can be estimated reliably with remained inelastic processes) or to the moment of time at which the system becomes practically not interacting (after expansion) because of large cumulative mean free path in compare to the system’s size. Though this criterion is likely to be too strict, and at some higher temperatures the approximation of conserved particle numbers should still work well. The main alternative criterion is based on comparison of collision rates of elastic and inelastic processes (as implemented in the Bleicher:2002dm ()). Using this criterion and some other ones, the chemical freeze-out line303030This is an approximation. In fact this should be a range in which particles of different particle species have their own freeze-out points. in the plane can be built for the hadron gas, see e. g. Cleymans:2005xv (); Andronic:2005yp (). At zero chemical potentials the chemical freeze-out temperature is approximately equal to . The remaining question is how good is the approximation of only the elastic collisions at . From the hydrodynamical description of the elliptic flow at RHIC it’s found that near the chemical freeze-out Dusling:2011fd (). The constant value provides a good description of the elliptic flow both at RHIC and LHC Bozek:2011ph (). It seems that the approximation of conserved particle numbers is not implemented in the bulk viscosity formula used in the nhngr (). The bulk viscosity obtained from it is very close to the one of this paper. These results support the choice of the approximation of only the elastic collisions at and show that the deviations are likely no more than in 2-3 times. Anyway, the numerical calculations by the Kubo formula through simulations of collisions are desirable along and around the chemical freeze-out line for more accurate calculations (though the procedure of collisions of particles introduce some errors itself Bass:1998ca (), which should be kept in mind).

Errors due to the Maxwell-Boltzmann statistics, used instead of the Bose-Einstein or the Fermi-Dirac ones, were found to be small for the vanishing chemical potentials313131It should be mentioned that if the particles of the -th particle species are bosons and if then there is a (local) Bose-Einstein condensation for them, which should be treated in a special way.. According to calculations for the pion gas in davesne (), the bulk viscosity becomes larger at and larger at for the vanishing chemical potential. Although the relative deviations of the thermodynamic quantities of the pion gas at the nonvanishing chemical potential are not more than 323232The relative deviations of the thermodynamic quantities grow with the temperature for some fixed value of the chemical potential and tend to some constant. the bulk viscosity becomes up to times more. The shear viscosity becomes less at and less at for the vanishing chemical potential and less at and less at for the . The corrections to the bulk viscosity of the fermion gas, according to calculations of the bulk viscosity source term, not presented in this paper, are of the opposite sign and approximately of the same magnitude. So that for the hadron gas the error due to the used classical statistics can be even smaller than for the pion gas.

The numerical calculations in Sec. VI of the viscosities with the total cross sections justify the choice of one constant cross section for all hadrons. It’s approximately equal to , corresponding to the effective radius (as given by the (7)), which is used in the estimations below.

The condition of applicability of the ideal gas equation of state is controlled by the dimensionless parameter which appears in the first correction from the binary collisions in the virial expansion and should be small. Here is the so called excluded volume parameter and is the mean volume per particle. One finds at , at and at for the vanishing chemical potentials. Along the chemical freeze-out line (its parametrization can be found in Gorenstein:2007mw ()) the grows from to with the temperature. From comparison with lattice calculations Borsanyi:2010cj () one can find that the corrections to the ideal gas equation of state are small at . One could suspect that even small corrections to the thermodynamic quantities can result in large corrections for the bulk viscosity, though this seems to be not the case. The errors in the bulk viscosity from the scale-violating contributions of the hadrons’ masses are less than the errors from the contributions to the trace of the energy-momentum tensor (for more details see Sec. VI).

One more important requirement, which one needs to justify the Boltzmann equation approach, is that the mean free time should be much larger than (the is the characteristic single-particle energy) danielewicz () or the de Broglie wavelength should be much smaller than the mean free path Arnold:2002zm () to distinguish independent acts of collisions and for particles to have well-defined on-shell energy and momentum. This condition gets badly satisfied for high temperatures or densities. The mean free path of the particle species is given by the formula (315) or the formula (310) if the inelastic processes can be neglected. The wavelength can be written as , where the averaged modulus of the momentum of the -th species is

 ⟨|→pk′|⟩=∫d3pk′|→pk′|f(0)k′(pk′)∫d3pk′f(0)k′(pk′)=2e−zk′T(3+3zk′+z2k′)z2k′K2(zk′)=√8mk′TπK5/2(zk′)K2(zk′), (11)

where , is the modified Bessel function of the second kind. As it follows from the (11) the largest wavelength is for the lightest particles, the -mesons. The elastic collision mean free paths are close to each other for all particle species. Hence, the smallest value of the ratio is for the -mesons. Its value is close to the value of the and is exponentially suppressed for small temperatures too. At the temperature () and the vanishing chemical potentials this ratio is equal to 0.18 (0.7). Along the chemical freeze-out line it grows from to with the temperature.

To go beyond these conditions one can use the Kubo (or Green-Kubo) formulas, for instance. In the jeon () the Kubo formulas were used to perform perturbative calculations of the viscosities in the leading order. Basing on this result, an example of effective weakly coupled kinetic theory of quasiparticle excitations with thermal masses and thermal scattering amplitudes was presented in the jeonyaffe (). There the function (appearing because of the temperature dependence of the mass) takes into account the next in the coupling constant correction to the energy-momentum tensor and the equation of state333333In the hadron gas it’s believed that the vacuum masses are large in compare to their thermal corrections for the most of the hadrons at temperatures or even higher ones. Then expanding over the thermal correction in the matrix elements, one would get even smaller corrections than the ones to the equation of state in coupling constants (because of coupling constants next to the matrix elements) in a perturbation theory, e. g. chiral perturbation theory.. For further developments see Arnold:2002zm (); Gagnon:2006hi (); Gagnon:2007qt (). For some other approaches see Blaizot:1992gn (); Calzetta:1986cq (); Calzetta:1999ps () and Arnold:1997gh () with references therein.

V Details of calculations

v.1 The Boltzmann equation and its solution

The calculations in this paper go close to the ones in groot () though with some differences and generalizations. Let’s start from some definitions. Multi-indices will be used to denote particle species with certain spin states. Indexes will be used to denote particle species without regard to their spin states (and run from 1 to the number of the particle species ) and to denote conserved quantum numbers343434In systems with only the elastic collisions each particle species have their own ”conserved quantum number”, equal to 1.. Quantifiers with respect to the indexes are omitted in the text where they may be needed which won’t result in a confusion. Because nothing depends on spin variables one has for every sum over the multi-indexes

 ∑k...=∑k′gk′..., (12)

where is the spin degeneracy factor. The following assignments will be used:

 n≡∑knk ≡ ∑k′nk′,na≡∑kqaknk,xk≡nkn,xk′≡nk′n,xa≡nan, ^μk ≡ μkT,^μa≡μaT,zk≡mkT,πμk≡pμkT,τk≡pμkUμT, (13)

where denotes values of conserved quantum numbers of the -th kind of the -th particle species. Everywhere the particle number densities are summed, the spin degeneracy factor appears and then gets absorbed into the or the by the definition. All other quantities with primed and unprimed indexes don’t differ, except for rates, the mean free times and the mean free paths defined in Appendix D, the commented below, the coefficients , and, of course, quantities, whose free indexes set the indexes of the particle number densities . Also the assignment will be used for compactness somewhere.

The particle number flows are353535The metric signature is used throughout the paper.

 Nμk=∫d3pk(2π)3p0kpμkfk, (14)

where the assignment is introduced. The energy-momentum tensor is

 Tμν=∑k∫d3pk(2π)3p0kpμkpνkfk. (15)

The local equilibrium distribution functions are

 f(0)k=e(μk−pμkUμ)/T, (16)

where is the chemical potential of the -th particle species, is the temperature and is the relativistic flow 4-velocity such that (with a frequently used consequence ). The local equilibrium is considered as perturbations of independent thermodynamic variables and the flow velocity over a global equilibrium such that they can depend on the space-time coordinate . Additional chemical perturbations could also be considered, but they don’t enter in the first order transport coefficients if they are small, as is discussed in Sec. II. The chemical equilibrium implies that the particle number densities are equal to their global equilibrium values. The global equilibrium is called the time-independent stationary state with the maximal entropy363636The kinetic equilibrium implies that the momentum distributions are the same as in the global equilibrium. Thus, a state of a system with both the pointwise (for the whole system) kinetic and the pointwise chemical equilibria is the global equilibrium.. The global equilibrium of an isolated system can be found by variation of the total nonequilibrium entropy functional landau5 () over the distribution functions with condition of the total energy and the total net charges conservation:

 U[f]=∑k∫d3pkd3x(2π)3fk(1−lnfk)−∑k∫d3pkd3x(2π)3βp0kfk−∑a,kλaqak∫d3pkd3x(2π)3fk, (17)

where are the Lagrange coefficients. Equating the first variation to zero, one easily gets the function (16) with , and

 μk=∑aqakμa, (18)

where are the independent chemical potentials coupled to the conserved net charges.

With , substituted in the (14) and the (15), one gets the leading contribution in the gradients expansion of the particle number flow and the energy-momentum tensor:

 N(0)μk=nkUμ, (19)
 T(0)μν=ϵUμUν−PΔμν, (20)

where the projector

 Δμν≡gμν−UμUν, (21)

is introduced. The is the ideal gas particle number density,

 nk=UμN(0)μk=12π2T3z2kK2(zk)e^μk, (22)

the is the ideal gas energy density,

 ϵ=UμUνT(0)μν=∑k∫d3pk(2π)3p0kf(0)k=∑knkek,ek≡mkK3(zk)K2(zk)−T, (23)

and the is the ideal gas pressure,

 P=−13T(0)μνΔμν=∑k13∫d3pk(2π)3p0k→pk2f(0)k=∑knkT=nT. (24)

Also the following assignments are used:

 e ≡ ϵn=∑kxkek,hk≡ek+T,h≡ϵ+Pn=∑kxkhk, ^ek ≡ ekT=zkK3(zk)K2(zk)−1,^e≡eT,^hk≡hkT=zkK3(zk)K2(zk),^h≡hT. (25)

Above is the enthalpy per particle, is the energy per particle and , are the enthalpy and the energy per particle of the -th particle species correspondingly, which are well defined in the ideal gas.

In the relativistic hydrodynamics the flow velocity needs somewhat extended definition. The most convenient condition which can be applied to the is the Landau-Lifshitz condition landau6 () (Section 136). This condition states that in the local rest frame (where the flow velocity is zero though its gradient can have a nonzero value) each imaginary infinitesimal cell of fluid should have zero momentum, and its energy density and the charge density should be related to other thermodynamic quantities through the equilibrium thermodynamic relations (without a contribution of nonequilibrium dissipations). Its covariant mathematical formulation is

 (Tμν−T(0)μν)Uμ=0,(Nμa−N(0)μa)Uμ=0. (26)

The next to leading correction over the gradients expansion to the can be written as an expansion over the 1-st order Lorentz covariant gradients, which are rotationally and space inversion invariant and satisfy the Landau-Lifshitz condition373737 Also this form of respects the second law of thermodynamics landau6 () (Section 136). (26):

 T(1)μν≡2η∘¯¯¯¯¯¯¯¯¯¯¯¯¯¯∇μUν+ξΔμν∇ρUρ=η(ΔμρΔντ+ΔνρΔμτ−23ΔμνΔρτ)∇ρUτ+ξΔμν∇ρUρ, (27)

where for any tensor the symmetrized traceless tensor assignment is introduced:

 ∘¯¯¯¯¯¯¯¯aμν≡(ΔμρΔντ+ΔνρΔμτ2−13ΔμνΔρτ)aρτ≡Δμνρτaρτ,Δμν  ρτΔρτ  σλ=Δμν  σλ. (28)

The equation (27) is the definition of the shear and the bulk viscosity coefficients. The term in the (27) can be considered as a nonequilibrium contribution to the pressure which enters in the (20).

By means of the projector (21) one can split the space-time derivative as

 ∂μ=UμUν∂ν+Δνμ∂ν=UμD+∇μ, (29)

where , . In the local rest frame (where ) the becomes the time derivative and the becomes the spacial derivative. Then the Boltzmann equations can be written in the form

 pμk∂μfk=(pμkUμD+pμk∇μ)fk=Celk[fk]+Cinelk[fk], (30)

where represents the inelastic or number-changing collision integrals (it is omitted in calculations in this paper if the opposite is not stated explicitly) and is the elastic collision integral. The collision integral has the form of the sum of positive gain terms and negative loss terms. Its explicit form is383838The factor cancels double counting in integration over momentums of identical particles. The factor comes from the relativistic normalization of the scattering amplitudes. (cf. jeon (); Arnold:2002zm ())

 Celk[fk] = ∑lγkl12∫d3p1l(2π)32p01ld3p′k(2π)32p′0kd3p′1l(2π)32p′01l(f′kf′1l−fkf1l) (31) × |Mkl|2(2π)4δ4(p′k+p′1l−pk−p1l),

where if and denote the same particle species without regard to the spin states and otherwise, is the square of the dimensionless elastic scattering amplitude averaged over the initial spin states and summed over the final ones. Index designates that and are different variables. Introducing as

 Wkl=|Mkl|264π2δ4(p′k+p′1l−pk−p1l), (32)

one can rewrite the collision integral (31) in the form as in groot () (Chap. I, Sec. 2)

 Celk[fk]=(2π)3∑lγkl∫p1l,p′k,p′1l(f′k(2π)3f′1l(2π)3−fk(2π)3f1l(2π)3)Wkl. (33)

The is related to the elastic differential cross section as groot () (Chap. I, Sec. 2)

 Wkl=sσklδ4(p′k+p′1l−pk−p1l), (34)

where is the usual Mandelstam variable. The has properties (due to time reversibility and a freedom of relabelling of order numbers of particles taking part in reaction). And e.g. in the general case. The elastic collision integrals have important properties which one can easily prove groot () (Chap. II, Sec. 1):

 ∫d3pk(2π)3p0kCelk[fk]=0, (35)
 ∑k∫d3pk(2π)3p0kpμkCelk[fk]=0. (36)

Also th