Charged Hadron Multiplicity Distribution at Relativistic Heavy Ion Colliders

The present article reviews facts and problems concerning charge hadron production in high energy collisions. Main emphasis is laid on the qualitative and quantitative description of general characteristics and properties observed for charged hadrons produced in such high energy collisions. Various features of available experimental data e.g., the variations of charged hadron multiplicity and pseudo-rapidity density with the mass number of colliding nuclei, center-of-mass energies and the collision centrality obtained from heavy-ion collider experiments are interpreted in the context of various theoretical concepts and their implications. Finally, several important scaling features observed in the measurements mainly at RHIC and LHC experiments are highlighted in the view of these models to draw some insight regarding the particle production mechanism in heavy-ion collisions.


Introduction
Quantum chromodynamics (QCD), the basic theory which describes the interactions of quarks and gluons is a firmly established microscopic theory in high energy collision physics.
Heavy-ion collision experiments provide us a unique opportunity to test the predictions of QCD and simultaneously to understand the two facets of high energy collision process: hard process (i.e, the small cross-section physics) and soft process (i.e, the large cross-section physics). [1][2][3] Nuclear collisions at very high energies such as collidier energies enable us to study the novel regime of QCD, where parton densities are high and the strong coupling constant between the partons is small which further decreases as the distance between the partons decreases. The parton densities in the initial stage of the collision can be related to the density of charged hadrons produced in the final state. With the increase in collision energy, the role of hard process (minijet and jet production) in final state particle production rapidly increases and offers a unique opportunity to investigate the interplay between various effects. In this scenario, the perturbative QCD (pQCD) lends a good basis for high energy dynamics and has achieved significant success in describing hard processes occuring in 1/50 high energy collisions such as scaling violation in deep inelastic scattering DIS, 4 hadronic-jet production in e + e − annihilation, 5,6 large-pt-jet production in hadron-hadron collisions, [7][8][9][10][11] etc. On the other hand, in soft processes such as hadron production with sufficiently small transverse momentum in hadronic and nuclear collisions, the interactions become so strong that the perturbative QCD (pQCD) does not remain applicable any more. Thus, there is no workable theory yet for non-perturbative QCD regime which can successfully describe these soft processes. Due to inapplicability of pQCD in this regime, experimental input-based phenomenological models are proven to be an alternative tool to increase our knowledge on the property of the basic dynamics involved in such collision processes. Furthermore, these soft hadrons which decouple from the collision zone in the late hadronic freezeout stage of the evolution are quite useful in providing the indirect information about the early stage of the collision. Several experimental information on the multiparticle production in lepton-hadron, hadron-hadron, hadron-nucleus and nucleus-nucleus have been accumulated in recent past over a wide range of energy. In this context, the bulk features of multiparticle production such as the average charged particle multiplicity and particle densities are of fundamental interest as their variations with the collision energy, impact parameter and the collision geometry are very sensitive to the underlying mechanism involved in the nuclear collisons. These can also throw some more light in providing insight on the partonic structure of the colliding nuclei.
In order to understand the available experimental data, a lot of efforts have been put forward in terms of theoretical and phenomenological models. However, the absence of any well established alternative, the existing problem of the production mechanism of charged hadrons continues to facilitate proliferation of various models. The development in this direction is still in a state of flux for describing the same physical phenomenon using different concepts and modes of operation. Most of these theoretical models are based on the geometrical, hydrodynamical and statistical approaches. However, the diverse nature of the experimental data poses a major challenge before physics community to uncover any systematics or scaling relations which are common to all type of reactions. Thus, a search for universal mechanism of charged hadron production common to hadron-hadron, hadron-nucleus and nucleus-nucleus interactions is still continued and needs a profound effort to draw any firm conclusion. Also, the complicated process of many-body interactions occuring in these collision processes is still quite difficult to make a clear understanding of the phenomena by analyzing the experimental data on the multiparticle production in the final state. In this regard, ongoing efforts for the extensive analysis of the experimental data available on charged hadron production in the view of some successsful phenomenological models can provide us a much needed insight 2/50 in developing a better understanding of the mechanism involved in the particle production.
Moreover, these can also be useful in revealing the properties of the nuclear matter formed at extreme conditions of energy and matter densities.
In this review, we attempt to give a succint description of most of the progress made in this field till date even it is not so easy for us. Furthermore, we believe that the references mentioned in this review will surely guide the readers but we can never claim that they are complete. We apologize to those authors whose valuable contributions in this field have not been properly mentioned.
The structure of this paper is framed in the following manner: At first, In section 2, we start with a brief description of different models used for the study of charged hadron productions in this review in a systematic manner. In section 3, the experimental results on charged hadron production at collider energies are presented along with the comparison of different model results. Further in section 4, we will provide some scaling relations for charged hadrons production and evaluate them on the basis of their universality in different collisions.

