Transverse Momentum Distribution and Elliptic Flow of Charged Hadrons in $U$+$U$ collisions at $\sqrt{s_{NN}}=193$ GeV using HYDJET++

Recent experimental observations of the charged hadron properties in $U+U$ collisions at $193$ GeV contradict many of the theoretical models of particle production including two-component Monte Carlo Glauber model. The experimental results show a small correlation between the charged hadron properties and the initial geometrical configurations (e.g. body-body, tip-tip etc.) of $U+U$ collisions. In this article, we have modified the Monte Carlo HYDJET++ model to study the charged hadron production in $U+U$ collisions at $193$ GeV center-of-mass energy in tip-tip and body-body initial configurations. We have modified the hard as well as soft production processes to make this model suitable for $U+U$ collisions. We have calculated the pseudorapidity distribution, transverse momentum distribution and elliptic flow distribution of charged hadrons with different control parameters in various geometrical configurations possible for $U+U$ collision. We find that HYDJET++ model supports a small correlation between the various properties of charged hadrons and the initial geometrical configurations of $U+U$ collision. Further, the results obtained in modified HYDJET++ model regarding $dn_{ch}/d\eta$ and elliptic flow ($v_{2}$) suitably matches with the experimental data of $U+U$ collisions in minimum bias configuration.


I. INTRODUCTION
The basic motivation of heavy ion collision experiments is to understand the properties and behaviour of quantum chromodynamics (QCD) at very high temperature and chemical potentials via analysing the data on multi-particle production and by matching experimental measurements to the simulation models for the entire evolution of the fireball. There are existing computational models which use the theoretical or phenomenological foundation of strong interactions to mimic the space-time evolution of collision experiments. One can broadly classify these models in two types: dynamical models [1][2][3][4][5][6][7][8] and semi dynamical models [9][10][11]. Dynamical models are those which consider the pre-equilibrium evolution as well as post equilibrium hydrodynamic evolution like IP-Glasma model etc [1][2][3][4][5][6][7][8]. However, most of the models are semi dynamical models which use a static initial condition at proper thermalization time and then evolve the system using viscous or ideal hydrodynamics like AMPT, MC-Glauber etc [9][10][11]. The particle production mechanism of both types of model are quite different. In dynamical models, the parton saturation is a viable mechanism for particle production e.g., IP-Glasma model is based on the ab-initio color glass condensate framework which combines the impact parameter dependent saturation model for parton distributions with an event-by-event classical Yang-Mills description of early-tile glasma fields [1]. Similarly EKRT model is based on the assumption of final state gluon saturation and thus the initial energy density and produced number of partons scales with atomic number and beam energy [2,3]. In KLN model, the inclusive production of partons is driven by the parton saturation in strong gluon fields [4,5]. In saturation regime, the multiplicity of produced partons should be proportional to atomic number [4,5]. On the other hand the particle production mechanism in semi-classical models are implemented via some phenomenological parameterization or using Monte Carlo event generator e.g., in MC-Glauber model, the particle production is based on static initial conditions and two-component parameterization in which first term is proportional to mean number of participants and second term is proportional to mean number of collisions [9]. In AMPT model initial conditions are obtained from HIJING event generator then ZPC for parton scatterings. After that Lund string model for hadronization and ART model to treat the hadronic scatterings [10]. UrQMD model describes the particle production at low and intermediate energies in terms of scatterings amongst hadrons and their resonances. At higher energies, the excitation of colour strings and their subsequent fragmentation is the particle production mechanism in this model [11].
Most of the simulation models are successful in providing the multiplicity of charged hadrons produced in various heavy ion collision experiments. Vast experimental data on multi-particle production and distributions with collision control parameters like centrality, rapidity and/or transverse momentum etc., put a stringent constraint on these models so that we can understand the production mechanism more deeply and make our models more realistic. To strengthen our understanding about quantum chromodynamics (QCD), these collider experiments collide various nuclei at different colliding energies. Recently RHIC experiment has collided uranium (U) nuclei at the center-of-mass energy √ s N N = 193 GeV [12]. As we know that uranium is a deformed nuclei (prolate in shape) so various kind of initial configurations are possible in U + U collision e.g., body-body, tip-tip, body-tip etc.
The various computational models previously predicted a large difference in multiplicity and elliptic flow between body-body and tip-tip configurations of U + U collisions [13,14].
However, the experimental data of multi-particle production in U + U collisions regarding multiplicity and elliptic flow (v 2 ) contradicts the earlier expectations of most of these computational and theoretical models and shows a small correlation between multiplicity (and/or v 2 ) and initial configurations of U + U collision [12]. This contradiction may have two possible reasons. Either the simulation models have something missing or experimentally we are not quite able to disentangle the events with different geometrical orientations. Thus we have to work on both the aspects since U + U collision in its various orientations is quite useful to understand wide range of physics. Quark gluon plasma (QGP) phase which is characterized by the observables like elliptic flow, jet quenching, charmonia suppression and multiplicity can be better understood in the collision of deformed uranium nuclei due to its initial geometry and specific orientation [14][15][16][17][18][19]. Further U + U collisions can provide a reliable tool to subtract the background elliptic flow effect from the signal so that one can detect the chiral magnetic effect (CME) [14]. In spherical nuclei, it is difficult to disentangle both these effect since the strength of both the signals generated from elliptic flow and CME is of similar strength in peripheral collisions. However, in U + U central collisions, the different geometrical orientations can provide a way to subtract the background signal from CME signal due to a measurable difference in their strength. Thus central collisions of U + U nuclei in tip-tip configuration can possibly be a good tool to characterize the signal of CME [13,20].
Very recently different methods have been proposed to modify some of the models to incorporate the experimental U + U observations in that particular simulation models [21][22][23]. The constituent quark model is also proposed to describe the experimental observation [24,25]. In this article we want to study the U + U collision at √ s N N = 193 GeV in body-body and tip-tip configurations by modifying HYDJET++ model which uses PYTHIA type initial condition for hard part and Glauber type initial condition for soft part. Further most of the existing models either consist of high p T particle production from jet fragmentation or involve low p T hadron production using thermal statistical processes. However, HYDJET++ model [26] consistently includes production of hard as well as soft p T hadrons, to calculate the charged hadron production in U + U collisions at The details on physics model and simulation procedure of HYDJET++ can be found in the corresponding manual [26,27]. The main features of HYDJET++ model are listed very briefly in this section.

