Introduction

Matter comprises quarks and gluons, which are the most fundamental constituents of nature, together with leptons and gauge bosons. The interactions between quarks and gluons are governed by quantum chromodynamics (QCD). At low temperatures or matter densities, quarks and gluons are confined in hadrons, whereas at high temperatures or matter densities, they are deconfined in an extended volume called the quark–gluon plasma (QGP). The phase transition at low temperature and high matter density is first-order, and at high temperature and low matter density, it is a smooth crossover, as predicted by lattice QCD [1]. It has been conjectured that a critical point (CP) exists in the nuclear matter phase diagram of temperature versus matter density between the first-order phase transition and smooth crossover [2-6]. The correlation length increases dramatically near the CP, causing large fluctuations in conserved quantities such as the net baryon number [7]. Searching for the subject of active research in heavy-ion collisions [8-17].

Light nuclei production is well-modeled by nucleon coalescence [18-25]. The coalescence model predicts that large baryon number fluctuations affect the production rate of light nuclei [26, 27]. For example, the production of tritons (t) is enhanced relative to that of deuterons (d) when there are extra fluctuations in the neutron density because a triton contains two neutrons (n), whereas a deuteron contains only one. As a result, the compound ratio $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ involving the multiplicities of protons ($N_{\text{p}}$), deuterons ($N_{\text{d}}$), and tritons ($N_{\text{t}}$) was enhanced. Similarly, the production of helium-3 (³He) and the ratio of $N_{{}^{3}he}{N}_{\text{n}}/N_{\text{d}}^2$ with respect to extra fluctuations in the proton density.

Following Ref. [26-28], the average deuteron multiplicity density in the coalescence model is given by $\begin{array}{c}{\overline{n}}_{\text{d}}/\left(\frac{3}{\sqrt{2}}A\right)=\langle \left({\overline{n}}_{\text{p}}+\delta n_{\text{p}}\right)\left({\overline{n}}_{\text{n}}+\delta n_{\text{n}}\right)\rangle \\ ={\overline{n}}_{\text{p}}{\overline{n}}_{\text{n}}+\langle \delta n_{\text{p}}\delta n_{\text{n}}\rangle \\ ={\overline{n}}_{\text{p}}{\overline{n}}_{\text{n}}\left(1+{\left(\alpha \text{Δ}{n}_{\text{n}}\right)}_{\text{d}}\right)\mathrm{,}\end{array}$ (1) where the protons and neutrons are assumed to be in thermal equilibrium with an effective temperature T_eff and $A={\left(\frac{2\pi }{m_N{T}_{\text{eff}}}\right)}^{3/2}$ is the shorthand notation (mN is the nucleon mass). Here, $\overline{n}$ denotes the average density. The neutron density fluctuations are denoted as $\text{Δ}{n}_{\text{n}}=\langle {\left(\delta n_{\text{n}}\right)}^2\rangle /{\overline{n}}_{\text{n}}^2\mathrm{,}$ (2) and $\alpha =\frac{\langle \delta n_{\text{p}}{n}_{\text{n}}\rangle }{{\overline{n}}_{\text{p}}{\overline{n}}_{\text{n}}}/\text{Δ}{n}_{\text{n}}=\frac{{\overline{n}}_{\text{n}}}{{\overline{n}}_{\text{p}}}\frac{\langle \delta n_{\text{p}}{n}_{\text{n}}\rangle }{\langle {\left(\delta n_{\text{n}}\right)}^2\rangle }$ denotes the correlations between proton and neutron number fluctuations. If the proton and neutron numbers fluctuate independently, α=0; if they fluctuate simultaneously, α=1. It is noteworthy that we have thus far neglected the details of deuteron formation, which is determined by the deuteron wavefunction and is often implemented by the Wigner function formalism (see below). The fluctuations affecting deuteron formation are those within the typical volumes of deuterons. We note this in Eq. (1) by the subscript ‘d’ in ${\left(\alpha \text{Δ}{n}_{\text{n}}\right)}_{\text{d}}$ to indicate that it is the average αΔn_n within the typical deuteron size that is relevant.

