Correlation femtoscopy study at NICA and STAR energies within a viscous hydrodynamic plus cascade model

TL;DR

This paper investigates particle emission characteristics in heavy-ion collisions at various energies using a viscous hydrodynamic plus cascade model, linking femtoscopy measurements to the equation of state of hot QCD matter.

Contribution

It introduces femtoscopic observable calculations within a viscous hydro + cascade model across a range of collision energies, exploring the impact of different equations of state.

Findings

Femtoscopic observables vary with collision energy and EoS.

The model reproduces key features of experimental femtoscopy data.

Insights into the QCD phase diagram from energy-dependent femtoscopic measurements.

Abstract

Correlation femtoscopy allows one to measure the space-time characteristics of particle production in relativistic heavy-ion collisions due to the effects of quantum statistics (QS) and final state interactions (FSI). The main features of the femtoscopy measurements at top RHIC and LHC energies are considered as a manifestation of strong collective flow and are well interpreted within hydrodynamic models employing equation of state (EoS) with a crossover type transition between Quark-Gluon Plasma (QGP) and hadron gas phases. The femtoscopy at lower energies was intensively studied at AGS and SPS accelerators and is being studied now in the Beam Energy Scan program (BES) at the BNL Relativistic Heavy Ion Collider in the context of exploration of the QCD phase diagram. In this article we present femtoscopic observables calculated for Au-Au collisions at $\sqrt{s_{NN}} = 7.7 - 62.4$ GeV in a viscous hydro + cascade model vHLLE+UrQMD and their dependence on the EoS of thermalized matter.

Tables

Table 1. Table 1: Values of hydrodynamic starting time τ 0 subscript 𝜏 0 \tau_{0} , initial state granularity R ⟂ subscript 𝑅 perpendicular-to R_{\perp} , R η subscript 𝑅 𝜂 R_{\eta} and shear viscosity over entropy ratio η / s 𝜂 𝑠 \eta/s adjusted for different collision energies in order to reproduce basic observables in the RHIC BES region. An asterisk marks the values of τ 0 subscript 𝜏 0 \tau_{0} which are adjusted instead of being set directly from Eq.  1 .

$\sqrt{s_{\mathrm{NN}}}$ [GeV]	${\tau }_0$ [fm/c]	$R_{⟂}$ [fm]	$R_{\eta }$ [fm]	$\eta /s$
7.7	3.2	1.4	0.5	0.2
8.8 (SPS)	2.83	1.4	0.5	0.2
11.5	2.1	1.4	0.5	0.2
17.3 (SPS)	1.42	1.4	0.5	0.15
19.6	1.22	1.4	0.5	0.15
27	1.0	1.2	0.5	0.12
39	0.9*	1.0	0.7	0.08
62.4	0.7*	1.0	0.7	0.08
200	0.4*	1.0	1.0	0.08

Table 2. Table 2: Extracted average pion emission times t ¯ ¯ 𝑡 \bar{t} as a function of s N ​ N subscript 𝑠 𝑁 𝑁 \sqrt{s_{NN}} in the center-of-mass system of colliding gold nuclei depending on the EoS used.

$\sqrt{s_{NN}}$	EoS	particlization surface		last interactions
$\left[GeV\right]$		$\overline{t}$ [fm/c]	RMS [fm/c]	$\overline{t}$ [fm/c]	RMS [fm/c]
7.7	1PT	7.24	2.84	13.15	6.56
	XPT	6.16	2.01	11.61	6.26
11.5	1PT	7.33	2.31	13.09	6.92
	XPT	6.36	1.91	11.57	6.41
19.6	1PT	6.88	2.16	13.18	7.56
	XPT	6.41	2.15	11.93	6.93
27	1PT	6.85	2.37	13.38	8.07
	XPT	6.40	2.39	12.62	7.57
39	1PT	7.17	2.75	13.98	8.30
	XPT	6.64	2.58	13.05	7.85
62.4	1PT	7.00	2.82	14.11	8.50
	XPT	6.60	2.63	12.72	7.81

Table 3. Table 3: Results of single (Eq. ( 5 )) and double Gaussian fits of model source emission functions shown in Fig.  7 . χ 2 / n ​ d ​ f superscript 𝜒 2 𝑛 𝑑 𝑓 \chi^{2}/ndf values in parenthesis correspond to the XPT EoS. All calculations are performed in PRF.

Single Gaussian	${\chi }^2/ndf$	Radius	1PT	XPT
	975 (1247)	$R_{out}$ [fm]	11.10 $\pm$ 0.01	10.03 $\pm$ 0.01
7.7 GeV	300 (251)	$R_{side}$ [fm]	4.19 $\pm$ 0.01	4.16 $\pm$ 0.01
	3259 (3878)	$R_{long}\left[fm\right]$	4.59 $\pm$ 0.01	4.20 $\pm$ 0.01
Double Gaussian	${\chi }^2/ndf$	Radius	1PT	XPT
	51.6 (73.0)	$R_{out}^{small}$ [fm]	8.91 $\pm$ 0.01 (${\lambda }_{out}^{small}$=0.66)	7.83 $\pm$ 0.01 (${\lambda }_{out}^{small}$=0.68)
