format: text

text:
# Untitled Document

COHERENT SCALAR FIELD OSCILLATIONS AND

THE EPOCH OF DECELERATION

by

J.C. Jackson**e-mail: john.jackson@unn.ac.uk

School of Computing and Mathematics

University of Northumbria at Newcastle

Ellison Building

Newcastle upon Tyne NE1 8ST, UK

ABSTRACT

The discovery of supernova SN 1997ff at $`z1.7`$ has confirmed the expected switch from cosmological acceleration to deceleration, as predicted by the concordance $`\mathrm{\Lambda }`$/CDM model. However its position in the SN Ia Hubble diagram suggests that the switch is too pronounced, which here is taken to mean that a cosmological constant is not an adequate description of the state of the vacuum. An ‘oscillessence’ model is invoked, with a scalar field $`\varphi `$ governed by a simple quadratic potential, which gives a better fit to the new data point. The field is undergoing coherent oscillations, and a key feature of the proposal is that we are towards the end of the second period of acceleration; a $`\mathrm{\Lambda }/\varphi `$ mix replaces $`\mathrm{\Lambda }`$/CDM, with $`\mathrm{\Omega }_\mathrm{\Lambda }0.4`$ and $`\mathrm{\Omega }_\varphi 0.6`$.

Key words: cosmology – observations – theory – dark matter.

1 INTRODUCTION