Similarly, the triton average multiplicity density is given by $\begin{array}{c}{\overline{n}}_{\text{t}}/\left(\frac{3\sqrt{3}}{4}{A}^2\right)=\langle \left({\overline{n}}_{\text{p}}+\delta n_{\text{p}}\right){\left({\overline{n}}_{\text{n}}+\delta n_{\text{n}}\right)}^2\rangle \\ ={\overline{n}}_{\text{p}}{\overline{n}}_{\text{n}}^2\left(1+{\left(\text{Δ}{n}_{\text{n}}\right)}_{\text{t}}+2{\left(\alpha \text{Δ}{n}_{\text{n}}\right)}_{\text{t}}+{\left(\beta \text{Δ}{n^{\prime }}_{\text{n}}\right)}_{\text{t}}\right)\mathrm{,}\end{array}$ (3) where $\text{Δ}{n^{\prime }}_{\text{n}}=\langle {\left(\delta n_{\text{n}}\right)}^3\rangle /{\overline{n}}_{\text{n}}^3$ and $\beta =\frac{\langle \delta n_{\text{p}}{\left(\delta n_{\text{n}}\right)}^2\rangle }{{\overline{n}}_{\text{p}}{\overline{n}}_{\text{n}}^2}/\text{Δ}{n^{\prime }}_{\text{n}}=\frac{{\overline{n}}_{\text{n}}}{{\overline{n}}_{\text{p}}}\frac{\langle \delta n_{\text{p}}{\left(\delta n_{\text{n}}\right)}^2\rangle }{\langle {\left(\delta n_{\text{n}}\right)}^3\rangle }$. Similar to α in Eq. (1), the β parameters denote the three-body correlations. When the proton and neutron numbers fluctuate independently, β=0, and when they fluctuate together, β=1. The subscript ‘t’ in Eq. (3) indicates that only fluctuations averaged within the typical volume of the triton size matter. Note that in Eq. (3), we simply write $n_{\text{n}}^2$ for the neutron pair density; however, for an identical particle pair, the multiplicity is N(N-1), so the pair fluctuations of Eq. (2) should be understood as those beyond Poisson fluctuations.

Thus, the compound ratio is given by $\frac{N_{\text{p}}{N}_{\text{t}}}{N_{\text{d}}^2}=\frac{{\overline{n}}_{\text{p}}{\overline{n}}_{\text{t}}}{{\overline{n}}_{\text{d}}^2}=\frac{1+{\left(\text{Δ}{n}_{\text{n}}\right)}_{\text{t}}+2{\left(\alpha \text{Δ}{n}_{\text{n}}\right)}_{\text{t}}+{\left(\beta \text{Δ}{n^{\prime }}_{\text{n}}\right)}_{\text{t}}}{2\sqrt{3}{\left(1+{\left(\alpha \text{Δ}{n}_{\text{n}}\right)}_{\text{d}}\right)}^2}\mathrm{.}$ (4) If one neglects α and β, then Eq. (4) is reduced to $\frac{N_{\text{p}}{N}_{\text{t}}}{N_{\text{d}}^2}\approx \frac{1}{2\sqrt{3}}\left(1+{\left(\text{Δ}{n}_{\text{n}}\right)}_{\text{t}}\right)\mathrm{,}$ (5) as in Ref. [26, 27]. Note that the factor $1/2\sqrt{3}$ originates from the thermal equilibrium assumption of nucleon abundances. Therefore, the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ may be a good measure of neutron density fluctuations, a large value of which can signal the CP.

A unique signature of the CP is the nonmonotonic behavior of the ratio $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ in the beam energy, where a peak of the ratio in a localized region of beam energy can signal large neutron fluctuations and the CP [26, 27]. The STAR experiment at RHIC recently observed non-monotonic behavior of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio in the top 10% of central Au+Au collisions as a function of the nucleon-nucleon center-of-mass energy ($\sqrt{s_{{}_{\text{NN}}}}$) in the Beam Energy Scan (BES) data [29].

A nonmonotonic bump was observed in the ratio localized in the energy region $\sqrt{s_{{}_{\text{NN}}}}=20-30$ GeV. It should be noted that the local bump is prominent in the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio of the extrapolated yields to all transverse momenta (p_T) but not as prominent in the measured fiducial range of $0.5<p_{\text{T}}/A<1.0$ GeV/c and $0.4<p_{\text{T}}/A<1.2$ GeV/c (where A is the mass number corresponding to each light nucleus in the ratio) [29].

