text:
# Secondary antiproton flux induced by cosmic ray interactions with the atmosphere

## I Introduction

The antiproton ($`\overline{p}`$) has a particular status in the spectrum of Cosmic Radiation mainly because of its particular production dynamics and kinematics. The main part of the Cosmic Ray (CR) $`\overline{p}`$spectrum measured in balloon and satellite experiments is well accounted for by assuming it to consist of a secondary flux, originating from the interaction between the nuclear CR flux and the interstellar matter in the galaxy (ISM) DO01 . It is expected to be dominant with respect to components from other possible origins. Such other contributions of primary origin and of major astrophysical interest, have been considered recently. In particular, the $`\overline{p}`$flux induced by annihilation of dark matter constituents ST88 ; BO99a ; BE99 , and by primordial black hole evaporation CA76 ; BA00  have been discussed. All these possible contributions are intimately entangled together and their phenomenological disentagling relies critically on the accuracy of the experimental data. The measurements of the $`\overline{p}`$flux thus provide a sensitive test of the production source and mechanism, and of the propagation conditions in the galaxy DO01 ; GA92 ; CH97 ; SI98 ; BI99 ; UL99 ; MO02 .

CR antiprotons have been experimentally studied for several decades by satellite or balloon borne experiments (see references in BO01 ). Several recent balloon experiments, like BESS BESS97 ; BESS98  and CAPRICE BO01 ; BO97  have collected new data samples whose analysis have provided determinations of the galactic $`\overline{p}`$flux over a kinetic energy range extending from about 0.2 GeV kinetic energy up to about 50 GeV. In these works, the values of the antiproton galactic flux were obtained by subtracting the calculated atmospheric $`\overline{p}`$flux from the values of the measured total flux.

Secondary galactic as well as atmospheric antiprotons are both produced in hadronic collisions by the same elementary reaction mechanism in nucleon-nucleon collisions between the incident CR flux and either ISM nuclei (mainly Hydrogen) in the galaxy, or atmospheric nuclei (mainly Nitrogen) in the atmosphere. The basic $`\overline{p}`$production reaction is the inclusive $`NN\overline{p}X`$ process, $`N`$ standing for nucleon and X for any final quantal hadronic state allowed in the process. The ratio of $`\overline{p}`$production in the galaxy and in the atmosphere scales with the ratio of matter thickness (in units of interaction length $`\lambda _I`$) crossed by protons in the two media. These thicknesses are known to be of the same order of magnitude. In addition, both flux are driven by similar transport equations (see GA92 ; SI98 ; PF96 ; GA99  for example), with however the escape term arising from convective and diffusive effects in the interstellar medium on the galactic flux GA92 , making a significant difference with the transport of flux in atmosphere, which tends to decrease the transported flux compared to the atmospheric transport conditions.

It can be shown using a Leaky Box Model (LBM) for the galactic transport and a simple slab model for the production in the upper atmosphere MA01 , that the ratio of the atmospheric $`\overline{p}`$flux at balloon altitude $`N_{atm}`$ to the galactic $`\overline{p}`$flux at TOA $`N_{gal}`$, is approximately:

$$\frac{N_{atm}(\overline{p})}{N_{gal}(\overline{p})}\frac{x_{Atm}}{\lambda _e}\frac{\sigma (pAtm\overline{p}X)}{\sigma (pp\overline{p}X)}\frac{m_H}{m_{Atm}}\frac{x_{Atm}}{\lambda _e}(\frac{m_H}{m_{Atm}})^{1/3}$$

Where $`x_{Atm}`$ is the thickness of atmosphere on top of balloon experiments, while $`\lambda _e`$ is the LBM escape length, $`\sigma (pp\overline{p}X)`$ and $`\sigma (pAtm\overline{p}X)`$ being the inclusive antiproton production cross sections on hydrogen and on the atmospheric nuclei respectively, while $`m_H`$ and $`m_{Atm}`$ are the hydrogen and the mean atmospheric nuclear mass respectively. Using $`x_{Atm}`$=3.9 g/cm2 (for a 38 km altitude) HE91  and $`\lambda _e`$=8 g/cm2 and 11.8 g/cm2 for particle rigidities of 3 GV and 10 GV respectively MATT , the above ratio is found to be of the order of 0.15 and 0.2, respectively. The contribution of the atmospheric antiproton production to the total flux measured in balloon experiments is thus not expected to be negligible with respect to the galactic component. The correction of the total flux from the atmospheric contribution therefore needs the latter to be calculated very carefully since the accuracy on the evaluation of this component sets a limit of accuracy on the final value of the measured galactic flux.

It must be emphasized that studying also the secondary proton flux in the atmosphere in this context is interesting since the latter is very sensitive to all the components of the simulation process, in particular to the secondary proton production cross section which contributes as well to the generation of the antiproton flux. The comparison of the calculation to the recently measured data provides a robust validation of the approach used and of the overall calculation, and firmly supports the reliability of the results reported here. This study is reported in a separate (companion) paper BA03 .

The present work is an extension of a resarch program aiming at the interpretation of satellite data and which first results on the flux of protons, leptons, and light ions, below the geomagnetic cutoff (GC), at satellite altitude, have been reported recently PAP1 ; PAP2 ; PAP3 .

The paper reports on the calculated $`\overline{p}`$atmospheric flux over the range from sea level up to satellite altitudes by Monte-Carlo simulation. The main features of the calculations are described in section II. The production cross sections used in the event generator are given in section III. The results are discussed in section IV. The work is concluded in section V.

## II Simulation conditions

As mentioned above, the flux of secondary atmospheric antiprotons has been investigated using the same simulation approach which has allowed to successfully account for the $`p`$, $`d`$, $`He`$, and $`e^\pm `$ experimental flux below the Geomagnetic Cutoff (GC) measured by the AMS experiment, as well as the experimental proton and muon flux in the atmosphere, the latter being studied together with the atmospheric neutrino flux LI03 .

The same computing environment has been used here for the charged particle propagation in the terrestrial environment including the atmosphere, as in the previous studies, with the event generator being dedicated however, based on the antiproton production cross section in nucleon-nucleon collisions.

The calculation proceeds by means of a full 3D-simulation program. Incident Cosmic Rays are generated on a virtual sphere chosen at a 2000 km altitude. Random events are generated uniformly on this sphere. The local zenith angle distribution of the particule momentum is proportional to $`\mathrm{cos}(\theta _z)d\mathrm{cos}(\theta _z)`$, $`\theta _z`$ being the zenithal angle of the particle, in order to get an isotropic flux at any point inside the volume of the virtual sphere. The geomagnetic cut-off is applied by back-tracing the particle trajectory in the geomagnetic field, and keeping in the sample only those particles reaching a backtracing distance of 10 Earth radii. Flux conservation along any allowed particle path in the geomagnetic field is ensured by application of Liouville’s theorem. The normal particle propagation as well as its back-tracing are performed using the adaptative Runge-Kutta integration method in the Geomagnetic field LI03 .

1. For the incident CR proton and helium flux, functional forms fitted to the 1998 AMS measurements were used PROT ; HELI  (see also BO99 ; ME00 ; SA00 ). The heavier components of the CR flux were not taken into account in the calculations (see LI03 ). For other periods of the solar cycle, the incident cosmic flux are corrected for the different solar modulation effects using a simple force law approximation SOLMOD .

2. Each particle is propagated in the geomagnetic field and interacts with nuclei of the local atmospheric density according to their total reaction cross section and producing secondary nucleons $`p`$, $`n`$, and antinucleons $`\overline{p}`$, $`\overline{n}`$, with cross sections and multiplicities. This important issue is discussed in section III below. The specific ionization energy loss is computed for each step along the trajectory.

3. In the following step, each secondary particle produced in a primary collision is propagated in the same conditions as incident CRs in the previous step, resulting in a more or less extended cascade of collisions through the atmosphere, which may include up to ten generations of secondaries for protons for the sample generated in this work BA03 ; PAP1 .

For the antinucleon inelastic collisions, only the annihilation reaction channel was taken into account. Non annihilating inelastic $`\overline{N}+A\overline{N}+X`$ ($`\overline{N}`$ standing for antinucleon) interactions whose contribution to the total reaction cross section $`\sigma _R`$, is small. It consists basically of the single diffractive dissociation cross section (for the proton target in individual $`\overline{N}p`$ collisions), and it would be of the same order of magnitude as for pp colisions, namely $``$10% of $`\sigma _R`$ or less at the energies considered here PBARINEL . It has been neglected at this stage. It will be included in the further developments of the calculation program.

