text:
# 𝐵_𝑑→ϕ𝐾_𝑆 CP asymmetries as an important probe of supersymmetry

## Abstract

The decay $`B_d\varphi K_S`$ is a special probe of physics beyond the Standard Model (SM), since it has no SM tree level contribution. Motivated by recent data suggesting a deviation from the SM for its time-dependent CP asymmetry, we examine supersymmetric explanations. Chirality preserving contributions are generically small, unless gluino is relatively light. Higgs contributions are also too small to explain a large asymmetry. Chirality flipping $`LR`$ and $`RL`$ gluino contributions actually can provide sizable effects without conflict with all related results. We discuss how various insertions can be distinguished, and argue the needed sizes of mass insertions are reasonable.

preprint:     December, 2002 KAIST-TH 2003/04  hep-ph/0212092 MCTP-03-17  MADPH-03-1327

$`B\varphi K`$ is a powerful testing ground for new physics. Because it is loop suppressed in the standard model (SM), this decay is very sensitive to possible new physics contributions to $`bss\overline{s}`$, a feature not shared by other charmless $`B`$ decays. Within the SM, it is dominated by the QCD penguin diagrams with a top quark in the loop. Therefore the time dependent CP asymmetries are essentially the same as those in $`BJ/\psi K_S`$: $`\mathrm{sin}2\beta _{\varphi K_S}\mathrm{sin}2\beta _{J/\psi K_S}+O(\lambda ^2)`$ worah .

Recently both BaBar and Belle reported the branching ratio and CP asymmetries in the $`B_d\varphi K_S`$ decay:

$`𝒜_{\varphi K}(t)`$ $``$ $`{\displaystyle \frac{\mathrm{\Gamma }(\overline{B}_{\mathrm{phys}}^0(t)\varphi K_S)\mathrm{\Gamma }(B_{\mathrm{phys}}^0(t)\varphi K_S)}{\mathrm{\Gamma }(\overline{B}_{\mathrm{phys}}^0(t)\varphi K_S)+\mathrm{\Gamma }(B_{\mathrm{phys}}^0(t)\varphi K_S)}}`$ (1)
$`=`$ $`C_{\varphi K}\mathrm{cos}(\mathrm{\Delta }m_Bt)+S_{\varphi K}\mathrm{sin}(\mathrm{\Delta }m_Bt),`$
where $`C_{\varphi K}`$ and $`S_{\varphi K}`$ are given by

$$C_{\varphi K}=\frac{1|\lambda _{\varphi K}|^2}{1+|\lambda _{\varphi K}|^2},\mathrm{and}S_{\varphi K}=\frac{2\mathrm{Im}\lambda _{\varphi K}}{1+|\lambda _{\varphi K}|^2},$$
(2)
with

$$\lambda _{\varphi K}e^{2i(\beta +\theta _d)}\frac{\overline{A}(\overline{B}^0\varphi K_S)}{A(B^0\varphi K_S)},$$
(3)
and the angle $`\theta _d`$ represents any new physics contributions to the $`B_d\overline{B_d}`$ mixing angle. They find $`\mathrm{sin}2\beta _{\varphi K}=S_{\varphi K}=0.39\pm 0.41`$ babar ; belle , which is a 2.7 $`\sigma `$ deviation from the SM prediction: $`\mathrm{sin}2\beta _{J/\psi K_S}=0.734\pm 0.054`$ nir . On the other hand, the measured branching ratio $`(89)\times 10^6`$ is not far from the SM prediction. The direct CP asymmetry in $`B_d\varphi K_S`$ is also reported by the Belle collaboration: $`C_{\varphi K_S}=0.56\pm 0.43`$.

In this paper, we wish to explain the deviation of $`\mathrm{sin}2\beta _{\varphi K_S}`$ from $`\mathrm{sin}2\beta _{J/\psi K_S}`$ within general SUSY models with $`R`$parity conservation. (See Refs. moroi  for related works.) There are basically two interesting classes of modifications: (i) gluino-mediated $`bsq\overline{q}`$ with $`q=u,d,s,c,b`$ induced by flavor mixings in the down-squark sector, and (ii) Higgs mediated $`bss\overline{s}`$ in the large $`\mathrm{tan}\beta `$ limit ($`\mathrm{tan}^3\beta `$ at the amplitude level). In the following, we analyze the effects of these two mechanisms on the branching ratio of $`B_d\varphi K_S`$, $`S_{\varphi K}`$, $`C_{\varphi K}`$, $`BX_s\gamma `$ and its direct CP asymmetry, $`B_s^0\overline{B_s^0}`$ mixing, and the correlation of the Higgs mediated $`bss\overline{s}`$ transition with $`B_s\mu ^+\mu ^{}`$. We will deduce that $`LL`$ and $`RR`$ insertions give effects too small to cause an observable deviation between $`S_{\varphi K}`$ and $`S_{J/\psi K}`$, unless the gluino and squarks have masses close to the current lower bounds. Once the existing CDF limit on $`B_s\mu ^+\mu ^{}`$ is imposed on the Higgs mediated $`bss\overline{s}`$, that also provides too small an effect. On the other hand, the down-sector $`LR`$ and $`RL`$ insertions can provide a sizable deviation in a robust way consistent with all other data.

In the general MSSM, the squark and the quark mass matrices are not diagonalized simultaneously. There will be gluino mediated FCNC of strong interaction strength, which may even exceed existing limits. It is customary to rotate the effects so they occur in squark propagators rather than in couplings, and to parametrize them in terms of dimensionless parameters. We work in the usual mass insertion approximation (MIA) mia , where the flavor mixing $`ji`$ in the down type squarks associated with $`\stackrel{~}{q}_B`$ and $`\stackrel{~}{q}_A`$ is parametrized by $`(\delta _{AB}^d)_{ij}`$. More explicitly,

$$(\delta _{LL}^d)_{ij}=\left(V_L^dM_Q^2V_L^d\right)_{ij}/\stackrel{~}{m}^2,\mathrm{etc}.$$

in the super CKM basis where the quark mass matrices are diagonalized by $`V_L^d`$ and $`V_R^d`$, and the squark mass matrices are rotated in the same way. Here $`M_Q^2`$, $`M_D^2`$ and $`M_{LR}^2`$ are squark mass matrices, and $`\stackrel{~}{m}`$ is the average squark mass. Because we are considering $`bss\overline{s}`$ transitions, we have $`(i,j)=(2,3)`$.

Since the decay $`B_d\varphi K_S`$ is dominated by SM QCD penguin operators, the gluino-mediated QCD penguins may be significant if the gluinos and squarks are not too heavy. For $`LL`$ and $`RR`$ mass insertions, the gluino-mediated diagrams can contribute to a number of operators already present in the SM, such as those with LH flavor-changing currents ($`O_{3,\mathrm{},6}`$) and the magnetic and chromomagnetic operators ($`O_{7\gamma }`$, $`O_{8g}`$). On the other hand $`LR`$ and $`RL`$ insertions do not generate $`O_{3,\mathrm{},6}`$. (All operators are explicitly defined in bbns .) But $`RR`$ and $`RL`$ insertions also generate operators with RH flavor-changing currents and so we must consider the operators $`\stackrel{~}{O}_{3,\mathrm{},6}`$, $`\stackrel{~}{O}_{7\gamma }`$ and $`\stackrel{~}{O}_{8g}`$ which are derived from the untilded operators by the exchange $`LR`$. In the SM, the $`\stackrel{~}{O}_{3,\mathrm{},6}`$ operators are absent and the $`\stackrel{~}{O}_{7\gamma ,8g}`$ operators are suppressed by light quark masses, $`m_s`$. But in SUSY, $`\stackrel{~}{O}_{7\gamma ,8g}`$ in particular can be enhanced by $`\stackrel{~}{m}/m_s`$ with respect to the SM, making them important phenomenologically. Of course, these gluino induced FCNC could also affect $`BJ/\psi K_S`$ in principle. But this gold-plated mode for $`\mathrm{sin}2\beta `$ is dominated by a tree level diagram in the SM, and the SUSY loop corrections are negligible.

We calculate the Wilson coefficients corresponding to each of these operators at the scale $`\mu \stackrel{~}{m}m_W`$; explicit expressions are in Ref. kkkpww . We evolve the Wilson coefficients to $`\mu m_b`$ using the appropriate renormalization group (RG) equations rg , and calculate the amplitude for $`B\varphi K`$ using the recent BBNS approach bbns  for estimating the hadronic matrix elements.