This is illustrated in the left panels of Figs. 1 where the STAR-measured $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratios from the total extrapolated yields and two measured p_T/A ranges are reproduced. The “extrapolation factor” (f), that is, the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio of the total light nuclei yields extrapolated to the entire p_T/A range [0–∞) divided by the ratio of those measured within a fiducial p_T/A range, is shown in the right panel of Fig. 1 for the two measured p_T/A ranges. It should be noted that the integrated yields of protons, deuterons, and tritons over the full momentum space in the STAR experiment were extrapolated from slightly different ranges of the scaled transverse momentum p_T/A. The aforementioned ratio for the measured fiducial range was obtained using particle yields within the same p_T/A range as in the STAR experiment. The statistical uncertainties between the fiducial yield of a given p_T spectrum and the extrapolated total yield are correlated, as are the systematic uncertainties. Thus, the statistical and systematic uncertainties are considered to be the quadratic difference of the corresponding uncertainties in the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ between the fiducial p_T/A and total ranges. As expected from the STAR results [29], the f values peaked in the $\sqrt{s_{{}_{\text{NN}}}}=20-30$ GeV range and were non-monotonic; the non-monotonic effect was of the order of 10%.

Fig. 1Full-screen View

(Color online) (Left) The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio in 0-10% central Au+Au collisions as functions of beam energy $\sqrt{s_{{}_{\text{NN}}}}$ measured by STAR [29]. The ratio of the extrapolated total yields and those from two measured p_T/A ranges are shown. The ratios from the measured p_T/A ranges are shifted in the horizontal axis for clarity. (Right) The extrapolation factor f, i.e., the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio from the extrapolated total yields divided by that from the fiducial yields in a given measured p_T/A range, is shown for the two measured p_T/A ranges as functions of $\sqrt{s_{{}_{\text{NN}}}}$. The statistical and systematic uncertainties on f are the quadratic difference of the corresponding uncertainties of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio between the given measured p_T/A range and the total range. The rightmost square is shifted in the horizontal axis for clarity

In this study, we use a toy model to generate nucleons and form d, t, and ³He by using a coalescence model. We study the ratios of the light nuclei yields as functions of the magnitude of neutron multiplicity fluctuations and examine the role of p_T acceptance in these ratios. First, we confined ourselves to a simple toy model to gain insight. We then attempted to conduct more realistic simulations mimicking STAR BES data to examine what STAR measurements may entail.

Coalescence model

The probability of forming a composite particle from particle 1 at position ${\overrightarrow{r}}_1$ with momentum ${\overrightarrow{p}}_1$ and particle 2 at position ${\overrightarrow{r}}_2$ with momentum ${\overrightarrow{p}}_2$ is estimated using the Wigner function: $W\left({\overrightarrow{r}}_1\mathrm{,}{\overrightarrow{p}}_1\mathrm{,}{\overrightarrow{r}}_2\mathrm{,}{\overrightarrow{p}}_2\right)=g\cdot 8\exp \left(-\frac{r_{12}^2}{{\sigma }_{\text{r}}^2}-\frac{p_{12}^2}{{\sigma }_{\text{p}}^2}\right)\mathrm{,}$ (6) where $\begin{array}{l}{r}_{12}=|{\overrightarrow{r}}_1-{\overrightarrow{r}}_2|\mathrm{,}\\ p_{12}=\mu \left|\frac{{\overrightarrow{p}}_1}{m_1}-\frac{{\overrightarrow{p}}_2}{m_2}\right|\mathrm{,}\end{array}$ and $\mu =\frac{m_1{m}_2}{m_1+m_2}$ denote the reduced mass of a two-body system. The parameter σ_r is the characteristic coalescence size in the configuration space and σ_p=1/σ_r (where $\hslash =c=1$) is that in the momentum space.

Based on the Wigner function expressed in Eq. (6), the root mean square (RMS) radius of the coalesced composite particle is calculated as $R=\frac{\sqrt{3m_1{m}_2/2}}{m_1+m_2}{\sigma }_{\text{r}}$. Therefore, the coalescence parameter σ_r can be determined from the particle size asfollows: ${\sigma }_r=\frac{m_1+m_2}{\sqrt{3m_1{m}_2/2}}R\mathrm{.}$ (7) Deuterons coalesce with protons and neutrons. The RMS size of the deuteron is R_d=1.96 fm [30], so for deuteron ${\sigma }_{\text{r}}=\sqrt{\frac{8}{3}}{R}_{\text{d}}=3.20$ fm. The triton (helium-3) is formed by the coalescence of a deuteron and neutron (proton), following the same prescription as Eq. (6). The Triton RMS size is R_t=1.59 fm [30]. Therefore, for triton ${\sigma }_{\text{r}}=\sqrt{3}{R}_{\text{t}}=2.75$ fm. For ³He we also assume σ_r=2.75 fm to be symmetric.

The gfactor represents the probability of the proper total spin of the composite particle. $g=(2S+1)/{\Pi }_{i=1}^2\left(2s_i+1\right)\mathrm{,}$ (8) where S is the spin of the composite particle (1 for d, 1/2 for both t and ³He) and si is the spin of each coalescing particle. For the d formed from p and n, the g-factor was 3/4. For t (³He), formed from d and n (p), it is 1/3.