The reaction products are counted whenever they cross, upwards or downwards, the virtual detection spheres (several can be defined in the program) at the altitude of the detectors: from sea level up to about 36 km for ground and balloon experiments (BESS, CAPRICE), 370 km for the AMS satellite experiment. Higher altitudes up to more than 10000 km were also investigated, with the purpose of understanding the dynamics of the population of quasi-trapped particles in the earth environment (see section IV.7). Each particle is propagated until it disappears by nuclear collision (annihilation for antinucleons), stopping in the atmosphere by energy loss, or escaping to the outer space beyond twice the generation altitude PAP1 ; PAP2 ; PAP3 ; LI03 .

In the terminology used in the following, one event is defined as the full cascade induced by an incoming CR particle interacting with one atmospheric nucleus. For each CR producing at least one secondary particle, the whole event is stored with all the relevant topological, dynamical, kinematical, and geographical informations. This includes the collision rank, geophysical location, altitude, momentum and particle type, and parent particle type, in form of event files. The collision rank is defined as the number of a given collision in the cascade initiated by the first CR interaction with atmosphere (rank 1).

The calculations do not include any adjustable parameter.

## III Cross sections

### III.1 Proton induced secondaries

#### III.1.1 Protons

The inclusive $`p+Ap+X`$ proton production cross sections used are described in ref PAP1 . They are based on the results of refs BA85  and AB92  for the two components corresponding basically to forward (or direct quasi elastic) and backward (or relaxed deep inelastic) productions respectively. The values obtained have been found in reasonable agreement with the results of the INCL model of intranuclear cascade calculations BO02 . This cross section allows to reproduce very successfully the atmospheric secondary proton flux down to the lowest altitude BA03 .

The cross section for secondary neutron production was taken equal to that of proton production. Similarly, the neutron induced cross sections of secondary nucleon production used was taken the same as for protons, the coulomb interaction making negligible differences over the considered energy range.

#### III.1.2 Antiprotons

The inclusive $`p+A\overline{p}+X`$ antiproton production cross section has been obtained by fitting a set of available experimental data between 12 GeV incident kinetic energy and 24 GeV/c incident momentum SU98 ; AB93 ; AL70 ; EI72 , using a modified version of the analytical formula proposed in ref KMN , the latter being referred to as KMN (for Kalinovski, Mokhov, Nikitin) in the following.

The invariant triple differential cross section is described by means of the formula used in KMN :

$$(E\frac{d^3\sigma }{d^3p})_{inv}=\sigma _RC_1A^{b(p_t)}(1x)^{C_2}\mathrm{exp}(C_3x)\mathrm{\Phi }(p_t)$$
(1)
In this relation the kinematical variables $`p_t`$ and $`x`$ are the transverse momentum and the fractional energy of the particle respectively as defined in KMN  (relation 3.26 in this reference), while $`\sigma _R`$ is the total $`p+A`$ reaction cross section. The function $`\mathrm{\Phi }(p_t)`$ was modified as:

$$\mathrm{\Phi }(p_t)=\mathrm{exp}(C_4p_t^2)+C_5\frac{\mathrm{exp}(C_6x_t)}{(p_t^2+\mu ^2)^4}\mathrm{exp}(\alpha \sqrt{s})$$
(2)
With $`b(p_t)=b_0p_t`$ within the range of 4-momentum transfer considered here. Figure 1 shows some of the results obtained HU01a  by fitting the set of available data studied for this work with relation 1. The energy dependent exponential factor in the second term of relation 2 was introduced for it was found that this term contributes only at low incident energy. The values of the parameters obtained are given in table 1. These values are significantly different from those given in the original work HU01a ; HU02  on a smaller number of data. This work is currently being extended up to 400 GeV incident proton energy DU03 .

### III.2 4He induced secondaries

#### III.2.1 Protons

The inclusive $`{}_{}{}^{4}He+Ap(n)+X`$ proton (neutron) production cross section was obtained as described in section III.1.1 above for the $`p+Ap+X`$ cross section (using the total reaction cross section from JA78  for this system), renormalized to the available experimental multiplicities measured for this reaction AR87 ; BA93 .

#### III.2.2 Antiprotons

The inclusive antiproton production cross section was evaluated by means of the wounded nucleon model BI76 ; FA84  using the experimental values of the total reaction cross sections for the $`{}_{}{}^{4}He+A`$ and $`p+A`$ systems, and of the $`\overline{p}`$production multiplicity in nucleon-nucleon collisions HU02  available.

In this model, the particle production multiplicity $`<n_{AB}>`$ in the collision between ions $`A`$ and $`B`$ is related to the multiplicity nucleon-nucleon $`(𝒩𝒩)`$ collision $`<n_{𝒩𝒩}>`$ by the relation:

$$<n_{AB}>=\frac{1}{2}\left(A\frac{\sigma _{pB}}{\sigma _{AB}}+B\frac{\sigma _{pA}}{\sigma _{AB}}\right)<n_{𝒩𝒩}>$$
(3)
with $`\sigma _{ij}`$ being the total reaction cross section between the $`i`$ and $`j`$ system. Using this model, the $`\overline{p}`$production multiplicity induced by the CR $`He`$ component on the nitrogen $`N`$ component of the atmosphere is found to be $`<n_{HeN}>1.55<n_{pN}>2.5<n_{PP}>`$.

### III.3 Total reaction cross sections

$``$Protons: The values of the proton total reaction cross section on nuclei used were obtained from the parametrization of LE83 , and checked on the carbon data from JA78 .
$``$4He: The 4He+A total reaction cross sections used were taken from JA78 .
$``$Antinucleons: The $`\overline{p}`$total reaction cross section was taken from CA79 , with the energy dependence from the data compilation of BA88 . The same production cross section and total reaction cross section have been assumed for $`\overline{n}`$production as for $`\overline{p}`$.

## IV RESULTS

A sample of about 35 106 CRs have been simulated, of which 20% were effectively propagated to the atmosphere (above GC), for detection altitudes going from sea level up to 104 km altitude, including the BESS/CAPRICE balloon altitude ($``$38 km), the AMS orbit altitude ($``$370 km), and the recent BESS measurement terrestrial altitude (2770 m). The flux at sea level was calculated to investigate the possibility of ground level measurement of atmospheric antiprotons with existing devices HU01b . This was achieved independently by BESS at mountain altitude and the results are discussed below.

### IV.1 Particle trajectories in the Earth magnetic field

The time of confinement of particles in the earth environment together with their particular trajectories, determine their status with respect to the three categories: trapped, semi-trapped, and non trapped (escape). Trapped particles are spiralling back and forth around and along the magnetic field lines long enough to drift many times around the earth (see for example SI62 ; ST55  and below). Trapped particles are practically not observed in the energy domain considered here. They are not dynamically forbiden however and a few trajectories with a few 102 bounces are observed, which corresponds to short-lived trapped particles. Quasi trapped particles are in similar kinematic conditions but accomplish only a limited number of bounces at mirror points before being absorbed or before escaping (see examples below). This concept appeared during the first years of radiation belts studies GA60  (see also the discussion in SI62 ). Escape particles do not match the kinematic conditions for being trapped at their production point and escape in a very short time to the outer space. All intermediate situations between the stereotypes of quasi-trapped and escape trajectories are in fact observed in the simulation results (see example in bottom left fig 2).

Figure 2 shows four examples of characteristic trajectories of antiprotons generated in this study. Each of the 4 panels gives a side view (projection on the meridian plane, top left), side view zoomed around the production point showing the spiralling trajectory of the particle (top right), top view (projection on the equatorial plane, bottom left), and 3-D representation (bottom right) of each of the selected trajectories.

The top left event (1.52 GeV kinetic energy) is an escape particle produced close to the North pole. Top right is a semi trapped single bounce event (0.54 GeV) annihilating in the atmosphere close to its production point. Bottom left is a longer lifetime, multi bounces, semi trapped event (2.48 GeV), drifting around the earth about three quarters of a turn before annihilation in the atmosphere in the SAA region. Bottom right is intermediate between semi-trapped (since it displays at least one clear bounce) and escape event (0.54 GeV). It is a type of event for which the first adiabatic invariant (magnetic moment conservation) is not conserved because of a large variation of the magnetic field along the radius of gyration SI62 .

### IV.2 General features of the simulated data

