Revisiting the centrality definition and observable centrality dependence
of relativistic heavy-ion collisions in PACIAE model

The ultra-relativistic heavy-ion collisions produce strongly interacting quark-gluon matter (sQGM) under extreme conditions of temperature and energy density at Relativistic Heavy Ion Collider (RHIC) [1, 2, 3, 4] and Large Hadron Collider (LHC) [5, 6, 7]. The centrality is one key physical character in studying the high energy heavy-ion collisions, because it is directly related to the interacting volume (overlap zone) of the collision system. This overlap zone depends on the impact parameter $b$ defined as the distance between the centers of two colliding nuclei in the plane transverse to the beam axis. The centrality of a nucleus-nucleus (AA) collision with impact parameter $b$ is usually defined as a percentile $c$ in the nucleus-nucleus total cross section $\sigma_{AA}$ [8]:

In the experiment, this centrality percentile $c$ of the nucleus-nucleus total cross section is usually assumed to be approximately equivalent to the fraction of charged particle multiplicity above a multiplicity cut of $N_{ch}^{cut}$, or to the energy deposited in the zero-degree calorimeter (ZDC) below a cut of $E_{ZDC}^{cut}$[8]:

The nucleus-nucleus total cross section in Eq. (1) is calculated by

i.e. by the nucleus-nucleus geometrical total cross section ($\pi b_{max}^{2}$) corrected with the fraction of events with at least one nucleon-nucleon collision. Meanwhile, the centrality percentile $c$ in a nucleus-nucleus collision is also assumed to be equivalent to the fraction of impact parameter distribution ($f(b)\propto bdb$). Therefore, a mapping relation of

is obtained. In above equation $c$ refers to the centrality percentile and $b_{max}$ is assumed to be

where the $R_{A}$, for instance, is radius of nucleus $A$ (it denotes the atomic number of nucleus either), the $d=0.546fm$ refers to the tail of the nuclear density profile, and the coefficient $f$ is a free parameter.

The PACIAE model [9, 10, 11] is a parton and hadron cascade model based on PYTHIA [12]. For nucleon-nucleon (NN) collisions, with respect to PYTHIA, the partonic and hadronic rescatterings are introduced before and after the hadronization, respectively. The final hadronic state is developed from the initial partonic hard scattering and parton showers, followed by parton rescattering, string fragmentation, and hadron rescattering stages. Thus, the PACIAE model provides a multi-stage transport description on the evolution of the collision system.

For a nucleus-nucleus (AA) collisions, the initial positions of nucleons in the colliding nuclei are sampled according to the Woods-Saxon distribution. Together with the initial momentum setup of $p_{x}=p_{y}=0$ and $p_{z}=p_{\rm beam}$ for each nucleon, a list containing the initial state of all nucleons in AA collision is constructed. A collision happens between two nucleons from different nucleus if their relative transverse distance is less than or equal to the minimum approaching distance: $D\leq\sqrt{\sigma_{\rm NN}^{\rm tot}/\pi}$. The collision time is calculated with the assumption of straight-line trajectories. All such nucleon pairs compose a nucleon-nucleon (NN) collision time list. The NN collision with least collision time is selected from the list and executed by PYTHIA (PYEVNW subroutine) with the hadronization temporarily turned-off, as well as the strings and diquarks broken-up. The nucleon list and NN collision time list are then updated. A new NN collision with least collision time is selected from the updated NN collision time list and executed by PYTHIA. With repeating the aforementioned steps till the NN collision list empty, the initial partonic state is constructed for a AA collision.

Then, the partonic rescatterings are performed, where the LO-pQCD parton-parton cross section [13, 14] is employed. After partonic rescattering, the string is recovered and then hadronized with the Lund string fragmentation scheme resulting in an intermediate hadronic state. Finally, the system proceeds into the hadronic rescattering stage and produces the final hadronic state observed in the experiments.

Thus PACIAE Monte-Carlo simulation provides a complete description of the NN and/or AA collisions, which includes the partonic initialization stage, partonic rescattering stage, hadronization stage, and the hadronic rescattering stage. Meanwhile, the PACIAE model simulation could be selected to stop at any stage desired conveniently. In this work, the simulations are stopped after the hadronic rescattering stage, i.e. at the final hadronic state.

In the PACIAE 2.0 and successors, we introduce the geometrical model of Eq.(4) and set the coefficient of $f$ to be equal to 2 and 1 for nucleus-nucleus and proton nucleus collisions [9], respectively. The results of this geometrical model were very consistent with the STAR/RHIC centrality definitions [9, 15], because of they have the same $b_{max}$ definition. Later on, the ALICE, ATLAS and CMS collaborations at LHC have observed that the $b_{max}$ should be extended to 20 $fm$, where the interaction is just approaching to zero [8]. This 20 $fm$ is larger than the $b_{max}$ calculated by Eq.(4) with above $f$ coefficient setting. A new set of $f$ coefficient is then required.

In the improved Monte Carlo Glauber (MC-Glauber, MCG) model simulation [16] the impact parameter $b$ bin, corresponding to a given centrality percentile bin, is sliced according to impact parameter distribution up to a $b_{cut}$, which is assumed to be at somewhere between $R_{A}+R_{B}$ and 20 $fm$. But the last $b$ bin is just assumed to be extended from $b_{cut}$ to the $20fm$. For the Pb-Pb, Xe-Xe, Au-Au and Cu-Cu collisions at the relativistic energies, the $b_{cut}$ is assumed to be approximately equal to 15.6, 13.8, 14.9, and 11.0 $fm$ [16], respectively. Substituting these $b_{cut}$ values into left side of Eq.(5) individually, the corresponding $f$ value of 4.21, 4.45, 3.42, and 3.74 is obtained sequentially. Therefore we assume the parameter $f$ in Eq. (5) equals to 4 for the nucleus-nucleus collision but equals to 2 for the proton-nucleus collision, and we just slice the last $b$ bin up to $b_{cut}$ in this work.