A. Hard multi-jet production
The model for the hard multi-parton production of HYDJET++ event is based on PYQUEN partonic loss model [28][29][30]. In brief the hard part of hadron production in HY-DJET++ uses PYQUEN [28] which includes generation of initial parton spectra according to PYTHIA and production vertices is measured at a given impact parameter. After that rescattering of partons is incorporated using an algorithm of the parton path in a dense medium along with their radiative and collisional energy loss. Finally hadronization takes place according to the Lund string model [31] for hard partons and in-medium emitted gluons. An important cold nuclear matter effect which is shadowing of parton's distribution function is included using Glauber-Gribov theory [32]. As a simplification to the model, the collisional energy loss due to scattering [33,34] with low momentum transfer is not con- The main modification which we have done in the present version of HYDJET++ is, to change the nuclear density profile function. However, this modification is not straightforward in HYDJET++ as done in AMPT by other authors [36] since HYDJET++ deals in cylindrical polar coordinates (ρ, z, ψ) instead of spherical polar coordinate system (r, θ, φ).
To make HYDJET++ work for U + U collisions, one has to transform the deformed Woods-Saxon nuclear density profile function from spherical polar to cylindrical polar coordinate system. In spherical polar coordinates the deformed Woods-Saxon for uranium nucleus is defined as follows [37]: where ρ 0 is calculated using a simple equation i.e., ρ 0 = ρ const 0 + correction and The correction term is calculated by using ρ const  [16,38]. Here the body-body and tip-tip configuration is mainly controlled by θ and all other coordinates integrated over same range.
However as shown in Ref. [36] one can change the range of φ to make other configurations as well but here we will stick to body-body and tip-tip configurations. In the conversion of nuclear density profile from spherical polar to cylindrical polar coordinate, we find a relation θ = tan −1 (r/z) and θ = tan −1 (z/r) for tip-tip and body-body configuration of U + U collision, respectively. Here r is basically ρ of cylindrical polar coordinate system and not spherical polar coordinate r. We follow this representation so that readers do not get confused it with nuclear density function (ρ). The values and range of ψ remains equal to φ during this coordinate conversion as far these two configurations are concerned. It is quite difficult to make conversion mapping between these two coordinate systems to incorporate random values of theta from its whole range i.e. 0 to π. Thus we reserve this topic for our future research work. To show the validity of our modification in deformed Woods-Saxon function and make the readers visualize, the nuclear density profiles (in cylindrical coordinate V is the jet production vertex and its coordinate will be (r cosψ, r sinψ) for tip-tip and body-body configuration both.
system) for non-deformed gold nucleus along with tip and body configuration of uranium nucleus are shown in Fig. 1, 2 and 3, respectively. Now the two quantities, nuclear thickness function (T A ) and nuclear overlap function (T AA ) can be calculated using this modified and deformed Woods-Saxon nuclear density profile function in cylindrical coordinates ρ(r, z, ψ) by following expressions [39](Please see Fig. 4 (a) and (b)): where r 1,2 (b, r, ψ) are the distances between the centers of colliding nuclei and the jet production vertex V (r cos ψ, r sin ψ), r is the distance from the nuclear collision axis to V , is the transverse distance from the nuclear collision axis to the effective boundary of nuclear overlapping area in the given azimuthal direction ψ.