7.7 GeV	58.7 (66.6)	$R_{side}^{small}$ [fm]	4.10 $\pm$ 0.01 (${\lambda }_{side}^{small}$=0.99)	4.08 $\pm$ 0.01 (${\lambda }_{side}^{small}$=0.99)
	130.3 (195.1)	$R_{long}^{small}$ [fm]	3.10 $\pm$ 0.01 (${\lambda }_{long}^{small}$=0.64)	2.93 $\pm$ 0.01 (${\lambda }_{long}^{small}$=0.74)
		$R_{out}^{large}$ [fm]	13.85 $\pm$ 0.01 (${\lambda }_{out}^{large}$=0.34)	12.88 $\pm$ 0.01 (${\lambda }_{out}^{large}$=0.32)
		$R_{side}^{large}$ [fm]	9.76 $\pm$ 0.01 (${\lambda }_{side}^{large}$=0.01)	9.45 $\pm$ 0.01 (${\lambda }_{side}^{large}$=0.01)
		$R_{long}^{large}$ [fm]	6.06 $\pm$ 0.01 (${\lambda }_{long}^{large}$=0.36)	6.24 $\pm$ 0.01 (${\lambda }_{long}^{large}$=0.26)

Equations

$${\tau }_0=2R/\sqrt{(\sqrt{s_{\mathrm{NN}}}/2m_N{)}^2-1},$$

$$\Delta P_{ijk}^{\alpha }$$

$$\Delta N_{ijk}^0$$

$$C({\mathbf{k}}^{\ast },\mathbf{P})=\int {\mathrm{d}}^3{\mathbf{r}}^{\ast }{S}^{\alpha }({\mathbf{r}}^{\ast },\mathbf{P})\overline{{\left|{\psi }_{-{\mathbf{k}}^{\ast }}^{\mathcal{S},{\alpha }^{\prime }\alpha }({\mathbf{r}}^{\ast })\right|}^2},$$

$$S(\mathbf{r})\sim \mathbf{exp}\left(-\frac{{{\mathbf{r}}^{\ast }}_{\mathbf{out}}^2}{4{{\mathbf{R}}^{\ast }}_{\mathbf{out}}^2}-\frac{{{\mathbf{r}}^{\ast }}_{\mathbf{side}}^2}{4{{\mathbf{R}}^{\ast }}_{\mathbf{side}}^2}-\frac{{{\mathbf{r}}^{\ast }}_{\mathbf{long}}^2}{4{{\mathbf{R}}^{\ast }}_{\mathbf{long}}^2}\right).$$

$${\left|{\psi }_{-{\mathbf{k}}^{\ast }}({\mathbf{r}}^{\ast })\right|}^2={\left|[\exp (-i{\mathbf{k}}^{\ast }{\mathbf{r}}^{\ast })+\exp (i{\mathbf{k}}^{\ast }{\mathbf{r}}^{\ast })]/\sqrt{2}\right|}^2=$$

$$=1+\cos (2{\mathbf{k}}^{\ast }{\mathbf{r}}^{\ast }),$$

$$C(\mathbf{q})=N\left(1+\lambda \exp (-R_{\mathrm{out}}^2{q}_{\mathrm{out}}^2-R_{\mathrm{side}}^2{q}_{\mathrm{side}}^2-R_{\mathrm{long}}^2{q}_{\mathrm{long}}^2)\right).$$

$$S(r^{\ast })=\frac{{\sum }_{i\ne j}{\delta }_{\Delta }({\mathbf{r}}^{\ast }-{\mathbf{r}}_{\mathbf{i}}^{\ast }+{\mathbf{r}}_{\mathbf{j}}^{\ast })}{N{\Delta }^3}$$

Full text

relativistic heavy-ion collisions, hydrodynamics, collective phenomena, Monte Carlo simulations, vHLLE, UrQMD

††preprint: draft

I Introduction

One of the main motivations for heavy-ion collision programs is to study a new state of matter, the Quark-Gluon Plasma (QGP) which is defined as a deconfined state of quarks and gluons Cabibbo and Parisi (1975); Shuryak (1980). The systematics of transverse momentum spectra, elliptic and higher order flow coefficients measured in heavy-ion collisions at BNL Relativistic Heavy Ion Collider (RHIC) and CERN Large Hadron Collider (LHC) energies confirmed the presence of strong collective motion and hydrodynamic behavior of the system Arsene et al. (2005); Back et al. (2005); Adams et al. (2005); Adcox et al. (2005); Aamodt et al. (2010); Abelev et al. (2015); Adam et al. (2015); Aamodt et al. (2011). While the hydrodynamic approach was successful in reproduction of elliptic flow measured at the top RHIC energies from the very beginning, it was unable to reproduce the femtoscopic correlation measurements. The problem was solved several years ago along with substantial improvements in hydrodynamic modeling Sinyukov et al. (2009); Broniowski et al. (2008); Pratt (2009). The improvements comprise a presence of pre-thermal transverse flow, an inclusion of shear viscous corrections to hydrodynamic evolution, an equation of state compatible with recent lattice QCD calculations, and a consistent treatment of hadronic stage (hadronic cascade phase). As a result, existing state-of-the-art hydrodynamic models can reproduce, besides the transverse momentum distributions and elliptic flow coefficients, also the pion femtoscopic measurements Werner et al. (2010); Bożek (2011); Karpenko et al. (2013).

