Theoretical developments on the initial state in relativistic particle collisions

. We discuss recent progress towards developing accurate initial state descriptions for heavy ion collisions focusing on weak coupling based approaches, that enable one to constrain the high-energy structure of nuclei from deep inelastic scattering or proton-nucleus collisions. We review recent developments to determine the event-by-event fluctuating nuclear geometry, to describe gluon saturation phenomena at next-to-leading order accuracy, and to include longitudinal dynamics to the initial state descriptions


Introduction
A crucial ingredient needed to simulate the space-time evolution in heavy ion collisions is the structure of the colliding nuclei at small momentum fraction x, which is the region probed in high-energy collisions.A consistent description of the heavy ion collision initial state together with e.g.deep inelastic scattering (DIS) data has been achieved in approaches based on collinear factorization and Color Glass Condensate (CGC).In the EKRT model based on collinear factorization [1] the partonic content of the nuclei is described in terms of nuclear parton distribution functions.In the CGC approach (implemented e.g. in the IP-Glasma [2] framework) the DIS and p+A cross sections and the time evolution of the color fields immediately after the heavy ion collision are described in terms of the universal Wilson line correlators.
There are also many other approaches to describe the initial state in heavy ion collisions including, for example, parametrization-based models such as T R ENTo [3] and different event generators (Pythia/Angantyr, EPOS, HIJING).In this contribution we, however, focus on weak coupling approaches with a direct connection to DIS and p+A collisions.

Probing nuclear geometry in photon-nucleus scattering
In heavy ion collisions the initial state eccentricities are transformed into momentum space anisotropies by the hydrodynamically evolving QGP.As such, a crucial input to the QGP simulations is the spatial distribution of nuclear matter at the initial condition and immediately after the collision before an approximatively thermalized QGP is formed.As such, there has been extensive activity in recent years to constrain the event-by-event fluctuating shape for of the proton and heavy nuclei.
Exclusive processes where the total momentum transfer is the Fourier conjugate to the impact parameter directly probe the nuclear geometry [4].Exclusive vector meson production and DVCS were extensively studied ad HERA.In recent years, vector-meson photoproduction has been studied in ultra peripheral photon-mediated collisions at RHIC and at the LHC where the photon-nucleus processes are available before the EIC era.Such measurements directly probe the nuclear geometry, down to very small x P ∼ 10 −5 , and are sensitive probes of non-linear QCD dynamics.Furthermore the potential to constrain nuclear PDFs in the poorly constrained small-x region has also been investigated recently [5].
When exclusive vector meson production has been measured as a function of momentum transfer [6,7], a spectrum that is more steeply falling compared to the one obtained as a Fourier transform of the Woods-Saxon density profile is obtained.This can be interpreted as a signature of saturation phenomena that effectively transform the nuclear density profile towards the black disc shape at small x P [8]. at a sub-nucleon scale.These results confirm the importance of sub-nucleon fluctuations to describe the measured incoherent J/y process at high energies, representing the first experimental step to use the quantum fluctuations of the gluon field to search for saturation effects in heavy nuclei.In addition, this measurement, when confronted to models, demonstrates that the contribution of the dissociative component to the total incoherent cross section depends on |t|.Thus, future analyses shall study the incoherent production of J/y as a function of rapidity and |t| [47].Finally, this analysis, together with recent measurements [18,20], indicate that new or improved theoretical models are needed to describe simultaneously the energy and |t|-dependence of both the coherent and the incoherent processes of J/y photoproduction, to gain a better understanding of saturation effects at a more fundamental level.
The ALICE Collaboration would like to thank all its engineers and technicians for their invaluable contributions to the construction of the experiment and the CERN accelerator teams for the outstanding performance of the LHC complex.In this conference, for the first time measurements of the incoherent J/ψ photoproduction cross section in ultra peripheral heavy ion collisions (i.e. in photon-nucleus collisions) as a function of the momentum transfer was reported by the ALICE and STAR collaborations [9].The incoherent production where the target dissociates is interesting, as it probes the event-by-event fluctuations in the target geometry [10].For example, the HERA data has been shown to prefer significant event-by-event geometry fluctuations for the proton [11].The new ALICE data is shown in Fig. 1, where it is compared to CGC calculations that either use spherical nucleons or include an event-by-event fluctuating nucleon geometry constrained by the HERA data [12].Although the overall cross section is overestimated by the theory calculation (i.e.not enough nuclear suppression is obtained), the tslope can be interpreted to prefer the calculation with nucleon substructure similar to that in protons at HERA kinematics.In addition to nucleon substructure fluctuations, in this conference recent progress towards probing the deformed structure of e.g.Uranium and Xenon in deep inelastic scattering was presented [13].