Figure 3 shows a few basic distributions of physics observables relevant to the dynamics of the process for two detection altitudes: 38 km (solid line) and 380 km (dashed line) corresponding to balloon and satellite (AMS) altitudes, inside and outside the atmosphere, respectively. The rank of the collision producing the antiproton (top left panel on the figure) appears to extend from 1 up to about 10 for the simulated sample, showing that $`\overline{p}`$s are produced up to the tenth generation in the collisions cascade. The distributions are a little different for the two altitudes, with a significantly larger number of $`\overline{p}`$occuring from first interaction at the lower altitude. The altitude distribution of the production point for the detection at 38 km (top right) shows a discontinuity at this altitude due to the incoming flux dominated by production from the upper layer of atmosphere. The mean production altitude is found around 46 km and 48 km for the lower and upper detection altitudes respectively. The particle momentum spectrum at the production point (bottom left) is found harder at the higher altitude. The number of bounces at the mirror points for particle trajectories spiralling around the magnetic field lines are found as expected very different for the two altitudes of detection: At 38 km, only a small population is seen to reach a number of bounces larger than a few units (5-6). This flux is significant however and must correspond to trajectories lying mostly outside of the atmosphere. At 380 km, the observed flux of $`\overline{p}`$ trajectories with more than one bounce is larger by about 2 orders of magnitude than at 38 km, corresponding to the population of quasi-trapped particles as discussed previously in PAP1 ; BU01 .

### IV.3 Antiproton flux at mountain altitude

The recent measurements of the $`\overline{p}`$flux at 2770 m of altitude by the BESS collaboration FU01  allows a sensitive test of the ability of the present simulation program to account for the observed flux since at this altitude the $`\overline{p}`$production occurs on the average after a casade of 4 collisions on the average (see previous section and fig 3). Another highly sensitive test of the overall calculation is provided in BA03  on the atmospheric proton flux.

Figure 5 shows the $`\overline{p}`$spectrum at 2770 m of altitude measured by BESS, compared to the simulation results. The latter has been run with the geometrical acceptance function of the BESS spectrometer given in FU01  (figures 4.37 to 4.39, see also AJ00 ), the overall acceptance angle being of the order of 25. The total CR $`p+He`$ and partial $`He`$ flux are shown on the same figure. The $`He`$ flux contribution is seen to produce a small fraction of about 5% of the full $`\overline{p}`$flux. Although the total flux calculated somewhat underestimates the experimental values, the overall agreement is good, the calculated values being on the average within one standard deviation from the experimental values. This gives confidence in the results of the calculations obtained for the other altitudes investigated and reported below.

Figure 5 compares the experimental zenith angle distributions of the $`\overline{p}`$flux to the calculated values (histograms) for the same kinetic energy bins as measured by BESS FU01 . On this figure, the overall agreement between data and calculations again appears to be good for all energy bins.

Note that no upward particle have been produced at this altitude in the simulated sample, as it could be expected HU02 .

### IV.4 Balloon data

In this section the atmospheric $`\overline{p}`$flux at ballon altitude is investigated for comparison with the atmospheric $`\overline{p}`$corrections made to the raw flux data in the BESS and CAPRICE experiments.

Figure 7 shows the values of the galactic $`\overline{p}`$flux obtained from the BESS and CAPRICE measurements respectively. These values have been obtained from the measured raw flux by subtraction of the atmospheric $`\overline{p}`$flux evaluated using an average of theoretical calculations for the BESS experiment BESS98 , and using the calculations of ref ST97  for CAPRICE. On the figure, the atmospheric flux calculated in ST97  is compared with the results from the present work (see also PF96 ). For the two sets of data, it appears that the the present calculations are in fairly good agreement with the atmospheric antiproton flux obtained from transport equation calculations and used to correct the raw flux measured. There is a slight trend however for the simulation results to be larger than those obtained from the differential equation approach by about 20% over the range 10-30 GeV. At low $`\overline{p}`$energies the opposite trend is observed and the simulation results are found significantly below the values obtained from the differential equation. One might say that the simulation results should be taken with care below 1 GeV because of the lack of experimental cross sections for low energy $`\overline{p}`$production, and thus of the large corresponding uncertainties on the results of the simulation over this range, however, it must be noted that the calculated cross sections for low $`\overline{p}`$momentum should be in principle reliable for the following reason. The $`\overline{p}`$distribution is naturally symmetric in the rapidity space. The fitting function has the same property since it depends only on variables matching this symmetry. Therefore a good fit to a set of experimental cross sections for $`\overline{p}`$rapidities above the center of mass rapidity $`Y_{cm}`$, automatically ensures the right behaviour of the calculated values for rapidities below $`Y_{cm}`$, i.e., for particle momenta in the laboratory frame, because of the symmetry law.

To conclude this section, the atmospheric $`\overline{p}`$flux calculated in this work confirm the corrections of the raw flux values measured in the BESS and CAPRICE experiments. This result updates and corrects a previous preliminary conclusion on the issue HU01b  recently quoted in MO03 .

The contributions of A$`>`$4 CR components were not included in the calculations, neither were those from the non-annihilating inelastic contributions in the $`\overline{p}`$propagation through the atmosphere. These contributions are small however and not likely to change the results by more than a few percents LI03 .

### IV.5 AMS altitude

Future satellite experiments in preparation, plan to measure the $`\overline{p}`$flux. A reliable knowledge of the atmospheric $`\overline{p}`$flux at satellite altitudes is therefore highly desirable for these experiments to step on explored grounds. The $`\overline{p}`$flux calculated for the altitude of the AMS orbit are presented in this section.

Figure 7 shows the expected downwards (Secondaries and reentrant Albedo, dashed histogram) and upwards (Splash Albedo particles, solid histogram) flux of atmospheric antiprotons at the altitude of AMS for two regions of lower geomagnetic latitude: equatorial ($`0<|\theta _M|<0.2`$ rad), intermediate ($`0.2<|\theta _M|<0.4`$ rad), and subpolar ($`0.8<|\theta _M|<0.9`$ rad). As expected, the flux is predicted larger for the lower latitudes than it is around the poles, because of the existence of quasi-trapped $`\overline{p}`$components at the low and intermediate latitudes. Note that the simulated flux is surprisingly predicted larger downwards than upwards (bottom left panel). This is in fact an effect of the spectrometer acceptance (taken to be 30 with respect to zenith), the mean angle for upward particle trajectories being 2 radians. The overall upward flux is larger than downward by a factor of about 2.5 HU01b ; HU02 . This shows that the future satellite measurements of antiproton flux at low latitudes will have to be corrected from the atmospheric contributions and will probably suffer more uncertainties than previously thought.

The lower right panel compares the $`\overline{p}`$data at TOA reported by AMS AMS1T , BESS BESS98 , and CAPRICE BO01 , to the flux calculated in the polar region where the AMS data have been measured. The calculated (downward) atmospheric $`\overline{p}`$component is at the percent level of the measured flux in the low energy $`\overline{p}`$range, and can be considered as negligible at all energies.

### IV.6 $`\overline{p}`$and $`\overline{n}`$flux at terrestrial altitudes

The $`\overline{p}`$and $`\overline{n}`$flux have been calculated also at sea level in order to provide a realistic order of magnitude of these flux for general purpose and for ground testing of embarked experiments.

$``$Antiprotons :

The flux of atmospheric antiprotons at sea level has been calculated with the same simulation program. Figure 8 shows the distributions obtained at sea level (left) and at 4000 m. The energy integrated flux is of the order of 0.4 10$`^3`$ $`\overline{p}`$$`s^1m^2sr^1`$ at all latitudes (see fig 9 below). At 4000 m (right panel on the figure), the flux raises to about 7 10$`^3`$ $`\overline{p}`$$`s^1m^2sr^1`$. These values are small but large enough for this flux to be measured by currently existing large acceptance detectors (BESS, CAPRICE), or in a near future by new detectors under construction like AMS and PAMELA.

$``$Antineutrons :

Atmospheric secondary antineutrons may also be of interest in ground or balloon measurements ME01 . Figure 8 shows the kinetic energy spectrum of the expected $`\overline{n}`$flux at sea level and at 4000 m (this latter altitude being that the Cerro-La-Negra observatory in Mexico where some experimental measurements of the antineutron flux are being considered CERRO ).

### IV.7 Flux dependence on the altitude

The $`\overline{p}`$flux has been calculated up to altitudes of 2 104 km with the aim of investigating the general features of the dynamics and kinematics of the particles in the more remote Earth environment than considered in the previous sections.

Figure 9 shows the altitude dependence of the energy integrated upward and downward $`\overline{p}`$flux in bins of latitude, assuming a geometrical acceptance of 30 for the detector. The calculated distributions display two main features:
1) In the atmospheric range of altitudes, a large peak of incoming flux centered around 20 km and corresponding to atmospheric secondaries, dominates the distribution ($``$ 50km, i.e., $``$TOA), with basically no associated outgoing flux (inside the quoted acceptance angle).
2) Above the atmosphere, surprisingly, the calculated upward and downward flux are found close to each other up to fairly high altitudes, namely $``$104 km, for the low and intermediate latitudes ($`|\theta _M|0.7`$ rad). This shows that a population of quasi-trapped, particles should be observed in this region of space, i.e., up to around five to ten thousands kilometers. This is confirmed by the lifetime of the particle between their production and their absorption, and by the number of bounces of the particles between the mirror points of their trajectories, which extend up to 100 seconds and several 102 bounces (see Fig. 3), respectively for the simulation sample produced. At higher latitudes with ($`|\theta _M|0.7`$ rad), the incoming flux progressively disappears, and the outgoing flux then corresponds to escape particles.