Recent lattice QCD calculations show that the transition from QGP to a hadron gas at high temperature and small $\mu_{B}$ is a crossover Brown et al. (1990); Aoki et al. (2006); Borsanyi et al. (2010); Cheng (2009), and there are no signs of critical behavior in the region of $\mu_{\rm B}/T<2$ Karsch et al. (2016). This is supported by a recent analysis of combined full RHIC ($\sqrt{s_{NN}}=200$ GeV) and LHC ($\sqrt{s_{NN}}=2.76$ TeV) data in viscous hydrodynamic + cascade model, where an elaborate model-to-data comparison using Bayesian framework suggests that the speed of sound at all temperatures cannot fall below the hadron resonance gas value of $\sim 0.15$, and that the resulting posterior distribution over possible equations of state of matter is compatible with the lattice QCD results Pratt et al. (2015). At the same time, there are predictions inspired by the lattice QCD calculations on a possible change of existing regime to a first-order phase transition occurring at lower energies and higher chemical potentials Masayuki and Koichi (1989); Barducci et al. (1994); Berges and Rajagopal (1999); Halasz et al. (1998); Scavenius et al. (2001); Antoniou and Kapoyannis (2003); Hatta and Ikeda (2003) and thus implying the existence of a critical point on the QCD phase diagram at a moderate value of chemical potential Stephanov (2005). These considerations motivated the BES program allowing one to study different parts of the QCD phase diagram at existing accelerators like SPS and RHIC, and, in future, at NICA and FAIR facilities.

It has been shown many years ago Pratt (1986); Bertsch et al. (1988); Rischke and Gyulassy (1996) that the long duration of particle emission related to a first order phase transition could reveal itself in the energy region of the onset of deconfinement as a strong increase of the Gaussian femtoscopic radius $R_{out}$, measured along the pair transverse momentum, compared with the nearly constant radius $R_{side}$, measured along the perpendicular direction in the transverse plane. As a result, one may expect a strong increase of the ratio $R_{out}/R_{side}$. In fact, a first order phase transition leads to a stalling of the mean expansion speed and a longer emission duration $\Delta\tau$ manifested as an increase of the radius measured along the beam direction $R_{long}$ and of the ratio of transverse femtoscopy radii $R_{out}/R_{side}$, respectively. A large data set of correlation functions of identical charged pions has been recently obtained by the STAR Collaboration within the RHIC BES at $\sqrt{s_{NN}}=$ 7.7, 11.5, 19.6, 27, 39, 62.4 GeV. The study of $R_{out}/R_{side}$ and $R_{out}^{2}-R_{side}^{2}$ behavior as a function of $\sqrt{s_{NN}}$ indicates a wide maximum near $\sqrt{s_{NN}}\sim$ 20 GeV, which is reported also by the PHENIX Collaboration Adare et al. . Could this wide maximum be related to the expected change of the type of phase transition?

To answer this and other questions, we study the sensitivity of Bose-Einstein correlations of identical pions to the EoS using a hybrid vHLLE+UrQMD model Karpenko et al. (2015). The model combines the UrQMD approach Bass et al. (1998); Bleicher et al. (1999) for the early and late stages of the evolution with numerical $(3+1)$-dimensional viscous hydrodynamical solution Karpenko et al. (2014) for the hot and dense expanding matter. A hydrodynamic approach has an essential advantage for the present analysis, since it allows one to simulate different scenarios of hadron / quark-gluon transition by changing the EoS and other transport coefficients input.

The paper is organized as follows: the details of the model are discussed in Section II; in Section III the femtoscopy formalism is described; in Section IV the results are presented and discussed; Section V is dedicated to conclusions.

II vHLLE+UrQMD model

The use of multi-component dynamical models for the description of dynamics of relativistic heavy-ion collisions at RHIC and higher energies is essential because hydrodynamics alone cannot describe the entire reaction. At the early stage of collision, thermalization of out-of-equilibrium quark-gluon system is assumed to occur, which allows one to describe the subsequent complex multi-particle dynamics using a relatively simple formalism of relativistic hydrodynamics. This requires an approach to calculate initial conditions for the hydrodynamic evolution111Note that for lower collision energies there exist one-, two-, and three-fluid hydrodynamical models (e.g. Ivanov et al. (2006)), which apply hydrodynamical description already for incoming nuclei.. As the matter expands, the characteristic mean free path of its constituents (quarks and gluons transforming into hadrons) becomes comparable to the system size. The interactions become less frequent, but they do not cease instantaneously, a process which can be simulated by switching from the hydrodynamical evolution to a hadronic cascade, usually with the help of Cooper-Frye formula Cooper and Frye (1974).