Toy-Model Simulation Details

A system of nucleons was generated using a toy model. The nucleons were assumed to be uniformly distributed in a cylinder of 10 fm length and 10 fm radius. The transverse momentum spectra of the nucleons were assumed to be $\text{d}N/\text{d}{p}_{\text{T}}\propto p_{\text{T}}\exp \left(-p_{\text{T}}/T\right)\mathrm{,}$ (9) where the “temperature” parameter is set to T=150 MeV. Note that Eq. (9) is not a thermal distribution for massive particles; we use the term “temperature” for convenience. The azimuthal angle of the momentum vector was uniformly distributed between 0 and 2π. The pseudorapidity η was assumed to be uniform between -1 and 1. The use of rapidity or pseudorapidity was insignificant.

In this study, we did not use a thermal model to predict the average multiplicity; rather, we set the average multiplicities of protons and neutrons to ${\overline{N}}_{\text{p}}={\overline{N}}_{\text{n}}=20$. This is simple, because our goal was to study only the effects of neutron multiplicity fluctuations. Therefore, the baseline of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio is not $1/2\sqrt{3}$ but is rather determined by the ratio of the degeneracy g factors as $\frac{1}{3}\cdot \frac{3}{4}/{\left(\frac{3}{4}\right)}^2=\frac{4}{9}$ [31].

We assigned Poisson fluctuations to the number of protons and only varied the fluctuation magnitude for the number of neutrons (we focused on the compound ratio $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$). The latter is achieved by using negative binomial distributions $P\left(N_{\text{n}}=k\right)=C_k^{k+r-1}{\left(1-p\right)}^k{p}^r$, where r and p are free parameters.

The mean and variance are ${\overline{N}}_{\text{n}}=\frac{r(1-p)}{p}$ and ${\sigma }_{N_{\text{n}}}^2=\frac{r(1-p)}{p^2}$, respectively. The fluctuation magnitude is given by ${\sigma }_{N_{\text{n}}}^2={\overline{N}}_{\text{n}}/p\equiv \theta {\overline{N}}_{\text{n}}$ (where $\theta \equiv 1/p\ge 1$), which is always larger than Poisson fluctuations unless the probability p=1 when the negative binomial distribution is reduced to Poisson. We used p to control the magnitude of the fluctuations in $N_{\text{n}}$ and selected the proper r value to obtain the desired average neutron multiplicity, ${\overline{N}}_{\text{n}}$.

It is worth noting that we utilized fluctuations in the total number of neutrons $N_{\text{n}}$ event-by-event to mimic fluctuations in the local neutron number density. The purpose of our study is to investigate the behavior of $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ as a function of the neutron density fluctuation magnitude but not to suggest that density fluctuations are caused by fluctuations in the total multiplicity. In our simulation, given a total $N_{\text{n}}$ in an event, the neutrons were randomly distributed in the cylinder volume. The fluctuations in $N_{\text{n}}$ determine the magnitude of the neutron density fluctuations, which is averaged over the entire volume. In other words, Δn_n in Eq. (2) can be quantified by fluctuations in $N_{\text{n}}$ beyond Poisson, $\text{Δ}{n}_{\text{n}}=\left({\sigma }_{N_{\text{n}}}^2/{\overline{N}}_{\text{n}}-1\right)/{\overline{N}}_{\text{n}}\equiv (\theta -1)/{\overline{N}}_{\text{n}}\mathrm{.}$ (10) Note that θ quantifies the fluctuations in $N_{\text{n}}$ in the unit of Poisson fluctuations, and $\theta -1$ quantifies those beyond Poisson fluctuations.

In total, 1.6 × 10⁸ events were simulated. Deuterons, tritons, and helium-3 were formed by coalescence. For each event, the deuteron formation via double loops over protons and neutrons was considered. For each deuteron, the formation of t and ³He was implemented by looping over the remaining nucleons. The nucleons were randomly reordered such that the probabilities of forming t and ³He were unbiased. A random number uniformly distributed between 0 and 1 was used to determine whether the two particles coalesced into a light nucleus. If the random number is smaller than the Wigner function value in consideration, the corresponding light nucleus is formed. Once a light nucleus formed, the coalescing particles were removed from further consideration.

Toy-Model Simulation Results and Discussions