The explicit expressions for $`\mathrm{\Gamma }(B\varphi K_S)`$ can be found in Refs. kkkpww ; bbns . In the numerical analysis presented here, we fix the SUSY parameters to be $`m_{\stackrel{~}{g}}=\stackrel{~}{m}=400`$ GeV. In each of the mass insertion scenarios to be discussed, we vary the mass insertions over the range $`\left|\delta _{AB}^d\right|1`$ to fully map the parameter space. We then impose two important experimental constraints. First, we demand that the predicted branching ratio for inclusive $`BX_s\gamma `$ fall within the range $`2.0\times 10^4<B(BX_s\gamma )<4.5\times 10^4`$, which is rather generous in order to allow for various theoretical uncertainties. Second, we impose the current lower limit on $`\mathrm{\Delta }M_s>14.9`$ ps<sup>-1</sup>.

A new CP-violating phase in $`(\delta _{AB}^d)_{23}`$ will also generate CP violation in $`BX_s\gamma `$. The current data on the direct CP asymmetry $`A_{\mathrm{CP}}^{bs\gamma }`$ from CLEO is cleo  $`A_{\mathrm{CP}}^{bs\gamma }=(0.079\pm 0.108\pm 0.022)(1.0\pm 0.030)`$, which is not particularly constraining. Within the SM, the predicted CP asymmetry is less than $`0.5\%`$, so a larger asymmetry would be a clear indication of new physics kn . Remarkably, despite the well-known importance of $`BX_s\gamma `$ for constraining $`(\delta _{AB}^d)_{23}`$, we will find in the following that the $`B\varphi K`$ branching ratio provides an independent constraint on the $`LR`$ and $`RL`$ insertions which is already as strong as that derived from $`BX_s\gamma `$. Where relevant, we will show our predictions for $`A_{\mathrm{CP}}^{bs\gamma }`$.

We begin by considering the case of a single $`LL`$ mass insertion: $`(\delta _{LL}^d)_{23}`$. Our results will hold equally well for a single $`RR`$ insertion. Though the $`LL`$ insertion generates all of the operators $`O_{3,\mathrm{},6,7\gamma ,8g}`$, qualitatively there are substantial cancellations between their SUSY contributions and so the effect of the $`LL`$ insertion is significantly diluted. Scanning over the parameter space as discussed previously, we find that $`S_{\varphi K}>0.5`$ for $`m_{\stackrel{~}{g}}=\stackrel{~}{m}=400`$GeV and for any value of $`|(\delta _{LL}^d)_{23}|1`$, the lowest values being achieved only for very large $`\mathrm{\Delta }M_s`$. If we lower the gluino mass down to $`250`$GeV, $`S_{\varphi K}`$ can shift down to $`0.05`$, but only in a small corner of parameter space. Similar results hold for a single $`RR`$ insertion. Thus we conclude that the effects of the $`LL`$ and $`RR`$ insertions on $`BX_s\gamma `$ and $`B\varphi K`$ are not sufficient to generate a (large) negative $`S_{\varphi K}`$, unless gluino and squarks are relatively light. Nonetheless, their effects on $`B_s\overline{B_s}`$ mixing could still be large and observable, providing a clear signature for $`LL`$ and $`RR`$ mass insertions (see Ref. kkkpww  for details).

Next we consider the case of a single $`LR`$ insertion. Scanning over the parameter space and imposing the constraints from $`BX_s\gamma `$ and $`\mathrm{\Delta }M_s`$, we find $`|(\delta _{LR}^d)_{23}|10^2`$. This is, however, large enough to significantly affect $`B\varphi K_s`$, both its branching ratio and CP asymmetries, through the contribution to $`C_{8g}`$.

In Fig. 1(a), we show the correlation between $`S_{\varphi K}`$ and $`C_{\varphi K}`$. Since the $`LR`$ insertion can have a large effect on the CP-averaged branching ratio for $`B\varphi K`$ we further impose that $`B(B\varphi K)<1.6\times 10^5`$ (which is twice the experimental value) in order to include theoretical uncertainties in the BBNS approach related to hadronic physics. Surprisingly, we see that the $`B\varphi K`$ branching ratio constrains $`(\delta _{LR}^d)_{23}`$ just as much as $`BX_s\gamma `$.

Studying Fig. 1, we can learn much about the viability and testability of scenarios with a single $`LR`$ insertion. For one thing, we can get negative $`S_{\varphi K}`$, but only if $`C_{\varphi K}`$ is also negative (this is excluded by the present data, but only at the $`1.4\sigma `$ level). Conversely, a positive $`C_{\varphi K}`$ implies that $`S_{\varphi K}>0.6`$. This correlation between $`S_{\varphi K}`$ and $`C_{\varphi K}`$ can be tested at B factories in the near future given better statistics. The correlation between $`S_{\varphi K}`$ and the direct CP asymmetry in $`BX_s\gamma `$ ($`A_{\mathrm{CP}}^{bs\gamma }`$) is shown in Fig. 1(b). We find $`A_{\mathrm{CP}}^{bs\gamma }`$ becomes positive for a negative $`S_{\varphi K}`$, while a negative $`A_{\mathrm{CP}}^{bs\gamma }`$ implies that $`S_{\varphi K}>0.6`$. It is also clear from the figure that the present CLEO bound on $`A_{\mathrm{CP}}^{bs\gamma }`$ does not significantly constrain the $`LR`$ model. Finally, we also find that the deviation of $`B_s\overline{B}_s`$ mixing from the SM prediction is small after imposing the $`B_dX_s\gamma `$ and $`B_d\varphi K_S`$ branching ratio constraints. Thus we conclude that a single $`LR`$ insertion can describe the data on $`S_{\varphi K}`$ but only if the measurement of $`C_{\varphi K}`$ shifts to negative values. This scenario can then be tested by measuring a positive direct CP asymmetry in $`BX_s\gamma `$ and $`B_d`$-$`\overline{B}_d`$ mixing consistent with the SM.

The last single mass insertion case is that of $`RL`$. We find that the $`RL`$ operator most conveniently describes the current data. For example, in Figs. 2(a) and (b), we show the correlation of $`S_{\varphi K}`$ with $`B(B\varphi K)`$ and $`C_{\varphi K}`$, respectively, in the $`RL`$ dominance scenario. We find that $`S_{\varphi K}`$ can take almost any value between $`1`$ and $`+1`$ without conflict with the observed branching ratio for $`B\varphi K`$. We also find that $`C_{\varphi K}`$ can be positive for a negative $`S_{\varphi K}`$, unlike in the $`LR`$ case. In particular, if we impose the $`B(B_d\varphi K^0)`$ constraint, then $`|C_{\varphi K}|>0.2`$ for $`S_{\varphi K}<0`$, except in very small regions of parameter space. This could be a useful check on the presence of the $`RL`$ operator. But the $`RL`$ operator leaves SM predictions intact in two important observables. First, within the range of allowed $`(\delta _{RL}^d)_{23}`$, the $`B_s`$-$`\overline{B}_s`$ mixing receives only very small contributions from SUSY. Second, because the $`RL`$ insertion generates the $`\stackrel{~}{O}_{3,\mathrm{},6,7\gamma ,8g}`$ operators while the SM generates only the untilded operators, there is no interference between the phases of the SUSY and SM contributions, and so $`A_{\mathrm{CP}}^{bs\gamma }=0`$.

The process $`BX_s\gamma `$ actually has a rather unconventional but well-motivated structure in the $`RL`$ scenario. In Ref. kane , it was found that one can have a strong cancellation between the SM and SUSY contributions to $`C_{7\gamma }`$ and $`C_{8g}`$, as is expected in the supersymmetric limit. Then the observed branching ratio is due to the $`\stackrel{~}{C}_{7\gamma }`$ operator, which implied an $`RL`$ insertion of about the same size as that needed here. Thus this earlier analysis motivates the $`RL`$ insertion for the present case.

We have considered separate insertions so the effects of each can be seen clearly. In many realistic SUSY models, combinations of various insertions will play a role phenomenologically, and of course we do not mean to imply that only single insertions should occur in nature. A particularly interesting case is one in which the $`RL`$ and $`RR`$ operators co-exist. While we do not study this model in detail here, we note that (as discussed in Ref. kane ) if both $`(\delta _{RR}^d)_{23}`$ and $`(\delta _{RL}^d)_{23}`$ are nonzero, then $`A_{\mathrm{CP}}^{bs\gamma }`$ can arise as an interference between $`\stackrel{~}{C}_{7\gamma }`$ and $`\stackrel{~}{C}_{8g}`$. Then the $`RR`$ insertion could shift $`\mathrm{\Delta }M_s`$ appreciably, while the $`RL`$ insertion gives the above result.

