Searching for the QCD critical point with net-proton number fluctuations

TL;DR

This paper investigates net-proton number fluctuations in heavy ion collisions to identify signs of the QCD critical point, using a phenomenological model that reproduces experimental trends in cumulant ratios.

Contribution

It introduces a model coupling critical mode fluctuations to protons and anti-protons, capturing both monotonic and non-monotonic behaviors observed experimentally.

Findings

The model reproduces the monotonic behavior of the first cumulant ratio.

The model captures the non-monotonic behavior of higher-order cumulant ratios.

Results depend on the coupling strength and the critical point location.

Abstract

Net-proton number fluctuations can be measured experimentally and hence provide a source of important information about the matter created during relativistic heavy ion collisions. Particularly, they may give us clues about the conjectured QCD critical point. In this work the beam-energy dependence of ratios of the first four cumulants of the net-proton number is discussed. These quantities are calculated using a phenomenologically motivated model in which critical mode fluctuations couple to protons and anti-protons. Our model qualitatively captures both the monotonic behavior of the lowest-order ratio as well as the non-monotonic behavior of higher-order ratios, as seen in the experimental data from the STAR Collaboration. We also discuss the dependence of our results on the coupling strength and the location of the critical point.

Tables

Table 1. Table 1. : Locations of the QCD critical point in the ( μ B , T subscript 𝜇 𝐵 𝑇 \mu_{B},T )-plane considered in this work. These locations in the QCD phase diagram are shown in Fig.  2 .

$C{P}_i$	${\mu }_{cp}$[MeV]	$T_{cp}$[MeV]
1	390	149
2	420	141
3	450	134

Equations

$$n_i(T,{\mu }_i)=d_i\int \frac{d^{3}k}{(2\pi {)}^3}{f}_i^0(T,{\mu }_i).$$

$$f_i^0=\frac{1}{(-1{)}^{B_i}+e^{(E_i-{\mu }_i)/T}}$$

$$m_i\sim m_0+g_i\sigma ,$$

$$\delta f_i=\frac{\partial f_i}{\partial m_i}\delta m_i=-\frac{g_i}{T}\frac{v_i^2}{{\gamma }_i}\delta \sigma ,$$

$$C^{i}_{n}=VT^{3}\frac{\partial^{n-1}(n_{i}/T^{3})}{\partial(\mu_{i}/T)^{n-1}}\bigg{|}_{T}\,,$$

$$C_n$$

$$J_i=\frac{gd}{T}\int \frac{d^{3}k}{(2\pi {)}^3}\frac{1}{E}{f}_i^0(1-f_i^0).$$

$$⟨(V\delta \sigma {)}^n{⟩}_c={\left(\frac{T}{V{H}_0}\right)}^{n-1}{\frac{{\partial }^{n-1}{M}_I}{\partial h^{n-1}}|}_r,$$

$$C_2^{\text{sing.}}\propto \frac{\partial M_I}{\partial h}.$$

$${\chi }_{\mu \mu }\simeq {\chi }_{\mu \mu }^{reg}+{\sigma }^2{\chi }_{\text{chiral}},$$

$$C_2=C_2^p+C_2^{\bar{p}}+g^2{\sigma }^2⟨(V\delta \sigma {)}^n⟩(J_p-J_{\bar{p}}{)}^2.$$

$$C_3=C_3^p-C_3^{\bar{p}}-g^3{\sigma }^3⟨(V\delta \sigma {)}^n⟩(J_p-J_{\bar{p}}{)}^3$$

$$C_4=C_4^p+C_4^{\bar{p}}+g^4{\sigma }^4⟨(V\delta \sigma {)}^n⟩(J_p-J_{\bar{p}}{)}^4.$$

$$\frac{C_2}{C_1}=\frac{{\sigma }^2}{M},\frac{C_3}{C_2}=S\sigma ,\frac{C_4}{C_2}=\kappa {\sigma }^2,$$

Full text

net-proton number fluctuations, QCD critical point, heavy-ion collisions

1 Introduction

