We have updated the parton and hadron cascade model PACIAE 2.0 (cf. Ben-Hao Sa, Dai-Mei Zhou, Yu-Liang Yan, Xiao-Mei Li, Sheng-Qin Feng, Bao-Guo Dong, Xu Cai, Comput. Phys. Comm. 183 (2012) 333.) to the new issue of PACIAE 2.1. The PACIAE model is based on PYTHIA. In the PYTHIA model, once the hadron transverse momentum pT is randomly sampled in the string fragmentation, the px and py components are originally put on the circle with radius pT randomly. Now it is put on the circumference of ellipse with half major and minor axes of pT(1+δp) and pT(1−δp), respectively, in order to better investigate the final state transverse momentum anisotropy.Manuscript title: PACIAE 2.1: An updated issue of the parton and hadron cascade model PACIAE 2.0Authors: Ben-Hao Sa, Dai-Mei Zhou, Yu-Liang Yan, Bao-Guo Dong, and Xu CaiProgram title: PACIAE version 2.1Journal reference:Catalogue identifier:Licensing provisions: noneProgramming language: FORTRAN 77 or GFORTRANComputer: DELL Studio XPS and others with a FORTRAN 77 or GFORTRAN compilerOperating system: Linux or Windows with FORTRAN 77 or GFORTRAN compilerRAM: ≈ 1GBNumber of processors used:Supplementary material:Keywords: relativistic nuclear collision; PYTHIA model; PACIAE modelClassification: 11.1, 17.8External routines/libraries:Subprograms used:Catalogue identifier of previous version: aeki_v1_0*Journal reference of previous version: Comput. Phys. Comm. 183(2012)333.Does the new version supersede the previous version?: Yes*Nature of problem: PACIAE is based on PYTHIA. In the PYTHIA model, once the hadron transverse momentum (pT) is randomly sampled in the string fragmentation, the px and py components are randomly placed on the circle with radius of pT. This strongly cancels the final state transverse momentum asymmetry developed dynamically.Solution method: The px and py component of hadron in the string fragmentation is now randomly placed on the circumference of an ellipse with half major axis of pT(1+δp) and the half minor axis of pT(1−δp) instead of the circle.Reasons for the new version: PACIAE is based on PYTHIA, where once the hadron transverse momentum (pT) is randomly sampled in the string fragmentation, the px and py components are randomly placed on the circle with radius of pT. This is not only strongly canceling the final state transverse momentum asymmetry developed dynamically, but also inconsistent with the ATLAS observation of the final state charged particle transverse sphericity being less than unity [8].Summary of revisions: The main revision is executed by randomly placing px and py components of the hadron transverse momentum pT, in the string fragmentation, on the circumference of an ellipse with half major axis of pT(1+δp) and half minor axis of pT(1−δp) instead of a circle.Restrictions: Depend on the problem studied.Unusual features:Additional comments: Email addresses: zhoudm@phy.ccnu.edu.cn (D.-M. Zhou), yanyl@ciae.ac.cn (Y.-L. Yan).Running time:•Using the attached input file of usux.dat (where the string fragmentation is selected and the elastic parton–parton interactions are considered only, the same later) to run 1000 events for the s=200GeV Non Single Diffractive pp collision by 21a.tar.gz takes 0.5 min.•Using the attached input file of usu.dat to run 10 events for the 10%–40% most central Au+Au collisions at sNN=200GeV by 21b.tar.gz takes 5 min.•Using the attached input file of usu.dat to run 10 events for the 10%–40% most central Au+Au collisions at sNN=200GeV by 21c.tar.gz takes 17 min.1.The large azimuthal anisotropy (the large second harmonic coefficient v2) of the emitted particle is an important feature of the hot and dense medium created in the ultra-relativistic nuclear collisions. This large v2 has contributed to the observation of a strongly coupled quark–gluon plasma (sQGP) in the nucleus–nucleus collisions at the RHIC energies [1–4].The nuclear overlap zone created in a nucleus–nucleus collisions at a given impact parameter possesses an almond-like spatial asymmetry. Because of the strong parton rescattering, the local thermal equilibrium and asymmetric pressure gradient may build up in this initial fireball. The asymmetric pressure gradient then drives a collective anisotropic expansion. The expansion along the almond minor axis (along the large pressure gradient) is faster than the one along the major axis. This results in a strong asymmetric transverse momentum azimuthal distribution and hence a large elliptic flow coefficient v2 of the final hadronic state.As mentioned in [6], PACIAE is a parton and hadron cascade model for the ultra-relativistic nuclear collisions and is based on PYTHIA [7]. In the PACIAE model, a nucleus–nucleus collision is decomposed into a sequence of nucleon–nucleon (NN) collisions according to the collision geometry and the NN total cross section. Each NN collision is performed, in turn, by the PYTHIA model with the string fragmentation switched-off temporarily and the diquark (anti-diquark) broken into quark pairs (anti-quark pairs) randomly. The parton rescattering then proceeds. This parton evolution stage is followed by the hadronization at the moment of partonic freeze-out (exhausting the partonic collisions). The Lund string fragmentation regime and/or phenomenological coalescence model is provided for the hadronization. Then the rescattering among produced hadrons is dealt with by the usual two body collision model [6].In the PYTHIA model [7] once the transverse momentum pT of a final state hadron generated from the string fragmentation and/or the unstable particle decay is randomly sampled, the px and py components are randomly placed on the circle with radius of pT. This px and py determination may strongly cancel the final state transverse momentum anisotropy developed dynamically. The charged particle transverse sphericity [8–10] may reach to unity (isotropic). This is inconsistent with the experimental observation that the charged particle transverse sphericity is less than unity [8]. Therefore we randomly placed the generated final state hadrons on the circumference of an ellipse with half major axis of pT(1+δp) and a half minor axis of pT(1−δp) instead of a circle. This change is also introduced in the particle/parton production process of hard scattering, multiple interactions, initial- and final-state parton showers, for the nucleus–nucleus collisions [7]. This change is even introduced in the deexcitation of energetic quark (anti-quark) when the phenomenological coalescence model [6] is selected for hadronization. Of course, a new pT should be recalculated by px and py after this change. Then the transverse momentum distribution of the final state hadron may be modified. However, if the deformation parameter δp is less than unity (a small perturbation) the change in transverse momentum distribution may be weak.From ideal hydrodynamic calculations [11] one knows that the integrated elliptic flow parameter is directly proportional to the initial spatial eccentricity of the nuclear overlap zone. Therefore, if the nuclear overlap zone is assumed to be an ellipse with major axis of b=2RA(1+δr) and minor axis of a=2RA(1−δr) for a symmetry nucleus–nucleus collision with nuclear radius of RA, we may assume (1)δp=Cδr where C is an extra model parameter. C=0 corresponds to the original case of px and py put on the circle randomly.In order to calculate δr we first calculate the reaction plane eccentricity [12](2)ϵrp=σy2−σx2σy2+σx2 according to the participant nucleons spatial distributions inside the nuclear overlap zone [6] in the PACIAE simulation. In the above equation, σx2=x2¯−x̄2 (the same for σy2) and x̄2 (x̄) denotes an average of x2 (x) over particles in a single event. The event average reaction plane eccentricity reads (3)〈ϵrp〉=〈σy2−σx2σy2+σx2〉. On the other hand, the geometric eccentricity [13] of the ellipse-like nuclear overlap zone is (4)ϵg=b2−a2b2. Letting ϵg=ϵrp, one approximately obtains (5)δr≃ϵrp24.The calculated charged particle v2(pT) at mid-rapidity (|η|<1) in the 10%–40% most central Au+Au collision at sNN=200GeV is compared with the corresponding STAR data [5] in Fig. 1. In this figure the STAR data are denoted by solid symbols: the black circles are measured with the event plane method (EP), red squares with Lee–Yang zero point method (L–YZ), and green triangles with four particle cumulant method (4 cumulant). The PACIAE results are given by open symbols: the black circles calculated with C=4, red squares with C=3, and green triangles with C=2. One sees in this figure that the STAR data [5] on the charged particle v2(pT) are able to be reproduced by the PACIAE calculations with C=3. The C=0 PACIAE results are too small compared to the STAR data. Fig. 1(Color online) The charged particle v2(pT) at mid-rapidity (|η|<1) in 10%–40% most central Au+Au collisions at sNN=200GeV. The solid symbols are the STAR data taken from [5] and open symbols are the PACIAE results calculated with C=2, 3, and 4.A similar comparison for charged particle v2(η) (0.15<pT<2GeV/c) in the 10%–40% most central Au+Au collision at sNN=200GeV is given in Fig. 2. Here one sees again that the STAR data on the charged particle v2(η) [5] are able to be reproduced by the PACIAE calculations with C=3. We have to mention here that in Fig. 3 of Ref. [5] the PHOBOS data [14] were introduced to compare with the STAR data and to complement the lack of the STAR data in 1.5<|η|<2.5 region. Because the PHOBOS data were measured for the 0%–40% most central Au+Au collisions at the same energy but in the full pT phase space, it is not suitable to compare the PHOBOS data with the STAR data. Therefore we do not include the PHOBOS data [14] in Fig. 2 here. Fig. 2(Color online) The charged particle v2(η) (0.15<pT<2) in 10%–40% most central Au+Au collisions at sNN=200GeV. The solid symbols are the STAR data taken from [5] and open symbols are the PACIAE results calculated with C=2, 3, and 4.We give the calculated charged particle transverse momentum distribution in 10%–40% most central Au+Au collisions at sNN=200GeV in Fig. 3. Fig. 3(a) and (b) are drawn for the full and partial (|η|<1) pseudo-rapidity phase space, respectively. In panels (a) and (b) the solid black circles, open red circles, open green triangles, and the blue line are calculated with C=4, 3, 2, and 0, respectively. We see in this figure that the charged particle transverse momentum distribution is really not sensitive to the parameter C both in the full and partial η phase space, provided the deformation parameter δp is less than unity (a small perturbation). Fig. 3(Color online) The charged particle transverse momentum distribution in 10%–40% most central Au+Au collisions at sNN=200GeV calculated by the PACIAE model with C=0, 2, 3, and 4: (a) in the full η phase space, (b) in |η|<1.In addition, an extra switching parameter of iparres is introduced in the new issue of PACIAE 2.1. We assume that the iparres=0 is for the elastic parton–parton rescattering only and iparres=1 for otherwise.