Gluon saturation at the precision level
Despite the fact that the leading order CGC calculations (that resum α s ln 1/x contributions to all orders) have been successful in describing large amount of small-x data, it is crucial to develop the theory to the next-to-leading order accuracy to enable precision level comparisons with the current and future measurements.Over the last couple of years, there has been an extensive effort in the community to bring the theory calculations describing the gluon saturation phenomena to the next-to-leading order accuracy.
At small-x cross sectios factorize to a convolution of Wilson lines and a hard impact factor.The energy dependence of the Wilson lines is described by perturbative Balitsky-Kovchegov or JIMWLK evolution equations.For a fully consistent NLO calculation, all these ingredients, including the non-perturbative initial condition for the high-energy evolution, need to be promoted to this order in α s .
The NLO BK evolution equation became available already in 2007 [14], and later the NLO JIMWLK equation has also been obtained [15,16] although for that there is currently no known method to solve it numerically.Typically the non-perturbative initial condition describing the proton structure at moderately small x has been extracted from fits to the proton structure function data [17].This became possible at NLO accuracy once the hard impact factor for DIS (the photon light front wave function at NLO) became available [18] (see also Ref. [19] for a complementary approach based on proton valence quark wave function).The first NLO fit has been reported in Ref. [20], and recently also a successful description of both the total and heavy quark production cross sections in DIS has been obtained [21].
In addition to total DIS cross sections, impact factors for many other scattering processes are currently known at NLO.These include, for example, exclusive vector meson production [22][23][24][25] and dijet/dihadron production [26][27][28]28] in DIS and inclusive hadron production in proton-nucleus collisions [29][30][31].Recently, first phenomenological applications at NLO accuracy have also become available.In particular consistent NLO calculations (with the caveat that the NLO BK evolution equation is approximated by a leading order equation into which dominant higher order corrections have been resummed) compared to available data exist for exclusive light and heavy vector meson production [22][23][24], inclusive π 0 production in proton-lead collisions [32] and dijet production in DIS [26] (although in that case the initial condition for the small-x evolution is not constrained by other collider data).Additionally, numerical results where leading order Wilson line correlators are used together with NLO impact factors exist for charged hadron production in proton-nucleus collisions [33].The obtained nuclear suppression factors for charged hadron and dijet production are shown in Figs. 2 and 3.
This rapid progress towards the NLO accuracy has brought the field to the point where precision level studies of saturation phenomena are becoming feasible.As saturation effects are typically expected to be only moderate [34], there is likely no smoking gun for gluon saturation even at the EIC.Instead, it will be crucial to perform global analyses where different DIS and p+A observables are simultaneously included.As there is always dependence on the non-perturbative input to the small-x evolution equation, in such analyses it will also be crucial to properly take into account uncertainties in this non-perturbative input (and in other non-perturbative ingredients such as in the vector meson wave function when calculating exclusive vector meson production).At the moment this non-perturbative input does not have any uncertainty estimates available at NLO [20], but first steps to include uncertainties in the extraction of the BK evolution initial condition [35] and the fluctuating proton geometry [36] at leading order have been taken recently.

