Revealing the collision energy dependence of $\eta/s$ in RHIC-BES Au+Au collisions using Bayesian statistics

TL;DR

This study uses Bayesian statistics to analyze how the shear viscosity to entropy ratio ($s$) varies with collision energy in RHIC Au+Au collisions, revealing a higher s at lower energies, influenced by baryon chemical potential.

Contribution

It introduces a Bayesian framework to quantify the energy dependence of s in heavy-ion collisions, highlighting the role of baryon chemical potential.

Findings

s increases at lower collision energies.

Bayesian analysis constrains s values with posterior distributions.

Dependence of s on baryon chemical potential  is suggested.

Abstract

We investigate the collision energy dependence of $\eta/s$ in a transport + viscous hydrodynamics hybrid model. A Bayesian analysis is performed on RHIC beam energy scan data for $Au+Au$ collisions at $\sqrt{s_{NN}}=19.6,39,$ and $62.4$ GeV. The resulting posterior probability distributions for the model parameters show a preference for a larger value of $\eta/s$ at $19.6$ GeV compared to $62.4$ GeV, indicating dependence on baryon chemical potential $\mu_B$.

Figures

Equations

$$\mathcal{L}(x)\propto \exp \left(-\frac{1}{2}(y(x)-y^{\text{exp}}{)}^T{\Sigma }^{-1}(x)(y(x)-y^{\text{exp}})\right)$$

Full text

I Introduction

Significant progress has been made in the past few years in determining QGP properties, such as the temperature dependence of shear viscosity over entropy density ratio $\eta/s$ Csernai:2006zz ; Niemi:2015qia ; Denicol:2015nhu ; Plumari:2016lul . However, in a recent RHIC beam energy scan study Karpenko:2015xea it was found that the best fit with the experimental data was reached using larger values of $\eta/s$ at lower collision energies, suggesting an additional dependency on the baryon chemical potential $\mu_{B}$.

It is generally difficult to determine the uncertainties associated with the extracted best-fit values of QGP properties, as the computational models used in the analysis typically have numerous interconnected parameters, which need to be tuned on large sets of experimental data. We tackle this issue in the present study by utilizing Bayesian statistics combined with Markov chain Monte Carlo (MCMC) to calibrate the computational model to data. The end result of such an analysis is a multidimensional probability distribution, which provides not only a set of data-calibrated parameter values, but a full uncertainty quantification as well. This approach has already been applied with great success to $Pb+Pb$ collisions at the LHC Bernhard:2016tnd .

In this article, we investigate the $\mu_{B}$ dependence of $\eta/s$, using the collision energy $\sqrt{s_{NN}}$ as the control parameter. We simulate $Au+Au$ collisions at $\sqrt{s_{NN}}=19.6,39,$ and $62.4$ GeV, using the same transport + (3+1)-D viscous hydrodynamics model as in Ref. Karpenko:2015xea . In addition of introducing the robust uncertainty quantification method described above, this revised analysis benefits also from the recently published experimental data on identified particle multiplicities and mean transverse momentum Adamczyk:2017iwn .

II Transport + hydrodynamics hybrid model

In the utilized hybrid model, the initial non-equilibrium evolution is handled by the UrQMD hadron-string cascade Bass:1998ca ; Bleicher:1999xi . The transition from hadron-string transport to hydrodynamics happens some time after the two colliding nuclei have passed through each other. Beyond this condition, the exact transition time is a free parameter of the model. At the transport-to-hydro transition, the microscopic particle properties are mapped to hydrodynamic densities using 3-D Gaussians with width parameters $R_{\textrm{trans}}$, $R_{\textrm{long}}$.

The hydrodynamic evolution is done using a (3+1)-D relativistic viscous hydrodynamics code Karpenko:2013wva with a constant value of $\eta/s$, which is given as an input. A chiral model equation of state Steinheimer:2010ib is used during the hydro evolution to accommodate for finite values of baryon chemical potential.

The transition from hydrodynamics back to hadron transport happens when the local rest frame energy density in hydro cells falls below the user-defined switching value $\epsilon_{SW}$. The ”Cornelius” routine Huovinen:2012is is used to construct the iso-energy density hypersurface, from which hadrons are sampled according to Cooper-Frye procedure. The rescatterings and resonance decays of the sampled hadrons are processed within UrQMD.

III Statistical analysis

According to Bayes’ theorem, the posterior probability for the model parameters is a product of our prior knowledge about the plausible range of values for each parameter, and the likelihood of any given combination of parameter values being the ”correct” one when compared to experimental data:

[TABLE]

where $\Sigma$ is the covariance matrix, representing the uncertainties related to the comparison of model output $y(x)$ and the data $y^{\textrm{\,exp}}$ for the input parameter combination $x$.