The participant nucleons number $N_{\rm part}$ and binary nucleon-nucleon collisions number $N_{\rm coll}$ are commonly used to represent the collision centrality. In MC-Glauber model, they are counted as the number of binary nucleon collisions happened and the number of nucleons which suffer at least one collision (wounded nucleons), respectively, in a nucleus-nucleus collision simulation within the boundary of $b_{cut}$. The $\langle N_{\rm{part}}\rangle$ and $\langle N_{\rm{coll}}\rangle$ denote their average value over events. And the nuclear overlap function $\langle T_{AA}\rangle$ is calculated by

where $\sigma_{NN}^{inel}$ is the inelastic nucleon-nucleon (NN) cross section.

In the optical Glauber model[17] used in PACIAE, the $N_{\rm part}$ and $T_{AB}$ are analytically calculated by

for the asymmetric AB collisions. In Eq.(8) the $\vec{s}$ refers to a vector of an area perpendicular to the beam axis $z$, $s=|\vec{s}|$, and $\rho(\vec{s},z)$ stands for the nuclear density in the volume element of $d^{2}sdz$ at point of $(\vec{s},z)$. The nuclear density distribution in a nucleus is assumed to be the spherically symmetric Woods-Saxon density distribution [15]

One assumes further:

[18]. The notarization relation of

is then required.

The relation of Eq. 5 with new set of $f$ coefficient, together with the Eqs. 8-10 are used in the PACIAE 2.2.2 to calculate the $b_{max}$ and $b$ bins as well as $N_{\rm coll}$ and $N_{\rm part}$ in the Pb-Pb and p-Pb collisions at $\sqrt{s_{NN}}$=5.02 TeV, Xe-Xe collisions at $\sqrt{s_{NN}}$=5.44 TeV, as well as Au-Au and the Cu-Cu collisions at $\sqrt{s_{NN}}$=0.2 TeV, respectively. The nucleon-nucleon inelastic cross section $\sigma_{NN}^{inel}$ are set as 41.6, 67.6 and 68.4 mb for $\sqrt{s_{NN}}$=0.2, 5.02 and 5.44 TeV[16], respectively. The results are given in the Tabs. 1-5. In these table, the corresponding improved MC-Glauber model results [16] are also given for the comparisons. We see in these tables that, the optical Glauber model results are well consistent with corresponding MC-Glauber ones within the error bar, except the most peripheral collisions.

The PACIAE model simulated charged-particle pseudorapidity distributions (open symbols) in Pb-Pb collision at $\sqrt{s_{NN}}$=5.02 TeV are compared with the corresponding ALICE data [19] (solid symbols) in Fig. 1 for ten centrality bins of 0-5%, 5-10%, 10-20%, …, 80-90%. This figure shows that the PACIAE 2.2.2 model results well reproduce the ALICE data from central to peripheral Pb-Pb collisions. In the $|\eta|>4$ region, the theoretical results are smaller than the experimental data. This has to by studied further in the next work by adjusting the parameters in Lund fragmentation function and/or trying other fragmentation functions.

We compare the PACIAE model simulated charged-particle transverse momentum distribution in Pb-Pb collision at $\sqrt{s_{NN}}$=5.02 TeV to the corresponding ALICE data [20] in Fig 2 for centrality bins of 0-5%, 5-10%, …, 70-80%. For better visibility, both the PACIAE results and the ALICE data, except 70-80% one, are rescaled by a factor of 10, $10^{2}$, …, $10^{8}$, respectively. In this figure we see the PACIAE model well reproduce the corresponding ALICE data in $p\rm{{}_{T}}<$6 GeV/c region. However in the $p\rm{{}_{T}}>$6 GeV/c region, the theoretical results are smaller than the experimental data. This might because the $p\rm{{}_{T}}$ distribution in heavy-ion collisions change from exponential-like at low $p\rm{{}_{T}}$ region to a power-law shape at high $p\rm{{}_{T}}$ region [21], and the $p\rm{{}_{T}}$ distribution is sampled by exponential-like distribution in this work. In the next study we shall try to sample particle $p\rm{{}_{T}}$ by power law distribution instead of exponential-like one in the $p\rm{{}_{T}}>$ 6 GeV/c region.

The $f$ coefficient in the impact parameter $b_{max}$ formula in PACIAE 2.2.2 model is reset to respond the ALICE, ATLAS, and CMS’s observation of the maximum impact parameter in heavy-ion collisions at relativistic energies should be extended to 20 $fm$. Consequently, the PACIAE model is updated to a new version of PACIAE 2.2.2 with convenience of studying the elementary nuclear collisions, proton-nucleus collisions, and the nucleus-nucleus collisions in an unified program version. The PACIAE model calculated impact parameter bin and the optical Glauber results of $N_{\rm part}$ and $N_{\rm coll}$ are well consistent with the corresponding improved MC-Glauber ones in the Pb-Pb, $p$-Pb, Xe-Xe, Au-Au, and Cu-Cu collisions at relativistic energies. PACIAE model simulated charged-particle pseudorapidity and transverse momentum distributions in Pb-Pb collisions at $\sqrt{s_{NN}}$=5.02 TeV well reproduce the corresponding ALICE data. These results prove the correction of the impact parameter centrality definition in the updated PACIAE model and the efficiency of PACIAE model itself.

Acknowledgments This work was supported by the National Natural Science Foundation of China under grant Nos.: 11775094, and 11905188 and by the 111 project of the foreign expert bureau of China. YLY acknowledge the financial support from the Continuous Basic Scientific Research Project (No, WDJC-2019-13) in CIAE.