Longitudinal dynamics
Many state-of-the-art descriptions of QGP evolution use 3+1D hydrodynamical simulations and hadronic afterburners.In order to fully describe longitudinal dynamics in heavy ion collisions, a realistic x-dependent initial condition is also necessary.This energy dependence has been included for example in the recent T R ENTo-3D [3] initial state parametrization.
In weak coupling approach one can again apply either collinear factorization or CGC to go beyond midrapidity.In the EKRT model [41] the input is an x-dependent nuclear PDF such as EPPS21 [42].In this conference, recent developments to the EKRT model were presented [38], including spatially dependent nuclear PDFs with event-by-event fluctuations, a dynamical event-by-event saturation criterium based on minijet production, minijet multiplicity fluctuations and global energy conservation.When this new 3D initial state description is R pPb , 3.0 < y < 3.5 Resummed Nuclear suppression factor for inclusive charged particle production in proton-lead collisions at the LHC compared to the LHCb data [37].Figure from Ref. [33].
unambiguously positive for the R range studied.Theory uncertainties in the NLO result can be divided into four classes; three of these are displayed in Fig. 2 (Bottom).We first show uncertainties from the unknown order N 2 LO contributions beyond the NLO impact factor; they are estimated by varying the running coupling scale c = 0.5 2 both in the NLO coefficient function where µ R = cP ?and in the Sudakov factor.Since they are parametrically of order ↵ 2 s ln 2 (P ?/µ 0 ), the band width grows with decreasing q ? .This illustrates the importance of controlling powers of ↵ s ln(P ?/µ 0 ) for future precision studies.
The second source of uncertainty are the missing contributions from the full NLO BK kernel.To gauge the sensitivity to these, we use two different formulations [85,86] of the kinematically constrained running coupling BK equation that differ by the additional resummation of single transverse logarithms [84].The blue area shows the corresponding sensitivity with relative variation of O(10%); including the full NLLx BK RGE can therefore significantly improve the overall precision of the computation.Thirdly, variations with respect to ↵ s,max are shown by the gray band in Fig. 2 (Bottom).Though as expected they grow at small q ?this sensitivity is mitigated, especially in large nuclei, because the scale (minimal transverse size) controlling the coupling is set instead by Q s .Lastly, power correction q 2 ?/P 2 ?, Q 2 s /P 2 ?uncertainties (not shown) previously discussed at LO [25,26] can be O(10%) for q ? .1.5 GeV and P ?= 4 GeV.Fig. 3 displays R eA , the ratio of the azimuthally averaged back-to-back dijet yield in e+A to e+p collisions.Such ratios minimize the aforementioned theory uncertainties as well as experimental ones.The top plot shows the q ?dependence of R eA for a large nucleus; for simplicity, we take A 1/3 = 6.At LO, it has a "Cronin" peak well-known from the corresponding ratio in protonnucleus (p+A) collisions [87]; in the CGC, it is generated by coherent multiple scattering that shifts the typical momentum imbalance to larger q ? in heavier nuclei [88].At NLO, we see that the Cronin enhancement is washed out by Sudakov corrections alone.A further strong effect is seen from the NLO contributions dominantly caused by the WW gluon TMD RGE which suppresses R eA analogously to the R pA case [89, 90].Qualitatively, Sudakov logs suppress configurations corresponding to small q ?(or large r bb 0 ) in the projectile.However since a fundamental consequence of gluon saturation is that even configurations with small r bb 0 are sensitive to nonlinear RG evolution with x, its precocious onset in large nuclei [91] leads to a suppression in R eA with A 1/3 .This is clearly demonstrated in the bottom plot.For fixed q ?= 1.5 GeV, one observes an increasing suppression with A 1/3 .The systematics of this suppression with A 1/3 and q ?are sensitive to the WW TMD RGE.Additional plots with different kinematic choices are provided in Figure 3. q? and A dependence (Top and Bottom respectively) of the nuclear modification factor ReA for the azimuthally averaged back-to-back dijet yield. the supplemental material.
While more detailed studies are necessary, our results are suggestive that inclusive back-to-back dijets in e+A collisions show strong potential to be a golden channel for gluon saturation at the EIC.Our conclusions can be strengthened by minimizing the stated theory uncertainties and by extending the comprehensive NLO study here to the di-hadron channel.Global analyses incorporating other e+A small x final states [85, 92-102] and analogous studies [103][104][105][106][107][108][109][110][111][112][113][114][115][116][117][118][119] in p+A collisions at RHIC and the LHC will further enable unambiguous determination of the dynamics of gluon saturation.
For simplicity we will focus sp region, where ⌘ s > 0 and large.moving nucleon has x 1 ∼ 1 wh cleon possesses the opposite be the characteristic scale, the sa distribution increases with decr pidities 1 will peak at small will be dominated by large mod in which p > q and expand a poses of this work, we can keep  [39].Figure adopted from Ref. [40].coupled to 3+1D relativistic hydrodynamics, a good description of key heavy ion observables away from midrapidity is obtained as illustrated in Fig. 4.
In CGC approach new developments towards a 3D initial condition were also presented.In the new McDipper initial condition [40] the initial quark and gluon production is calculated in a similar manner as inclusive particle production discussed in Sec. 3, with the proton smallx structure being described by a parametrization fitted to HERA data.Although currently the parton production is computed at leading order accuracy, extensions to higher order accuracy are in principle possible.Using initial state estimators, a good description of key heavy ion observables is obtained, although the flow decorrelation is underestimated as shown in Fig. 5.
The x-dependence can also be calculated perturbatively by solving the JIMWLK equation.In this conference, the initial state geometry and momentum correlations obtained from the JIMWLK evolution in Ref. [43] were presented.These results suggest that the initial momentum correlations are short-range in rapidity, unlike the event geometry for which cor-relations vanish much more slowly when the rapidity separation increases.The longitudinal structure of nuclei governed by the JIMWLK evolution can also be coupled to 3D classical Yang-Mills simulations to determine the time evolution after the collision before QGP is formed.First such implementation was recently shown in Ref. [44].When coupled to 3+1D hydrodynamical simulations, a good description of particle spectra, mean p T and flow harmonics at midrapidity is obtained, but again not enough longitudinal decorrelation for flow is obtained.This can be seen to suggest a need for an additional source of fluctuations.