Finally we consider a completely different class of contributions to $`bss\overline{s}`$. At large $`\mathrm{tan}\beta `$, FCNC’s can also be mediated by exchanges of neutral Higgs bosons, as found in $`BX_s\mu ^+\mu ^{}`$ huang , $`B_{s,d}\mu ^+\mu ^{}`$ kolda , $`\tau 3\mu `$, etc. Babu:2002et . Likewise, $`bss\overline{s}`$ could be enhanced by neutral Higgs exchange, which is flavor dependent since the Higgs coupling is proportional to the Yukawa couplings. The effective coupling for $`bss\overline{s}`$ is basically the same with that for $`B_s\mu ^+\mu ^{}`$ up to a small difference between the muon and the strange quark masses (or their Yukawa couplings). Imposing the upper limit on this decay from CDF cdf  during the Tevatron’s Run I, $`B(B_s\mu ^+\mu ^{})<2.6\times 10^6`$, we can derive a model-independent upper limit on this effective coupling. We find that $`S_{\varphi K_S}`$ cannot be smaller than 0.71 for such couplings, so that we cannot explain the large deviation in $`S_{\varphi K_S}`$ with Higgs-mediated $`bss\overline{s}`$ alone.

Now let us provide possible motivation for values of $`(\delta _{LR,RL}^d)_{23}10^2`$ that we find phenomenologically are needed in order to generate $`S_{\varphi K}`$. In particular, at large $`\mathrm{tan}\beta `$ it is possible to have double mass insertions which give sizable contributions to $`(\delta _{LR,RL}^d)`$. First a $`(\delta _{LL}^d)_{23}`$ or $`(\delta _{RR}^d)_{23}10^2`$ is generated. The former can be obtained from renormalization group running even if its initial value is negligible at the high scale. The latter may be implicit in SUSY GUT models with large mixing in the neutrino sector moroi . Alternatively, in models in which the SUSY flavor problem is resolved by an alignment mechanism using spontaneously broken flavor symmetries, or by decoupling, the resulting $`LL`$ or $`RR`$ mixings in the $`23`$ sector could easily be of order $`\lambda ^2`$ align ; decoupling . However as discussed above, this size of the $`LL`$ and/or $`RR`$ insertions cannot explain the measured CP asymmetry in $`B_d\varphi K_S`$. But at large $`\mathrm{tan}\beta `$, the $`LL`$ and $`RR`$ insertions can induce the $`RL`$ and $`LR`$ insertions needed for $`S_{\varphi K}`$ through a double mass insertion bjkp :

$$(\delta _{LR,RL}^d)_{23}^{\mathrm{ind}}=(\delta _{LL,RR}^d)_{23}\times \frac{m_b(A_b\mu \mathrm{tan}\beta )}{\stackrel{~}{m}^2}.$$

One can achieve $`(\delta _{LR,RL}^d)_{23}^{\mathrm{ind}}10^2`$ if $`\mu \mathrm{tan}\beta 10^4`$GeV, which could be natural if $`\mathrm{tan}\beta `$ is large (for which $`A_b`$ becomes irrelevant). Note that in this scenario both the $`LL(RR)`$ and $`LR(RL)`$ insertions would have the same CP violating phase, since the phase of $`\mu `$ here is constrained by electron and down-quark electric dipole moments. Lastly, one can also construct string-motivated $`D`$-brane scenarios in which $`LR`$ or $`RL`$ insertions are $`10^2`$ kkkpww .

In this letter, we considered several classes of potentially important SUSY contributions to $`B\varphi K_S`$ in order to see if a significant deviation in its time-dependent CP asymmetry $`S_{\varphi K_S}`$ could arise from SUSY effects. The Higgs-mediated FCNC effects and models based on the gluino-mediated $`LL`$ and $`RR`$ insertions both give contributions too small to alter $`S_{\varphi K}`$ significantly. On the other hand, the gluino-mediated contribution with $`LR`$ and/or $`RL`$ insertions can lead to sizable negative in $`S_{\varphi K_S}`$ (as reported experimentally) as long as $`|(\delta _{LR,RL}^d)_{23}|10^310^2`$. As a byproduct, we found that nonleptonic $`B`$ decays such as $`B\varphi K`$ begin to constrain $`|(\delta _{LR,RL}^d)_{23}|`$ as strongly as $`BX_s\gamma `$. Besides producing no measurable deviation in $`B^0`$-$`\overline{B}^0`$ mixing, the $`RL`$ and $`LR`$ operators generate definite correlations among $`S_{\varphi K}`$, $`C_{\varphi K}`$ and $`A_{\mathrm{CP}}^{bs\gamma }`$, and our explanation for the negative $`S_{\varphi K}`$ can be easily tested by measuring these other observables. In particular $`C_{\varphi K}`$ can be positive only for an $`RL`$ insertion. A pure $`RL`$ scenario also predicts a vanishing direct CP asymmetry in $`BX_s\gamma `$. However, if the $`RL`$ insertion is accompanied by an $`RR`$ insertion, then the resulting direct CP asymmetry in $`BX_s\gamma `$ can be large kane  and the $`RR`$ contribution to the $`B_s\overline{B_s}`$ mixing can induce significant shifts in $`\mathrm{\Delta }M_s`$ and the phase of $`B_s\overline{B}_s`$ mixing. Finally, we also point out that the $`|(\delta _{LR,RL}^d)_{23}|10^2`$ can be naturally obtained in SUSY flavor models with double mass insertion at large $`\mathrm{tan}\beta `$, and in string-motivated models kkkpww .

###### Acknowledgements.

We thank S. Baek, H. Davoudiasl, B. Nelson, K. Tobe and J. Wells for helpful comments, and A. Kagan for pointing out the sign error in the corrleation between $`S_{\varphi K}`$ and $`A_{\mathrm{CP}}^{bs\gamma }`$. PK and JP are grateful to MCTP at Univ. of Michigan for the hospitality. This work is supported in part by KRF grant KRF-2002-070-C00022 and KOSEF through CHEP at Kyungpook National University (PK and JP), by National Science Foundation under grant PHY00–98791 (CK), by Department of Energy (GK and HW), by a DOE grant No. DE-FG02-95ER40896 and in part by the Wisconsin Alumni Research Foundation (LW).


format: text

text:
# References

A Market of Inhomogeneous Threshold Cellular Automata

Dietrich Stauffer\* and Gérard Weisbuch\**

\*Department of Physics, Indian Institute of Technology, Kanpur 208016, India; visiting from Institute for Theoretical Physics, Cologne University, D-50923 Köln, Euroland

\**Laboratoire de Physique Statistique, ENS, 24 rue Lhomond, F-75231 Paris Cedex 5, Euroland

This article summarizes some physics aspects of the market model of Weisbuch and Stauffer, Physica A (2003): How do demand and quality expectation adjust to each other in a buyer dominated market like cinema visits? And how can business cycles be modelled?

“Social percolation”  applied percolation theory to a marketing problem: How to quality expectations of potential buyers and the quality of a product adjust to each other ? The model lead to self-organized criticality at the percolation threshold. Thus fluctuations were large, some products were a big hit, most were flops and lots of products were somewhere in the middle. Fig.1a  shows one of the simulations, and Fig.1b a corresponding result from a real market (cinema films; see also ). There is qualitative agreement; note that the simulation rank plot extends over four decades, the one with real data only over two decades, in vertical direction.

For this “threshold” conference, instead of this percolation model we apply threshold cellular automata to this marketing problem: How many potential buyers actually buy ? Different people have different opinions, and thus our $`N`$ automata have different thresholds $`\theta _i,i=1,2,\mathrm{}N`$. This is what is meant with “inhomogeneous”. Each site $`i`$ on a large square lattice carries besides the threshold $`\theta _i`$ a spin $`S_i`$ which is +1 for buying and –1 for not buying. The spin orientation at the next time step is fully determined by the sum $`Q=_kS_k`$ of the four neighbour spins $`k`$ of $`i`$: If this “quality” $`Q`$ is larger than the threshold $`\theta _i`$, site $`i`$ buys; otherwise it does not buy. The magnetization $`M=_iS_i/N`$ is the normalized difference between the numbers of buyers and non-buyers. We start with all spins up and then use random sequential updating of the spins, with one time step meaning that on average every site was updated once.