From these calculations it can be concluded that there should exist a significant flux of quasi-trapped particles extending approximately over a decade of altitudes, from about 50 km (TOA) up to $``$104 km, depending on the particle energy and latitude. This flux has been observed already by the AMS experiment in the lower part of the altitude range (380 km) PAP1 . Note that the issue was discussed long ago in a pioneering paper about the electron flux GR77 . See also BA03  (companion paper) for a complementary discussion. The energy spectrum of this flux extends up to around 10 GeV, which is about the upper momentum limit (8.5 GeV) for which particles can match the simple geometrical condition that the gyration radius is smaller than the mean trajectory radius to the upper atmosphere (for equatorial latitude trajectories, at the limit of large, close to $`\pi `$/2, pitch angles).

Future embarked experiments should take these features into account, even though, in principle, the accurate knowledge of the kinematical conditions of the particle at the detection point allows to know whether it is or not of atmospheric origin.

## V Summary and conclusion

In summary, the secondary antiproton flux produced by the Cosmic Ray proton and Helium flux on the atmosphere has been calculated by Monte-Carlo simulation. The flux calculated for the altitude of 2770 m is in fair agreement with the recent BESS measurements. At sea level, it is small but measurable and could provide a natural facility for testing the identification capability of existing experiments or of future devices. For balloon altitudes, the calculated flux has been found in agreement with the values calculated in previous works. At satellite altitudes (380 km) it appears to be negligible compared to the CR flux for polar latitudes, and of the same order of magnitude as for the high balloon altitudes for equatorial and intermediate latitudes (below the geomagnetic cutoff), indicating that it will have to be taken into account in future measurements of the galactic antiproton flux at similar altitudes.

Acknowledgements: The authors are indebted to V. Mikhailov for helpful discussions on the issue and for pointing ref GR77  to them. They are also grateful to M. Nozaki, M. Fujikawa, and the BESS collaboration, for providing their measurements of $`\overline{p}`$flux at 2770 m prior to publication.


format: text

text:
# Comment on “Conformally flat stationary axisymmetric metrics”

## Abstract

In , the authors claim to have found a previously overlooked family of stationary and axisymmetric conformally flat spacetimes, contradicting an old theorem of Collinson . In both  and  it is tacitly assumed that the isometry group is orthogonally transitive. Under the same assumption, we point out here that Collinson’s result still holds if one demands the existence of an axis of symmetry on which the axial Killing vector vanishes. On the other hand if the assumption of orthogonal transitivity is dropped, a wider class of metrics is allowed and it is possible to find explicit counterexamples to Collinson’s result.

PACS Numbers: 04.20.Jb , 02.40.Hw

In  the general form of conformally flat spacetimes admitting an Abelian 2-parameter orthogonally transitive group $`G_2`$ of isometries acting on timelike orbits is found. In that paper the assumption of orthogonal transitivity is not mentioned explicitly, but from the form of the line element given in equation (1) which García and Campuzano take as their starting point, it is clear that orthogonal transitivity has been assumed. The authors then prove that not all those spacetimes can be diagonalized or, in other words, that the there are line-elements with no hypersurface-orthogonal timelike Killing vector field. Thus, these cases are properly stationary (non-static) solutions. This result seems to contradict an important theorem due to Collinson , see also .

However, Collinson’s result still holds in an appropriate sense: if we require that the spacetime be truly axially symmetric — meaning that it contains a non-empty axis of symmetry —, and not merely cyclically symmetric , then none of the new stationary metrics survive. Of course, every spacetime with a $`G_2`$ acting on timelike orbits can be made cyclically symmetric by simply closing the orbits of a spacelike Killing vector, that is, identifying two values of the appropriate coordinate. Nevertheless, this does not mean that the spacetime contains an axis of symmetry (the 2-dimensional cylinder is cyclically symmetric, but there is no axis: the axis is “outside” the manifold). For a metric to have an axis it is necessary that the cyclic Killing vector vanishes on it, see e.g.  and references therein. As a matter of fact, this is independent of whether or not the axis is regular . If the axis is regular, then the elementary flatness condition must also be satisfied .

The cyclic Killing vector of the new stationary solutions found in , formula (51), is given by $`_\varphi `$. A necessary condition so that $`_\varphi `$ vanishes on a would-be axis is that its scalar products with every other vector vanish there. This would require that — using the notation of  — either $`e^{G(x,y)}=0`$ somewhere, which is clearly not possible as the whole metric would be zero, or that $`x=0`$ and $`1+x^2=0`$, which is also impossible. Thus, none of the stationary solutions found in  is axially symmetric, neither with a regular nor with a singular axis.

The conclusion is that Collinson theorem remains valid in the orthogonally transitive case if one requires that the spacetime contains a non-empty axis. Combining the results in  and , we have

###### Theorem 1

Every conformally flat stationary cyclically symmetric spacetime in which the isometry group is orthogonally transitive, is given by the solutions in . If in addition axial symmetry is required, then the spacetime is necessarily static.

If the assumption of orthogonal transitivity is dropped then, as pointed some years ago , Collinson’s result is not valid as the following counterexample shows. The Stephani metric  with zero fluid expansion is:

$$ds^2=(1+\frac{1}{4}Kr^2)^2(dx^2+dy^2+dz^2V^2dt^2)$$

where $`K=0,\pm 1`$ and $`V`$ is given by

$$V=a(t)+b(t)r^2+𝐫𝐜(t)$$

with $`a`$, $`b`$, $`c_1`$, $`c_2`$ and $`c_3`$ arbitrary functions of $`t`$. This metric is a conformally flat perfect fluid solution (with fluid velocity proportional to $`_t`$). The fluid energy density and pressure are given by

$$\mu =3Kp=3K+(aK+4b)(1+\frac{1}{4}Kr^2)/V$$

so that the density is constant, but the pressure normally depends on all four coordinates. In general the Stephani metric admits no isometries, however, as shown in , if we take $`K=+1`$, $`2a=12A\mathrm{sin}(Bt)`$, $`8b=1+2A\mathrm{sin}(Bt)`$, $`c_1=c_2=0`$ and $`c_3=A\mathrm{cos}(Bt)`$ where $`A`$ and $`B`$ are non-zero constants, the metric admits a (complete) two-dimensional isometry group with Killing vectors

$$\stackrel{}{X}=_t+B\left[\frac{1}{2}Kxz_x+\frac{1}{2}Kyz_y+(1\frac{1}{4}Kr^2+\frac{1}{2}Kz^2)_z\right]$$

$$\stackrel{}{Y}=x_yy_x$$

This has an axis $`x=y=0`$ where the spacelike Killing vector $`\stackrel{}{Y}`$ vanishes and the ‘tilted’ Killing vector $`\stackrel{}{X}`$ is timelike at least in a region around the origin $`x=y=z=0`$ if $`B`$ is sufficiently small. The vector $`\stackrel{}{X}`$ is not hypersurface orthogonal as may be checked by direct calculation or by noting that if it were hypersurface orthogonal it would be a Ricci eigenvector  which it clearly is not (provided that $`\mu +p0`$). Furthermore note that the fluid velocity vector does not lie in the two-plane spanned by $`\stackrel{}{X}`$ and $`\stackrel{}{Y}`$ so that the fluid flow is convective; this is necessarily the case since, otherwise, the isometry group would be orthogonally transitive e.g. . Note also that there are similar examples when $`K=1,0`$; see  for full details.

The properly axially symmetric spacetimes include the physically relevant cases in which a fluid is to be matched to an exterior spacetime across a spatially compact boundary, such as in the case of rotating stars. Observe also that the following classic statement  is only true if the proviso of orthogonal transitivity is added: The only conformally flat stationary and axially symmetric non-convective perfect fluid spacetime is the spherically symmetric interior Schwarzschild solution with constant density.


text:
# Particle Dispersion on Rapidly Folding Random Hetero-Polymers

## Abstract

We investigate the dynamics of a particle moving randomly along a disordered hetero-polymer subjected to rapid conformational changes which induce superdiffusive motion in chemical coordinates. We study the antagonistic interplay between the enhanced diffusion and the quenched disorder. The dispersion speed exhibits universal behavior independent of the folding statistics. On the other hand it is strongly affected by the structure of the disordered potential. The results may serve as a reference point for a number of translocation phenomena observed in biological cells, such as protein dynamics on DNA strands.