Conclusions
The initial state of heavy ion collisions contains a vast amount of interesting fundamental physics, and a realistic initial state description is also crucial to probe in detail the QGP properties.At the moment there is a rapid development to include a realistic description for longitudinal dynamics which enables one to compare simulations of heavy ion collisions to observables away from midrapidity, and to understand the saturation effects at precision level.
The initial state can be inferred directly from heavy ion collisions, or probed in other scattering processes such as deep inelastic scattering or proton-nucleus collisions.A simultaneous description of heavy ion initial state and other collider data can be obtained in collinear factorization based approaches such as EKRT, or in Color Glass Condensate based implementations such as IP-Glasma.For EKRT, recent developments to include longitudinal dynamics was presented in this conference.In CGC based initial state descriptions, there is currently rapid progress in the field to include higher order corrections.These developments are crucial to both enable accurate studies of gluon saturation phenomena especially in the next decade when the Electron-Ion Collider [45] becomes operational, and to develop precise initial state descriptions for heavy ion collisions including a perturbatively calculated x-dependence.

Figure 2 :
Figure 2: Cross section for the incoherent photoproduction of J/y vector mesons in ultra-peripheral Pb-Pb collisions at p sNN = 5.02 TeV measured at midrapidity.The uncorrelated uncertainty (statistical and systematic added in quadrature) is indicated with the vertical bar, while the correlated uncertainty by the grey band.The width of each |t| range is given by the horizontal bars.The lines show the predictions of the different models described in the text.The bottom panel presents the ratio of the integral of the predicted to that of the measured cross section in each |t| range.The relative uncertainties on the ratios calculated from GSZ are 45%.

7 Figure 1 .
Figure 1.Incoherent J/ψ production as a function of squared momentum transfer |t| measured by ALICE compared to theory calculations with (MS-hs) and without (MS-p) nucleon substructure fluctuations.Figure from Ref. [9]

Figure 3 .
Figure 3. Nuclear suppression factor for dijet production in deep inelastic scattering in EIC kinematics as a function of dijet momentum q T .Figure from Ref. [26].

Figure 4 .
Figure 4. Pseudorapidity distribution of charged hadron multiplicity calculated from the 3D EKRT model[38] compared to the ALICE data.

Figure 5 .
Figure 5. decorrelation estimated form the initial spatial eccentricity in the McDipper model compared to the STAR data[39].Figure adopted from Ref.[40].