As long as the thresholds do not vary with time, the above problem is similar to a zero-temperature simulation of a random-field Ising model. However, after going to a movie you may not go immediately to another one except when it is very good; in contrast, when you have not seen a good movie since a long time you may visit even a mediocre one. Thus $`Q`$ is moved up if $`S_i`$ now is +1, and moves down if $`S_i=1`$. Since the quality $`Q`$ as reported from the neighbours can take only the values $`4,2,0,2,4`$, a threshold 2.1 is equivalent to a threshold 3 or 3.9. Thus for simplicity we take our thresholds as odd integers $`5,3,1,1,2,3`$, initially distributed randomly, and change them at each time step by $`2S_i=\pm 2`$. When the threshold has reached +5 the next updating of this spin will change $`S_i`$ to –1 and $`\theta _i`$ to 3, since $`Q`$ cannot reach 5. (If one allows the thresholds the change by only a small amount, the dynamics gets slow and is described by a scaling law ; we ignore this complication here.)

This model thus describes a market dominated by demand of buyers and assuming sufficient supply of sellers. The judgment of the buyers is determined by what their neighbours do; if a sufficient number of neighbours buy, then also the central person buys. Which number of buying neighbours is sufficient in this sense depends on the time-dependent threshold $`\theta _i(t)`$. This threshold increases by 2 (i.e. by one more neighbour) whenever the person buys, and decreases by the same amount when the person does not buy.

Results:

Fig.2 shows that in a stationary equilibrium the buyers are not randomly distributed but form small clusters. The cluster size is much smaller than the lattice size (usually one million sites); in contrast to social percolation , we do not observe self-organized criticality. The fluctuations in the stationary magnetization decrease as $`1/\sqrt{N}`$ as usual: Fig.3. The magnetization ($`M=1`$ initially) decays to zero with damped oscillations, in one (Fig.4), two (Fig.5), three (Fig.6) and four (Fig.7) dimensionsions (hypercubic lattices). (Also the triangular lattice behaves similarly; not shown).

When you have drunken one bottle of wine, you should not immediately buy and drink the next one. Thus a site which just bought should wait $`\tau `$ time steps before considering to buy again. During this waiting time its spin stays at its value and counts as buyer. Analogously, after a decision not to buy, the spin stays at $`1`$ for $`\tau `$ iterations and counts as not buying, to preserve up-down symmetry. The previous description thus is the special case $`\tau =0`$. Then the oscillations are damped much less, as can be seen from Figs. 4-6. The amplitude of the oscillations decays exponentially with time at intermediate time scales (not shown).

Durable goods like houses are often bought only once in life. Setting in some sense $`\tau =\mathrm{}`$ but sacrificing the up-down symmetry kept until now, we assume that once a site has bought in never buys again. We still start with all $`S_i=1`$ and $`\theta _i`$ random. Fig.8 then shows the numbers of actually still buying people; this number initially may vary non-monotonically with time but later decays exponentially with time until the whole market is saturated forever.

When we go back to the above reversible model with $`\tau =0`$ and start with randomly oriented spins instead of all spins parallel, then the various clusters oscillate with different phases, showing oscillations for the overall magnetization $`M(t)`$ only for small lattices. In that case it is better to observe the autocorrelation function $`<M(T+t)M(T)>`$ versus $`t`$ for large times $`T`$ such that equilibrium is established; this function then looks similar to the $`M(t)`$ of Fig.3. More information, including Fourier transforms, on this way of initialization are given in .

At the conference, many suggestions were made by the audiance, like putting this model on a Barabási-Albert or small-world network, allowing for quenched or annealed solution, tune the model to the percolation threshold. J. Kertész suggested to put in a fraction $`p`$ of antiferromagnetic bonds, i.e. mistrust of neighbours: If this neighbour likes the movie it must be bad. We found that for $`\tau =0`$ the oscillations vanish if at least half of the bonds are antiferromagnetic, while for $`\tau 2`$ one overshooting takes place even if all bonds are antiferromagnetic.

DS thanks a DFG/BMZ grant for Indian-German collaboration for supporting this visit.


text:
# 1 Introduction

## 1 Introduction

If Quantum chromodynamics (QCD) is correct there is no doubt that under conditions that exist in Lattice Gauge Theory (LGT) calculations, a large volume quark-gluon plasma (QGP) is expected to be created. However one can question whether high energy heavy ion collisions in the RHIC collider reproduce the LGT calculation conditions well enough, that production of a detectable large volume of QGP occurs. This has more or less been assumed in many theoretical calculations. However, to our knowledge no-one has shown that the actual dynamics existing at RHIC would allow this to occur with detectable cross-sections. Two Lorentz contracted heavy ions pass through each other in the short times available. There are turbulent ever changing dynamics of the environment resulting from their interaction. This does not give one assurance that even the very basic LGT requirements of thermal and chemical equilibrium can be met. Certainly this has not been demonstrated by any theoretical work we are aware of. Therefore from the early eighties onward, as documented in the first four RHIC experimental workshops from 85-90, we have concluded this is not likely to happen.

Instead our view has been that if QGP is created in a high energy heavy ion collider such as RHIC or LHC, it is more probable that local fluctuations would produce one to many droplets (clusters, chunks or bubbles) of QGP. These would possibly be detectable, especially if they were localized in phase space. Leon van Hove had this view, and did string model calculations, which resulted in small droplets of QGP being formed by the breaking of stretched strings. These QGP bubbles were localized in rapidity, and gave rise to rapidity bumps or peaks in the rapidity or psuedo-rapidity distribution.

We have previously published a treatment of the single (spherical-no longitudinal expansion) bubble case (similar to the Van Hove type), which serves as an introductory paper for the general subject. It is possible that with enough statistics one could in principle find single to a few bubble events. We have concluded that multiple bubble formation (which we mentioned but did not do calculations for in Ref. ) is the more probable general case. Therefore this case should be treated for realistic attainable statistics. Of course a large number of bubbles in the multiple bubble case will obscure the resolution and observation of single bubble phenomena such as rapidity bumps etc. Furthermore the overall result would likely appear similar to a thermal model. However, this likely occurs because particles from different parts of space go into the same phase space. Assume all the bubbles were Van Hove spherical bubbles at rest (i.e. located at mid rapidity). Even though they are at different space points, they would add up to be equivalent to one big Van Hove spherical bubble. If we now give motion to each bubble along the beam direction, they will smear out in momentum space. Therefore to detect the effects of the multiple bubble case, we must carefully choose a region where phase space focusing will allow the addition of multiple bubble effects. Thus we must find a part of momentum space that is highly correlated to position space. If we could force bubbles to remain at rest, then their particles would add together forming a rapidity bump at $`\eta =0`$. Even if this happened it would be difficult to get a clean signal over background in this region, because background particles from soft fragmentation end up at central rapidity. Thus it is important to find a phase space region free of such background.

At RHIC the pre-hadronic matter is being pushed in the transverse direction building up transverse momentum $`p_t`$. Bubbles that are pushed along with this flow will hadronize into particles focused over a limited range of angles. In this paper we will model such bubbles, and address backgrounds which mimic bubble effects. These are mainly jets and to some extent resonance’s.

The above states the goals of this paper. We will treat the case of an approximately maximum number of multiple bubbles in one outer ring around the blast region. These bubbles contribute to the final hadronization of particles coming from the QGP. This may not occur in the actual case. Therefore only future data analysis, can shed light on this. We will find that a reasonable theoretical treatment can lead to methods for detecting QGP bubbles.

In this paper we are essentially limiting ourselves to simulations and analysis suitable for comparison with forthcoming RHIC data. We will utilize relevant parts of the considerable body of data that has already been published. It is expected that in the near future data from RHIC for 65GeV/n $`Au`$ colliding with 65 GeV/n $`Au`$ could test these ideas. However, similar methods could be applied to higher energies at RHIC and LHC.

## 2 General Considerations

At RHIC the pre-hadronic matter is being pushed outward in the transverse direction. Particles with higher transverse momentum ($`p_t`$) are pushed more than particles with lower $`p_t`$. Analysis of pions by Hanbury Brown-Twiss (HBT) methods have shown that at low $`p_t`$ the source size is about 6 fm, while above 0.8 GeV/c the source size becomes about 2 fm. This allows phase space focusing to form a reasonable signal in the 2 fm source size region. These measurements imply that the viewed region of the initial position space becomes smaller as one selects higher $`p_t`$ particles. Pions at a $`p_t`$ $``$ 1-2 GeV/c will be coming from the outer regions of the expanding fireball in regions where the HBT original source size radial width is of the order of 2 fm. Softer pions will mainly have a radius of within about 6 fm. This supports a rough estimate of $``$ 6+2 fm to be the transverse radius of the fireball that is emitting hadronization particles. One should keep in mind that this is a quantum mechanical system with dynamical and turbulent changes. Therefore the previous arithmetic, and subsequent arithmetic is to be considered in the sense of very crude estimates.