For the early stage of collision different approaches (or models) have been used in literature, such as (MC-)KLN, IP-Glasma, EPOS, HIJING, Glauber model etc. Those approaches have been developed for full RHIC or LHC energies, and lose their applicability as the collision energy decreases. As for the Glauber model, it can estimate initial energy density profiles in transverse direction only. Therefore we choose to use UrQMD to simulate the initial stage of the collision. We enforce a transition to hydrodynamical description at a hyper-surface of constant longitudinal proper time $\tau_{0}=\sqrt{t^{2}-z^{2}}$. The minimal value of the starting time $\tau_{0}$ is taken to be equal to the average time for the two colliding nuclei to completely pass through each other:

[TABLE]

where R is average radius of the nucleus and $m_{N}$ is nucleon mass.

At $\tau=\tau_{0}$ energy, momentum and baryon/electric charges of hadrons are distributed to fluid cells $ijk$ around each hadron’s position according to Gaussian profiles:

[TABLE]

where $P^{\alpha}$ and $N^{0}$ are 4-momentum and charge of a hadron, $\{\Delta x_{i},\Delta y_{j},\Delta\eta_{k}\}$ are the distances between hadron’s position and center of a hydro cell $ijk$ in each direction, $\gamma_{\eta}={\rm cosh}(y_{p}-\eta)$ is the longitudinal Lorentz factor of the hadron as seen in a frame moving with the rapidity $\eta$, and $C$ is a normalization constant. The normalization constant $C$ is calculated so that the discrete sum of energy depositions to the hydrodynamic cells equals to the energy of the hadron. The width parameters $R_{\perp}$ and $R_{\eta}$ control granularity of the produced initial state.

For all collision energies in consideration, the resulting initial energy density is large enough for the dense parts of the system to reside in the QGP phase.

The following 3-dimensional viscous hydrodynamic expansion is simulated with the vHLLE code Karpenko et al. (2014). Another input to the hydrodynamic part is the EoS, for which we use chiral model EoS Steinheimer et al. (2011) or bag model EoS Kolb et al. (2000). Whereas the present version of the chiral model EoS has a crossover type transition between QGP and hadronic phases for all baryon densities, the bag model EoS has a first order phase transition between the phases also for all baryon densities. Therefore below we dub chiral model EoS as “XPT EoS”, and bag model EoS as “1PT EoS”. For both EoS there are publicly available tables, computed in full physically allowed $T-\mu_{B}$ region, which makes them particularly useful for hydrodynamic computations with fluctuating initial conditions. Pressure as a function of energy density from both EoS is demonstrated on Fig. 1.

Fluid to particle transition, or particlization, is set to happen at a hypersurface of constant (hydrodynamic) energy density $\epsilon_{\rm sw}=0.5$ GeV/fm3, when the hydrodynamic EoS corresponds to hadronic phase. The particlization hypersurface is reconstructed with the CORNELIUS subroutine Huovinen and Petersen (2012). At this hypersurface, individual hadrons are sampled using the Cooper-Frye formula including shear viscous corrections to the distribution functions. The hadronic rescatterings and decays are treated with the UrQMD cascade.

The initial state parameters $R_{\perp}$, $R_{\eta}$, hydrodynamic starting time $\tau_{0}$ and shear viscosity over entropy ratio $\eta/s$ in fluid phase are tuned for different collision energies in order to approach basic experimental observables in the RHIC Beam Energy Scan region: (pseudo)rapidity distributions, transverse momentum spectra and elliptic flow coefficient Karpenko et al. (2015). The resulting values of the parameters are presented in Table 1. The tuning has been made with the XPT (chiral model) EoS, and in present work we use the same set of parameter values (i.e. do not re-tune them) for the simulations with the 1PT (bag model) EoS.

III Femtoscopy formalism

Since the first demonstration of the sensitivity of the Bose-Einstein correlations to the spatial scale of the emitting source done almost 60 years ago by G. Goldhaber, S. Goldhaber, W. Lee and A. Pais Goldhaber et al. (1960), the momentum correlation technique was successfully developed and is known presently as a “correlation femtoscopy”. It was successfully applied to the measurement of the space-time characteristics of particle production processes in high energy collisions, especially in heavy-ion collisions Podgoretsky (1989a, b); Lednicky (2004); Pratt (1984); Lisa et al. (2005). Femtoscopy correlations are studied by means of a two-particle correlation function. In a production process of a small enough phase space density, the correlations of two particles emitted with a small momentum $k^{*}=|{\bf k}^{*}|$ in the pair rest frame (PRF)222Calculations made in PRF are denoted by asterisk. are dominated by the effects of their mutual final state interaction (FSI) and quantum statistics (QS), depending on the PRF temporal ($t^{*}=t_{1}^{*}-t_{2}^{*}$) and spatial (${\bf r}^{*}={\bf r}_{1}^{*}-{\bf r}_{2}^{*}$) separation of particle emission points. Usually, one can neglect the temporal separation Lednický and Lyuboshits (1982); Lednický (2009) and in such equal-time approximation describe these effects by properly symmetrized wave functions at a given total pair spin $\cal{S}$, $\left[\psi_{-{\bf k}^{*}}^{\cal{S},\alpha^{\prime}\alpha}({\bf r}^{*})\right]^{*}$, representing solutions of the scattering problem viewed in the opposite time direction. So the complex conjugate, negative sign of the vector ${\bf k}^{*}={\bf p}_{1}^{*}=-{\bf p}_{2}^{*}$ and the detected channel $\alpha$ being the entrance one. Since the FSI factorization requires the FSI duration much larger than the particle production time, the relative momentum should be small also in the intermediate channels $\alpha^{\prime}$, so that one may consider the particles in these channels belonging to the same isospin multiplets as those in the detected channel $\alpha$. Particularly, for identical particles, the multi-channel problem reduces to the single elastic transition $\alpha\to\alpha$ only. Assuming further sufficiently smooth behavior of single-particle spectra in a narrow correlation region (smoothness assumption) Podgoretsky (1989a, b), one can neglect the space-time coherence and write the correlation function at a given ${\bf k}^{*}$ and pair three-momentum ${\bf P}$ as

