Hydrodynamic modeling of heavy-ion collisions

This contribution presents a theoretical overview of hydrodynamic modelling of heavy-ion collisions, with highlights on some recent developments. In particular, the formulation of anisotropic hydrodynamics, the role of hydrodynamic fluctuations, and the non-linear coupling of flow coefficients will be discussed.


Introduction
Heavy-ion collisions carried out at the Relativistic Heavy-Ion Collider and the Large Hadron Collider create a Quark-Gluon Plasma (QGP), a hot medium that behaves like a perfect fluid.The fluidity of QGP is essential for understanding many of the observed phenomena in high energy nucleus-nucleus collisions, in particular for soft probes such as long-range multi-particle correlations whose dynamical properties are dominated by the bulk medium evolution.
To quantify long-range multi-particle correlations, the harmonic flow V n are often introduced.Following a Fourier decomposition of the emitted single-particle spectrum in the azimuthal direction, harmonic flow characterizes momentum asymmetry in azimuth order-by-order.The curly braces in Eq. ( 1) notate average over the emitted particles of an individual collision event.Note that the V n are defined as complex quantities in Eq. ( 1), with magnitude and phase indicated by v n and Ψ n , respectively.For n = 2 for instance, the so-called elliptic flow [1] V 2 describes asymmetric emission of particles from directions in-and out-of reaction plane.In experiments, extensive evidence of QGP fluidity has been obtained involving V n in various aspects in nucleus-nucleus collisions.In addition to the simplest observables related to the rms average of flow magnitude v n {2} = |v n | 2 (double brackets . . .stand for an event average), more precised details of harmonic flow have been investigated in recent measurements.These include fluctuations and correlations of the flow magnitude v n and correlations of phase Ψ n , corresponding respectively to cumulant TeV.Figures from Ref. [12] and Ref. [13].
of harmonic flow v n {m} (m = 4, 6, etc.), event-by-event flow distributions P(v n ) [2], symmetric cumulants sc(n, m) [3] and event-plane correlators [4].Remarkably, an appropriate modeling of the medium evolution with hydrodynamics successfully describes all of the observed flow signatures (see [5,6] for recent reviews).Fig. 1 presents results of hydro simulations for the most recent measurements of the rms flow harmonics, for the Pb+Pb collisions at √ s NN = 5.02 TeV, where agreement between hydro and experiments is clearly seen.One may refer to Refs.[7,8,9,10,11] for hydro calculations of some other types of flow observables, regarding recent measurements of the nucleus-nucleus collisions.

Ingredients of hydrodynamic modeling
Relativistic hydrodynamics plays a central role in the modeling of heavy-ion collisions, in which the evolution of QGP medium is described by solving an appropriate equation of motion.Hydrodynamics is an effective theory for the low-energy (long wavelength) degrees of freedom of an evolving thermal system, consisting of a set of conservation laws.The conservation of energy-momentum, for instance, leads to hydrodynamic equation of motion where energy-momentum tensor T µν is constructed in terms of hydrodynamic variables: Flow velocity u µ , energy density e, and pressure P, ∆ µν = g µν − u µ u ν is a projection tensor, which helps to define the covariant form of spatial gradient, ∇ µ ≡ ∆ µν ∂ ν .The terms π µν and Π in Eq. ( 3) are shear and bulk viscous corrections to the stress tensor, respectively, containing a gradient expansion of hydrodynamic variables.To the first order in the expansion, the shear and bulk stress tensors are given by the Navier-Stokes hydrodynamics, with coefficients shear viscosity η and bulk viscosity ζ determined by the underlying microscopic dynamics.
For the case of QGP, it has been realized that hydrodynamic description of the medium demands a small shear viscosity over entropy ratio, η/s ∼ O(1/4π) [14], which is close to the theoretical expectation for a strongly coupled system [15].To avoid acausal modes in practical simulations, second order viscous corrections are taken into account as well, which depend on second order transport coefficients, τ π , τ Π , etc. Eq. ( 2) is coupled and closed by an equation of state, for which results from lattice QCD are generally incorporated [16].With respect to the system evolution in heavy-ion collisions, hydrodynamics cannot be applied to the very early stages after collisions, during which the system is far from local equilibrium, or very late stages when the system starts to decouple kinematically.Despite a large local pressure anisotropy at the initial stages of medium evolution (see Section 3.1 for more details), in hydro modeling an abrupt switch-on of hydro is taken into account around time τ 0 ∼ O(1) fm/c, given an appropriate initial distribution of energy (or entropy) density as input for the equation of motion, Eq. ( 2).In hydro modeling, initial density profile can be realized by a variety of effectively models, such as MC-Glauber [17], MC-KLN [18], EKRT [10] and IP-Glasma [19].In these models, the energy (or entropy) deposition from collisions is realized through different prescriptions, but initial state fluctuations of nucleon-nucleon collisions are typically implemented on an event-by-event basis.For instance, the IP-Glasma model is inspired by the idea of gluon saturation [20].Of particular importance are the so-called initial state anisotropies, which are defined through the energy density profile (n ≥ 2) [21], E n quantify geometric deformations of the initial density profile, owing to the overlapped collision geometry and fluctuations.Upon medium collective expansion, E n provide the geometric origin of momentum anisotropies of the final state particle spectrum, namely, V n .Regarding the calculations of flow harmonics, it is now widely acknowledged that a large source of uncertainty of hydro modeling comes from initial state.At late stages of the system evolution, hydrodynamics breaks down along with particle generation from the decoupled medium (Cooper-Frye freeze-out [22]), which is further followed in the hydro modeling by hadron-cascade dominated evolution and resonance decay [23,24].