B. Soft 'thermal' hadron production
The soft part of HYDJET++ is the thermal hadronic state generated on the chemical and thermal freeze-out hypersurface obtained from the parameterization of relativistic hydrodynamics with a given freeze-out condition [40,41]. The first and foremost modification which we have done in soft part is to change the nuclear density profile function for deformed uranium nucleus as discussed in above subsection. After that we have to modify the freeze-out hypersurface to properly include the effect of nuclear deformation via change in number of participants.
There are various ways to generate the initial conditions for chemical and thermal freezeout hypersurface [42][43][44]. However we first want to start here with the hydrodynamic evolution of this freeze-out hypersurface i.e., the hydrodynamic evolution laws for QCD medium.
In HYDJET++, the QCD medium is assumed to evolve according to the Bjorken boostinvariant hydrodynamics. Therefore the cooling laws for energy density and temperature are as follows [39]: respectively. In above equations, ǫ 0 , and T 0 are the initial energy density, and temperature at initial proper time τ 0 at which the local thermal equilibrium has been established. The initial energy density at τ 0 and at impact parameter b = 0 is calculated by estimating the energy density inside the co-moving volume of longitudinal size i.e., ∆z for tip-tip and ∆r for body-body configuration. The expression of total initial transverse energy deposition in the mid-rapidity region is as follows [39]: where T AA can be calculated by using Eq. (2) and (4) for tip-tip and body-body, respectively.
. p T is the first p T moment of the inclusive differential minijet cross-section which is determined by the dynamics of the nucleon-nucleon interactions at the corresponding c.m.s. energy. The initial energy density at a given impact parameter can be calculated from the following expression [39]: where S AA (b) is effective transverse area of the nuclear overlapping zone at impact parameter b [39] and is calculated as: Now to calculate the initial temperature in our calculations we have used a parameterization based on ideal thermal gas approximation [43,45] where T 0 (b = 0, τ 0 ) and baryon chemical potential µ 0 (b = 0, τ 0 ) can be calculated from the collision energy using the following relations: Here the parameters a, b, c, d, and e have been determined from the best fit of the par- using the following relation so that one can convert the fixed freeze-out hypersurface into a centrality(or N part ) dependent hypersurface which is much needed modification in soft particle production in HYDJET++: We have treated the µ B as centrality independent since the value of baryon chemical potential is small at highest RHIC energies and thus the effect of change due to centrality dependence should not affect the multiplicity by more than 5% [46]. Further hadron multiplicities are calculated using the effective thermal volume approximation and Poisson multiplicity distribution around its mean value, which is supposed to be proportional to the number of participating nucleons at a given impact parameter of A-A collision. We have shown the change in effective thermal volume between body-body and tip-tip configuration with respect to b/R A in Fig. 6. We have also plotted the variation of chemical freeze-out temperature with respect to b/R A on the same plot (Fig. 6) for body-body and tip-tip configuration of U + U collision at √ s N N = 193 GeV. Feed-down corrections from two-and three-body decays of the resonances with branching ratios are taken according to SHARE particle decay table [47] when calculating the final multiplicity of the particles.

C. Elliptic flow
Non-central collisions generate an initial spatial asymmetry of almond shape in the plane transverse to the reaction plane. The re-interactions among the reaction products in the initial state converts this spatial anisotropy into particle momentum anisotropy. In other words the spatial anisotropy in the collision zone results in anisotropic pressure gradients that generate stronger (weaker) collective flow in the direction of the major (minor) axis of the almond-shaped reaction zone. This phenomenon is called elliptic flow and is measured by v 2 . The elliptic flow coefficient v 2 is determined as the second-order Fourier coefficient in the hadron distribution over the azimuthal angle ψ relative to the reaction plane ψ R [27], Here, In HYDJET++ framework, the reaction plane of order two is zero for all the events. The above Eq. (12) can be rewritten in a simpler form as follows [48] v As we know that most of the elliptic flow arises due to the contribution of soft hadrons having lower transverse momentum and the role of hadrons having large transverse momentum is rather subdued. In HYDJET++ model, soft particle emission takes from a freeze-out hypersurface at the time of freezeout. Consequently, the elliptic flow arises in HYDJET++ model is not directly related to the initial spatial anisotropy (ǫ 0 ) of the participating nucleons as it is in other models like AMPT etc. In HYDJET++, we create a fireball which have geometrical irregularities in different directions of phase space at the time of freezeout and we assume that these irregularities are somewhat related with the initial spatial distribution of the participating nucleons in the collision region but in an involved manner. The shape of the fireball in the transverse region x − y at the freezeout can be approximated by an ellipse in non-central collision. Radii R x and R y of the ellipse at a given impact parameter b can be parameterized [49][50][51][52] in terms of spatial anisotropy at freezeout The transverse radius R ell (b, φ) of the fireball in the given azimuthal direction φ is related to spatial anisotropy at the time of freezeout as: where R 0 denotes the freeze-out transverse radius in central collision.