One should note that particles above 2 GeV/c will likely have jets as a sizeable source of contamination. Thus we believe we can work with a window in $`p_t`$ of 0.8 to 2 GeV/c, to search for multiple bubble effects. We chose the upper end of the window to avoid the hard physics region above 2 GeV/c, which would increase the background for the effects we are looking for. The lower end of the range is chosen to maintain the space momentum correlation of our signals, and therefore enhance them.

Sometimes we vary the window somewhat to investigate a particular point. From Ref.  we found that the single plasma bubbles have a mean $`p_t`$ of about 0.5 GeV/c. Thus the azimuthal angular range in $`\mathrm{\Delta }\mathrm{\Phi }`$ can be crudely estimated for the bubbles in our window. With an average $`p_t`$ of about 1 GeV/c, and the above right angle momentum, we form an $``$ 30$`{}_{}{}^{c}irc`$ angle. Thus we can assume an approximation that spherical bubbles have an angular range in $`\mathrm{\Delta }\mathrm{\Phi }`$ of about 30 degrees. This however is very similar to the angular spread of jet fragmentation, making it difficult to separate the spherical bubbles $`\mathrm{\Delta }\mathrm{\Phi }`$ distribution from the $`\mathrm{\Delta }\mathrm{\Phi }`$ distribution of jet fragmentation. It might be noted that there are arguments for jet quenching which would improve our signal compared to the assumed background. However we will ignore jet quenching in our simulations in order to be very conservative in drawing our conclusions.

We know there is a longitudinal expansion-the Landau fireball effect. The value of this longitudinal expansion has to be determined from analyzing the data. However the Landau longitudinal expansion of our bubbles used in our previous paper was probably too large. Therefore we will choose a reasonable value intermediate between that and the value for a spherical bubble (which has zero longitudinal expansion) for our simulations.

Figure 1 contains a sketch of the assumed bubble geometry, and details of how bubbles are embedded are given in the caption. In the language of Van Hove, the string stopping after breaking is not complete, so that longitudinal expansion is left in the strings or bubbles. The longitudinal expansion will increase the angular spread of $`\mathrm{\Delta }\eta `$ due to bubbles, which will distinguish bubbles from jet fragmentation. However it will spread out the bubble signal, but this may not matter. Increasing the energy in the bubble, which also reasonably could occur, would enhance the bubble signal. This leads to our using an about 50 range for the bubbles in the psuedo-rapidity direction for our simulations. The angular range of the bubbles in the azimuthal direction is about 30, as previously discussed.

From HBT work, previously referred to, we can estimate the bubble would have a radius of about 2 fm, and thus has a diameter of about 4 fm. The rapid transverse expansion in the blast region pushes high density pre-hadronic matter from the central regions outward to where it hadronizes. We assume a single outer ring of bubbles at the outer circumference of the blast region would be filled with bubbles at hadronization. The bubbles provide the hadronization coming from the QGP. An angular range of 30 is one twelfth of 360. Since each bubble has a diameter of about 4 fm, the circumference to cover the entire azimuth, would be approximately 48 fm. A circumference of 48 fm implies a radius of about 8 fm for the ring of bubbles. The number of 8 fm is consistent with our HBT picture presented previously. Inside this outward shell of bubbles there may be other overlapping bubbles, but the ring we have chosen will predominantly contribute to the mid rapidity region.

Using the work in Ref.  that bubbles of a 2 fm radius (one of the choices) bubble size for a RHIC event would have about 40 domains (which we call bubbles) with energies of about 3 GeV per bubble going into charged pions. We calculated the energy per domain (bubble) using the information in their paper. For our simulations we felt it was reasonable to use the 2 fm bubble size in employing the methods in  since it was the choice we arrived at from the HBT work. We also decided that in this first simulation it was reasonable to use our calculation (based on) of about 3 GeV per bubble hadronizing into charged pions, in order to avoid being arbitrary. Obviously this value has to be considered a parameter which could be determined in conjunction with data analysis when it becomes available. However it should be noted we are not using any of the other detailed work in Ref. .

We used an average of 13 bubbles in a ring at approximately mid-rapidity in each central 65 GeV/n $`Au`$ on 65GeV/n $`Au`$ event. This fills the ring of bubbles, whose contributions would dominate what is seen at central rapidities at the RHIC STAR Detector. This detector would be the most likely near future source of experimental results to check these ideas. We used an average of 3.25 GeV per bubble. The energy was increased from 3.0 GeV per bubble, since we are producing more than pions going into charged particles. This led to an average of 1.95 charged particles going into the cuts we will use.

In our bubble scenario each QGP bubble is an uncharged, gluon dominated, color singlet system. Thus when the bubble hadronizes the total charge of the particles coming from the bubble is zero. Since we are selecting a $`p_t`$ range where we expect the bubble concentration to be rich, we should therefore see a suppression of charge fluctuations. This occurs because the charge fluctuations coming from a localized QGP bubble should be less by a factor $``$4, than charge fluctuations coming from an ideal pion gas. This subject will be addressed in the next section.

## 3 Charge Suppression Effects Due to QGP Formation

We now address recent papers on charge fluctuation suppression calculations, of the ratio of positively charged and negatively charged pions as a signal for QGP formation. In the letters of S. Jeon, and V. Koch, and M. Asakawa, U. Heinz, and B. Mueller, they concluded that a parameter $`D=\frac{4<\mathrm{\Delta }Q^2>}{<N_{ch}>}`$ for mesons evaluated for event by event charge fluctuations, is $``$4 for an ideal pion gas not produced by QGP. They estimate it would be approximately 1 for pions originating from a QGP both from LGT or ideal gas calculations. Correcting for resonance effects in actual observations, they concluded the observed $`D`$ would become 3-2 respectively. Therefore they argued observing these reduced charge fluctuations would serve as a distinct signal for QGP production. We pointed out that their method of treatment allowed color charge fluctuations, and kinematic mixing effects, which they overlooked. These color charge fluctuations, and kinematic mixing effects were shown by us to be important, and model dependent. They could modify the result by washing out or even entirely eliminating the charge suppression effects predicted. In regard to the kinematic mixing effects, it should be noted that in their treatments, even if locally in space one has a charge fluctuation suppression of $``$ 70% to 50%, kinematic mixing of position space into momentum space causes electric charge fluctuations to increase over a wide region, thus reducing suppression. Cuts on the particles one measures will also reduce the calculated suppression. If one gains particles from other parts of position space, that reduces the calculated suppression. Thus we concluded that it is unlikely their arguments could support their predictions for observation of electric charge fluctuation suppression, even if large volumes of QGP were produced as they assumed.

## 4 Electric Charge Fluctuation Suppression in the Bubble Scenario

In our bubble scenario each QGP bubble is a very localized color singlet and uncharged system. Thus when the bubble hadronizes the total charge of the particles coming from the bubble has to be approximately zero. Therefore due to the localization, color charge fluctuations, and most of the kinematic mixing effects which can drastically change the suppression predicted in , are negligible. Since we are selecting a $`p_t`$ range where we expect the bubble effects to be substantial, we should see a suppression of charge fluctuations in the charged particles coming from the QGP bubbles . The large predicted reduction by a factor of $``$4 of charge fluctuations which is reduced to $``$2-3 by resonance’s in, will not be achieved by our measurements because we take account of the presence of background particles, and the loss of some of the plasma particles out of our cuts. However there will be predicted a measurable reduction, and it might be considerable.

In a study of resonance effects, we will find that conclusions from our simulations for our chosen signals are not significantly affected by resonance background. Thus we expect under the bubble scenario that we are following, there will be predicted observable charge suppression. This will be good evidence for QGP formation if observed. One needs to await analysis and availability of the relevant data to test, and if it appears relevant, to optimize these ideas.

## 5 HIJING Based Model

Our first objective is to construct a model which will hopefully take account of the most important effects due to QGP bubbles, background, and resonance effects etc. This would allow a direct comparison of the model predictions and the future RHIC data in a reasonably quantitative manner. We now make a model based on the HIJING event generator. For $`Au`$ on $`Au`$ at 65 GeV/c per nucleon. HIJING is a good choice to base a simulation on since HIJING has been successfully used to fit, and help with the analysis of RHIC data in numerous instances. However HIJING has an important missing part for our $`Au`$ on $`Au`$ simulation. That is elliptic flow which has been measured at RHIC. Not taking account of elliptic flow would unrealistically modify our simulation results. Therefore we have modified HIJING to include relevant elliptic flow effects. HIJING has two relevant sources of particle production: Jets which fragment into particles which are referred to as jet particles, and the soft particles which come from beam jet fragmentation. The jet particles are not flat in azimuth but bunch around the jet axis. The beam jets fragmentation’s are very flat in azimuth.