The theoretical and experimental investigation of the phase diagram of strongly interacting matter is an important subject of modern high energy physics. One of the unresolved questions concerns the existence and location of the QCD critical point (CP) in the $T$ and $\mu$ plane. Strong fluctuations of the critical mode, $\sigma$, in the vicinity of CP, although not directly observable, are expected to couple to physically measurable quantities, such as fluctuations of conserved charges [1,​ 2].

Fluctuations of the net-proton number serve as an experimental probe of baryon number fluctuations. Recent, but still preliminary results of the STAR Collaboration [3,​ 4,​ 5] show a non-monotonic beam energy dependence of the ratios of higher order net-proton number cumulants. However, the interpretation of the data is still unclear [6,​ 7,​ 8,​ 9] and therefore effective models are needed to improve our understanding of these quantities.

One of such models was developed in [10], where the impact of resonance decays on net-proton number cumulant ratios was studied. This model could qualitatively describe the non-monotonic behavior of the $C_{3}/C_{2}$ and $C_{4}/C_{2}$ ratios. However, it also showed a strong non-monotonic behavior of the $C_{2}/C_{1}$ ratio which is not observed experimentally. Recently, this model was re-examined [11] to take into account the scaling properties of the baryon number and chiral susceptibilities obtained within effective models [12,​ 13]. This reduces the effect of critical fluctuations in the net-proton number variance and thus allows for a better description of the STAR data.

Here we discuss the beam energy dependence of the ratios of net-proton number cumulants obtained using the refined model from Ref. [11] and study their dependence on the coupling strength between critical mode and (anti)protons as well as their dependence on the location of the critical point.

2 Model setup

As a baseline model to calculate the net-proton number cumulants we choose the hadron resonance gas (HRG) model in which the number density of each particle species is given by the ideal gas formula,

[TABLE]

Here $d_{i}$ is the degeneracy factor and

[TABLE]

is the equilibrium distribution function, where $E_{i}=\sqrt{\vec{p}^{2}+m_{i}^{2}}$ is the dispersion relation and $\mu_{i}=B_{i}\mu_{B}+S_{i}\mu_{S}+Q_{i}\mu_{Q}$ the chemical potential of a particle of mass $m_{i}$, baryon number $B_{i}$, strangeness $S_{i}$ and electric charge $Q_{i}$. $\mu_{B}$, $\mu_{Q}$ and $\mu_{S}$ denote baryon, strangeness and charge chemical potentials.

Since in HRG model the QCD pressure is approximated by a sum of partial ideal gas pressures corresponding to different particles, there are only thermal fluctuations in this approximation. To include critical fluctuations on top of the thermal ones we follow the phenomenological approach employed in Ref. [10]. In this approach the particle mass is assumed to be composed of critical and non-critical parts as suggested in linear sigma models,

[TABLE]

where $m_{0}$ is a non-critical contribution and $g_{i}$ is the coupling strength between the critical mode and the particle of type $i$. Critical mode fluctuations modify the distribution function into $f_{i}=f_{i}^{0}+\delta f_{i}$, where the change of the distribution function due to critical mode fluctuations reads

[TABLE]

with $v_{i}^{2}=f_{i}^{0}((-1)^{B_{i}}f_{i}^{0}+1)$ and $\gamma_{i}=E_{i}/m_{i}$.

Fluctuations of the particle number in the thermal medium can be quantified in terms of cumulants. The n-th order cumulant of i-th particle species reads

[TABLE]

where temperature $T$ is kept constant. In this work we consider the first four cumulants of the net-proton number, $N_{p-\bar{p}}=N_{p}-N_{\bar{p}}\,$, which are given by [10]

[TABLE]

where $C_{n}^{p}$ and $C_{n}^{\bar{p}}$ are n-th order proton and anti-proton cumulants obtained within the baseline model, respectively, $\langle(V\delta\sigma)^{n}\rangle_{c}$ is n-th critical mode cumulant and

[TABLE]

Moreover, the contributions of other particles and resonance decays are neglected.