Anisotropic hydrodynamics
Applicability of viscous hydrodynamics is constrained by the convergence of gradient expansion, which however is conceptually challenged by the issue of local pressure anisotropy.In the initial stages of heavyion collisions, evolution towards local thermal equilibrium of a quark-gluon system experiences strong longitudinal expansion along the beam axis, which is reflected by the ratio of longitudinal to transverse pressures, P L /P T .These pressures are defined in kinematics as (assuming z-axis as the beam axis) so that for a locally thermalized system, when the phase-space distribution function gets isotropized, P L /P T = 1.Since deviation from the thermalized distribution function in hydrodynamics is reflected in dissipation, pressure anisotropy is related to viscous corrections. 1 It is only with small viscous corrections, namely when P L /P T is close to unity, can one expect a good convergence behavior of the gradient expansion.Illustrated in Fig. 2 is a qualitative picture of the time evolution of pressure anisotropy in the initial stages of heavyion collisions, based on theoretical analyses of weakly coupled QCD systems (cf.Ref. [25]).Pressure anisotropy is potentially large around τ 0 ∼ 1 fm/c, although at which hydro starts to be applied in present modelings of heavy-ion collisions.Distinguished from a canonical formulation of viscous hydrodynamics (such as Chapman-Enskog expansion), in which the gradient expansion is realized with respect to an isotropized distribution function, in the so-called anisotropic hydrodynamics (ahydro) [26,27], pressure anisotropy is absorbed into an anisotropic background.Such a modification greatly improves the applicability of hydrodynamics, especially accounting for a large pressure anisotropy.An example is shown in Fig. 3 for the calculation of number density, for a special case of Bjorken flow (medium expansion dominated in the longitudinal direction) and relaxation time approximation.For a Bjorken flow and relaxation time approximation, the background distribution in ahydro can be formulated as with momentum anisotropy captured by a new variable ξ(τ).A gradient expansion around Eq. ( 7) gives rise to the viscous anisotropic hydrodynamics (vahydro) in Fig. 3. Comparing to viscous hydrodynamics with 2nd or 3rd order viscous corrections, convergence of ahydro and vahydro to the exact solution is impressive.Especially, one sees that vahydro applies to cases when 4πη/s approaches to order of 10 3 , which is far beyond local thermal equilibrium.One may refer to Ref. [28] for a detailed review of ahydro.
Recent progress in ahydro formulation has been achieved in various aspects [29,30,31,32,33].For instance, it is realized that one can close ahydro equations of motion by an extra equation matching the longitudinal pressure P L [34,35].