To take account of the observed elliptic flow we modify the distribution in azimuth of the soft particles (beam jets) so that we develop a $`\mathrm{cos}2\mathrm{\Phi }`$ component about a fixed axis for each $`Au`$ on $`Au`$ simulated event. Into each central $`Au`$ on $`Au`$ event we have added ,on the average, 13 adjoining bubbles in a single ring in the central rapidity region to replace the mini-jets. We are assuming in essence that region is perhaps the source of bubble production. Each bubble contributes from 1 to 4 charged particles to the $`\eta `$ range of +0.75 to -0.75, with $`p_t`$ greater than 0.8 GeV/c and less than 2.0 GeV/c. This $`p_t`$ cut has its lower bound chosen to maintain the space momentum correlation which enhances our signal, while the upper bound is chosen to avoid contamination from the hard physics region above 2 Gev/c, and it is our most relevant $`p_t`$ cut. However some of the time we use 1.2 GeV/c for the lower bound to investigate and separate various effects as indicated on some figures and in the text.

The total charge of each bubble was set to zero, which is appropriate for a QGP bubble. We then generated 100,000 bubble events with impact parameter ranging from 0.0 to 4.0 fm. We also generated 100,000 events of our modified HIJING taking into account relevant elliptic flow effects.

From this point onward, we have made our assumptions on bubble geometry, and how to embed them. Therefore our efforts in the remainder of the paper are devoted to how to separate background and resonance effects etc. Thus allowing us to detect the effects of the bubbles in the RHIC experiment cited.

In Fig. 2 we show the $`\mathrm{\Delta }\mathrm{\Phi }`$ correlation generated by the above simulations, including the embedded bubbles and relevant elliptic flow in our modified HIJING. For comparison we make use of our modified HIJING which has the beam jets modified by elliptic flow, and the expected mini-jets predicted by HIJING.

For $`Au`$ on $`Au`$ with an energy of 65 GeV per nucleon we used the standard 2 GeV/c QCD cutoff. One obtains an average of 17.6 jets per event from which 13.3 charged particles fall into our cuts (Fig. 2). It appears the STAR detector at RHIC is the best bet for experimentally checking in detail the theoretical calculations in this paper. Therefore we have added two track merging that one sees in the STAR detector. This inefficiency only effects small $`\mathrm{\Delta }\mathrm{\Phi }`$’s. In Fig. 2, we see that the plasma bubbles have a stronger correlation than the standard mini-jets. In order to make quantitative comparisons for angles less than 60 where we expect the bubbles to contribute, we calculate $`\chi ^2`$ in all angular correlation calculations with this cut. In Fig. 2 the $`\chi ^2`$ is 53 for 8 bins (Degrees of Freedom). This is an order of $``$8$`\sigma `$ effect. The theory predicting the number of mini-jets is not exact, and we wish to be conservative. Therefore we ask the question, by how much of a factor do we have to increase the HIJING predicted jets to make the bubble effect difference in $`\mathrm{\Delta }\mathrm{\Phi }`$ become indistinguishable from HIJING with the added jets.

In Fig. 3 we show that from an 100,000 event simulation, that arbitrarily doubling the number of mini-jets in HIJING causes the $`\chi ^2`$ to drop to 6 resulting in no difference between bubbles, and arbitrarily increasing the number of mini-jets generated by HIJING to double. This represents a very conservative, and probably an excessively overdone approach, especially since reasonable arguments exist that actually jet quenching rather than enhancement occurs. The agreement between HIJING plus bubbles and HIJING plus double jets is also good if we choose a tighter cut with a $`p_t`$ range 1.2 to 2.0 GeV/c (see Fig. 4, for which $`\chi ^2=10`$).

If we bin in $`\mathrm{\Delta }\eta `$ and plot HIJING plus bubbles, and HIJING plus double jets, we see the correlation for the bubbles compared to the jets in Figs. 5-8 for various $`\mathrm{\Delta }\eta `$ ranges. The $`\chi ^2`$ for the 4 Figs. are 1, 4, 14, 1 for 3DF. Only Fig. 7, shows a difference, with some statistical significance of over 3$`\sigma `$ for particles which are separated by $`\mathrm{\Delta }\eta `$ of 1.05. The difference is that the correlation is wider for the bubbles compared to the jets. This difference in width along the $`\eta `$ direction, is expected from the Landau longitudinal expansion of the bubbles. With more statistics in the simulation, we can expect that the differences in width along the $`\eta `$ direction would become more evident, and have better statistical significance due to the nature of the effects longitudinal expansion produces.

Let us look at an angular correlation of the angle between particles (opening angle $`\mathrm{cos}\mathrm{\Theta }`$). Figure 9 indicates that the bubbles have a wider correlation than the double jets. The $`\chi ^2`$ for Fig. 9 is 26 for 8 DF which is also $`3\sigma `$. Next let us do $`\mathrm{\Delta }\mathrm{\Phi }`$ correlations for like and unlike charges separately. In Fig. 10 we show this correlation for the $`p_t`$ range 0.8 to 2.0 GeV/c. We see that the difference between the unlike and like charges $`\mathrm{\Delta }\mathrm{\Phi }`$ correlation is larger for the bubbles than for double jets. The $`\chi ^2`$ is 70 for DF=24 (3X8=24) which is a 5$`\sigma `$ effect. This difference is due to the zero charge of the bubbles while jets only have a reduced charge. This represents the charge fluctuation suppression effects we are looking for. We can form a measure of charge fluctuations by looking at the charge difference, event by event, for particles which lie in our $`p_t`$ range (0.8 to 2.0 GeV/c), and $`\eta `$ range ($`|\eta |<0.75`$). For the mean of the charge difference, with our average of 108 particles per event, we get 4.0 positive charges per event. The width for the double jets is 10.4 particles. The square root of 108 is 10.4. Thus the width for the double jets is consistent with a purely random charge fluctuation result. The width for the bubbles is 9.7 particles thus being consistent with 95% of a random charge fluctuation result. When pairs of charges are created and go into the $`p_t`$ window (0.8 to 2.0 GeV/c) then the charge fluctuations are reduced. However when one charge goes into the window and one outside there is a random addition.

Since we see a net positive charge this implies baryon transport to the central region. We are summing over impact parameters of 0.0 to 4.0 fm, and thus have a varying fluctuation of baryon transport. For more central events the net positive charge would be larger, and for the less central events this positive charge would be smaller. This effect causes a larger than random charge fluctuation. It also appears that these effects cancel out for our HIJING simulation with the double jets.

The bubble fluctuations appear somewhat smaller. This is because not all the bubble are contained in the above cuts. Pions from the bubble end up having $`p_t`$ near the lower edge of the $`p_t`$ range. The kaons are boosted to the mid-range, and protons are near the upper range. As stated previously. The upper end of the range is chosen to avoid the hard physics region above 2.0 GeV/c, and the lower end of the range is chosen to maintain the space momentum correlation.

In the future the planned Time of flight system which surrounds the central TPC at RHIC is expected to become available. We will then be able to do much more detailed treatments similar to the above. We can analyze $`\pi ^+`$ and $`\pi ^{}`$, $`K^+`$ and $`K^{}`$, and $`p`$ and $`\overline{p}`$, instead of just using, all positive and all negative charge pairs. Eventually when there is enough experimental statistics available to compare with, we can do multiple particle correlation’s for larger numbers of particles.

## 6 Estimating Resonance Effects with a Resonance Gas Model

In order to estimate resonance effects, we use a model based on simple thermal ideas, plus resonance’s expected in a hot hadronic system. For thermal particle production we will use a simple factorized form for $`p_t`$ and $`y`$ of our particles and resonance’s.

$$\frac{d^2N}{dp_tdy}=Ap_te^{\frac{M_t}{Te}}e^{\frac{y^2}{2\sigma _y^2}}.$$
(1)
The inputs are the mass of the particle or resonance ($`M`$)(GeV), the calculated transverse mass($`M_t`$)(GeV),and the thermal temperature of the particle ($`Te`$)(GeV). A Gaussian width in rapidity $`\sigma _y`$ is also used. The resonance mass is smeared by a Breit-Wigner form. A typical resonance is the $`K^{}`$ which has a form:

