3D Multi-system Bayesian Calibration with Energy Conservation to Study Rapidity-dependent Dynamics of Nuclear Collisions

Considerable information about the early-stage dynamics of heavy-ion collisions is encoded in the rapidity dependence of measurements. To leverage the large amount of experimental data, we perform a systematic analysis using three-dimensional hydrodynamic simulations of multiple collision systems -- large and small, symmetric and asymmetric. Specifically, we perform fully 3D multi-stage hydrodynamic simulations initialized by a parameterized model for rapidity-dependent energy deposition, which we calibrate on the hadron multiplicity and anisotropic flow coefficients. We utilize Bayesian inference to constrain properties of the early- and late- time dynamics of the system, and highlight the impact of enforcing global energy conservation in our 3D model.


Introduction
High-energy nucleus collisions at the Relativistic Heavy Ion Collider (RHIC) and the Large Hadron Collider (LHC) produce quark-gluon plasma (QGP).As it cools, the QGP undergoes rapid 3-dimensional hydrodynamic expansion.A popular and partly successful theoretical description of QGP evolution in heavy-ion collisions has been the longitudinal boost-invariant fluid, an approximate symmetry present in a significant region around mid-rapidity.Multistage models utilizing longitudinal boost invariance made use of 2D models of energy deposition in the initial state to compare with mid-rapidity measurements [1][2][3].In recent years 3D models of the initial state, coupled to three dimensional viscous hydrodynamics, have been developed and tested against experimental data extending beyond mid-rapidity [4][5][6].In this study, we use such a 3D model in conjunction with experimental data of high-energy collisions from four RHIC experiments in a Bayesian framework to quantify constraints on the initial state and the properties of the QGP.
In particular, we utilize a 3D model in a Bayesian inference framework to study the effect of the temperature dependence of the specific shear and bulk viscosity coefficients of the QGP on dynamics at forward and backward rapidities.We do this in a multi-system calibration, using data from Au-Au and d-Au collisions at the same energy (200 GeV).We highlight the impact of using an initial-state model featuring energy-momentum conservation, which propagates to the hydrodynamics the incoming energy in the collision.The present study serves as an exploratory first step toward a Bayesian system scan of all high-energy data probing the 3D dynamics of heavy-ion collisions at RHIC and eventually also the LHC.

The model
A variety of 3D initial state models have been developed in recent years, including extensions of the 2D T R ENTo and Monte Carlo Glauber models to three dimensions [5,7].In this study we use the 3D Glauber Model from Ref. [4,5].The positions of collisions in the transverse plane are determined using the same prescription as in the 2D Glauber model.At these positions, strings then extend in the longitudinal direction and are decelerated by a classical string tension [4,8].The energy of the string, which is then transferred to the medium, is proportional to the energy lost during the collision.A rapidity loss is defined as a function of the incoming rapidity of the participants and piece-wise parameterized at three incoming rapidity points.A larger rapidity loss corresponds to a greater amount of energy deposited in the midrapidity region.The initial state model is used to initialize the 3D hydrodynamic evolution in MUSIC, after which the UrQMD hadronic afterburner simulates decays and scatterings.In past analyses using different initial condition models, a normalization parameter was used to scale the deposited energy density to fit final charged hadron multiplicities for different energies or collision systems.In this analysis, for a given event, our model is constrained to a precise energy density profile in 3D due to the conservation of energy, and therefore, such a scaling parameter is not applicable.In this way, the same set of parameters is varied for both d-Au and Au-Au.

Bayesian Analysis and Model Emulation
A Bayesian analysis enables us to obtain probability distributions for model parameters that are constrained by a particular set of experimental measurements and their corresponding uncertainties-that is, obtain posterior probability distributions from mildly informed priors.Our analysis includes 20 model parameters.In order to make model predictions at such a large number of parameter sets, we make use of computationally-cheaper Gaussian Process Emulators (GPs) at the expense of introducing emulator uncertainty.After dimensional reduction of the observable space, GP surrogates are trained to reproduce the model predictions for a set of 10 principal components.The emulators are systematically validated against a set of model predictions for parameter sets not used in emulator training.Additionally, the analysis as a whole is validated by means of closure tests, in which the model predictions at given parameter sets are used as pseudo-experimental data.We look for posterior distributions of particularly sensitive parameters to converge on the "truth" values of the parameter set used for closure.
The experimental data used in this first iteration of a 3D analysis are the charged hadron multiplicity and elliptic flow measurements as a function of pseudorapidity for various centrality bins in both Au-Au and d-Au collisions from the PHOBOS, BRAHMS, STAR, and PHENIX experiments [9][10][11][12].The different data sets extend to varying pseudorapidiy, though all go well beyond the mid-rapidity region.We use the 20-70% and 0-5% centrality bins for the v 2 in Au-Au and d-Au, respectively.Experimental uncertainties are taken to be uncorrelated across centrality and pseudorapidity bins.