Hydrodynamic fluctuations
In addition to the initial state fluctuations which have quantum origins, there are also thermal fluctuations during the whole evolution stages of the quark-gluon medium.In the absence of these thermal fluctuations in hydrodynamic simulations, hydro modeling is defective.Inclusion of thermal fluctuations in fluid dynamics, namely, hydrodynamic fluctuations, has become a novel focus in the heavy-ion community, considering  especially its influence on the studies of small colliding systems and the QCD critical point searching in the beam energy scan program.
Hydrodynamic fluctuations can be formulated stochastically through the fluctuation-dissipation relations [37].With respect to the energy-momentum tensor T µν , one is allowed to add a stochastic tensor S µν , whose two-point auto-correlation is related to the corresponding dissipations.For the Navier-Stokes hydrodynamics [38], which clearly demonstrates that the strength of hydrodynamic fluctuations is determined by shear viscosity and bulk viscosity.Generally, for a system with larger dissipation, one would also expect a stronger effect of hydrodynamic fluctuations.Despite the complexity of implementing hydrodynamic fluctuations in a realistic hydro modeling of heavy-ion collisions, there have been several attempts towards an estimate of the effect of hydrodynamic fluctuations in heavy-ion collisions.For a 0+1D expanding conformal system, inclusion of hydrodynamic fluctuations can be reformulated to form a set of effective kinetic equations, which couple to the background Bjorken flow [40].It is then noticed that a fractional order, O(∇ 3/2 ), in the gradient expansion arises from hydrodynamic fluctuations, due to the non-linearities of the hydro equation of motion [41].For a 1+1D expanding conformal system, solution of the noisy viscous hydrodynamics can be achieved semianalytically [39], as a consequence of the analytically solvable background Gubser flow [42,43].In the context of heavy-ion collisions, it is first interesting to notice that the two-point auto-correlation of hydrodynamic fluctuations, Eq. ( 9), boils down to a form which inversely depends on the event-averaged multiplicity.Such a statement supports application of hydro modeling to small colliding systems, provide that multiplicity productions in these collisions are sufficiently high.Taking into account simultaneously contributions from initial state fluctuations, one can further explore the effect of hydrodynamic fluctuations on the flow harmonics.In Fig. 4, the additional increase of flow velocity anisotropies due to hydrodynamic fluctuations is plotted for the ultra-central proton-proton, proton-lead and lead-lead systems.Although the overall increase is not significant (smaller in lead-lead, larger in proton-proton), a clear trend of enhancement from the second order harmonics to the fourth order harmonics exhibits.
The statement that hydrodynamic fluctuations are more significant for higher order harmonics gets also examined in realistic numerical simulations [44].Fig. 5 presents the results of factorization ratio of flow harmonics of order 2 and 3, as a function of pseudo-rapidity, from one of the numerical hydro simulations

Nonlinear mode couplings in harmonic flow
Hydro modeling of the medium evolution in heavy-ion collisions is remarkably successful in calculating harmonic flow.A flow paradigm is established accordingly based on event-by-event simulations for nucleusnucleus collisions, from which relations among the harmonic flow V n , initial state geometric anisotropies E n , and the dynamical properties of the fluid medium, is proposed.Driven by the recent progress towards a precision measurement of flow harmonics, this flow paradigm needs to be examined at a more detailed level in the application of hydrodynamics.In particular, by doing so, one would expect to reduce uncertainties on the extracted transport coefficients of QGP.
Initial anisotropies in Eq. ( 5) are small quantities bounded by one, by definition.A linear eccentricity scaling relation can thereby be introduced between V n and E n , by cutting the expansion over E n at the linear order.Indeed, hydrodynamic simulations have verified that, to a good approximation, in nucleus-nucleus collisions [45,46], when collision centrality percentile is not very large.Note that the information of the initial state geometry is contained separately in E n , the coefficient κ n which we refer to as the linear response coefficient depends only on medium dynamical properties.Given this linear relation, one is allowed to relate, for instance, the event-by-event fluctuations of flow P(v n ), directly to that of initial anisotropies P(ε n ) [47,48,49].There are two cases, however, in which one has to consider beyond the linear eccentricity scaling relation in hydro modelings.The first involves lower order flow harmonics (n = 2, 3) in nucleus-nucleus collisions of large centralities.In these collision events, geometry of the overlapped region is so deformed that ε 2 cannot be treated as small.Constrained by rotational symmetry and analyticity of the relation, the correction to the linear eccentricity scaling relation of V n is of cubic order, where κ n is the cubic order response coefficient.Event-by-event hydro simulations have been carried out which found supports to the relation in Eq. ( 11) [50].The role of cubic response is not obvious in the rms magnitudes of elliptic flow V 2 .However, in the more differential analysis of the event-by-event v 2 distributions, it is noticed that cubic response is crucial in explaining the observed cumulant ratios.In Fig. 6, the ratio between the fourth order moment and square of the rms value of v 2 is plotted as a function of centrality percentile, for the Pb+Pb collisions at the LHC energy.Linear eccentricity scaling implies that the ratio is identical in v 2 and ε 2 .Comparing to the results calculated from the initial state modelings, the apparent difference observed when centrality becomes larger than 30% indicates the significance of cubic response.
In the other case concerning higher order flow harmonics (n ≥ 4), non-linearities are more essential.It is the consequence that the linear part of higher order flow get strongly suppressed by viscous corrections [21].Following the similar strategy mentioned before, one writes flow harmonics in an expansion of the couplings of E 2 and E 3 , with respect to rotational symmetry in the azimuthal direction.Further re-expressing E 2 by V 2 , and E 3 by V 3 , one has for V 4 and V 5 a separation between the linear part and the nonlinear couplings [51], where we have defined nonlinear flow response coefficients χ 422 and χ 523 .Similar relations can be obtained for V 6 and V 7 , where multiple nonlinear coupling terms and nonlinear flow response coefficients present.One notices that χ 422 and χ 523 are measurables in experiments which rely only on medium properties, with the assumption that the linear part and the nonlinear part of V 4 and V 5 are uncorrelated [51,52].These new observables have been measured in experiments at the LHC energy [53,54], in comparison with results from recent hydro simulations.Eq. ( 12) can also be applied to more complicated flow signatures involving higher order flow.For instance, it is realized that the symmetric cumulant, sc(n, m) which measures mixings of flow magnitudes, can be factorized into event-by-event flow fluctuations and event-plane correlations [7].