In general, cumulants of the critical mode cannot be calculated analytically. Following the approach introduced in Ref. [10], we model them using universality class arguments which state that different physical systems belonging to the same universality class exhibit the same critical behavior close to the critical point [14]. Under the assumption that QCD belongs to the same universality class as the three-dimensional Ising model [16,​ 17,​ 18] we can identify the QCD order parameter, $\sigma$, with the magnetization, $M_{I}$, the order parameter of the spin model. Hence, the critical mode cumulants can be written as [10]

[TABLE]

where $r=(T-T_{c})/T_{c}$ is the reduced temperature and $h=H/H_{0}$ the reduced magnetic field. The critical point is located at $r=h=0$.

In the net-proton number cumulants, the singular part of the second cumulant receives a contribution from the first derivative of the order parameter with respect to the reduced magnetic field,

[TABLE]

The right-hand side of this equation is the magnetic susceptibility of the Ising model which, due to universality, can be identified with the chiral susceptibility of QCD. The $C_{2}$, however, is related to the baryon number susceptibility which is known to diverge weaker than the chiral one [19,​ 12,​ 13]. Therefore, the model introduced in Ref. [10] requires some modifications [11]. This can be done using the following relation obtained within the effective model calculations on the mean field level [12,​ 13],

[TABLE]

in which the singular contribution to the baryon number susceptibility is proportional to the chiral susceptibility times the squared order parameter and $\chi_{\mu\mu}^{reg}$ is the regular part of the baryon number susceptibility. To obtain such a form of the second cumulant, the proton mass in Eq. (2) should be replaced by the order parameter, $\sigma$, such that the new $C_{2}$ reads

[TABLE]

The modified higher order cumulants are

[TABLE]

and

[TABLE]

Since the cumulants are volume-dependent it is convenient to consider their ratios in which this dependence cancels out,

[TABLE]

where $M=C_{1}$ is the mean, $\sigma^{2}=C_{2}$ the variance, $\kappa=C_{4}/C_{2}^{2}$ the kurtosis and $S=C_{3}/C_{2}^{3/2}$ the skewness.

To use universality class arguments discussed above, a mapping between the QCD phase diagram and reduced temperature and magnetic field of the spin model is needed. Such a mapping is non-universal and has to be modeled for each system separately. In this work we use a linear mapping [20,​ 21] in which the critical point is located at $r=h=0$, the $r$ axis is tangential to the QCD first order phase transition line and the positive direction of the $h$ axis points towards the hadronic phase. Schematically, this is shown in Fig. 1, where the green line denotes the first order phase transition and the filled band shows lattice QCD constraints on the location of the chiral crossover region.

To calculate the order parameter as well as its cumulants we use the parametric representation of the magnetic equation of state [22]. For a more detailed discussion of the mapping, lattice limits as well as the magnetic equation of state we refer the reader to the papers [10,​ 11].

Finally, assuming that the matter created during a heavy ion collision forms a thermal medium characterized by temperature and chemical potentials, experimental data on event-by-event multiplicity fluctuations can be compared with model results. To this end, we calculate the net-proton number cumulants at the chemical freeze-out. The chemical freeze-out conditions used in this work were obtained by the analysis of hadron yields [23,​ 24,​ 25,​ 26,​ 27,​ 28]. The blue line in Fig. 1 shows the recently obtained parametrization [15].

3 Numerical results

In this section we discuss numerical results on net-proton number cumulant ratios obtained within the current model. The set of model parameters includes the coupling strength $g$ between (anti)protons and the critical mode, the parameters of the magnetic equation of state as well as the size of the critical region in the $(T,\mu)$ plane. Their values as well as a detailed discussion can be found in Refs. [10,​ 11]. Moreover, the location of the QCD critical point is unknown. To study the effect of its position in the QCD phase diagram on the refined model results, we consider three different locations of the CP listed in Tab. 1 and shown in Fig. 2, where the distance to the freeze-out curve is the farthest for CP1 and closest for $CP_{3}$.