The study of random motion on complex structures is essential to the understanding of dispersion phenomena observed in numerous physical systems, ranging from epidemics spreading in complex networks and information transport in modern communication networks such as the internet Albert and Barabási (2002); Volchenkov et al. (2002). In biological cells, the transport of macromolecules is accomplished by a variety of translocation processes in which carrier molecules move along complex fibrous polymer networks, e.g. myosin translocation on actin fibers Amblard et al. (1996) or transport on microtubules Caspi et al. (2000). If the involved topologies are scale-free, diffusion is often anomalous, i.e. the mean square displacement of a particle violates the linear dependence on time $`X^2(t)t^\gamma `$ with $`0<\gamma 1`$ Bouchaud and Georges (1990). Depending on the underlying microscopic dynamics, subdiffusive ($`\gamma <1`$) as well as superdiffusive ($`\gamma >1`$) behavior is observed. For instance, when a particle moves along a polymer in a complex folding state, it may jump to a neighboring location in Euclidean space which is distant in chemical coordinates. Effectively, the particle moves superdiffusively along the chain Sokolov et al. (1997); Manna et al. (1989) and performs a random walk known as a Lévy flight. This mechanism may explain fast target localization of regulatory proteins moving along DNA strands Berg et al. (1981). Lévy flights have been observed in a variety of systems, ranging from chaotic systems Geisel et al. (1985) and particle dispersion in turbulent flows Porta et al. (2001) to foraging animals Viswanathan et al. (1996); Levandowsky et al. (1997) and climate changes Ditlevsen (1999). Lévy flights are characterized by an exponent $`0<\mu <2`$ which quantifies the degree of superdiffusion and is related to the heuristic dispersion relation $`X(t)t^{1/\mu }`$. When Lévy flights successfully mimic single trajectories, the associated stochastic evolution equations are no longer of the Fokker-Planck type but rather generalizations thereof which involve fractional differential operators. Fractional models have contributed considerably to the understanding of these systems, in fact the terms fractional kinetics and fractional dynamics have been coined to classify them Fogedby (1994); Hilfer (2000); Metzler and Klafter (2000); Brockmann and Sokolov (2002); Brockmann and Geisel (2003); Zaslavsky (2002). Of particular interest are systems in which the cause for superdiffusive dispersion and the heterogeneity of the environment interact antagonistically.

In this Letter we introduce and investigate a model for superdiffusive particle dispersion on flexibly folding random hetero-polymers. We focus on the interplay between long range Lévy type transitions due to folding and the quenched random disorder caused by the heterogeneity of monomers of the chain. Based on simple assumptions on the hopping rate and configurational dynamics, we derive a fractional Fokker-Planck equation (FFPE) describing the motion of the particle along the polymer. We compute the relaxation properties as a function of the effective potential strength and the Lévy exponent $`\mu `$. We find that the dispersion speed depends considerably on $`\mu `$, but becomes universal on larger spatial scales apart from a discontinuous change at $`\mu =2`$ (i.e. for ordinary diffusion). Furthermore, the relative concentration of monomers and thus the particular shape of the potential does not affect the ordinary diffusion process ($`\mu =2`$), but strongly affects all superdiffusive processes, a result we believe to be crucial for the understanding of transport phenomena in living cells.

Consider the scenario depicted in Fig. 1. A particle is attached to a hetero-polymer and performs a random walk along the chain. Let $`x`$ denote the chemical coordinate with an inter-monomer spacing of $`a`$. The chain is flexible, and rapidly changing its conformational state defined by the Euclidean coordinate $`𝐑(x)`$ of each monomer. The heterogeneity of the chain is modeled by a potential $`V(x)`$ which specifies the probability of the particle being attached to site $`x`$. In a thermally equilibrated system this probability is proportional to the Boltzmann factor $`\mathrm{exp}\left[\beta V(x)\right]`$. The dynamics of the particle is governed by the rate $`w(x|y,t)`$ of making a transition $`yx`$ at time $`t`$. We assume that transitions occur only between monomer sites which are close in Euclidean space, i.e. when $`|𝐑(y,t)𝐑(x,t)|a`$. We make the simplest possible ansatz for this rate to take into account the requirements of Gibbs-Boltzmann statistics, the potential heterogeneity of the chain and the complexity of conformational states,

$$w(x|y,t)=\frac{1}{\tau }e^{\beta [V(x)V(y)]/2}\mathrm{\Gamma }(x,y;t).$$
(1)
where the parameter $`\tau `$ is the typical microscopic time constant and the function $`\mathrm{\Gamma }(x,y;t)`$ is defined by

