Particle production in Ultra-relativistic Heavy-Ion Collisions : A Statistical-Thermal Model Review

The current status of various thermal and statistical descriptions of particle production in the ultra-relativistic heavy-ion collisions experiments is presented in detail. We discuss the formulation of various types of thermal models of a hot and dense hadron gas (HG) and the methods incorporated in implementing the interactions between hadrons. We first obtain the parameterization of center-of-mass energy ($\sqrt{s_{NN}}$) in terms of temperature ($T$) and baryon chemical potential ($\mu_B$) obtained by analyzing the particle ratios at the freeze-out over a broad energy range from the lowest Alternating Gradient Synchrotron (AGS) energy to the highest Relativistic Heavy-Ion Collider (RHIC) energies. The results of various thermal models together with the experimental results for the various ratios of yields of produced hadrons are then compared. We have derived some new universal conditions emerging at the chemical freeze-out of HG fireball which demonstrate the independence with respect to the energy as well as the structure of the nuclei used in the collision. Further, we perform the calculation of various transport properties of HG such as shear viscosity-to-entropy density ratio ($\eta/s$) etc. using thermal model and compare with the results of other models. We also present the calculation of the rapidity as well as transverse mass spectra of various hadrons in the thermal HG model. The purpose of this review article is to organize and summarize the experimental data obtained in various experiments with heavy-ion collisions and then to examine and analyze them using thermal models so that a firm conclusion regarding the formation of quark-gluon plasma (QGP) can be obtained.