Figure 2 shows the p_T spectra of the generated neutrons and the coalesced d, t, and ³He in the left panel, and the right panel those spectra in the scaled p_T/A. No extra fluctuations are included beyond Poisson in this figure (i.e., θ=1); therefore, the proton spectrum is identical to that of neutrons. In the right panel, the product of the proton and neutron spectra and the product of the proton and squared neutron spectra are shown as smooth histograms. These products are the corresponding d and t spectra if the coalescence parameters σp=0 and σr=∞ in Eq. (6). The coalesced d and t spectra are steeper than those of the products because of the finite σ_p and σ_r implemented in our coalescence model.

Fig. 2Full-screen View

(Color online) (Left) Transverse momentum p_T spectra of deuteron, ³He and triton calculated by the coalescence model from a system of average ${\overline{N}}_{\text{p}}=20$ protons and ${\overline{N}}_{\text{n}}=20$ neutrons at “temperature” T=150 MeV (Eq. (9)) randomly distributed within a cylinder of 10 fm radius and 10 fm length. No extra fluctuations are included beyond Poisson, i.e., θ=1. (Right) The same spectra plotted as functions of p_T/A (A is the corresponding mass number). Superimposed in curves are the products of $\left(\text{d}{N}_{\text{p}}/\text{d}{p}_{\text{T}}\right)\times \left(\text{d}{N}_{\text{n}}/\text{d}{p}_{\text{T}}\right)$ and $\left(\text{d}{N}_{\text{p}}/\text{d}{p}_{\text{T}}\right)\times {\left(\text{d}{N}_{\text{n}}/\text{d}{p}_{\text{T}}\right)}^2$ at the same p_T/A value, arbitrarily scaled to compare to the shapes of the deuteron and triton spectra, respectively

Figure 3 shows in the left panel the yield ratios of $N_{\text{d}}/N_{\text{p}}$ and $N_{\text{t}}/N_{\text{d}}$ as functions of the neutron fluctuation magnitude Δn_n in Eq. (10). The $N_{\text{t}}/N_{\text{d}}$ ratio increases with Δn_n, which is consistent with the expected stronger effect of the neutron density fluctuations on $\text{t}$ than on d production. The $N_{\text{d}}/N_{\text{p}}$ ratio decreases slightly, but is statistically significant. This is counterintuitive, as one would expect no dependence because ${\overline{N}}_{\text{n}}$ is the same for all values of Δn_n; events with more neutrons would be balanced out by those with fewer neutrons in terms of deuteron production. However, ${\overline{N}}_{\text{p}}=20$ is fixed and the fluctuations in $N_{\text{p}}$ are treated as uncorrelated with those in $N_{\text{n}}$ in our study. The production of deuterons will be “saturated” in events with large $N_{\text{n}}$–those neutrons would not get “equal” share of protons to form deuterons. Consequently, $N_{\text{d}}/N_{\text{p}}$ is smaller for a larger Δn_n.

Fig. 3Full-screen View

(Color online) Yield ratios of $N_{\text{d}}/N_{\text{p}}$ and $N_{\text{t}}/N_{\text{d}}$ (left) and $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ (right) as functions of Δn_n of Eq. (10), the magnitude of neutron multiplicity fluctuations beyond Poisson. Numbers inside the parentheses indicate the uncertainty to the corresponding last digit. The light nuclei are formed by coalescence from a system of average ${\overline{N}}_{\text{p}}=20$ protons and ${\overline{N}}_{\text{n}}=20$ neutrons at T=150 MeV (Eq. (9)) randomly distributed within a cylinder of 10 fm radius and 10 fm length

The right panel of Fig. 3 shows the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as a function of Δn_n. The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ clearly increases with Δn_n. The slope was approximately 0.477. were significantly nonzero, and the intercept was approximately 0.439. These values are approximately equal to the expected degeneracy factor of 4/9. The slight deviations may be due to the different size parameters of the deuterons and triton as the Wigner function parameters differ.