$$\mathrm{\Gamma }(x,y;t)=\{\begin{array}{cc}1\hfill & \text{if}\left|𝐑(x,t)𝐑(y,t)\right|a,\hfill \\ 0\hfill & \text{otherwise}.\hfill \end{array}$$
(2)
$`\mathrm{\Gamma }(x,y;t)`$ reflects the dependence of transitions on the time dependent conformational state of the chain and is symmetric, i.e. $`\mathrm{\Gamma }(x,y;t)=\mathrm{\Gamma }(y,x;t)`$. The function $`\mathrm{\Gamma }(x,y;t)`$ can be interpreted as a time dependent connectivity matrix (Fig. 1B). The propagator $`p(x,t)`$ of a particle initially ($`t=0`$) at the origin evolves according to the master equation,

$$_tp(x,t)=\text{d}y\left[w(x|y,t)p(y,t)w(y|x,t)p(x,t)\right]$$
(3)
in which the rate is given by Eq. (1). The geometrical factor $`\mathrm{\Gamma }(x,y;t)`$ varies erratically and can be regarded as a stochastic process evolving on a time scale $`\tau _g`$, which is generally different from the hopping time constant $`\tau `$. Averaging (denoted by $`[]`$) over conformational states the dynamics reads $`_t\left[p\right]=\left[p\right]`$ where the operator $``$ is defined by the rhs of Eq. (3). If conformational changes occur on smaller time scales than the hopping ($`\tau _g\tau `$) we may substitute$`\left[p\right]\left[\right]\left[p\right]`$, which represents a mean field approximation. In mean field, Eq. (1) is given by

$$\left[w(x|y,t)\right]=\frac{1}{\tau }e^{\beta \left[V(x)V(y)\right]/2}\rho (|xy|),$$
(4)
where $`\rho (|xy|)=\left[|𝐑(x,t)𝐑(y,t)|a\right]`$ is the probability that two given sites $`x`$ and $`y`$ are neighbors in Euclidean space. If the folding process is stationary, this probability is time independent, and due to translation invariance along the chain it is a decreasing function of distance in chemical space. The specific functional form of $`\rho (x)`$ determines the asymptotics of Eq. (3). Consider the situation depicted in Fig. 1C., where the chain is knotted such that non-local transitions occur on a typical scale $`\sigma >a`$. On larger scales $`\rho (x)`$ vanishes. In this case, a Kramers-Moyal expansion of the rhs of Eq. (3) yields the FPE $`_tp=V^{}p+D\mathrm{\Delta }p`$, in which the diffusion coefficient is given by $`D(\sigma /a)^2/\tau `$ and the gradient force is determined by the potential $`V(x)`$ along the chain. The situation changes drastically for the type of chain sketched in Fig. 1D. For a freely flexible chain the quantity $`\rho (x)`$ follows an inverse power law with increasing chemical distance, i.e. $`\rho (x)1/|x|^{1+\mu }`$. Typically $`\mu <2`$ de Gennes (1979) and thus $`\rho (x)`$ lacks a well defined variance and consequently a typical scale in long range transitions. A particle moving along such a chain will behave superdiffusively and perform a Lévy flight in chemical coordinates. Inserting $`\rho (x)1/|x|^{1+\mu }`$ with $`0<\mu <2`$ into Eq. (4) and subsequently into Eq. (3) the asymptotics is governed by a fractional Fokker-Planck equation (FFPE),

$$D^1_tp=e^{\beta V/2}\mathrm{\Delta }^{\mu /2}e^{\beta V/2}ppe^{\beta V/2}\mathrm{\Delta }^{\mu /2}e^{\beta V/2}.$$
(5)
A detailed derivation is given in Ref. Brockmann and Sokolov (2002). Here, $`D`$ is the generalized diffusion coefficient and the operator $`\mathrm{\Delta }^{\mu /2}`$ is a generalization of the ordinary Laplacian,

$$\left(\mathrm{\Delta }^{\mu /2}f\right)(x)=C_\mu \text{d}y\frac{f(y)f(x)}{|xy|^{1+\mu }}$$
(6)
with $`C_\mu =\pi ^1\mathrm{\Gamma }(1+\mu )\mathrm{sin}(\pi \mu /2)`$. In contrast to the ordinary Laplacian, $`\mathrm{\Delta }^{\mu /2}`$ is a non-local, singular integral operator, reflecting the superdiffusive behavior of the process. The boundary case $`\mu =2`$ represents the limit of ordinary diffusion, i.e. Eq. (5) reduces to an ordinary FPE. When the potential vanishes, $`V0`$, Eq. (5) becomes $`_tp=D\mathrm{\Delta }^{\mu /2}p`$ and is solved by the propagator of the symmetric Lévy stable process of index $`\mu `$, i.e. $`p(x,t)=(Dt)^{1/\mu }L_\mu \left(x/(Dt)^{1/\mu }\right)`$ with $`L_\mu (z)=(2\pi )^1\text{d}k\mathrm{exp}(ikz|k|^\mu )`$ Brockmann and Sokolov (2002). In its general form Eq. (5) describes the dynamics of Lévy flights in external potentials obeying Gibbs-Boltzmann thermodynamics.

In the following we investigate the relaxation properties of Eq. (5) in random potentials $`V`$. Since the shape of the potential is determined by the ordering of different types of monomers along the chain, $`V(x)`$ will be bounded and fluctuate about some average. Furthermore, it will generally possess a typical correlation length $`\xi `$. Without loss of generality we let $`V=0`$ and $`V^2=V_0^2`$. The correlation length $`\xi `$ is defined by $`\xi =V_0^2_0^{\mathrm{}}\text{d}yC(y)`$ where $`C(y)=V(y)V(0)`$ is the correlation function. The most straightforward way to incorporate these attributes into a model is by using Gaussian random phase potentials,

$$V(x)=\frac{1}{2\pi }\text{d}k\varphi (k)e^{i(kx+\vartheta (k))},$$
(7)
which are defined by a set of random uncorrelated phases $`\vartheta (k)`$ (with $`\vartheta (k)=\vartheta (k)`$) and the power spectrum $`S(k)=\text{d}xe^{ikx}C(x)`$ (with $`\varphi (k)\varphi ^{}(k^{})=2\pi S(k)\delta (kk^{})`$). The pdf associated with this choice of $`V(x)`$ is Gaussian with zero mean and variance $`V_0^2=(2\pi )^1\text{d}kS(k)`$. Fig. 2A. and 2B. show two realizations of random phase potentials, each one with a different power-spectrum (and correlation function).

The relaxation properties are determined by the eigenvalue spectrum $`E(k)`$ of the evolution operator $``$ defined by the rhs of Eq. (5). In order to compute the spectrum, the FFPE can be transformed by means of $`p(x,t)=e^{\beta V(x)/2}\psi (x,t)`$. This yields a fractional Schrödinger equation with identical spectral properties,

$`_t\psi `$ $`=`$ $`\psi \text{with}`$ (8)
$``$ $`=`$ $`D\left(\mathrm{\Delta }^{\mu /2}+U\right).`$ (9)
The operator $``$ is symmetric and the effective potential $`U`$ is related to the original potential $`V`$ by $`U(x)=e^{ϵv}\mathrm{\Delta }^{\mu /2}e^{ϵv}`$, where $`v(x)=V/V_0`$ is a rescaled potential of unit variance and $`ϵ=\beta V_0/2=V_0/2k_BT`$ is the potential strength in units of $`k_BT`$.

For vanishing potential $`ϵ=0`$, we have $`U0`$ and $`_0=D\mathrm{\Delta }^{\mu /2}`$ which describes free superdiffusion when $`\mu <2`$. The spectrum of $`_0`$ is given by $`E_0(k)=Dk^\mu `$. The wave number $`k>0`$ defines the spatial scale of the corresponding mode. When a potential is present, the spectrum can be written as $`E(k)=D_\mu (k;ϵ)k^\mu `$ where $`D_\mu (k;ϵ)`$ quantifies the relaxation properties on scales $`k^1`$ with the unperturbed $`k^\mu `$-behavior as a reference. If $`D_\mu (k;ϵ)/D<1`$ the process relaxes more slowly compared to free superdiffusion. The spectrum $`E(k)`$ can be obtained for weak potentials by perturbation theory. If $`ϵ1`$, the effective potential $`U`$ can be treated as a small perturbation, for $`U=𝒪(ϵ)`$. Up to second order in $`ϵ`$ the quantity $`D_\mu (k;ϵ)`$ reads

$$D_\mu (k;ϵ)/D=\left(14ϵ^2G_\mu (k)\right),$$
(10)
where the effect on relaxation is provided by the function

$`G_\mu (k)`$ $`={\displaystyle \frac{1}{8\pi }}{\displaystyle \text{d}qS(q)g_\mu (k/q)}\text{with}`$ (11)
$`g_\mu (z)`$ $`={\displaystyle \frac{1}{z^\mu }}\left({\displaystyle \frac{1}{(1+z)^\mu z^\mu }}+{\displaystyle \frac{1}{|1z|^\mu z^\mu }}2\right).`$ (12)
Figs. 3A and 3B depict $`G_\mu (k)`$ as a function of $`k`$ in units of the inverse correlation length $`\xi ^1`$ for the two types of random phase potentials defined in Fig. 2A and 2B. The solid line depicts the limiting case of ordinary diffusion ($`\mu =2`$). The potential slows down the ordinary diffusion process ($`G_\mu (k)>0`$) on scales larger than the correlation length, and speeds it up ($`G_\mu (k)<0`$) on scales smaller than $`\xi `$. The function $`G_\mu (k)`$ has a pronounced minimum at $`k\xi ^1`$. Moderately superdiffusive processes ($`\mu 1`$) behave in a similar fashion, exhibiting the highest variation for $`k\xi ^1`$. On the other hand, $`G_\mu (k)`$ differs strongly for different $`\mu `$ in the asymptotic regime $`k\xi ^1`$. Note also that on small scales ($`k>\xi ^1`$) almost all processes relax faster than without the potential. In fact, $`G_\mu (k)<0`$ for $`2\mu \mu _c`$ where $`\mu _c=2\mathrm{ln}3/\mathrm{ln}20.415`$. Surprisingly, this is no longer valid for strongly superdiffusive processes with $`\mu <\mu _c`$. For instance, in the case $`\mu =0.2`$ (dashed line in Figs. 3A and 3B) $`G_\mu (k)`$ is positive for $`k>\xi ^1`$, implying that strongly superdiffusive processes are slowed down even on small scales. Comparing potential types, we see that the relaxation is different for each potential, but these differences become less important in the asymptotic regime, which is governed by $`G_\mu (k)`$ as $`k0`$. This limit can be computed from Eqs. (11, 12), observing that $`_0^{\mathrm{}}\text{d}qS(q)=\pi `$ and $`g_2(z0)=1`$ and $`g_{\mu <2}(z0)=1/4`$,

$$\underset{k0}{lim}G_\mu (k)=\{\begin{array}{cc}1/4\hfill & \mu <2\hfill \\ 1\hfill & \mu =2,\hfill \end{array}$$
(13)
Hence, the asymptotic behavior is universal, with the exception of the limiting case of ordinary diffusion, and is independent of properties of the potential. The range of validity of the limit (13), however, strongly depends on $`\mu `$. The limit is not attained for marginal exponents (e.g. $`\mu =0.2`$ and $`1.9`$) even on scales several orders of magnitude larger than the correlation length (Fig. 3A and 3B). The above results are valid for small potential strengths $`ϵ`$. For higher effective potential strengths we investigate the asymptotics numerically. The quantity of interest is the normalized generalized diffusion coefficient $`D_\mu (ϵ)`$ defined by

$$D_\mu (ϵ)=\underset{k0}{lim}D_\mu (k;ϵ)/D.$$
(14)
In the perturbative regime Eqs. (10) and (14) yield the universal relation $`D_\mu (ϵ)=1ϵ^2`$ for $`\mu <2`$ and $`D_\mu (ϵ)=14ϵ^2`$ for ordinary diffusion. Figs. 3C and 3D compare these results to those obtained numerically. Although the numerics deviates from perturbation theoretic predictions for greater potential strengths $`ϵ`$, the universality still holds, i.e. the asymptotics ($`k0`$) is independent of $`\mu `$ and of the statistical properties of the potential. The crucial property is the non-locality of the process (i.e. $`\mu =2`$ vs. $`\mu 2`$). Thus, as soon as the folding properties of the chain permit scale free transitions ($`\mu 2`$), the behavior of $`D_\mu (ϵ)`$ changes abruptly.

The pdf of random phase potentials is symmetric with respect to the mean, i.e. a value $`V`$ is as likely to occur as $`V`$ along the chain. For a number of hetero-polymers this assumption is inadequate. Consider the simple model copolymer depicted in Fig. 2C. The chain consists of a random arrangement of monomers, each one equipped with an intrinsic local potential parity $`v_{}`$ and $`v_+`$, with $`v_{}=v_+>0`$ and an interaction range which we assume to be a Gaussian $`f(xx_n)`$ centered at the monomer site $`x_n`$,

$$V(x)=v_nf(xx_n),$$
(15)
with $`f(x)=\mathrm{exp}[x^2/2\sigma ^2]`$ and $`\sigma a`$. The $`v_n`$ are randomly drawn from a pdf $`p(v)=c\delta (vv_{})+(1c)\delta (vv_+)`$. The relative concentration of low and high energy monomers is given by $`c`$ and $`1c`$, respectively. The parameter $`c`$ determines the shape of the overall potential. When $`c<1/2`$ ($`c>1/2`$) the potential consists of a series of localized peaks (troughs). Mean and variance of the potential are $`V=(1c)v_{}+cv_+`$ and $`(VV)^2=(1c)c\delta v^2`$ with $`\delta v=v_+v_{}`$. Figs. 4A and 4B depict the results obtained for the generalized diffusion coefficient $`D_\mu (ϵ)`$ on these types of copolymers for three values of $`c`$, each one representing one of the situations depicted in Fig. 2C. The parameter $`\delta v`$ was chosen such that the variance is identical in all potentials. Although the value of $`D_\mu (ϵ)`$ in the weak potential regime ($`ϵ1`$) is consistent with the one observed in random phase potentials, for greater values of $`ϵ`$ a striking deviation occurs. On one hand, the ordinary diffusion process ($`\mu =2`$) is nearly insensitive to the shape of the potential, all functions $`D_2(ϵ)`$ coincide. On the other hand, the superdiffusive process exhibits a more (less) pronounced decrease with increasing $`ϵ`$ when $`c=0.95`$ ($`c=0.05`$) as compared to an unbiased concentration of monomer types.

The results reported in this Letter predict for superdiffusive behavior on folding polymers that the dispersion speed depends strongly on the specific arrangement of various types of monomers. This is in sharp contrast to the case of ordinary diffusion, which solely depends on magnitude variations of the potential. Therefore we expect that potential heterogeneity is an essential ingredient in superdiffusive translocation phenomena of proteins along biopolymers.

###### Acknowledgements.

D. B. thanks W. Noyes for interesting comments and discussion.


text:
# Are we missing baryons in Galaxy Clusters?

## 1 INTRODUCTION

The recent analysis of the angular power spectrum of the Cosmic Microwave Background (CMB) obtained from WMAP (Bennett et al. 2003) has provided constraints on the cosmological parameters (Spergel et al. 2003) that confirms with greater accuracy the current energy density of the Universe to be comprised by about 73 per cent of dark energy and 27 per cent of matter, mostly non-baryonic and dark. In particular, the quoted constraint on the baryon density, $`\mathrm{\Omega }_\mathrm{b}`$, is $`0.0224\pm 0.0009h_{100}^2`$, and on the total matter density, $`\mathrm{\Omega }_\mathrm{m}`$, is $`0.135_{0.009}^{+0.008}h_{100}^2`$. Consequently, the cosmic baryon fraction, $`\mathrm{\Omega }_\mathrm{b}/\mathrm{\Omega }_\mathrm{m}`$, is equal to $`0.166_{0.013}^{+0.012}`$ and the ratio between baryon and cold dark matter density, $`\mathrm{\Omega }_\mathrm{c}=\mathrm{\Omega }_\mathrm{m}\mathrm{\Omega }_\mathrm{b}`$, is equal to $`0.199_{0.019}^{+0.017}`$. These values are expected to be maintained in regions at high overdensities that collapse to form galaxy clusters.

The clusters baryon budget is composed mainly from the X-ray luminous baryons, $`M_{\mathrm{gas}}`$, of the intracluster medium (ICM) that becomes hotter upon falling into the cluster dark matter halo by gravitational collapse. Other contributions come from the baryonic stellar mass in galaxies, $`M_{\mathrm{gal}}`$, and from other “exotic” sources, like intergalactic stars and a still poorly defined baryonic dark matter. Given the large uncertainties on the relative contribution from baryons that are not accounted for in either $`M_{\mathrm{gas}}`$ or $`M_{\mathrm{gal}}`$, I qualify these as “other baryons”, $`M_{\mathrm{ob}}`$, as already done in a previous work (Ettori 2001) in which I discussed the constraints on cluster baryon budget from BOOMERANG and MAXIMA-I data. The tighter constraints on the cosmological parameters provided from WMAP allow now more firm conclusions.

Therefore, one can put the following relation between the relative amount of baryons in the Universe and in clusters with total gravitating mass, $`M_{\mathrm{tot}}`$:

$$\frac{\mathrm{\Omega }_\mathrm{b}}{\mathrm{\Omega }_\mathrm{m}}=\frac{M_\mathrm{b}}{M_{\mathrm{tot}}}=\frac{f_{\mathrm{gas}}}{YBC}+\frac{f_{\mathrm{gal}}}{B}+\frac{f_{\mathrm{ob}}}{B},$$
(1)
where $`f_{\mathrm{gas}}=M_{\mathrm{gas}}/M_{\mathrm{tot}}`$, $`f_{\mathrm{gal}}=M_{\mathrm{gal}}/M_{\mathrm{tot}}`$ ($`0.010_{0.004}^{+0.005}h_{100}^1`$ in Fukugita et al. 1998), $`f_{\mathrm{ob}}=M_{\mathrm{ob}}/M_{\mathrm{tot}}`$, and $`Y`$ is the parameter representing the cosmic depletion of baryons at the virial radius with respect to the global value ($`0.92\pm 0.06`$ from the hydrodynamical simulations of the Santa Barbara Project; Frenk et al. 1999). I parametrize the uncertainties on the measurements of the total gravitating mass and gas mass through the factors $`B`$ and $`C`$, respectively. These factors act to increase the total mass estimates (i.e. $`B>1`$) if corrections to the hydrostatic equilibrium equation are required for bulk motions of the ICM or non-thermal pressure support, and to lower the true gas mass (i.e. $`C>1`$) if clumpiness is present in the ICM that is assumed to be smoothly distributed (e.g. Mathiesen et al. 1999).

In this Letter, I will analyze the equation 1 (i) to assess the consistency between the cosmic and the cluster baryon budget and (ii) to put significant constraints on $`f_{\mathrm{ob}}`$. In a consistent way, I adopt the WMAP results on the Hubble constant, $`H_0`$, of $`71_{0.03}^{+0.04}`$ km s-1 Mpc-1 to rescale all the measured quantities.

## 2 The cosmic baryon budget from WMAP

Figure 1 shows the allowed $`1\sigma `$ region from the WMAP results on the cosmic baryon budget with respect to the observed gas and baryon fraction for a sample of relaxed galaxy clusters both at low and high redshift. All these values are estimated at the overdensity of 200 with respect to the critical density (i.e. within the cluster region that numerical simulations show to be virialized in a $`\mathrm{\Omega }_\mathrm{m}`$ independent way; see, e.g., Evrard et al. 2002). Using a Bayesian approach (Press 1996), I measure the average gas (baryon) fraction and confirm that is consistent between the two samples: $`0.107_{0.019}^{+0.028}(0.133_{0.026}^{+0.029})`$ and $`0.111_{0.063}^{+0.069}(0.133_{0.067}^{+0.063})`$ (the error-weighted means of the gas fraction distribution are $`0.116\pm 0.005`$ and $`0.111\pm 0.010`$ in the low$`z`$ and high$`z`$ sample, respectively. Note that these values do not include any correction by the baryonic depletion that is expected to raise the gas/baryon fraction by about 8/6 per cent at this overdensity. See also Fig. 1 and text that follows). These estimates are consistent with other, independent, recent determinations (e.g. Allen et al. 2001, Pratt & Arnaud 2002) and consistently lower than the baryon budget required from WMAP results. It is worth noticing that only the highest estimates of the gas (baryon) mass fraction (e.g. A426, A2142, RXJ1350) are perfectly consistent with WMAP results, whereas the other clusters are systematically below them.

To investigate the systematics that could affect this estimate, I change the values of the factors $`B`$ and $`C`$ and study their influence on the baryon fraction. The factor $`B`$, which parametrizes the uncertainties on $`M_{\mathrm{tot}}`$ is expected to be between 1 and 1.15 from the cluster mass profiles recovered from both X-ray and lensing data (e.g. Allen et al. 2001). The factor $`C`$ represents the level of clumpiness that affects the estimate of $`M_{\mathrm{gas}}`$ in X-ray analysis and that simulations show to be lower than 1.2 (Mathiesen et al. 1999).

Figure 2 shows that higher values of $`B`$ and $`C`$ require a more relevant role to be played by $`f_{\mathrm{ob}}`$, giving a $`2\sigma `$ positive detection for typical value of $`Y`$ when $`B`$ and $`C`$ are 15 per cent larger than the null hypothesis of reliable estimates of both $`M_{\mathrm{gas}}`$ and $`M_{\mathrm{tot}}`$ from X-ray analysis.

Moreover, I can estimate from the observables the ratio $`R_{\mathrm{gas}}=(\mathrm{\Omega }_\mathrm{b}/\mathrm{\Omega }_\mathrm{m}f_{\mathrm{gal}}/B)/f_{\mathrm{gas}}\times (YBC)(f_{\mathrm{ob}}/f_{\mathrm{gas}}+1)`$ and to evaluate the probability that $`R_{\mathrm{gas}}>1`$ and a no-zero value of $`f_{\mathrm{ob}}`$ is required from the data (see Fig. 3). I obtain values of $`R_{\mathrm{gas}}`$ between 1.3 ($`B=1,C=1`$) and 1.7 ($`B=1.15,C=1.15`$), with an interval accepted at the 95 per cent confidence level of 0.5–2.7. More significantly, this ratio has to be larger than 1 at 76.6 ($`B=1,C=1`$) and 92.5 ($`B=1.15,C=1.15`$) per cent confidence level. This result gives a high confidence to the conclusion that a significant amount of baryons has to be present apart from what is observed both at X-ray and optical wavelength.

## 3 Conclusions

By comparing the recent cosmological constraints from measurements of the angular power spectrum of the temperature anisotropy in the CMB done with WMAP with the observed distribution of the gas mass fraction in clusters of galaxies, I conclude that

1. galaxy clusters with the highest observed gas/baryon fraction are well in agreement with the WMAP estimate of the cosmic baryon budget. On average, however, a disagreement in the order of 15–20 per cent is present, the observed cluster baryon fraction being a lower estimate of the cosmic one. This implies that estimates of the cold dark matter density, $`\mathrm{\Omega }_\mathrm{c}`$, done by applying mean results from large sample of objects tend to overestimate it. In this perspective, a more “realistic” result is provided from the highest estimate in the distribution of the cluster gas (baryon) mass fraction;
2. this dark baryonic component appears to be between a fraction, and up to 2 times, the measured gas fraction. Values of $`R_{\mathrm{gas}}=(f_{\mathrm{WMAP}}f_{\mathrm{gal}})/f_{\mathrm{gas}}`$ larger than 1 are required from the data with a level of confidence of about 80 per cent and more, if we are underestimating (overestimating) the total (gas) mass. The most probable values are $`f_{\mathrm{ob}}/f_{\mathrm{gas}}=R_{\mathrm{gas}}1=`$ 0.3 ($`B=C=1`$) and 0.7 ($`B=C=1.15`$), and lower than 1.7 at 95 per cent confidence level.

It is very unlikely that galactic objects, such as “colored” (red, brown, white, beige) dwarfs, stellar remnants and other species of MACHOs (see review in Gilmore 1999 and Evans 2003), or intergalactic ones formed from tidal disruption of cluster dwarfs, like planetary nebulae (e.g. Ciardullo et al. 2002), red-giant-branch stars (e.g. Ferguson et al. 1998) and supernovae (Gal-Yam et al. 2003), can be responsible for such amount of baryons. It is reasonable to believe that they can contribute by about 0.2 times $`f_{\mathrm{gal}}`$, or $`0.04f_{\mathrm{gas}}`$.

The most plausible suspect to give so large contribution is then a X-ray warm ($`10^5`$ K $`<T<`$ $`10^7`$ K) intracluster medium (W-ICM). Large scale cosmological and hydrodynamical simulations by Cen & Ostriker (1999) and Davé et al. (2001) show that the mass fraction at redshift 0 is largely dominated from a warm medium, with a relative contribution in mass that is about 2 times larger than the amount of hot ($`T>10^7`$ K) baryons. However, less than 30 per cent of it falls in overdensities $`60`$ that are typical for bound structures in a $`\mathrm{\Lambda }`$CDM universe. Furthermore, Bonamente et al. (2002) present evidence of excess in the soft X-ray emission between 0.2–0.4 keV in 50 per cent of the 38 clusters in their sample of Rosat PSPC observations. They list several suggestions on how to explain this excess, originally observed in extreme ultraviolet (Lieu et al. 1996), both as thermal and non-thermal component. If we assign this emission to the baryons that are lacking in our budget, we interpret it as thermal emission due either to the diffuse/halo component of unresolved X-ray faint cluster galaxies or W-ICM. In the first case, we are forced to consider an inexplicable large amount of X-ray emitting member galaxies. More plausible is then the hypothesis of W-ICM, even though its cooling time tends to be very short with the bulk of the radiation in emission line if this gas is not primordial. Fabian (1997) suggested that it can be located in turbulent mixing layers lying between embedded cold clouds and the ICM. However, the traditional picture on the efficiency of cooling processes in the cluster cores is not supported anymore after that XMM and Chandra observations did not report evidence of gas cooler than 1–2 keV (e.g. Peterson et al. 2003) and showed a strong interplay between ICM, the central active galaxy (e.g. Fabian 2002) and merging cool clumps (e.g. Markevitch et al. 2000, Mazzotta et al. 2003). On the other side, the production of thermal energy per particle due to supernovae related to the star formation activity is in the order of $`0.4(\eta /0.1)(N_{\mathrm{SNII}}/10^9)(10^{13}M_{}/M_{\mathrm{gas}})`$ keV for a given efficiency $`\eta `$ in converting the kinetic energy of the explosion into thermal energy through galactic winds, and adopting typical values of the cluster gas mass and number of type II supernovae as required from the observed ICM metallicity. While this energy per gas particle is not enough to stop cooling the hot ICM, it can easily accommodate for the survival of the warm component.

Intriguingly, the stronger soft excess detected in Rosat data from Bonamente et al. (2002) is measured in objects like A85 and A1795, that are the ones lacking most of the baryons with respect to the cosmic budget as plotted in Fig. 1. Significant detection is also present in A2029, A2199 and A3571, whereas a marginal detection is associated to A2142. Nevalainen et al. (2003; see also Kaastra et al. 2003) confirm with XMM the excess in the soft X-ray emission in A1795 and that this excess is best fitted by a thermal component with a characteristic temperature of 0.8 keV, which is about an order of magnitude higher than what required from Rosat data but still consistent with our energetic arguments. Using their estimation of the atom density of the W-ICM in the core of A1795, and assuming that it is broadly distributed like the ICM, one can infer a $`f_{\mathrm{W}\mathrm{ICM}}/f_{\mathrm{gas}}=R_{\mathrm{gas}}10.43`$. In general, values of $`f_{\mathrm{W}\mathrm{ICM}}/f_{\mathrm{gas}}`$ between 0.1 and 0.5 are expected.

To summarize, clusters seem to have similar behavior in accumulating the same relative amount of baryons. It is then their peculiar thermal history due to the interplay of merging actions and/or activity of the central active galaxy that provides the baryonic ingredients and cook the baryonic pie that we taste and show in Figure 4. To prepare it, I have considered only the eight nearby clusters which provide a more reliable estimate of $`f_{\mathrm{gal}}`$ and are less affected from systematics in the determination of $`f_{\mathrm{gas}}`$ (see discussion in Ettori et al. 2003). I have also corrected the gas fraction by the depletion factor $`Y0.92`$. The baryonic pie is then made of 70 per cent of hot ICM, with $`1\sigma `$ range between 56 and 89 per cent and a distribution of the calculated $`f_{\mathrm{gas}}/f_{\mathrm{WMAP}}`$ that spans between 28 and 143 per cent at $`2\sigma `$ confidence level (higher upper limit observed in A426 that has a most probable $`f_{\mathrm{gas}}/f_{\mathrm{WMAP}}`$ value of 115 per cent). The cold, stellar component is responsible for 13 ($`1\sigma `$: 8–19) per cent, with an observed distribution in the sample between 2 and 37 per cent ($`2\sigma `$ lower and upper limit, with the latter reached in A2199, which has a central value of 21 per cent). Finally, a third ingredient, probably a warm ICM, contributes by about 17 (0–33) per cent (and a probability to be larger than 0 of 73 per cent) with a distribution that goes from –29$`\pm `$15 per cent in A426 to 40$`\pm `$12 per cent in A1795, one of the objects with the largest detected soft excess (Bonamente et al. 2002, Kaastra et al. 2003).

## ACKNOWLEDGEMENTS

I am grateful to Hans Böhringer and Andy Fabian that pointed out a relevant weakness in the assumptions done in the original manuscript.


text:
# Inflation model with lower multipoles of the CMB suppressed

## ACKNOWLEDGMENTS

We acknowledge the use of CMBFAST cmbfast  package for our numerical calculations. F.T. thanks the Japan Society for the Promotion of Science for financial support.