Wounded nucleon model
In 1958, Glauber presented his first collection of various papers and unpublished work. 12 Before that, there were no systematic calculations for treating the many-body nuclear system as either a projectile or target. Glauber's work put the quantum theory of collisions of composite objects on a firm basis. In 1969 Czyz and Maximon 13 applied the Glauber's theory in its most complete form for proton-nucleus and nucleus-nucleus collisions. By using Glauber's theory, finally A. Bialas et al. 14 first proposed the wounded nucleon model which was based on the basic assumption that the inelastic collisions of two nuclei can be described as an incoherent composition of the collisions of individual nucleons of the colliding nuclei. In this approach, the collective effects which may occur in nuclei were neglected. According to their assumption, in nucleon-nucleus collisions, a fundamental role is played by the mean number of collisions ν suffered by the incident nucleon with the nucleons in the target nucleus. 15 Similarly, nucleus-nucleus collision is also described in terms of the number of wounded nucleons (w). For the nucleon-nucleus collisions, there is a simple relation between ν and w w = ν + 1 (1) but no such type of relation exist for the nucleus-nucleus collisions. Motivated by the data available on nucleon-nucleus (N − A) interactions, 16 the average multiplicity follows approxi-3/50 mately the formulan wheren H is the average multiplicity in nucleon-nucleon collisions. For nucleus-nucleus (A-B) collisions generalization of this picture implies that the average multiplicity is where the number of wounded nucleons (participants) in the collision of A and B is the sum of wounded nucleons in the nucleus A and the nucleus B i.e.,w AB =w A +w The extension of Glauber model was used to describe elastic, quasi-elastic, the total cross sections. 13,[17][18][19][20][21] In the wounded nucleon model, the basic entity is the nucleonnucleon collision profile P (b) defined as the probability of inelastic nucleon-nucleon collision at impact parameter b. Once the probability of a given nucleon-nucleon interaction is known, the probability of having n such interactions in collision of nuclei A and B is given as: 22 where T AB (b) is thickness function and the σ N N inel is the inelastic nucleon-nucleon cross-section. The total cross-section is given by: The total number of nucleon-nucleon collisions is given by: The number of participants ( wounded nucleons) at a given impact parameter is given by: with is the probability per unit transverse area of a given nucleon being located in the target flux tube of A or B and ρ A(B) (s, z A(B) ) is the probability per unit volume, normalized to uniy for finding the nucleon at location (s, z A(B) ).
The rapidity density of particles in nucleus-nucleus (A − B) collision is given by: 23 where y is the rapidity in c.m. system, F A (y) and F B (y) are the contribution from a single 4/50 wounded nucleon in A and B. Thus: The model gives a good description of the data, with the condition (9) being well satisfied, except at rapidities close to the maximal values. 23 It can be seen from Eq.(8) that for nucleon-nucleon collision, one have: and thus for the ratio one obtains consequently, one have at y=0, which implies that the value of the ratio R AB at mid-rapidity is fully determined by the number of wounded nucleons and independent of the function F (y). Unlike the N part scaling observed in charged hadron multiplicity, pseudo-rapidity density at mid-rapidity does not scale linearly with N part /2. 22 It was conjectured 24-28 that at sufficiently high energy the particle production in nucleus-nucleus collisions will be dominated by hard processes. However, the gross features of particle production at CERN SPS energies were found to be approximately consistent 29 with the N part scaling as accomodated by the wounded nucleon model. A better agreement with the data is found in two-component model for estimating the pseudo-rapidity density in wounded nucleon model as shown by Ref.. 30 In A-A collisions hadron production from the two processes scales as N part /2 (number of nucleon participant pairs) and N bin (number of binary N-N collisions), respectively. 30 According to this assumption the pseudo-rapidity density of charged hadrons is given as: where x quantifies the relative contributions of two components arising from hard and soft processes. The fraction x corresponds to the contribution from hard processes and the remaining fraction (1 − x) describes the contributions arising from the soft processes.

Wounded quark Model
Charged hadron production by using the concept of the constituent or wounded quarks has been widely used for many years. [31][32][33][34] In the wounded quark picture, nucleus-nucleus collisions are effectively described in terms of the effective number of constituent quarks participating in the collision process along with the effective number of collisions suffered by each of them. Recently, the idea of wounded quarks was resurrected by S. Eremin et al. 35 in which they modified the overlap function by increasing the nucleon density three times and introduced one more parameter the quark-quark interaction cross-section, which reproduced the data well. They have further shown that the charged hadron density at mid-rapidity can be described well by the wounded quark model. This problem was further investigated in several papers [36][37][38][39][40] representing the analysis of various spectra, SPS data, total multiplicities and the energy deposition. De et al. 37 have shown that the data on K ± and p/p favors the N q−part scaling over the N part scaling whereas the pions do not agree well with such scaling law. Recently, we proposed a wounded quark model 41 which is primarily based on the previous work by C. P. Singh et al.. [42][43][44] In this picuture, during the collision, a gluon is exchanged between a quark of projectile or first nucleus and a quark belonging to target or other collliding nucleus. The resulting colour force is then somewhat stretched between them and other constituent quarks because they try to restore the colour singlet behaviour.
When two quarks separate, the colour force builds up between them and the energy in the colour-field increases, the colour tubes thus formed finally break-up into new hadrons and/or quark-antiquark pairs. We consider a multiple collision scheme in which a valence quark of the incident nucleon suffers one or more inelastic collisions with the quarks of target nucleons. In a nucleon-nucleon collision only one valence quark of each nucleon (i.e, target and projectile) interacts while other quarks remain as spectators. 40 Only a part of the entire nucleon energy is spent for secondary production at mid-rapidity. The other spectator quarks are responsible for forming hadrons in the nucleon fragmentation region. In the case of nucleus-nucleus collisions, more than one quark per nucleon interacts and each quark suffers more than one collision due to a large nuclear size since large travel path inside the nucleus becomes available. If we search a universal mechanism of charged particle production in the hadron-hadron, hadron-nucleus and nucleus-nucleus collisions, it must be driven by the available amount of energy required for the secondary production, and it also depends on the mean number of participant quarks. The main ingredients of our model are taken from the paper by Singh et. al.. [42][43][44] Charged hadrons produced from A − A collisions are assumed to result from a somewhat unified production mechanism common to p − p collisions at various energies.