[TABLE]

where the overline describes the averaging over the total pair spin $\cal{S}$ and summing over the intermediate channels $\alpha^{\prime}$. It is implied that particles are produced in a complex process with equilibrated spin and isospin projections and so the separation distribution (source function) $S^{\alpha}({\bf r}^{*},{\bf P})$ is independent of $\cal{S}$ and $\alpha^{\prime}$.

Experimentally, a two-particle correlation function is defined as a ratio $C({\bf q})=A({\bf q})/B({\bf q})$. $A({\bf q})$ is a measured distribution of the difference ${\bf q}={\bf p}_{1}-{\bf p}_{2}$, where ${\bf p}_{1}$ and ${\bf p}_{2}$ are three-momenta of two considered particles taken from the same event, while $B({\bf q})$ is a reference distribution of pairs of particles taken from different events. The momentum difference is usually calculated in the longitudinally co-moving system (LCMS), where the longitudinal pair momentum vanishes. The vector ${\bf q}$ is usually expressed in terms of $q_{\rm out}$ , $q_{\rm side}$ , $q_{\rm long}$, where the “long” axis is directed along the beam, “out” - along the pair transverse momentum and “side” - perpendicular to the latter one in the transverse plane.

To perform a quantitative analysis of femtoscopic correlations, an analytical form of $S$ is often used so that the result of the integration procedure in Eq. (4) can be compared with a correlation function $C$ obtained from an experiment. The source function is usually considered independent of the relative momentum ${\bf q}$ and its Gaussian shape is assumed:

[TABLE]

The widths in three directions (out, side and long) are called “Gaussian femtoscopy radii”. In LCMS $q_{\rm out}=\gamma_{t}(q_{out}^{*}+\beta_{t}(m_{1}^{2}-m_{2}^{2})/m_{12})$, $q_{side}=q_{side}^{*}$, $q_{long}=q_{long}^{*}$, where $\gamma_{t}$ and $\beta_{t}$ are the LCMS Lorentz factor and velocity of the pair. Since the space-time separation in PRF and LCMS is related by the Lorentz boost in the out direction: $r_{\rm out}^{*}=\gamma_{t}(r_{\rm out}-\beta_{t}t)$, $r_{\rm side}^{*}=r_{\rm side}$, $r_{\rm long}^{*}=r_{\rm long}$, the PRF and LCMS Gaussian radii coincide except for $R_{\rm out}^{*}=\gamma_{t}R_{\rm out}$. In present paper we consider the correlation function of two identical pions neglecting their FSI, so that

[TABLE]

and Eqs. (4, 5, 6) yield the 3-dim Gaussian form of the correlation function. This form is usually used to fit the LCMS Gaussian radii according to:

[TABLE]

where $N$ is the normalization factor and $\lambda$ is the correlation strength parameter, which can differ from unity due to the contribution of long-lived emitters and a non-Gaussian shape of the correlation function; $R_{out}$, $R_{side}$, $R_{long}$ are the Gaussian femtoscopy radii in the LCMS frame. Eq. (7) assumes azimuthal symmetry of the production process, which forbids the presence of the cross-terms except for $q_{out}q_{long}$. We neglect the latter assuming further the invariance under longitudinal boosts. Generally, in case of a correlation analysis with respect to the reaction plane, all three cross-terms $q_{i}q_{j}$ contribute.

The described fitting procedure allows one to compare extracted femtoscopy radii from the model with existing experimental data. This can be considered as a standard approach. A disadvantage of this approach is due to the fact that the Gaussian parametrization can suppress important information that could be derived from the long non-Gaussian tails of source functions. The PHENIX and STAR collaborations have recently started to apply a new “imaging technique” in order to extract directly the source function Adamczyk et al. (2013). In contrast to the standard approach, the source imaging allows one to extract a real non-Gaussian source function, being in this sense a model-independent one.

