Searching for the QCD critical point via the rapidity dependence of cumulants

The search for a possible critical point in the QCD phase diagram is ongoing in heavy ion collision experiments at RHIC which scan the phase diagram by scanning the beam energy; a coming upgrade will increase the luminosity and extend the rapidity acceptance of the STAR detector. In fireballs produced in RHIC collisions, the baryon density depends on rapidity. By employing Ising universality together with a phenomenologically motivated freezeout prescription, we show that the resulting rapidity dependence of cumulant observables sensitive to critical fluctuations is distinctive. The dependence of the kurtosis (of the event-by-event distribution of the number of protons) on rapidity near mid-rapidity will change qualitatively if a critical point is passed in the scan. Hence, measuring the rapidity dependence of cumulant observables can enhance the prospect of discovering a critical point, in particular if it lies between two energies in the beam energy scan.

The search for a possible critical point in the QCD phase diagram is ongoing in heavy ion collision experiments at RHIC which scan the phase diagram by scanning the beam energy; a coming upgrade will increase the luminosity and extend the rapidity acceptance of the STAR detector. In fireballs produced in RHIC collisions, the baryon density depends on rapidity. By employing Ising universality together with a phenomenologically motivated freezeout prescription, we show that the resulting rapidity dependence of cumulant observables sensitive to critical fluctuations is distinctive. The dependence of the kurtosis (of the event-by-event distribution of the number of protons) on rapidity near mid-rapidity will change qualitatively if a critical point is passed in the scan. Hence, measuring the rapidity dependence of cumulant observables can enhance the prospect of discovering a critical point, in particular if it lies between two energies in the beam energy scan.
A central goal of heavy ion collision experiments is to map the QCD phase diagram as a function of temperature T and baryon chemical potential µ B [1][2][3]. At zero µ B , the phase diagram features a continuous crossover from quark-gluon plasma (QGP) to ordinary hadronic matter as a function of decreasing T [4][5][6][7][8]. Increasing µ B corresponds to doping the QGP with an excess of quarks over antiquarks, and it is an open question whether the crossover becomes a sharp first order phase transition beyond some critical point [3,9]. At nonzero µ B where lattice calculations become extremely difficult [10,11], there are no first-principles theoretical calculations which provide reliable guidance as to the whether there is a critical point in the phase diagram of QCD, or its location if it does exist [12][13][14][15]. Model calculations suggest the existence of a critical point, but disagree wildly on its location in the (µ B , T ) plane [14,15]. Reducing the beam energy increases the µ B of the QGP produced in a heavy ion collision [3,[16][17][18] (principally because lower energy collisions make less entropy but also because they deposit more of their baryon number in the plasma) but it also reduces the temperatures achieved. So, these experiments can scan the crossover (and potentially critical) regime of the phase diagram out to some value of µ B corresponding to the lowest energy collisions that reach the crossover (critical) temperature [1,2]. If a critical point is located in the regime that is within reach, it may be detected experimentally. The search for a critical point in the phase diagram of QCD at the Relativistic Heavy Ion Collider (RHIC) is currently underway, with collisions at energies ranging from √ s = 200 AGeV down to √ s = 7.7 AGeV, producing fireballs that freeze out with chemical potentials in the range 25 MeV < ∼ µ B < ∼ 400 MeV [17,18]. Phase I of the RHIC beam energy scan (BES-I) was completed in 2014, with no signs of a critical point for µ B < 200 MeV and with tantalizing but inconclusive results at larger µ B , in collisions with 19.6 AGeV ≥ √ s ≥ 7.7 AGeV [1,2,9,[19][20][21]. Phase 2 of the scan (BES-II), to begin in 2019 [9,21], will focus on this regime with increased luminosity and consequently much higher statistics. One of the improvements planned before BES-II is an upgrade of the inner Time Projection Chamber (iTPC) at STAR, which will extend its rapidity acceptance for protons from |y| < 0.5 in BES-I to |y| < 0.8 in BES-II [22].
The energy of a heavy ion collision sets the initial T and µ B of the QGP which is created, with lower energy collisions being more baryon-rich. The QGP then follows a trajectory in the (T, µ B ) plane as it expands and cools. If there is a critical point in the QCD phase diagram within the range of µ B which is accessible in the Beam Energy Scan (BES), then at some collision energies the fireball produced may pass through or near the critical region, while at higher (lower) collision energies the fireball produced will pass the critical point on the low (high) µ B side. The theoretical challenge is to describe the unique signatures of this scenario which would be observable in data from the BES.
A critical point in a thermodynamic system is characterized by an enhanced correlation length. Although the correlation length itself is not observable, because the critical order parameter σ couples to all hadrons the n'th cumulant moments κ n of the event-by-event distribution of the measured multiplicity N of various particle species, for example κ 4 [N ] = (δN ) 4 − (δN ) 2 2 , scale with powers of the correlation length ξ near a critical point [23]. Protons couple more strongly to σ than pions or kaons, making cumulants of the proton multiplicity good observables with which to look for critical fluctuations [24,25]. The non-Gaussian cumulants κ 3 [N ] = (δN ) 3 ∼ ξ 9/2 and κ 4 [N ] ∼ ξ 7 scale with higher powers of the correlation length than the Gaussian cumulants and are therefore more sensitive to critical behavior [23,25]. Furthermore, an analysis that is valid for any critical point in the same (3d Ising) universality class has shown that κ 4 will also change sign near a QCD critical point [26,27]. The dependence of fluctuation measures on the total rapidity acceptance has been studied before [28,29] upon assuming that µ B , and hence ξ and κ 4 , do not depend on rapidity. However, the baryon density does depend non-trivially on rapidity at RHIC BES energies (see e.g. Refs. [30,31]). Since the correlation length and fluctuations become very large near a critical point, the rapidity dependence of the baryon density gives rise to a strong, non-trivial rapidity dependence of the cumulants near a critical point which was not incorporated in previous work. Furthermore, since each rapidity is associated with a different value of µ B and therefore probes a different part of the critical region on the phase diagram, we shall see that integrating over the full rapidity acceptance averages out interesting features in the cumulants which are characteristic of critical behavior. Instead, we propose that binning the cumulants in rapidity gives a more crisp picture of the critical regime, and demonstrate that the rapidity dependence of these binned cumulants near mid-rapidity will change qualitatively if a critical point is passed in the BES. We therefore propose this observable as a complementary means by which to observe the presence of a critical point at the BES. In particular, it provides a new, and distinctive, signature by which to determine whether downward steps in the collision energy take us past a critical point in the phase diagram.
Rapidity-dependence of µ B : In this Letter, we illustrate the effect that the rapidity-dependence of µ B at RHIC BES energies has on the rapidity-dependence of cumulants, and propose using this toward discovering (or ruling out) a critical point using RHIC BES data.
Since the baryon density at freezeout is symmetric in the spacetime rapidity y s for symmetric heavy-ion collisions, for small |y s | the deviation from boost invariance takes the form with µ B,0 and α constants that depend on the beam energy √ s. We shall use this form for illustrative purposes, noting of course that it cannot be relied upon at large |y s |. As we have discussed, the basis of the BES is that downward steps in √ s yield upward steps in µ B,0 . For illustrative purposes, we shall pick three values of µ B,0 within the BES range, and see what happens if these steps were to happen to take us past a possible critical point. The value of α has been measured in SPS collisions with √ s = 17.3 AGeV, where α = 50 MeV [30]. At this (and all higher, and some lower) beam energies, α > 0 because the baryon number density is peaked at roughly two units of rapidity below the beam rapidity, meaning that it is less at y s = 0 than at larger |y s | [3,16]. In AGS collisions with √ s = 5.5 AGeV, though, the beam rapidity is low enough that the baryon number density is peaked at y s = 0, and α < 0 [16]. The √ s at which α changes from positive to negative is not known, but is likely near the lower end of the BES range. Ultimately, measurements of the ratios of the mean particle number distribution for different species from the RHIC BES should be used to measure how µ B at freezeout depends on y s at each BES collision energy, and hence to determine the value of α at each energy. For illustrative purposes here, we shall investigate the consequences of choosing α = 50 MeV at each of our three values of µ B,0 as well as checking how things change if we choose α = −50 MeV instead at our largest value of µ B,0 . Cumulants in the critical regime: Order parameter fluctuations near a critical point induce fluctuations in the event-by-event particle multiplicities. Throughout this work, we will consider the cumulants of protons, as these are expected to be most sensitive to critical fluctuations [25]. From previous work [23,25,26,28], the contribution to the fourth cumulant of the proton multiplicity distribution coming from critical fluctuations (denoted by the subscript σ) takes the form where the x-integral is a spacetime integral over the freezeout hypersurface, where T , µ B and consequently ξ(µ B , T ) and K 4 (µ B , T ) (proportional to the kurtosis of the event-by-event distribution of the fluctuating order parameter, see below) can take on different values at different points on the freezeout hypersurface, where g is the σ-proton-proton coupling which we set to the same benchmark value g = 7 used in Ref. [25], where the pintegral is a momentum-space integral over the protons at the point x, where χ p = f p (1 − f p )/T if we assume local equilibrium, with f p the Fermi-Dirac distribution boosted by the radial flow velocity at the point x, and where γ p = p 2 + m 2 /m, with m the proton mass. We follow Ref. [28] and use a blast wave model to obtain the radial flow velocity and freezeout hypersurface, taking the freezeout curve in the (T, µ B ) plane from the fit to experimental data found in Ref. [17,18]. Following Ref. [28], we shall make the approximation χ p ≈ f p /T and use the Boltzmann distribution for f p , allowing us to do some of the integrals analytically. The shape of the dependence of K 4 and ξ on µ B and T are governed by universal properties of critical fluctuations. A critical point in the QCD phase diagram, if it exists, is known to be in the same universality class as the 3d Ising model [32][33][34][35][36]. The mapping of the Ising variables (r, h) onto the QCD variables (µ B , T ) is not universal, but for illustrative purposes we employ the widely-used assumption [37] that the Ising r axis (and hence the line of first order transitions) is parallel to the QCD µ B axis, and the Ising h axis is parallel to the QCD T axis. For illustrative purposes, we shall place a hypothetical QCD critical point at µ B = 260 MeV, T = 160 MeV. The 3d Ising universality then determines K 4 at some point away from the critical point in terms of the direction in which that point lies in the (µ B , T ) plane, and ξ in terms of this angle, the distance away from the critical point, and one non-universal parameter whose choice determines the contour on the phase diagram where ξ = 1 fm, as illustrated in Fig. 1. Because once ξ is less than 1 fm the magnitude of κ 4 ∝ ξ 7 is negligible, for simplicity we set ξ = 0 outside the critical regime. (For details, see Refs. [25,26,38,39].) Following Ref. [28], we cast the momentum integration in terms of the momentum-space rapidity y and transverse momentum p ⊥ , which are measured in experiment: yc+∆y/2 yc−∆y/2 dy pmax pmin p ⊥ dp ⊥ 2π 0 dψ .
(3) We have introduced a finite acceptance in both rapidity and transverse momentum. We will keep p min = 0.4 GeV and p max = 2 GeV throughout. We shall compute κ 4 using two different kinds of rapidity cuts, either varying ∆y with y c = 0, in which case |y| < y max ≡ ∆y/2, or varying y c with fixed bin width ∆y.
To simplify the interpretation of our results, we shall show the critical contribution to the cumulants normalized by the average number of protons, ω 4,σ ≡ κ 4 [N ] σ / N . This cumulant ratio has the advantage that if the background (noncritical) contribution were Poisson-distributed it would contribute ω 4,σ = 1, meaning that our results in Figures 1 and 2 should be interpreted as critical contributions to be added to a background of order 1.
Results and conclusions: In this Section, we demonstrate that the rapidity-dependence of µ B makes the rapidity dependence of cumulants sensitive to critical fluctuations in a way that yields distinctive, qualitative, observable consequences. In Fig. 1 we first compute the dependence of the cumulant ratio ω 4,σ on the total rapidity acceptance y max . This dependence was studied previously in Ref. [28] upon assuming that µ B itself is constant in rapidity; we find striking consequences of the rapidity-dependence of µ B . Next, motivated by the expanded rapidity coverage that the STAR iTPC upgrade will bring, we compute ω 4,σ for bins in rapidity, something that has not been considered previously. We find that the rapidity dependence of the cumulant ratio is a sensitive and interesting probe of critical behavior.
In Fig. 1 we consider a hypothetical set of scenarios motivated by the possibility that there may be a QCD critical point within the energy range to be explored by the RHIC BES. We imagine a critical point at µ c B = 260 MeV, and in the first three rows of the figure we consider heavy ion collisions with three decreasing values of the beam energy such that freezeout at midrapidity occurs at µ B,0 = 200, 230, 240 MeV. In all three rows, we choose α = 50 MeV, corresponding to the measured value from SPS collisions with √ s = 17.3 GeV and µ B,0 = 237 MeV. Because a real critical point may lie at larger µ c B than this, where α may become negative, in the fourth row we flip the sign of α. The right column of Figure 1 shows ω 4,σ binned in rapidity bins of width ∆y = 0.4 centered around y = ±y c , an observable which to our knowledge has not been considered before. This is a more sensitive observable to the unique features of critical behavior than the dependence on the total rapidity acceptance in the center column because it isolates contributions coming from more similar values of µ B , and the correlation length and other features of the critical regime are sensitive to µ B near µ c B . We see many interesting qualitative features in the rapidity dependence of ω 4,σ . For example, if µ B,0 is in the red region, where ω 4,σ is negative and relatively small in magnitude, larger and positive contributions to ω 4,σ can be found at larger rapidity. This can be seen in the middle panels of the first and second rows, but it is much more striking in the right panels, indicating the value of binning in rapidity.
On the other hand, if µ B,0 lies in the blue region, in the right column the largest value of ω 4,σ is obtained for the bin centered at y = 0, with ω 4,σ decreasing with increasing rapidity while staying positive if α > 0 as in the third row or decreasing with increasing rapidity while becoming negative if α < 0 as in the fourth row. Both the sign change and the non-monotonic behavior in ω 4,σ , as a function of the rapidity acceptance in the center panel of Fig. 1 and even more so as a function of the rapidity bin in the right panel of Fig. 1 are new results of this work. They arise from the rapidity dependence of µ B at freezeout in collisions at RHIC BES energies, and provide distinctive signatures if decreasing the beam energy in this scan takes µ B,0 past a critical point.
To complement Fig. 1, Fig. 2 shows the cumulant ratio ω 4,σ binned in rapidity for a fixed beam energy (fixed µ B,0 ) as the location of the critical point µ c B is changed. There are several features of binning the cumulants in rapidity which we believe will make doing so an important way to probe the critical region, if a critical point is discovered in the RHIC BES. First, ω 4,σ increases with |y c | if freezeout at mid-rapidity occurs at a µ B,0 that is well below µ c B , in the red region, whereas it will decrease with |y c | if µ B,0 is closer to or larger than µ c B , in the blue region. This remains true even if α changes sign, as demonstrated in the bottom row of Fig. 1. Furthermore, a sign change in ω 4,σ as a function of y will be easier to see upon binning in |y c | since not doing so, as in the middle panels, can obscure it by mixing data from different regions in rapidity. Even in cases where the sign change in ω 4,σ is visible in the middle column, as in the first and second rows, it happens at a lower rapidity in the right panel than in the middle panel, making it more feasible to observe at STAR via binning in |y c |.
We conclude that the rapidity-dependence of µ B at  FIG. 1. In the left column we see that we have assumed the existence of a critical point (red dot) at (µ c B , T ) = (260, 160) MeV whose critical region, bounded by the contour where ξ = 1 fm, is colored red and blue. The colors denote the sign of ω4,σ, with ω4,σ > 0 in blue and ω4,σ < 0 in red. Different rows correspond to different assumptions for where on the phase diagram a heavy ion collision freezes out, cf. collisions with varying beam energy. The black circles show where freezeout occurs at mid-rapidity, from top to bottom with µB,0 = 200, 230, 240, 240 MeV. The black dashed curves show how the freezeout conditions change with increasing spacetime rapidity, with the circle, square, and triangle indicating freezeout at ys = 0, 0.6, and 1.2, respectively. In the top three rows, we have chosen α = 50 MeV (see eq. 1) while for the bottom row we have chosen α = −50 MeV. The middle column shows how ω4,σ computed for a rapidity acceptance |y| < ymax depends on ymax. The right column shows how ω4,σ computed in a pair of bins with width ∆y = 0.4 centered at ±yc depends on yc. The results in the middle column sum over a wide range of rapidities (with |y| between 0 and ymax) which freezeout with a range of µB, meaning that features from the left column are more directly visible in the right column than in the middle. In both the center and right columns, the black dotted lines show ω4,σ with α = 0, i.e. what would have been obtained if µB = µB,0, denoted by the black circles in the left column, everywhere. The results shown in the right and middle column should not be taken as quantitative predictions since they depend on the many assumptions that we made for illustrative purposes; they are illustrative of qualitative features to be expected in the rapidity-dependence of cumulants if steps in beam energy take us past a critical point. RHIC energies may result in qualitative signatures of critical fluctuations manifest in the rapidity-dependence of the cumulant ratio ω 4,σ . Complementary to scanning the phase diagram by taking steps in beam energy, the rapidity dependence of µ B provides additional scans of small regions of the phase diagram. We have seen that non-monotonicity and a sign change of the critical contribution to ω 4,σ as a function of rapidity will arise if the BES includes energies on both sides of a critical point. Binning the cumulants in rapidity provides a sensitive probe of these effects. Signatures of critical behavior in the √ s-dependence of ω 4,σ can therefore be crosschecked by looking for qualitative changes in the rapiditydependence of ω 4,σ between nearby beam energies. We have made arbitrary choices at many points, for illustrative purposes. A future quantitative study should include investigation of changes to these choices, as well as an investigation of consequences of variations in the baryon density at a given rapidity, for example along the lines of Ref. [31]. We note however that if the value of α is determined from experimental data as in Ref. [30] this will incorporate the most important such consequence. Future studies should also include an analysis of the quantitative effects of non-equilibrium dynamics on the values of √ s or y c at which the qualitative features we have found occur [37,[40][41][42].