The first step of our discussion is the comparison between the $C_{2}/C_{1}$ ratio obtained using the model from Ref. [10], where the n-th net-proton number cumulant is given by Eq. (2), and the refined model [11] for the critical point location CP1. This is shown in Fig. 3. Results obtained using the original model exhibit clear non-monotonic behavior and deviate strongly from the non-critical baseline (the black-dotted line) even for small coupling, $g=3$, which becomes more pronounced for $g=5$ (the red solid and dashed lines, respectively). Using the current approach we find a substantial reduction of criticality in $C_{2}/C_{1}$ ratio, even for larger values of $g$ (see the blue curves in Fig. 3). The refined model results for the $C_{2}/C_{1}$ ratio agree with the experimental data from the STAR Collaboration [5]. On the other hand, the original model would require exceptionally small coupling strength in order to capture the experimentally observed behavior.

The net-proton cumulant ratios obtained in the refined model for different locations of the critical point (as listed in Tab. 1) and with a fixed value of coupling, $g=5$, are shown in Fig. 4. We find that a non-monotonic behavior of cumulant ratios becomes more pronounced when the critical point is closer to the freeze-out line. Moreover, the deviation from the non-critical HRG baseline becomes larger for higher order cumulant ratios.

Finally, Fig. 5 shows the coupling strength dependence of net-proton number cumulant ratios obtained for CP3. We find a strong $g$ dependence of all ratios which is expected since in our refined model the n-th cumulant scales as $g^{2n}$, according to Eqs. (7) and (11)-(13). When our model results are compared to the STAR data [5], we find a qualitative agreement with the $C_{2}/C_{1}$ and $C_{4}/C_{2}$ ratios. On the other hand, the $C_{3}/C_{2}$ ratio does not follow the systematics seen in the data, i.e. our model results overshoot the HRG baseline while the data stay below.

Our results suggest that the appropriate choice of model parameters as well as the location of the QCD critical point allows us to describe some of the experimentally observed cumulant ratios. Especially, the smooth dependence of $C_{2}/C_{1}$ and strong increase of $C_{4}/C_{2}$ at low beam energies, $\sqrt{s}<20\,$GeV, seen by the STAR Collaboration, suggest that the QCD critical point may be located close to the freeze-out curve. However, in this case the $C_{3}/C_{2}$ ratio should increase beyond the non-critical baseline, which is not seen in the experimental data. Therefore it seems unlikely that the QCD critical point is close to the freeze-out curve. This conclusion, however, requires additional theoretical and experimental justifications due to uncertainties in the model parameters as well as in the experimental data.

4 Conclusions

We presented ratios of net-proton number cumulants obtained within an effective model in which the coupling between (anti)protons and critical mode fluctuations is introduced by connecting particle masses to the order parameter. We modified the existing approach [10] to take into account the correct scaling properties of the baryon number susceptibility as dictated by the universality hypothesis.

Model results were compared with the recent experimental data on net-proton number fluctuations from the STAR Collaboration. We find a substantial reduction of the signal coming from the presence of the QCD critical point in the $C_{2}/C_{1}$ ratio which stays in agreement with the experimental data. Moreover, we find that the model discussed in the present work allows us to describe some of the experimentally observed features in the net-proton number cumulant ratios. Particularly, smooth dependence of $C_{2}/C_{1}$ and increase of $C_{4}/C_{2}$ at lower beam energies ($\sqrt{s}<20\,$GeV) suggest that the critical point may be located close to the freeze-out curve. However, the experimentally observed $C_{3}/C_{2}$ ratio does not follow the behavior expected from such a scenario.

Therefore, it seems unlikely that the QCD critical point is located close to the phenomenological freeze-out curve. However, because of uncertainties on both theoretical as well as experimental sides, this statement requires further investigation.

*M. S. and M. B. acknowledge support of the Polish National Science Center (NCN) under grant Polonez UMO-2016/21/P/ST2/04035. M. B. was partly supported by the program "Etoiles montantes en Pays de la Loire 2017". K. R. and C. S. were supported in part by the Polish National Science Center (NCN) under Maestro grant DEC-2013/10/A/ST2/00106. K. R. also acknowledges partial support of the Polish Ministry of Science and Higher Education. We thank F. Geurts and J. Thäder for providing the preliminary STAR data [5] shown in Figures 3 and 5. *