As previously mentioned, the measured nonmonotonic feature of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as a function of the beam energy $\sqrt{s_{{}_{\text{NN}}}}$ by STAR appears to depend on the considered p_T/A range [29], which is stronger for the extrapolated yield ratio than for those with limited p_T/A acceptance. Therefore, it is important to investigate whether a nonmonotonic extrapolation factor can result from trivial physics, such as nucleon spectral shape changes as a function of the beam energy. To this end, we first examined the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ as a function of p_T/A in our toy model simulation. This is shown in Fig. 4 left panel shows the various θ values. In the extreme case of coalescence with an identical momentum, the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ was independent of p_T. The falling and rising characteristics of the shape are determined by the p_T-distribution in Eq. (9), and the coalescence parameters of the light nuclei. The shapes are similar for various fluctuation magnitudes, and the overall ratio increases with θ as expected. The right panel of Fig. 4 shows extrapolation factor f as a function of Δn_n. (Note that the total yields in simulation are known of course, not from extrapolation, but we keep use of the term “extrapolation factor”.) Two fiducial ranges are depicted: p_T/A=0.1–0.2 and 0.3–0.4 GeV/c. The extrapolation factor f varies with Δn_n; however, the variation is relatively small, less than 1% for the p_T/A=0.3–0.4 GeV/c range for θ=2, which is the fluctuation magnitude of a factor of 2 of that of Poisson. This implies that even when there are enhanced fluctuations in a local region in the beam energy $\sqrt{s_{{}_{\text{NN}}}}$, the extrapolation factor remains approximately the same and does not cause a change in the shape of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio vs. $\sqrt{s_{{}_{\text{NN}}}}$ from the limited p_T range to the extrapolated full p_T range. However, we examined only relatively narrow p_T/A ranges at small p_T/A values because of statistical considerations. The fiducial p_T/A ranges in the STAR analysis are relatively large and wide. Therefore, no direct comparisons can be made between our study and the STAR results thus far.

Fig. 4Full-screen View

(Color online)(Left) The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as functions of p_T/A. (Right) The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ extrapolation factor f from p_T/A=0.1–0.2 GeV/c and 0.3–0.4 GeV/c as a function of Δn_n. The light nuclei are formed by coalescence from a system of average ${\overline{N}}_{\text{p}}=20$ protons and ${\overline{N}}_{\text{n}}=20$ neutrons at T=150 MeV (Eq. (9)) randomly distributed within a cylinder of 10 fm radius and 10 fm length

More Realistic Simulations

In heavy-ion collisions, p_T distributions can be significantly altered by the collective radial flow. The shape of the p_T distribution affects the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as a function of p_T/A. To investigate the effect of the nucleon p_T spectral shape on the extrapolation factor of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio, we first vary the temperature parameter T in Eq. (9). Figure 5 shows the T-dependence of the extrapolation factotr f for p_T/A=0.3–0.4 GeV/c, with a given neutron fluctuation of Δn_n=0.1 (i.e., θ=3). As expected, the f value varies with T as expected, and the variation is monotonic. Thus, the T-dependent p_T spectral shape does not cause an artificial non-monotonic $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ enhancement from limited to full p_T acceptance in any given localized T range, and correspondingly, the beam energy range.

Fig. 5Full-screen View

(Color online) The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ “extrapolation factor” for the p_T/A=0.3–0.4 GeV/c acceptance as a function of T, with θ=3 or Δn_n=0.1. The light nuclei are formed by coalescence from a system of average ${\overline{N}}_{\text{p}}=20$ protons and ${\overline{N}}_{\text{n}}=20$ neutrons randomly distributed within a cylinder of 10 fm radius and 10 fm length

The collective radial flow in heavy ion collisions creates a correlation between the momentum of a particle and its freeze-out radial position. This correlation is not present in the above simulations, but can be important for coalescence. To fully comprehend the implications of STAR data [29], we performed simulations using more realistic kinematic distributions. We obtained freeze-out information from the measured data [32, 33] as described below. The other simulation details are the same as those described in Sect. 3.

The charged hadron mulitplicity N_ch, the inclusive proton multiplicity density $\text{d}{N}_{\text{p}}/\text{d}y$, the chemical freeze-out temperature T_chem and volume ($V_{\text{chem}}\equiv \frac{4\pi }{3}{R}_{\text{sphere}}^3$), the kinematic freeze-out temperature T_kin and the average collective radial flow velocity $\langle \beta \rangle$ in the Blast-Wave parameterization have been reported for the STAR experiment using various beam energies [32, 33]. We concentrate only on the central 0-5% collisions and parameterize the freeze-out quantities as functions of N_ch. The parameters have fit uncertainties; however, because we were only interested in the input values in our simulation to study light nuclei coalescence, we simply used the parameterized values.

• The measured $\text{d}{N}_{\text{p}}/\text{d}y$ is shown in Fig. 6 (upper left) as a function of N_ch. The measured protons originate from two sources: transport protons whose abundance decreases with energy, and produced protons whose abundance increases with energy. We thus parameterize $\text{d}{N}_{\text{p}}/\text{d}y$ by combination of two functions, one decreasing and the other increasing with N_ch. The parameterization is given by $\text{d}{N}_{\text{p}}/\text{d}y=318e^{-N_{\text{ch}}/148}+68N_{\text{ch}}^{1.4}$, as superimposed.

Fig. 6Full-screen View