A. Pseudorapidity distributions
We have generated one million events for each centrality class for each of the configuration (tip-tip and body-body) separately using HYDJET++. Probability distribution curves for body-body and tip-tip events are shown in Fig. 7. We start our analysis with pseudorapidity distribution of charged hadrons. Pseudorapidity distribution of charged hadrons is a useful collisions in most central events. We have also plotted the corresponding experimental data [53,54] for comparison. observable which can help us to understand various properties of the fireball formed and the particle production process.
In Fig. 8, we have plotted the pseudorapidity distributions of charged hadrons produced  GeV for most central events [55]. One can see that dn ch /dη at midrapidity in most central U + U collisions is larger than the most central Au + Au collisions. Moreover it can be observed from the plot that the particle multiplicity in 5 − 10% tip-tip configuration of U + U collision is also larger than most central Au + Au collision. The shape of distribution at larger rapidities is somewhat different in U + U collision than Au + Au collisions. However as we already mentioned in the model formulation section that HYDJET++ uses Bjorken boost invariant hydrodynamics which is not very much applicable at larger rapidities. Thus the observations at large rapidities may change if a proper hydrodynamical treatment is incorporated in HYDJET++ at large rapidities. In Fig. 10, we have presented the variation of pseudorapidty distribution with η in body-body collisions between uranium nuclei. We have again presented the experimental multiplicity in most central Au + Au collisions at √ s N N = 200 GeV on this plot. Here again we found that the dn ch /dη at midrapidity in most central U + U collisions is greater than dn ch /dη of most central Au + Au collisions.
In Fig. 11, dn ch /dη with respect to η is shown for most central tip-tip collision. Further we have presented the jet (hard) part and hydro (soft) part separately to show their relative contribution in the total multiplicity. From Fig. 11, one can see that the hard part has relatively low contribution than the soft part and hydro part is almost 3 times larger than the jet part. One can also see that the jet part is almost flat in central rapidity region and the dip at η = 0 is mainly due to soft part of particle production. Further we have compared these most central tip-tip results with the most central body-body results. One can see that the combined multiplicity (soft plus hard) is larger in most central tip-tip configuration than the most central body-body configuration. Jet part has also the same behaviour. However, soft part shows an opposite behaviour. Here the body-body soft multiplicity is larger than tip-tip results. Similarly Fig. 12 presents the variation of dn ch /dη with respect to η for most peripheral tip-tip configuration along with separate soft and jet part. Further we have compared these results with most peripheral body-body configuration. Here we found that the combined multiplicity is larger in body-body configuration than corresponding tiptip result. Furthermore both jet as well as hydro part is larger in comparison to tip-tip configuration. Even the hydro part in body-body configuration is larger than the overall multiplicity in tip-tip configuration in most peripheral events.
In Fig. 13 Thus the difference in multiplicity between body-body and tip-tip is smaller in HYDJET++ as compared to AMPT. Another difference between HYDJET++ and AMPT is the sharp decrease in dn ch /dη by increasing η in AMPT as compared to HYDJET++ results. In Fig. 14, we have calculated the dn ch /dη at midrapidity in minimum bias configuration using HYDJET++. We have used a pseudorapidity cut as |η| < 0.5. Further we have compared HYDJET++ results with the experimental results obtained by PHENIX collaboration [53].
We found that the minimum bias data is successfully reproduced by HYDJET++ in the case of dn ch /dη at midrapidity.  5 − 10% centrality class at intermediate p T range. However, at low and high p T range, the multiplicity is larger in Au+Au collision than the body-body configuration of U +U collision.
In body-body configuration (see Fig. 18), the slope of distribution is more in comparison to tip-tip configuration for given centrality class due to the effect of transverse flow. The difference of p T distribution for both the configurations can be seen at intermediate and large p T region for central collision. As shown in Fig. 19, for most-peripheral collisions there is small difference between tip-tip and body-body configurations in low p T region only.
As we know that most of the low-p T particles are due to thermal production and high-p T  particles are due to jet fragmentation. Thus, in peripheral collision the initial configuration of nuclei affects the thermal part mostly (as shown in Fig. 12) and very small difference in jet-part but in central collisions, initial configuration mostly affect the jet-fragmentation part at higher p T (as shown in Fig. 11).
In with the experimental data. We have observed a suitable match between data and the model results. Fig. 21 demonstrates the variation of elliptic flow with respect to transverse momentum (p T ) for various centrality class in body-body configuration of U + U collisions. We have shown a comparison of v 2 for tip-tip and body-body in central collisions in Fig.   23. We find that the elliptic flow of body-body configuration is slightly larger than the elliptic flow in tip-tip configuration and as we move towards larger p T , this difference in v 2 between two configurations increases with the increase in p T . We have also plotted the STAR experimental data of U +U collisions in 0−5% and 0−0.5% centrality class with |η| ≤ 1 [12]. The thought behind calculating elliptic flow for 0 − 0.5% centrality class in STAR is that they should consists mainly tip-tip events of U + U collisions. STAR collaboration has done the calculation of v 2 for charged hadrons. By comparison we observe that STAR data of 0 − 0.5% centrality class have lower v 2 in comparison to our tip-tip as well as body-body results. However 0 − 0.5% centrality class data matches with our tip-tip results when p T < 1 GeV (see inset of Fig. 23). Further the experimental data from 0−5% is greater than our tiptip result when p T < 1 GeV but it matches with our tip-tip results for intermediate and large p T . Another important observation is that the v 2 in body-body configuration is higher than both these data sets along with tip-tip results from HYDJET++. In peripheral collisions (see Fig. 24), the qualitative difference between v 2 in tip-tip and body-body configuration is same. However the magnitude of difference in v 2 of charged hadrons is quite visible even at body and tip-tip configuration. Green stars are experimental data taken from Ref. [12] and AMPT model data is taken from Ref. [36].
from central to peripheral which is actually due to a increase in eccentricity going from central to peripheral collisions. However from here it is clear that in central collisions the difference in magnitude between body-body and tip-tip collisions is small. However, in semiperipheral as well as in peripheral collisions, one can distinguish between body-body and tip-tip events by observing the v 2 magnitude of charged hadrons. We have also shown the results obtained in Ref. [36] using AMPT model in two different modes (string melting mode and default mode). We found that the qualitative behaviour of variation of v 2 with centrality in HYDJET++ is quite opposite to AMPT model and shows a small difference in v 2 for tip-tip and body-body configuration in central events and a large difference in peripheral events. On the other side AMPT has shown opposite behaviour. We have also plotted the STAR experimental data [12] of v 2 as a function of centrality for minimum bias events. We found that the experiment results are nearly in between the HYDJET++ model results for tip-tip and body-body configurations. However, the experimental data is between the AMPT-Default mode results for tip-tip and body-body configurations in central and mid-central events but not in peripheral events. In AMPT-SM mode, the experimental data is very close to tip-tip configuration results.
In summary, we have calculated and shown the pseudorapidity density and transverse momentum distributions of charged hadrons produced in U + U collisions at √ s N N = 193 GeV in various initial geometrical configurations. In present study, it has been shown that the correlation between multiplicity and initial geometrical configurations of U +U collisions is small which is in accordance with the recent experimental observation. However, the experimental results are quite preliminary due to complexity in disentangling the tip-tip and body-body events. We have shown the midrapidity charged-particle multiplicity distribution from HYDJET++ model which is in good agreement with the experimental results for minimum bias events. Further, we have shown the evolution of elliptic flow with p T and centrality in different configurations of U + U collisions. It has been observed that elliptic flow generated in body-body collisions is larger than tip-tip collisions but the difference in magnitude of v 2 is small in central collisions and large in peripheral events. Further, we found that our tip-tip results of elliptic flow matches with STAR experiment result of 0-0.5% centrality class when p T < 1. At last, we have observed that the experimental results of v 2 as a function of centrality for minimum bias events are nearly in between the tip-tip and body-body configuration results of our model. However, this is not the case for AMPT results. Finally we may conclude that our present study will shed some light on the particle production mechanism and the evolution of the fireball created in various geometrical configurations of U + U collisions specially the entanglement of hard (jet) and soft (hydro) part in body-body and tip-tip configurations.