Multiplicity and pseudorapidity density distributions of charged particles produced in pp, pA and AA collisions at RHIC \& LHC energies
Sumit Basu, Sanchari Thakur, Tapan K. Nayak, Claude A. Pruneau
MMultiplicity and pseudo-rapidity densitydistributions of charged particles produced in pp,pA and AA collisions at RHIC & LHC energies
Sumit Basu , , Sanchari Thakur , Tapan K. Nayak , , andClaude A. Pruneau Department of Physics & Astronomy, Wayne State University, MI 48201 USA Lund University, Department of Physics, Division of Particle Physics, Box 118,SE-221 00, Lund, Sweden Variable Energy Cyclotron Centre, HBNI, Kolkata-700064, India National Institute of Science Education and Research, HBNI, Jatni 752050, India CERN, CH 1211, Geneva 23, SwitzerlandE-mail: [email protected]
E-mail: [email protected]
E-mail: [email protected]
E-mail: [email protected]
Abstract.
Multiplicity and pseudorapidity ( η ) density ( dN ch /dη ) distributions of chargedhadrons provide key information towards understanding the particle productionmechanisms and initial conditions of high-energy heavy-ion collisions. However,detector constraints limit the η -range across which charged particle measurementscan be carried out. Extrapolating the measured distributions to large η -range byparameterizing measured distributions and by using calculations from event generators,we characterize the production of charged particles over the full kinematic range. In thepresent study, we use three different ans¨atze to obtain quantitative descriptions of theshape of pseudorapidity distributions of charged hadrons produced in pp, p-A, and A-Acollisions for beam energies ( √ s NN ) ranging from a few GeV to a few TeV correspondingto RHIC and LHC energies. We study the limiting fragmentation behavior in thesecollisions and report evidence for participant-scaling violations in high-energy collisionsat the TeV scale. We additionally examine measured pseudorapidity distributionsto constrain models describing initial conditions of particle production. We predictthe centrality dependence of charged particle multiplicity distributions at FAIR andNICA energies and give an estimation of charged particle multiplicity at η = 0 for theproposed HE-LHC and FCC energies. a r X i v : . [ nu c l - e x ] N ov ultiplicity and pseudorapidity density distributions
1. Introduction
Multiplicity and pseudorapidity ( η ) density distributions of charged particles alongwith the transverse momentum ( p T ) spectra constitute some of the basic observablesfor understanding the particle production mechanisms in high-energy elementaryparticle and heavy-ion collisions [1–4]. The dependence of these distributions on thecolliding particle species, collision energy, and collision centrality have been extensivelydiscussed in the literature [5–10]. In proton-proton (pp) collisions, these distributionsprovide precise calibration of particle production models such as PYTHIA [11] andHERWIG [12], which are used to make predictions of various searches including thoseof physics beyond the standard model. These measurements play an important rolein the study of heavy-ion collisions at ultra-relativistic energies in which short-livedsystems consisting of nuclear matter at extreme conditions of temperature and energydensity are created. There is evidence that this matter undergoes a phase transitionfrom a confined state to a de-confined state of quark-gluon plasma (QGP) [13, 14]. Thecomparison of the charged particle distributions in pp, p–A, and A–A collisions areessential to characterize the formation of QGP and understand the particle productionmechanisms.The measured charged-particle multiplicity and p T distributions are dominatedby final state interactions and the state of matter at freeze-out. Nonetheless, thesedistributions are also sensitive to the initial stages of the collision. At small Bjorken- x (expressed as x = p T √ s .e − y ∼ p T √ s .e − η , y being the rapidity), the gluon density of theparton distribution functions (PDF) of the proton grows and is expected to reach asaturation domain [15–17]. So the particle productions at large collision energies andforward rapidities are characterized by a large number of gluons [18–23]. With presenttheoretical understandings and available deep inelastic scattering (DIS) experimentsat HERA [24, 25], it becomes possible to study the expected growth and saturationof gluon density at high collision energies. The measurements of produced charged-particle multiplicity, p T , and η distributions are typically restricted to the mid-rapidityregion. Such restrictions arise in part because of favorable kinematic conditions atmid-rapidity and largely because of experimental limitations at forward rapidities.Understanding the particle production dynamics, including effects of nuclear stopping,color transparency, jet quenching, and long range correlations, requires the measurementof particle production in full pseudorapidity ranges. However, energy-momentumconservation dictates that particle production must vanish at or beyond the beamrapidity. It is thus of interest to consider pseudorapidity ansatz that assumes vanishingdensity at such large rapidities. These may then be compared to production modelsand may, in principle, be used to estimate the total charged-particle production in ppand A–A collisions. The precision achievable with such extrapolations is obviouslylimited by the quality of the ansatz but it can be tested with existing pp and A–Acollision models. However, no specific or widely accepted pseudorapidity rapidity modelis currently available in the literature to carry out such extrapolation. In this work, ultiplicity and pseudorapidity density distributions N part ) [29, 30]. But at higher energies of the CERN Large Hadron Collider (LHC),the N part scaling has been observed to be broken [31, 32]. One of the main reasonsfor this scale breaking is the enhancement in gluon productions at high energies (lowBjorken- x ) [16, 17, 33].The paper is organised as follows. In sec. 2, we examine the measured charged-particle multiplicity, and pseudorapidity distributions, dN ch /dη , observed in pp, p–A,and A–A collisions across a wide range of beam energies and compare these withresults from selected event generators. In sec. 3, we parameterize these dN ch /dη distributions using three different ans¨atze to obtain a satisfactory model one canintegrate over the full pseudorapidity range spanned by particle production. Suchan extrapolation requires that we examine, in sec. 4, the measured distributions inthe vicinity of the beam rapidity and study the applicability of the notion of limitingfragmentation. In sec. 5, we use the favored ansatz to estimate the total number ofcharged particles produced per N part as a function of collision energy and centrality.We inspect whether the charged particle production scales with N part irrespective of thecollision energy. Using the parameterization of the pseudorapidity density distributions,we give predictions for these distributions as well as total charged particle multiplicitiesfor lower collision energies of the future experiments at FAIR (Facility for Anti-protonand Ion Research) at GSI, Germany and NICA (Nuclotron-based Ion Collider fAcility)at JINR, Russia. Additionally, we extend the charged particle multiplicity density at η = 0 for the proposed HE-LHC (High-Energy LHC) and FCC (Future Circular Collider)at CERN. Finally, in sec. 6, we explore whether selected initial condition scenarios canbe meaningfully constrained by measured particle multiplicity distributions. The paperis summarized in sec. 7.
2. Charged-particle multiplicity distributions
In this section, we present the charged particle multiplicity density at mid-rapidityand pseudorapidity distributions from available experimental data for pp, p–¯p, d–Au,p–Pb, Au–Au, and Pb–Pb collisions. These data are compared to calculations fromevent generators. For pp and p–¯p collisions, the multiplicities are calculated with ultiplicity and pseudorapidity density distributions p T distributions,and anisotropic flow [41–45]. THERMINATOR (THERMal heavy IoN generator) isa statistical hadronization model commonly used to estimate the relative abundancesof particles species produced in relativistic heavy-ion collisions. It enables arbitraryimplementations of the shape of the freeze-out hyper surface and the expansion velocityfield. The multiplicities were computed including HBT effects and 3+1 dimensionalprofiles [46]. The charged-particle multiplicity density at mid-rapidity dN ch /dη | η =0 has been reportedfor different colliding systems, collision centrality and collision energies. The averagenumber of participants ( (cid:104) N part (cid:105) ) characterizes the collision centrality and collidingsystem. In Fig. 1, we present a compilation of the measurements of scaled charged-particle multiplicity density at mid-rapidity, (cid:104) N part (cid:105) dN ch /dη | η =0 , as a function of collisionenergy in pp [47–49], p–¯p [50, 51], Au–Au [5, 52–57], Pb–Pb [58–60], d–Au [61], and p–Pb [62–64] collisions observed at Fermilab, RHIC, and LHC energies. Results frompp and p–¯p collisions are for non-single diffractive (NSD) as well as inelastic (INEL)collisions, whereas those from Au–Au and Pb–Pb collisions correspond to most centralcollisions.The multiplicity densities measured in pp (p–¯p) and A–A collisions exhibit ratherdifferent dependence as a function of collision energy. These dependencies can becharacterised with power-law fits performed separately for A–A , NSD pp, and INEL pp(p–¯p) collisions. We find that the √ s dependence of the multiplicity density of pp (p–¯p)collisions are well matched by power laws of the form ( s NN ) α with exponent α = 0 . α = 0 .
11 for INEL and NSD collisions, respectively. In contrast, the multiplicity ultiplicity and pseudorapidity density distributions (GeV) NN S 10 = η ) η / d c h ( d N 〉 pa r t N 〈 ALICE 5.02 TeVALICECMSATLASSTARPHENIXPHOBOSBRAHMSNA50
Central collision (AA)
EPOSUrQMDTHERMINATOR /pp pp (NSD)
ALICECDFUA5 s ∝ s ∝ /pppp (INEL) ALICENA49PHOBOSPYTHIAEPOSd+AUp+Pb AMPT SM s ∝ (Pb+Pb) (PHOBOS)(ALICE) Figure 1.
Compilation of measurements of the beam-energy dependence of charged-particle multiplicity density at mid-rapidity, scaled by the average number ofparticipating nucleon pair ( (cid:104) N part (cid:105) / densities observed in A–A collisions exhibit a steeper dependence on the beam energyand are best described with a power law exponent α = 0 . √ s NN ≥
100 GeV. In the case of A–A systems, onefinds that AMPT SM and UrQMD predictions are in good agreement with data overthe entire √ s NN range considered in this work. We additionally find that EPOSpredictions are also in reasonable agreement with data from both p–p and A–A systemsover a wide range of beam energies. However, the single THERMINATOR predictionconsidered at √ s NN = 200 GeV is found to considerably underestimate the measuredcharged-particle density. Overall, PYTHIA , EPOS, AMPT, and UrQMD are foundto reproduce reasonably well the observed √ s NN power law behavior even though theyare based on rather different interaction and transport models. The √ s NN evolutionof the (cid:104) N part (cid:105) dN ch /dη | η =0 , an integrated observable, is not a strong discriminant of ultiplicity and pseudorapidity density distributions - - h / d c h d N ALICE p+p - - h h / d c h d N PHOBOS 2760 GeVALICE 19.6 GeV62.4 GeV130 GeV200 GeV Pb+Pb 0-5% EPOSUrQMDAMPT(0-5%)Au+Au 0-6%THERMINATOR(0-5%)
Figure 2.
Comparison of selected experimental dN ch /dη distributions of measuredin pp and p–¯p collisions (upper panel) and Au–Au and Pb–Pb collisions (lowerpanel) with calculations performed with the PYTHIA, AMPT, UrQMD, EPOS andTHERMINATOR models. these models and their underlying particle production mechanisms. Indeed, althoughhistorically important, measurements of scaled charge particle density at central rapidityprovide only a rather limited amount of information about the specific particle creationand transport mechanisms involved in pp and A–A collisions [10]. Additional insightinto these mechanisms, however, may be gained from charged-particle pseudorapiditydistributions. Figure 2 presents a compilation of σ dσ/dη distributions measured in pp,Au–Au, and Pb–Pb collisions at the SPS, RHIC and LHC. The upper panel of Fig. 2 displays distributions measured in pp collisions at energiesranging from 0.9 to 13 TeV [47–50], and p–¯p collisions at 0.2 TeV [51] (open symbols), ultiplicity and pseudorapidity density distributions . ≤ √ s NN ≤
200 GeV [5,54–57,61],and top 5% central Pb–Pb collisions at 2.76 TeV [58, 60].First focusing our attention to the upper panel of Fig. 2, we note that only theUA5 data cover a wide enough pseudorapidity range capable of revealing the full shapeof the η distribution measured in pp collisions while the measurements reported bythe ALICE collaboration are limited to the central rapidity region in pp collisions.One nonetheless observes that the particle density produced in pp (p–¯p) collisionsrises monotonically, as expected, with beam energy. One also notes that the measuredpseudorapidity distributions all feature a dip, centered at η = 0, whose depth and widthincreases with rising √ s . The presence of this dip is in part associated with partialtransparency and limited stopping power of the colliding protons (or anti-protons) [73].The dip may also result, in part, from the measurement being reported as a functionof pseudorapidity. A boost invariant rapidity distribution would indeed yield a broaddip in pseudorapidity owing to the presence of mass term in the denominator of y → η transformation Jacobian.The pseudorapidity distributions measured in pp (p–¯p) collisions are comparedwith Monte Carlo calculations performed with PYTHIA 6.4 [11] (dash lines) andEPOS [44] (solid lines) event generators. One observes that both PYTHIA and EPOSapproximately reproduce the magnitude and η dependence of the distributions: bothmodels indeed capture the essential features of the measured distributions, includingthe presence of the widening and deepening dip, centered at η = 0, with increasing √ s . However, PYTHIA appears to be in better agreement with the data than EPOS at √ s = 0.2, 2.76, and 8 TeV. Observe, in particular, that EPOS does not reproduce thedip structure seen in p–¯p data at 0.2 TeV.Let us next examine the pseudorapidity distributions reported by the PHOBOSand ALICE collaborations shown in the lower panel of Fig 2. The PHOBOS datacover the range − . ≤ η ≤ . − . ≤ η ≤ .
5. One finds that the pseudorapidity distribution observedat √ s NN = 19.6 GeV features an approximate Gaussian shape peaked at η = 0,while distributions observed at higher beam-energy are much broader and feature dipstructures qualitatively similar to those observed in pp collisions. However, the dipstructures observed in A–A collisions are typically shallower and wider than thoseobserved in pp collisions.The experimental data are compared to calculations based on the UrQMD, AMPT,EPOS, and THERMINATOR models shown as solid lines in Fig. 2. The calculationswere performed within p T ranges reproducing the experimental conditions of thePHOBOS and ALICE experiments. Overall, we note that all four models qualitativelyreproduce the observed distributions. However, we also note that best agreement withthe measured pseudorapidity density distributions is obtained with the EPOS modelat beam energies √ s NN ≤
200 GeV, while at 2.76 TeV, none of the models reproducethe observed distributions quantitatively. Overall, all four models considered manage to ultiplicity and pseudorapidity density distributions η = 0,albeit with insufficient depth. EPOS’ predictions are ∼
5% low at RHIC energies andapproximately ∼
2% high at 2.76 TeV. AMPT and UrQMD, on the other hand, seemto reproduce the measured densities rather well at mid-rapidity, across the entire √ s NN range presented in Fig. 2, but fail to reproduce the overall shape at higher η values.THERMINATOR, on the other hand, is doing rather poorly at √ s = 200 GeV.
3. Parameterization of the pseudorapidity distributions
Although the PHOBOS and ALICE data presented in the lower panel of Fig. 2 coverlarge ranges of pseudorapidity, they do not capture the full range of particle productioninvolved at top RHIC energy and at the LHC. In fact, most reported measurementsof charged-particle pseudorapidity distributions are limited to somewhat narrow rangesof pseudorapidity and do not account for the full particle production. For instance,the measured distributions reported by the ALICE collaboration for Pb-Pb collisionscover a fairly wide range, | η | < .
5, but this range is quite narrow relative to the beamrapidity ( y beam ∼ . √ s NN dependence of the multiplicity densities observed in pp and A–A collisions result froman overall increase in the produced multiplicity per participant pair or does it resultsimply from a shift in the particle production towards central rapidity, due possibly tolarger stopping in A–A collisions?Very few experiments are equipped to cover the entire pseudorapidity range requiredto answer these questions unambiguously. In the absence of such measurements, weseek to parameterize the measured η distributions to obtain sensible extrapolations atforward/backward rapidities that may be used to estimate the total charged-particleproduction.In the transverse direction, extrapolation of measured particle densities, p T dNdp T , tozero and infinite p T is relatively straightforward because the cross-section must vanishat these limits [72]. Evidently, models used to integrate p T spectra are not constrainedoutside of the fiducial p T acceptance, but the fact that the cross-section vanishes at nulland infinite p T significantly constrains the shape of p T distributions. Uncertainties asto the exact shape of the p T distribution outside of the fiducial acceptance thus yieldsystematic uncertainties on the integral of the distributions.We seek to use the same concept towards extrapolating at forward and backwardrapidities beyond the fiducial pseudorapidity acceptance. The pseudorapidity densitymust obviously vanish at suitably large values of | η | but extrapolation beyond themeasurement acceptance does require one makes assumptions about the overall shape ofthe distributions. In this work, we explore three fitting functions to fit and extrapolate ultiplicity and pseudorapidity density distributions ) η / d c h ( d N pp data pA data AA data ) η / d c h ( d N pp data pA data AA data η − − ) η / d c h ( d N pp data η − − pA data η − AA data
A Ap p A p
Eq. (1)Eq. (2)Eq. (3) (Symmetric Trapezoidal)(Sum of two Gaussians)(Difference of Two Gaussians)
Figure 3.
Schematic representation of dN ch /dη distributions of pp, p–A, and A–A collisions using three different ans¨atze. The distributions shown as dotted linesare parameterized with three fit functions: symmetic trapezoidal (upper panels),sum of two Gaussian distributions (middle panels), and difference of two Gaussiandistributions (lower panels). The analysis of the shapes of the pseudorapidity distributions is based on thedistributions presented in Fig. 2. We first note that the pseudorapidity densitydistributions produced in symmetric collisions (e.g., pp and A–A ) are symmetric about η = 0, but distributions of the pA collisions feature a pronounced forward/backwardasymmetry, with an excess of particles being produced in the nucleus direction. Wefurther note that the shape of the pseudorapidity distributions may be characterizedby one broad peak with approximate Gaussian shape or two peaks of approximatelyGaussian shape separated by a trough. For simplicity, we thus consider three basicshapes defined according to: f T ( η ) = c (cid:112) − / ( α cosh η ) e ( | η |− β ) /a , (1) f G + G ( η ) = A e − ( η − µ σ + A e − ( η − µ σ , (2) ultiplicity and pseudorapidity density distributions f G − G ( η ) = A e − η σ − A e − ( η − µ )22 σ . (3) − − − − η d c h d N η − − − − D a t a / F i t (1) p+p
13 TeV8 TeV7 TeV − − − − η − − − − (2) pp+ − − − − η − − − − (3)Eq. Figure 4.
Beam energy dependence of charged particle pseudorapidity densitydistributions measured in minimum bias pp collisions by the ALICE collaboration [47–49] and in p–¯p collisions by the CDF collaboration [50]. Dashed lines show best fitsobtained with Eqs. (1-3).
Equation (1) is motivated by observed symmetric trapezoidal functions witha plateau around mid-rapidity [55]. It was used by the PHOBOS collaborationto model their measurements and extract the produced total particle multiplicity.Although Eq. (1) adequately reproduces some of the measured distributions, its built-insymmetry about η = 0 limits its applicability to symmetric collision systems exclusively.Equations (2) and (3) involve sum and differences of two Gaussian distributions,respectively. Equation (1) features four parameters ( c, α, β, a ), while Equation (2)involves six parameters ( A , µ , σ , A , µ , σ ). Equation (3) features five parameters( A , σ , A , µ, σ ) but reduces to four for symmetric collisions, which are characterizedby µ = 0.Figure 3 shows schematic diagrams of dN ch /dη distributions obtained for pp, p–A,and A–A collisions obtained with the ans¨atze embodied by Eqs. 1-3. The red and bluelines and associated shaded areas schematically represent the contributions from nucleonparticipants originating from either incoming nuclei. In the middle panel row, the bluedash line corresponds to the sum of both contributions. The shape of the distributionis here determined by the relative contributions of the incoming nuclei as well as the ultiplicity and pseudorapidity density distributions − − − η / d c h d N η − − − D a t a / F i t Au+Au 0 6 %
200 (GeV)13062.419.6
Au+Au 0 6%Eq.(1) − − − η − − − Pb+Pb 0 5 %Eq.(2) − − − η − − − Cu+Cu 0 6 %200 (GeV)62.422.4Eq.(3)
Figure 5.
Beam energy dependence of charged particle pseudorapidity densitydistributions measured in central Cu–Cu and Au–Au collisions by the PHOBOScollaboration [54–57] and in Pb–Pb collisions by the ALICE collaboration [58]. Dashedlines show best fits obtained with Eqs. (1-3) from left to right panels respectively. degree of nuclear stopping achieved in the collisions. In the bottom panel row, therelative contributions and stopping are modeled with a difference of two Gaussians asper Eq. (3) and illustrated with the black dash line. In the upper panel, the trapezoidalansatz sums contributions from both incoming nucleai, and is thus not easy to visualize.In each case, the overlap region can be visualized as a measure of the dip at η = 0 foroverall distribution. If the overlap region decreases, then the dip is pronounced and ifthe overlap region increases, the overall distribution becomes flat.The three functions are used to fit the available experimental data displayed inFigs. 4 – 6. The parameter µ of Eq. (3) is set to zero for symmetric collisions but leftunconstrained for asymmetric systems. Fits to pseudorapidity distributions measuredin pp and p–¯p collisions are displayed in Fig. 4; those to Cu–Cu, Au–Au, and Pb–Pbcollisions data are shown in Fig. 5; whereas those to asymmetric systems, d–Au and p–Pbcollisions, are presented in Fig. 6. The upper panels of each figure display experimentaldata and fits obtained with the three functions in distinct panels horizontally, while thelower panels of the figure show ratios of the measured data to the fits. The fits werecarried out with the ROOT least square minimization function [65]. Data uncertaintiesincluded in the fits were set as quadratic sums of statistical and systematic errorsreported by the experiments. The goodness of fit is reported in terms of χ per degreesof freedom ( χ / NDF) in Tab. 1.We find that the three functions fit the pp data relatively well with χ / NDFtypically smaller than 3. However, best fits are achieved with Eqs. (1) and (3). Similarly, ultiplicity and pseudorapidity density distributions − − η / d c h d N η − − D a t a / F i t d+Au : 0 20%p+Pb : 0 5% Eq. (1) − − η − − (2) p+Pb 5020 GeVd+Au 200 GeV Eq. − − η − − Eq.(3)
Figure 6.
Beam energy dependence of charged particle pseudorapidity densitydistributions measured in minimum bias d–Au and p–Pb collisions measured by thePHOBOS collaboration [61] and the ATLAS collaboration [62]. Dashed lines showbest fits obtained with Eqs. (1-3). we find that all three functions provide reasonably sensible parameterizations of the Au–Au, Pb–Pb, and Cu–Cu data presented in Fig. 5. We note, however, that Eq. (1) yieldsfits with smaller deviations from the data, on average, in the central rapidity region.We also find that Eq. (2) does not fully reproduce the dip structure observed at centralrapidity in high-energy datasets.As anticipated, fits with Eq. (1) fail to describe pseudorapidity density distributionsmeasured for asymmetric systems but Eqs. (2,3) produce reasonable fits. Note, however,that deviations in excess of 5% are observed at | η | > . η rangesup to the beam rapidity (where experimental measurement has severe limitations) toestimate the number of produced charge particles. Altogether, we conclude that Eq. (3)provides the best descriptions of the pseudorapidity density distributions, regardless ofcollision system size and energy considered in this work. ultiplicity and pseudorapidity density distributions √ s NN χ /NDF(GeV) Eq. (1) Eq. (2) Eq. (3)pp MB 900 1.056 0.552 0.826pp 2360 0.691 1.367 0.742pp 2760 2.670 14.805 1.495pp 7000 0.458 14.805 1.495pp 8000 1.103 9.320 0.157pp. 13000 0.416 1.312 0.0145p–¯p 630 2.355 2.636 2.144p–¯p 1800 0.986 0.184 0.155CuCu (0-6%) 22.4 1.1806 1.503 1.026CuCu 62.4 0.802 0.778 0.766CuCu 200 0.877 1.095 1.185AuAu 19.6 0.596 0.725 0.592AuAu 62.4 2.184 2.074 0.412AuAu 130 1.889 2.176 0.179AuAu 200 1.103 0.262 0.419PbPb (0-5%) 2760 1.813 1.562 0.943PbPb 5020 1.319 4.216 1.462dAu top 5% 200 No Fit 2.149 3.57pPb mult class 5020 No Fit 3.299 2.118 Table 1. χ /NDF of the fits of data presented in Figs. 4-6 with Eqs. (1-3). MBdenotes minimum bias distribution.
4. Multiplicity distributions at large η : limiting fragmentation Fits with Eqs. (1-3) of pseudorapidity distributions measured in Pb–Pb collisions at2.76 and 5.02 TeV, discussed in the previous section, successfully model the databut are under constrained at large rapidities. In this section, we use the notion oflimiting fragmentation to provide constraints on the shape of these distributions atlarge rapidity. In high-energy hadronic collisions, the limiting fragmentation [26–28, 67]concept stipulates that pseudorapidity densities reach a fixed or universal curve close tothe beam rapidity ( y beam ). This implies that the particle production in the rest frameof one of the colliding hadrons is (approximately) independent of the collision center-of-mass energy. Many explanations have been suggested to interpret this behavior,including gluon saturation in the color glass condensate (CGC) framework [18–20, 66].Indeed, parton distribution functions measured in deep inelastic scattering show that,at very high collision energies, gluons densities largely dominate those of quarks. Thissuggests that the medium produced in these collisions mostly consists of gluons. Withincreasing collision energy, the gluon density increases, eventually leading to saturation.In the previous section, we found that Eq. (3) provides the best fit to the ultiplicity and pseudorapidity density distributions | η | > .
5. In this context, we investigatewhether the notion of limiting fragmentation can further constrain our modeling of theparticle density distributions. - - - beam - y h = h ) h / d c h ( d N - - - beam - y h > ) pa r t / ( . * < N = h ) h / d c h ( d N Eq. (1)Eq. (2)Eq. (3) ALICE 0-5%Pb+Pb2760 GeV5020 GeVXe+Xe5440 GeV
Figure 7.
Limiting fragmentation behavior for Au–Au, Cu–Cu, Xe-Xe and Pb–Pbcollisions at large η -ranges, plotted as a function of η − y beam . The y-axis in the lowerpanel is scaled by the number of participating nucleons pair. Recent studies of limiting fragmentation have shown that Glauber-inspired modelsof particle production in heavy-ion collisions generally fail to reproduce limitingfragmentation [68, 69] behaviour, especially at LHC energies. These studies indicatethat the particle production is a function of the combination of N part and number ofcollisions ( N coll ), as the ratio between the two depends non-trivially on the collisionenergy. Hence, if the nuclei-sized domains are uncorrelated, one generically expectslimiting fragmentation to be broken and which is also true in Color Glass type models.In Ref. [68], the authors have argued that wounded parton models, provided the nucleon ultiplicity and pseudorapidity density distributions x , could in principle reproduceboth multiplicity dependence with energy and limiting fragmentation. The differentcalculations can be verified by studying the limiting fragmentation behaviour of particleproduction by re-plotting the the pseudorapidity density distributions measured incentral Cu–Cu, Au–Au, Xe-Xe and Pb–Pb collisions at RHIC and LHC energies asa function of shifted rapidities, η − y beam .The upper panel of Fig. 7 shows the pseudorapidity distributions for centralcollisions at different colliding energies as a function of η − y beam for Xe–Xe [31] and Pb–Pb [32] systems at the LHC, and Au–Au collisions at the RHIC energies. We observethat the distributions tend to converge towards a single curve close to η − y beam ∼ (cid:104) N part (cid:105) /
2, we obtain the distributions shown in the lower panel ofFig. 7. We observe that the scaled distributions obtained from Xe–Xe, Au–Au, andPb–Pb collisions at several energies closely overlap and more or less follow a universallimiting fragmentation behavior near η − y beam = 0.We further test the notion of limiting fragmentation with fits of the data presentedin Fig. 7 with Eqs. (1-3). Fits of the different data sets, presented in the figure, indeedmerge together near the beam rapidity. In order to further validate the different ans¨atze,the fits were performed by restricting the fit regions and then extrapolating to higher η . This is verified for Xe-Xe collision (at √ s NN = 5.44 TeV) and Pb-Pb collisions (atboth √ s NN =2.76 TeV and √ s NN =5.02 TeV), by fitting the experimental data in theranges (i) | η | ≤ . − < ( η − y beam ) <
3. We find that the extrapolations ofthese fits in the beam rapidity are in near perfect agreement, with a maximum mutualdifference of 1%. We also verified that integrals of the fits, yielding total charged-particlemultiplicity, differ by less than 5%. Additionally, we further checked the validity of thelimiting fragmentation hypothesis by considering fits to two hybrid datasets. Thesehybrid datesets were constructed by joining data points from LHC energies in the range − . < ( η − y beam ) < − . (cid:104) N part (cid:105) scaledvalues from the STAR 200 GeV data points in the range − . < ( η − y beam ) < . N totalch ) production at LHC energies and compare with values obtained at RHICenergies. We discuss the extraction of N totalch in detail in the next section. ultiplicity and pseudorapidity density distributions
5. Total charged-particle multiplicity
We proceed to determine the total charged-particle multiplicity, N ch , produced in Cu–Cu, Au–Au, and Pb–Pb collisions by integration of the fitted pseudorapidity densities,constrained by limiting fragmentation, over the full range of particle production.Figure 8 presents values of N ch scaled by (cid:104) N part (cid:105) / (cid:104) N part (cid:105) for Pb–Pb collisions at 2.76 TeV and Au–Au collisions at 200 GeV. Experimental data pointsreported by the ALICE [70] and PHOBOS [71] collaborations are shown with red andblue dash curves, respectively. Total charged-particle production values are obtained byintegration of the fitted Eq. (1-3) in the range − y beam ≤ η ≤ y beam . Values obtainedwith Eqs. (1), (2), and (3) are shown with solid red, blue and green points, respectively. 〉 part N 〈 〉 pa r t N 〈 / c h N Eq. 3Eq. 2Eq. 1data
Pb+Pb 2760GeV Au+Au 200GeV
Eq. 3Eq. 2Eq. 1data
Figure 8.
Total charge particle multiplicity scaled by the number of participantpair, (cid:104) N part (cid:105) /
2, as a function of (cid:104) N part (cid:105) based on Eqs. (1-3). Red and blue dashlines correspond to data reported by the ALICE and PHOBOS collaboration based onmeasured charge particle multiplicity measured in the range | η | ≤ . − y beam ≤ η ≤ + y beam . The shaded bands represent error bars corresponding tothe correlated systematic uncertainties reported by the experiments [70, 71]. We find that the scaled values of N ch (red triangles and red circles) obtained byintegration of Eq. (3) are consistent, within uncertainties (represented by shaded bands),with those reported by the PHOBOS and ALICE collaborations. Only the N ch valuesobtained at the lowest (cid:104) N part (cid:105) shown appreciably underestimate the PHOBOS data. ultiplicity and pseudorapidity density distributions N ch obtained by integration of Eq. (1) follow a similar trend whilethose obtained with Eq. (2) tend to systematically underestimate the values reportedby PHOBOS. Overall, we find the best agreement with PHOBOS data is achieved usingEq. (3), with deviations of the order of 0.5% compared to 1% with the other twoequations. Hereafter, we use the differences of the three fit extractions as an estimateof the systematic errors associated with the extrapolation procedure based on fits ofEq. (3) to obtain the total charged-particle multiplicities. æ part N Æ c h t o t a l N + b part Fit 1: a N ) (1+bN2 part
NFit 2: a
Au+Au Cu+CuPb+Pb
Figure 9.
Centrality dependence of the total charged-particle multiplicity, estimatedfrom integrals of Eq. (3) across the range − y beam ≤ η ≤ + y beam , in pp, d–Au, p–Pb,Cu–Cu, Au–Au, and Pb–Pb collisions at RHIC and LHC energies. We next proceed to use fits of the measured pseudorapidity distributions withEq. (3) to extract values of N ch for several colliding systems, collision energies, andcollision centralities. Results are shown in Fig. 9 as a function of the number ofparticipating nucleons in pp collisions at 19.6 GeV, 200 GeV and 2.76 TeV, Au–Aucollisions at 19.6, 62.4, 130, and 200 GeV, Cu–Cu collisions at 22.4, 62.4 and 200 GeV,d–Au collisions at 200 GeV, Pb–Pb collisions at 2.76 and 5.02 TeV, and p–Pb collisions ultiplicity and pseudorapidity density distributions N part .We further examine the N ch dependence on N part by considering parameterizationsof this dependence with (a) a linear function aN part + b , and (b) a power law a N part (1+ bN part ), shown in Fig. 9 with black dashed and red solid lines, respectively. Wefind that the power-law parameterization provides a better description of the evolutionof N totalch with N part . Notably, the linear fit fails to describe the evolution of N totalch with N part at LHC energies. Deviations are observed for peripheral collisions with bothparameterizations. Moreover, both the linear and power law functions provide a ratherpoor description of the computed multiplicities in the case of p–Pb collisions. æ part N Æ æ pa r t N Æ / c h t o t a l N p+Pb 5020 GeVd+Au 200 GeV19.6 GeV62.4 GeV130 GeV200 GeV 22.4 GeV62.4 GeV200 GeV Linear extrapolation2760 GeV5020 GeV p+p 2760 GeVp+p 200 GeVp+p 19.6 GeVAu+Au Cu+Cu Pb+Pb 5440 GeVXe+Xe Figure 10.
Centrality dependence of total charged multiplicity per participantpair in pp, d–Au, p–Pb, Cu–Cu, Au–Au, and Pb–Pb collisions at RHIC and LHCenergies [5, 58, 59]. ultiplicity and pseudorapidity density distributions æ part N Æ æ pa r t N Æ / = h ) h / d c h ( d N Figure 11.
Centrality dependence of charged-particle multiplicity density at mid-rapidity in Cu–Cu, Au–Au, Pb–Pb and Xe-Xe collisions at RHIC and LHC energies.
In order to further examine the evolution of N totalch with N part , we plot thecentrality dependence of total charged particle multiplicities scaled by the numberof participant pair in Fig. 10. We observe that for Cu–Cu and Au–Au collisions atRHIC energies, scaled values of N totalch / ( (cid:104) N part (cid:105) /
2) are essentially independent of thecollision centrality, whereas ( dN ch /dη ) η =0 / ( (cid:104) N part (cid:105) / N part in these collision systems. This implies that the shape of the η density distribution changes with centrality and becomes more peaked with increasingcentrality. In contrast, we find that, at LHC energies, both N totalch / ( (cid:104) N part (cid:105) /
2) (Fig. 10)and ( dN ch /dη ) η =0 / ( (cid:104) N part (cid:105) /
2) (Fig. 11) display monotonic increase with N part . ForLHC collisions, the ratio N ch / ( (cid:104) N part (cid:105) /
2) shows a growth, compatible with a power-law behavior. A similar behavior is observed for p–Pb collisions at 5.02 TeV (Fig. 10).The observed violation of participant scaling at LHC energies is in sharp contrastto the near perfect scaling observed at RHIC energies. Furthermore, a scaling violationis observed for both charged-particle multiplicity density at mid-rapidity as well as thetotal number of charged particles. The causes of these violations can be manifolded. ultiplicity and pseudorapidity density distributions x at LHC much lower compared to that atRHIC. At RHIC energies, a transverse mass, m T , of 1 GeV corresponds to x ∼ − at y = 0, whereas at LHC it corresponds to x ∼ − . Bjorken- x values are even lowerat large η . The gluon density is expected to grow and reach saturation with lowering x [74]. At the LHC, one gets to the small x domain where gluon productions dominatesthereby producing large number of additional particles with no relation to the numberof participants. This is consistent with the CGC formalism of the initial state of thecolliding nuclei.Alternatively, particle production at high energy may be described in terms of atwo component model involving soft and hard components, σ total = σ soft + σ hard , in which σ soft represents the cross-section for soft particle production and is proportional to N part ,whereas σ hard , the cross-section for high- p T particle production, is proportional to thenumber of inelastic nucleon-nucleon collisions ( N coll ). A significant increase of σ hard fromRHIC to LHC, relative to σ soft could then possibly explain the observed departure from N part scaling. We use the power law obtained in the previous section to “predict” the total charged-particle production as a function of the number of participants at the FAIR and NICAfacilities, expected to become online in 2022 and 2025, respectively. To calculate thesepredictions, we first remark that the shape of the (cid:104) N part (cid:105) dependence of the centralrapidity particle density for RHIC energies is essentially invariant with respect to √ s NN .To illustrate this approximate invariance, we plot central multiplicity densities scaled tothe corresponding multiplicity density at √ s NN = 200 GeV as a function of (cid:104) N part (cid:105) forseveral collision systems and energies in Fig. 12. The scaling factors were determined asthe ratio of multiplicity density at central rapidity measured at different beam energies √ s NN to the multiplicity density observed at central rapidity in √ s NN = 200 GeV Au–Au collisions. These are listed for each collision system and energy in the upper panelof the figure. The scaled densities are compared to the CGC initial condition motivatedfit (discussed in the next section) to the data at √ s NN = 200 GeV, shown as a blue dashline. We observe from Fig. 12 that the overlap of the data points is reasonable atenergies lower than √ s NN = 200 GeV, which makes it possible to predict the particledensities at lower collision energies. The scaling factors are plotted as a function of √ s NN in the lower panel of the figure and fitted with a first order polynomial shownby the red dash line. We extract the coefficients a and b , and use these to obtainscaling factors for NICA and FAIR energies. These scaling factors are used to obtainpredictions of collision centrality evolution of the central particle density per participant, dN ch /dη | η =0 / (cid:104) N part (cid:105) /
2. This is shown in Fig. 13 as a function of (cid:104) N part (cid:105) . ultiplicity and pseudorapidity density distributions 〉 part N 〈 D a t a sc a l ed b y G e V × × ×
130 GeV 0.922 ×
200 GeV 1.278 × × × × × p Pb 1 × d Au 200 GeV in GeV NN s10 S c a li ng F a c t o r +c b f= axa= 4.205b= 0.302c= 0.139/ndf = 1.067 χ Figure 12. (Upper) Centrality dependence of charged-particle multiplicity densityscaled to that of AuAu collisions at √ s NN = 200 GeV. (Lower) Scaling factors forcharged-particle multiplicity density to the data at 200 GeV. The High-Energy Large Hadron Collider (HE-LHC) [75] and the FCC [76] acceleratorsproposed at CERN will achieve unprecedented large collision energies for pp as well asheavy-ions. The expected energies for Pb-Pb collisions are 11 TeV and 39 TeV for HE-LHC and FCC, respectively. It is thus imperative to make predictions for the numberof produced particles at such large energies. The scaling technique used to extrapolatethe particle multiplicities for collision energies lower than √ s NN = 200 GeV is notappropriate for extending to higher energies as the approximate N part scaling is broken(as per figure 11). The indication of the scale breaking for multiplicity density at mid-rapidity is also evident by a closer look to upper panel of Fig. 12 at √ s NN = 2760 GeVand √ s NN = 5020 GeV for N part > . s . ± . NN ), we canpredict the charged particle multiplicity densities at mid-rapidity for Pb–Pb collisionsat 11 TeV and 39 TeV. The extrapolation gives (cid:104) N part (cid:105) dN ch /dη | η =0 as 13 . ± .
504 and19 . ± . η = 0 are ≈ ±
93 and ≈ ± ultiplicity and pseudorapidity density distributions 〉 part N 〈 〉 pa r t N 〈 / = η ) η / d c h ( d N
20A ( 6.27 ) GeV 30A ( 7.62 ) GeV45A ( 8.96 ) GeV
Expectation of CBM (FAIR) : Fixed Target Pb+Pb ) NN s ( Lab E NN s& NICA : collider Au+Au Figure 13.
Expected evolution of charged-particle multiplicity density with centralityfor CBM (FAIR) energies. one should be careful by considering the present knowledge of gluon saturation picture,which could limit the particle production in these energies and push the multiplicitytowards a lower value than expected from these extrapolation.
6. Multiplicity density from initial condition motivated models
The collision centrality dependence of the ratio N ch / (cid:104) N part (cid:105) is expected to be somewhatsensitive to the initial state conditions of heavy-ion collisions [79, 80]. The measuredevolution of charged-particle multiplicity distributions vs. collision centrality, presentedin Fig. 14 for selected collision systems, may thus be used to contrast predictionsobtained with different models. We focus our discussion on the Glauber [52, 53] andcolor glass condensate [79, 80] models.Within the Glauber model, a soft/hard two-component model is used toparameterize the particle production as a function of collision centrality according to dN ch dη (cid:12)(cid:12)(cid:12) AA = n pp (cid:20) (1 − x ) N part xN coll (cid:21) , (4)where N part and N coll represent the number of soft and hard scatterings, respectively, and n pp denotes the average number of produced charged particles per unit pseudorapidity ultiplicity and pseudorapidity density distributions æ pa r t N Æ / = h ) h / d c h ( d N æ part N Æ D a t a / F i t Glauber(a) æ part N Æ
200 GeV130 GeV62.4 GeV19.6 GeV2760 GeV p+p
CGC (b) æ part N Æ
200 GeV130 GeV62.4 GeV19.6 GeV
Au+Au
EKRT (c)
Figure 14.
Parameterization of the N part dependence of charged-particle multiplicitydensity per participant pair for symmetric collision systems fitted with initial conditionsaccording to (a) Glauber, (b) CGC, and (c) EKRT models. in pp collisions. The variable x , representing the fraction of hard collisions, is hereconsidered a fit parameter. The fit results of hard scattering component x , is within therange of 0.10 to 0.16 and in agreement with previous measurements. Panel (a) of Fig. 14displays fits (green dash lines) of data from Cu–Cu, Au–Au, and Pb–Pb collisions acrossa wide span of beam energies. To carry out the fits, we evaluated values of n pp vs. √ s based on the parameterization, n pp ∝ s . NN , presented for (NSD) p–p collisions in Fig. 1.In the context of the Color Glass Condensate model, one expects that small x gluons overlap and recombine thereby reducing the overall number of gluons and thenumber of hadrons they hadronize into. The charged-particle density is hence modeledaccording to [19]: dN ch dη ≈ N αpart ( √ s NN ) γ , (5)where α and γ are free parameters. Fits based on this model are shown in Fig. 14(b). By contrast, models based on final state gluon saturation, e.g., EKRT [77],predict a decreasing trend in charged-particle multiplicities per participant nucleon withincreasing collision centrality according to dN ch dη = C
23 1 . (cid:16) N part (cid:17) . ( √ s NN ) . , (6)where C is the only free parameter. While the Glauber and CGC initial conditionsparameterizations shown in panels (a) and (b) provide excellent agreement with ultiplicity and pseudorapidity density distributions 〉 part N 〈 〉 part N 〈 〉 part N 〈 D a t a /I C GlauberCGCEKRTGlauberCGCEKRT 〉 pa r t N 〈 / = η ) η / d c h ( d N Initial Condition (IC) d+Au 200 GeVp+Pb 5020 GeV
Figure 15. d+Au and p+Pb asymmetric collisions fitted with different initialconditions according to Glauber, CGC and EKRT models. measured data, one finds fits based on Eq. (6), presented in Fig. 14 (c) are in starkdisagreement with the data, owing evidently to the fixed N part power smaller thanunity.We extend this study to d–Au and p–Pb collision systems in Fig. 15 using theparameterizations (4-6). We find that, in these two systems, the soft/hard two-component model and the EKRT Eq. (6) provide a relatively poor representation of thedata. Overall then, we conclude the CGC inspired parameterization, Eq. (5), providesa suitable description of the evolution of the charged-particle multiplicity density with N part in both symmetric and asymmetric collision systems.However, we note that recent event-by-event calculations carried out using next-to-leading order EKRT model [78], with saturation for soft particle production andviscous hydrodynamics for the space-time evolution of the produced matter, can welldescribe the multiplicity density discussed above. In addition, the recent theoreticaldevelopment on initial conditions known as TRENTO [81, 82] initial conditions alsoprovides a successful description of the densities (as well as several other observables)in p–p, p–Pb, Au–Au, and Pb–Pb collisions both at RHIC and LHC energies. ultiplicity and pseudorapidity density distributions
7. Summary
We have presented a comprehensive study of the multiplicity and pseudorapiditydistributions of charged particles produced in p–p, p–Pb, d–Au, Cu–Cu, Au–Au, andPb–Pb collisions at energies ranging from a few GeV to several TeV, correspondingto the available experimental data at RHIC and LHC. The experimental data havebeen compared to calculations of selected event generators, including PYTHIA,EPOS, AMPT, UrQMD, and THERMINATOR, which feature different physics modelassumptions. We find these event generators qualitatively reproduce the observedparticle densities at | η | = 0. However, none are able to satisfactorily explain measureddistributions over a broad range of pseudorapidities. With the goal of extrapolating themeasured data to forward rapidities, to estimate the total charged particle productionin various collision systems, and to obtain the dependence on the collision energy,we have studied three different functional forms to describe the experimental data onthe pseudorapidity distributions. Among these functional forms, the difference of twoGaussian distributions, Eq. (3), is found to best reproduce the measured multiplicitydensities observed in different collision systems and collision energies.Furthermore, we used Eq. (3) to estimate the total charged-particle production andstudy the evolution of multiplicity density at central rapidity ( dN ch /dη/ (cid:104) N part (cid:105) / | η =0 )as a function of collision centrality and collision energy. At beam energies √ s NN ≤
200 GeV, the charged-particle rapidity density exhibits a modest increase with (cid:104) N part (cid:105) while the total charge production is approximately independent of collision centrality.In contrast, at LHC energies, both the particle density at mid-rapidity and the totalcharge particle production exhibit a rapid increase with (cid:104) N part (cid:105) . We thus conclude thatthere is a qualitative change in the particle production at LHC relative to RHIC. AtRHIC energies, the multiplicity density at mid-rapidity increases with (cid:104) N part (cid:105) while thetotal particle production remains fixed. That implies the pseudorapidity distributionnarrows with increasing (cid:104) N part (cid:105) thereby yielding a larger central rapidity density albeitwith a fixed integral. At the LHC, by contrast, both the central rapidity density and thetotal charged-particle production increase with (cid:104) N part (cid:105) . One then has entered a differentregime of particle production in which both the central rapidity and total multiplicitiesper participant monotonically increase with (cid:104) N part (cid:105) .We found that the limiting fragmentation hypothesis holds at the TeV energy scaleand thus can be used to approximately constrain the shape of dN/dη distributionsand their integrals over the full range of particle production. In addition, we havestudied charged-particle multiplicity productions considering different initial conditions.We observe that CGC like initial condition is best suited to describe the publisheddata for both symmetric and asymmetric type of collisions. We have extended theparticle production studies to lower collision energies corresponding to those of upcomingaccelerator facilities of FAIR at GSI, Darmstadt and NICA at JINR Dubna. We haveextrapolated the charged particle multiplicity densities at η = 0 for expected heavy-ioncollisions at the proposed HE-LHC and FCC at CERN. ultiplicity and pseudorapidity density distributions Acknowledgments:
SB & CP acknowledge the financial support by the U.S. Department of Energy Officeof Science, Office of Nuclear Physics under Award Number DE-FG02-92ER-40713.
References [1] A. Biallas, M. Bleszynski, and W. Czyz, Nucl. Phys. B , 461 (1976).[2] J. D. Bjorken, Phys. Rev. D , 140 (1983).[3] D. Kharzeev and M. Nardi, Phys. Lett. B , 121 (2001).[4] J.F. Grosse-Oetringhaus and K. Reygers, J. Phys. G: Nucl. Part. Phys. , 083001 (2010).[5] B. Alver et al. [PHOBOS Collaboration],Phys. Rev. C 83 (2011) 024913.[6] L. Adamczyk et al. (STAR Collaboration) Phys. Rev. C , 044904 (2017).[7] A. Toia et al. (ALICE Collaboration) J. Phys. G , 124007 (2011).[8] R. Sahoo, A. N. Mishra, N. K. Behera and B. K. Nandi, Adv. high-energy Phys. , 612390(2015).[9] E. K. G. Sarkisyan, A. N. Mishra, R. Sahoo and A. S. Sakharov, Phys. Rev. D , 054046 (2016).[10] S. Basu, T.K. Nayak, K. Datta, Phys. Rev. C , 064902 (2016).[11] Torbjorn Sj¨ostrand, Stephen Mrenna, and Peter Skands, J. High Energy Phys. , 026 (2006).[12] G.Corcella et. al. Jour. High Ener. Physics. (2001)010.[13] J. Adams et al. (STAR Collaboration) Nucl. Phys. A , 102 (2005).[14] U. W. Heinz and M. Jacob, arXiv:nucl-th/0002042 [nucl-th].[15] L. McLerran, J. Dunlop, D. Morrison, R. Venugopalan, Nucl. Phys. A (2011) 1.[16] E. Iancu, K. Itakura, L. McLerran, Nucl. Phys. A (2002) 327.[17] E. Iancu, L. McLerran, Phys. Lett. B (2001) 145.[18] Jamal Jalilian-Marian Phys. Rev. C 70 , 027902 (2004).[19] D. Kharzeev and E. Levin, Phys. Lett. B , 79 (2001).[20] K. Dusling and R. Venugopalan, Phys. Rev. Lett. (2012) 262001.[21] F.O. Duraes, A.V. Giannini, V.p. Goucalves and F.S. Navarra, Phys. Rev. C 94, 024917 (2016).[22] T. Lappi, H. Mantysaari, Nucl. Phys. A 926, 186 (2014).[23] A. H. Rezaeian, Phys. Lett. B 727, 218 (2013).[24] M. Derrick et al. (ZEUS Collaboration) Phys. Lett. B (1995) 576.[25] C. Adlof et al. (H1 Collaboration), Phys. Lett. B (2001) 183.[26] J. Benecke, T. T. Chou, C. N. Yang and E. Yen, Phys. Rev. , 2159 (1969).[27] R. Beckmann, S. Raha, N. Stelte, R.M. Weiner, Phys. Lett.
B105 , 411 (1981).[28] F. Gelis, A.N. Stasto, R. Venugopalan, Eur. Phys. J. C , 489 (2006).[29] B. Alver et al. (PHOBOS Collaboration) Phys. Rev. C , 024903 (2016).[30] G. Torrieri, EPJ Web of Conferences, 04002 (2011).[31] S. Acharya et al. (ALICE Collaboration) arXiv:1805.04432 [nucl-ex].[32] J. Adam et al. (ALICE Collaboration) Phys. Lett B , 567 (2017).[33] L. D. McLerran, Raju Venugopalan, Phys. Rev. D (1994) 2233.[34] M. Bleicher et al., , Phys. Lett. B 435 , 9 (1998).[35] M. Bleicher, S. Jeon, V. Koch, Phys. Rev. C , 061902(R) (2000)[36] S. Haussler, H. Stocker and M. Bleicher, Phys. Rev. C , 021901(R) (2006).[37] A. Chatterjee, S. Chatterjee, T. K. Nayak, N. R. Sahoo, J. Phys. G: Nucl. Phys. J. Phys. ,125103 (2016)[38] Z.-W. Lin et al. Phys. Rev. C , 064901 (2005).[39] Z.-W. Lin et al. , Phys. Rev. C , 011902 (2001).[40] B. Zhang et al. , Phys. Rev. C , 067901 (2000).[41] A.G. Knospe et al. Phys. Rev. C , 014911 (2016).[42] K. Werner et al. , Phys. Rev. C , 064903 (2014).[43] K. Werner, Phys. Rev. Lett. , 152301 (2007). ultiplicity and pseudorapidity density distributions [44] K. Werner et al. , Phys. Rev. C , 044904 (2010).[45] M. Nahrgang et al. , Phys. Rev. C , 024907 (2014).[46] M. Chojnacki, A. Kisiel, W. Florkowski and W. Broniowski, Comput. Phys. Commun. (2012)746.[47] K. Aamodt et al. [ALICE Collaboration], Eur. Phys. J. C (2010) 89.[48] J. Adam et al. [ALICE Collaboration], Eur. Phys. J. C (2017) no.1, 33.[49] J. Adam et al. [ALICE Collaboration], Phys. Lett. B (2016) 319.[50] F. Abe et al. [CDF Collaboration], Phys. Rev. D (1990) 2330.[51] K. Alpgard et al. [UA5 Collaboration] Phys. Lett.B (1982) 183, G. J. Alner et al. [UA5Collaboration] Phys. Rept. (1987) 247.[52] B. I. Abelev et al. [STAR Collaboration], Phys. Rev. C (2009) 034909.[53] B. Alver et al. [PHOBOS Collaboration], Phys. Rev. Lett. (2009) 142301.[54] B. B. Back et al. [PHOBOS Collaboration], Phys. Rev. C (2006) 021901.[55] B. Alver et al. [PHOBOS Collaboration], Phys. Rev. C (2011) 024913[56] B. B. Back et al. [PHOBOS Collaboration], Phys. Rev. Lett. (2001) 102303[57] B. B. Back et al. [PHOBOS Collaboration], Phys. Rev. Lett. (2003) 052303[58] E. Abbas et al. [ALICE Collaboration], Phys. Lett. B (2013) 610.[59] J. Adam et al. [ALICE Collaboration], Phys. Lett. B 772 (2017) 567.[60] J. Adam et al. [ALICE Collaboration], Phys. Rev. Lett. (2016) 222302.[61] B. B. Back et al. [PHOBOS Collaboration], Phys. Rev. C (2005) 031901[62] G. Aad et al. [ATLAS Collaboration], Eur. Phys. J. C (2016) no.4, 199[63] B. Abelev et al. [ALICE Collaboration], Phys. Rev. Lett. (2013) no.3, 032301[64] J. Adam et al. [ALICE Collaboration], Phys. Rev. C (2015) no.6, 064905.[65] R. Brun and F. Rademakers, Nucl. Inst. and Meth. in Phys. Res. A (1997) 81.[66] A. Krasnitz and R. Venugopalan, Nucl. Phys. A (2002) 209.[67] B. Kellers, G. Wolschin, Prog. Theor. Exp. Phys. (2019) 053D03.[68] K. J. Gon¸calves, A. V. Giannini, D. D. Chinellato and G. Torrieri, Phys. Rev. C (2019)054901.[69] P. Sahoo, P. Pareek, S. K. Tiwari and R. Sahoo, Phys. Rev. C (2019) 044906.[70] J. Adam et al. [ALICE Collaboration], ‘Phys. Lett. B (2016), 373.[71] B. B. Back et al. [PHOBOS], Phys. Rev. C (2006), 021902.[72] V. Khachatryan et al. [CMS Collaboration] Jour. High Ene. Phys. (2010) 041.[73] G. Wolschin, Phys. Rev. C , 014905 (2015).[74] L. A. Harland-Lang, A. D. Martin, P. Motylinski, and R.S. Thorne, Eur.Phys.J. C (2015) no.5,204[75] A. Abada et al. [FCC Collaboration], Eur. Phys. J. Special Topics (2019) 1109.[76] A. Dainese et al. PoS HardProbes2018 (2019) 005 (e-Print: 1901.10952 [hep-ph]).[77] K. J. Eskola, K. Kajantie, P. V. Ruuskanen, and K.Tuominen, Nucl. Phys.
B 570 , 379 (2000).[78] K. J. Eskola, H. Niemi and R. Paatelainen, Nucl. and Part. Phys. Proc., (2016) 161.[79] J.L. Albacete, C. Marquet, Progress in Particle and Nuclear Physics 76, 1 (2014).[80] J. L. Albacete, A. Dumitru, Y. Nara, J. Phys.: Conf. Ser. , 012011 (2011).[81] J. Scott Moreland, Jonah E. Bernhard, and Steffen A. Bass, Phys. Rev. C (2015) 011901.[82] C. Shen et al. Comput. Phys. Commun.199