6/50
Based on the experimental findings by PHOBOS Collaboration, recently Jeon and collaborators have shown that the total multiplicity obtained at RHIC can be bounded by a cubic logarithmic term in energy. 45 Therefore, we propose here a new parameterization involving a cubic logarithmic term so that the entire p − p experimental data starting from low energies (i.e. from 6.15 GeV) upto the highest LHC energy (i.e. 7 TeV) can suitably be described as: 41 In Eq. (15), α is the leading particle effect and √ s a is the available center-of-mass energy (i.e., The generalization of the above picture for the case of nucleus-nucleus collisions is finally achieved as follows: 41 The parametrization in Eq. (17) thus relates nucleus-nucleus collisions to hadron-nucleus and further to hadron-proton collisions and the values of the parameters a ′ , b ′ , c ′ , and d ′ remain unaltered which shows its universality for all these processess.
In creating quark gluon plasma (QGP), greater emphasis is laid on the central or head-on collisions of two nuclei. The mean multiplicity in central collisions can straight forwardly be 7/50 generalized 41 from Eq. (17) as : The pseudo-rapidity distribution of charged particles is another important quantity in the studies of particle production mechanism from high energy h − h and A − B collisions, which, however, is not yet understood properly. It has been pointed out that (dn ch /dη) can be used to get the information on the temperature (T ) as well as energy density (ρ) of the QGP . [51][52][53] For the pseudo-rapidity density of charged hadrons, we first fit the experimental data of (dn ch /dη) pp η=0 for collision-energy ranging from a low energy to very high energy. One should use the parameterization upto squared logarithmic term in accordance with Ref.. 45 Hence, using a parameterization for central rapidity density as: we obtain the values of the parameters a ′ 1 = 1.24, b ′ 1 = 0.18, and c ′ 1 = 0.044 from the reasonable fit to the p-p data. 41 Earlier many authors have attempted to calculate the pseudo-rapidity density of charged hadrons in a two component model of parton fragmentation. 54,55 Its physical interpretation is based on a simple model of hadron production: longitudinal projectile nucleon dissociation (soft) and transverse large-angle scattered parton fragmentation (hard). However, this assumption which is based on a nucleon-nucleon collision in the Glauber model is crude and it looks unrealistic to relate participating nucleons and nucleon-nucleon binary collisions to soft and hard components at the partonic level. Here we modify the two component model of pseudo-rapidity distributions in A-A collisions as given in Eq. (14) for the wounded quark scenario and assume that the hard component which basically arises due to multiple parton interactions 56 scales with the number of quark-quark collisions (i.e. N AB q ν AB q ) and soft component scales with the number of participating quarks (i.e. N AB q ). Thus, the expression for (dn ch /dη) AB η=0 in A-B collisions can be parameterized in terms of p-p rapidity density as follows: 41 In order to incorporate η dependence in central A-B collisions, we further extend the 8/50 model by using the functional form: where β, γ and σ are fitting parameters and (dn ch /dη) AB η=0 is the central pseudo-rapidity density in A-B collisions obtained from Eq. (20).  62 The inclusive spectra for charged particle multiplicity in p − p collisions is given as follows: 63

Dual Parton Model
when all string contributions are identical means all individual scattering are same then this expression reduces to the simple expression here k is the average number of inelastic collisions and dN pp 0 dη the charged multiplicity per unit pseudo-rapidity in an individual p − p collision.
To calculate the weights σ k for the occurence of k inelastic collisions, a quasi-eikonal model has been used. The σ k is given as follows: 62 Here ξ = ln(s/s o ), σ p = 8πγ p exp(∆ξ) and Z = 2C E γ p exp(∆ξ)/(R 2 + α ′ p ξ). In Eq. (23), σ p is the Born term given by Pomeron exchange with intercept α p (0) = 1 + ∆. According to a well known identity which is known as AGK cancellation, 65 all multiple scattering contributions vanishes identically in the single particle inclusive distribution dσ/dη 9/50 and only the Born term (σ p ) contribution is left. The parameter R 2 and α ′ p control the tdependence of the elastic peak and C E contain the contribution of diffractive intermediate states. The total p − p cross-section in this prescription is: 64 Thus one can calculate the pseudo-rapidity distribution of charged hadron from Eq.(22) and using expectation value of k as follows: Further, the cross-section for ν inelastic collisions in pA scattering after considering the AGK cancellation is as follows: 62 Here T A (b) represents the nucleon profile function as defined in the wounded nucleon model. nucleons. On the contrary the multiplicity for p − A scattering in DPM is given as follows: withν = Aσ pp /σ pA . Thus dN pA dη scales with the number of binary collisions. Further, the AGK cancellation implies that at mid-rapidity dσ AA /dη = A 2 dσ pp /dη, which implies As a function of the impact parameter, the average number of binary nucleon-nucleon collisions n coll (b) can be expressed as follows: 62 Capella and Ferreiro 62 authors have also incorporated the corrections arising due to shadowing effects by including the contribution of triple Pomeron graphs. The suppression in multiplicity from shadowing in A − A collision for a particle produced at mid-rapidity can be obtained by replacing the nuclear profile function where can be obtained as follows:

The Color Glass Condensate Approach
Color Glass Condensate 66-82 (CGC) is an effective theory that describes the gluon content of a high energy hadron or nucleus in the saturation regime. At high energies, the nuclei gets contracted and the gluon density increases inside the hadron wave functions and at small x, the gluon density is very large in comparison to all other parton species (the valence quarks) and the sea quarks are suppressed by the coupling α s , since they can be produced from the gluons by the splitting g → qq. System that evolve slowly compared to natural time scales are generally glasses. The word Color is because CGC is composed of colored gluons. The word Glass is because the classical gluon field is produced by fast moving static sources. The distribution of these sources is real. System that evolves slowly compared to natural time scales are generally glasses. The word condensate means the gluon distribution has maximal phase space density for momentum modes and the strong gluon fields are self generated by the hadron. This effective theory approximates the description of the fast partons in the wave function of a hadron. This framework has been applied in a range of experiments e.g., from DIS to proton-proton, proton-nucleus and nucleus-nucleus collisions. One of the early success of the CGC was the description of multiplicity distributions in DIS experiments. 53 When the nucleus is boosted to a large momentum-then due to Lorentz contraction in the transverse plane of nuclei, partons have to live on thin sheet in the transverse plane. Each parton occupies the transverse area π/Q 2 and can be probed with the cross section σ ∼ α s (Q 2 )π/Q 2 . On the other side the total transverse area of the nucleus is S A ∼ πR 2 A . Therefore if the number of partons exceeds then they start to overlap in the transverse plane and start interacting with each other which prevents further growth of parton densities. At this situation the transverse momenta of the partons are of the order of Q s The multiplicity of the produced partons should be proportional to: 30 In the first approximation, the multiplicity in this high density regime scales with the number of participants. However, there is an important logarithmic correction to this from the 11/50 evolution of parton structure functions with Q 2 s . The coefficient of proportionality is given as follows: 30 where N c =3 is the number of colour, xG(x, Q 2 s ) is the gluon structure function in nucleus and ρ part is the density of participants in the transverse plane.
Number of produced partons is given as follows: 30 where c is the "parton liberation" coefficient accounting for the transformation of virtual partons in the initial state to the on-shell partons in the final state. Integrating over transverse coordinate and using Eq.(35), one can obtain: 30 Here c is assumed to be close to unity in the context of local "parton hadron duality" hypothesis 84 according to which at p ⊥ ∼ Q s , the distribution of produced hadrons mirrors that of the produced gluons. Using Eq. (35) and (37), one can evaluate the centrality dependence as follows: 30 Here Λ QCD is QCD scale parameter.
The rapidity density can be evaluated in CGC Effective Field Theory (EFT) by using the following expression: 83 where σ AA is the inelastic cross-section of nucleus-nucleus interaction. Further, the differential cross-section of gluon production in a A − A collision is written as: 83 In the saturation regime, where S A is the nuclear overlap area.
Since the rapidity (y) and Bjorken variable x are related by ln 1 x = y, the x-dependence of the gluon structure function translates into the following rapidity dependence of the saturation Integrating over transverse momentum in (39) and (40), the rapidity distribution comes as follows: 83 where const is energy-independent, Q 2 s ≡ Q 2 s (s; y = 0) and Q s,min(max) are defined as the smaller (larger) value of Eq. (43); at

Model based on the percolation of strings
The physical situation described by the Glasma, the system of purely longitudinal fields in the region between the parting hadrons, is analogous to that underlying in the string percolation model (SPM). In the SPM one considers Schwinger strings, which can fuse and percolate, [85][86][87][88] as the fundamental degrees of freedom. The effective number of strings, including percolation effects, is directly related to the produced particle's rapidity density. The SPM 89 for the distribution of rapidity extended objects created in heavy-ion collisions combines the generation of lower center-of-mass rapidity objects from higher rapidity ones with asymptotic saturation in the form of the well-known logistic equation for population dynamics: where ρ ≡ ρ(∆, Y ) is the particle density, Y is the beam rapidity, and ∆ ≡ |y| − Y . The variable -∆ plays the role of evolution time, parameter δ controls the low density ρ evolution and parameter A responsible for saturation. The Y-dependent limiting value of ρ, determined by the saturation condition ∂ρ ∂(−∆) =0, is given by ρ Y = 1 A . Whereas the separation between the 13/50 region ∆ > ∆ 0 (i.e, low density and positive curvature) and the region ∆ < ∆ 0 (i.e, high density and negative curvature) is defined by Integrating Eq. (46) between ρ 0 and some ρ(∆), one can obtain: 89 The particle density in SPM is proportional to N A (i.e, ρ ∝ N A ), therefore we can write To be more specific regarding the quantities δ, ρ Y and ∆ 0 , δ does not strongly depend on Y , the parameter ρ Y is the normalized particle density at mid-rapidity and is related to gluon distribution at small x i.e, ρ ∝ e λY . As ρ Y increases with rapidity Y, energy conservation arguments gives ∆ 0 = −αY with 0 < α < 1. Now Eq. (48) can be rewritten as The particle density ρ is a function of both ∆ and the beam rapidity Y whereas according to the hypothesis of limiting fragmentation, for ∆ larger than some Y -dependent threshold, the density ρ remains a function of ∆ only. In the Glauber model, energy-momentum conservation constrains the combinatorial factors at low energy but in the framework of SPM energymomentum conservation is accounted by reducing the effective number of sea strings via where which yields In the string percolation model (SPM), the charged particle multiplicity at mid-rapidity in proton-proton interactions is given by: 133 such that for symmetric nucleus-nucleus case, one has: 133 In the SPM, the power law dependence of the multiplicity on the collision energy is the same in pp and AA collisions.
The pseudo-rapidity dependence of AA collisions is given by: 90 where J is the Jacobian J = coshη

Charged Hadron Multiplicity Distributions
Here, we discuss various observed features of the measured charged hadron multiplicity in hadron-hadron, hadron-nucleus and nucleus-nucleus interactions at different energies and their comparison in the view of above mentioned models. Feynman was first to point out that the multiplicity spectrum observed in protonproton collisions at asymptotically large energies becomes independent of the center-of-mass energy( √ s) as √ s → ∞. [105][106][107] This assumption is known as F eynman scaling. He concluded this fact primarily based on the phenomenological arguments about the exchange of quantum numbers between the colliding particles. This energy independent behaviour of the height of rapidity distribution around mid-rapidity naturally implies that the total multiplicity after integration over rapidity involves ln √ s behaviour since y max = ln √ s m N , where m N is the nucleon-mass. However up to √ s = 1800 GeV the experimental data does not indicate that the height of rapidity distribution around mid-rapidity (i.e, ( dn dη ) η=0 ) gets saturated. Thus, Feynman scaling is violated by continuous increase in the central rapidity density. Later on, it was realized that the gluons arising from gluon-bremsstrahlung processes give QCD-radiative corrections. 108 Recently, it was noticed by PHOBOS Collaboration from the p − p and/or p −p data that central plateau height i.e., ( dn dη ) η=0 grows as ln 2 √ s which in turn will give ln 3 √ s type behaviour in the total multiplicity. 45 clearly indicates that this scaling type of behaviour is not supported by experiments. Based on the above findings, we have included ln 3 √ s behaviour as in Eq. (15) and have shown its applicability in describing the p − p and p −p data quite successfully for the available entire energy range till date.
In Fig. 1, the inelastic (filled symbols) and non-single diffractive (NSD) data (open symbols) of charged hadron multiplicity in full phase space for p − p collisions at various centerof-mass energies from different experiments e.g., ISR, UA5 and E735 is shown. Inelastic data at very low energies (filled symbols) are used because NSD data are not available for these energies and also the trend shows that the difference between inelastic and NSD data is very small at lower energies. Further, this data set is fitted with three different functional forms in order to make a simultaneous comparison with our parameterization as given in Eq. (15).
The short-dashed line has the functional form as : a 2 + b 2 s 1/4 which is actually inspired by the Fermi-Landau model. It provides a resonable fit to the data at higher √ s N N with a 2 = 5.774 and b 2 = 0.948. 111 However, since a 2 summarizes the leading particle effect also, its value should not be much larger than two. The dotted line represents the functional form as : In Table 1 and Table 2, the wounded quark model 41 results are organized in tabular form for ten centrality classes and are compared accordingly with the RHIC multiplicity data. 121

Charged hadron pseudo-rapidity distributions
The pseudo-rapidity distribution of charged hadrons is another quantity in the studies of the particle production mechanism in high energy hadron-hadron and nucleus-nucleus collisions. Pseudo-rapidity density is a well defined quantity which is sensitive to the initial con-  GeV, which can be well explained as: 117 However, the measurements in Pb-Pb collisions at 2.76 TeV observed at LHC exceed the above logarithmic dependence. This logarithmic behaviour describe Au-Au and Pb-Pb data quite successfully with a little disagreement for Cu-Cu data. A better description of the data for energy dependence in A-A collision is given by the following relation: 117, 118 which is quite successful in describing the data at RHIC and LHC, but fails to explain the lower energy data as it overestimates in lower energy region below 17.3 GeV. 117  the experimental pseudo-rapidity density data at mid-rapidity normalized by the participating quarks calculated in wounded quark model. Open points are wounded quark model result which is well described by the following functional form: A little deviation from this functional from can be clearly seen for Cu − Cu collisions at 19.6 and 62.4 GeV , which arises due to the difference in number of quark collisions from Au − Au collisions at the same energy and signifies the role of hard processes involved.

KNO Scaling
The energy dependence of the multiplicity distribution observed in high energy collision using a variety of colliding systems is one of the primary issue in search for any systematics/scaling in multiparticle production. Based on the assumption of the Feynman scaling for the inclusive particle production cross section at asymptotic energies Koba, Nielsen and Olesen (KNO) in 1972 proposed scaling property for the multiplicity distribution at asymptotic energies. 138 According to which the normalized multiplicity distributions of charged particles should become independent at asymptotic energies. According to KNO scaling, the probability P(n) to produce n charged particles in the final states is related to a scaling function ψ(z) where ψ(z) is universal function independent of energy and the variable z = n n stands for normalized multiplicity. σ n and σ tot corresponds to n particle production and the total cross section respectively. Thus, the rescaling of P n measured at different energies via stretching or shrinking the axes by the average multiplicity n leads the rescaled curves coinciding with one another. However, the experimental data show only an approximate KNO scaling behaviour for multiplicity istributions in the different energy regions and hence Buras et. al 139 proposed Modified KNO scaling in which variable z ′ = n−α n−α . Here parameter α depends solely on the reaction. In order to describe the properties of multiplicity distributions of final state particle at different energies, it is more convenient to study their moments. KNO scaling also implies that the multiplicity moments C q = n n q P n ( n nP n ) q = n q n q 29/50 are energy independent since where q is the order of the moment. Later on, a comparative study of charged particle multiplicity distributions 143 arising from non single diffractive inelastic hadronic collisions at √ s=30 GeV to 1800 GeV including the data in UA5(SPS) and E735 (Tevatron) experiments again confirmed the deviations from so called KNO scaling behaviour in the form of a shoulder like structure clearly appeared in the collider data. This shoulder like structure observed in the collider data arises due to the superposition of the distribution of the particles from some other process different from KNO producing process incoherently superimposed on the top of the KNO producing process. 141

30/50
A strong linear increase of the C q moments with energy and strong KNO scaling violation at √ s=7 TeV in form of the observed change of slope in P n confirm these earlier measurements.
Recently, A. Capella et al. 64 have computed pp multiplicity distributions at LHC in the framework of a multiple-scattering model(DPM) as shown in Fig.13 and Fig. 14

Intermittency
A remarkably intense experimental and theoretical activity has been performed in search of scale invariance and fractality in soft hadron production processes. In this search, a lot of efforts have been performed by investigating all types of reactions ranging from e + e − annihilation to nucleus-nucleus collisions at different energies. In this investigation, Self-similarity (or intermittency) is a power-law behaviour of a function of the form f (x) ∼ x a and thus it also gives the scaling property f (λx) ∼ λ a f (x). Since the scaling is the underlying property of critical phenomena, it is natural to think of a phase transition wherever scaling (or selfsimilarity) is found. Self-similarity is also connected with a fractal pattern of the structure. Similarly the scaled factorial moments are defined as: where n q = ∞ n=1 n q P n . Intermittency is defined as the scale invariance of factorial moments with respect to changes in the size of phase-space cells or bins (say δy), for small enough δy: or This give the wide, non-statistical fluctuations of the unaveraged rapidity distribution at all scales. Such a kind of behaviour is a manifestation of the scale invariance of the physical process involved. In Fig. 15 and 16, the intermittent behaviour of charged hadrons produced in the nucleus-nucleus interactions are shown in rapidity-and azimuthal-space. 153 The slope (φ q )in a plot of lnF q vs -lnδy at a given positive integer q is called intermittency index or intermittency slope which provide us an opportunity to characterize the apparently irregular PS) predictions were also found to be consistent with the data. The CELLO collaboration analysed the 3D intermittency signal in e + e − annihilation and found good agreement with the LUND model. The hadron-proton (π + − p and K + − p) collisions were analyzed by the NA22 collaboration. The data for p−p collisions at 630 GeV were analyzed by the UA1 Collaboration and an indication for the increase of the intermittency signal for the low-multiplicity sample was found. A multi chain version of DPM including mini-jet production was compared to NA22 and UA1 data but the slopes were found to be too small. [164][165][166] The NA22 group performed the analysis for hadron-nucleus interactions and they found a weaker intermittency signal for large targets. The hadron-nucleus and nucleus-nucleus interactions were also studied by KLM collaboration in 1d and 2d. Decrease in the intermittency signal was observed for larger projectile nuclei and this decrease is smaller than expected from the increase of the mean multiplicity in the collision. The EMU-01 collaboration performed the intermittency analysis for different colliding systems and found a similar dependence. Further, the results for the SFMs in 3D spectra for e + e − , hadron-proton, hadron-nucleus and nucleus-nucleus collisions have shown that the dependence of the SFM on the resolution is stronger than a power law.
The failure of Monte Carlo calculations to replicate the observed strength of the intermittency signal is quite evident in the nucleus-nucleus interactions. 167 A deeper investigations on 1D, 2D and 3D phase space in different kinematic variables (η, φ, p T ) was performed by OPAL collaboration for high statistics study at LEP upto fifth order and down to very small bin sizes. 160, 161 Later on, it was shown that the Negative Binomial distribution (NBD) faces 33/50 difficulties to describe the high statistics genuine multiparticle correlations and no one of the conventional multiplicity distributions, including NBD can describe the high statistics data on intermittency of OPAL measurements. 162,163 Here, we discuss some important systematics and regularities observed in recent years in this intermittency study. The phase-space density of multiparticle production is anisotropic, and the upward bending in the SFM plot is a direct consequence of this anisotropy. 168,169 Thus, it is suggested that in the (2d) SFM analysis, the phase-space should be partitioned asymmetrically 170 for taking into the account the anisotropy of the phase-space by introducing a 'roughness' parameter called as the Hurst exponent (H). The Hurst exponent is calculated by fitting the (1d) SFM with the Och's formula. [171][172][173] For H < 1.0, the φ-direction is partitioned into finer intervals than the η-direction, for H > 1.0, the vice-versa is true while, H = 1 means that the phase-space is divided similarly in both the directions. Ochs and Wosiek 174 observed that 1D factorial moments even if they do not strictly obey the intermittency power law over the full rapidity range, still obey the generalized power law where the function g depends on the energy and the bin width δy. Eliminating g yields the linear relation This relation is found to be well satisfied by the experimental data for various reactions. The dependence of φ q on q can be examined by establishing a connection between intermittency and multifractality. 34 The generalized Renyi dimensions (D q ) are related to φ q as where D being the topological dimension of the supported space and the anamolous dimensions d q are defined as If we plot the ratio of anamolous dimensions as a function of the order q for various reactions, the q-dependence is claimed to be indicative of the mechanism causing intermittent behaviour and all points fall on a universal curve 34/50 parameterized as Above parameterization has been derived by Brax and Peschanski 175 for a self-similar partonbranching process. µ (0 ≤ µ ≤ 2) is called Levy stable index. 176 For µ = 2 Eq. (69) reduces to Eq. (68). The multifractal behaviour characterized by Eq. (68) and (69) reduces to a monofractal behaviour. 177-179 µ = 0 will imply an order-independent anamolous dimension which will happen if intermittency were due to a second-order phase transition. 165 Consequently, monofractal behaviour may be a signal for a quark-gluon plasma phase transition. A noticeable experimental fact that the factorial moments of different order follow simple relation (65), means that correlation function of different orders are not completely independent but are interconnected in some way. Intermittency can be further studied in the framework of Ginzberg-Landau (GL) model as used to describe the confinement of magnetic fields into fluxoids in type II superconductor. According to which, the ratio d q /d 2 should respect the with ν being a universal quantity independent of the underlying dimension. The value of ν was first derived analytically by R. C. Hwa et al. 180 in Ginzburg-Landau model and was confirmed later in a laser experiment. 181 If there is any phase transition, it may not necessarily always be a thermal one as the new phase formed is not essentially characterized by the thermodynamical parameters. The possibility of the simultaneous existence of two nonthermal phases (in analogy to different phases of the spin-glass) can be investigated by the intermittency parameter 176 In case of the existence of such different phases in a self-similar cascade mechanism, λ q should have a minimum at some value q = q c . The regions q < q c resembles liquid phase with a large nuber of small fluctuations and q > q c resembles dust phase with a small number of large fluctuations. Another way to measure the intermittency strength in terms of the intermittency index in the framework of a random cascading model like the α-model. 182 According to which the strength parameter α q is related to the Renyi dimensions as 183 A thermodynamic interpretation of multifractality can also be given in terms of a constant specific heat C, provided the transition from monofractal to multifractal is governed by a