Here we will use the vHLLE+UrQMD model to study the effect of a non-Gaussian shape of the source function and its dependence on the nature of the phase transition. This model provides the information on particle four-momenta and four-coordinates of the emission points allowing one to calculate the correlation function with the help of the weight procedure. For non-interacting identical pions, the weight is given in Eq. (6).

IV Results and Discussion

As mentioned above, the parameters of the model were adjusted to approach experimental data for (pseudo-)rapidity distributions, transverse momentum spectra and elliptic flow coefficients with the XPT EoS, corresponding to a crossover type transition Karpenko et al. (2015). Then, switching the EoS from XPT to 1PT leads to only small differences in multiplicities and rapidity distributions of the produced hadrons in the model. At the same time, hydrodynamic evolution with the 1PT EoS leads to somewhat decreased mean $p_{\perp}$, elliptic and triangular flow coefficients Karpenko et al. (2016). Those trends are explained by a less violent transverse expansion with the 1PT EoS.

In this section, we study the space-time characteristics of the hadron emission in the model and present the results for the Bose-Einstein correlations of identical pions, obtained with the two aforementioned EoS’s in a wide collision energy range $\sqrt{s_{NN}}$ = 7.7 - 62.4 GeV covered by the BES program at RHIC.

IV.1 Pion emission time distributions

In the model one has access to the space-time points of particle production in last collisions and resonance decays, in addition to their momenta. In Fig. 2 we visualize the averaged time distributions of last interaction points of pions from the model simulations at the lowest, middle and highest collision energies ($\sqrt{s_{NN}}=7.7$, $19.6$ and $62.4$ GeV, respectively) using the 1PT and XPT EoS’s.

A detailed information on the pion emission times as a function of $\sqrt{s_{NN}}$ for all simulated collision energies and EoS’s is given in Table 2. The time distributions of midrapidity pions have been obtained at the particlization surface (points of their creation) and at the points of last interactions.

From the Table 2 one can conclude an apparently weak dependence of the average pion creation time $\bar{t}$ at the particlization surface on collision energy. It is an interplay of longer pre-thermal and shorter hydrodynamic stage at lower collision energies: at $\sqrt{s_{NN}}<39$ GeV, the hydro stage starts at $\tau_{0}=2R/(\gamma v_{z})$ when the two colliding nuclei have completely passed through each other, and the value of $\tau_{0}$ is as large as 3.2 fm/c at $\sqrt{s_{NN}}=7.7$ GeV. On the other hand, the duration of hydro stage becomes shorter as a collision energy decreases because of lower initial energy density at the hydro starting time. For the average pion creation times, the differences between 1PT EoS and XPT EoS are largest at the lowest collision energy in consideration. In addition, the 1PT EoS leads to larger root-mean-square (RMS) values of the time distributions, and the difference is again largest at the lowest collision energy. Because of re-scatterings and, more importantly, resonance decays in the final stage of hadronic cascade, the points of last interactions correspond to larger values of $\bar{t}$, which also depend weakly on collision energy. The cascade somewhat smears the relative difference between the 1PT and XPT scenarios, both for $\bar{t}$ and $RMS$.

IV.2 Three-dimensional correlation radii in LCMS.

An example of one-dimensional projections of three-dimensional correlation function obtained with the vHLLE+UrQMD model using the 1PT EoS and XPT EoS is shown in Fig. 3. The analysis involved the simulations performed for gold-gold collisions at $\sqrt{s_{NN}}=7.7$ GeV with applied cuts on event centrality of $0-5\%$ and the pair transverse momentum $k_{T}$. The second one pertains to the range of (0.15 - 0.25) GeV/c.

One can see that the pion correlation functions at small $q_{i}$ (i = “out”, “side”, “long”) are not well described by the Gaussian function in Eq. (7). The observed difference between correlation functions calculated with the 1PT and XPT EoS’s is noticeable in the “out” and “long” directions. In Fig. 4 this fact is demonstrated by the ratios of individual projections. The ratios in the “out” and “long” directions reach values up to 1.03 at small $q_{out}$ and $q_{long}$. A percent deviation from unity at small $q_{side}$ values appears due to the finite cuts on $q_{out}$ and $q_{long}$.

In Fig. 5 we present the $m_{T}$-dependence of the three-dimensional femtoscopy LCMS radii calculated at $\sqrt{s_{NN}}$ = 7.7, 11.5, 19.6, 27, 39, 62.4 GeV using the 1PT and XPT EoS’s, and a comparison of the obtained results with those ones obtained by the STAR collaboration Adamczyk et al. (2015).