(Color online) The measured proton multiplicity density $\text{d}{N}_{\text{p}}/\text{d}y$ (upper left), the chemical freeze-out temperature T_chem (upper right), the kinetic freeze-out temperature T_kin (lower left), and the average radial velocity $\langle \beta \rangle$ (lower right) as functions of the charged hadron multiplicity density $\text{d}{N}_{\text{ch}}/\text{d}y$ [32, 33]. The parameterizations to the data are superimposed and used as inputs to our simulations

• The chemical freeze-out temperature T_chem is shown in Fig. 6 (upper right) as a function of N_ch. The T_chem increases with N_ch and then saturates, so we parameterize it by $T_{\text{chem}}=167\left(1-e^{-N_{\text{ch}}/144}\right)$ MeV.

• The effective sphere radius for the chemical freeze-out volume appears linear in $N_{\text{ch}}^{\text{1/3}}$, so we parameterize it as $R_{\text{sphere}}=3.1+0.42N_{\text{ch}}^{1/3}$. We allow the intercept to be a free parameter because the fit is otherwise not as good.

• The kinetic freeze-out temperature T_kin is shown in Fig. 6 (lower left) as a function of N_ch. It is found to decrease with N_ch, and parameterized by $T_{\text{kin}}=142e^{-N_{\text{ch}}/1827}$ MeV.

• The average radial flow velocity is shown in Fig. 6 (lower right) as a function of N_ch. It increases with N_ch and appears to saturate at large N_ch, so we parameterize it as $\langle \beta \rangle =0.6\left(1-e^{-N_{\text{ch}}/241}\right)$.

The information on the nucleons to be input into the coalescence model is that of kinetic freeze-out. The kinetic freeze-out volume V_kin is obtained by assuming the system expands adiabatically, so $V_{\text{kin}}/V_{\text{chem}}={\left(T_{\text{chem}}/T_{\text{kin}}\right)}^{3/2}$. In our simulation, we assumed that the collision zone was a cylinder with a radius R=7 fm, and the length of the cylinder was determined by $V_{\text{kin}}/\left(\pi R^2\right)$. Protons and neutrons were positioned uniformly and randomly inside the cylinder. The mean numbers of protons and neutrons were assumed to be ${\overline{N}}_{\text{p}}={\overline{N}}_{\text{n}}=2\times \text{d}{N}_{\text{p}}/\text{d}y$ over two units of speed. Protons and neutrons were first generated according to the thermal distribution at T_kin. $\text{d}N/\text{d}{p}_{\text{T}}\propto p_{\text{T}}{m}_{\text{T}}{e}^{-m_{\text{T}}/T}\mathrm{.}$ (11) The generated protons and neutrons are boosted radially with a boost velocity dependent on the radial position of the nucleon (r) within the cylinder: $\beta ={\beta }_{\text{s}}{\left(r/R\right)}^n\mathrm{,}$ (12) where the surface velocity is ${\beta }_{\text{s}}=(1+n/2)\langle \beta \rangle$. In this study, the parameter n for all the collision energies was set to 1.

Figure 7 shows the calculated p_T/A spectra of triton as an example for selected beam energies for θ=1 and θ=10. The p_T/A spectra had the same shape for θ=1 and θ=10 for a given beam energy. The right panel shows the yield fractions of the protons (red), deuterons (blue), triton (green) in the limited p_T range of $0.5<p_{\text{T}}/A<1.0$ GeV/c. It is interesting to note that while the proton fiducial yield fraction steadily decreases with the beam energy, as expected from the flattening of the proton p_T spectra, the deuteron and triton fiducial yield fractions first increase with the energy and then decrease. This is because the yields in the low p_T/A<0.5 GeV/c region are disproportionately more dominant at lower energies for the coalesced light nuclei, and more so for triton. However, no differences were observed in the yield fractions for θ=1 and θ=10. This implies that the yield extrapolations are the same for different fluctuation magnitudes. This, in turn, suggests that the effect of extrapolation on the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio is the same, independent of the fluctuation magnitude, as illustrated below.

Fig. 7Full-screen View

(Color Online) (Left panel) The p_T/A spectra of triton for different beam energies. Two values of θ are shown: θ=1 for Poisson (hollow markers) and θ=10 for enhanced neutron fluctuations (filled markers). (Right panel) The yield fractions of proton (red), deuteron (blue), triton (green) within $0.5<p_{\text{T}}/A<1.0$ GeV/c as functions of $\sqrt{s_{{}_{\text{NN}}}}$ with θ=1 (open marker) and θ=10 (filled marker)

