Hydrodynamic Collectivity in Proton--Proton Collisions at 13 TeV

In this paper, we investigate the hydrodynamic collectivity in proton--proton (pp) collisions at 13 TeV, using iEBE-VISHNU hybrid model with HIJING initial conditions. With properly tuned parameters, our model simulations can remarkably describe all the measured 2-particle correlations, including integrated and differential elliptic flow coefficients for all charged and identified hadrons ($K_S^0$, $\Lambda$). However, our model calculations show positive 4-particle cumulant $c_{2}\{4\}$ in high multiplicity pp collisions, and can not reproduce the negative $c_{2}\{4\}$ measured in experiment. Further investigations on the HIJING initial conditions show that the fluctuations of the second order anisotropy coefficient $\varepsilon_{2}$ increases with the increase of its mean value, which leads to a similar trend of the flow fluctuations. For a simultaneous description of the 2- and 4- particle cumulants within the hydrodynamic framework, it is required to have significant improvements on initial condition for pp collisions, which is still lacking of knowledge at the moment.

The high energy proton-lead (p-Pb) and protonproton (p-p) collisions at the LHC were originally aimed to provide the reference data for the high energy nucleusnucleus collisions. However, various unexpected phenomena have been observed in these small systems, especially in the high multiplicity region. One surprising discovery is the long-range "ridge" structures in two-particle azimuthal correlations with a large pseudo-rapidity separation in high multiplicity p-Pb and p-p collisions [33][34][35][36][37]. Such long-range correlation structures were firstly discovered in Au-Au and Pb-Pb collisions and interpreted as a signature of the collective expansion. In general, the theoretical interpolations for the long-range "ridge" structure in the small systems can be classified into three big categories: final state interactions, such as hydrodynamic expansion [38][39][40][41][42][43], parton cascade [44][45][46][47], hadronic rescattering [48], rope and shoving mechanism [49], initial state effects related to the gluon saturation [50][51][52][53][54][55][56][57] and combinations of both initial and final state effects [58]. For recent theoretical progresses, please refer to [7,59].
In experiments, one of the crucial questions on the "ridge" structure is whether it arises from correlations of all particles related to the collective flow or it only involves with the correlations from few particles, e.g. from resonance decays or jets, which is defined as non-flow. In small pp and p-Pb systems, the non-flow contributions are always significant, even in case that collective expansion has been developed. It is thus necessary to remove such non-flow effects before comparing the data with the model calculations. Based on different assumptions, various non-flow subtraction methods, e.g. template fit [36,60,61] and peripheral subtraction [37] have been applied to the measurements of 2-particle correlations in pp collisions at 13 TeV, which yield different non-flow subtracted results. Currently, it is still unclear which one is a better approach to remove the non-flow effects.
Compared with the 2-particle correlations, multiparticle cumulants are less influenced by the non-flow effects, which are expected as one of the key observables to evaluate the anisotropic collectivity of the small systems. Besides, the multi-particle and 2-particle cumulants show different sensitivities to event-by-event flow arXiv:1801.00271v2 [nucl-th] 21 Apr 2018 fluctuations [1,62,63]. An extensive measurement of these different cumulants could provide tight constraints on the initial state fluctuations.
In order to extract real values of the flow coefficients, the 2-, 4-, 6-and 8-particle cumulants are expected to carry positive, negative, positive and negative signs, respectively. Such "changing sign pattern" has been observed in the measured 2-and multi-particle cumulants in Pb-Pb collisions at the LHC, where the created QGP fireballs undergo fast collective expansion [8,10]. Based on the similar idea, it is proposed to measure the 2and multi-particle cumulants to evaluate the collectivity in high multiplicity pp collisions at 13 TeV. However, it was found that the standard multi-particle cumulants in the small systems are still largely affected by the residual non-flow and its fluctuations, which even presents fake flow signals with the "right sign". Recently, Ref. [64,65] developed 2-and 3-subevent methods for the multi-particle cumulants, which further remove the residual non-flow from jet. The related measurements from ATLAS have confirmed the observations of positive 2-particle cumulants and negative 4-particle cumulants in high multiplicity pp collisions at 13 TeV [66]. It is thus the proper time to study and evaluate these possible collective flow signal, using the hydrodynamic calculations. In this paper, we will implement iEBE-VISHNU hybrid model with HIJING initial conditions to study 2-and 4particle correlations in pp collisions at √ s NN = 13 TeV, together with a detailed examinations of the initial state fluctuations from HIJING.