One can see that the model reasonably describes the $m_{T}$-dependence of radii for all beam energies with both EoS’s. As for the radii, they show different trends in the “out”, “side” and “long” directions. Whereas the $R_{side}$ using both EoS’s practically coincide, the $R_{\rm out}$ with the 1PT EoS is generally larger (however, not more than 0.5 fm at any collision energy) than for the XPT EoS. This also leads to larger values of the $R_{out}/R_{side}$ ratio using the 1PT EoS. The difference comes from a weaker transverse flow developed in the fluid phase with the 1PT EoS as compared with the XPT EoS333A similar influence of the transverse flow on $R_{out}$ and $R_{side}$ has been observed for the RHIC and LHC energies Karpenko et al. (2013).. A longer lifetime of the fluid phase in the 1PT scenario also results in a larger values of $R_{long}$ as compared with the XPT scenario. Whereas one could expect that at lower collision energies in the 1PT EoS a larger fraction of the fluid phase evolution occurs in the mixed phase with zero speed of sound leading to an increase of evolution time and $R_{long}$, we did not observe such a trend in the model. The reason is that at lower collision energies in the model a sizable amount of radial flow is developed at pre-hydro stage. At the same time, the $R_{out}/R_{side}$ ratio at lowest collision energies shows a clear EoS dependence.

The $R_{out}/R_{side}$ and $R_{out}^{2}-R_{side}^{2}$ as a function of $\sqrt{s_{NN}}$ were studied at fixed $m_{T}$ by the STAR collaboration Adamczyk et al. (2015). A wide maximum near $\sqrt{s_{NN}}$ $\sim$ 20 GeV/c in both excitation functions was observed. This observation is however accompanied by rather large systematic error bars. We have calculated the very same quantities in the model and compared them with experimental data. The result of comparison is shown in Fig. 6.

One can see that due to large experimental error bars the model calculations involving the XPT EoS are in a strong agreement with the data within the error bars at all energies, whereas the 1PT EoS overestimates the data. However, in the model taking into account the XPT EoS we observe a monotonic increase in excitation functions of both quantities, meanwhile the 1PT EoS results in a non-decreasing behavior of the quantities. The XPT EoS “works” better for lowest collision energies that might be seen earlier from a better description of individual radii in that energy region shown in Fig. 5. A study of the $R_{out}/R_{side}$ ratio looks traditional in the modern femtoscopy since the $R_{out}$ and $R_{side}$ radii are both reduced by flow, thus their ratio is a more robust against the flow effects.

As mentioned above, the parameters of the model were adjusted to approach the basic hadronic observables: rapidity, transverse momentum distributions and elliptic flow coefficients within the BES region, but not femtoscopic ones. No model tuning has been made for the femtoscopy, therefore the obtained radii may be considered as a free “prediction” even though the experimental data already exists.

IV.3 Source emission functions

In a Monte Carlo model one has an access to the space-time characteristics of produced particles, which allows one to avoid the complicated procedure of solving the integral equation (see Eq. 4) as it is done in experiment Adamczyk et al. (2013). The source emission function can be calculated directly as:

[TABLE]

Here $r^{∗}_{i}$ and $r^{∗}_{j}$ are the particles space positions, $r^{∗}$ is the particles separation in PRF; $\delta_{\Delta}$ = 1 if $|x|<p/2$ and 0 otherwise, $p$ is a size of the histogram bin. The denominator in Eq. (8) takes care for the normalization by a product of the number of pairs $N=\sum_{i\neq j}1$ and a bin size $\Delta^{3}$.

Fig. 7 demonstrates an example of one-dimensional projections of source emission function $S(r^{*})$ derived from the model directly.

One can see that calculations involving the 1PT EoS lead to longer visible tails in the projections as compared with the XPT EoS, especially for “out” direction. It is related to a weaker transverse flow developed in the fluid phase and a longer lifetime of the fluid phase taking place in the 1PT EoS. A similar observation has been reported for the “out” femtoscopic radii in the previous section.

A set of functions consisting of a single Gaussian, double Gaussian, Gaussian + Exponential, Gaussian + Lorentzian, and Hump function Alt et al. (2010) was tested for a description of one-dimensional projections of source emission functions. The best description was obtained with the Hump-function and the double Gaussian function. The last one gives only slightly worse $\chi^{2}$ than the Hump-function, but allows one for a clear interpretation of parameters and a more stable fit. The single Gaussian and double Gaussian fit functions are shown in Fig. 7. The parameters extracted from these fits are presented in Table 3.

The fit of projections of source emission function with a single Gaussian gives a large value of $\chi^{2}/ndf$. The fit with a double Gaussian allows one to get much better values of $\chi^{2}/ndf$ and obtain a better description of the tails of projections of source emission functions until $\sim$ 60 fm in “out” and $\sim$ 25 fm in “side” and “long” directions, respectively. The radii extracted from the double Gaussian fit have a small component $R_{i}^{small}$ of 4-12 fm and a large component $R_{i}^{large}$ of 8-20 fm (as usual, $i$ denotes “out”, “side” and “long” directions). It reflects the fact that pion source consists of direct particles (described by the first component) and re-scatterings (the second one). The difference between radii extracted from the source emission functions obtained with the two EoS’s is seen for both components - $R_{i}^{small}$ and $R_{i}^{large}$, but it is rather small, less than 0.5 fm. The radii are larger for the 1PT scenario being consistent with the three-dimensional femtoscopic radii reported above.

It is interesting to note that in case of the single Gaussian fit the values of the radii are approximately equal to the ones derived from the double Gaussian fit and averaged quadratically over relative contributions of small and large radii. It means that the one-dimensional Gaussian radii roughly reflect the main features of double Gaussian fits.