In practice, we produce samples of the Bayesian posterior distribution using Markov chain Monte Carlo, where the initial positions of the random walkers are based on the prior (in this case, a uniform distribution in a 5-dimensional hypercube), and the likelihood function determines the probability for a walker to accept a proposed step direction in the parameter space. We use $\mathcal{O}(1000)$ random walkers to produce a sufficient number of samples of the posterior distribution.

The MCMC method necessitates fast evaluations of the likelihood function, but running the actual hybrid model simulations to determine $y(x)$ for arbitrary $x$ requires an infeasible amount of computational effort. We circumvent this problem by utilizing Gaussian process (GP) emulators, which provide an efficient method for predicting $y(x)$ with quantitative uncertainty.

In order to use GP emulators, they need to be conditioned on training data. For this analysis, approximately 100 training data points were produced for each collision energy. A Latin hypercube method was used to sample the training points, to ensure a representative sample of the whole input parameter space.

IV Results

The analysis was performed using data for charged particle multiplicity $N_{ch}$ Abelev:2008ab , pseudorapidity distribution $dN_{ch}/d\eta$ Alver:2010ck , and elliptic flow $v_{2}\{\textrm{EP}\}$ ($\sqrt{s_{NN}}=62.4$ GeV) Adamczyk:2012ku ; Magdy:2017kji or $v_{2}\{2\}$ ($\sqrt{s_{NN}}=19.6$ and $39$ GeV) Adamczyk:2012ku . Identified particle observables include charged pion HBT radii $R_{\textrm{out}}$, $R_{\textrm{side}}$ and $R_{\textrm{long}}$ Adamczyk:2014mxp , $K^{+}/\pi^{+}$ ratio Abelev:2008ab , and multiplicities and mean transverse momentum of $\pi^{+}$, $\pi^{-}$, $K^{+}$, $K^{-}$ Adamczyk:2017iwn and $\Omega$ Aggarwal:2010ig ; Chatterjee:2015fua ; Adamczyk:2015lvo . Proton data was not included in this analysis, as proton yields suffer from additional uncertainties related to feed-down corrections.

We have increased the weight of $v_{2}$ in the analysis by a factor of 5 and $N(\Omega)$ by a factor of 4, as these observables have been found to improve constraints on $\eta/s$ and $\epsilon_{SW}$, respectively. As an example, figure 1 shows the effect of weighting on $v_{2}$ and $N(\Omega)$ at $\sqrt{s_{NN}}=39$ GeV. The range of emulator predictions for model outputs, based on 100 random draws from the posterior distributions, is clearly in better agreement with the measured values when extra weighting is introduced. The emulator quality is demonstrated in the figure by showing both the GP predictions and the real model output, when the median values of posterior probability distributions were used as input.

Figure 2 shows one-dimensional projections of posterior probability distributions for $\eta/s$ for the three collision energies, with and without weighting. If no weighting is introduced, shear viscosity remains largely unrestricted, especially at $\sqrt{s_{NN}}=39$ GeV. At $\sqrt{s_{NN}}=19.6$ GeV, there is a visible preference towards values over 0.2, while for $\sqrt{s_{NN}}=62.4$ GeV the probability distribution peaks at $\eta/s\approx 0$. These tendencies are more emphasized when the importance of $v_{2}$ and $\Omega$ observables are increased in the analysis, while the posterior distribution at $\sqrt{s_{NN}}=39$ GeV develops a peak at intermediate value $\eta/s\approx 0.1$.

We present a rough illustration of the collision energy dependencies of the posterior distributions of all five model parameters in figure 3. Even with the introduction of additional weighting, the 90% confidence limits remain large for all parameters. Especially the width parameters, controlling the creation of the initial density profile for the hydrodynamics, are poorly constrained by the data. This uncertainty about the initial state is naturally reflected also on the other parameters.

The differences between peak and median values (open and solid symbols in Fig. 3) indicate that the distributions are skewed, as already illustrated for $\eta/s$ in figure 2. The consistent right-skewness of transverse smearing $R_{\textrm{trans}}$ suggests that the optimal value is more likely to be 1.0 fm or less. The strongest statement can be made on the value of hydro-to-transport switching energy density, which has both the peak and the median on the range $\approx 0.3-0.4$ GeV/fm3 for all investigated collision energies, making it very likely that the optimal value will be found within this range.

V Summary

Utilizing Bayesian statistics and Gaussian process emulators, we have performed a state-of-the-art model-to-data comparison on a transport + hydrodynamics hybrid model describing $Au+Au$ collisions at RHIC beam energy scan. While the collision energy dependence of the posterior probability distribution strongly suggests that shear viscosity over entropy density ratio $\eta/s$ depends on baryon chemical potential $\mu_{B}$, the present uncertainties prohibit precise statements about the ”correct” parameter values. We expect the situation to improve once the dynamics of the initial fluidization are better understood.

VI Acknowledgements

This work has been supported by NSF grant no. NSF-ACI-1550225 and by DOE grant no. DE-FG02-05ER41367. CPU time was provided by the Open Science Grid, supported by DOE and NSF.