I. INTRODUCTION
One of the main purpose of various heavy-ion collision programmes running at various places such as Relativistic Heavy-Ion Collider (RHIC) at Brookhaven National Laboratory (BNL) and Large Hadron Collider (LHC) at CERN is to understand the properties of strongly interacting matter and to study the possibility of a phase transition from a confined hot, dense hadron gas (HG) phase to a deconfined and/or chiral symmetric phase of quark matter called as quark-gluon plasma (QGP) [1][2][3][4][5][6][7][8]. By colliding heavy-ions at ultrarelativistic energies such a phase transition is expected to materialize and QGP can be formed in the laboratory. Unfortunately, the detection of the QGP phase is still regarded as an uphill task. However, the existence of a new form of a matter called as strongly interacting QGP (sQGP) has been demonstrated experimentally [9]. There are many compelling evidences e. g., elliptic flow, high energy densities and very low viscosity etc. [10]. However, we do not have supportive evidence that this fluid is associated with the properties quark deconfinement and/or chiral symmetry restoration which are considered as direct indication for QGP formation [10]. Although various experimental probes have been devised, but a clean unambiguous signal has not yet been outlined in the laboratory. So our prime need at present is to propose some signals to be used for the detection of QGP. However for this purpose, understanding the behaviour and the properties of the background HG is quite essential because if QGP is not formed, matter will continue to exist in the hot and dense HG phase. In the ultra-relativistic nucleus-nucleus collisions, a hot and dense matter is formed over an extended region for a very brief time and it is often called a 'fireball'. The quark matter in the fireball after subsequent expansion and cooling will be finally converted into HG phase. Recently, the studies of the transverse momentum spectra of dileptons [11][12][13][14][15][16][17] and hadrons [18][19][20][21] are used to deduce valuable information regarding temperature, energy density of the fireball. The schematic diagram for the conjectured space-time evolution of the fireball formed in the heavy-ion collisions is shown in Fig. 1 [22]. The space-time evolution consists of four different stages as follows : (i) In the initial stage of collisions, labeled as "Pre-equilibrium" in Fig. 1, processes of parton-parton hard scatterings may predominantly occur in the overlap region of two colliding nuclei, thus depositing a large amount of energy in the medium. The matter is still not in thermal equilibrium and perturbative QCD models can describe the underlying dynamics as a cascade of freely colliding indicate equilibration. The level of equilibration of particles produced in these experiments is tested by analyzing the particle abundances [18, or their momentum spectra [18-21, 30, 37, 38, 49] using thermal models. Now, in the first case one establishes the chemical composition of the system, while in the second case, additional information on dynamical evolution and collective flow can be extracted. Furthermore, study of the variations of multiplicity of produced particles with respect to collision energy, the momentum spectra of particles, and ratios of various particles have led to perhaps one of the most remarkable results corresponding to high energy strong interaction physics [5].  [50]. But in order to account for the realistic behaviour of HG, a short range repulsion must also be introduced. The importance of such correction is more obvious when we calculate the phase transition using IHG picture which shows the reappearance of hadronic phase as a stable configuration in the simple Gibbs construction of the phase equilibrium between the HG and QGP phases at very high baryon densities or chemical potentials. This anomalous situation [51][52][53] cannot be explained because we know that the stable phase at any given (T, µ B ) is the phase which has a larger pressure. Once the system makes a transition to the QGP phase, it is expected to remain in that phase even at extremely large T and µ B due to the property of asymptotic freedom of QCD.
Moreover, it is expected that the hadronic interactions become significant when hadrons are densely packed in a hot and dense hadron gas. One significant way to handle the repulsive interaction between hadrons is based on a macroscopic description in which the hadrons are given a geometrical size and hence they will experience a hard-core repulsive interaction when they touch each other and consequently a van-der Waals excluded-volume effect is visible. As a result, the particle yields are essentially reduced in comparison to that of IHG model and also the anomalous situation in the phase transition mentioned above disappears.
Recently, many phenomenological models incorporating excluded-volume effect have been widely used to account for the properties of hot and dense HG [54][55][56][57][58][59][60][61][62][63][64][65]. However, these descriptions usually suffer from some serious shortcomings. First, mostly the descriptions are thermodynamically inconsistent because one does not have a well defined partition function or thermodynamical potential (Ω) from which other thermodynamical quantities can be derived e. g. the baryon density (n B ) = ∂Ω/∂µ B . Secondly, for the dense hadron gas, the excluded-volume model violates causality (i.e., velocity of sound c s in the medium is greater than the velocity of light). So, although some of the models explain the data very well but such shortcomings make the models mostly unattractive. Sun et al. [65] have incorporated the effect of relativistic correction in the formulation of an EOS for HG. However, such effect is expected to be very small because the masses of almost all the hadrons present in the hot and dense HG are larger than the temperature of the system, so they are usually treated as non-relativistic particles except pions whose mass is comparable to the temperature but most of the pions come from resonances whose masses are again larger than the temperature of HG [66]. In Ref. [67], two-source thermal model of an ideal hadron gas is used to analyze the experimental data on hadron yields and ratios. In this model, the two sources, a central core and a surrounding halo, are in local equilibrium at chemical freeze-out. It has been claimed that the excluded-volume effect becomes less important in the formulation of EOS for hadron gas in the two-source model.
Another important approach used in the formulation of an EOS for the HG phase is meanfield theoretical models [68,69] and their phenomenological generalizations [70][71][72]. These models use the local renormalizable Lagrangian densities with baryonic and mesonic degrees of freedom for the description of HG phase. These models rigorously respect causality.
Most importantly they also reproduce the ground state properties of the nuclear matter in the low-density limit. The short-range repulsive interaction in these models arises mainly due to ω-exchange between a pair of baryons. It leads to the Yukawa potential V (r) = (G 2 /4 π r) exp(−m ω r) , which further gives mean potential energy as U B = G 2 n B /m ω . It means that U B is proportional to the net baryon density n B . Thus U B vanishes in the n B → 0 limit. In the baryonless limit, hadrons (mesons) can still approach point-like behaviour due to the vanishing of the repulsive interactions between them. It means that in principle one can excite a large number of hadronic resonances at large T . This will again make the pressure in the HG phase larger than the pressure in the QGP phase and the hadronic phase would again become stable at sufficiently high temperature and the Gibbs construction can again yield HG phase at large T . In some recent approaches this problem has been cured by considering another temperature dependent mean-field U V DW (n, T ), where n is the sum of particle and anti-particle number densities. Here U V DW (n, T ) represents van-der Waals hardcore repulsive interaction between two particles and depends on the total number density n and is non-zero even when net baryon density n B is zero in the large temperature limit [62,64]. However, in the high-density limit, the presence of a large number of hyperons and their unknown couplings to the mesonic fields generates a large amount of uncertainty in the EOS of HG in the mean-field approach. Moreover, the assumption of how many particles and their resonances should be incorporated in the system is a crucial one in the formulation of EOS in this approach. The mean-field models can usually handle very few resonances only in the description of HG and hence are not as such reliable [64].
In this review, we discuss the formulation of various thermal models existing in literature quite in detail and their applications to the analysis of particle production in ultra-relativistic heavy-ion collisions. We show that it is important to incorporate the interactions between hadrons in a consistent manner while formulating the EOS for hot, dense HG. For repulsive interactions, van-der Waals type of excluded-volume effect is often used in thermal models while resonances are included in the system to account for the attractive interactions. We precisely demand that such interactions must be incorporated in the models in a thermodynamically consistent way. There are still some thermal models in the literature which lack thermodynamical self-consistency. We have proposed a new excluded-volume model where an equal hard-core size is assigned to each type of baryons in the HG while the mesons are treated as pointlike particles. We have successfully used this model in calculating various properties of HG such as number density, energy density etc.. We have compared our results with those of the other models. Further, we have extracted chemical freeze-out parameters in various thermal models by analyzing the particle-ratios over a broad energy range and parameterized them with the center-of-mass energy. We use these parameterizations to calculate the particle-ratios at various center-of-mass energies and compare them with the experimental data. We further provide a proposal in the form of freeze-out conditions for a unified description of chemical freeze-out of hadrons in various thermal models. An extension of the thermal model for the study of the various transport properties of hadrons will also be discussed. We also analyze the rapidity as well as transverse mass spectra of hadrons using our thermal model and examine the role of any flow existing in the medium by matching the theoretical results with the experimental data. Thus the thermal approach indeed provides a very satisfactory description of various features of HG by reproducing a large number of experimental results covering wide energy range from Alternating Gradient Synchrotron (AGS) to Large Hadron Collider (LHC) energy.

II. FORMULATION OF THERMAL MODELS
Various types of thermal models for HG using excluded-volume correction based on vander Waals type effect have been proposed. Thermal models have often used the grand canonical ensemble description to write the partition function for the system because it suits well for systems with large number of produced hadrons [24] and/or large volume. However, for nonrelativistic statistical mechanics, the use of a grand canonical ensemble is usually just a matter of convenience [73]. Furthermore, the canonical ensemble can be used in case of small systems (e. g. peripheral nucleus-nucleus collisions) and for low energies in case of strangeness production [74] due to canonical suppression of the phase space. Similarly some descriptions also employ isobaric partition function in the derivation of their HG model. We succinctly summarize the features of some models as follows :

A. Hagedorn Model
In the Hagedorn model [58], it is assumed that the excluded-volume correction is proportional to the pointlike energy density ǫ 0 . It is also assumed that the density of states of the finite-size particles in total volume V can be taken as precisely the same as that of pointlike particles in the available volume ∆ where ∆ = V − Σ i V 0 i where V 0 i is the eigen volume of the i th particle in the HG. Thus, the grand canonical partition function satisfies the relation They further modified the canonical partition function Z C by introducing a step-function in the volume so as to incorporate excluded-volume correction into the formalism. Therefore, the grand canonical partition function (11) finally takes the form : The above ansatz is motivated by considering N particles with eigen-volume V 0 in a volume V as N pointlike particles in the available volume V − V 0 N [59]. But, the problem in the calculation of above Eq. (12) is the dependence of the available volume on the varying number of particles N [59]. To overcome this difficulty one should use Laplace transformation of Eq. (12). Using the Laplace transform, one gets the isobaric partition function as : or, : where x = V − V 0 N andμ = µ − T V 0 ξ. Finally, we get a transcendental type of equation as follows : where,μ The expressions for the number density, entropy density, and energy density in this model can thus take a familiar form like : These equations resemble with equations (5), and (7) as given in Cleymans-Suhonen model [54] with µ replaced byμ. The above model can be extended for a hadron gas involving several baryonic species as follows : where,μ with i = 1, · · · , h. Particle number density for the i th species can be calculated from following equation : Unfortunately, the above model involves cumbersome, transcendental expressions which are usually not easy to calculate. Furthermore, this model fails in the limit of µ B = 0 becausẽ µ B becomes negative in this limit.

D. New Excluded-Volume Model
Singh et al., [60] have proposed a thermodynamically consistent excluded-volume model in which the grand canonical partition function using Boltzmann approximation can be written as follows : where g i and λ i = exp(µ i /T ) are the degeneracy factor and the fugacity of i th species of baryons, respectively. k is the magnitude of the momentum of baryons. V 0 i is the eigenvolume assigned to each baryon of i th species and hence j N j V 0 j becomes the total occupied volume where N j represents the total number of baryons of j th species. We can rewrite Eq.
(23) as : where integral I i is Thus we have obtained the grand canonical partition function as given by Eq. (24) by incorporating the excluded volume effect explicitly in the partition function. The number density of baryons after excluded-volume correction (n ex i ) can be obtained as : So our prescription is thermodynamically consistent and it leads to a transcendental equation : Here R = i n ex i V 0 i is the fractional occupied volume. It is clear that if we put the factor ∂R/∂λ i = 0 and consider only one type of baryons in the system then the above Eq.
For a multi-component hadron gas, Eq. (27) takes the following form : Using the method of parametric space [75], we write : where a i is the parameter and t gives the space. We finally get the solution of Eq. (29) as follows : where t is a parameter such that : and, If λ i ' s and t are known, one can determine a i . The quantity t is fixed by setting a 1 = 0 and one obtains t = 1/I 1 V 1 λ 1 , here the subscript 1 denotes the nucleon degree of freedom and h is the total number of baryonic species. Hence by using R and ∂R/∂λ i one can calculate n i . It is obvious that the above solution is not unique. Since it contains some parameters such as a i , one of which has been fixed to zero arbitrarily. Alternatively, one can assume that [60] : Here an assumption is made that the number density of i th baryon will only depend on the fugacity of same baryon. Then the Eq. (29) reduces to a simple form as : The solution of Eq. (35) can be obtained in a straight forward manner as [60] : where Now, R can be obtained by using the following relation : where Here X is the ratio of the occupied to the available volume. Finally, n ex i can be written as: The solution obtained in this model is very simple and easy. There is no arbitrary parameter in this formalism so it can be regarded as a unique solution. However, this theory still depends crucially on the assumption that the number density of i th species is a function of the In order to remove above discrepancies, Mishra et al. [26] have proposed a thermodynamically consistent EOS for a hot and dense HG using Boltzmann's statistics. They have proposed an alternative way to solve the transcendental Eq. (27). We have extended this model by using quantum statistics into the grand canonical partition function so that our model works even for the extreme values of temperature and baryon chemical potential.
Thus the Eq. (25) can be rewritten as follows [25]: and the Eq. (27) takes the following form after using the quantum statistics in the partition function : where I ′ i is the partial derivative of I i with respect to λ i . We can write R in an operator equation form as follows [26,76] : i is the density of pointlike baryons of ith species and the operatorΩ has the form : Using Neumann iteration method and retaining the series uptoΩ 2 term, we get After solving Eq. (45), we finally get the expression for total pressure [60] of the hadron gas as : Here P meson j is the pressure due to jth type of meson.
Here we emphasize that we consider the repulsion arising only between a pair of baryons and/or antibaryons because we assign each of them exclusively a hard-core volume. In order to make the calculation simple, we have taken an equal volume V 0 = 4πr 3 /3 for each type of baryons with a hard-core radius r = 0.
where S i is the strangeness quantum number of the ith hadron, and n s i (n s i ) is the strange (anti-strange) hadron density, respectively. Using this constraint equation, we get the value of strange chemical potential in terms of µ B . Having done all these things, we proceed to calculate the energy density of each baryon species i by using the following formula : Similarly entropy density is : It is evident that this approach is more simple in comparison to other thermodynamically consistent excluded-volume approach which often possesses transcendental final expressions [59,73]. Our approach does not involve any arbitrary parameter in the calculation. Moreover, this approach can be used for extremely low as well as extremely large values of T and µ B where all other approaches fail to give a satisfying result since we do not use Boltzmann's approximation in our calculation.

III. STATISTICAL AND THERMODYNAMICAL CONSISTENCY
Recently, question of statistical and thermodynamical consistency in excluded-volume models for HG has widely been discussed [66]. In this section, we re-examine the issue of thermodynamical and statistical consistency of the excluded-volume models. In RGSG model [59], the single particle grand canonical partition function (12) can be rewritten as follows : where, Here, in this model N in the available volume (V −V 0 N) is independent of µ. Differentiating Eq. (49) with respect to µ, we get the following Eq. : Multiplying both sides of above Eq. (51) by T /Z ex G , we get : We know that the expressions for statistical and thermodynamical averages of number of baryons are : and respectively. Using Equations (53) and (54) in Eq. (52), we get [66] : Now, we calculate the statistical and thermodynamical averages of number of baryons in our excluded-volume model. The grand canonical partition function in our model (i.e. Eq. (23)) can take the following form : where z is given by Eq. (50). We use Boltzmann's statistics for the sake of convenience and consider only one species of baryons. In our model, N, present in the available volume However, for multicomponent system one cannot use "fixed N", because in this case, van-der Waals approximation is not uniquely defined [77]. So, we use average N in our multicomponent grand partition function. Although, at high temperatures it is not possible to use one component van-der Waals description for a system of various species with different masses [77]. Now differentiating Eq. (57) with respect to µ, we get : Multiplying both sides of above Eq. (58) by T /Z ex , we get : Using the definitions (53) and (54), Eq. (59) can take the following form : or, Here n is the thermal average of number density of baryons, n 0 is the number density of pointlike baryons and : The second term in Eq. (61) is the redundant one and arises because N , present in the , is a function of µ. We call this term as 'correction term'. In etc. are approximately equal to each other in our model also. Thus, our excluded-volume model is not exactly thermodynamically consistent but it can safely be taken as consistent because the correction term in the averaging procedure appears as negligibly small.

DATA
In this section, we review various features of hadron gas and present our comparisons between the results of various HG models and the experimental data.

A. Chemical Freeze-out Criteria
Thermal models provide a systematic study of many important properties of hot and dense HG at chemical freeze-out (where inelastic interactions cease). To establish a relation between chemical freeze-out parameters (T, µ B ) and √ s N N , a common method used is to fit the experimental hadron ratios. Many papers [24,40,41,[78][79][80] have appeared in which T and µ B are extracted in terms of √ s N N . In Ref. [24], various hadron ratios are analyzed from √ s N N =2.7 GeV to 200 GeV and chemical freeze-out parameters are parameterized in terms of √ s N N by using following expressions : and, Here, a, b, c, d, e, and T lim (the limiting temperature), are fitting parameters. Various authors [78,79] have included the strangeness suppression factor (γ s ) in their model while extracting the chemical freeze-out parameters. In the thermal model, γ s is used to account for partial equilibration of the strange particles. Such situation may arise in the elementary p − p collisions and/or peripheral A − A collisions and mostly the use of γ s is warranted in such cases [78,81]. Moreover, γ s ≈ 1 has been found in the central collisions at RHIC [82,83]. We do not require γ s in our model as an additional fitting parameter because we assume strangeness is also fully equilibrated in the HG. Also, it has been pointed out that inclusion of γ s in thermal model analysis does not affect the values of fitting parameters T and µ B much [24]. Dumitru et al. [80] have used inhomogeneous freeze-out scenario for the extraction of T and µ B at various √ s N N . In a recent paper [84], condition of vanishing value of κ σ 2 or equivalently m 4 = 3χ 2 are used to describe the chemical freeze-out line where κ, σ, m 4 , and χ are kurtosis, the standard deviation, fourth order moment, and susceptibility, respectively. In Ref. [85], first time experimental data on κ σ 2 and Sσ, here S is skewness, has been compared with the lattice QCD calculations and hadron resonance gas model to determine the critical temperature (T c ) for the QCD phase transition. Recently, it is shown that the freeze-out parameters in heavy-ion collisions can be determined by comparing the lattice QCD results for the first three cumulants of net electric charge fluctuations with the experimental data [86]. In Fig. 3, we have shown the energy dependence of thermal parameters T and µ B extracted by various authors. In all the studies, similar behaviour is found except in the Letessier et al. [79], which may be due to usage of many additional free parameters such as light quark occupancy factor (γ q ), an isospin fugacity etc.. We have also extracted freeze-out parameters by fitting the experimental particle-ratios from the lowest SIS energy to the highest RHIC energy using our model [25]. For comparison, we have shown the values obtained in other models, e. g., IHG model, Cleymans-Suhonen model and RGSG model in Table I. We then parameterize the variables T and µ B in terms of √ s N N as follows [87] : where the parameters a,b,c,d and e have been determined from the best fit : a = 1.482 ± and e = 0.015 ± 0.01 GeV −3 . The systematic error of the fits can be estimated via quadratic deviation, δ 2 [24] defined as follows :     [25,54,59]. Points are the experimental data [88][89][90][91][92]. RHIC data are at midrapidity. out by C. Alt et al. for K + /π + and K − /π − ratios [88]. In Fig. 4, we have shown the parametrization of the freeze-out values of baryon chemical potential with respect to √ s N N and similarly in Fig. 5, we have shown the chemical freeze-out curve between temperature and baryon chemical potential [25].

B. Hadron ratios
In an experimental measurement of various particle ratios at various centre-of-mass energies [93][94][95][96][97], it is found that there is an unusual sharp variation in the Λ/π − ratio increasing upto the peak value. This strong variation of the ratio with energy indicates the critical temperature of QCD phase transition [98] between HG and QGP [99,100] and a nontrivial information about the critical temperature T C ≈ 176 MeV has been extracted [100]. Fig. 6 shows the variation of Λ/π − with √ s N N . We compare the experimental data with various thermal models [25,54,59] and find that our model calculation gives much better fit to the experimental data in comparison to other models. We get a sharp peak around centre-of-mass energy of 5 GeV and our results thus almost reproduce all the features of the experimental data.
In Fig. 7, we have shown the variations of φ/π ratio with √ s N N . The φ yields in the The energy dependence of φ/π ratio. Lines are the results of various thermal models [25,54,59]. Points are the experimental data [88][89][90][91][92]. RHIC data are at midrapidity. thermal models are often much higher in comparison to data. We notice that no thermal model can suitably account for the multiplicity-ratio of multi strange particle since φ is ss hidden-strange quark combination. However, quark coalescence model assuming a QGP formation has claimed to explain the results [101] successfully. In the thermal models, the results for the multistrange particles raises doubt over the degree of chemical equilibration for strangeness reached in the HG fireball. We can use an arbitrary parameter γ s as used in several models. The failures of thermal models in these cases may indicate the presence of QGP formation but it is still not clear. In Fig. 8, we have shown the energy dependence of antiparticle to particle ratios e. g. K − /K + ,p/p,Λ/Λ, andΞ + /Ξ − . These ratios increase sharply with respect to √ s N N and then almost saturate at higher energies reaching the value equal to 1.0 at LHC energy. On comparison with the experimental data we find that almost all the thermal models describe these data successfully over all center-of-mass energies. However, RGSG model [59] fails to describe the data at SPS and RHIC energies in comparison to other models [25].

C. Thermodynamical Properties
We present the thermal model calculations of various thermodynamical properties of HG such as entropy per baryon (s/n B ) and energy density etc. and compare the results with the predictions of a microscopic model URASiMA event generator developed by Sasaki [102]. URASiMA (ultra-relativistic AA collision simulator based on multiple scattering algorithm) is a microscopic model which includes the realistic interactions between hadrons.
In URASiMA event generator, molecular-dynamical simulations for a system of a HG are It is certainly encouraging to find an excellent agreement between the results obtained with two widely different approaches.   causality or not. In Fig. 11, we have plotted the variations of the total hadronic pressure P as a function of the energy density ǫ of the HG at a fixed entropy per particle using our model calculation [25]. We find for a fixed s/n, the pressure varies linearly with respect to energy density. In Fig. 12, we have shown the variation of c s (c s 2 = ∂P/∂ǫ at fixed s/n) with respect to s/n. We find that c s ≤ 0.58 in our model with interacting particles. We get c s = 0.58 (i.e. 1/ √ 3) for an ideal gas consisting of ultra-relativistic particles. This feature endorses our viewpoint that our model is not only thermodynamically consistent but it does not involve any violation of causality even at large density. Similarly in RGSG model [59], we do not notice that the value of c s exceeds 1 as shown in Fig. 12. It should be mentioned that we are using full quantum statistics in all the models taken for comparisons here. However, we find that the values in the RGSG model cannot be extracted when temperature of the HG exceeds 250 MeV . No such restriction applies for our model.

E. Universal freeze-out criteria
One of the most remarkable success of thermal models is in explaining the multiplicities and the particle ratios of various particles produced in heavy-ion experiments from the lowest SIS energy to maximum LHC energy. Some properties of thermal fireball are found Further, Cleymans et al. [87] have found that all the above conditions separately give a satisfactory description of the chemical freeze-out parameters T and µ B in an IHG picture only. Moreover, it was also found that these conditions are independent of collision energy and the geometry of colliding nuclei. Furthermore, Cleymans et al. [87] have hinted that incorporation of excluded-volume correction leads to wild as well as disastrous effects on these conditions. The purpose in this section is to reinvestigate the validity of these freezeout criteria in excluded-volume models. Along with these conditions, a condition, formulated by using percolation theory is also proposed as a chemical freeze-out condition [104]. An assumption is made that in the baryonless region the hadronic matter freezes-out due to hadron resonances and vacuum percolation, while in the baryon rich region the freeze-out takes place due to baryon percolation. Thus the condition which describes the chemical freeze-out line is formulated by following equation [104] : where V h is the volume of a hadron. The numbers 1.24 and 0.34 are obtained within percolation theory [106].
In Fig. 13, we have shown the variation of E/N with respect to √ s N N at the chemical freeze-out point of the fireball. The ratio E/N shows a constant value of 1.0 in our model and it shows also a remarkable energy independence. Similarly the curve in IHG model shows that the value for E/N is slightly larger than one as reported in [87]. However, results support the finding that E/N is almost independent of energy and also of the geometry of the nuclei. Most importantly we notice that the inclusion of the excluded-volume correction does not change the result much which is contrary to the claim of Cleymans et al. [87].
The condition E/N ≈ 1.0GeV was successfully used in the literature to make predictions [36] of freeze-out parameters at SPS energies of 40 and 80 A GeV for Pb-Pb collisions long before the data were taken [25]. Moreover, we have also shown in Fig. 13, the curves in the Cleymans-Suhonen model [54] and the RGSG model [59] and we notice a small variation with √ s N N particularly at lower energies. In Fig. 14 Our finding lends support to the crucial assumption of HG fireball achieving chemical equilibrium in the heavy-ion collisions from the lowest SIS to RHIC energy and the EOS of the HG developed by us indeed gives a proper description of the hot and dense fireball and its subsequent expansion. However, we still do not get any information regarding QGP formation from these studies. The chemical equilibrium once attained by the hot and dense HG removes any memory regarding QGP existing in the fireball before HG phase.
Furthermore, in a heavy-ion collision, a large amount of kinetic energy becomes available and part of it is always lost during the collision due to dissipative processes. In thermal description of the fireball, we ignore the effect of such processes and we assume that all available kinetic energy (or momentum) is globally thermalized at the freeze-out density.
Experimental configuration of the collective flow developed in the hot, dense matter reveals the unsatisfactory nature of the above assumption. There were several analytic calculations for η and η/s of simple hadronic systems [113][114][115][116][117][118][119] along with some sophisticated microscopic transport model calculations [120][121][122] in the literature. Furthermore, some calculations predict that the minimum of shear viscosity to entropy density is related with the QCD phase transition [123][124][125][126][127]. Similarly sound velocity is an important property of the matter created in heavy ion collision experiments because the hydrodynamic evolution of this matter strongly depends on it. A minimum in the sound-velocity has also been interpreted in terms of a phase transition [23,124,[128][129][130][131][132][133] and further, the presence of a shallow minimum corresponds to a cross-over transition [134].
In view of the above, it is worthwhile to study in detail the transport properties of the HG in order to fully comprehend the nature of the matter created in the heavy-ion collisions as well as the involved phase transition phenomenon. In this section, we have used thermal models to calculate the transport properties of HG such as shear viscosity to entropy ratio etc. [25].
We calculate the shear viscosity in our thermal model as was done previously by Gorenstein et al. [135] using RGSG model. According to molecular kinetic theory, we can write the dependence of the shear viscosity as follows [136] : where n is the particle density, l is the mean free path, and p is the average thermal momentum of the baryons or antibaryons. For the mixture of particle species with different masses and with the same hard-core radius r, the shear viscosity can be calculated by the following equation [135] : where n i is the number density of the ith species of baryons (anti-baryons) and n is the total baryon density. In Fig. 15, we have shown the variation of η/s with respect to temperature as obtained in our model for HG having a baryonic hard-core size r = 0.5 fm, and compared our results with those of Gorenstein et al. [135]. We find that near the expected QCD phase In Fig. 16, we have shown the variation of η/s with respect to µ B at a very low temperature (≈ 10 MeV) [25]. Here we find that the η/s is constant as µ B increases upto 700 MeV and then sharply decreases. This kind of valley-structure at low temperature and at µ B around 950 MeV was also obtained by J. W. Chen et al. [116] and K. Itakura et al. [118].
They have related this structure to the liquid-gas phase transition of the nuclear matter.
As we increase the temperature above 20 MeV, this valley-like structure disappears. They further suspect that the observation of a discontinuity in the bottom of η/s valley may correspond to the location of the critical point. Our HG model yields a curve in complete agreement with these results. Fig. 17 represents the variation of η and η/s with respect to temperature at a fixed µ B (= 300 MeV), for HG having a baryonic hard-core size r = 0.8 fm. We have compared our result with the result obtained in Ref. [118]. Here we find that η increases with temperature in our HG model as well as in the simple phenomenological calculation [118], but decreases with increasing temperature in low temperature effective field theory (EFT) calculations [116,118]. However, η/s decreases with increasing temperature in all three calculations and η/s in our model gives the lowest value at all the temperatures in comparison to other models. In Fig. 18, we have shown a comparison between η calculated in our HG model with the results obtained in a microscopic pion gas model used in Ref. predicts lower values of η in comparison to our results. The most probable reason may be that the calculations have been done only for pion gas in the microscopic model while at low temperatures the inclusion of baryons in the HG is very important in order to extract a correct value for the shear viscosity [25]. Figure 19 shows the variation of η/s with respect to √ s N N in our model calculation. We have compared our results with that calculated in [135]. There is similarity in our results at lower energies while our results significantly differ at higher energies. Although both the calculations show that η/s decreases with increasing √ s N N .
The study of the transport properties of non-equilibrium systems which are not far from an equilibrium state has yielded valuable results in the recent past. Large values of the elliptic flow observed at RHIC indicates that the matter in the fireball behaves as a nearly perfect liquid with a small value of the η/s ratio. After evaluating η/s in strongly coupled theories using AdS/CFT duality conjecture, a lower bound was reported as η/s = 1 4π . We surprisingly notice that the fireball with hot, dense HG as described in our model gives transport coefficient which agree with those given in different approaches. Temperature and baryon chemical potential dependence of the η/s are analyzed and compared with the results obtained in other models. Our results lend support to the claim that knowledge of the EOS and the transport coefficients of HG is essential for a better understanding of the dynamics of the medium formed in the heavy-ion collisions [25].  [25] and other model calculation [135].

G. Rapidity and Transverse Mass Spectra
In order to suggest any unambiguous signal for QGP, the dynamics of the collisions should be understood properly. Such information can be obtained by analyzing the properties of various particles which are emitted from various stages of the collisions. Hadrons are produced at the end of the hot and dense QGP phase, but they subsequently scatter in the confined hadronic phase prior to decoupling (or "freeze-out") from the collision system and finally a collective evolution of the hot and dense matter occurs in the form of transverse, radial or elliptic flow which are instrumental in shaping the important features of the particle spectra. The global properties and dynamics of freeze-out can be at best studied via hadronic observables such as rapidity distributions and transverse mass spectra etc. [139]. There are various approaches for the study of rapidity as well as transverse mass spectra of HG [22,30,[140][141][142][143][144][145][146][147][148][149][150][151][152][153][154][155][156][157]. Hadronic spectra from purely thermal models usually reveal an isotropic distribution of particles [158] and hence the rapidity spectra obtained with the purely thermal models do not reproduce the features of the experimental data satisfactorily. Similarly the transverse mass spectra from the thermal models reveal a more steeper curve than that observed experimentally. The comparisons illustrate that the fireball formed in heavy-ion collisions does not expand isotropically in nature and there is a prominent input of collective flow in the longitudinal and transverse directions which finally causes anisotropy in the rapidity and transverse mass distributions of the hadrons after the freeze-out. Here we mention some kinds of models of thermal and collective flow used in the literature. Hydrodynamical properties of the expanding fireball have been initially discussed by Bjorken and Landau for the central-rapidity and stopping regimes, respectively [22,140].
However, collisions even at RHIC energies reveal that they are neither fully stopped, nor fully transparent. As the collision energy increases, the longitudinal flow grows stronger and leads to a cylindrical geometry as postulated in Ref. [30,141,142]. They assume that the fireballs are distributed uniformally in the longitudinal direction and demonstrate that the available data can consistently be described in a thermal model with inputs of chemical equilibrium and flow, although they have also used the experimental data for small systems only. They use two simple parameters : transverse flow velocity (β r ) and temperature (T ) in their models. In Ref. [143], non-uniform flow model is used to analyze the spectra specially to reproduce the dip at midrapidity in the rapidity spectra of baryons by assuming that the fireballs are distributed non-uniformly in the longitudinal phase space. In Ref.
[ [149][150][151], the rapidity-dependent baryon chemical potential has been invoked to study the rapidity spectra of hadrons. In certain hydrodynamical models [152], In this section, we study the rapidity and transverse mass spectra of hadrons using thermal approach. We can rewrite Eq. (42) in the following manner [49] : It means that the invariant distributions are [30,141] : If we use Boltzmann's approximation, Eq. (71) differs from the one used in the paper of Schnedermann et al. [141] by the presence of a prefactor (1 − R) − λ i ∂R ∂λ i . However, we measure all these quantities precisely at the chemical freeze-out using our model and hence quantitatively we do not require any normalizing factor as is required in Ref. [141]. We use the following expression to calculate the rapidity distributions of baryons in the thermal model [49] : Here y is the rapidity variable and m T is the transverse mass (m T = m 2 + p T 2 ). Also V is the total volume of the fireball formed at chemical freeze-out and N i is the total number of ith baryons. We assume that the freeze-out volume of the fireball for all types of hadrons at the time of the homogeneous emissions of hadrons remains the same. It can be mentioned here that in the above equation, there occurs no free parameter because all the quantities g, V , λ, R etc. are determined in the model. Although, the above Eq. (72) describes the experimental data only at midrapidity while it fails at forward and backward rapidities, so we need to modify it by incorporating a flow factor in the longitudinal direction. Thus the resulting rapidity spectra of ith hadron is [30,49,141]: where dN i dy th can be calculated by using Eq. (72). The expression for average longitudinal velocity is [143,161]: Here η max is a free parameter which provides the upper rapidity limit for the longitudinal flow velocity at a particular √ s N N and it's value is determined by the best experimental fit. The value of η max increases with the increasing √ s N N and hence β L also increases.
Cleymans et al. [151] have extended the thermal model [149], in which the chemical freezeout parameters are rapidity-dependent, to calculate the rapidity spectra of hadrons. They use the following expression for rapidity spectra : where dN i 1 dy is the thermal rapidity distribution of particles calculated by using Eq. (72) and ρ(y F B ) is a Gaussian distribution of fireballs centered at zero and given by : Similarly we calculate the transverse mass spectra of hadrons by using following expression [49] : where K 1 m T T is the modified Bessel's function : Above expression for transverse mass spectra arises from a stationary thermal source alone which is not capable of describing the experimental data successfully. So, we incorporate flow velocity in both the directions in the Eq. (77), longitudinal as well as transverse, in order to describe the experimental data satisfactorily. After defining the flow velocity field, we can calculate the invariant momentum spectrum by using the following formula [141,162] : While deriving Eq.(79), we assume that the local fluid velocity u µ give a boost to an isotropic thermal distribution of hadrons. Now the final expression of transverse mass spectra of hadrons after incorporation of flow velocity in our model is [49] : Here I 0 p T sinhρ T is the modified Bessel's function : where ρ is given by ρ = tanh −1 β r , with the velocity profile chosen as β r = β s ξ n [30,141].
β s is the maximum surface velocity and is treated as a free parameter and ξ = r/R 0 . The average of the transverse velocity can be evaluated as [154] : < β r >= β s ξ n ξ dξ ξ dξ = 2 2 + n β s .
In our calculation we use a linear velocity profile, (n = 1) and R 0 is the maximum radius of the expanding source at freeze-out (0 < ξ < 1) [154]. The symbols are the HBT data for freeze-out volume V HBT for the π + [163]. A ′ ,B ′ and C ′ are the total freeze-out volume and A ,B and C depict the dV /dy as found in our model for π + , K + and K − , respectively. D represents the total freeze-out volume for π + calculated in the Ideal HG model. D ′ is the the total freeze-out volume for π + in our model calculation using Boltzmann's statistics [49].
In Fig. 20, we have shown the variations of V and dV /dy with the √ s N N calculated in our excluded-volume model and compared with the results of various thermal models.
We show the total freeze-out volume for π + calculated in our model using Boltzmann's statistics. We see that there is a significant difference between the results arising from quantum statistics and Boltzmann's statistics [49]. We also show the total freeze-out volume for π + in IHG model calculation by dash-dotted line D. We clearly notice a remarkable difference between the results of our excluded-volume model and that of IHG model also.
We have also compared predictions from our model with the data obtained from the pion interferometry (HBT) [163] which in fact reveals thermal (kinetic) freeze-out volumes. The results of thermal models support the finding that the decoupling of strange mesons from the fireball takes place earlier than the π-mesons. Moreover, a flat minimum occurs in the curves around the center-of-mass energy ≈ 8 GeV and this feature is well supported by HBT  GeV . Dashed and dotted lines are the transverse mass spectra due to purely thermal source for π + and proton, respectively. Solid and dash-dotted lines are the results for π + and proton, respectively obtained after incorporation of flow in thermal model [49]. Symbols are the experimental data [165].
data. In Fig. 21, we present the rapidity distribution of π + for central Au + Au collisions at √ s N N = 200 GeV over full rapidity range. Dotted line shows the distribution of π + due to purely thermal model. Solid line shows the rapidity distributions of π + after the incorporation of longitudinal flow in our thermal model and results give a good agreement with the experimental data [164]. In fitting the experimental data, we use the value of We find that the results slightly differ and hence it shows a small dependence on η max [49]. Fig. 22 represents the rapidity distributions of pion at various √ s N N calculated by using Eq. (75) [151]. There is a good agreement between the model results and experimental data at all √ s N N .
In Fig. 23, we show the transverse mass spectra for π + and proton for the most central collisions of Au + Au at √ s N N = 200 GeV . We have neglected the contributions from close agreement with the experimental data [165]. In Fig. 25, we show the p T spectra of π − , K − andp for the P b−P b collisions at √ s N N = 2.76 T eV at the LHC. Our calculations again give a good fit to the experimental results [166]. We also compare our results forp spectrum with the hydrodynamical model of Shen et al. [167], which successfully explains π − , and K − spectra but strongly fails in the case ofp spectrum [49]. In comparison, our model results show closer agreement with the experimental data. Shen et al. [167] have employed (2 + 1)dimensional viscous hydrodynamics with the lattice QCD-based EOS. They use Cooper-Frye prescription to implement kinetic freeze-out in converting the hydrodynamic output into the particle spectra. Due to lack of a proper theoretical and phenomenological knowledge, they use the same parameters for P b − P b collisions at LHC energy, which was used for Au − Au collisions at √ s N N = 200 GeV . Furthermore, they use the temperature independent η/s ratio in their calculation. After fitting the experimental data, we get β s = 0.80 (β r = 0.53) at this energy which indicates the collective flow becoming stronger at LHC energy than that observed at RHIC energies. In this plot, we also attempt to show how the spectra for π − will change at a slightly different value of the parameter i.e., β s = 0.88 [49].

V. SUMMARY AND CONCLUSIONS
The main aim in this article is to emphasize the use of the thermal approach in describing various yields of different particle species that have been measured in various experiments running at various places. We have discussed various types of thermal approaches for the formulation of EOS for HG. We have argued that, incorporation of interactions between hadrons in a thermodynamically consistent way is important for the realistic formulation of HG from both qualitatively and quantitatively point of view. We have presented the systematic study of the particle production in heavy-ion collisions from AGS to LHC energy.
We have observed from this analysis that the production of the particles seems to be oc- In conclusion, the net outcome of this review is indeed a surprising one. The excludedvolume HG models are really successfull in describing all kinds of features of the HG formed in ultra-relativistic heavy-ion collisions. Most important property indicated by such description is the chemical equilibrium reached in such collisions. However, the description is still a geometrical one and does not involve any microscopic picture of interactions. Moreover, it's relativistic and field theoretic generalizations are still needed in order to make the picture a more realistic description. But it is amazing to find that these models still work much better than expected. Notwithstanding these remarks, we should add that Lattice QCD results are now available for the pressure, entropy density and energy density etc. for the entire temperature range from T = 0 to higher values at µ = 0. Here low temperature phase of QCD is the HG and recently our excluded-volume model reproduces these properties quite in agreement with Lattice results [168]. We have also used these calculations in the precise determination of the QCD critical end point [169]. Thus we conclude that the excluded-volume models are successfull in reproducing numerical results obtained in various experiments and, therefore, further research is required to show how these descriptions are connected with the microscopic interactions.