Fig. 8 shows a $\sqrt{s_{NN}}$-dependence of the small and large radii and their relative contributions extracted from the double Gaussian fit depending on EoS.

The radii increase with increasing $\sqrt{s_{NN}}$ for both types of calculations. The visible difference between small and large radii in “out” direction (see Fig.8 (a)) decreases with increase of $\sqrt{s_{NN}}$. The relative contributions of small and large radii to “out” direction are equal to $\sim 0.65$ and $\sim 0.35$ (see Fig.8 (d)) and practically do not depend on either $\sqrt{s_{NN}}$ or type of EoS. The radii in “side” direction seem to be independent of $\sqrt{s_{NN}}$ (see Fig.8 (b), (e)). The radii in “long” projection almost coincide for both types of EoS (see Fig.8 (c)), but the relative contributions of them as a function of $\sqrt{s_{NN}}$ demonstrate a difference depending on EoS. The relative contribution of the large radii has a tendency to increase with $\sqrt{s_{NN}}$ and is larger in case of the 1PT scenario (see Fig.8 (f)).

Of course, the best comparison with experiment is a direct comparison of source emission functions from the model with the extracted ones experimentally. Nevertheless, this study shows us that the use of a double Gaussian fit also reflects a lot of interesting features of source emission functions, while a single Gaussian fit used in many experiments can be a rather risky procedure due to poor description of the source by this function. However, as it was demonstrated above, the radii extracted from the single Gaussian fit are equal to the properly averaged double Gaussian radii, giving, in principle, a realistic information on the source. This result is quite encouraging since it is much easier to study the three-dimensional radii than the source emission functions.

V Summary

We have presented the first study of pion femtoscopy in viscous hydro + cascade model vHLLE+UrQMD in the energy range of the BES program at RHIC. It is shown that the chiral model EoS Steinheimer et al. (2011) (XPT EoS), which has a crossover-type transition between QGP and hadron gas phases, in the fluid phase results in a quite reasonable reproduction of three-dimensional pion femtoscopic radii measured by the STAR collaboration.

The “out” Gaussian femtoscopic radii obtained with the bag model EoS (1PT EoS) are systematically larger as compared with the XPT EoS; the “side” radii coincide for both types of EoS; the “long” radii are also somewhat larger for the 1PT EoS.

The 1PT EoS results in a systematically worse reproduction of the data, however the differences between two EoS’s are not so large. The $R_{out}/R_{side}$ ratio and $R_{out}^{2}-R_{side}^{2}$ are in agreement with the STAR results within the error bars at all collision energies using the XPT EoS, but their energy dependences observed in the model are quite monotonic as opposed to the broad maximum around $\sqrt{s_{NN}}=20$ GeV reported by STAR. At the same time, the 1PT EoS overestimates experimental data points for both $R_{out}/R_{side}$ and $R_{out}^{2}-R_{side}^{2}$. In particular, the latter EoS does not reproduce the femtoscopic radii even at the lowest energy considered, $\sqrt{s_{NN}}=7.7$ GeV.

Parameters of the model were adjusted in Karpenko et al. (2015) based on rapidity, transverse momentum spectra and elliptic flow data in the BES region for the XPT EoS scenario. No readjustment for the 1PT EoS has been made, which poses an open question whether the differences in femtoscopic radii between the two EoS will be even smaller if the readjustment is made for each EoS scenario individually. Also, no additional parameter tuning has been made for the femtoscopic observables, therefore the results may be considered as “free model predictions” even though the experimental data already exists.

We find that a better overall description of the femtoscopic radii would require about 1 fm/c shorter duration of pion emission with the present setup of the model. If it is realized, then at lower energies the 1PT scenario will be closer to the data. This in turn may indicate a change of the nature of phase transition at energies less than $\sqrt{s_{\rm NN}}=20$ GeV and this should be an incentive for future experiments at NICA and FAIR facilities. It is an open question whether a new set of parameters more suitable for the femtoscopic radii description can be found.

In addition to traditional femtoscopic radii, we have calculated source emission functions of pion pairs. We show that it is possible to distinguish calculations with the two different EoS. The projections of source emission functions onto “out” direction are wider for the use of the 1PT EoS. For “side” direction these projections coincide for both scenarios; for “long” direction the projections obtained with the 1PT EoS are also wider in comparison with calculations using the XPT EoS. This observation is related to a weaker transverse flow developed in the fluid phase and a longer lifetime of the phase in case of the 1PT EoS used.

In order to describe the source emission functions quantitatively a set of different fitting functions has been tested. It is shown that the use of a double Gaussian fit to the source emission function gives a reasonable description and allows one for a simple interpretation of the obtained small and large radii.

So far we have performed femtoscopic analysis with vHLLE+UrQMD model only. As a next step we plan to extend the analysis using 3-fluid hydrodynamics-based event generator THESEUS Batyuk et al. (2016). In THESEUS the hydrodynamical description of heavy ion reaction starts earlier, which results in different sensitivity to hydrodynamic EoS especially in the NICA energy range.