Constraints on Initial State and Transport Properties
Constraints on the rapidity loss parameters are shown on the left side of Fig. 1.As expected, the strongest constraints are in the incoming rapidity region corresponding to the beam rapidity at 200 GeV.We verified that most observables are strongly sensitive to the value of rapidity loss parameters at 4 and 6 units, consistent with the recent Bayesian analysis in [13].
The analysis also shows little constraining power in regions of lower incoming beam rapidity.
The data also indicate a preference for a specified large bulk viscosity, as shown on the right panel of Fig. 1.The preferred viscosity is so large that it extends beyond the 90% prior for ζ/s.We verified (not shown) that a large bulk viscosity is found whether one uses on only Au-Au or only d-Au data.
For completeness, we repeated the Bayesian analysis using only mid-rapidity measurements (|η| < 1).The result is shown as the green dashed line in the right panel of Fig. 1.We find that one can describe mid-rapidity data with small values of ζ/s, indicating that the large values of the specific bulk viscosity found in the 3D analysis are truly from the rapiditydependent measurements.In fact, the constraints on the bulk viscosity found with the midrapidity measurements are generally consistent with the ones found in a previous JETSCAPE calibration [2].This is not obvious, since the model and data sets used in the present work and in Ref [2] are rather different; further verification will be necessary to ascertain that this conclusion remains the same after more measurements are added to the analysis.

Normalization of Initial State Energy Deposition
To successfully describe charged hadron multiplicity measurements, previous analyses scaled the initial state energy density by a normalization parameter that may be fit to each system or energy independently.In the absence of such an explicit normalization in a Bayesian analysis, the remaining parameters sensitive to the multiplicity compensate for a good fit to the multiplicity data in different systems.The most sensitive are the rapidity loss parameters, which define a 3D energy profile.Another one is the bulk viscosity, which affects the final multiplicities through viscous heating of the system during the hydrodynamic evolution.When the model is constrained only at mid-rapidity, the rapidity loss acts in effect as a normalization parameter to match with the multiplicity at mid-rapidity.Indeed, as we noted, when we elect to use only mid-rapidity data in our analysis, we obtain constraints on the bulk viscosity similar to those of previous 2D analyses.When forward/backward rapidity data is included in the analysis, however, the rapidity loss is constrained by global conservation of energy, and an absolute scale degree of freedom is no longer allowed.To optimize the absolute scale to match the magnitude of the multiplicity, the bulk viscosity at low temperature needs to be large.The results we see in this work are consistent with this interpretation.By contrast, we did not find strong constraints on the shear viscosity (not shown) even when rapidity-dependent observables are included.

Outlook
In this work, we have performed a multi-system Bayesian comparison of rapidity-dependent measurements of charged hadron multiplicities and elliptic flow with a 3D model of the initial state and hydrodynamic evolution of high-energy nuclear collisions.We have noted the strong constraints on the rapidity loss and the system dependence of the calibration.We have, in particular, seen indications of a preference for a large bulk viscosity and discussed the role of energy-momentum conservation in a 3D model on this result.The JETSCAPE collaboration is actively working on 3D calibrations using a more extensive set of available experimental data.
In addition, a 3D calibration of the properties of the QGP fluid throughout its space-time evolution utilizing all available soft hadron data is a crucial baseline for further studies that couple hard physics to a hydrodynamic simulation for a complete description of heavy ion collisions in realistic jet+hydro simulations, a long-standing goal of the community.

Figure 1 .
Figure 1.The 90% confidence intervals (C.I.) for the prior and posterior of the rapidity loss as a function of the incoming rapidity (left) and specific bulk viscosity as a function of temperature (right).

Table 1 .
System and observable breakdown of the experimental data used in this analysis.