There is now a strong concensus that the basic cosmological parameters are known, and that we are living in a spatially flat accelerating Universe, with $`\mathrm{\Omega }_{\text{CM}}0.3`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }1\mathrm{\Omega }_{\text{CM}}0.7`$ (CM $``$ Cold Matter, that is Cold Dark Matter plus the baryonic component). This concencus is based primarily upon observations of Type Ia supernovae (SNe Ia) (Schmidt et al. 1998; Riess et al. 1998; Perlmutter et al. 1999), coupled with observations of Cosmic Microwave Background (CMB) anisotropies (Efstathiou et al. 1999) and Large Scale Structure (LSS) information (Bridle et al. 1999; Lasenby, Bridle & Hobson 2000; Efstathiou et al. 2001). Separately these observations place different constraints on the parameter values, but in combination degeneracies are removed. However, the CMB measurements have been refined through a number of ground-based and balloon-borne experiments (Halverson et al. 2001; Lee et al. 2001; Netterfield et al. 2001), particularly significant being the detection of multiple peaks in the CMB angular power spectrum, and the data sets are now beginning to tell the same story when considered separately (Balbi et al. 2000; de Bernardis et al. 2001; Jaffe et al. 2001; Pryke et al. 2001; Stompor et al. 2001). Support for flatness comes from measures of the angular size of the acoustic horizon at decoupling, via detection of the first Doppler peak in the CMB angular spectrum (de Bernardis et al. 2000; Hanany et al. 2000); for this reason I shall restrict my attention to spatially flat models.

This impressive convergence of independent measures of the various cosmological parameters is coming to be known as ‘concordance’, and the corresponding cosmological model as the concordance model. The concordance model predicts a transition from acceleration to deceleration at a redshift $`z=(2\mathrm{\Omega }_\mathrm{\Lambda }/\mathrm{\Omega }_{\text{CM}})^{1/3}10.67`$, beyond which the matter density is dominant; the discovery of SN 1997ff (Gilliand & Phillips 1998), together with fortuitous retrospective multicolour photometric observations which have established its type, a redshift $`z1.7`$ and a distance modulus, provide strong evidence that the epoch of deceleration has been observed (Riess et al. 2001). I shall address here what appears to be a discordant note arising from these observations, which is a strong suggestion that the switch to deceleration is too violent. Figure 1 is a Hubble diagram for SNe Ia, reproducing the results presented in Figure 11 of Riess et al. (2001), in which the points with error bars are redshift-binned data from Riess et al. (1998) and Perlmutter et al. (1999); here the magnitude $`\mathrm{\Delta }m`$ is relative to an empty non-accelerating universe. The cross shows SN 1997ff at $`z=1.7`$, $`\mathrm{\Delta }m=0.5`$, and the corresponding 68% and 95% confidence regions are indicated. The continuous curve is an unweighted least-squares fit to the $`z<1`$ points, ($`\mathrm{\Omega }_{\text{CM}}=0.28`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.72`$), which clearly illustrates the point I am making here. The dashed curve shows an unweighted least-squares fit to all the points ($`\mathrm{\Omega }_{\text{CM}}=0.44`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.56`$); the extra matter produces more deceleration, but destroys the good fit to the $`z<1`$ points. The possibility discussed here is that this behaviour is evidence that a simple cosmological constant is not an adequate description of the present state of the vacuum, and the purpose of this letter is to present a simple alternative which is compatible with the SN Ia Hubble data, and might at the same time do least damage to the concordance picture.

2 $`\mathrm{\Lambda }/\varphi `$ MODELS

Some time ago (before the case for acceleration was as strong as it now appears to be) (Jackson 1998a,b; Jackson & Dodgson 1998; see also Frieman et al. 1995) I considered models dominated by a homogeneous dynamical scalar field $`\varphi `$ governed by an inflationary potential corresponding to an ultra-light inflaton:

$$V(\varphi )=\frac{\mathrm{\Lambda }}{8\pi G}+\frac{1}{2}\omega _c^2\varphi ^2$$
$`(1)`$
where $`\omega _c`$ is the associated Compton frequency, here taken to be somewhat larger than the inverse of the present Hubble time. The scalar field has effective density $`\rho _\varphi =(\dot{\varphi }^2+\omega _c^2\varphi ^2)/2`$ and pressure $`p_\varphi =(\dot{\varphi }^2\omega _c^2\varphi ^2)/2`$. The most convenient formulation of the Friedmann equations governing the scale factor $`R(t)`$ is in terms of the instantaneous values of the deceleration parameter $`q(t)`$, and the various density parameters $`\mathrm{\Omega }=8\pi G\rho /3H^2`$; in the spatially flat case these are

$$q=\mathrm{\Omega }_{\text{CM}}/2+(1+3w)\mathrm{\Omega }_\varphi /2\mathrm{\Omega }_\mathrm{\Lambda }$$
$`(2)`$
and

$$\mathrm{\Omega }_{\text{CM}}+\mathrm{\Omega }_\mathrm{\Lambda }+\mathrm{\Omega }_\varphi =1$$
$`(3)`$
where $`w`$ is the ratio $`p_\varphi /\rho _\varphi `$.

The scalar field is governed by the equation

$$\ddot{\varphi }+3H\dot{\varphi }+\omega _c^2\varphi =0,$$
$`(4)`$
where $`H=\dot{R}/R`$ is Hubble’s ‘constant’. Initially (when $`H>>\omega _c^1`$) the field $`\varphi `$ undergoes a slow roll towards the minimum in $`V(\varphi )`$ at $`\varphi =0`$, during which phase the effective cosmological constant is $`\mathrm{\Lambda }+4\pi G\omega _c^2\varphi ^2`$. This is the standard inflationary picture, but in this context the pressure-free matter is assumed to be dynamically dominant during the slow-roll phase. Thereafter, when $`H<\omega _c^1`$, $`\varphi `$ undergoes coherent oscillations, and the speculation here (and in Jackson 1998a,b; Jackson & Dogson 1998) is that the Universe has entered a phase in which the scalar field is dynamically dominant and is executing such oscillations. During this phase periods of violent deceleration (when $`\rho _\varphi +3p_\varphi >0`$) alternate with periods of not-so-violent acceleration (when $`\rho _\varphi +3p_\varphi <0`$), which behaviour can account for the SN 1997ff observations. Figure 2 is a typical example, showing the evolution of the deceleration parameter $`q(t)`$ and the age/Hubble time ratio $`t/t_H`$. In what follows I shall assume generically that we are currently at point P, at the end of the second period of acceleration, to maximise the effect of the latter. In the example shown the parameter values at P are $`\mathrm{\Omega }_{\text{CM}}=0.012`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.433`$, $`\mathrm{\Omega }_\varphi =0.555`$; their sum is 1 as dictated by equation (3), and $`q=0`$ at P fixes $`w`$ (or equivalently $`\varphi `$ and $`\dot{\varphi }`$) and hence the model according to equation (2). Locally these values give the best fit to the SN Ia data, including SN 1997ff; the corresponding curve is shown in Figure 3. Observation points located during the first period of acceleration do no better in this respect than straightforward $`\mathrm{\Lambda }/\text{CM}`$ models, as in Figure 1. Even later periods of acceleration might be considered, but the corresponding Hubble curves generally have features for which there is no evidence.

A interesting possibility is that just the quadratic term in the potential (1) might suffice, with no true cosmological constant, but this is discounted by the dashed curve in Figure 3, with $`\mathrm{\Omega }_{\text{CM}}=0`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$, $`\mathrm{\Omega }_\varphi =1`$; deceleration sets in too soon in this case, and addition of some Cold Matter would clearly make the situation worse. The picture which emerges here is that of a mixed $`\mathrm{\Lambda }/\varphi `$ model, with a low Cold Matter content, of the order of the baryon content of the Universe (see Coc et al. (2001) for a recent review). If SN Ia data were the only consideration then this model would be perfectly viable, but almost certainly not compatible with CMB and LSS. However, there are acceptable models (i.e. when SN 1977ff is included) with more Cold Matter (but not more than $`\mathrm{\Omega }_{\text{CM}}0.1`$), which nevertheless probably are compatible. This universe is heading for heat death, with a final oscillatory flourish before oblivion; this is illustrated in Figure 4.

I have allowed blatant fine-tuning here, with no attempt to account for the balance between the two components of dark energy. Following the seminal work of Ratra & Peebles (1988) and Wetterich (1988), who considered decaying potentials of the form $`V(\varphi )\varphi ^\alpha `$ and $`V(\varphi )\mathrm{exp}(\alpha \varphi )`$, the manufacture of scalar-field models has evolved into a large industrial concern, particularly post 1998 and particularly with view to accounting for the balance between dark energy and dark matter in a natural way (see Bean and Melchiorri (2001) for a recent review). Some of the corresponding potentials also engender oscillatory behaviour (for example Frieman et al. 1995; Skordis & Albrecht 2000), and might serve as a alternatives to equation (1) in the present context. Alternatively, if the simplicity of potential (1) is to be preferred, then an anthropic argument might be found to account for the putative balance. The generic term ‘oscillessence’ might be reserved for dark energy which allows multiple periods of acceleration/deceleration (cf. Caldwell, Dave & Steinhardt 1998).

FIGURE CAPTIONS

Figure 1. Type 1a supernovae and $`\mathrm{\Lambda }`$/CM models. Hubble diagram relative to to an empty non-accelerating universe. The cross indicates SN 1997ff, with 68% and 95% confidence regions. The continuous curve shows a concordance model ($`\mathrm{\Omega }_{\text{CM}}=0.28`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.72`$); the dashed curve is the best fit to all the data points ($`\mathrm{\Omega }_{\text{CM}}=0.44`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.56`$).

Figure 2. An oscillating model, showing the evolution of the deceleration parameter $`q(t)`$ (continuous) and the ratio $`t/t_H`$ (dashed). The current observation point is located at P, where $`\mathrm{\Omega }_{\text{CM}}=0.012`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.433`$, $`\mathrm{\Omega }_\varphi =0.555`$.

Figure 3. Type 1a supernovae and $`\mathrm{\Lambda }/\varphi `$ models. Hubble diagram relative to to an empty non-accelerating universe. The cross indicates SN 1997ff, with 68% and 95% confidence regions. The continuous curve shows the best fitting $`\mathrm{\Lambda }/\varphi `$ model, with parameters as in Figure 2. The dashed curve has parameters $`\mathrm{\Omega }_{\text{CM}}=0`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$, $`\mathrm{\Omega }_\varphi =1`$. The dash-dot curve is the concordance model as in Figure 1.

Figure 4. Oscillating model showing the evolution of various density parameters $`\mathrm{\Omega }=8\pi G\rho /3H^2`$; $`\mathrm{\Omega }_{\text{CM}}`$ (dotted), $`\mathrm{\Omega }_\mathrm{\Lambda }`$ (dashed), $`\mathrm{\Omega }_\varphi `$ (continuous). The dash-dot curve is a composite parameter based upon the active gravitational density $`\rho _\varphi +3p_\varphi `$.

REFERENCES

Skordis C., Albrecht A., 2000, astro-ph astro-ph/0012195

Balbi et al., 2001, ApJ, 2000, 545, L1; astro-ph/0005124

Bean R., Melchiorri A., 2001, astro-ph/0110472

Bridle S.L., Eke, V.R., Lahav O., Lasenby A.N., Hobson M.P., Cole S., Frenk C.S., Henry J.P., 1999, MNRAS, 310, 565

Caldwell R.R., Dave R., Steinhardt P.J., 1998, Phys. Rev. Lett., 80, 1582

Coc A., Vangioni-Flam E., Cassé M., Rabiet M., 2001, astro-ph/0111077

de Bernardis P. et al., 2000, Nature, 404, 955

de Bernardis P. et al., 2001, astro-ph/0105296

Efstathiou G., Bridle S.L., Lasenby A.N., Hobson M.P., Ellis R.S., 1999, MNRAS, 303, 47

Efstathiou et al., 2001, astro-ph/0109152

Frieman J.A., Hill C.T., Stebbins A., Waga I., 1995, Phys. Rev. Lett., 75, 2077

Gilliland R.L., Phillips M.M., 1998, IAU Circ. 6810

Hanany S. et al., 2000, ApJ Letters, 545, 5

Halverson N.W. et al., 2001, astro-ph/0104489

Jackson J.C., 1998a, MNRAS, 296, 619

Jackson J.C., 1998b, Mod. Phys. Lett. A, 13, 1737

Jackson J.C., Dodgson M., 1998, MNRAS, 297, 923

Jaffe et al., 2001, Phys. Rev. Lett., 86, 3475; astro-ph/0007333

Lee A. T. et al., 2001, astro-ph/0104459

Lasenby A.N., Bridle S.L., Hobson, M.P., 2000, Astrophys. Lett. & Communications, 37, 327

Netterfield C.B. et al., 2001, astro-ph/0104460

Ratra B., Peebles P.J.E., 1988, Phys. Rev., D37, 3406

Perlmutter S. et al., 1999, ApJ, 517, 565

Pryke, C., 2001, astro-ph/0104490

Riess A.G. et al., 1998, AJ, 116, 1009

Riess A.G. et al., 2001, astro-ph/0104455

Schmidt B.P. et al., 1998, ApJ, 507, 46

Stompor et al., 2001, astro-ph/0106451

Wetterich, C., 1988, Nucl. Phys. B, 302, 668


text:
# Tight-binding study of high-pressure phase transitions in titanium: alpha to omega and beyond

\[

## Abstract

We use a tight-binding total energy method, with parameters determined from a fit to first-principles calculations, to examine the newly discovered $`\gamma `$ phase of titanium. Our parameters were adjusted to accurately describe the $`\alpha `$Ti-$`\omega `$Ti phase transition, which is misplaced by density functional calculations. We find a transition from $`\omega `$Ti to $`\gamma `$Ti at 102 GPa, in good agreement with the experimental value of 116 GPa. Our results suggest that current density functional calculations will not reproduce the $`\omega `$Ti-$`\gamma `$Ti phase transition, but will instead predict a transition from $`\omega `$Ti to the bcc $`\beta `$Ti phase.

\]

Structural transformations in titanium have received a great deal of experimental  and theoretical attention. This Letter is motivated by a recent experimental study which revealed a previously unsuspected phase transition in titanium at 116 GPa from the $`\omega `$Ti phase to a new $`\gamma `$Ti phase. We have been able to confirm these experiments by performing highly accurate tight-binding calculations of the phase diagram of Ti.

At room temperature the group-IV metals zirconium and hafnium transform under pressure from the ground state hexagonal close packed (hcp) phase to the intermediate pressure $`\omega `$ phase (space group $`P6/mmmD_{6h}^1`$, Pearson Symbol hP3, Strukturbericht Designation: C32) at 2.2 GPa and 38 GPa, respectively. At 35 GPa and 71 GPa, respectively, the metals transform from the $`\omega `$ phase to a body centered cubic (bcc) structure. One would logically assume that titanium also follows this transformation sequence, and indeed the transition from hcp $`\alpha `$Ti to $`\omega `$Ti takes place at a pressure variously given as 9 GPa or 20-90 GPa. However, no room temperature pressure driven transition from $`\omega `$Ti to $`\beta `$Ti has been observed, although first-principles calculations predict a transition at 98 GPa, and $`\beta `$Ti exists at room pressure and temperatures above 1155 K. Recently, however, Vohra and Spencer found that at 116 GPa titanium transforms from the $`\omega `$ phase to a previously unsuspected $`\gamma `$Ti phase. The new phase has a two-atom body-centered orthorhombic unit cell, space group $`Cmcm`$$`D_{2h}^{17}`$, Pearson symbol oC4, with the atoms at the points $`(0,\pm yb,\pm c/4)`$, where $`a`$, $`b`$, and $`c`$ are the lengths of the primitive vectors in the full orthorhombic unit cell, and $`y`$ is an internal parameter.

This crystal structure is observed in various materials, including $`\alpha `$U, the random alloy AgCd, and a metastable form of gallium. (The $`\alpha `$U phase has the Strukturbericht designation A20.) With appropriate choices of the parameters this structure can reproduce several higher symmetry phases. In particular, when $`b/a=\sqrt{3}`$ and $`y=1/6`$, it becomes the hcp structure, while when $`b/a=c/a=\sqrt{2}`$ and $`y=1/4`$ the atoms are on the sites of a bcc cell. Wentzcovitch and Cohen used this pathway to describe a possible theoretical model for the hcp $``$ bcc transition in magnesium.

Examination of the $`\gamma `$Ti structure by first-principles techniques requires a minimization of the total energy with respect to three parameters (e.g., $`b/a`$, $`c/a`$, and $`y`$) at several volumes. This is impractical because of the high computational demand of first principles methods. We have instead chosen to study the $`\alpha `$-$`\omega `$-$`\gamma `$ transformation sequence using the much faster NRL tight-binding method. This method has been shown to reproduce the ground-state phase, elastic constants, surface energies, and vacancy formation energies of the transition metals. The tight-binding parameters in Ref.  were found by fitting to a Local Density Approximation (LDA) database of total energies and eigenvalues for the fcc and bcc structures. The parameters correctly predicted the ground state structures of all of the transition and noble metals, including the hcp metals and manganese. However, upon examination, we found that the titanium parameters from Ref.  do not predict the correct position for the $`\omega `$Ti phase. In fact, no $`\alpha `$Ti-$`\omega `$Ti phase transition is seen.

We therefore developed a new set of tight-binding parameters according to the procedures of Ref. , fit to an expanded database of first-principles calculations. In particular, our database includes the fcc, bcc, simple cubic, hcp, and $`\omega `$ structures. The eigenvalues and total energies were generated using the general-potential Linearized Augmented Plane Wave (LAPW) method, using the Perdew-Wang 1991 Generalized Gradient Approximation (GGA) to density functional theory. We fit our tight-binding parameters to both total energies and band structures, using the parametrization described by equations (7), (8), (9), and (11) of Ref. . The RMS error in fitting the energies for the lowest energy phases (hcp, $`\omega `$, fcc, and bcc) was less than 1 mRy/atom. The band structure RMS error is about 10 mRy for the occupied bands of the hcp and $`\omega `$ structures.

In agreement with previous calculations, our first-principles results show that the $`\omega `$Ti phase is slightly lower in energy (about 0.5 mRy/atom) than $`\alpha `$Ti. This implies a $`5`$ GPa $`\omega `$Ti-$`\alpha `$Ti phase transition. We have adjusted our tight binding parametrization to shift the $`\omega `$Ti phase equilibrium upwards by 0.8 mRy/atom. This produces an $`\alpha `$Ti-$`\omega `$Ti phase transition at 6 GPa, in good agreement with the 9 GPa transition found in Ref.  As we shall see, this has important consequences for the $`\omega `$Ti-$`\gamma `$Ti phase transition.

We have tested the fit in a variety of ways. In Table I we show the equilibrium lattice constants for $`\alpha `$-, $`\beta `$\- and $`\omega `$Ti, as determined by our tight-binding parameters, our first-principles LAPW calculations, and experiment. The TB agreement with experiment is comparable to that achieved by the first-principles calculations.

We also examined the behavior of a variety of crystal structures using our tight-binding parameters. Fig. 1 shows the energy/volume behavior of a number of low energy structures. As expected, the observed phases ($`\alpha `$Ti and $`\omega `$Ti) are followed by close-packed stacking fault phases (9R, 4H, and fcc), and then other common metallic phases.

We further checked the behavior of our tight-binding Hamiltonian by determining the elastic constants and phonon frequencies in $`\alpha `$Ti, as shown in Tables II and III, respectively. Compared to experiment, we find an RMS error of 22 GPa for the elastic constants, and 32 cm-1 for the phonon frequencies. This is typical of the predictive capability of the tight-binding method for hcp metals.

We studied the $`\alpha `$-$`\omega `$-$`\gamma `$ transition path in titanium by fixing the volume of a given phase, and then minimizing the total energy as a function of the other parameters ($`c/a`$ for $`\alpha `$ (hcp) and $`\omega `$; $`b/a`$, $`c/a`$, and $`y`$ for $`\gamma `$). The pressure was calculated in one of two ways: by differentiating an extended Birch fit, and by calculating the pressure by numerical differentiation of the total energy with respect to volume. The enthalpy of each phase, $`H(P)=E+PV`$ is then calculated by both methods. In Fig. 2 we show the enthalpy of the $`\omega `$Ti, $`\gamma `$Ti, and bcc ($`\beta `$Ti) phases in the transition region. From the plot we see that the $`\omega `$Ti-$`\gamma `$Ti phase transition takes place at about 102 GPa, compared to the experimental result of 116 GPa. We also see a $`\gamma `$Ti-$`\beta `$Ti phase transition at about 115 GPa. This is not seen experimentally, but it suggests that we may expect a higher pressure $`\gamma `$Ti-$`\beta `$Ti phase transition, which would complete the $`\alpha `$-$`\omega `$-$`\beta `$ transition sequence seen in Zr and Hf, albeit with an interloping $`\gamma `$Ti phase between $`\omega `$Ti and $`\beta `$Ti. More details of the phase transitions predicted by our Hamiltonian are shown in Table IV.

In the absence of the $`\gamma `$Ti phase, Fig. 2 shows that there would be an $`\omega `$Ti-$`\beta `$Ti phase transition at 110 GPa. This is in good agreement with the prediction made from the LAPW/GGA calculations in our database, 105 GPa, and with the LMTO/GGA prediction of 98 GPa found in Ref. .

The behavior of titanium in the $`\alpha `$-, $`\beta `$-, and $`\gamma `$Ti phases is explored further in Fig. 3, which shows the lattice parameters $`b/a`$, $`c/a`$, and $`y`$ as a function of the volume. We see that at a volume of about 85 Bohr3/atom there is an abrupt change from hcp $`\alpha `$Ti into the lower symmetry $`\gamma `$Ti phase. From this point the structure merges more or less continuously into bcc $`\beta `$Ti at about 70 Bohr3/atom. Note, however, that none of these phases is observed in the volume range 74-108 Bohr3, as this is the region where the $`\omega `$Ti phase is stable.

Finally, we note that our LAPW calculations and other first-principles calculations place the $`\omega `$Ti phase slightly lower in energy than the $`\alpha `$Ti phase, leading to a direct transition from $`\omega `$Ti to $`\beta `$Ti at 105 GPa. Hence, the essential difference between the first-principles calculations and our TB model is the ordering of the $`\alpha `$Ti and $`\omega `$Ti phases.

In conclusion, our tight-binding Hamiltonian provides a good description of the low pressure behavior of titanium, and shows the correct $`\alpha `$-$`\omega `$-$`\gamma `$ transition sequence as reported in recent experiments. Our work suggests that current first-principles density functional calculations, which place the $`\omega `$Ti phase below the $`\alpha `$Ti phase, will also fail to predict the stability of the $`\gamma `$Ti phase under pressure.

We thank I. I. Mazin for reminding us of Ref. , and T. A. Adler for pointing us to the $`\alpha `$U structure. This work was supported by the U. S. Office of Naval Research. The development of the tight-binding codes was supported in part by the U. S. Department of Defense Common HPC Software Support Initiative (CHSSI).


text:
# Measurement of the Cosmic-Ray Antiproton to Proton Abundance Ratio between 4 and 50 GeV

thanks: Now at School of Physics and Astronomy, University of Minnesota, MN 55455thanks: Now at California Institute of Technology, Pasadena, CA 91125thanks: Now at Physics and Geology Department, Northern Kentucky University, Highland Heights, KY 41099thanks: corresponding author: schubnel@umich.edu

## Abstract

We present a new measurement of the antiproton to proton abundance ratio, $`\overline{p}/p`$, in the cosmic radiation. The HEAT-pbar instrument, a balloon borne magnet spectrometer with precise rigidity and multiple energy loss measurement capability, was flown successfully in Spring 2000, at an average atmospheric depth of 7.2 g/cm2. A total of 71 antiprotons were identified above the vertical geomagnetic cut-off rigidity of 4.2 GV. The highest measured proton energy was 81 GeV. We find that the $`\overline{p}/p`$ abundance ratio agrees with that expected from a purely secondary origin of antiprotons produced by primary protons with a standard soft energy spectrum.

preprint: HEP/123-qed

Antiprotons constitute a rare but interesting component of the cosmic radiation. They are secondary cosmic-ray particles, generated in nuclear interactions of high-energy cosmic rays with the interstellar medium (ISM). It remains an open question whether there are significant additional contributions that have a different and perhaps more exotic origin. The kinematic threshold energy for $`\overline{p}`$ production in p-p collisions of primary cosmic-ray protons causes a $`\overline{p}`$ energy spectrum, and a $`\overline{p}/p`$ intensity ratio, that decline rapidly from a few GeV towards lower energies. Solar modulation inside the heliosphere softens this ‘cutoff’ and leads to uncertainties in the predicted flux at low energy. Additional sources of antiprotons might be evaporating primordial black holes (PBH), or annihilating supersymmetric particles. The PBH contribution would be expected to be significant at energies well below 1 GeV, while supersymmetric particle annihilations, for instance neutralinos, could also affect the antiproton intensity at higher energy, above several GeV. At those higher energies, the $`\overline{p}`$ energy spectrum will be essentially unaffected by uncertainties due to solar modulation.

The energy spectrum of antiprotons measured near Earth carries the imprint of losses during propagation through the Galaxy and thus, is a sensitive probe of the confinement environment of protons. For instance, if the propagation path length $`\lambda `$ and the diffusion coefficient for protons depended on energy E in the same way as has been observed for the heavy cosmic ray nuclei (i.e. $`\lambda E^{0.6}`$), one would expect the antiproton fraction $`\overline{p}/p`$ to gradually decrease with increasing energy above a few GeV. This behavior would be equivalent to that of the intensities of secondary spallation nuclei such as Li, Be, and B, relative to those of their heavier primary parents C and O. In order to obtain a self-consistent model of the propagation of protons in the Galaxy the observed $`\overline{p}`$ spectrum must also be compared with measurements of positrons and gamma rays, which also result from nuclear interactions in the ISM (mostly via $`\pi ^+`$ and $`\pi ^0`$ decay).

It would be difficult to explain if the antiproton fraction $`\overline{p}/p`$ were found to be constant or if it increased at higher energy. Models of extragalactic origin are unlikely because of constraints in intergalactic transport Adams et al. (1997). Closed Galaxy models have been suggested Peters and Wetergaard (1977) which would boost the $`\overline{p}/p`$ ratio at high energy, but these would also predict a higher abundance of He3 than has been observed Beatty et al. (1993); Reimer et al. (1998). It has been proposed that the primary proton spectrum in distant regions of the Galaxy is harder than locally observed Moskalenko et al. (1998). This would also lead to enhanced $`\overline{p}`$ production in these regions, although not necessarily to an enhanced $`\overline{p}`$ fraction near Earth. Thus it appears that a constant or rising $`\overline{p}/p`$ fraction might indeed require the presence of a primary and possibly exotic source of antiprotons.

Observationally, the situation has been unclear. Measurements at low energies, in particular the series of observations with the BESS instrument Orito et al. (2000); Maeno et al. (2000), have provided the $`\overline{p}`$ energy spectrum with good statistical accuracy from $``$ 0.2 to 3 GeV. These results are in good agreement with interstellar secondary production models. PBH contributions, if they exist, are hidden by the uncertainties due to solar modulations. Above 5 GeV, the results reported in the three previous measurements Mitchell et al. (1996); Hof et al. (1996); Boezio et al. (1994, 2001) are statistics limited, and no solid conclusion about the shape of the energy dependence of the antiproton fraction can be drawn. The current HEAT-pbar experiment has been developed to clarify this situation with a series of balloon flights.

A schematic view of the HEAT-pbar instrument is shown in Fig. 1. It consists of a superconducting magnet with a drift tube hodoscope at its center, combined with stacks of multi-wire proportional chambers above and below the hodoscope. Two layers of scintillators, one at the very top and one at the very bottom of the instrument, provide time-of-flight (ToF) measurements and, together with a scintillator just above the hodoscope, form the event trigger.

The spectrometer has been described in detail Barwick et al. (1997). It consists of 479 drift tubes in twenty four layers. It is mounted in the room-temperature bore of a superconducting magnet with a 10 kG central field, and measures the particle trajectory within the magnetic field, providing both the particle rigidity (momentum/charge) and charge sign. It has a single-point tracking accuracy of better than 70 $`\mu `$m, an average track length of 58 cm, and typically measures 15 points along a particle trajectory. The continuous-tracking approximation then yields a mean Maximum Detectable Rigidity (MDR) of 170 GV, and the rigidities of particles up to $``$ 60 GV can be reliably measured.

Antiproton flux measurements require excellent particle identification for background discrimination. The primary sources of background to the $`\overline{p}`$ flux are electrons, and negatively charged muons, pions and kaons produced in the atmosphere as well as in the material above and in the detector. To provide mass discrimination between antiprotons and these particles we measure multiple samples of the ionization energy loss. The logarithmic rise in the mean rate of energy loss for a relativistic charged particle is used to determine the Lorentz factor of the particle from which, together with the rigidity measurement, the mass is obtained.

The multiple dE/dx detector consists of a stack of 140 segmented multi-wire proportional chambers, each providing a measurement of the specific ionization loss. In order to maximize the particle identification power, the proportional chambers are filled with Xenon (and 5% CH4), which exhibits an increase in ionization loss rate of 70% between minimum ionization and relativistic saturation Labrador et al. (1999).

The HEAT-pbar instrument was flown from Ft. Sumner, NM on June 3, 2000. The detector was at float altitude for 22 hours, at an average atmospheric overburden of 7.2 g/cm2. More than 1.9 million events were recorded over an integrated live time of 16.2 hours at a residual pressure between 4.5 and 8.6 mbar. The average vertical geomagnetic cut-off rigidity along the flight path is 4.2 GV. The instrument performed flawlessly during the flight.

Albedo particles, which mimic antiparticles in the spectrometer, can be efficiently rejected with the ToF measurement. The flux of upward-going relativistic particles is roughly 10-3 that of relativistic downward-going particles, and thus a rejection power of 1000 is required to keep the contamination of these particles below 1% in the final antiproton data set. The standard deviation in the velocity distribution for relativistic protons is 0.093c. This results in a rejection power against upward moving particles which is several orders of magnitude better than required.

In analyzing events recorded by the dE/dx system, first suitable selections on tracking quality are made and events are selected for which the incident particle transits the entire dE/dx chamber system. The effect of high-energy tails in the Landau distribution is minimized by computing, for each event, a restricted average specific ionization $`dE/dx_{res}`$ which is the average ionization signal measured by 50% of the dE/dx chambers with signals smaller than the median for this event. Note that the precise value of the selected fraction is not very critical, but that 50% is close to optimal. Fig. 2a and b show histograms of the $`dE/dx_{res}`$ response for negative and positive particles, respectively, having rigidities in the range 4.5–6 GV. In these distributions, the peak at smallest $`dE/dx_{res}`$ corresponds to protons and antiprotons, the next peak to $`\pi ^+`$/$`\mu ^+`$ and $`\pi ^{}`$/$`\mu ^{}`$, and the peak at large $`dE/dx_{res}`$ to $`e^+`$ and $`e^{}`$. Compared to the antiproton flux, the kaon production is small, although not negligible. Our Monte Carlo simulations show that our event selection criteria reduce this contribution to negligible levels, since kaons result from interactions in or near the instrument. A somewhat better mass resolution can be obtained by properly accounting for the dependence of the $`dE/dx_{res}`$ signal on the rigidity within a given rigidity interval, but these histograms, which are representative of the equivalent histograms at higher energy, demonstrate clearly that particle identification is achieved for both antiprotons and positrons. The large sampling for each energy loss measurement produces highly gaussian distributions, and thus for each energy interval, the $`\overline{p}/p`$ ratio can be obtained from fits to the restricted average dE/dx distributions, such as those shown in Fig. 2a and b for the rigidity interval from 4.5 – 6 GV.

The results are summarized in Table 1. In order to obtain the number of protons and antiprotons at the top of the atmosphere, we correct for particle production in the atmosphere (total average column density 7.2 g/cm2) above the instrument. In addition, corrections for interaction and annihilation losses of protons and antiprotons in the atmosphere and in the instrument (maximum column density 5.7 g/cm2 above the lower ToF counter for a vertically traversing particle) are applied. The corrections assume that all particles that interact inside the gondola are rejected by our selection criteria.

The calculated background of secondary antiprotons and protons produced in the atmosphere was based on Pfeifer et al. Pfeifer et al. (1996) for the antiprotons and Papini et al. Papini et al. (1996) for the protons. Interaction and annihilation losses are based on the measured cross sections quoted in Kuzichev et al. Kuzichev et al. (1994) and Denisov et al. Denisov et al. (1973), accounting in detail for the total material traversed by a particle in passing through the atmosphere, aluminum shell, and detector material. The number of antiprotons and protons in each energy bin obtained after applying all of these corrections are shown in Table 1, along with the resulting $`\overline{p}/p`$ ratios. The errors quoted in this table are purely statistical. Systematic errors resulting from uncertainties in correcting the particle numbers in the instrument to the top of the atmosphere and in the background due to particle mis-identification are estimated to be less than 4% of the $`\overline{p}/p`$ ratio.

Our results are shown in Fig. 3, along with previous measurements by others, and a number of recent calculations of the $`\overline{p}/p`$ ratio resulting from secondary production of antiprotons in the interstellar medium. Only recent measurements have been included in this figure Orito et al. (2000); Maeno et al. (2000); Mitchell et al. (1996); Hof et al. (1996); Boezio et al. (1994, 2001).

Many predictions for the $`\overline{p}/p`$ ratio have been published over the years. We show here theoretical curves that are consistent with the now well measured flux ratio in the low energy region around 1 GeV. The result of calculations by Simon et al. Simon et al. (1998) are shown in the figure as a shaded band. The calculations are based on the leaky box model and the uncertainties in the flux prediction, reflected by the band in the figure, are primarily uncertainties in the galactic path length distribution. The dashed and solid line in Fig 3 show the results of calculations by Moskalenko et al. Moskalenko et al. (1998) within a self-consistent CR propagation model. The dashed line represents the case of a proton injection spectrum that is much harder than locally observed, which has been proposed to explain the observed high continuum gamma-ray emission above $``$ 1 GeV Hunter et al. (1997). A standard proton injection spectrum, consistent with the locally observed one, is reflected in the solid line. The sensitivity of the $`\overline{p}/p`$ ratio to the nucleon injection spectrum above a few GeV makes antiproton measurements at energies above a few GeV an important test for CR models. Our data are in good agreement with the ‘standard spectrum’ calculations Moskalenko et al. (1998) at high energy, and do not support an antiproton to proton ratio approaching 10-3 at energies above 20 GeV, in contrast to recent CAPRICE measurements Boezio et al. (2001). Our result does not support models which are based on hard nucleon injection spectra. At energies covered by the measurements presented here, secondary $`\overline{p}`$ production with a nucleon injection spectrum consistent with the locally observed one describes the data well.

The HEAT-pbar instrument is scheduled for additional balloon flights and we expect to statistically improve the data and to further clarify the experimental situation.

###### Acknowledgements.

This work was supported by NASA grants No. NAG 5-5058, No. NAG 5-5220, No. NAG 5-5223, and No. NAG 5-5230, and by financial assistance from our universities. We wish to thank the National Scientific Balloon Facility and the NSBF launch crews for their excellent support of balloon missions and we acknowledge contributions from D. Kouba, M. Gebhard, S. Ahmed, and P. Allison.


text:
# Statistically Preserved Structures and Anomalous Scaling in Turbulent Active Scalar Advection

## Abstract

The anomalous scaling of correlation functions in the turbulent statistics of active scalars (like temperature in turbulent convection) is understood in terms of an auxiliary passive scalar which is advected by the same turbulent velocity field. While the odd-order correlation functions of the active and passive fields differ, we propose that the even-order correlation functions are the same to leading order (up to a trivial multiplicative factor). The leading correlation functions are statistically preserved structures of the passive scalar decaying problem, and therefore universality of the scaling exponents of the even-order correlations of the active scalar is demonstrated.

The old riddle of anomalous scaling in generic turbulent statistics had been solved recently, albeit in the modest context of passive scalar and passive vector advection 01CV ; 01ABCPV ; 01CGP . Motivated by some recent numerical indications 01CMMV , we propose in this Letter that under generic conditions a similar understanding can be extended to active scalar (or vector) advection by a generic turbulent velocity field. The difference between the two problems is that passive fields leave the velocity statistics intact, satisfying a purely linear problem into which the velocity field is injected as an advection term, while active fields influence the statistics of the velocity field. This development brings us therefore closer to understanding anomalous scaling in the nonlinear velocity problem itself.

Passive fields are advected in turbulence with typical dynamics of the form

$$\frac{\varphi (𝒓,t)}{t}+𝒖(𝒓,t)\mathbf{}\varphi (𝒓,t)=\kappa ^2\varphi (𝒓,t)+f(𝒓,t).$$
(1)
Here the field $`\varphi (𝒓,t)`$ can be a passive scalar $`c(𝒓,t)`$ or passive vector $`𝒗(𝒓,t)`$. $`\kappa `$ is the diffusivity and $`f(𝒓,t)`$ is a white random force of zero mean with compact support in $`𝒌`$-space, acting on the largest scales of the order of the outer scale $`L`$ only. The advecting velocity field is a generic turbulent field at high Reynolds number with an extended scaling range, typically exhibiting anomalous scaling. We use the word “generic” here to distinguish from model velocity fields with $`\delta `$-function temporal correlations 01FGV . “Anomalous scaling” means that multi-point correlation functions are homogeneous functions of their arguments, with exponents that cannot be guessed from dimensional analysis. Thus for example the field $`\varphi (𝒓,t)`$ has simultaneous multi-point correlation functions

$$F^{(n)}(𝒓_1,𝒓_2,\mathrm{},𝒓_n)\varphi (𝒓_1,t)\varphi (𝒓_2,t)\mathrm{}\varphi (𝒓_n,t)_f,$$
(2)
where pointed brackets with subscript $`f`$ refer to averaging over the statistics of the advecting velocity field and of the forcing. Anomalous scaling means that

$$F^{(n)}(\lambda 𝒓_1,\mathrm{},\lambda 𝒓_n)=\lambda ^{\zeta _n}F^{(n)}(𝒓_1,\mathrm{},𝒓_n),$$
(3)
with $`\zeta _n`$ having a non-trivial dependence on $`n`$.

In recent work 01ABCPV ; 01CGP  it was clarified why and how passive fields exhibit anomalous scaling. The key is to consider a problem associated with Eq. (1) which is the decaying problem in which the forcing $`f(𝒓,t)`$ is put to zero. The problem becomes then a linear initial value problem, $`\frac{\varphi }{t}=\varphi `$, with a formal solution

$$\varphi (𝒓,t)=𝑑𝒓^{}𝑹(𝒓,𝒓^{},t)\varphi (𝒓^{},0),$$
(4)
with the operator $`𝑹T^+\mathrm{exp}\left(_0^t𝑑s(s)\right)`$, and $`T^+`$ being the time ordering operator. Define next the now time dependent correlation functions of the decaying problem:

$$C^{(n)}(𝒓_1,\mathrm{},𝒓_n,t)\varphi (𝒓_1,t)\mathrm{}\varphi (𝒓_n,t).$$
(5)
Pointed brackets without subscript $`f`$ refer to the decaying object in which averaging is with respect to realizations of the velocity field only. As a result of Eq. (4) the decaying correlation functions are developed by a propagator $`𝒫_{\underset{¯}{𝒓}|\underset{¯}{𝝆}}^{(n)}`$, (with $`\underset{¯}{𝒓}𝒓_1,𝒓_2,\mathrm{},𝒓_n`$):

$$C^{(n)}(𝒓_1,\mathrm{},𝒓_n,t)=𝑑\underset{¯}{𝝆}𝒫_{\underset{¯}{𝒓}|\underset{¯}{𝝆}}^{(n)}(t)C^{(n)}(𝝆_1,\mathrm{},𝝆_n,0).$$
(6)
In writing this equation we made explicit use of the fact that the initial distribution of the passive field $`\varphi (𝒓,0)`$ is statistically independent of the advecting velocity field. Thus the operator $`𝒫_{\underset{¯}{𝒓}|\underset{¯}{𝝆}}^{(n)}`$ can be written explicitly

$$𝒫_{\underset{¯}{𝒓}|\underset{¯}{𝝆}}^{(n)}(t)𝑹(𝒓_1,𝝆_1,t)𝑹(𝒓_2,𝝆_2,t)\mathrm{}𝑹(𝒓_n,𝝆_n,t).$$
(7)
The key finding 01ABCPV ; 01CGP  is that the operator $`𝒫_{\underset{¯}{𝒓}|\underset{¯}{𝝆}}^{(n)}`$ possesses left eigenfunctions of eigenvalue 1, i.e. there exist functions $`Z^{(n)}(𝒓_1,𝒓_2\mathrm{}𝒓_n)`$ satisfying

$$Z^{(n)}(𝒓_1,\mathrm{},𝒓_n)=𝑑\underset{¯}{𝝆}𝒫_{\underset{¯}{𝒓}|\underset{¯}{𝝆}}^{(n)}(t)Z^{(n)}(𝝆_1,\mathrm{},𝝆_n).$$
(8)
The functions $`Z^{(n)}`$ are referred to as “statistically preserved structures”, being invariant to the dynamics, even though the operator is strongly time dependent and decaying. How to form from these functions infinitely many conserved variables in the decaying problem was shown in 01ABCPV . The functions $`Z^{(n)}(\underset{¯}{𝒓})`$ are homogeneous functions of their arguments, with anomalous scaling exponents $`\zeta _n`$. More importantly, it was shown that the correlation functions of the forced case, $`F^{(n)}(\underset{¯}{𝒓})`$, have exactly the same scaling exponents as $`Z^{(n)}(\underset{¯}{𝒓})`$ 01CGP . In the scaling sense

$$F^{(n)}(\underset{¯}{𝒓})Z^{(n)}(\underset{¯}{𝒓}).$$
(9)
This is how anomalous scaling in passive fields is understood.

Consider next the dynamics of an active field. For concreteness we will think about temperature in turbulent Boussinesq convection, but the ideas are immediately generalized to any other type of active field. The equation of motion is formally identical to Eq. (1):

$$\frac{T(𝒓,t)}{t}+𝒖(𝒓,t)\mathbf{}T(𝒓,t)=\kappa ^2T(𝒓,t)+f(𝒓,t),$$
(10)
but now the velocity field is affected by the temperature. For an incompressible fluid of unit density 71MY ,

$$\frac{𝒖}{t}+𝒖\mathbf{}𝒖=\mathbf{}p+\nu ^2𝒖+\alpha gT\widehat{z}.$$
(11)
Here $`p`$, $`\nu `$, $`\alpha `$, $`g`$ and $`\widehat{z}`$ are the pressure, kinematic viscosity, thermal expansion coefficient, acceleration due to gravity and a unit vector in the upward direction respectively. The appearance of $`T`$ in the equation for $`𝒖`$ is crucial, and changes the scaling exponents of $`𝒖`$ from Kolmogorov (approximately) to Bolgiano (approximately) 71MY . It makes no sense now to consider the decaying problem for $`T`$; as this field decays, the statistics of the velocity field changes, and there is very little to say. On the other hand, we can learn a great deal from considering the forced solutions. We note that both the passive and the active forced equations can be solved in the same way:

$`T(𝒓,t)=`$ $`{\displaystyle 𝑑𝒓^{}𝑹(𝒓,𝒓^{},t)T(𝒓^{},0)}`$
$`+`$ $`{\displaystyle 𝑑𝒓^{}_0^t𝑑\tau 𝑹(𝒓,𝒓^{},t\tau )f(𝒓^{},\tau )}.`$ (12)
with a similar expression for the passive scalar. For $`t\mathrm{}`$, the first term decays to zero, but the second term exists. Next we encounter the difference between the active and passive case. Computing the average of $`T`$, or any other odd moment, we get a finite result:

$$T(𝒓,t)=𝑑𝒓^{}_0^t𝑑\tau 𝑹(𝒓,𝒓^{},t\tau )f(𝒓^{},\tau ).$$
(13)
We cannot decorrelate the random forcing $`f`$ from the operator $`𝑹`$ which involves the velocity field. The forcing $`f`$ is correlated with $`T`$, which is itself correlated, via Eq.(11), with $`𝒖`$, and therefore the field $`T`$ has a first and higher odd moments: the probability density function (pdf) of $`T`$ will be asymmetric. For the passive case, the equivalent Eq. (13) can be simplified by decorrelating $`f`$ from $`𝑹`$. As the odd moments of $`f`$ vanish, the passive scalar has zero odd moments, and its pdf is symmetric. This difference is fundamental.

We now argue that all the even moments, and all even correlation functions of the active and passive fields may be identical to leading scaling order. Consider for example the second order correlation function $`F^{(2)}(𝒓_1,𝒓_2)`$ of the active scalar:

$`F^{(2)}(𝒓_1,𝒓_2)={\displaystyle 𝑑𝒓^{}𝑑𝒓^{\prime \prime }_0^t𝑑\tau ^{}𝑑\tau ^{\prime \prime }}`$
$`𝑹(𝒓,𝒓^{},t\tau ^{})f(𝒓^{},\tau ^{})𝑹(𝒓,𝒓^{\prime \prime },t\tau ^{\prime \prime })f(𝒓^{\prime \prime },\tau ^{\prime \prime }).`$ (14)
As before, in the case of the active scalar the forces cannot be decorrelated from the operators. For the passive scalar this can be done with impunity, leading to the correct form of the correlation function, i.e.

$`\varphi (𝒓_1)\varphi (𝒓_2)_f`$ $`=`$ $`{\displaystyle 𝑑𝝆_1𝑑𝝆_2_0^t𝑑\tau 𝒫_{𝒓_1,𝒓_2|𝝆_1𝝆_2}^{(2)}(t\tau )}`$ (15)
$`\times `$ $`f(𝝆_1)f(𝝆_2).`$
It was shown in 01CGP  how this expression indeed yields the forced correlation function in agreement with Eq.(9). We observe however that any correlation function $`AB`$ can be written in terms of its Gaussian decomposition and the “connected part”: $`AB=AB+\stackrel{~}{A}\stackrel{~}{B}`$, where $`\stackrel{~}{A}AA`$. Identify now $`A`$ with the product of operators $`𝑹`$ in Eq. (14), and $`B`$ with the product of forcing, and realize that Eq.(15) is the Gaussian decomposition contribution to Eq. (14). So the solution (15) exists as one contribution to (14), the other being what remains of the correlation, i.e. the part in which the forcing cannot be decorrelated from the operator. A similar phenomenon occurs for all the even order correlation functions. But we know from the theory of the passive scalar that the solution (15) and the equivalent ones for all the even-order correlation functions are identified with statistically preserved structures whose scaling is anomalous 01CGP . We propose that the contributions to the even-order correlation of the active scalar which are identical to the even-order correlations of the passive scalar are leading, and are therefore expected to be realized in experiments and simulations. The expectation is that the connected contributions, that are purely “forced” solutions, are either subleading or have the same exponents. This is a crucial conjecture that needs to be tested against examples. In 01CMMV  it was shown that this is correct for a particular 4th order correlation function of the active and passive fields.

To exemplify the generality and plausibility of these ideas we examine in detail a model of active and passive scalar for which the statistical object can be computed to high accuracy. We consider a variant of the shell model studied in ref. Brand :

$`{\displaystyle \frac{T_n}{t}}`$ $`=`$ $`\stackrel{~}{a}k_n(u_{n1}T_{n1}hu_nT_{n+1})+\stackrel{~}{b}k_n(u_nT_{n1}hu_{n+1}T_{n+1})\kappa k_n^2T_n+f\delta _{n,0}`$ (16)
$`{\displaystyle \frac{u_n}{t}}`$ $`=`$ $`ak_n(u_{n1}^2hu_nu_{n+1})+bk_n(u_nu_{n1}hu_{n+1}^2)\nu k_n^2u_n+T_n`$ (17)
$`{\displaystyle \frac{C_n}{t}}`$ $`=`$ $`\stackrel{~}{a}k_n(u_{n1}C_{n1}hu_nC_{n+1})+\stackrel{~}{b}k_n(u_nC_{n1}hu_{n+1}C_{n+1})\kappa k_n^2C_n+f\delta _{n,0}.`$ (18)
In this model $`n`$ stands for index of a shell of k-vector $`k_n=k_0h^n`$, with $`n=0,1,\mathrm{},N1`$. We take $`h=2`$, and the parameters used in the simulation are $`a=0.01`$, $`\stackrel{~}{a}=\stackrel{~}{b}=b=1`$, $`k_0=1`$, $`\kappa =\nu =5\times 10^4`$. The number of shells is $`N=30`$.

Without the coupling to $`T_n`$, the velocity equation has an inviscid unstable Kolmogorov fixed point, $`u_nk_n^{1/3}`$. This is changed by the coupling Brand , and the system of equations for $`T_n`$ and $`u_n`$ exhibits an inviscid unstable Bolgiano fixed point, $`u_nk_n^{3/5}`$, $`T_nk_n^{1/5}`$. The chaotic dynamics renders the statistics of the velocity strongly non-Gaussian, cf. inset in Fig 1. The exponents for the active scalar are markedly anomalous, whereas for the velocity appear closer to normal, see Fig. 1.

To demonstrate the first point of the fundamental difference between the active and passive fields, we show in Fig. 2 the pdf’s of $`x=\stackrel{~}{\varphi }_n/\stackrel{~}{\varphi }_n^2^{1/2}`$ where $`\stackrel{~}{\varphi }`$ is $`\stackrel{~}{T}_n`$ or $`C_n`$, for $`n=14`$. One clearly sees the symmetry of the pdf of the passive scalar, in contradistinction to the asymmetry of the pdf of the active scalar. This is typical to all $`n`$ in the inertial range. This is a demonstration of the discussion after Eq.(13). For the passive scalar the odd moments vanish, whereas for the active scalar they all exist.

The situation is altogether different for the statistics of even moments. To demonstrate the difference we plot in Fig. 3 the (typical) pdf of $`\stackrel{~}{T}_n^2`$ and $`C_n^2`$ for $`n=9`$ and 14. In plotting we realize that the passive scalar is defined up to a constant, so for the passive scalar the pdf is plotted for the rescaled variable $`\beta C_n^2`$, where $`\beta =\stackrel{~}{T}_n^2/C_n^20.6327`$. Note that there is only one numerical freedom $`\beta `$, constant for all $`n`$ in the inertial range. We find very close agreement of all the pdf’s in the inertial range.

The identity of the pdf’s of $`\stackrel{~}{T}_n^2`$ and $`C_n^2`$ translates automatically to the identity (up to a constant $`\beta ^m`$) of the even-order structure functions $`S^{(2m)}(k_n)\stackrel{~}{\varphi }_n^{2m}`$, where $`\stackrel{~}{\varphi }_n=\stackrel{~}{T}_n`$ or $`C_n`$. This is demonstrated in Fig. 4. We see that the 2nd, 4th and 6th-order structure functions are barely distinguishable, with the same scaling exponents in the inertial range.

Finally, the identity of the statistics of the squares of the passive and active scalars transcends structure functions. We consider next multi-point correlation functions, and in Fig. 5 compare the correlation functions $`\stackrel{~}{T}_n^2\stackrel{~}{T}_{n+5}^2`$ and $`\stackrel{~}{T}_n^2\stackrel{~}{T}_{n+5}^2\stackrel{~}{T}_{2n}^2`$ to their passive counterparts.

The conclusion is that again the multi-point correlation functions are indistinguishable once the passive ones are rescaled by $`\beta ^q`$ where $`q`$ is the over-all order of the correlation function. In summary, we offered a discussion and detailed numerical evidence to support the conjecture that under generic conditions the even-order correlation functions of active scalars can be understood via the emerging theory of Statistically Preserved Structures of the passive scalar counterpart. We have by now two examples, together with 01CMMV , which share however a similar coupling of the active scalar to the velocity field. The importance of this conjecture warrants considerable further work to uncover its provisos and delineate its generality.

###### Acknowledgements.

This work had been supported in part by the Research Grants Council of Hong Kong SAR, China (CUHK 4119/98P and CUHK 4286/00P), the European Commission under a TMR grant, the German Israeli Foundation, and the Naftali and Anna Backenroth-Bronicki Fund for Research in Chaos and Complexity. TG thanks the Israeli Council for Higher Education and the Feinberg postdoctoral Fellowships program at the WIS for financial support.


text:
# Numerical simulation of interference experiments in a local hidden variables model.

## I Introduction

The problem of interference has haunted modern physics practically from the beginning of the twentieth century. Even though experiments, where a time ordered sequence of single impacts has been recorded, have only recently been developed , the thought experiments predicting the very same outcome ran into considerable conceptual problems much earlier. For a comprehensive review of the problem see Franco Selleri’s book .

The problem can be stated in the following way. Consider a double slit interferometer. Only one particle is thought to propagate from the source to the screen. Therefore the single particle passes a single slit without any interaction with other particles, and impinges on the screen. However, even though the particle cannot be physically influenced by the existence of the second slit, its impact probability at certain spots of the screen depends explicitly on the existence of this second slit. How then, can the particle avoid these spots even though there seems no physical means to that end? For this reason, we do not have a single theoretical model at present, which describes the impact of this particle by a consistent, time ordered sequence of events, along the particle’s path from the source, through the slits and onto the screen, where all interactions are local. The only model, which is capable of such an analysis, is non-local. It is based on David Bohm’s theory of hidden variables . However, this theoretical model depends explicitly on potentials, changing the trajectory of the particle between the slits and the screen, it is therefore unsuitable to describe photon interference measurements.

We wish to close this gap in the present paper by a model, which describes the location of every single photon in the interferometer at every single moment. Moreover, it allows to recover all the results of the experiments while remaining perfectly local. The model has been implemented in a numerical simulation scheme. We perform these simulations mainly for photons, even though a slightly modified model could also account for the recent C60 single-slit interference experiments. The aim of our presentation is to demystify the problem by stating very clearly and in unambiguous terms how the statistical picture arises from limited control in the experiments, and how the physical interactions, which are thought to be deterministic, come into play. It will be seen that the statistical picture underpins the usual treatment in classical field theory, e.g. in the scalar Kirchhoff model, and that it is the very same picture which has been adopted for quantum mechanical purposes. We make one modification to the current understanding of this process, i.e. that the phase of a photon’s field does change discontinuously during the scattering at the slit. And we generally limit the model to the experimentally verified situations, that is to the high intensity regime.

However, the low intensity regime, as in electron interference measurements or recent measurements with molecules, is also discussed extensively. It is shown that in one case one obtains identical results, if the single impacts follow a Gaussian distribution around the actual slits. In the second case, the molecule interference experiments, our investigation is still incomplete but we can make it plausible that a similar model to the one suggested here also accounts for these experiments.

## II Theoretical model

In our theoretical description of photon scattering we use essentially Richard Feynman’s path integral method, with a few but decisive modifications. To incorporate the phase of a photon along a given path we describe photons as electromagnetic fields of limited extensions. The compatibility of such a photon model with quantum mechanics has been shown quite some time ago . We omit higher Fourier coefficients, which are due to the confinement of the field, and consider every photon as essentially monochromatic. The field of the photon possesses an arbitrary phase $`\varphi _0`$ at the point of emission. This arbitrariness is due to our lack of knowledge about the individual photon. To simplify the model we assume that the fields of all photons are polarized in the same direction, arbitrary polarization only has the effect that the intensity of the beam will have to be higher to obtain the same effects. The phase of the electromagnetic field of the single photon varies along its optical path from the source to the slit, a difference of $`\mathrm{\Delta }\varphi =2\pi `$ describes two separate points with a distance $`\lambda `$, the wavelength of the photon. The phase at the slit of our interferometer is therefore $`\varphi _1=\varphi _0+2\pi L_1/\lambda `$, where $`L_1`$ is the distance between the source and the slit.

The scattering of single photons at the slit environment is described by a statistical model. Since currently no theoretical treatment of the problem exists, which would allow us to relate the interaction of the photon with the slit unambiguously to a direction, into which the photon is scattered, we formalize the process in the following way:

The phase $`\varphi _1`$ of the photon together with a random phase $`\overline{\varphi }_1`$ of the slit environment shall trigger a process which changes the direction of the photon to a scattering angle $`\theta _1`$. The angles of all photons measured will comply with a Gaussian distribution centered around zero. In practice we have simulated two photons simultaneously, using the phase of the second photon $`\varphi _1(2)`$ as the input from the slit environment $`\overline{\varphi }_1(1)`$, and vice versa. Since both phases are random, this procedure does not couple the photons in any way, it merely guarantees that the angular distribution after the slit is Gaussian. The halfwidth of the distribution is one of our input parameters. The relevant lines in our code are shown in Appendix A.

Clearly, if the photon changes its direction, then it has interacted with the slit environment. Since we have no deterministic model of this process, we also do not know, what the phase of the photon after the interaction will be. To understand in more detail the nature of this problem, consider the classical scalar model of interference experiments. Here, the scalar field $`u(𝐫)`$ is thought to be emitted from a point like source. After scattering at the slits it is adsorbed at the screen. Throughout the interferometer the field is at every point in phase with the field at the source. It is thus only by a rigid connection of field amplitudes via the phase that interference effects occur at all within this model.

While such a model works perfectly well for laser beams with their high coherence length and where all photons can be thought to be in phase, it is less appropriate for ordinary monochromatic light. These emissions originate from thermal processes, therefore their phases must generally be random, even though their frequency, e. g. from the decay of an electron excitation in an atom, can be sharply defined. Therefore we explore a different possibility in this paper, and we shall show, that also for a statistical distribution of initial phases interference effects can be obtained, provided the scattering process of a single photon at one slit changes the phase of the photon discontinuously. It goes without saying that the results are the same, if we use monochromatic laser light instead of monochromatic light from emissions of thermally excited atom.

In our model, as indeed in every model which constructs interferences from uncorrelated photons, the single photon will have an arbitrary phase at the source. Moreover, the phase of the photon is not connected to the phase of any other photon. Therefore, if we adopt the view that the field of the single photon remains in-phase with the field at its source, we will not describe any interference effects. We conclude from this fact that the field cannot remain in phase. Since the photon interacts only at one point of our system, at the slit of the interferometer, it can only be this point, where the phase changes.

Here, we make the only truly novel conjecture in this paper: the phase angle $`\varphi _1^{}`$, after the scattering process, shall be proportional to the angle, into which the photon is scattered:

$$\varphi _1^{}=\alpha \theta _1$$
(1)
For reasons of consistency we choose the constant $`\alpha `$ equal to one. The conjecture can then be linked to the description of field deflection in classical electrodynamics, where a reflection of a field at a surface leads to a change of its phase by $`\pi `$ . Since in classical field theories the phase is continuous throughout the system, the phase shift is related to changes throughout the system. In case of arbitrary initial phases we can only recover this feature, if we set the phase itself equal to the angle of scattering. In our view this model is a generalization of classical concepts to the case, where we have given up the rigid connection of fields throughout the system.

The photons are scattered at a single slit and propagate in a straight line towards the detector screen. Here, the fields of individual photons are superimposed and the total field at a point P(x,y) is given by the expression:

$$E(x,y)=\underset{i=1}{\overset{n}{}}cos(\varphi _2(i))$$
(2)
where n is the total number of photons impinging at this point, and $`\varphi _2`$ is the phase at their moment of impact. The phase is related to the scattering angle $`\theta _1`$ by the expression:

$$\varphi _2=\theta _1+2\pi \frac{\sqrt{L^2+y^2}}{\lambda }\frac{y}{L}=tan\theta _1$$
(3)
Here $`L`$ is the distance between the plane of the slits and the plane of the detector screen. The details of the model are displayed in Fig. 1. The intensity at the point P(x,y) of the detector screen is consequently:

$$I(x,y)=E(x,y)^2$$
(4)
The intensity we observe is thus the intensity due to the total field after superposition of all contributions. It is for this reason that the model is generally limited to the high intensity regime. However, we will indicate in the discussion, how a slightly modified version of the same model can also account for single molecule interferences. Furthermore, we shall demonstrate that single electron experiments, where a periodic image is gradually evolving, can be understood without the occurrence of any superpositions from separate slits and thus without interferences proper.

## III Numerical implementation

The theoretical model requires a random input at two points of the interferometer: (i) at the emission of a photon from the source; and (ii) at the scattering of the photon at the slit environment. We have used a random number generator for both processes , mapping the random number  onto the initial phase \[$`0,2\pi `$\]. The dimensions of the apparatus were chosen similar to an experimental one. We also chose monochromatic visible light of $`\lambda =500`$nm. The distance $`d`$ between two slits of the interferometer equals 0.5 $`\mu `$m. The distance $`L`$ between the source and the slit and the slit and the screen in our simulations was 10mm (see Fig. 1).

The simulation of an experimental run is straightforward. For a pair of photons, we create the initial phases $`\varphi _0(1,2)`$, and compute the phase after covering the distance to the slit. Here, we use the phase of the second photon $`\varphi _1(2)`$ as the random phase of the slit environment for photon one and vice versa. The two random inputs are used to create two random outputs, which now comply with a Gaussian distribution (see Appendix A). This Gaussian distribution describes the scattering angle of the photons after the slit. The phase angle of the photons at this point is set equal to the scattering angle, and both photons propagate until they hit the screen. We record the points of impact and the phases at the moment of impact, before simulating the next pair. In all our simulations we simulated 100000 photons per slit. We have assumed that the slit is very narrow, in practice we have always used the same point as the slit position. The simulation method has been tested for slits of variable width, but apart from a blurred image due to the variable position of the photons in the plane of the slits we did not obtain any new feature compared to the results obtained with a very narrow slit. Since the slit width is generally unknown we have omitted this part of the numerical analysis.

To evaluate the photon impacts we have assumed that photons interact if their distance is less than 10 nm on our detector screen. This is no limitation of the generality of our model, since we could always increase the number of photons and thus allow for higher resolution. The total field at each of these separate points was calculated and then smoothed with a Gaussian of 0.1 $`\mu `$m halfwidth. The plots of the number of photons and the intensity are shown in units of $`I_0`$, the maximum intensity, and $`N_0`$, the maximum number of impacts.

We have also accounted for a numerical problem in the transition from single slit diffraction to multiple slit interferometry. In a single slit, and assuming that the scattering angle is equal to the phase angle, the intensity distribution of the photons on the screen is essentially due to the intensity distribution of the photon fields. All photons, at a given point, will have the same phase. This is not the case in a multiple slit interferometer. Here we have to sum up contributions from different slits, and the total intensity is due to the total electromagnetic field. The critical points at the screen, in a numerical model, are the points where the phase difference between separate slits yields constructive interference, while the phase from one single slit, or the effect of photon diffraction, yields a minimum of the field. For a double slit interferometer, the $`\mathrm{cos}\varphi _2`$ at the median point between two slits is zero (because $`\varphi _2\pi /2`$), while the interference between the two separate slits is constructive (because both phases are equal). Since it is largely a matter of convention in field physics, whether a constant phase is added in all calculations, we have used the same approach in this numerical model. For multiple slit simulations we therefore have increased the final phase $`\varphi _2`$ by a value of $`0.2\pi `$ and for every single photon.

## IV Results

The main results of our simulations are shown in Figs. 2 to 6. In all plots we show the normalized impact numbers (dash dotted grey curve) as well as the normalized intensity distribution (black curve). The position of the slits is indicated by a dash dotted vertical line. We also show a simulated greyscale image of the detector screen, this filmstrip describes how the image would look in a real interferometer.

The height of the peaks depends, for a given number of slits, essentially on the halfwidth $`\sigma `$ of the scattering angles. Since this parameter reflects the actual experimental situation, or rather the unknown scattering process in the slit environment, we cannot set it to a theoretically required value. For this reason we can only change it so that the results of our simulation are in accordance with experimental data. In consequence the absolute value of a specific peak reflects to a certain extent our choice of parameters. But the change of the peak height, e.g. if a double slit and a four slit image are compared, does not depend on our settings, it is a genuine consequence of the theoretical assumption that the interference pattern is built up by superposition of the electromagnetic fields of independent photons.

Analyzing the difference between a single slit diffraction and double slit interference (Fig. 2 and Fig. 3), it becomes immediately clear why, in the first case, the peak is at the position of the slit, while in the second case it indicates the median position between the two slits. This is simply due to the interference of photons passing through different slits, which leads to destructive interference at both slit positions, thus creating the distinct intensity minima. In Fig. 4 we show, how the distance between the slits and the halfwidth effect the intensity distribution. For this image we have increased the distance between the slits to $`2\lambda `$, the halfwidth $`\sigma `$ is now 1.5. The number of observed peaks in this case is substantially increased. The image is equal to the experimental images obtained in a Young interferometer .

The effect of destructive interference leads to an interesting phenomenon, if the fringes for two and four slits are compared. Because in the latter case (Fig. 5) the peak at the center of the screen is reduced, it seems thus that the opening of two new slits effectively decreases the number of impinging photons. As the number distribution reveals, this is clearly not the case. The intensity is reduced due to destructive interference of the electromagnetic fields carried by each individual photon. Generally speaking the intensity peaks become more and more similar with an increasing number of slits, so that in the limit of infinite slits, or an infinite chain of slits, the intensity distribution would be purely periodic while the distance between the peaks is equal to the distance between the slits. This is, what we actually observe in diffraction experiments with either photons (x-ray diffraction), or electrons (electron diffraction).

Comparing the distribution of photon impacts with the intensity distribution, it is clear that the main feature, which creates the interference pattern, is the loss of intensity at certain points of the detector screen. These points are due to destructive interference, in the multiple slit interferometer, or to minima of the photon fields, in the single slit interferometer. We take this as a general result of our simulation: periodic intensity distributions in a process, or interference fringes, are always due to a loss of intensity at specific points.

Our records of the simulation contain the information about each individual photon: its trajectory and its phase at the moment of impact. Information, which according to the current understanding prohibits the interference of photons from different slits (this is the principle of the so called which-path or welcher Weg experiments). Here, we retain all the information and still find an interference pattern.

This is due to the fact that the vanishing interference pattern in a Welcher Weg experiment is no consequence of the information about the path of the photon, but, in the way the experiment is often implemented, due to the vanishing superposition of the photon’s electromagnetic fields. In practice, such an experiment can be carried out with polarized photons in a two slit interferometer, if photons of horizontal (vertical) polarization pass only through the first (second) slit. If we simulate this situation, then we have to account for the vector characteristics of photon polarization. This, in turn, means that the electromagnetic field of the horizontally polarized photons is perpendicular to the electromagnetic field of the vertically polarized photons. And then the superposition of the fields occurs independently, and in two dimension, the two directions perpendicular to the photon’s vector of motion. The result of this simulation is shown in Fig. 7. The two single slits are still clearly visible, the maxima, furthermore, are at the position of the slits, clearly photons from one slit have not interacted with photons from the other slit. The typical two slit interference pattern can be made to reappear by inserting a diagonal polarizer into the path of all photons: since all photons then have the same vector of polarization, the superposition occurs again only in one dimension and the interference pattern is recovered.

## V Discussion

The local model has two key ingredients: a local scattering event at a single slit, where the photon is scattered into a new direction, and a local superposition at the focal plane of the interferometer, where the electromagnetic fields of single photons are superimposed. The second feature, the local superposition of fields, is also the main feature of classical models of interference. It is, as David Deutsch remarked in his book , also the main issue in quantum mechanical models, even though not every photon can be thought to be real in every situation. One has to conclude, therefore, that no model exists, which predicts interference effects from the scattering events alone. This would imply a non-local connection between the individual slits. Theoretically speaking we have no reason to believe that such a model is physically meaningful. Therefore, even given the present results, the conceptual problem in the low intensity regime still seems to prevail. To our knowledge these low intensity experiments so far have only been performed with electrons, and C60 molecules . We are currently investigating the experiments with C60 molecules and the emerging picture is roughly the following: due to its high energy the molecule possesses a considerable dipole moment, which oscillates along its path. This oscillating moment plays the same role as the periodic electromagnetic field of a photon. In the slit environment the molecule is scattered into a new direction due to interaction between the molecule and the atoms of the slit. The scattering angle depends on the phase of the molecular dipole moment. The molecule can only be detected after ionization. Since the ionization process will not be successful in every case, not all molecules will actually be detected. And if the success rate of ionization depends on the dipole moment of the molecule, then the angle of scattering is linked to the probability of detection. We will present the exact numerical results for this experiments in a future publication.

The impact of single electrons gradually building up the interference pattern of the statistical theory, as in Tonomura’s 1989 experiments, seems to be the strongest proof, that a single electron undergoes interference on the screen, even though there is no other electron present, which could do the trick. In the view of the Many Worlds interpretation of quantum mechanics , the electron in this case interacts with an electron from a parallel universe, a ’shadow’ electron. As already mentioned, there is so far no unambiguous proof that the same behavior can be seen for single photons: photon interference measurements have only be performed in the high intensity regime. This makes electron interferences all the more important. And, if these experiments are correct, then it could be seen as a demonstration of the validity of the Many Worlds interpretation. However, an identical result can be obtained, if the single electrons do not interact with ’shadow’ electrons, but if the angular distribution of scattered electrons after the slit system is very narrow. So narrow, in fact, that the overlap between the impact of electrons from separate slits is close to zero. We have simulated such a measurement with photons, an identical result could be obtained for electrons if the lengthscale is reduced by three orders of magnitude. As seen in Fig. 8, the peaks of the number distribution coincide with the peaks of the intensity distribution. This is fundamentally different from the previous cases, where the intensity peaks were always between the positions of single slit. It is, in short, the proof that there is no interference at work under these conditions. The periodic distribution at the screen is merely an image of the periodic arrangement of slits of the interferometer.

Following this line of research, it is also possible that the experiments with single C60 molecules  are no indication of any phase effects at all, but simply due to a specific choice of experimental parameters. If we assume that the initial molecular beam before the slit system is wider than the distance between the single slits, so that molecules will pass through the central slit and two side-slits, then the image on the screen will show one main peak between two smaller peaks. Apart from determining the absolute position of the image in relation to the central slit, there is no way to decide, whether we deal with interference, or just a statistical distribution of molecule impacts due to their passing through the slit system. We have simulated this situation by counting the impacts of single photons, if 40% of the photons pass through a side-slit. The Gaussian in this case was assumed to be narrow enough to retain the side peaks as distinct features. This can essentially be obtained in an experiment with molecular beams by making the molecular velocity high enough. The result of this simulation is shown in Fig. 9. It can be seen that the essential features of the experimental image - one main peak, two side peaks, the visibility of the minima very low - are recovered if the parameters are suitably chosen. Also in this situation we do not deal with interference, but merely with the impact of particles scattered at the slit system.

This simple model can easily be checked experimentally. Since the occurrence of ”interference fringes” depends on the statistical spread and thus the velocity of the C60 molecules, the fringes should disappear, if the velocity decreases below a certain threshold. In fact, the predictions of the model can be tested by varying the velocity. We predict that the visibility depends on molecular velocity, while in interference experiments proper the position of the maxima and minima depends on the velocity. If the intensity distribution is only due to kinetic effects, then the positions will not be altered by a variation of velocity.

Let us go back to our analysis of the double-slit interferometer in the introduction and to the problem, how photons can avoid regions of low intensity on the detector screen. It becomes clear now, that this is actually not the case. The number density of photons on the screen in all cases of physical interference is more or less Gaussian. The periodicity of the intensity, the observed interference pattern, is only due to interference between the fields of individual photons. The number density, therefore, is generally not related to the intensity distribution. It appears thus that quantum mechanics, which arrives at an oscillating number density (the probability distribution $`|\psi |^2`$) to account for the intensity in interference experiments, gives a misleading and even wrong account of the actual physical situation.

## VI Experimental tests

The conjecture used to construct the local model also for a random distribution of photon phases has interesting consequences, which can in principle be checked experimentally. Consider an experiment with polarized monochromatic light. Then after the slit the direction into which a photon is scattered is related to the phase after the interaction with the slit environment. Therefore photons at different points of the detector screen will have different phases. The difference between e.g. the first and second maximum in a double slit interferometer (see Fig. 4) corresponds to a phase difference of

$$\mathrm{\Delta }\varphi \pi $$
(5)
For polarized photons, this difference can be measured with an optical modulator, since the phase at the front end of the modulator determines the rotation of the field vector . Subjecting photons after the modulator to a polarization measurement should yield a different result for photons impinging at different points of the detector screen.

## Acknowledgements

A number of people have shared my interest in interference phenomena during the last year, most notably R. Stadler, who contributed more to this paper then he is willing to take credit for. I am also obliged to J. Gavartin, he provided the code for the numerical scattering model. And finally I have to thank the UCL HiPerSPACE center for providing the computer infrastructure, which made the simulations possible.

## Appendix A

The actual code fragment, which transforms a uniform distribution  into a Gaussian distribution $`\mathrm{exp}(x^2/2)`$ includes only 15 lines. The main features are: first create a random number PHASE1, add the phase from the propagation to the slit, and throw away multiples of one (or of $`2\pi `$): this random number is now transformed from the interval  (PHI1) to the interval \[-1,1\] (V1):

10 PHASE1 = RANDOM()

PHI1 = PHASE1 + Z

K = PHI1

PHI1 = PHI1 - FLOAT(K)

V1 = 2. * PHI1 - 1.

Now do the same for a second number PHASE2, which results in V2:

PHASE2 = RANDOM()

PHI2 = PHASE2 + Z

K = PHI2

PHI2 = PHI2 - FLOAT(K)

V2 = 2. * PHI2 - 1.

Calculate the radius and see, whether it is inside the halfwidth of the desired Gaussian: if it is, then compute the final values, if it is not, then choose two new random numbers:

RADIUS = V1\*V1+V2\*V2

IF (RADIUS.GT.1.0) GOTO 10

FACTOR = SQRT(-2.\*LOG(RADIUS)/RADIUS)

GAUSS1 = V1\*FACTOR

GAUSS2 = V2\*FACTOR

These two variables, GAUSS1 and GAUSS2, which comply with a Gaussian of halfwidth one, are then multiplied by our input parameter $`\sigma `$, the chosen Gaussian halfwidth of our angular distribution. In the simulations the halfwidth was set to $`6\times 10^5`$.


text:
# THE LUMINOSITY DEPENDENCE OF CLUSTERING

## 1 Introduction

Ever since the first galaxy surveys it has been clear that the galaxy distribution is not uniform and that the measured clustering is dependent on the properties of the sample in question. Differential clustering as a function of morphological type, colour, luminosity and wavelength probed have all been measured and imply the existence of biases between the distribution of galaxies and mass. It is this complex distribution of galaxies which galaxy formation models have tried to reproduce. In the era of new large redshift surveys it is now possible to test to much higher accuracy the large scale structure predictions of these models. This talk concentrates on the luminosity dependence of clustering which is perhaps the most controversial of these measurements.

Typical galaxy formation models start with a smoothed random density field, which is left to evolve. It is assumed that galaxies will form at the peaks of this density field and that the mass of the galaxy is related to the peak height. If we then assume that galaxy luminosity is related to galaxy mass we come to the conclusion that brighter galaxies are found at the higher peaks. As higher peaks are more strongly clustered , more luminous galaxies should also be more strongly clustered. Hence these galaxy formation models predict the amount of variation in galaxy clustering as a function of absolute magnitude that we should observe. They all show a strong variation over a wide range of magnitudes with a large increase of clustering strength for the very bright galaxies . Although the models do disagree on the exact amount of variation, this effect should be large enough to be observable and the observations may allow us to select the more realistic models.

Many previous studies have claimed to have observed this effect (e.g. ) but there have also been plenty of non-detections (e.g. ). This controversy is likely due to the limited sizes and small number of galaxies in past surveys which have restricted the accuracy of the measurements. I use two recent, large, but very different redshift surveys to try and resolve this controversy.

## 2 Using the Infra-Red: The PSCz survey

The publicly available PSCz survey  was selected from the IRAS Point Source Catalogue and contains nearly 15000 galaxy redshifts to a limiting $`60\mu `$m flux ($`f_{60}`$) of 0.6Jy. The median redshift is $`z0.03`$ and the survey covers about $`84\%`$ of the sky (Figure 1). The main incompleteness is at low galactic latitudes and is due to galactic extinction which is masked out in the analysis. The redshift-space correlation length, $`s_0`$, is measured using standard methods for the whole survey and for various sub-samples of the catalogue.

Galaxies with a high star formation rate (SFR) tend to be bright in the far infra-red (FIR) because of the thermal emission from dust heated by young stars. They will hence have warmer FIR colours and so splitting the survey by colour will test the dependence of clustering on star formation rate.

Figure 2a shows the results for three sub-samples selected on their $`f_{100}/f_{60}`$ colour. The dashed line is a best fit to these results and it is clear that the clustering amplitude is less for the galaxies with a lower $`f_{100}/f_{60}`$ ratio, and hence warmer colour. This implies that these warmer, more star-forming galaxies are found outside the dense regions of clusters which is confirmed by the observations that later type optical galaxies are found less clustered than earlier galaxy types.

The results for two luminosity selected samples and a range of volume-limited samples are plotted in Figure 2b. The best fit to the colour sample results is converted to a luminosity prediction via the observed PSCz colour-luminosity relation. This is the dashed line in Figure 2b and it can be seen that it is consistent with the data points. The dotted line is the prediction of a typical galaxy formation model described above and is clearly inconsistent with the data. The conclusion is that the assumption made in the galaxy formation models that galaxy luminosity is related to galaxy mass is not valid in the far infra-red, and that a galaxy’s far infra-red luminosity is dominated by its SFR .

## 3 Using the Optical: The 2dF Galaxy Redshift Survey

The 2dF Galaxy Redshift Survey was selected in the $`b_J`$ band from the APM survey and has a limiting apparent magnitude of $`b_J=19.45`$. The survey covers approximately 2000 square degrees and is split into two broad strips, one in the South Galactic Pole region (SGP) and the other in the direction of the North Galactic Pole (NGP). The aim of the survey is to observe a total of 250,000 galaxies and as the survey is still on-going, the results that follow are from an interim study. To date over 170,000 unique galaxy redshifts have been measured, with a median redshift of $`z0.1`$. The survey is described in detail in , and the preliminary data release is available at http://www.mso.anu.edu.au/2dFGRS/. The large number of galaxies, over a large area of the sky, makes the survey ideal for measuring the large scale structure of the Universe. It is also large enough to allow accurate estimates of the clustering to be measured for sub-samples of the data, adding extra dimensions to the output of the survey. It is this fact which makes it ideal for looking for the dependence of clustering on galaxy properties.

The clustering is measured  using the projected correlation function, $`\mathrm{\Xi }(\sigma )`$, which is obtained by projecting $`\xi (\sigma ,\pi )`$ onto the $`\sigma `$ axis, giving a correlation function measurement free from redshift-space distortions. The real-space correlation length, $`r_0`$, can then be measured from $`\mathrm{\Xi }(\sigma )`$.

Figure 3 shows the projected correlation function results for a single volume limited sample defined to include all galaxies with $`21.5<M_{b_J}5\mathrm{log}h<18.5`$. This sample is further split into three disjoint one magnitude slices. The brightest bin (solid line) clearly shows stronger clustering than the fainter two bins over a wide range of scales. The error bars come from the scatter between mock catalogues. Another feature of this plot is that the slopes of the measured correlation functions are very similar. Using a single volume is an important robustness test as the volume probed is the same for each sample and hence the cosmic variance is the same for each line. The disadvantage is that the brighter bins contain a lot fewer galaxies.

To overcome the lack of galaxies in the brighter bins of a single volume, a series of half-magnitude wide volume limited samples are used. They are created by changing the absolute magnitude limits to maximise the number of the very bright galaxies in the brightest bins. Thus the clustering can be measured as a function of median absolute magnitude. The disadvantage here is that different volumes are then being analysed and so could introduce cosmic variance.

Figure 4a shows the estimates of $`r_0`$ measured independently in the NGP and the SGP for the range of volume limited samples. The results are consistent within the errors and is another important test of the robustness of the clustering effect. In Figure 4b the NGP and SGP measurements have been combined to give the overall 2dFGRS results and these are compared to the results from other surveys. The results confirm the rapid increase of clustering for samples with $`M`$ brighter than $`M^{}`$. The error bars are smaller than previous measurements even though the scatter between mock catalogues is used. The Poisson and bootstrap re-sampling methods of previous surveys will tend to underestimate the errors.

Although the increase in clustering appears sharp in Figure 4, Figure 5 shows the relative bias, $`b/b^{}=(r_0/r_0^{})^{0.5}`$, which is measured as a function of luminosity and compared to the results of . The 2dFGRS results are actually well fit by the linear relation,

$$b/b^{}=0.85+0.15L/L^{}.$$
(1)
## 4 Summary

$``$ Cold IRAS galaxies are found to be more clustered than hot IRAS galaxies.
$``$ Bright and faint IRAS galaxies are found similarly clustered, consistent with the colour relation.
$``$ Bright optical galaxies are found to be significantly more clustered than faint optical galaxies and the measurement presented is the most accurate yet made.

Acknowledgements.  I thank the rest of the PSCz and 2dFGRS teams for all their dedicated work.