II. THE MODEL AND SET-UPS
iEBE-VISHNU [67] is an event-by-event simulation version of the early developed hybrid model VISHNU [68] that combines (2+1)-d viscous hydrodynamics VISH2+1 [69,70] to describe the QGP expansion with a hadron cascades model UrQMD [71,72] to simulate the evolution of hadronic matter. In the viscous hydrodynamics part, VISH2+1 solves the transport equations for the energymomentum tensor T µν and shear stress tensor π µν with a state-of-art equation of state (EoS) s95-PCE [29,73] as an input. For simplicity, we neglect the bulk viscosity, net baryon density and heat conductivity and assume that specific shear viscosity η/s is a constant. The hydrodynamic evolution matches the hadron cascade simulations at a switching temperature T sw , where various hadrons are emitted from the switching hyper-surface for the succeeding UrQMD evolution.
In our calculations, we use the modified HIJING model [74] to generate fluctuating initial profiles for the succeeding iEBE-VISHNU simulations in high energy pp collisions. In HIJING [75][76][77], the produced jet pairs and excited nucleus are treated as independent strings, where the hard jet productions are calculated by pQCD, and the soft interactions are treated as gluon exchange within Lund string model. Here, we assume these strings break , where σ R is the Gaussian smearing factor of the string.
Following [79], the initial energy density profiles in the transverse plane are constructed from the energy decompositions of emitted partons of HIJING together with an additional Gaussian smearing [79] ( where σ 0 is the Gaussian smearing factor, p i is the momentum of the produced parton i and K is an additional normalization factor, U 0 is the initial flow velocity of the corresponding fluid cell. Here, we assume zero transverse initial flow and only consider the partons within the spacial mid-rapidity |η s | < 1 (for related details, please refer to [79]).
In iEBE-VISHNU simulations with HIJING initial conditions, the hydrodynamic starting time τ 0 , switching temperature T sw , specific shear viscosity η/s, the Gaussian smearing width σ R and σ 0 , and normalization factor K are free parameters, which need to be fixed by experimental data. In general, one uses total multiplicity, p T spectra of identified hadrons and integrated flow harmonics of all charged hadrons to tune the related parameters in the hydrodynamic calculations [31,79,80]. However, not all these needed data are available in pp collisions at 13 TeV. For example, the p T spectra has not been measured and released. Here, we assume that slope of p T spectra do not significantly change from 7 TeV to 13 TeV, and use these slopes of pions and protons at 7 TeV [81], the total multiplicity [37] and v 2 {2} [37, 61] at 13 TeV to partially constrain the parameters in HIJING and iEBE-VISHNU. In fact, these available data are not enough to fully fix all these parameters in our model calculations, especially considering that the measured v 2 {2} from CMS [37] and ATLAS [60,61] differ by 20%. We thus select four possible parameter sets as listed in Tab. I, that roughly fit slope of p T spectra at 7 TeV and fit the measured v 2 {2} from either ATLAS or CMS collaborations, for the following calculations and investigations in Sec. III. The predicted p T spectra of pions, kaons and protons in 0-0.1% pp collisions at 13 TeV are shown in Fig. 1, where the centralities are cut by the multiplicity of all charged hadrons at mid-rapidity |η| < 2.4 [86]. The spectra at 13 TeV have similar slopes as the ones at 7 TeV, which can be measured in experiments in the near future.

III. RESULTS AND DISCUSSION
The 2-and multi-particle correlations are common measurements to study anisotropic azimuthal correlations, which can be calculated with the Q-cumulant method [82] and the Generic Framework [48]. In our calculations, these two methods are identical, since it does not involve any inefficiency in azimuthal acceptance or tracking efficiency as happened in experiments. We thus follow the same procedure as in our early study [48] to calculate the 2-and multi-particle correlations as well as the related flow harmonics.

A. 2-particle cumulant
With the four sets of parameters listed in Table I, we calculate the integrated flow harmonics v 2 , v 3 and v 4 in pp collisions at 13 TeV, using iEBE-VISHNU model with HIJING initial conditions. The 2-particle correlation method with a pseudorapidity gap |∆η| > 0, kinematic cuts 0.3 < p T < 3.0 GeV/c and |η| < 2.4 is applied in our calculations [87]. The results and the comparisons to the experimental data are shown in Fig. 2 TeV, calculated by iEBE-VISHNU with HIJING initial conditions. The CMS and ATLAS data are taken from [37] and [61], respectively. eral subtraction" method [37], the solid crosses and stars are the ATLAS v 2 {2} measurements with "template fit" and "peripheral subtraction" methods [60]. As shown in Fig. 2, our model calculations reproduce the multiplicity dependence of the integrated v 2 {2} from low multiplicity (∼ 30) to high multiplicity (∼ 160). More specifically, our calculations with parameter sets I, II and III fit the CMS and ATLAS measurements with the "peripheral subtraction" method, and the one with para-IV describes the ATLAS data from the "template fit" method. For low multiplicity N ch ∼ 30, the iEBE-VISHNU calculations fail to describe the CMS or ATLAS measurements with the "peripheral subtraction" method, which do not significantly decrease as the data. Figure 2 also compares v 3 {2} from iEBE-VISHNU and from CMS measured with the "peripheral subtraction" method [37]. Our calculations with III roughly reproduce the CMS data, while the results from para-I, para-II and IV obviously over-estimate the data. It thus shows that the v 3 {2} measurements could further constrain the initial conditions. We also predict v 4 {2} as a function of multiplicity. Currently, the related experimental data is not publicly available, but can be compared with our model calculations in the near future [83].
In Fig. 3, we calculate the p T -differential v 2 (p T ) of all charged and identified hadrons, using iEBE-VISHNU simulations and with the 2-particle cumulant method with a pseudorapidity gap |∆η| > 0. Panel (a) shows a comparison between model and data for all charged hadrons, where the data from CMS and ATLAS are measured by the "peripheral subtraction" method [37] and "template fit method" [60], respectively. For para-I, para-II and III, our calculations roughly describe the CMS and ATLAS measurements within p T < 2.0 GeV/c. In contrast, the calculations with para-IV slightly over predict the data , and for pions, kaons and protons (c) in high multiplicity pp collisions at 13 TeV, calculated by iEBE-VISHNU with HIJING initial condition. The CMS and ATLAS data are taken from [37] and [60], respectively. above 1.0 GeV/c. Fig. 3 (b) and (c) show the v 2 (p T ) for K 0 S and Λ and for π, K and p for the multiplicity range 80 < N sel ch < 120. Clear mass orderings between K 0 S and Λ and among π, K and p are seen in our iEBE-VISHNU calculations. In the hydrodynamic language, the radial flow blue-shifts the lower-p T to higher p T with the mass-dependent effects, which leads to the observed mass ordering among various hadron species. Panel (b) also shows that the v 2 mass splitting between K 0 S and Λ is more significant for the calculations with para-III, which indicates a stronger radial flow. This consists with the results of p T -spectra in Fig. 1, which shows that p T -spectra from para-III are flatter than other ones due to the larger radial flow.

B. 4-particle cumulant
In experiments, the observed negative c 2 {4}, together with the positive c 2 {6} and negative c 2 {8}, is interpreted as a signature of collective expansion in the small systems. However, although iEBE-VISHNU with HIJING initial conditions can describe the measured 2-particle correlations for both charged and identified hadrons, it fails to reproduce the negative c 2 {4} as measured by CMS and ATLAS with the standard cumulant method [37] and three-subevent method [61] (0.3 < p T < 3.0 GeV). Fig. 4 shows that, for four parameter sets of HIJING initial conditions, iEBE-VISHNU always predicts positive values of c 2 {4} in the high multiplicity regime. We have also checked that these positive values are not caused by the specific cumulant method, possible non-flow contributions or multiplicity fluctuations in our model calculations (please also refer to the Appendix for details).
In pure hydrodynamics, 2 v 2 2 2 , which is influenced by both flow fluctuation and the mean value, and can be evaluated by the related v 2 distribution P (v 2 ). Due to the approximate linear relationship between v 2 and ε 2 , P (v 2 ) almost follows P (ε 2 ) of the initial condition model [28,84,85]. We thus further check P (ε 2 ) distributions of HIJING in Fig. 5, which shows that the fluctuations σ ε2 increases with the mean values ε 2 . In other words, the narrower P (ε 2 ) distribution with smaller σ ε2 has a smaller mean value ε 2 , and vice versa. In  Fig. 4. For para-III, c 2 {4} ε has a small negative value. However, the pp fireballs in our hydrodynamic simulations do not evolve enough long time to translate that negative c 2 {4} ε into an definite negative c 2 {4} as measured in experiment. For a simultaneous description of the 2-and 4-particle cumulants within the framework of hydrodynamics, other initial condition models for pp collisions should be further developed and investigated.

IV. SUMMARY
In this paper, we studied the 2-and 4-particle cumulants in proton-proton (pp) collisions at 13 TeV, using iEBE-VISHNU hybrid model with HIJING initial condition. With properly tuned parameters, our model calculations quantitatively describe the measured 2-particle cumulants, including the integrated second and third order flow coefficient v n {2} (n = 2, 3) for all charged hadrons. In addition, iEBE-VISHNU also reproduces the p T -differential elliptic flow v 2 (p T ) for all charged and identified hadrons (K 0 S and Λ) in the high multiplicity pp collisions. We also predicted the v 2 (p T ) of pions, kaons and protons, which shows similar characteristic mass ordering feature as the cases in Pb-Pb and p-Pb collisions, and can be further examined in experiments.
However, our iEBE-VISHNU calculations with HIJING initial conditions always give positive values of 4-particle cumulants c 2 {4} for various parameter sets, which can not reproduce the negative c 2 {4} measured by CMS and ATLAS. Further investigations showed that this positive c 2 {4} are not caused by possible non-flow contributions, multiplicity fluctuations or the multi-particle cumulant method applied in our calculations. In fact, it is originated from the imprinted fluctuation pattern of the HIJING initial conditions, where the fluctuations of the eccentricity ε 2 increase with the increase of the mean value. Due to the approximate linear response between v 2 and ε 2 , the mean value and fluctuations of the second order flow coefficient v 2 present a similar trend, which fails to simultaneously fit the measured 2-and 4-particle cumulants with various possible parameters. In order to simultaneously fit the flow-like data in high multiplicity pp collisions at 13 TeV within the framework of hydrodynamics, other initial condition model for the small systems should be further developed. Besides hydrodynamics, it is also necessary to investigate these 2-and 4-particle cumulants in high energy pp collisions within other theoretical approaches, to better understand the physics in high energy pp collisions. In Ref. [64,65], it was argued that 4-particle cumulants with 2-and 3-subevent methods could further remove the residual non-flow, e.g. from the contributions of jets. Considering that non-flow effects in iEBE-VISHNU simulations are not influenced by jet, but mainly contributed by resonance decays, we implement the standard method to calculate the 4-particle cumulant c 2 {4} in Sec. III. In Fig. 6 (a), we further compare the c 2 {4} from the standard method and from the 2-subevent and 3-subevent methods, using iEBE-VISHNU simulations with HIJING initial conditions (para-I). It shows a good agreement for these three methods, which indicates that the nonflow in our calculations have been cleanly removed by the standard 4-particle cumulant, which are not necessary to further implement the 2-and 3-subevent methods that require larger statistical runs. and using the standard 4-particle cumulant method for cases (b) and (c). As shown in Fig. 6 (b), these two results from cases (a) and (b) almost overlap within error bars, which indicates that the 4-particle cumulant c 2 {4} can properly describe the flow and flow fluctuations in a small collision system with N ch ∼ 100. In Fig. 6 (b), we also study the effects of hadronic evolution on c 2 {4} through comparing case (b) the pure hydrodynamic calculations and case (c) the full iEBE-VISHNU simulations, which shows that the hadronic scatterings and decays slightly decrease c 2 {4}.