35/50
Bernoulli type of fluctuation. Bershadskii proposed a phenomenological relation among the D q and C as 183 C =0 in the monofractal phase which becomes non-zero finite in the multifractal phase.
The experimental data on intermittency do not satisfy the criterion for the phase transition and thus one can say the observed intermittency patterns are not suitable probe for QGP formation. 167

Negative binomial multiplicity distribution
Another interesting and well known feature observed in multiplicity distribution of produced charged particles in high energy collisions is the occurrence of a negative binomial distribution (NBD) which holds good for a variety of collision processes. It has been observed experimentally for hadron-hadron, hadron-nucleus and nucleus-nucleus collisions over a wide range of energy in finite rapidity as well as in the full phase space. It has further been found that the NBD is valid not only for hadronic, but for semileptonic and leptonic processes as well which emphasize that it is a general property of multiparticle production process regardless of type of colliding particles. To be more specific, the NBD qualitatively describes well the multiplicity distribution almost in all inelastic high energy collision processes except for the data particularly at the highest available collider energies where some deviations seem to be apparent. The negative binomial probability distribution for obtaining n charged particles in the final state is given as follows: wheren and k are two parameters varying with energy and determined by the experimental data. The parameter n has the interpretation of average multiplicity and k is the parameter influencing the width of the distribution . These two parameters are related to the dispersion D = n2 −n 2 as follows: It is observed that average multiplicity (n) increases with energy ( NBD is quite interesting quantity as for increasing k(k → ∞), the probability distribution (P N BD ) gets narrower tending to be the Poisson distribution while for the decreasing values of parameter k, the probability distribution becomes broader and broader than the Poisson distribution. For k = 1, the P N BD is given by the Bose-Einstein (or geometric distribution).
Under the limit of large multiplicity (n ch → large), the P N BD goes over to a gamma distribution. The negative binomial distribution with fixed parameter k exhibits F scaling and asymptotic KNO scaling. 3 Interpretation of the empirical relation (73) in terms of underlying production mechanism 37/50 common to hadronic, leptonic and semi leptonic processes is still a challenging problem till date. Some attempts have already been made to derive the NBD from general principles using a stochastic model, 185 a string model, 186 a cluster model, 187 a stationary branching process, 188 and a two-step model of binomial cluster production and decay. 189 However, it is still difficult to understand why the same distribution fits such diverse reactions. Moreover, we are still lacking a microscopic model explaining the behaviour of parameter k of the NBD distribution which according to fits decreases with the energy. Recently the NBD form has theoretically been derived in a simplified description of the QCD parton shower with or without hadronization.
In particular, the NBD qualitatively describes the multiplicity distributions in almost all inelastic high energy collisions except for the data recently in collider experiment at much higher energies showing the deviation from the NBD distribution. For the charged particle multiplicity distribution at √ s=900 GeV SPS energy, a single NBD function could suitably describe the data only for small pseudo-rapidity intervals at the mid-rapidity region whereas for large intervals, shoulder-like structure appeared in the multiplicity distribution. Appearance of sub-structures in multiplicity distributions at higher energies and in larger pseudo-rapidity intervals has been attributed to weighted superposition or convolution of more than one functions due to the contribution of more than one source or process of particle production. The two-component model of A. Giovannini and R. Ugoccioni 141,190 is quite useful in explaining the multiplicity distribution data at the cost of increased number of adjustable parameters, in which they used weighted superposition of two NBDs representing two class of events, one arising from semi-hard events with minijets or jets and another arising from soft-events without minijets or jets.
where α sof t signifies the contribution of soft events and is a function of √ s, the other parameters are functions of both the √ s and the η c , have their usual meanings as in Eq. (74) indicating respective components. For small pseudo-rapidity intervals, a single NBD functions explains the distribution data at √ s= 0.9 and 2.36 TeV reasonably well but at √ s= 7 TeV, a single NBD function appears inadequate while weighted superposition of two NBDs explain the data satisfactorily as shown in Fig. 17. Whereas, for large pseudo-rapidity intervals (η c < 2.4), the weighted superposition of two NBDs agree better than a single NBD function with the distribution data for all available LHC energies clearly indicating the role of soft and semi-hard processes in TeV range.

Limiting Fragmentation or Longitudinal Scaling observed for charged hadrons
In order to understand the collision dynamics completely, the study of particle production away from midrapidity region provides some more insight. The properties of charged hadron production in the fragmentation region can be studied well by viewing the data at different energies in the rest frame of one of the colliding nuclei as the particles near beam and target rapidity are supposed to be governed by the hypothesis of limiting fragmentation. 106 The hypothesis of limiting fragmentation for inclusive particle distributions in high energy hadron-hadron collisions was first proposed by J. Benecke et al. 106 in 1969. According to this hypothesis, the produced particles in the rest frame of projectile or target will approach a limiting distribution e.g., the pseudo-rapidity density ( dn dη ) as function of η ′ (η − Y beam ) approach a limiting curve in the fragmentation region. This picture is based on the geometrical picture of scattering as considered by Yang and co-workers. [191][192][193] According to their assumption, in the laboratory frame the projectile nuclei undergoes a Lorentz contraction in the form of a thin disk during the collision with the target nuclei which gets further and further compressed with the increasing energy. However, the momentum and quantum number transfer process between the projectile and target does not appreciably change during this compression. This behaviour leads to limiting distribution of the produced particles in the fragmentation region independent of the center-of-mass energy. Central to this limiting fragmentation hypothesis was the assumption of the constancy of the total hadronic cross sections at asymptotically large center-of-mass energies. In other words, the probability of the interaction does not change rapidly with the further increase of energy of the incident colliding nuclei. It is expected that the excitation and breakup of a hadron would become independent of the center-of-mass energy and the produced particles in the fragmentation region would approach a limiting distribution. Later on, it was realized that the total hadronic cross sections are not constant, they grow at a slow rate with the increasing center-of-mass energy ( √ s). 194,195 This slow increase in the cross sections is still continued upto the LHC energy. In spite of this fact limiting fragmentation is found to remain valid over a wide range of energy. This asymptotic property has been observed experimentally in a variety of collision processes such as hadron-hadron, 29,196 hadron-nucleus, 29,197 and nucleus-nucleus 120,124,125,153,[198][199][200][201][202][203][204][205][206] for produced charged hadrons in fragmentation region at different energies. It can be clearly observed that this scaling fea- ture covers a more extended range of η than expected according to the original hypothesis of limiting fragmentation due to which the term extended longitudinal scaling is often used to describe this phenomenon. For most central collision, the distributions in fragmentation region are indeed observed to be independent of the collision energy over a substantial shifted pseudorapidity (i.e., η ′ = η − Y beam ) range. In Fig. 18 between the calculations and the experimental data in the mid-rapidity region is seen which arises mainly due to the violation of k ⊥ formalism in that regime as k ⊥ formalism becomes less reliable the further one is from the dilute-dense kinematics of the fragmentation regions. [208][209][210][211] In Fig. 19, the pseudo-rapidity distributions at different energies with different colliding nuclei lie on a common curve in fragmentation region over a broad range of η ′ indicating that limiting fragmentation holds good for A-A collision over a wide range of energy and is independent of the size of the colliding nuclei. This gives a clear indication that this universal curve is an important feature of the overall interaction and not simply a breakup effect. 201 In the CGC approach, the physical picture behind the limiting fragmentation is based on the black disk limit. According to which, in the rest frame of the target nucleus, the projectile nucleus is highly contracted and due to its large number of gluons looks black to the partons in the target nucleus which interact with the projectile nucleus with unit probability.. 213 The of the two colliding hadrons. This mechanism provides a natural explanation of the observed rapidity spectra, in particular their linear increase with increasing rapidity distance from the maximal rapidity and the short plateau in the central rapidity region.
In order to understand the nature of hadronic interactions that led to limiting fragmentation, a new mechanism is required which can give a better understanding of the subject.
Although, the limiting fragmentation behavior in small x region is quite successfully described by Color Glass Condensate (CGC) approach. However, the procedure employed in this approach needs further improvements. One of them is the impact parameter dependence of the unintegrated gluon distribution functions. Whereas, in large x region, the phenomenological extrapolation is used which also need to be better constrained and consistent

Summary and Conclusions
In this review, we have attempted to draw certain systematics and scaling laws for hadronhadron, hadron-nucleus and nucleus-nucleus interactions. Any distinct deviation from these relations observed in ultra-relativistic nuclear collisions will be the indicator of a new and exotic phenomenon occuring there. However, there are other effect which causes deviations in the scaling features such as mixing of contributions from soft and hard processes, etc. Present experimental results show that the scaling features of charged hadrons previously observed at lower energies do not hold good at higher energies. Thus, we are still lacking a scaling law which is universal to all types of reactions and can give some basic understanding of the mechanism involved in the particle production. The observed deviations from the well established scaling laws are of great interest as these clearly hint towards the increasing role of hard processes contributing in the production of final state particles at higher energies. These hard processes of the primordial stage of the collision contribute significantly to the hadronization process at the late stage of the evolution. The applicability of the two-component approach in nucleonnucleon interaction picture and in wounded quark picture to evaluate mid-rapidity density of charged hadrons is quite promising as the weighted superposition of hard and soft processes comes into picture and successfully describe the available experimental data over a large range of energy. Our study of the charged hadron production using phenomenological models along with the experimental data essentially highlights some kinds of scaling relations for these high energy collision processes. N part scaling of pseudo-rapidity density at mid-rapidity holds only for limited range of energy, while for the entire range of energy (including LHC) it does not seem to hold good which clearly hints towards the contributions of some other processes. On the other hand, the wounded quark model provides a more realistic explanation of charged 42/50 hadron production in terms of the basic parton structure involved in these processes. However, the charged hadron pseudo-rapidity density at mid-rapidity in wounded quark model also does not show an exact N q scaling type behaviour and the deviations are clearly observed which arises due to difference in the number of quark collisions for different colliding nuclei as involved In conclusion, the multiparticle production in the ultra-relativistic hadron-hadron, hadronnucleus, and nucleus-nucleus collisions is still a burning topic in high energy physics because it throws light on the production mechanism. Since we still lack explanation of these processes in QCD, we hope that our investigations regarding universal systematics and scaling relation on charged hadron production will certainly lead us to a correct phenomenological model for the production mechanism. In addition to this, recent studies regarding fluctuations and correlations is expected to provide deeper insight into different stages of the system evolution with the participant nucleon eccentricity (ǫ part ) . 215,216 Further, elliptic flow scaled with ǫ part for the same number of participating nucleons in Cu − Cu and Au − Au system show a scaling pattern upto 3.5 GeV /c in p T and across ±5 units of pseudo-rapidity. 215,216 Thus, it will definitely be interesting to investigate these fluctuations and correlations in future to draw a more consistent and unified picture of particle production and their evolution in heavy-ion collisions.