Summary
The success of hydrodynamic modeling of heavy-ion collisions is perhaps the most convincing evidence of the QGP fluidity in nucleus-nucleus collisions.In recent measurements at different colliding energies at RHIC and the LHC, soft probes of various types associated with the evolving medium were explored at an unprecedented level of resolution, which further pushes the limits of hydro modeling.Anisotropic hydrodynamics is developed, which abandons local pressure isotropization at the initial stages of the medium evolution.Effect of hydrodynamic fluctuations has been analyzed as well.At a qualitative level, a new term of the order of O(∇ 3/2 ) is discovered in the gradient expansion due to hydrodynamic fluctuations.Quantitative predictions are achieved regarding two-particle correlations in heavy-ion collisions, with more realistic simulations containing hydrodynamic fluctuations carried out.In addition to these theoretical developments, a flow paradigm which systematically analyzes the nonlinear couplings of harmonic flow is established.The fluid dynamical modelling of nuclear collisions continues to garner empirical success.The developments discussed herein should take this framework to its next level of sophistication.

Fig. 1 .
Fig. 1. Results of rms flow magnitude for harmonics n = 2, 3 and 4 from two separate groups of hydro simulations with initial state generated by EKRT (left) and IP-Glasma (right), with respect to the most recent measurements of Pb+Pb collisions at √ s NN = 5.02

Fig. 2 .
Fig.2.Pre-equilibrium evolution of a quark-gluon system in heavy-ion collisions characterized in terms of the pressure anisotropy P L /P T .

Fig. 3 .
Fig. 3. Convergence of ahydro and vahydro to the exact solution of Boltzmann equation with relaxation time approximation for a 0+1D expanding system.Figure from Ref. [36].

Fig. 4 .
Fig.4.Additional increase of flow velocity correlations induced by hydrodynamic fluctuations in a 1+1D expanding conformal system, with respect to harmonic order n = 2, 3 and 4, for the ultra-central proton-proton, proton-lead and lead-lead systems.Figures from Ref.[39].
Fig.4.Additional increase of flow velocity correlations induced by hydrodynamic fluctuations in a 1+1D expanding conformal system, with respect to harmonic order n = 2, 3 and 4, for the ultra-central proton-proton, proton-lead and lead-lead systems.Figures from Ref.[39].

Fig. 5 .
Fig. 5. Factorization breaking of the second and the third order harmonics, calculated from event-by-event hydrodynamics with also hydrodynamic fluctuations.Figures from Ref.[44].
Fig. 5. Factorization breaking of the second and the third order harmonics, calculated from event-by-event hydrodynamics with also hydrodynamic fluctuations.Figures from Ref.[44].
[50]6.Ratio of fourth order moment to the square of the rms value of v 2 magnitude.Experimental results are compared with event-by-event hydrodynamic simulations, and initial models where cubic response is absent.Figure from Ref.[50].