Figure 8 left panel shows the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ calculated using our coalescence model as a function of $\sqrt{s_{{}_{\text{NN}}}}$. The ratios are shown for the entire p_T range and for the two limited p_T/A ranges measured by STAR: $0.4<p_{\text{T}}/A<1.2$ GeV/c and $0.5<p_{\text{T}}/A<1.0$ GeV/c. The results for two values of θ are depicted: θ=1 corresponds to Poisson fluctuations in the neutron multiplicity, and θ=10 corresponds to enhanced fluctuations with a magnitude corresponding to ~10% increase in the measured $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio by STAR (cf Eq. (10), where a typical ${\overline{N}}_{\text{n}}\sim 60-100$, as shown in Fig. 6). Similar to the STAR measurements [29], the calculated $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratios show weak energy dependence. The calculated $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ values were lower than the measured values.

Fig. 8Full-screen View

(Color online) (Left panel) The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as functions of $\sqrt{s_{{}_{\text{NN}}}}$ for the entire p_T range (red) and for two limited p_T/A ranges (blue and green). Two values of θ are shown: θ=1 for Poisson (hollow markers) and θ=10 for enhanced neutron fluctuations (filled markers). (Right panel) The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ extrapolation factor for the p_T/A=0.4–1.2 GeV/c (blue) and 0.5–1.0 GeV/c (green) acceptance as functions of $\sqrt{s_{{}_{\text{NN}}}}$, with θ=1 (hollow markers) and 10 (solid markers)

The right panel of Fig. 8 show the extrapolation factors f of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ from p_T/A=0.4–1.2 GeV/c and 0.5–1.0 GeV/c, respectively. The extrapolation increases monotonically with an increase $\sqrt{s_{{}_{\text{NN}}}}$ and is smooth. For a given p_T/A range, the f values were almost identical for the two θ values; thus, we is so for all θ values we have simulated. The chosen θ=10 corresponds to Δn_n = 0.1 for ${\overline{N}}_{\text{n}}=90$ only, which corresponds to different values of Δn_n at various beam energies. We checked our results for f using a fixed Δn_n = 0.1 value at all beam energies studied. No qualitative differences were observed between groups. Our results suggest that p_T-extrapolation alone does not create, enhance, or reduce a local peak in the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as a function of $\sqrt{s_{{}_{\text{NN}}}}$. In other words, the fluctuation increase in the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as a function of $\sqrt{s_{{}_{\text{NN}}}}$~20–30 GeV region implied by the STAR-extrapolated $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio [29] (such as a jump from the hollow red points to the solid red points at $\sqrt{s_{{}_{\text{NN}}}}=20-30$ GeV in our simulation) should also be present in the ratios in the limited fiducial p_T/A regions (a jump from the hollow blue/green points to the solid blue/green points). Therefore, we postulate that the local peak in the extrapolated $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ in the STAR measurements [29] is unlikely to be caused by the physics of coalescence.

Summary

We simulated light nuclei production from a system of nucleons with varying magnitudes of neutron multiplicity fluctuations, Δn_n. It was found that the light nuclei yield ratio $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ increased linearly with Δn_n, confirming the findings of Ref. [26, 27]. We have further investigated the effect of finite acceptance by studying the “extrapolation factor” f for the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio as a function of Δn_n and the p_T spectra parameter T. The f value is found to be monotonic as a function of T and Δn_n; the variations in the latter were relatively small.

We also conducted a coalescence model study with realistic kinematic distributions of nucleons using measured freeze-out parameters, kinetic freeze-out temperature, and collective radial flow velocity. The $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio was calculated using the coalescence model as a function of beam energy, and a weak beam energy dependence was found, similar to the experimental data. The extrapolation factor f was found to be smooth and monotonic in beam energy and independent of the magnitude of the neutron density fluctuations. We conclude that the extrapolation of the $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ ratio in p_T does not create, enhance, or reduce a local peak in the ratio of the beam energy. Therefore, our study suggests that the measured enhancement of the p_T-extrapolated $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$ [29] is unlikely to be caused by coalescence.

Our study provides a necessary benchmark for light nuclei ratios as a probe for nucleon fluctuations, an important observation in the search for the critical point of nuclear matter. In the present study, we assumed that the fluctuations in the number of protons and neutrons were independent. One may implement varying degrees of correlation between these fluctuations and study their effects on $N_{\text{t}}{N}_{\text{p}}/N_{\text{d}}^2$. In addition, the effects of clustering or clumping [34-36], out-of-equilibrium effects [34] and feedback from excited states [37, 38] have not been considered. We leave these studies for future work.