Charge Conservation and Higher Moments of Charge Fluctuations
CCharge Conservation and Higher Moments of Charge Fluctuations
Scott Pratt
Department of Physics and Astronomy and National Superconducting Cyclotron LaboratoryMichigan State University, East Lansing, MI 48824 USA
Rachel Steinhorst
Department of Physics and AstronomyMichigan State University, East Lansing, MI 48824 USA (Dated: August 21, 2020)Higher moments of distributions of net charge and baryon number in heavy-ion collisions havebeen proposed as signals of fundamental QCD phase transitions. In order to better understandbackground processes for these observables, models are presented which enable one to gauge theeffects of local charge conservation, decays of resonances and clusters, Bose symmetrization, andvolume fluctuations. Monte Carlo methods for generating samplings of particles consistent withlocal charge conservation are presented, and are followed by a review of simple analytic modelsinvolving a single type of charge with a constant experimental efficiency. The main model consistsof thermal emission superimposed onto a simple parameterization of collective flow, known as ablast-wave, with emission being consistent with individual canonical ensembles. The spatial extentof local charge conservation is parameterized by the size and extent over which charge is conserved.The sensitivity of third and fourth order moments, skewness and kurtosis, to these parameters,and to beam energy and baryon density is explored. Comparisons with STAR data show that asignificant part of the observed non-Poissonian fluctuations in net-proton fluctuations are explainedby charge and baryon-number conservation, but that measurements of the STAR collaboration forfluctuations of net electric charge significantly differ from expectations of the models presented here. a r X i v : . [ nu c l - t h ] A ug I. INTRODUCTION
The fluctuation of conserved charges is a standard means by which to investigate and classify phase transitions.At the critical point correlation lengths diverge, which results in peaks in charge fluctuations as one approaches thecritical point. For systems with first-order phase transitions, fluctuations turn into phase separation and fluctuationmeasures are no longer intensive quantities. The growth of fluctuations becomes increasingly dramatic as one considersprogressively higher-order fluctuations. In a volume V , fluctuations of a charge Q can be defined as M N ≡ V (cid:104) ( Q − Q ) N (cid:105) = 1 V (cid:88) n P n ( n − n ) N , (1)when particles have unit charge. The measure is increasingly sensitive to the tails of the multiplicity distribution, P n , as n increases. The free energy, F ∼ a ( n − n ) + a ( n − n ) · · · , is minimized for n = n , but the quadratic partvanishes at the critical point, a →
0, which allows the fluctuations to grow and be dominated by higher-order terms.The shapes of the tails of the distribution are then profoundly altered.The properties of the QCD transition, deconfinement and the restoration of chiral symmetry, are not well understoodat finite baryon density. There exists the possibility that this transition is a true phase transition, with a criticalpoint at several times nuclear density and with a critical temperature close to the pion mass. If this is the case, itbegs the question as to whether the conditions for phase separation or for critical phenomena can be reproduced inthe laboratory. Heavy ion collisions at high energy, measured at the Relativistic Heavy Ion Collider (RHIC) or at theLHC, can produce mesoscopic regions at temperatures of a few hundred MeV, which is well above the expectationsfor a critical temperature, and densities of several times nuclear matter density. These densities might or might notbe sufficiently high to investigate the phase transition.High-energy heavy ion collisions are characterized by strong explosive collective flow. Measurement is confined tothe outgoing asymptotic momenta, but because of strong flow, correlations in coordinate space manifest themselvesas correlations in relative momentum. Thus, measurements of correlations binned by relative momentum or chargefluctuations within some defined region of momentum space serve as surrogates for the corresponding observables incoordinate space. Indeed, measurements of charge and baryon number fluctuations have been performed at RHIC.By adjusting the beam energy of the colliding nuclei, experiments at RHIC have explored conditions at which novelphase phenomena might occur. Fluctuations of electric charge and baryon number have been especially popular. Aninitial scan of beam energies [1–21] was rather inconclusive, but measurements with greatly improved statistics arecurrently being undertaken and analyzed. The main thrust of these studies is to search for evidence of a QCD phasetransition with a critical point at finite baryon density, with the hope that the phenomena one expects for an idealizedequilibrated system [22, 23] becomes manifest in the measured debris of heavy-ion collisions. Even if no true phasetransition exists, charge fluctuations, or equivalently susceptibilities, are fundamental properties of the quark gluonplasma and can be investigated with lattice gauge theory [24].In addition to the finite system size (event multiplicities might number in the thousands), the novel states of mattercreated in heavy-ion collisions persist for (cid:46)
10 fm/ c . This severely limits the degree to which phases can separate orto which critical fluctuations can grow. This also limits the extent to which conserved charges can separate from oneanother. For example, if a strange and an anti-strange quanta are produced together their separation is limited bydiffusion, which is difficult to calculate in lattice gauge theory [25, 26], but can be roughly extracted experimentallyand is gaining attention theoretically [27–30]. To justify thermal models of fluctuations of conserved charge based onequilibration, sufficient time is required for particles to enter and exit some defining volume.The first goal of this paper is to gauge the degree to which charge-balance correlations affect higher-order cor-relations. Charge balance functions, which are two-particle correlations related to charge conservation, have beenmeasured extensively [31–41], and modeled theoretically [42–49]. In addition to making it difficult for phases toseparate or for critical correlations to emerge, local charge conservation also represents its own source of correlation,which needs to be understood as a potential source of background before making firm arguments to have observedphenomena related to phase transitions. It is well known that charge-balance correlations are readily measurable atthe two particle level, N = 2 in Eq. (1), and that they explain the bulk of the N = 2 fluctuation measurement.However, their impact on N = 3 , C = M − M , (2)which is based on cumulants, subtracts much of the contribution to M coming from purely two-particle correlations.However, charge conservation can involve multiple particles, and the degree to which a cumulant-based measure, likethe kurtosis, is affected by charge conservation is not fully understood. Relations based on a uniform acceptanceprobability and for a single type of charge were worked out in [50], and provide significant insight into how higher-order correlations are affected by charge conservation. The goal of this study is to extend such ideas to a more realisticpicture, which takes into account the conservation of all three types of charge (baryon number, electric charge andstrangeness) and applies a more realistic model of experimental acceptance and efficiency. The interplay of chargeconservation with chemical equilibrium, decays, and Bose statistics are all considered.To understand the role of chemical equilibrium and decays, a model is presented which creates small volumes inwhich the net charges B, Q and S are each fixed at some value. Even if the net charges are all zero, charged particlesexist in combinations that conserve the net charge. Theoretical methods for exact calculation of the canonical ensembleand a method for Monte Carlo generation of statistically independent events are presented here. The method lendsitself to including decays and accounting for experimental acceptance and efficiency. The physical picture of treatingsmall volumes as independent patches was previously done for calculation of charge balance functions in [51–53], andwas also applied in [54]. Following the terminology in [54], we sometimes refer to these sub-volumes as patches. Themethod presented here creates perfectly independent samplings, and can generate billions of such patches within afew hours. This enables highly accurate calculations of higher moments with minimal numeric cost.After a brief review of cumulants and the definitions of skewness and kurtosis in Sec. II, the method for exactcalculation of the canonical ensemble describing a multi-component, multi-charge hadron gas is presented in Sec.III. These techniques extend those used for canonical ensembles used to study isospin fluctuations of a hadron gas[55], nuclear fragmentation [56], the level density of a Fermi gas [57], and the effect of restricting a quark-gluonplasma (QGP) to having fixed charge, including being in an overall color singlet [58]. Exact methods for calculatingcorrelations up to fourth-order are presented. Unfortunately, when including realistic acceptance effects and complexdecays, the exact expressions are no longer viable. However, as shown in Sec. IV, the exact expressions show howsample events can be generated. Each event, defined by a set of particles and momenta, is generated with perfectindependence from the others, and perfectly reproduces the canonical expressions from Sec. III. Section V extendsthe previous sections to show how Bose correlations can be included. Section VI considers the case of a single type ofcharge and a uniform efficiency, i.e. the probability of any particle being recorded is set to some fixed value. Much ofthis discussion repeats what is said in [50], and is included for completeness. This provides for a physical discussion ofhow charge conservation, volume fluctuations, chemical equilibrium, decays, clustering, and Bose corrections shouldaffect higher moments.The heart of this study is presented in Sec. VII. Here, the patches are assigned collective velocities consistent withthe collective flow deduced from heavy-ion collisions. A canonical sampling of particles is generated from each patch,followed by a simulation of their decays. The particles are then overlaid onto the acceptance of the STAR detector atRHIC. Each patch is uncorrelated with any other patch, so moments can be calculated by averaging over the indepen-dent contributions from the patches. Results are displayed alongside results from the STAR Collaboration. The sizeand sign of the fluctuations of the net-proton distributions are consistent with observations, but the calculations ofthe net charge distributions differs qualitatively from STAR observations. A detailed discussion of the lessons derivedfrom this study is presented in Sec. VIII. II. CUMULANTS, SKEWNESS AND KURTOSIS
For the manuscript to be more self-contained, a brief review of cumulants and the definition of skewness and kurtosisis presented here. Cumulants of a charge distribution are defined by C = (cid:104) Q (cid:105) , (3) C = (cid:104) ( Q − Q ) (cid:105) ,C = (cid:104) ( Q − Q ) (cid:105) ,C = (cid:104) ( Q − Q ) (cid:105) − C , where Q is the net charge. Here, Q might refer to baryon number, to strangeness, or to the electric charge measuredin units of e . Rather than showing the cumulants C n , ratios are presented to help minimize trivial dependences onsystem size. The skewness, S , is a measure of the third moment, S = C C / . (4)This definition has the advantage in being dimensionless, but it does not become independent of volume in the limitof large volumes. Thus, it is more common to consider the ratio Sσ = S (cid:112) C = C C , (5)which becomes an intensive measure in the limit of larger volumes. However, in this study we consider the ratio Sσ C = C C , (6)which is also intensive, and approaches unity for a uncorrelated emission, i.e. a Skellam distribution.The kurtosis is a measure of four-particle correlations, K = C C , (7)but instead of K , one typically chooses Kσ = C C , (8)where σ = C , to find an intensive measure of the fluctuation. For a measure to be intensive, it should be independentof volume in the large-volume limit. For small volumes, charge conservation alters the average densities of variousspecies, which is known as canonical suppression. Canonical suppression also distorts the higher moments for smallervolumes.The ratios C /C and C /C approach simple values in the limit that the distributions would be Poissonian. ForPoissonian emission the observation of a charge in one region of momentum space is uncorrelated with the emissioninto any other space. Thus, particles are correlated only with themselves. If charges appear only in integral positiveunits, one can apply the usual expression for the Poissonian moments where the mean is η , C = n = η (9) C = (cid:104) ( n − n ) (cid:105) = η,C = (cid:104) ( n − n ) (cid:105) = η = C ,C = (cid:104) ( n − n ) (cid:105) − C = η. If there exist both positive and negative charges, the distribution of the net charge can be derived by convoluting thetwo distributions. Convoluting two Poissonians results in a Skellam distribution. If the mean number of positives is η + and the mean number of negatives is η − , the distribution of net charge for a Skellam distribution, Q = n + − n − ,yields the following cumulants C = η + − η − , (10) C = η + + η − ,C = η + − η − = C ,C = η + + η − = C . Thus, if charges are produced in an uncorrelated fashion in increments of ±
1, the skewness and kurtosis become
S σ C = C C = 1 , (11) Kσ = C σ = 1 , where σ ≡ (cid:104) ( Q − Q ) (cid:105) = η + + η − . Even though most of the literature focuses on Sσ = C /C , this study presentsresults for C /C so that one can better understand deviations from the uncorrelated baseline.Moments depend on the efficiency α with which particles are measured. In the limit of vanishing efficiency alldistributions of positives or of negatives tend to become Poissonian [59], and the distribution of the net charge willthus become Skellam. Then Eq. (10) shows that as α → C /C = C /C = 1. This can be understood by seeingthat as α →
0, the moments are dominated by the probability of observing either zero charges or a single charge. Theprobability of observing a single charge is α →
0, while the probability of observing two charges is proportional to α , which is negligible. This assumption would fall through if multitple charges were observed on individual particle,but for final-state hadrons the charges are only ± /
2, and the odd moments will have an odd symmetry, C n (1 / δα ) = (cid:26) C n (1 / − δα ) , n = 2 , , · · ·− C n (1 / − δα ) , n = 3 , , · · · . (12) III. RECURSIVE TECHNIQUES FOR GENERATING CANONICAL PARTITION FUNCTIONS
For non-interacting particles the canonical partition function can be calculated exactly, or at least to the level thatall partitions of A ≤ A max hadrons are taken into account, with the exact solution being reached at A max = ∞ . Forour case, we conserve three quantities: the electric charge Q , the baryon number B and the strangeness S . For states i with energies E i , the partition function, Z ( Q, B, S ) = (cid:88) i,Q i = Q,B i = B,S i = S e − βE i , (13)where Q i , B i and S i are the discrete values of the conserved quantities for the state i , can be calculated recursively.The function Z A ( Q, B, S ) refers to the subset of states with A hadrons, Z ( Q, B, S ) = (cid:88) A ≥ Z A ( Q, B, S ) . (14)The recursive procedure begins with Z A =0 (0 , ,
0) = 1 , (15)the canonical partition function of the vacuum. The contribution for a given A , Z A ( Q, B, S ), can be written as Z A ( Q, B, S ) = 1 A (cid:88) h z h Z A − ( Q − q h , B − b h , S − s h ) , (16)where z h is the single-particle partition function for hadron species h , which has charges q h , b h and s h . This wasproved in [56], and can be understood by realizing that one can count all the ways to arrange A hadrons with a givencharge by considering all the ways to arrange one hadron multiplied by all the ways to arrange the remaining hadrons.To avoid double counting, a factor of 1 /A is applied. For a fixed charge the probability to have A hadrons is, P ( A ) = Z A ( Q, B, S ) (cid:80) A Z A ( Q, B, S ) = Z A ( Q, B, S ) Z ( Q, B, S ) . (17)In practice, the sum over A is cut off at some A max , but in our studies here that cutoff is made large enough thatcontributions to Z for A > A max are negligible. Thus, once one builds the partition function from A = 0 to A max onehas the partition function for all Q, B, S .Once the partition function is calculated one can also calculate the multiplicities and moments of observing specificspecies. For example, the multiplicity of species h in a system with charge Q, B, S is (cid:104) N h (cid:105) = z h Z ( Q − q h , B − b h , S − s h ) Z ( Q, B, S ) . (18)This also provides expressions for the various charges, e.g., (cid:104) Q (cid:105) = (cid:88) h q h z h Z ( Q − q h , B − b h , S − s h ) Z ( Q, B, S ) . (19)Spectra can also be calculated. For species h with spin j h , dN h d p = (2 j h + 1)Ω(2 π (cid:126) ) e − E h ( p ) Z ( Q − q h , B − b h , S − s h ) Z ( Q, B, S ) . (20)Second-order moments can also be calculated exactly, (cid:104) N h N h (cid:48) (cid:105) = δ hh (cid:48) z h Z ( Q − q h , B − b h , S − s h ) Z ( Q, B, S ) + z h z h (cid:48) Z ( Q − q h − q h (cid:48) , B − b h − b h (cid:48) , S − s h − s h (cid:48) ) Z ( Q, B, S ) . (21)It is straightforward to extend this expression to higher-order fluctuations.These expressions can also be extended to consider non-additive conservation laws. Net isospin conservation of ahadron gas was invoked in [55], i.e. restricting the states to being in an iso-singlet. Quark-gluon states restricted tobeing in both an iso-singlet and a color singlet were addressed in [58]. Bose and Fermi corrections are discussed inSec. V.For the case of a single kind of charge, one can see how the the recursive method above yields the same result aswhat one would expect by writing down the partition function for a system with ( A − Q ) / A + Q ) / Z A,Q = (cid:40) z A [( A − Q ) / A + Q ) / , A − Q is even , , A − Q is odd . , (22)where z is the partition function of a single charge. One can readily see that this is consistent with the recurrencerelations, Z A,Q = zA { Z A − ,Q − + Z A − ,Q +1 } (23)= zA (cid:26) z A − [( A − Q − / A + Q ) / z A − [( A − Q ) / A + Q − / (cid:27) = zA (cid:26) z A − ( A − Q ) / A − Q ) / A + Q ) / z A − ( A + Q ) / A − Q ) / A + Q ) / (cid:27) = z A [( A − Q ) / A + Q ) / . (24)This result is also equivalent to expectations based on setting reaction rates equal. If one assumes that pairs arecreated with some rate β , and that they are destroyed with some rate αN + N − , where N + and N − are the number ofpositive and negative charges, N + + N − = A . Setting the rates equal, α ( A − Q )( A + Q )4 Z A,Q = βZ A − ,Q . (25)One can see that if one chooses β/α = z that Eq.s (25) and (22) are consistent. If the net charge is zero, the result iseven simpler, P ( A | Q ) = Z A,Q (cid:80) A Z A,Q , (26) P ( A | Q = 0) = z A [( A/ (cid:40) (cid:88) A =even z A [( A/ (cid:41) − . Aside from the assumptions that Q is fixed and that there exists only one kind of charge, Eq.(25) also requiresthat Bose and Fermi quantum statistical corrections are negligible, and that only unit charges exist. Despite theseshortcomings, this picture is useful in that it allows one to see how multiplicity fluctuations are affected by chargeconservation in a simple model. IV. GENERATION OF UNCORRELATED SAMPLE EVENTS
Complicated experimental acceptances are difficult to incorporate into expressions for the moments. It is theneasiest to generate entire events via Monte Carlo, and filter the events through the acceptance. The Monte Carloprocedure involves choosing a hadron proportional to the number of ways the system might have such hadron, i.e. aproduct of the partition function of the individual hadron multiplied by the partition function of the remainder. Theprocedure becomes:1. Calculate and store the partition function, Z A ( Q, B, S ), up to some size A ≤ A max for values of Q, B, S thatmight ultimately couple back to a given A = A max for the given total values Q, B, S .2. For total charge
Q, B, S , choose the number of hadrons A proportional to Z A ( Q, B, S ) /Z ( Q, B, S ), where Z ( Q, B, S ) = (cid:80) A ≤ A max Z A ( Q, B, S ).3. Choose a hadron h proportional to the probability z h Z A − ( Q − q h , B − b h , S − s h ) /Z A ( Q, B, S ). If Bose degeneracyis to be taken into account this procedure is slightly modified as described in Sec. V.4. Choose the momentum proportional to the thermal weight e − E p /T .5. Repeat (3,4) but with A, Q, B, S being replaced by A − , Q − q h , B − b h , S − s h . The procedure is finished when A = 0.Bose effects can be included by altering the second and third steps above. This is addressed in Sec. V.Storing the partition function can require substantial memory for large A max because the indices Q, B and S mustalso vary over a range of order ± A max , so memory usage roughly scales with A . Because one is usually interestedin calculations with total charge near zero, one can ignore partition functions for charges that cannot couple back tothe fixed overall charge at A max . Once A exceeds A max /
2, the calculations here cut off values of
Q, B and S thatcould not ultimately affect the Q = B = S = 0 partition function for A = A max . Even with this savings, partitionfunctions with A max = 250 could require approximately 12 GB of memory, and need on the order of 10 minutes tocalculate on a single processor. For A max = 125, less than a GB of memory was needed and partition functions couldbe calculated in less than a minute. For hadron gases at temperatures of 150 MeV, A max = 250 was sufficient forpatch volumes (cid:46)
700 fm . If multiple patch volumes are to be explored for the same temperature, computationaltime can also be saved by realizing that the partition functions scale as Ω A . Thus, if one performs a calculation forsome initial volume Ω , scaling can provide results for new volumes with minimal computation.Once the partition function is calculated, event generation is remarkably fast. The time to generate an event scaleslinearly with the volume, or equivalently, linearly with the average number of particles generated. Running sufficientevents to generate a million individual particles can be accomplished within a few seconds on a single CPU. UnlikeMetropolis methods where events are modified by considering small changes to existing events, such as in [54], eachevent in this method is perfectly independent of previous events. V. BOSE AND FERMI STATISTICS
Including Bose and Fermi statistics into the recursive relations for partition functions is straight-forward, and wasshown in [55, 57]. The method is related to that used for calculating the effects of multi-boson interference for pioninterferometry [60]. In a fixed volume the partition function can be first treated as the usual procedure of accountingfor n identical particles being in different single-particle states. This includes the 1 /n ! term to account for the factthat the particles are indistinguishable, i.e. the Gibbs paradox. If m (cid:96) indistinguishable particles are in the samesingle-particle state (cid:96) , one must correct the weight by a factor of m (cid:96) ! for each level, which can also be thought of asthe analog of the symmetrized relative wave function with all the momenta being equal. For fermions, the weightbecomes ( − (cid:96) − m (cid:96) . As demonstrated in [57], the recurrence relation to the partition function then becomes Z A ( Q, B, S ) = A (cid:80) h (cid:80) n Z A − n ( Q − nq h , B − nb h , Q − nq h ) z h,n ( ± n − , (27)where z h,n is the partition function for n particles in some level, z h,n = (cid:80) (cid:96) e − nβ(cid:15) (cid:96) , (28)and (cid:96) refers to single-particle levels of energy (cid:15) (cid:96) . The ± n for any level (cid:96) is of the order e − β(cid:15) (cid:96) lower than the previous term. Thisfactor is largest for zero momentum, and for pions becomes e − βm , where m is the pion mass. For the zero-momentumlevel at a temperature of 150 MeV, the factor is e − m/T ≈ .
4, and as the system cools the factor falls slightly [61].For a more characteristic thermal momentum the factor is ≈ .
1. For heavier particles the factor is always small inthe context of relativistic heavy ion collisions. For example, for a zero-temperature ρ meson the factor is a fraction ofa percent.Given that symmetrization is only being applied to pions, which are bosons, one can incorporate these correctionsinto the Monte Carlo procedure outlined in Sec. IV. For fermions this might be problematic because of the negativeweights coming from the ( − n − factors in Eq (27), but fortunately this is unneccessary because the degeneracy offermions is negligible in the systems considered here. For pions, the algorithm is adjusted by treating each value of n as being a different species, with charges nq h and with the partition function calculated with a reduced temperature T → T /n . If one picks such a species in step 3 of the algorithm, n pions are generated, all with the same momentum.For a finite system, the n pions would be assigned small relative momenta on the order of the inverse system size.It is well known that bosonic effects can broaden multiplicity distributions, consistent with negative binomialdistributions [62, 63], making them super-Poissonian. One of the goals of this study is to discern how bosonicstatistics alter the skewness and kurtosis. VI. MODELS WITH UNIFORM EFFICIENCY AND FIXED CHARGE
Even for a volume of fixed charge, finite efficiency and acceptance leads to non-zero fluctuations. The degree towhich these fluctuations affect the skewness and kurtosis was worked out in [50] for emission of a fixed charge wherethe probability of any charge being observed was a constant α . If α were zero or unity, there would be no fluctuations,and because the charge on those particles that are not observed must fluctuate exactly opposite to the charge thatis observed, the even moments must be symmetric about α = 1 /
2, and the odd moments must be anti-symmetric.One of the most important results of [50] is that for fixed Q and α the ratios of cumulants depend only on α and thevariance and mean of the underlying multiplicity distribution. Even though Q is fixed, the net number of chargedparticles M can fluctuate.From [50], the probability that M charged particles with total charge Q will result in a measured charge q = n + − n − due to a uniform efficiency α is the convolution of two binomial distributions P ( q | M, Q ) = ( M + Q ) / (cid:88) n + =0 ( M − Q ) / (cid:88) n − =0 [( M + Q ) / M − Q ) / M + Q ) / − n + ]![( M − Q ) / − n − ]! n + ! n − ! α n + + n − (1 − α ) M − n + − n − . (29)After some tedious algebra, one can find the cumulants for fixed multiplicity, C (cid:48) = q = αQ, (30) C (cid:48) = (cid:104) ( q − q ) (cid:105) M = α (1 − α ) M,C (cid:48) = (cid:104) ( q − q ) (cid:105) M = α (1 − α )(1 − α ) Q,C (cid:48) = (cid:104) ( q − q ) (cid:105) M − (cid:104) ( q − q ) (cid:105) M = α (1 − α ) − α (1 − α ) M. Here, the primes emphasize that the averages (cid:104)· · · (cid:105) M . denote that they consider only those events with fixed basemultiplicities M . The even moments are all linear in M , while the odd moments are linear in Q . It might appearthat due to the linearity one might replace M with (cid:104) M (cid:105) for a fluctuating base multiplicity, but the second term inthe fourth cumulant − (cid:88) M P ( M ) (cid:104) ( q − q ) (cid:105) M (cid:54) = − (cid:104) ( q − q ) (cid:105) . (31)Instead, (cid:104) ( q − q ) (cid:105) = α (1 − α ) M . (32)This provides a contribution to C , and the cumulants and their ratios become [50] C = αQ, (33) C = α (1 − α ) M ,C = α (1 − α )(1 − α ) Q,C = α (1 − α ) M − α (1 − α ) M +3 α (1 − α ) (cid:104) ( M − M ) (cid:105) ,C C = (1 − α ) MQ ,C C = (1 − α )(1 − α ) ,C C = (1 − α ) QM ,C C = 1 + 3 α (1 − α )( ω M − ,ω M ≡ (cid:104) ( M − M ) (cid:105) M .
The relative variance of the multiplicity of the base distribution, ω M , is unity for a Poissonian distribution. In thatcase C /C and C /C fall below unity for non-zero α . The assumptions for this relation are that the charge andvolume are fixed, that the efficiency is uniform, and that Bose and Fermi symmetrization is ignored. This relationis important as it implies that if one understands the efficiency and the second moment of the base multiplicitydistribution, one can construct a baseline of cumulant ratios, and attribute any deviation from the baseline as due tofluctuations in Q , which is precisely the goal.A reasonable value of α can be taken from balance function analysis. If a charge is observed, there should be abalancing charge emitted nearby, and detected with probability α . This should correspond to the integrated strengthof the charge balance function, which for electric charge is in the neighborhood of 0.35 assuming the full acceptanceof the STAR TPC. Of course, α is not a constant. If a charge is observed in the center of the detector, its balancingcharge has a better chance of being observed than for one observed near the periphery of the acceptance. Nonetheless,for the purposes of roughly setting expectations, this suggests that C /C and C /C can significantly differ fromunity. In the next five subsections, VI A thrugh VI E, we consider how various effects, aside from fluctuating the basecharge Q , might push ω M to be either super-Poissonian, ω M > ω M <
1. We discuss the effectsof volume fluctuations, charge conservation with a single type of charge, decays, Bose condensation, and finally thesensitivity to considering a realistic collection of resonances accounting for all three types of conserved charge.
A. Volume Fluctuations
Heavy ion experiments measure collisions spanning a range of impact parameters. Even for a fixed impact parameter,energy deposition might significantly vary depending on how many nucleons actually collided or how many jets wereproduced. If more energy is deposited into a fixed volume, it might expand further before it hadronizes, resultingin larger volumes when the system hadronizes. Experimental analyses attempt to minimize these fluctuations byconstraining a given fluctuation measurement to a specific centrality bin, where “centrality” might be defined bymultiplicity, transverse energy, or energy deposition in a forward calorimeter. To reduce auto-correlation, centralitymeasurements are usually constrained to particles other than those used to construct the moments. Nonetheless, it isinevitable that a range of initial conditions is explored within any centrality bin, and one might thus expect the basemultiplicity distribution to broaden. If the patch volumes are fixed at some constant value, volume fluctuations canbe thought of as a fluctuation in the number of patches.It is tempting to expect that the ratio of cumulants would be independent of the number of patches, and thusimpervious to volume fluctuations. Independently, each cumulant C n scales linearly with the number of patchesbecause there are no cross correlations between patches. Thus, if P p ( N ) is the probability of having N patches, andif c n is the cumulant for a single patch, the cumulant for the overall system is C n = (cid:88) N P p ( N ) N c n , (34)and the ratio of cumulants is C n C m = (cid:80) N P p ( N ) N c n (cid:80) M P p ( M ) M c m (35)= c n c m . Thus, the ratio of cumulants is independent of the number of patches, and so independent of the overall volume,though the ratio can still depend on the volume of an individual patch.Unfortunately, although the ratios of cumulants are independent of the overall volume, they are not impervious tofluctuations of the overall volume, or equivalently, to fluctuations of the number of patches [64–66]. To illustrate this,one can consider a system with a distribution P p ( N ) as in Eq. (34), where the average number of patches is N . For n > δQ = Q − Q . Because Q varies with the number of patches we rewriteit as δQ = δQ N + δN Q , (36) δQ N ≡ Q − N Q , where Q = c is the average charge in a single patch. Here δQ N is the fluctuation of the charge for a specific number0 N of patches, and δN = N − N . Inserting these definitions into the expressions for the cumulants, C = κ c , (37) C = (cid:88) N P p ( N ) (cid:104) ( δQ N + δN Q ) (cid:105) = κ c + κ c ,C = (cid:88) N P p ( N ) (cid:104) ( δQ N + δN Q ) (cid:105) = κ c + κ c ,C = (cid:88) N P p ( N ) (cid:104) ( δQ N + δN Q ) (cid:105) − C = κ c + 3 κ c + 6 κ c c + κ c , where κ n are the cumulants of P p ( N ), κ = (cid:88) N N P p ( N ) = N , (38) κ = (cid:88) N ( N − N ) P p ( N ) ,κ = (cid:88) N ( N − N ) P p ( N ) ,κ = (cid:88) N ( N − N ) P p ( N ) − κ . If the number of patches is fixed, i.e. the overall volume does not fluctuate, then κ n> = 0 and each cumulant satisfies C n = N c n . The ratio of cumulants then cancels the factor N . Once the volume fluctuates, the ratios of cumulants ofthe charge distribution depend on the ratios of cumulants of P p ( N ). Methods have been constructed to reduce thedependence on volume fluctuations [65–68], and have been applied to STAR data [69, 70]. These methods are built onthe assumption that the observable used to identify the volume is correlated to the fluctuating charge only throughthe fact that they both scale with the volume. If the phase space used for the centrality measure is clearly distinct ofthe phase space over which charge is measured, this should be a good approximation.Equations (37) and (38) are built on the assumption that the average charge scales linearly with the number ofpatches. This differs from the assumptions going into Eq. (33), where it was assumed that charge was fixed within theoverall volume. If the total charge in the volume is absolutely fixed despite fluctuations in the number of patches, thenvolume fluctuations affect the answer in that the relative variance of the multiplicity distribution, ω M , is increasedby volume fluctuations as described in Eq. (33). B. Chemically Equilibrated Canonical Distribution
One can calculate the relative variance, ω M , for an equilibrated system of a single type of charge, where themeasurement directly samples the equilibrated canonical distribution. This ignores decays, which in later stagesproceed without the regeneration needed to maintain equilibrium. The results of the equilibrated canonical ensemblein Eq. (26) can be used to find ω M , the relative variance of the multiplicity distribution, and thus generate themoments using Eq. (33). For large average multiplicities M , the value of ω M approaches unity, the Poissonian limit.For small M it approaches two. This is expected, because for a fixed charge, you can only add particles pair-wise. Infact, as will be shown in the next section, if one generates a Poissonian number of pairs, the relative variance will be ω M = 2 for all multiplicities. Because ω M varies from from two to unity in a canonical ensemble, as the system sizeincreases from zero to infinity, the values of C /C and C /C stay below unity as shown in Fig. 1.The fact that ω M is above unity for small systems (due to only being able to sample even numbers of charges) isknown as canonical suppression. This lowers the mean multiplicity and raises the relative variance, ω M . As seen inFig. 1 it also raises the ratio C /C relative to its value for a large system. C. Decays
A system of uncorrelated neutral particles that decay to charged particles also leads to fixed zero net charge. Butunlike the results for the equilibrated canonical ensemble above, this can lead to ratios of C /C >
1. If each neutral1
Multiplicity M C / C canonical ensemble, α =0 . canonical ensemble, α =0 . decays, N bodies =2 decays, N bodies =4 ,α =0 . decays, N bodies =4 ,α =0 . ω M = σ M / M canonical ensembledecays, N bodies =2 decays, N bodies =4 FIG. 1. Upper panel: The relative variance, ω M , of the charged multiplicity distribution is shown for three cases for a systemcarrying only a single type of ± unit charge. For a neutral equilibrated system (red dashed line) in a canonical ensemble, ω M approaches unity, the Poissonian limit for higher mean multiplicities M . In the low multiplicity limit the multiplicity will beeither M = 0 or M = 2, which gives ω M = 2. This is in contrast to a system where one has a Poissonian distribution of neutralparticles which all decay into pairs (green dashed line) which gives ω M = 2 for all multiplicities, or if the neutral particles alldecay to four charges, which gives ω M = 4. Lower panel: According to Eq. (33), which was derived in [50], the ratio C /C is determined by ω M and the acceptance probability α . It is unity for ω M = 2, below unity for ω M <
2, and above unity for ω M >
2. This shows that charge conservation in an equilibrated system pushes C /C below unity, whereas if an equilibratedsystem of neutral particles decays, and if the decay products do not reform into the resonances, the resulting ratio of C /C isunity for two-body decays, or above unity for four-body decays. particle decays into N bodies charged particles, where the charges are ±
1, the charged particle multiplicity distributionsbecome (cid:104) M ch (cid:105) = N bodies M , (39) (cid:104) ( M ch − M ch ) (cid:105) = N (cid:104) ( M − M ) (cid:105) ,ω M = (cid:104) ( M ch − M ch ) (cid:105) M ch = N bodies (cid:104) ( M − M ) (cid:105) M . If the emission of neutrals is Poissonian, the result is simple, ω M = N bodies . Fig. 1 illustrates how this pictureaffects C /C . In this case, because ω depends only on N bodies and does not change with multiplicity or system size,the cumulants are also independent of multiplicity. For N bodies = 2 the cumulant ratios do not even depend on theefficiency. For N bodies > C /C exceeds unity, so if charge creation proceeded through the creation anddecay of neutral clusters it would be easy to generate large values of C /C . This has been discussed at length in [71].In an equilibrated system, decays and recombination have equal rates. However, as the system decouples recombi-nation stops and decays proceed until only stable hadrons remain. During the hadronic phase the number of charged2 µ (MeV) ω M = σ M / M , C / C ω M C /C FIG. 2. Fluctuations for a canonical ensemble of pions at fixed charge, Q = 0, with an effective chemical potential, µ , appliedto adjust the net pion number in a volume of 500 fm at a temperature of 100 MeV. The relative variance of the multiplicitydistribution (dashed red line) and the ratio C /C of the net-charge distribution (solid green line) grow dramatically as µ approaches the pion mass. Super-radiant effects can take place once µ reaches m π = 139 .
57 MeV. Heavy-ion collisions areexpected to decouple with effective chemical potentials near 75 MeV, which is well below the onset of large fluctuations. particles nearly doubles due to these decays. Thus, if the systems do equilibrate, then decay after chemical freeze-out,one would expect the ratio C /C to lie somewhere between the value for the canonical ensemble in Fig. 1 andthe value for pure decays with N bodies = 2. STAR’s experimental results for net proton fluctuations indeed satisfythis expectation, but their measured C /C for net charge fluctuations exceed unity, which can only be attained if N bodies >
2. Very few hadronic decays proceed via more than one charged pair, but one could have decays of clus-ters or hot spots. Other possible explanations for having C /C > D. Bose correlations
It has long been understood that Bose correlations induce super-Poissonian fluctuations [62, 63]. Here, we illustratehow Bose effects combine with charge conservation to determine C /C . Partition functions for a gas of positive andnegative pions, m π = 139 .
57 MeV/ c , were considered to be kinetically equilibrated but with a chemical potentialenforcing a fixed average density. If a gas is created at chemical equilibrium one expects µ = 0, but if it cools whilemaintaining a fixed number of pions, and if the pion number is fed by decays, an effective chemical potential, µ , isrequired. This changes the average number of pions by adjusting the phase space density as f ( (cid:126)p ) = e − ( E p − µ ) /T − e − ( E p − µ ) /T . (40)Because the system is out of equilibrium, and because the number of positives and negatives are nearly the same, theeffective chemical potential has the same sign for both π + and π − . At decoupling one expects kinetic temperaturesto fall near 100 MeV and the pion chemical potential to grow to perhaps as high as 75 MeV [61]. This estimatecan be understood by the fact that the phase space occupancies should stay roughly constant for fixed entropy perparticle for an expansion at fixed entropy, which suggests a chemical potential of approximately 50 MeV by the timethe system cools to 100 MeV. Decay products further feed the phase space occupancy. Pion condensation occurswhen this chemical potential reaches the pion mass. In the absence of Bose effects, this requires roughly doublingthe phase space density of pions as compared to the expectations above. If such conditions were realized Bose effectscould result in super-radiance [60], which should be accompanied by large multiplicity fluctuations.As a function of the effective chemical potential µ , the partition function for a pion gas in a volume of 500 fm andat a temperature of 100 MeV was calculated from Eq. (27). For µ in the range of of 75 MeV the relative varianceof the charged multiplicity distribution is ω M ≈ .
2, which modestly increases C /C as shown in Fig. 2 for a fixedefficiency of α = 0 .
3. More dramatic results for C /C are not expected unless the chemical potential reaches within afew MeV of the pion mass. As discussed above, this is not expected. However, if the number of pion sources fluctuatedwildly from one event to another, and if there were some events with twice the number of sources emitting into the3same phase space, super-radiance might occur in some small fraction of the events. Such behavior would stronglycontradict expectations based on chemical equilibrium. E. Hadron Gas
Thus far, all the simple examples presented in this section considered a system with one type of conserved charge,but in actuality a hadron gas obeys the conservation of three types of charge: baryon number, strangeness and electriccharge. This invalidates the use of the recurrence relation of Eq. (33) and requires the application of Eq. (16), or ifBose corrections are included, Eq. (27). Here, we calculate the canonical ensemble, Z A ( Q, B, S ), using the recurrencerelation, then apply the Monte Carlo techniques of Sec. IV to generate sample sets of particles. The calculation isbased on a large number, ∼ α is assumed, independent of species, momenta, or whether theproducts came from a weak decay (charged pions and kaons were not decayed).Emission is assumed to come from sub-volumes, or patches, of fixed size. Because each sub-volume is independent,there are no correlations between the various sub-volumes. Also, because cumulants scale linearly with the numberof sub-volumes, ratios of cumulants can depend only on the patch volume, not on the overall volume. However, if thenumber of patches fluctuates, the result is modified by volume fluctuations as shown in Sec. VI A.Here, the sensitivity to three parameters is investigated. The three varied parameters are the patch volume Ω, thebaryon density ρ B , and the fixed efficiency α . When one parameter is varied, the other two are set at default values.The defaults are set at T = 150 MeV, ρ B = 0, α = 0 . . The electriccharge density is set to half of the given baryon density. After these sensitivities are studied, calculations at finitebaryon density are modified so that the net baryon and electric charges responsible for the non-zero charge densityare allowed to fluctuate across sub-volumes according to a Poissonian distribution. This is motivated by the fact thatcharges transported far away from mid-rapidity at early times cause this non-zero baryon density, so local chargeconservation should play little role in the distribution of net charge across patches.The dependence on patch size is exhibited in Fig. 3. For small patches the emission of baryons is discouragedbecause each baryon must be accompanied by the existence of an anti-baryon. Thus, the Boltzmann factor, e − M/T ,for the mass is compounded by the necessity of having a second accompanying massive particle. Once the volumesbecome larger, and the mean number of baryon pairs exceeds unity, the observation of a baryon might be accompaniedby seeing one fewer baryon in the remainder of the system, rather than seeing an extra anti-baryon. In fact, thisis exactly what happens in a canonical ensemble in the large volume limit. In that limit, if a baryon is observed insome small amount of phase space, the mean number of baryons observed in the remainder of the system is decreasedby one half, while the mean number of anti-baryons is increased by one half. It is expected that the dependence ofcumulant ratios on patch size should disappear as the patch size increases beyond the threshold for the mean baryondensity to be unity. Figure 3 displays this behavior by displaying the ratio C /C as a function of patch volume.Because the mean multiplicity for charged particles exceeds unity at a smaller volume, the ratio levels off faster fornet charge than for net protons. Qualitatively, the same behavior was seen for a system with a single type of chargein Fig. 1.The sensitivity to the acceptance probability, α , is displayed in Fig. 4. As α approaches zero, distributions tend tobecome Poissonian, because the number of particles measured is dominantly either zero or unity. Eq. (33) describessimilar behavior for a system of a single type of charge. For net charge, the ratio is symmetric about α = 0 . α = 0 . α → C /C and C /C as a function of baryon density. The patch volume is fixed at 200 fm and thedensities chosen correspond to a fixed number, 4 , , , · · · , of baryons. The choice is made to display C /C ratherthan C /C because C /C would be unity for uncorrelated emission. The ratios all decrease, moving further fromunity, as baryon density increases. The densities displayed cover the range of what might be reached in heavy-ionsystems at the point when temperatures fall below 150 MeV.For the calculations illustrated in Fig. 5, the net baryon number, electric charge, and strangeness were all fixed.Because of local charge conservation and the limits of diffusion, balancing charges created after the collision areconstrained to stay within the same neighborhood. This constraint is adjusted by setting the size of the patchvolume. However, even for measurements at mid-rapidity, non-zero charges can arise due to the transfer of baryonnumber and electric charge from the projectile and target rapidities. These intruder charges are not typically balancedby charges within the acceptance of the detector, and thus their fluctuation should be considered separately. Figure6 shows how C /C and C /C behave as a function of baryon density, if the non-zero baryon number and electric4 Patch Volume (fm ) C / C ChargeProtons
FIG. 3. The ratio of cumulants is shown for a ρ B = 0 system as a function of the patch volume. For small patches a baryonis always accompanied by an extra anti-baryon, but for larger systems the observation of a baryon might also enhance theprobability that one fewer baryons is present in the remainder. C / C ChargeProtons
FIG. 4. As the fixed acceptance probability α approaches zero, distributions become Poissonian, and the ratios C /C approach unity. Even for perfect acceptance, the fluctuations are non-zero for a multi-charge system because the conservedcharges carried by a proton can be balanced by an array of other species, and charge within the proton and anti-proton sectoris not conserved. Thus, the ratio is not symmetric about α = 0 . charge are chosen to fluctuate within the patch according to a Poissonian distribution. As expected, the ratios rise ascompared to the fixed case shown in Fig. 5, but they do not exceed unity. Discerning these fluctuations of the actualnet charges is of particular interest. These results suggest that the starting point for the ratios C /C and C /C isbelow unity for random uncorrelated net charges accompanied by an ensemble of baryon charges. VII. BLAST WAVE MODEL WITH A FULL HADRON GAS AND COMPARISON TO STAR RESULTS
The calculations of the previous section were based on a simple picture, where each sub-volume emitted particleswhose probability of being observed was uniform, denoted by α . In practice, this probability depends on where thesub-volume is located within the overall reaction volume. A sub-volume in the region with spatial rapidity, | η s | > ± ± . y is less than5 B (fm ) C / C , C / C C / C (Charge) C / C (Charge) C / C (Protons) C / C (Protons) FIG. 5. The ratios C /C and C /C are shown for net charge and for net baryon number as a function of the baryon density.The ratios fall increasingly below the Poissonian limit of unity as the baryon number is increased. For this calculation the netbaryon charge was fixed at B = ρ B V and the net electric charge was fixed at Q = B/ B (fm ) C / C , C / C C / C (Charge) C / C (Charge) C / C (Protons) C / C (Protons) FIG. 6. The ratios C /C and C /C are shown for net charge and for net baryon number as a function of the baryon density.Calculations differ from those in Fig. 5 in that the net charge fluctuates according to a Poissonian distribution. The ratios riserelative to those in the fixed-charge case, but they remain below unity. the pseudo-rapidity η . This difference is magnified for more massive or slower particles. In fact, STAR analyses ofidentified particles often enforce cuts that only consider particles with true rapidities − . < y < .
5. Of course,even particles within the rapidity acceptance must exceed some minimum transverse momentum, and the efficiencyfor being detected is imperfect. For particles identified only by charge, the efficiencies are typically (cid:38) η s moves with a rapidity η s , and particles emitted from a sub-volumewith rapidity η s have rapidities y ≈ η s . This simple mapping is smeared by thermal motion. Collective radial flowand cooling combine to better align the spatial rapidity η s and the measured rapidity y . The thermal spread for pionsis ≈ . ≈ . η s , should affect the result. If theextent in η s is small, there is an enhanced probability that a charge and its balancing charge will both be identifiedand cancel one another when assigning the net charge for an event. If the region extends over a large rapidity range,the effects of charge balance are minimized because there is a better chance that one charge will be observed while theother is outside the acceptance. The extent of the region in the transverse direction is less important. For emittingregions at the edge of the fireball, which have more collective velocity, there is a modest increase of having balancing6charges both pushed into the acceptance. In addition to the size of the region over which particles from a givensub-volume are emitted, the overall size of the sub-volume at the point where chemical freeze-out occurs also plays arole as it sets the degree to which canonical suppression affects the results. As discussed in the previous section, thismatters only for for smaller sub-volumes.For the model in this section, canonical sub-volumes are overlaid onto a blast-wave parametrization of collectiveflow. A filter is applied to the calculations, representing the experimental STAR acceptance and efficiency. This shouldbe sufficiently realistic to make meaningful comparisons to STAR data. The four main parameters for the blast-wavemodel describe the radial flow u ⊥ , the kinetic freeze-out temperature, T k , and the baryon chemical potential andtemperature at chemical freeze-out, µ c and T c . The chemical freeze-out temperature µ c and chemical potential, T c ,are chosen to fit relative particle yields, while T k and u ⊥ are determined by simultaneously fitting the spectra ofspecies with varying mass, typically π, K, p . A variety of parameterizations exist, such as having the velocity increaselinearly from the origin, or the transverse rapidity, or having a sharp cutoff in radius vs. assuming a Gaussian profile.Depending on the choice, the value of T k and u ⊥ vary, but are typically in the neighborhood of T k ≈
100 MeV and u ⊥ ≈ .
6. For increasing multiplicities, the reaction volumes can expand and cool further, which leads to modestlyincreased values of u ⊥ and modestly decreased values of T k for either more central or for more energetic collisions.For this study the chemical freeze-out parameters were taken from [73, 74], which extracted T c and µ c for a variety ofbeam energies. For this section, Bose statistics were ignored. Over 300 hadron species were included in the analysis.The collective flow velocities were chosen to reproduce the mean transverse momenta of both pions and protons.The distribution of spatial rapidities η s over which all particles are emitted is chosen to be Gaussian, with a width σ that depends on beam energy as dNdη s ∼ e − η s / σ , (41) σ = 0 . y beam , where ± y beam are the rapidities of the incoming beams. This choice reproduces the rapidity widths measured atRHIC to roughly 10% accuracy [75]. The distribution of spatial rapidities for particles of a given sub-volume wasalso spread according to a Gaussian with a width σ η . The parameter σ η describes how far particles created from thesame sub-volume may have separated by the time of emission. For increasingly larger values of σ η , the chance thatany two observed particles are correlated decreases. The rapidity of a sub-volume was then distributed according toa Gaussian with width (cid:113) σ − σ η , so that the emission overall was characterized by the width σ . Whereas the otherblast-wave parameters were chosen to describe spectra and yields, the parameter σ η is related to charge conservation.If charge is created early, and especially if the diffusion constant is large, the width might be close to one unit ofspatial rapidity, whereas if all the charge were to come after hadronization, the width might more likely be a fewtenths. Blast-wave models of charge balance functions, which are also driven by charge conservation and diffusion,suggest widths should be of the order of a third, but variations of a factor of two were not ruled out [51]. One goalof this section is to investigate the sensitivity to σ η .For this blast-wave model, each sub-volume was also assigned a radial transverse velocity according to a Gaussian, dNd u ⊥ ≈ e − ( u x + u y ) / u ⊥ , (42)and a small spread in (cid:126)u ⊥ . However, given that for this study, cuts are not being considered in transverse momenta orazimuthal angle, varying (cid:126)u ⊥ has little effect. Particles from each differential volume element were assigned momentausing the Monte Carlo algorithm described in [76].The calculations of this section are filtered through a simplified model of the STAR detector’s acceptance andefficiency. For unidentified particles, pseudo-rapidities are required to be between ± . c and below 2 GeV/ c . For identified particles, rapidities were restricted tobeing between ± .
5, and transverse moment were required to be between 200 MeV/ c and 1.6 GeV/ c for pions orkaons, and between 400 MeV/ c and 2 GeV/ c for protons and antiprotons. Because the STAR data are correctedfor efficiency, the efficiencies were set to unity. In order to compare results to measurement from STAR, the baryondensities and chemical freeze-out temperatures were mapped to beam energy according to the analysis of [74], whichextracted chemical freeze-out temperatures and chemical potentials by considering ratios of particle yields. Theparameter u ⊥ , which controls the transverse radial flow, was chosen along with the kinetic breakup temperature, T k ,to simultaneously roughly fit the mean transverse momenta of both the protons and pions reported [77]. The kinetic7freeze-out parameters were chosen according to the parameterization, T k = 120 − f MeV , (43) u ⊥ = 0 . . − . f, (44) f = 1ln(2) ln (cid:34) √ s − √ s √ s f − √ s (cid:35) , (45)where √ s = 7 . √ s f = 200 GeV. Decays were simulated, and the products of weak decays (aside frompions or charged kaons) were included in the analysis. Undoubtedly, a more realistic model of the acceptance mightchange the ratios, but given that these are ratios, and that the overall efficiency and acceptances are not wildly off, amore rigorous model of the acceptance is unlikely to change the result by more than a few percent.The summation over sub-volumes was performed by summing over 400 values of u ⊥ and η s consistent with thedistributions described above. For each value, emission of N sample = 2 × sub-volumes was simulated so thatcumulants could be calculated for each value of η s and u ⊥ . Because emission from different sub-volumes is uncorrelated,the cumulants for the total emission are simply the sum of the cumulants of each sub-volume. Further, because thenet cumulants behave linearly in N sample , the ratios of net cumulants is independent of N sample . s NN (GeV) (c) net protons C / C (b) net kaons (a) net charge = 25= 50= 100= 200= 400 STAR(preliminary) s NN (GeV) (f) net protons C / C (e) net kaons (d) net charge = 25= 50= 100= 200= 400 STAR(preliminary) s NN (GeV) (i) net protons = 25= 50= 100= 200= 400 STAR(preliminary) C / C (h) net kaons (g) net charge FIG. 7. (color online) Ratios of cumulants from blast wave models, begining with C /C (panels a - c), then skewness (panels d - f ) and finally, kurtosis (panels g - i ), are plotted for for different values of the sub-volume Ω. For smaller sub-volumes theratios C /C , C /C and C /C increase. For net charge, the ratios approach an asymptotic value once Ω begins to pass ∼ , whereas for net protons the ratios appear to approach the limit at somewhat higher values of Ω. STAR measurementsfor net protons are not wholly dissimilar to the blast-wave calculations here, but those for net charge differ greatly. Additionalphysics from volume fluctuations might explain how C /C might exceed unity, but it is difficult to explain how this mighthappen for net-charge distributions while leaving C /C of the net-proton distribution unchanged. Figure 7 illustrates how results are sensitive to the size of the sub-volume, Ω. This is the effective volume at whichcharge conservation is enforced at chemical freeze-out. For these calculations σ η was fixed at 0.3. For a volume of 100fm the average number of hadrons is several dozen. For small sub-volumes, where in a grand canonical ensemble thetypical number of charges would be zero or one, the thermodynamic cost of having a second charge to balance the first8charge reduces the probability of having any charges. As mentioned earlier, this is known as canonical suppression,and it lowers both the multiplicities and the moments. For larger sub-volumes, this thermodynamic cost vanishes asthe system might have had fewer balancing charges of the same sign, rather than an extra balancing charge of theopposite sign. The characteristic volume for the disappearance of canonical suppression is the volume where the meannumber of pairs exceeds unity. In general, for charged particles this sets in at around 50 fm , but for baryons thecharacteristic volume is closer to 100 fm because of their being heavier and fewer. Thus, the proton moments forvolumes of 50 fm and 200 fm differ more noticeably. The time for a fluid element to expand and cool to the pointwhere it reaches chemical freeze-out tends to be on the order of 5 fm/ c . These times are shorter for matter on theperiphery, at a lower beam energy or centrality. The times are longer for fluid elements at the center, at higher beamenergy, or in a more central collision. The maximum transverse distance a charge can travel before chemical freeze-outis ≈
10 fm if they move in opposite directions, and if charge moves diffusively or if the charge is created later in thereaction, the separation should be significantly less. As discussed below, charge can separate further along the beamaxis, and because that is not well understood, the size of the sub-volume carries a large uncertainty. Anywhere from50 fm to a few hundred fm might be reasonable.Balancing charges separate from one another, and depending on when they are created, they can diffuse apartfrom one another. This separation, represented by the parameter σ η , accounts for the separation of balancing chargesboth before and after hadronization, or both before and after chemical freezeout. After chemical freeze-out chargesmight separate and mix between the sub-volumes. Because some of the separation comes after chemical freeze-out,this distance might exceed the scales representing the size of the sub-volume Ω. Charge can spread further in thelongitudinal direction due to the strong initial longitudinal collective flow at early times. This enhancement to theseparation depends sensitively on when charge pairs are created. Matter thermalizes at an early time, where largecollective velocity gradients along the beam axis are expected. If the motion is diffusive, the separation along thebeam axis depends logarithmically on the ratio of the final time to the initial creation time [42]. Thus, if a pair iscreated at 0.2 fm/ c , the separation for times, 0 . < τ < . c is as important as the additional separation theygain during the times 1 . < τ < . c . For this reason the size of charge spread, σ η , in spatial rapidity might beanywhere between a few tenths of a unit of rapidity to a full unit. Figure 8 shows the sensitivity of the moments tothis parameter. For large σ η the observation of a charge is less likely to influence the observation of a second charge,which is similar to having a lower efficiency. For low efficiencies one expects the behavior to be more Poissonian, and C /C and C /C to be closer to unity. Indeed, this is the case, but the dependence on σ η , as shown in Fig. 8, isnegligible for net charge and modest for net protons.The beam energy dependence mainly derives from the fact that the baryon density is higher for lower beam energies.Given that C falls with increasing beam energy, while C increases because multiplicity increases, it is no surprisethat the ratio C /C falls with increasing beam energy. The skewness, C /C , should approach unity if emission israndom, i.e. a Skellam distribution. The effects of local charge conservation keep C /C <
1. In the limit that thereare no anti-baryons, which is the limit of high net baryon density and equivalently low beam energy, the assumptionthat baryons are deposited amongst the sub-volumes according to a Poissonian distribution also drives C /C for netprotons closer to unity. Combined with higher canonical suppression for heavier particles, the sensitivity to Ω, shownin Fig. 7, or to σ η , illustrated in Fig. 8, is more pronounced for blast-wave calculations of net protons than of netkaons or of net charge.The measure of kurtosis, C /C , varies only modestly with beam energy. For net charge, the ratio does not varyfar from 0.8 in model calculations. In contrast, C /C from model calculations for net protons varies from 0.6 to 0.9depending on the values of Ω and σ η . The models also exhibit a modest dependence on beam energy for C /C . Theratio rises approximately 20 percent as the beam energy increases to 20 GeV, then plateaus and falls 10% until itbecomes flat. Similar to the C /C ratio, this ratio stays below the Skellam limit.Figures 7 and 8 also display results from the STAR Collaboration for data recorded for central collisions, 0-5%centrality [10, 20, 21]. For net proton fluctuations, the blast-wave model results are similar to STAR measurementsfor all three ratios. It is difficult to conclude the meaning of the solid agreement shown in the C /C ratios. Thisagreement covers net charge, net kaons and net protons, as long as the sub-volume Ω is of the order of 50 fm or greater.Interpreting the meaning of this agreement is difficult because it is mainly driven by having the model correctly matchthe particle yields with multiplicity. Unlike the C /C and C /C ratios there is no simple baseline for a Skellamdistribution. The Skellam baseline would depend on knowing the charged particle multiplicity distributions, i.e. themoments of the multiplicity distributions for protons plus antiprotons, kaons plus antikoans, or all charged particles.For the net-proton distribution, the behavior of C /C as a function of beam energy seems consistent with theexperimental uncertainties, but the statistical errors in the experimental data are so large, that little can be concludedaside from the fact that the models predict that C /C should be in the neighborhood of 0.75. These ratios from themodels can show a modest sensitivity to beam energy, but any such trends are overwhelmed by the statistical errorsof the experimental results at the moment.For both net kaons and for net charge, C /C and C /C lie above the range of model predictions, but the large9 s NN (GeV) (c) net protons C / C (b) net kaons (a) net charge = 0.1= 0.3= 0.5= 0.7 STAR(preliminary) s NN (GeV) (f) net protons C / C (e) net kaons (d) net charge = 0.1= 0.3= 0.5= 0.7 STAR(preliminary) s NN (GeV) (i) net protons = 0.7= 0.5= 0.3= 0.1 STAR(preliminary) C / C (h) net kaons (g) net charge FIG. 8. (color online) Ratios of moments are displayed for different values σ η , which sets the longitudinal size of over whichsub-volumes emit charge. A modest sensitivity is found for net protons, while the moments for net charge were fairly insensitive. experimental error bars forbid one from stating this with great confidence. Should the experimental results withimproved statistics confirm this discrepancy, it will be difficult to explain unless a significant number of charges areemitted from large clusters. In contrast, the model should always give ratios below unity for C /C , regardless of thechoice of parameters. These failures of the blast-wave model for net kaons and for net charge are discussed in theupcoming summary, Sec. VIII. VIII. SUMMARY
The principal goals of this study were to clarify background contributions for higher moments of charge distributionsmeasured by the STAR Collaboration at RHIC, and to state the degree to which current experimental results areeither consistent or inconsistent with these contributions. By background, this refers to sources of fluctuations besidesthose that arise from baryon number or charge clustering due to processes such as phase separation. The list ofsuch sources includes charge conservation, Bose corrections, volume fluctuations, and the decays of resonances. Inorder to gain better insight both simple semi-analytic models with a single type of conserved charge, similar to thework performed in [50], and a more realistic blast-wave model which includes a more realistic accounting of theSTAR acceptance, similar to what was applied in [54], were investigated. By using a highly efficient algorithm forMonte Carlo generation of particles according to the canonical ensemble, results were produced with small statisticaluncertainties. This enabled the exploration of sensitivities to critical parameters of the model.Of the various background correlations, charge conservation provided the strongest non-Poissonian contributions.Because fourth-order cumulants were defined to subtract contributions from second-order correlations, one mighthave expected a small contribution. Consistent with the result from [50] for a single type of charge, it was found thatgenerating sets of particles from a sub-volume equilibrated according to the canonical ensemble produced values of C /C and C /C which were significantly lower than the Skellam value of unity, which is what one would expectfor uncorrelated emission. In contrast, the contribution from two-particle decays does not change either C /C or0 C /C . This conclusion persisted for the more realistic blast-wave model which incorporated the conservation of allthree charges and included a filter of the STAR acceptance. The correlation varied modestly according to the size ofthe canonical sub-volume for small sub-volumes due to canonical suppression. A modest sensitivity was also foundto the spatial extent of this volume along the beam axis, as the overlay of collective longitudinal flow onto the finiteacceptance in rapidity effectively lowers the probability for two balancing charges to both be observed. Due to the factthat baryon density falls with increasing beam energy, the strength of such correlations do depend on beam energy.But this sensitivity was not dramatic. Thus, even though this background contribution is rather large, it is quitesmooth with respect to beam energy, so if sharp non-monotonic structures are observed experimentally with respectto beam energy, such structures are not likely to be driven by charge conservation.A second source of background arises from multi-particle symmetrization of the outgoing pions. By extending therecursive techniques applied for the canonical ensemble to include symmetrization, it was found that such effectsshould not affect the skewness or kurtosis unless the pion phase space density were to become surprisingly large. Inorder for the symmetrization effects to become large, the pion phase space density in the absence of symmetrizationwould have to double. If this were the case, the pion spectra would be more dramatically altered by Bose effects andthe measured HBT radii would have to be significantly altered.Volume fluctuations are somewhat of a wildcard for background processes. As shown in Sec. VI A, such fluctuationscan significantly increase the C /C ratio. The STAR Collaboration invested great effort in minimizing their impact,but it is difficult to gauge the degree to which such effects might have persisted. Volume fluctuations also increasesimilar moments of the multiplicity distribution, which is constructed by counting charged particles rather than netcharge. It is critical for the experiments to simultaneously present moments of the multiplicity distribution alongsidethose of the net-charge distributions. This sensitivity was also illustrated with the simple results from a system withone charge and uniform acceptance, restated in Eq. (33) from the previous work of [50]. Because volume fluctuationsshould similarly increase the C /C ratio for net charge, net strangeness and net baryon number, behavior of suchratios for one type of charge that are not seen in other types of charge can be considered as originating from someother effect.Due to the way in which cumulants are constructed, the third and fourth-order cumulants should be imperviousto the effects of two-particle decays. However, decays of clusters that produce four or more charged particles docontribute to C . Here, it was found that if a significant fraction of charged particles come from such decays, thehigher moments can be profoundly altered, as was already known from the work in [71]. However, the decays of suchclusters would also affect the multiplicity distribution, which again underscores the importance of the experimentsto simultaneously analyze fluctuation of net charge and of the the multiplicity distribution. If such effects wereimportant, it would suggest novel contributions to the dynamics of charge production, outside of the usual paradigmof creating equilibrated distributions of hadrons.The results of the blast-wave calculations were displayed alongside STAR results in Sec. VII. The experimentalresults had much larger statistical errors than the calculations, which limits the conclusions that can be drawn. Thefluctuations of net protons were not far from the range of those calculated here. This is consistent with chargeconservation being the dominant source of non-Poissonian behavior, i.e. C /C (cid:54) = 1 and C /C (cid:54) = 1. By no meansdoes this suggest that this is evidence for a lack of more novel sources of correlation in the baryons, such as that arisingfrom phase separation. The experimental uncertainties may be currently too large to unmask such phenomena.The observed moments of the net-charge distribution are perplexing. Although there are large statistical uncer-tainties in the experimental data, it appears that C /C exceeds expectations of the blast-wave calculations, andeven lies above unity. Given that the net-proton distributions are roughly in line with the model predictions, thisdiscrepancy cannot be explained by volume fluctuations, or equivalently by the systematics of event binning. Becausephase transition phenomena are expected to manifest themselves mainly in the net-proton distributions, this wouldsuggest that the decays of larger clusters into four or more charged particles might be present. If such processes werepresent, it would motivate a rethinking of models of chemical evolution and charge production in heavy-ion collisions.In order to confirm or dismiss this hypothesis, it is imperative that a simultaneous analysis of charge (not net charge)multiplicity distributions be undertaken from the same data sets with the same cuts on centrality. Given the greatlyimproved data sets currently being analyzed by STAR in the Beam Energy Scan II program at RHIC [78], this puzzleshould be clarified in the next year or two.Finally, the studies here help point the way to future improvements in modeling. The picture of independentcanonical sub-volumes is crude. It does incorporate the truth that charge conservation is enforced locally, over somelength scale, and is sufficient to provide the understanding of how large such effects might be. However, in realitybaryon number, electric charge and strangeness are created and evolve in different ways. Strangeness tends to becreated early in the collision, thus allowing the balancing strange and anti-strange quarks to separate before theemission of the hadrons to which they are asymptotically assigned. In contrast, up and down quarks are more likelyto be produced later in the reaction. Thus, the characteristic canonical volumes should have different sizes anddifferent longitudinal extents depending on whether one is considering up, down or strange quarks. This is also true1for off-diagonal correlations, e.g. correlations between strange and up. These only appear in the hadronic phase. Adiagrammatic formalism has been developed for evolving such two-, three-, and four-body correlations as a functionof the positions of each charge. These equations are based on knowing the chemical evolution of the system andthe diffusion constants for each type of charge. Such calculations have been performed for two-body correlationsand compared to experimentally measured charge-balance functions [79, 80] or fluctuations [27]. Unfortunately, themethod for three- and four-body correlations is challenging to implement numerically [81]. It is tractable, but wouldrequire significant effort. The study presented here suggests that such calculations would be warranted if one wantedto reproduce the moments of these distributions for each type of charge, and especially if one wants to consider crossterms [82], such as moments involving powers of both charge and net baryon number. Otherwise, given that chargeconservation effects are expected to evolve smoothly with beam energy, one could simply see whether measurementsof the ongoing Beam Energy Scan at RHIC unveil any sharp features, and assign such features to more novel typesof physics. ACKNOWLEDGMENTS
This work was supported by the Department of Energy Office of Science through grant number DE-FG02-03ER41259and through grant number DE-FG02-87ER40328. R. Steinhorst was additionally supported by the MSU ProfessorialAssistantship program and by the Director’s Research Scholars program at the National Superconducting CyclotronLaboratory. The work benefited by discussions within the BEST Collaboration, which is a topical collaborationfunded by the Department of Energy Office of Science, and especially from discussions with Dmytro Oliinychenko. [1] B. Abelev et al. [STAR], Phys. Rev. C , 024906 (2009).[2] A. Sarkar [STAR], Nucl. Phys. A , 796-801 (2014).[3] J. Xu [STAR], J. Phys. Conf. Ser. , no.1, 012073 (2017).[4] J. Thder [STAR], Nucl. Phys. A , 320-323 (2016).[5] X. Luo [STAR], Nucl. Phys. A , 911c-914c (2013).[6] T. J. Tarnowsky [STAR], J. Phys. G , 124054 (2011).[7] M. Aggarwal et al. [STAR], Phys. Rev. Lett. , 022302 (2010).[8] A. Chatterjee [STAR], PoS CORFU2018 , 164 (2019).[9] J. Adam et al. [STAR], Phys. Rev. C , no.1, 014902 (2019).[10] L. Adamczyk et al. [STAR], Phys. Rev. Lett. , 092301 (2014).[11] A. Sarkar [STAR], EPJ Web Conf. , 13013 (2013).[12] X. Luo [STAR], PoS CPOD2013 , 019 (2013).[13] N. R. Sahoo [STAR], DAE Symp. Nucl. Phys. , 766-767 (2012).[14] X. F. Luo [STAR], J. Phys. Conf. Ser. , 012003 (2011).[15] N. R. Sahoo [STAR], [arXiv:1101.5125 [nucl-ex]] (2011).[16] N. R. Sahoo [STAR], DAE Symp. Nucl. Phys. , 1028-1029 (2011).[17] T. K. Nayak [STAR], Nucl. Phys. A , 555C-558C (2009).[18] D. McDonald [STAR], Nucl. Phys. A , 907c-910c (2013).[19] D. McDonald, Ph.D. thesis, Rice University, https://drupal.star.bnl.gov/STAR/files/mcdonaldthesis.pdf (2013).[20] J. Adam et al. [STAR], [arXiv:2001.02852 [nucl-ex]] (2020).[21] L. Adamczyk et al. [STAR], Phys. Lett. B , 551-560 (2018) [arXiv:1709.00773 [nucl-ex]].[22] M. A. Stephanov, K. Rajagopal and E. V. Shuryak, Phys. Rev. Lett. , 4816-4819 (1998).[23] M. A. Stephanov, K. Rajagopal and E. V. Shuryak, Phys. Rev. D , 114028 (1999).[24] S. Borsanyi, Z. Fodor, S. D. Katz, S. Krieg, C. Ratti and K. Szabo, JHEP , 138 (2012).R. Bellwied, S. Borsanyi, Z. Fodor, S. D. Katz, A. Pasztor, C. Ratti and K. K. Szabo, Phys. Rev. D , no. 11, 114505(2015).[25] G. Aarts, C. Allton, A. Amato, P. Giudice, S. Hands and J. I. Skullerud, JHEP , 186 (2015).[26] A. Amato, G. Aarts, C. Allton, P. Giudice, S. Hands and J. I. Skullerud, Phys. Rev. Lett. , no. 17, 172001 (2013).[27] M. A. Aziz and S. Gavin, Phys. Rev. C , 034905 (2004).[28] J. I. Kapusta and C. Young, Phys. Rev. C , no. 4, 044902 (2014).[29] J. I. Kapusta and C. Plumberg, Phys. Rev. C , no. 1, 014906 (2018).[30] J. Hammelmann, J. M. Torres-Rincon, J. B. Rose, M. Greif and H. Elfner, Phys. Rev. D , no. 7, 076015 (2019).[31] H. Wang, Ph.D. Thesis, arXiv:1304.2073 [nucl-ex].[32] B. I. Abelev et al. [STAR Collaboration], Phys. Lett. B , 239 (2010).[33][33] N. Li et al. [STAR Collaboration], Indian J. Phys. , 923 (2011). [34] C. Alt et al. [NA49 Collaboration], Phys. Rev. C , 024914 (2007).[35] J. Adams et al. [STAR Collaboration], Phys. Rev. Lett. , 172301 (2003).[36] M. M. Aggarwal et al. [STAR Collaboration], Phys. Rev. C , 024905 (2010).[37][37] C. Alt et al. [NA49 Collaboration], Phys. Rev. C , 034903 (2005).[38] B. Abelev et al. [ALICE Collaboration], Phys. Lett. B , 267 (2013).[39] L. Adamczyk et al. [STAR Collaboration], Phys. Rev. C , no. 2, 024909 (2016).[40] L. Adamczyk et al. [STAR Collaboration], Phys. Rev. C , no. 6, 064911 (2013).[41] B. I. Abelev et al. [STAR Collaboration], Phys. Rev. Lett. , 251601 (2009).[42] S. A. Bass, P. Danielewicz and S. Pratt, Phys. Rev. Lett. , 2689 (2000).[43] P. Bozek, Phys. Lett. B , 247 (2005).[44] P. Bozek, W. Broniowski and W. Florkowski, Acta Phys. Hung. A , 149 (2005).[45] S. Pratt, W. P. McCormack and C. Ratti, Phys. Rev. C , 064905 (2015).[46] S. Pratt, J. Kim and C. Plumberg, Phys. Rev. C , no. 1, 014904 (2018).[47] B. Ling, T. Springer and M. Stephanov, Phys. Rev. C , no. 6, 064901 (2014).[48] S. Cheng, S. Petriconi, S. Pratt, M. Skoby, C. Gale, S. Jeon, V. Topor Pop and Q. H. Zhang, Phys. Rev. C , 054906(2004).[49] S. Pratt, Phys. Rev. Lett. , 212301 (2012).[50] O. Savchuk, R. V. Poberezhnyuk, V. Vovchenko and M. I. Gorenstein, Phys. Rev. C , no.2, 024917 (2020).[51] S. Schlichting and S. Pratt, Phys. Rev. C , 014913 (2011).[52] S. Schlichting and S. Pratt, [arXiv:1005.5341 [nucl-th]] (2010).[53] S. Pratt, S. Schlichting and S. Gavin, Phys. Rev. C , 024909 (2011).[54] D. Oliinychenko, S. Shi and V. Koch, [arXiv:2001.08176 [hep-ph]].[55] S. Cheng and S. Pratt, Phys. Rev. C , 044904 (2003).[56] S. Pratt and S. Das Gupta, Phys. Rev. C , 044603 (2000).[57] S. Pratt, Phys. Rev. Lett. , 4255-4259 (2000).[58] S. Pratt and J. Ruppert, Phys. Rev. C , 024904 (2003).[59] A. Bzdak and V. Koch, Phys. Rev. C , 044904 (2012).[60] S. Pratt, Phys. Lett. B , 159-164 (1993).[61] C. Greiner, C. Gong and B. Muller, Phys. Lett. B , 226-230 (1993).[62] P. Carruthers and C. C. Shih, Phys. Lett. B , 242-250 (1983).[63] P. Carruthers, E. M. Friedlander, C. C. Shih and R. M. Weiner, Phys. Lett. B , 487-492 (1989).[64] M. Gazdzicki and S. Mrowczynski, Z. Phys. C , 127-132 (1992).[65] M. I. Gorenstein and M. Gazdzicki, Phys. Rev. C , 014904 (2011).[66] M. Gazdzicki, M. I. Gorenstein and M. Mackowiak-Pawlowska, Phys. Rev. C , no.2, 024907 (2013).[67] V. V. Begun, M. I. Gorenstein and K. Grebieszkow, J. Phys. G , no.7, 075101 (2015).[68] E. Sangaline, [arXiv:1505.00261 [nucl-th]].[69] X. Luo and N. Xu, Nucl. Sci. Tech. , no.8, 112 (2017).[70] S. Esumi and T. Nonaka, [arXiv:2002.11253 [physics.data-an]].[71] A. Bzdak, V. Koch, D. Oliinychenko and J. Steinheimer, Phys. Rev. C , no.5, 054901 (2018).[72] M. Tanabashi et al. [Particle Data Group], Phys. Rev. D , no.3, 030001 (2018).[73] S. Das [STAR], Nucl. Phys. A , 891c-894c (2013).[74] L. Kumar [STAR], Central Eur. J. Phys. , 1274-1277 (2012).[75] C. E. Flores [STAR], Nucl. Phys. A , 280-283 (2016).[76] S. Pratt and G. Torrieri, Phys. Rev. C , 044901 (2010).[77] B. I. Abelev et al. [STAR], Phys. Rev. C , 024911 (2010).[78] C. Yang [STAR], Nucl. Phys. A , 800-803 (2017).[79] S. Pratt and C. Plumberg, Phys. Rev. C , no.4, 044916 (2019).[80] S. Pratt and C. Plumberg, [arXiv:1904.11459 [nucl-th]] (2019).[81] S. Pratt, Phys. Rev. C , no.1, 014914 (2020).[82] N. Abdelwahab et al. [STAR], Phys. Rev. C92