$`W(M)`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Gamma }_{(q)}^2M_k^{}^2}{(M^2M_k^{}^2)^2+\mathrm{\Gamma }_{(q)}^2M_k^{}^2}}`$
$`\mathrm{\Gamma }(q)`$ $`=`$ $`{\displaystyle \frac{2\mathrm{\Gamma }_k^{}(q/q_0)^3}{1+(q/q_0)^2}}.`$ (2)
where $`W(M)`$ is the mass weighting of the $`K^{}`$, and $`\mathrm{\Gamma }(q)`$ is the total width (GeV). $`M_k^{}`$ is the mass (GeV) of the $`K^{}`$ and $`q`$(GeV/c) is the $`\pi ^{}K`$ C.M. momentum, with $`q0`$ (GeV/c) being the saturation momentum. The powers 3 and 2 are derived from $`2l+1`$ and $`2l`$, where $`l`$ is the angular momentum of the $`\pi ^{}K`$ system. We also need to add elliptic flow to our particle production. A very simple $`p_t`$ dependent $`V2`$ parameter is used for our resonance gas as we used for the beam jets in HIJING. We used the form

$$E\frac{d^3N}{dp^3}=\frac{1}{2\pi p_t}\frac{d^2N}{dp_tdy}[1+2V_2p_t\mathrm{cos}2(\mathrm{\Phi }\mathrm{\Phi }_R)]$$
(3)
| Table I | | | | |
| --- | --- | --- | --- | --- |
| particle | number | temp | width | $`V_2`$ |
| $`\pi ^+`$ | 115 | .250 | 3.00 | .12 |
| $`\pi ^{}`$ | 115 | .250 | 3.00 | .12 |
| $`\rho ^0`$ | 350 | .290 | 2.90 | .12 |
| $`\omega `$ | 252 | .290 | 2.90 | .12 |
| $`\eta `$ | 505 | .290 | 2.90 | .12 |
| $`k^+`$ | 133 | .260 | 2.20 | .09 |
| $`k^{}`$ | 116 | .260 | 2.20 | .09 |
| $`k_s`$ | 125 | .260 | 2.20 | .09 |
| $`k^+`$ | 19 | .260 | 2.20 | .06 |
| $`k^{}`$ | 17 | .260 | 2.20 | .06 |
| $`k^0`$ | 35 | .260 | 2.20 | .06 |
| $`p`$ | 45 | .320 | 1.80 | .03 |
| $`\overline{p}`$ | 33 | .320 | 1.80 | .03 |
| $`\mathrm{\Lambda }`$ | 71 | .360 | 1.80 | .03 |
| $`\overline{\mathrm{\Lambda }}`$ | 50 | .360 | 1.80 | .03 |
| $`\mathrm{\Sigma }^+`$ | 33 | .360 | 1.80 | .03 |
| $`\mathrm{\Sigma }^{}`$ | 33 | .360 | 1.80 | .03 |
| $`\overline{\mathrm{\Sigma }}^+`$ | 21 | .360 | 1.80 | .03 |
| $`\overline{\mathrm{\Sigma }}^{}`$ | 21 | .360 | 1.80 | .03 |
| $`\mathrm{\Xi }^0`$ | 17 | .330 | 1.80 | .03 |
| $`\mathrm{\Xi }^{}`$ | 17 | .330 | 1.80 | .03 |
| $`\overline{\mathrm{\Xi }}^0`$ | 13 | .330 | 1.80 | .03 |
| $`\overline{\mathrm{\Xi }}^+`$ | 13 | .330 | 1.80 | .03 |
| $`\mathrm{\Omega }^{}`$ | 4 | .300 | 1.80 | .03 |
| $`\mathrm{\Omega }^+`$ | 4 | .300 | 1.80 | .03 |

Here there is only one parameter which is $`V2`$. We can also add jets from HIJING as we did above. In Refs.  one expects a thermal resonance gas system should have a suppression of charge fluctuation of the order of 20% (if $`D`$ is $`3`$ due to resonance effects whereas it would be $`4`$ without resonance effects). Using the resonance to pion ratio of Ref.  we can generate central 65 GeV/n $`Au`$ on 65 Gev/n $`Au`$ events that are very close to the HIJING simulation, using a Boltzman temperature of 0.180 GeV, and a Gaussian width rapidity of 2 units. It is important to note that in order to make comparisons between the different models, we need to reproduce the single particle ($`p_t`$) and psuedo-rapidity distributions.

Then for a pseudo-rapidity of $`0.75<\eta <0.75`$, we can look at the distribution of net charge. From this net charge distribution we can determine the charge suppression. Keeping the final yield of particles the same, we increase the resonance’s until we achieve a 20% reduction in charge fluctuation which represents the resonance effect used in \[9-11\]. This is roughly consistent with observations. The net charge mean is 9.5 with a width of 20. For a random charge fluctuation system the width should be 25. Adding jets from HIJING to our resonance gas does not appreciably change the above results. In table 1, we give the number of particles and resonance’s used and the temperatures and rapidity widths for each plus the $`V2`$ parameter. Figure 11 shows the $`\mathrm{\Delta }\mathrm{\Phi }`$ comparison with resonance plus jets to HIJING plus bubbles. We have a $`\chi ^2`$ of 22 for 8 $`DF`$ which is a $`2.5\sigma `$ difference, and thus the two are statistically equivalent to being the same, since we consider at least $`3\sigma `$ required for minimal statistical significance. In Fig. 12 we make a $`p_t`$ cut ($`1.2<p_t<2.0`$GeV/c) and see that HIJING plus bubbles has a considerably increased correlation as a function of $`\mathrm{\Delta }\mathrm{\Phi }`$ in our usual 0-60 cut region. There is a $`\chi ^2`$ of 206 which corresponds to over 20$`\sigma `$. We see that resonance’s do not decay most of the time into particles with enough $`p_t`$ values to satisfy this cut. Thus by making a tighter $`p_t`$ cut the correlation in $`\mathrm{\Delta }\mathrm{\Phi }`$ due to resonance’s is greatly reduced, and we can use this tightening of cuts to reveal the correlation in $`\mathrm{\Delta }\mathrm{\Phi }`$ without significant contributions from resonance effects. We now can make a charge fluctuation analysis of the two models inside our standard cuts. The results for resonance plus jets is a 1% reduction in the width of the net charge (10.05 compared to 10.15), whereas for HIJING plus bubbles we see a reduction of 5% in the width (9.73 compared to 10.29).

## 7 Simulation of Bubbles with 3 Times the Energy per Bubble

Up to this point in our model building, we have calculated the energy in a bubble using the work in Ref. , and replacing one of their domains with a bubble of about 2 fm radius, which is one of their choices they mention. Our reason for selecting a 2 fm bubble is based on using the previously mentioned HBT work on source size, and adapting it to our bubble scenario. That seems somewhat reasonable as a starting point, in order to avoid an arbitrary selection by us. However, one has to admit that the energy per bubble is really a parameter and that their value could be off by a considerable factor. If our use of calculated energy per bubble based on Ref.  estimates is too large, we feel it would become increasingly difficult as the energy per bubble decreased to observe multiple bubbles effects using the methods explored above. There are many unknown uncertainties in the energy per bubble used from Ref. , including the turbulent fluctuation phenomena which exist. Therefore it is reasonable to speculate what would happen if the bubble energies were bigger. Therefore we have rerun our bubble simulation, using 3 times the energy per bubble, that was used originally. Figure 13 shows the $`\mathrm{\Delta }\mathrm{\Phi }`$ correlation with the bigger energy bubbles compared to the regular modified HIJING and jets of Fig. 2. This is a very large effect with a $`\chi ^2`$ of 1840, which is well over 20 $`\sigma `$. Next we show the unlike and like charge sign comparison in Fig. 14. Here the $`\chi ^2`$ is 4399, thus we see again a very large difference well over 20 $`\sigma `$ between like and unlike charge pairs. If we look at the charge difference distribution of the 3 times more energy per bubble simulation, we find a mean of 4 with a width of 7.6. A pure random charge fluctuation case has a width of 10.4 thus leading to about a 27% reduction of charge fluctuation, with 3/4 of the particles in our cuts coming from the plasma bubbles.

## 8 Most Extreme Pure Neutral Resonance Case

Finally, let us explore the most extreme totally unjustified case, which should give a maximum charge suppression by generating a pure neutral resonance system. This possibility has never been seen in heavy ions collisions. It violates isospin symmetry along with other reasonable considerations. We generated this case merely to demonstrate that our bubbles focus much more of their decay particles into the cuts than a resonance case can. We will consider charged particles generated by only resonance’s decaying into the region of our standard cuts ($`0.8<p_t<2.00.75<\eta <0.75`$). For our neutral resonance’s we will use the $`\rho `$ and the $`\sigma `$ mesons. The $`\rho `$ is given by the $`I=1`$, P-wave $`\pi \pi `$ phase shifts, while the $`\sigma `$ is given by the $`I=0`$, S-wave $`\pi \pi `$ phase shift(see Ref. ).

In the heavy ion final state, we expect many resonance’s will be formed by pion re-scattering in the final state. Therefore the $`\sigma `$ meson will have a mass shift due to re-scattering like the $`a1`$ meson does in diffractive production. This re-scattering will create a threshold peak in the di-pion effective mass spectrum. We can show the unlike charge di-pion correlation as a function of the effective mass if we plot the ratio of the unlike charge pions over the like charge pions. In Fig. 15 we show this ratio versus the di-pion effective mass. There is a threshold bump given by the $`\sigma `$ and a second bump given by the $`\rho `$. In Fig. 16 we show the $`\mathrm{\Delta }\mathrm{\Phi }`$ correlation that we have generated using our two resonance’s. Again we must make sure that both the single particle $`p_t`$ and pseudo-rapidity spectra are the same as in our other simulations. The values used are shown in Table 2. In the table we see that we need two different temperature sources of neutral resonance’s in order to obtain the $`p_t`$ spectrum.

This time we compare our neutral resonance system with the resonance gas system plus bubbles. The $`\chi ^2`$ for these two in Fig. 16 is 11 which is within $`0.5\sigma `$ the same. We have compared this extreme pure neutral resonance case to our standard bubbles case where $`\frac{1}{4}`$ of the particles come from the bubbles, while 100% of the particles come from the neutral resonance’s. However again if we make a narrower $`p_t`$ cut ($`1.2<p_t<2.0`$) the bubble correlation is larger by a $`\chi ^2`$ of 36 which is a $`5\sigma `$ effect (Fig. 17). Also the width in $`\eta `$ is different than that of the bubbles. This can be seen by looking at the opening angle correlation (Fig. 18) which should be compared to Fig. 9. The $`\chi ^2`$ in Fig. 18 is 50 which is a $`7\sigma `$ effect. The reduction of charge fluctuation for the pure neutral resonance system is an interesting result. This system has a net charge of zero, and would have no charge fluctuation if there were no decays. It has only a 6% reduced fluctuation in our case. The cause of this effect is decay particles leaking out of the cuts. This 6% is the total amount possible from neutral resonance’s decaying into our cuts even in this extreme unrealistic case. If we increase our bubble energy by 3 times we achieve a 27% reduction in charge fluctuations. For the resonance gas plus bubbles case which we use in Figs. 16-18, we find a 7% reduction of charge fluctuations. This is a greater reduction than the pure neutral resonance case. In our bubble model the decay particles are well contained and focused into our cut region.

Thus considering the foregoing analyses, we conclude that background resonance effects will not to any significant degree affect our conclusions based on our simulations.

| Table II | | | | |
| --- | --- | --- | --- | --- |
| particle | number | temp | width | $`V_2`$ |
| $`\sigma `$ | 69 | .700 | 1.80 | .04 |
| $`\rho ^0`$ | 20 | .700 | 1.80 | .04 |
| $`\sigma `$ | 1246 | .335 | 1.80 | .07 |
| $`\rho ^0`$ | 380 | .335 | 1.80 | .08 |

## 9 Summary and Conclusions

From the early days of QGP theories we have shared the viewpoint of Leon Van Hove, that localized in phase space bubbles of QGP are more likely to be the origin of the hadronization products which originate from a QGP. This is especially the case for a QGP produced at high energy heavy ion colliders such as RHIC and LHC.

In a previous publication we have considered the case of one to at most a few separated QGP bubbles being produced. This could conceivably occur, and be observed with sufficient statistics. We referred to multiple bubble formation in Ref. , which is the most likely high cross-section case. However we did not calculate that case in Ref. . In this paper we conservatively calculated the multiple bubble case. The bubbles were embedded in HIJING, suitably modified for relevant elliptic flow effects.

We demonstrated that resonance effects will not significantly affect our conclusions. The first test of our ideas will come from forthcoming STAR data at RHIC. Therefore we made every reasonable effort to cast our predictions in a form which could be directly and quantitatively compared with the STAR data for 65GeV/n $`Au`$ colliding with 65 GeV/n $`Au`$. Obviously the methods and models we used could also be used for higher energy RHIC data, and eventually for LHC data. It is expected that high statistics data on 100GeV/n $`Au`$ colliding with 100 GeV/n $`Au`$ will subsequently become available in the near future and allow a more critical comparison with these ideas. We used a 2 fm radius bubble size based on HBT work. This work shows that above $`p_t`$ of 0.8 GeV/c the source size has a radius of about 2 fm. This is consistent with a 2 fm radius bubble size in Ref.  (one of their choices). We showed that a reasonable estimate of the outer shell of a ring of adjoining bubbles around the blast region, at central rapidity, when hadronization takes place, is located at a radius of about 8 fm. The width of this outer shell is roughly $`\pm `$ 2 fm about the 8 fm radius. We added a reasonable estimate of longitudinal expansion to the bubbles in order to take account of the expected Landau effect, but comparisons with the data can determine this value.

The lower bound in $`p_t`$ is chosen to maintain the space momentum correlation corresponding to the 2 fm source size, which is consistent with the HBT work and enhances our signal. The upper bound in $`p_t`$ is chosen to avoid the contamination from the hard physics region above 2 GeV/c. However at times we also raised the lower bound of the $`p_t`$ cut to 1.2 GeV/c, for various reasons as explained in the text. Since we estimated the $`\mathrm{\Phi }`$ angular range to which each bubble contributes is about 30, we used an average of 13 adjoining bubbles around a ring in the central rapidity region. This fills the $`\mathrm{\Phi }`$ coverage. Our standard cut in $`\eta `$ was $`0.75`$ to $`+0.75`$ to cover the central region dominated by our bubble ring hadronization products. Each bubble is considered to be a localized gluon dominated, uncharged color singlet system. Thus when the bubble hadronizes the total charge coming from the bubble was set equal to zero.

We have demonstrated by simulations signals for the QGP bubbles. Below we list signals which had a statistical significance of $`5\sigma `$ or more:

* The angular correlation of charged particles as a function of $`\mathrm{\Delta }\mathrm{\Phi }`$, and also opening angles in the 0-60 $`\mathrm{\Delta }\mathrm{\Phi }`$ range, where our bubbles contribute.
* Predicted suppression of electric-charge fluctuations.
* For correlation vs $`\mathrm{\Delta }\mathrm{\Phi }`$: Fig. 2 compares modified HIJING and modified HIJING with the addition of bubbles, which results in an $`8\sigma `$ difference effect (see text) for our 0-60 cut where the effect of bubbles is expected to occur. This cut is used on all angular distributions. To remove this difference requires the very conservative extreme measure of doubling the number of jets (as shown in Fig. 3 and text), when there are arguments that jet quenching occurs rather than enhancement.
* In Fig. 10 we see that when we compare the $`\mathrm{\Delta }\mathrm{\Phi }`$ correlations for like and unlike charge particles there is a $`5\sigma `$ difference in the case of the bubbles plus modified HIJING, and modified HIJING plus double jets. This is caused by charge suppression due to bubbles.
* We pointed out that earlier work on charge fluctuation suppression overlooked color charge effects, and kinematic mixing effects, which could drastically reduce, or even more or less eliminate the calculated suppression. In our localized bubble scenario these effects are negligible. In our analysis the predicted observed charge suppression in the 65 GeV/n $`Au`$ on 65 GeV/n $`Au`$ simulations will range from about 5% for the Ref.  based value of about 3.25 GeV per bubble going into charged particles, to about 27% for about 10 GeV per bubble going into charged particles. In the case of 10 GeV per bubble going into charged particles, we find the difference in the $`\mathrm{\Delta }\mathrm{\Phi }`$ correlation between our modified HIJING in Fig. 2 and modified HIJING plus bubbles is huge (see Fig. 13 and text), and well over $`20\sigma `$. Also in Fig. 14 using the same modified HIJING plus bubbles model the difference in the $`\mathrm{\Delta }\mathrm{\Phi }`$ correlation between like and unlike charge correlation’s is huge, and well over $`20\sigma `$.