text:
# Ultra-High Energy Neutrino Fluxes and Their Constraints

## I Introduction

The farthest source so far observed in neutrinos was supernova SN1987A from which about 20 neutrinos in the 10-40 MeV range where detected sn1987A . Extraterrestrial neutrinos of much higher energy are usually expected to be produced as secondaries of cosmic rays interacting with ambient matter and photon fields and should thus be associated with cosmic ray sources ranging from our Galaxy to powerful active galactic nuclei (AGN) review\_nu . Traditional neutrino telescopes now under construction aim to detect such neutrinos up to $`10^{16}`$eV by looking for showers and/or tracks from charged leptons produced by charged current reactions of neutrinos in ice, in case of AMANDA amanda ; amanda\_limit  and its next generation version ICECUBE icecube , in water, in case of BAIKAL baikal ; baikal\_limit , ANTARES antares , and NESTOR nestor , or underground, in case of MACRO MACRO  (for recent reviews of neutrino telescopes see Ref. nu\_tele ).

On the other hand, the problem of the origin of the cosmic rays themselves is still unsolved, especially at ultra-high energies (UHE) above $`4\times 10^{19}`$eV, where they lose energy rapidly by pion production and pair production (protons only) on the cosmic microwave background (CMB) reviews ; bs . For sources further away than a few dozen Mpc this would predict a break in the cosmic ray flux known as Greisen-Zatsepin-Kuzmin (GZK) cutoff gzk , around $`50`$EeV. This break has not been observed by experiments such as Fly’s Eye Fly , Haverah Park Haverah , Yakutsk Yakutsk  and AGASA agasa , which instead show an extension beyond the expected GZK cutoff and events above $`100`$EeV (however, the new experiment HiRes Hires  currently seems to see a cutoff in the monocular data Hires2001 ). This has lead to the current construction of the southern site of the Pierre Auger Observatory auger , a combination of an array of charged particle detectors with fluorescence telescopes for air showers produced by cosmic rays above $`10^{19}`$eV which will lead to an about hundred-fold increase of data. The telescope array, another planned project based on the fluorescence technique, may serve as the optical component of the northern Pierre Auger site planned for the future ta . There are also plans for space based observatories such as EUSO euso  and OWL owl  of even bigger acceptance. These instruments will also have considerable sensitivity to neutrinos above $`10^{19}`$eV, typically from the near-horizontal air-showers that are produced by them auger\_nu . Furthermore, the old Fly’s Eye experiment baltrusaitis  and the AGASA experiment agasa\_nu  have established upper limits on neutrino fluxes based on the non-observation of horizontal air showers.

In addition, there are plans to construct telescopes to detect fluorescence and Čerenkov light from near-horizontal showers produced in mountain targets by neutrinos in the intermediate window of energies between $`10^{15}`$eV and $`10^{19}`$eV fargion ; mount . The alternative of detecting neutrinos by triggering onto the radio pulses from neutrino-induced air showers is also investigated currently radhep . Two implementations of this technique, RICE, a small array of radio antennas in the South pole ice rice , and the Goldstone Lunar Ultra-high energy neutrino Experiment (GLUE), based on monitoring of the moons rim with the NASA Goldstone radio telescope for radio pulses from neutrino-induced showers glue , have so far produced neutrino flux upper limits. Acoustic detection of neutrino induced interactions is also being considered acoustic .

The neutrino detection rates for all future instruments will crucially depend on the fluxes expected in various scenarios. The flux of ”cosmogenic” neutrinos created by primary protons above the GZK cutoff in interactions with CMB photons depends both on primary proton spectrum and on the location of the sources. The cosmogenic neutrino flux is the only one that is guaranteed to exist just by the observations of ultra-high energy cosmic rays (UHECRs) and was studied soon after the discovery of the CMB cosmogenic . Note, however, that there is no firm lower bound on the cosmogenic neutrino flux if the UHECR sources are much closer than the GZK distance.

If sources are located beyond the GZK distance and the proton flux extends beyond the GZK cutoff, the neutrino fluxes can be significant. This possibility is favored by the lack of nearby sources and by the hardening of the cosmic ray spectrum above the “ankle” at $`5\times 10^{18}`$eV. It is also suggested by possible correlations of UHECR arrival directions with compact radio quasars corr\_radio  and more significant correlations with BL Lacertae objects corr\_bllac , some of which possibly also luminous in GeV $`\gamma `$rays egret\_bllac  and at distances too large to be consistent with the absence of the GZK cutoff.

Whereas roughly homogeneously distributed proton sources can naturally explain the UHECR flux below the GZK cutoff (see, e.g., Refs. berezinsky2002  and photons ), the highest energy events may represent a new component. They may be new messenger particles, which propagate through the Universe without interacting with the CMB messengers , or may have originated as extremely high energy ($`E10^{23}`$eV) photons, which can propagate several hundred Mpc (constantly loosing energy) and can create secondary photons inside the GZK volume photons . Decaying super-heavy relics from the early Universe (see Ref. bs ; kt  for reviews) can also explain UHECRs and predict UHE neutrino fluxes detectable by future experiments.

Another speculative possibility is to explain UHECRs beyond the GZK cutoff by the UHE protons and photons from decaying Z-bosons produced by UHE neutrinos interacting with the relic neutrino background zburst1 . The big drawback of this scenario is the need of enormous primary neutrino fluxes that cannot be produced by known astrophysical acceleration sources without overproducing the GeV photon background kkss , and thus most likely requires a more exotic top-down type source such as X particles exclusively decaying into neutrinos gk . As will be shown in Sect. VI, even this possibility is also significantly constrained by existing measurements.

Active galactic nuclei (AGN) can be sources of neutrinos if protons are accelerated in them review\_nu . In the present paper we consider only the two representative limits of low and high optical depth for pion (and neutrino) production in the source. In the first case the protons accelerated in the AGN freely escape and neutrinos are produced only in interactions with the CMB (cosmogenic neutrinos). For the second case we discuss an example of possible high neutrino fluxes from a non-shock acceleration AGN model neronov , in which primary protons lose all their energy and produce neutrinos directly in the AGN core.

Motivated by the increased experimental prospects for ultra-high energy neutrino detection, in the present paper we reconsider flux predictions in the above scenarios with our recently combined propagation codes code1 ; code2 ; photons ; kkss . Our main emphasize is thereby on model independent flux ranges consistent with all present data on cosmic and $`\gamma `$rays. For any scenario involving pion production the fluences of the latter are comparable to the neutrino fluences. However, electromagnetic (EM) energy injected above $`10^{15}`$eV cascades down to below the pair production threshold for photons on the CMB and EM energy above 100 GeV also cascades down due to the pair production on infrared/optical background. The resulting intensity and spectrum of $`\gamma `$rays below 100 GeV is rather insensitive to these backgrounds berez75 ; ca  (for review see bs ). The cascade thus gives rise to a diffuse photon flux in the GeV range which is constrained by the flux observed by the EGRET instrument on board the Compton $`\gamma `$ray observatory egret . For all neutrino flux scenarios the related $`\gamma `$ray and cosmic ray fluxes have to be consistent with the EGRET and cosmic ray data, respectively.

Sect. II summarizes the numerical technique used in this paper. In Sect. III we discuss the cosmogenic neutrino flux and its dependence on various source characteristics. We specifically find an upper limit considerably higher than typical fluxes in the literature and remark why it is higher than the Waxman-Bahcall (WB) wb  and even the Mannheim-Protheroe-Rachen (MPR) mpr  bounds for sources transparent to cosmic and $`\gamma `$rays. In Sect. IV we review neutrino flux predictions in top down scenarios where UHECRs are produced in decays of super-massive particles continuously released from topological defect relics from the early Universe. Sect. V discusses neutrino fluxes in scenarios where the cosmic rays observed at the highest energies are produced as secondaries of these neutrinos from interactions with the relic cosmological neutrino background, often called Z-burst scenario. In Sect. VI we focus on a combination of top-down and Z-burst scenarios considered by Gelmini and Kusenko gk , namely super-heavy particles mono-energetically decaying exclusively into neutrinos”(see, however, berkaros ). In Sect. VII we discuss possible high neutrino fluxes from a non-shock acceleration AGN model neronov . Finally, in Sect. VIII we conclude.

## II Numerical Technique

Our simulations are based on two independent codes that have extensively been compared down to the level of individual interactions. Both of them are implicit transport codes that evolve the spectra of nucleons, $`\gamma `$rays, electrons, electron-, muon-, and tau-neutrinos, and their antiparticles along straight lines. Arbitrary injection spectra and redshift distributions can be specified for the sources and all relevant strong, electromagnetic, and weak interactions have been implemented. For details see Refs. code1 ; code2 ; kkss .

Relevant neutrino interactions for the Z-burst scenario are both the s-channel production of Z bosons and the t-channel production of W bosons. The decay products of the Z boson were taken from simulations with the PYTHIA  Monte Carlo event generator using the tuned parameter set of the OPAL Collaboration OPAL .

The main ambiguities in propagation of photons concern the unknown rms magnetic field strength $`B`$ which can influence the predicted $`\gamma `$ray spectra via synchrotron cooling of the electrons in the EM cascade, and the strength of the universal radio background (URB) which influences pair production by UHE $`\gamma `$rays pb . Photon interactions in the GeV to TeV range are dominated by infrared and optical universal photon backgrounds (IR/O), for which we took the results of Ref. Primack . The resulting photon flux in GeV range in not sensitive to details of the IR/O backgrounds.

Concerning the cosmology parameters we chose the Hubble parameter $`H_0=70\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$ and a cosmological constant $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, as favored today. These values will be used in all cases unless otherwise indicated.

For the neutrinos we assume for simplicity that all three flavors are completely mixed as suggested by experiments superK  and thus have equal fluxes. For each flavor we sum fluxes of particles and antiparticles.

Predictions for the all fluxes in both codes agree within tens of percents. Only for the Z-burst scenarios the photon fluxes agree only within a factor $`2`$ between the two original codes due to the different implementation of Z-decay spectra and the ambiguities in photon propagation mentioned above. This difference has no influence on the conclusions of this paper.

In the present investigation we parameterize power law injection spectra of either protons (for UHECR sources) or neutrinos (for Z-burst models) per co-moving volume in the following way:

$`\varphi (E,z)`$ $`=`$ $`f(1+z)^mE^\alpha \mathrm{\Theta }(E_{\mathrm{max}}E)`$ (1)
$`z_{\mathrm{min}}zz_{\mathrm{max}},`$
where $`f`$ is the normalization that has to be fitted to the data. The free parameters are the spectral index $`\alpha `$, the maximal energy $`E_{\mathrm{max}}`$, the minimal and maximal redshifts $`z_{\mathrm{min}}`$, $`z_{\mathrm{max}}`$, and the redshift evolution index $`m`$. The resulting neutrino spectra depend insignificantly on $`z_{\mathrm{min}}`$ in the range $`0z_{\mathrm{min}}0.1`$ where local effects could play a role, and thus we will set $`z_{\mathrm{min}}=0`$ in the following.

To obtain the maximal neutrino fluxes for a given set of values for all these parameters , we determine the maximal normalization $`f`$ in Eq. (1) by demanding that both the accompanying nucleon and $`\gamma `$ray fluxes are below the observed cosmic ray spectrum and the diffuse $`\gamma `$ray background observed by EGRET, respectively.

## III The Cosmogenic Neutrino Flux

### III.1 Dependence on unknown parameters

In this section we discuss the case when primary UHECRs produce cosmogenic neutrinos as well as $`\gamma `$rays during propagation. For EM propagation we use $`B=10^9`$ G and the intermediate URB strength estimate of Ref. pb . These parameters only influence the $`\gamma `$ray flux at UHEs, but not in the GeV range where the flux only depends on the total injected EM energy. Therefore, in this scenario the resulting neutrino fluxes are insensitive to the poorly known UHE $`\gamma `$ray absorption because the “visible” UHE flux is always dominated by the primary cosmic rays and not by the secondary $`\gamma `$ray flux, as can be seen in Figs. 9 and 10 below.

Figs. 1, 2, 3, 4, and 5 show the dependencies of the cosmogenic neutrino flux (average per flavor) for which associated cosmic and $`\gamma `$ray fluxes are consistent with the data, on the maximal source redshift $`z_{\mathrm{max}}`$, maximal injection energy $`E_{\mathrm{max}}`$, redshift evolution index $`m`$, spectral power law index $`\alpha `$, and cosmological parameters, respectively. In each figure the line for the highest neutrino flux corresponds to a significant contribution of the accompanying $`\gamma `$-ray flux to the observed flux at the EGRET region, whereas for the other lines the $`\gamma `$-ray flux gives negligible contributions to the EGRET flux.

Fig. 2 shows that a change of the primary proton maximum energy changes the secondary neutrino flux only in the high energy region. The reason for this behavior is the shape of the pion production cross section which at the lowest energies is dominated by the single pion $`\mathrm{\Delta }`$resonance. Figs. 1 and 3 show that the cosmogenic neutrino flux at the lowest energies mostly depends on the maximum redshift $`z_{\mathrm{max}}`$ and the evolution index $`m`$. This is especially relevant for experiments with their main sensitivity below $`10^{19}`$eV such as ICECUBE (see Fig. 10 below). Fig. 4 shows that the maximal neutrino flux significantly increases with decreasing proton injection power law index. Fig. 5 shows that the variation of the maximal neutrino flux with cosmology parameters $`\mathrm{\Omega }_\mathrm{\Lambda }`$ and $`H_0`$ is rather modest, about 50% for values discussed in recent years.

### III.2 Active galactic nuclei as UHECR sources

Here we consider AGN sources for the primary UHECR flux, with the typical evolution parameters $`m=3.4`$ for $`z<1.9`$ and $`m=0`$ for $`1.9<z<2.9`$ agn\_evol . The only free remaining parameters are the power law index $`\alpha `$ and the maximum energy $`E_{\mathrm{max}}`$ for the proton injection spectrum, and the flux normalization $`f`$ in Eq. (1). We first checked that for the case $`\alpha =2`$ we agree with the WB bound. Our cosmogenic neutrino fluxes agree reasonably well with the corresponding ones shown in Fig. 4 of Ref. ess  when taking into account maximal mixing assumed in the present paper.

Fig. 6 shows the maximal cosmogenic neutrino flux for sources with hard spectra, $`\alpha =1`$, as a function of the maximal proton energy $`E_{\mathrm{max}}`$ and a mono-energetic injection spectrum at $`E=3\times 10^{21}`$eV. For the case $`\alpha =1`$, MPR mpr  computed the sum of the cosmogenic neutrino flux and the neutrino flux directly emitted by the sources which are assumed to be transparent to neutrons. Our fluxes shown here only include the cosmogenic flux. Nevertheless our fluxes overshoot the comparable MPR fluxes by up to a factor 5 at energies below the peak flux. This is most likely due to a combination of the different cosmology (see Fig. 5) and a different implementation of multi-pion production which influences interactions of nucleons at high energies, thus at high redshift and in turn the low energy tail of cosmogenic neutrinos.

Fig. 6 and also Figs. 2-4 demonstrate in a general way that it is easy to exceed the WB bound and even the MPR bound for injection spectra harder than about $`E^2`$. This is because Waxman-Bahcall restricted themselves to nucleon injection spectra softer than $`E^2`$ and sources smaller than nucleon interaction lengths wb . Thus, their bound does not directly apply to the cosmogenic neutrino flux. In addition, in our opinion, their assumptions on the injection spectra are too narrow: Possible scenarios with hard injection spectra and the AGN redshift evolution assumed here (which is the same as the one used by Waxman-Bahcall) include cases where UHECRs are accelerated by the electromotive force produced by magnetic fields threading the horizons of spinning super-massive black holes in the centers of galaxies bg  or by reconnection events around forming galaxies colgate .

In any scenario involving pion production for the creation of $`\gamma `$rays and neutrinos, the fluxes per flavor are approximately related by $`F_\nu (E)F_\gamma (E)/3`$. Assuming smooth spectra and comparing with the EGRET $`\gamma `$ray fluence, energy conservation implies

$$E^2F_\nu (E)6\times 10^2\mathrm{eVcm}^2\mathrm{s}^1\mathrm{sr}^1.$$
(2)
This ultimate bound is also shown in Fig. 6 and Fig. 7 below. It corresponds to the MPR limit for optically thick sources. The maximal $`E^2j(E)`$ of the fluxes in Figs. 6 and 7 below indeed reach this $`\gamma `$ray bound Eq. (2).

Note that the two theoretical bounds shown in Fig. 6 represent fluxes per neutrino flavor under the assumption of maximal neutrino mixing. They are thus about a factor 2 lower than in Refs. wb ; mpr  which show muon neutrino fluxes in the absence of mixing where the tau-neutrinos fluxes are negligible and electron neutrino fluxes are about a factor 2 smaller. The WB and MPR bounds represent upper neutrino flux limits for compact sources such as AGN and $`\gamma `$ray bursts in case of small optical depth for nucleons. We discuss a specific non-shock acceleration scenario in Sect. VII, for which both of the above bounds are not valid because the source optical depth for nucleons is large (for reviews on AGN neutrino fluxes see, e.g., Ref. review\_nu ).

### III.3 General case of arbitrary source evolution

In this section we consider more general UHECR sources and relax the restrictions on their redshift evolution. Fig. 7 shows that cosmogenic neutrino fluxes higher than both the WB and MPR limits are possible even for relatively soft $`E^2`$ proton injection spectra, if the redshift evolution is stronger than for AGNs: The curve for $`E^2`$ in Fig. 7 corresponds to the evolution parameters $`m=5`$, $`z_{\mathrm{max}}=3`$ and $`E_{\mathrm{max}}=10^{22}`$eV, the curve for $`E^{1.75}`$ to $`m=4.5`$, $`z_{\mathrm{max}}=3`$ and $`E_{\mathrm{max}}=10^{23}`$eV, and the curve for mono-energetic injection to $`m=4`$, $`z_{\mathrm{max}}=3`$, and $`E_{\mathrm{max}}=10^{21}`$eV.

Fig. 8 shows that between $`10^{18}`$eV and $`10^{20}`$eV the energy loss rate of protons on the CMB is dominated by pair production instead of pion production. The former does not contribute to neutrino production but the EM cascades initiated by the pairs lead to contributions to the diffuse $`\gamma `$ray background in the GeV range. Thus, the cosmogenic neutrino flux is the more severely constrained the bigger the fraction of cosmic ray power is in the range $`10^{18}\mathrm{eV}E10^{20}`$eV. This is mostly important for soft injection spectra and explains why the total neutrino energy fluence decreases with increasing $`\alpha `$ in Fig. 7.

### III.4 Comparison with experimental limits and future sensitivities

Figs. 9 and 10 shows a scenario maximized over all 4 parameters in comparison to existing neutrino flux upper limits and expected sensitivities of future projects, respectively. Both of these fall into two groups, detection in water or ice or underground, typically sensitive below $`10^{16}`$eV, and air shower detection, usually sensitive at higher energies. Existing upper limits come from the underground MACRO experiment MACRO  at Gran Sasso, AMANDA amanda\_limit  in the South Pole ice, and the Lake BAIKAL neutrino telescope baikal\_limit  in the first category, and the AGASA ground array agasa\_nu , the former fluorescence experiment Fly’s Eye baltrusaitis , the Radio Ice Čerenkov Experiment RICE rice , and the Goldstone Lunar Ultra-high energy neutrino experiment GLUE glue  in the second category. Future experiments in the first category include NT200+ at Lake Baikal baikal\_limit , ANTARES in the Mediterranean antares , NESTOR in Greece nestor , AMANDA-II amandaII , and ICECUBE icecube , the next-generation version of AMANDA, whereas the air shower based category includes the Auger project auger\_nu , the Japanese telescope array ta\_nu , the fluorescence/Čerenkov detector MOUNT mount , and the space based OWL owl\_nu  and EUSO euso  experiments. The vertical bars for the MOUNT sensitivity characterize the uncertainties due to the not yet determined zenith angle range and sensitivity to the fluorescence component. The OWL sensitivity estimate is based on deeply penetrating atmospheric showers induced by electron or muon-neutrinos only owl\_nu  and may thus be considerably better if tau neutrinos, Čerenkov events, and Earth skimming events are taken into account skim . The same applies to the EUSO project euso . The AMANDA-II sensitivity will lie somewhere between the ANTARES and ICECUBE sensitivities amandaII . The maximized fluxes shown in Figs. 9 and 10 are considerably higher than the ones discussed in Refs. cosmogenic ; stecker ; pj ; ydjs ; ess , and should be easily detectable by at least some of these future instruments, as is obvious from Fig. 10.

## IV Neutrino Fluxes in Top-Down Scenarios

Historically, top-down scenarios were proposed as an alternative to acceleration scenarios to explain the huge energies up to $`3\times 10^{20}`$eV observed in the cosmic ray spectrum bhs . In these top-down scenarios UHECRs are the decay products of some super-massive “X” particles of mass $`m_X10^{20}`$eV close to the grand unified scale, and have energies all the way up to $`m_X`$. Thus, The massive X particles could be metastable relics of the early Universe with lifetimes of the order the current age of the Universe or could be released from topological defects that were produced in the early Universe during symmetry-breaking phase transitions envisaged in Grand Unified Theories (GUTs). The X particles typically decay into leptons and quarks. The quarks hadronize, producing jets of hadrons which, together with the decay products of the unstable leptons, result in a large cascade of energetic photons, neutrinos and light leptons with a small fraction of protons and neutrons, some of which contribute to the observed UHECR flux. The resulting injection spectra have been calculated from QCD in various approximations, see Ref. bs  for a review and Ref. frag  for more recent work. In the present work we will use the spectra shown in Fig. 11 which are obtained from solving the DGLAP equations in modified leading logarithmic approximation (MLLA) without supersymmetry for X particles decaying into two quarks, assuming 10% nucleons in the fragmentation spectrum.

For dimensional reasons the spatially averaged X particle injection rate can only depend on the mass scale $`m_X`$ and on cosmic time $`t`$ in the combination

$$\dot{n}_X(t)=\kappa m_X^pt^{4+p},$$
(3)
where $`\kappa `$ and $`p`$ are dimensionless constants whose value depend on the specific top-down scenario bhs . Extragalactic topological defect sources usually predict $`p=1`$, whereas decaying super-heavy dark matter of lifetime much larger than the age of the Universe corresponds to $`p=2`$ bs . In the latter case the observable flux will be dominated by decaying particles in the galactic halo and thus at distances smaller than all relevant interaction lengths. Composition and spectra will thus be directly given by the injection spectra which are most likely inconsistent with upper limits on the UHE photon fraction above $`10^{19}`$eV gamma\_comp , see Fig. 11. In addition, decaying dark matter scenarios suffer in general from a more severe fine tuning problem and predict a smaller neutrino flux than extragalactic topological defect model scenarios. We will therefore here focus on the latter, with $`p=1`$.

Fig. 12 shows the results for $`m_X=2\times 10^{14}`$GeV, with $`B=10^{12}`$G, and again the minimal URB consistent with data Clark ; pb . These parameters lead to optimistic neutrino fluxes for the maximal normalization consistent with all data. For more detailed recent discussions of top-down fluxes see Refs. slsc ; slby .

## V The Z-Burst Scenario with Acceleration Sources

In the Z-burst scenario UHECRs are produced by Z-bosons decaying within the distance relevant for the GZK effect. These Z-bosons are in turn produced by UHE neutrinos interacting with the relic neutrino background zburst1 . If the relic neutrinos have a mass $`m_\nu `$, Z-bosons can be resonantly produced by UHE neutrinos of energy $`E_\nu M_Z^2/(2m_\nu )4.2\times 10^{21}\mathrm{eV}(\mathrm{eV}/m_\nu )`$. The required neutrino beams could be produced as secondaries of protons accelerated in high-redshift sources. The fluxes predicted in these scenarios have recently been discussed in detail in Refs. fkr ; kkss . In Fig. 13 we show an optimistic example taken from Ref. kkss  for comparison with the other scenarios. As in Refs. fkr ; kkss  no local neutrino over-density was assumed. The sources are are assumed to not emit any $`\gamma `$rays, otherwise the Z-burst model with acceleration scenarios is ruled out, as demonstrated in Ref. kkss . We note that no known astrophysical accelerator exists that meets the requirements of the Z-burst model kkss . Also note that even exclusively emitting neutrino sources appear close to be ruled out already by the GLUE bound.

## VI Mono-energetic Super-heavy Relic Neutrino Sources

Since no known astrophysical accelerator exists that produces sufficiently strong neutrino beams up to sufficiently high energies for the Z-burst scenario to work kkss , one can speculate about more exotic sources. Gelmini and Kusenko gk  have considered a top-down type source for the Z-burst scenario, namely super-heavy particles mono-energetically decaying exclusively into neutrinos. In Fig. 14 we show predictions for one example of this kind of model with the same maximal energy $`E_{\mathrm{max}}=10^{23}\mathrm{eV}=m_X/2`$ as for case of Z-burst model. Again, no local neutrino over-density was assumed and the calculation assumed a minimal URB flux and a magnetic field $`B=10^{12}`$G. Note that, for the same UHECR flux, in this scenario the photon flux in the EGRET region is larger than in the Z-burst model with power-law neutrino sources, compare Fig. 13, because all secondary photons from Z-boson decays at redshifts $`z>3`$ contribute only to the EGRET energy range. Thus, the normalization of the photon flux to the EGRET measurements leads to a decrease of the UHE proton and photon fluxes, see Fig. 14. The EGRET measurement therefore considerably constrains the parameter space of this model. Neutrino masses $`m_\nu 0.1`$ eV are required, which allow X-particle masses $`m_X10^{14}`$GeV implying secondary photon fluxes that are below the measured level in the EGRET energy range.

Let us also note though the neutrino flux in UHECR region is much smaller in Gelmini-Kusenko model in compare with Z-burst model, the GLUE bound constrains both models in the peak of neutrino flux in similar way.

## VII Neutrino fluxes in a non-shock acceleration model for AGN

Although the recent exciting discoveries by the Chandra X-ray space observatory added much to our knowledge of structures of the jets of powerful radiogalaxies and quasars, they didn’t solve the old problems and, in fact, brought new puzzles. If the observed X-ray emission is due to synchrotron energy losses of electrons, the energy of such electrons should be of the order of 100 TeV. Electrons of such energies loose all their energy on a typical scale of only 0.1 kpc. In order to explain observed X-ray data, such 100 TeV electrons should be created uniformly over the jet length of order of 100 kpc. The conventional shock-wave acceleration mechanism in this case would require a jet uniformly filled with 1000 identical shocks, following one another along the jet axis.

In Ref. neronov  a “non-shock acceleration” version of the electron synchrotron model was proposed, namely assuming that electrons are not accelerated in the jet, but are instead the result of pair production by very high energy (VHE) $`\gamma `$rays interacting with the CMB. The typical attenuation length of $`10^{1618}`$ eV photons is of order 100 kpc, comparable with the lengths of large scale extragalactic jets.

In this model the VHE $`\gamma `$-rays are produced by accelerated protons interacting with the ambient photon fields (supplied, for example by the accretion disk around the massive black hole) through photo-meson processes. At the same time those protons produce neutrinos which are emitted in the direction of the jet. Therefore, this model predicts a high neutrino flux comparable in power with the $`\gamma `$-ray flux. The detailed numerical simulations of proton acceleration in the central engine of the AGN new  show that the collimated jet of almost mono-energetic VHE protons (linear accelerator) can be created in the electro-magnetic field around the black hole and the energy of those protons can be converted into photons and neutrinos, while protons can be captured inside of the source. The nucleon flux leaving the AGN is well below observed cosmic ray flux in this scenario. Furthermore, since all nucleons leaving the source are well below the GZK cutoff, there is no cosmogenic contribution to the neutrino flux from these sources.

Fig. 15 shows a typical prediction of the diffuse neutrino flux in this scenario. This flux is beyond the WB limit which is not applicable in this case because the sources are optically thick for nucleons with respect to pion photo-production. The flux is consistent with MPR bound for optically thick sources.

In the AGN model discussed above, blazars would be seen by neutrino telescopes as point-like sources with neutrino fluxes which are smaller or of the same order as the photon flux emitted by these same sources and which are detectable by $`\gamma `$-ray telescopes.

## VIII Conclusions

Based on our transport code we reconsidered neutrino flux predictions and especially their maxima consistent with all cosmic and $`\gamma `$ray data, for cosmogenic neutrinos produced through pion production of UHECRs during propagation, and for the more speculative Z-burst scenario and top-down scenarios. We pointed out that one can easily exceed the WB bound and, in the most optimistic cases, even the MPR bound for cosmogenic neutrinos in scenarios with cosmic ray injection spectra harder than $`E^{1.5}`$, maximal energies $`E_{\mathrm{max}}10^{22}`$eV, and redshift evolution typical for quasars, or stronger. Given our poor knowledge on the origin of UHECRs, in our opinion these are possibilities that should not be discarded at present, especially since they would lead to considerably increased prospects of ultra-high energy neutrino detection in the near future. We also show that for non-shock AGN acceleration models the AGN neutrino fluxes can reach the $`\gamma `$ray bound Eq. (2) around $`10^{16}`$eV which represents the ultimate limit for all scenarios of $`\gamma `$ray and neutrino production involving pion production.

Finally, we note that fluxes as high as for the optimistic scenarios discussed here would also lead to much stronger constraints on the neutrino-nucleon cross section at energies beyond the electroweak scale. This is important, for example, in the context of theories with extra dimensions and a fundamental gravity scale in the TeV range tev .

## Acknowledgments

We would like to thank Andrey Neronov and Igor Tkachev for fruitful discussions and comments. We thank Cyrille Barbot, Zoltan Fodor, and Andreas Ringwald for useful comments on the first version of the manuscript. We also thank David Seckel for information on the RICE experiment.


format: text

text:
# Violence in the Dark Ages

## 1 Introduction

An epoch of galaxy outflows is a virtually inevitable consequence of hierarchical structure formation. In such models, subunits merge and accrete mass to form objects of ever-increasing size with a rich substructure. Yet while the sizes of galaxies increase strongly over time, the sizes of the structures within them remain largely fixed. Thus such properties as the scale of associations of O and B stars (McKee & Ostriker (1997)) and the energy deposition due to the resulting supernovae (SNe) are likely to be weak functions of the mass of the host galaxy, depending mainly on environmental effects such as the distribution of stellar material, gas, and dust (Omukai (2000)). Hence, while the ejecta of conglomerations of Type-II SNe are primarily confined to the interstellar medium in Milky Way size galaxies, they are easily able to escape from the dwarf galaxies whose formation at high-redshift ushered in the first period of intense star formation in the universe. Thus the end of the dark ages, before galaxies formed, is likely to have been accompanied by a potentially drastic change of both the properties of the Intergalactic Medium (IGM) and the evolution of collapsing protogalactic dwarf systems (Scannapieco, Thacker & Davis 2001).

From an observational point of view, the evidence for such an epoch is equally compelling. Analyses of C$`\mathrm{IV}`$lines in QSO absorption spectra uncover an IGM that has been widely and inhomogeneously enriched by stellar material (Cowie, et al. (1995)), and heating from outflows is necessary to avoid overproduction of the X-ray background (Pen (1999)). Such mechanisms are believed to be needed to explain the X-ray luminosity-temperature relation ($`L_x:T`$) of the hot gas in galaxy clusters, which is widely discrepant with predictions that consider only heating from gravitational collapse (Kaiser (1991)). Note that a recent argument by Voit & Bryan (2001) shows that the precise details of outflow heating are comparatively unimportant to the shape of the $`L_x:T`$ relation (excessively heated gas will float to the edges of a cluster where it is undetectable, further, simulations by Pearce et al. (2000) demonstrate that cooling alone will produce a break in the $`L_x:T`$ relation). Finally, the starbursting galaxies that form the low-redshift end of this epoch have been detected directly, both in optical and infrared observations (Pettini et al. (2001)) at $`3<z<4`$, and optical observations (Frye, Broadhurst & Benitez (2002); Hu et al. (2002)) of lensed galaxies at $`4<z<6.56`$.

Motivated by this wide range of observational evidence and theoretical pointers, we have carried out a detailed numerical study of galaxy outflows, their impact on the IGM, and the formation of galaxies within it. While previously studies of outflows and cosmological metal enrichment have been attempted, these have either employed approximate semi-analytical methods (Nath & Chiba (1995); Madau, Ferrara & Rees 2001; Scannapieco & Broadhurst (2001)), low-resolution simulations (Aguirre et al. 2001a ; Aguirre et al. 2001b ; Cen & Ostriker (1999)) that are unable to identify the dwarf galaxies that are likely to be the primary source of metals (Madau et al. (2001); Scannapieco & Broadhurst (2001); Scannapieco, Ferrara & Madau 2002), or a high-resolution scheme in which supernovae explosions were not adequately accounted for (Gnedin & Ostriker (1997)). We also conducted an earlier study (Scannapieco et al. (2001)) in which we demonstrated that outflows can reduce the number density of dwarf galaxies at high redshift. This effect was seen to be a combination of ‘baryonic stripping’, whereby outflows strip gas from low over-density halos, and a delay in formation times caused by outflows slowing accretion on to the host halo. However, due to the significant wall clock time required to conduct the simulations we were unable to integrate them beyond a redshift of $`z=6`$. In this paper we build on our previous work, firstly by improving the physical model through the inclusion of variable metallicity and, secondly, we incorporate a new numerical technique to integrate to a directly observable epoch. By adopting a simple model for outflow generation and enrichment with stellar material in high-resolution cosmological simulations, we demonstrate that the of inclusion winds is able to explain a diverse set of poorly understood properties of objects ranging from moderately-sized galaxies to the absorption clouds that trace the most tenuous cosmological structures.

## 2 Simulations of Galaxy Outflows

Based on the latest measurements of the Cosmic Microwave Background and the number abundance of galaxy clusters (Balbi et al. (2000); Vianna & Liddle (1996)) we focus our attention on a Cold Dark Matter cosmological model with parameters $`h=0.65`$, $`\mathrm{\Omega }_0`$ = 0.35, $`\mathrm{\Omega }_\mathrm{\Lambda }`$ = 0.65, $`\mathrm{\Omega }_b=0.05`$, $`\sigma _8=0.87`$, and $`n=1`$, where $`h`$ is the Hubble constant in units of 100 km s-1 Mpc-1, $`\mathrm{\Omega }_0`$, $`\mathrm{\Omega }_\mathrm{\Lambda }`$, and $`\mathrm{\Omega }_b`$ are the total matter, vacuum, and baryonic densities in units of the critical density, $`\sigma _8^2`$ is the variance of linear fluctuations on the $`8h^1\mathrm{Mpc}`$ scale, and $`n`$ is the “tilt” of the primordial power spectrum.

As in our earlier work (Scannapieco et al. (2001)), simulations were conducted with a parallel OpenMP based implementation of the “HYDRA” code (Thacker & Couchman (2000)) that uses the Adaptive Particle-Particle, Particle-Mesh algorithm (Couchman (1991)) to calculate gravitational forces, and the Smoothed Particle Hydrodynamic (SPH) method (Lucy (1977); Gingold & Monaghan (1977)) to calculate gas forces. We simulate a periodic box of size 5.2 $`h^1`$ comoving Mpc containing $`192^3`$ dark matter and $`192^3`$ gas particles with masses of $`2.5\times 10^6M_{}`$ and $`5.0\times 10^5M_{}`$ respectively. This is sufficient to resolve the smallest objects that are able to cool efficiently without molecular hydrogen, with masses $`>10^8M_{}`$. We chose this limit as primordial molecular hydrogen is extremely fragile, and is likely to be quickly dissociated by the first stars that form at very high redshift (Haiman, Rees & Loeb (1997)). Note that while alternative schemes for creating cosmological H2 have been suggested, star formation and outflow generation is likely to be inefficient in objects that rely on molecular cooling, causing them to have a small overall impact on IGM heating and enrichment (Madau et al. (2001); Scannapieco et al. (2002)). Forces are softened using a fixed physical Plummer softening length of 1.54 kpc, yielding a minimum hydrodynamic resolution of 1.8 kpc; gas densities and energies are calculated using the standard SPH smoothing kernel method (Scannapieco et al. (2001)), with the kernel tuned to smooth over 52 particles; and radiative cooling is calculated using standard tables (Sutherland & Dopita (1993)). Because the epoch of reionization is unknown, and also for reasons of simplicity, we do not include a fiducial photoionization background in the simulation.

Under heavy particle clustering, and when implemented with a minimum resolution limit, the SPH method can degenerate to an order $`N^2`$ algorithm. While some authors resort to merging particles to avoid this problem, this breaks the one-to-one relationship between mass and particle number and leads to a degradation of accuracy in the hydrodynamic calculation. Our solution to this problem (Thacker (2002)) is to calculate an aggregate solution in regions below the resolution limit, grouping particles in sets of 4-30 that may change from step to step, while never actually merging particles together. The errors introduced by this method are far smaller than those created by merging particles, and the algorithm remains less than order $`N^{3/2}`$ under heavy clustering, which translates to a three-fold performance improvement over our previous code (Scannapieco et al. (2001)).

To study the propagation of galaxy outflows, it is first necessary to identify where galaxies are forming. Typically, this is done using the friends-of-friends group finding algorithm which relies upon calculating inter-neighbor distances (Efstathiou et al. (1985)), however this technique is very computationally expensive. We therefore adopt an alternative approach in which we search for regions where gas particles exceed a density threshold $`\delta _c`$ times the mean gas density. Provided more than $`2.6\times 10^7M_{}`$ (our effective mass resolution limit) is contained in such a region, the object is counted as a galaxy. This method performs well in comparison to the friends-of-friends technique (Scannapieco et al. (2001)), and throughout this paper we set $`\delta _c=500.`$

Star formation is modeled by converting a fixed fraction, $`ϵ_{\mathrm{sf}}=0.1`$, of the gas in a collapsed object into stars in a single starburst. Subsequent starbursts within a galaxy are likely to be smaller and less efficient, so our model is a conservative lower limit on star formation activity. The gas involved in a starburst is tagged to prevent further bursts until a sufficiently large reservoir of ‘fresh’ gas is accumulated. Star formation can only occur when over one-third of the galaxy, and at least 52 new particles, are untagged, although testing has shown that these thresholds can be varied over a wide range with only a small effect on our results (Scannapieco et al. (2001)). Finally, we also impose an upper limit of 1000 particles per burst, corresponding to $`5\times 10^8M_{}`$ of stars, which prevents unphysically large groups of supernovae from exploding simultaneously and generating enormous outflows (Ferrara, Pettini & Shchekinov (2000)).

The structure of galaxy outflows has been studied in detail previously (Mac Low & Ferrara (199); Martel & Shapiro (2000); Mori, Ferrara & Madau 2001), yielding the conclusion that in starbursting disk galaxies with masses $`<10^{10}M_{}`$ “blow-out” occurs. In these galaxies, the superbubbles around groups of Type-II SNe punch out of the plane of the disk and shock the surrounding IGM although they fail to excavate the interstellar medium itself. Guided by this scenario, we model the sub-resolution physics of outflow generation by rearranging the IGM surrounding a galaxy while leaving the central object intact.

Hence we construct our outflows from gas particles within a radius twice that of the central galaxy ($`r2r_g`$), but outside the galaxy itself ($`r>r_g`$). This local IGM is arranged randomly into two concentric spherical shells of radius $`2.0r_g`$ and $`1.8r_g`$ of equal mass. This multi-shell structure assures that the outflows will be sufficiently well-resolved radially to be reasonably treated by the SPH algorithm. Sufficient resolution in each of the shells is assured by arranging the particles in an anti-correlated fashion such that all particles have an inter-neighbor distance in excess of half the average inter-particle spacing. We have tested this approach in detail against one-dimensional analytical models of expanding shells, and have found that it reproduces these results almost exactly (Scannapieco et al. (2001)).

The particles in these shells are assigned velocities that conserve their overall center of mass velocity and angular momentum, and include a velocity boost in the radial direction. The strength of this boost is determined by assuming that a fixed fraction of the mechanical energy from supernovae from each starburst is channeled into the outflow, $`ϵ_{\mathrm{wind}}=0.1`$. We estimate this kinetic input to be $`2\times 10^{51}`$ ergs per supernova, to take into account the contribution from stellar winds, and assume that one supernova occurs for every 100 $`M_{}`$ of stars formed (Gibson (1997)).

Our previous simulations side-stepped the issue of metal enrichment, taking a fixed value of $`0.05Z_{}`$ for all gas particles. In this work we adopt a more realistic model in which the primordial metallicity is taken to be zero, and gas particles are enriched by ejecta from Type-II SNe. The metal content of each gas particle is therefore allowed to evolve in a similar fashion to hydrodynamic variables. We assume a standard yield of 2 $`M_{}`$ of metals per SNe, corresponding to an effective metal yield $`p=0.02`$, and distribute this material equally between the galaxy and the outflow. For these parameters a galaxy will be enriched following the first outflow event to a metallicity $`Z_{gal}`$ given by, $`\mathrm{log}_{10}(Z_{gal}/Z_{})\mathrm{log}_{10}(ϵ_{\mathrm{sf}}\times f\times p/Z_{})=1.3`$, if we take $`f=0.5`$ as the fraction of metals retained by the galaxy. Following an outflow event, the total metal yield, by mass, for the galaxy is distributed evenly among gas particles that determined to be in the galaxy ($`\delta M_{particle}^Z=(1/N_{gal})fM_{total}^Z`$, where $`M_{total}^Z`$ is the total metal yield, by mass, of the starburst, $`N_{gal}`$ is the number of particles in the galaxy and $`\delta M_{particle}^Z`$ is the increase in the metal mass of a particle in the galaxy). Similarly, the metal yield deposited in the outflow is distributed evenly among gas particles determined to be in the outflow ($`\delta M_{particle}^Z=(1/N_{outflow})(1f)M_{total}^Z`$). Individual gas particles can undergo multiple episodes of enrichment in response to a series of outflows.

Enriched particles retain their metals for the remainder of the simulation, and the spatial metal distribution is calculated by convolving with the SPH smoothing kernel. The smoothed field is then used self-consistently in computing gas cooling. This may be viewed as assuming the metals within a given particle’s smoothing length are well mixed, a reasonable approximation given that the Rayleigh-Taylor instability will mix gases in a time roughly proportional to the size of the region divided by the sound speed.

## 3 Outflows and Structure Formation

We have conducted two simulations: a fiducial run in which we fix the star-formation and outflow efficiencies to be $`ϵ_{\mathrm{sf}}=ϵ_{\mathrm{wind}}=0.1`$, and a comparison “no-outflows run” in which we fix $`ϵ_{\mathrm{sf}}=0.1`$ but do not include metal enrichment or outflows. In the fiducial run the gas is initially free of metals, while in the no-outflows run the metallicity is fixed at $`0.05Z_{}`$ at all times. While particle-particle interactions halted our previous simulations at $`z=6`$, the pseudo-merging scheme discussed above allowed us to reach $`z=4`$ in this study. A total of 4515 outflows were generated in this case, with initial radial velocities ranging from 80 to 300 km/s, consistent with direct observations (Pettini et al. (2001); Frye, Broadhurst & Benitez (2002)).

The cumulative number density of galaxies within the simulations as a function of baryonic mass and redshift is shown in Figure 1. Below a limiting mass of $`5\times 10^9M_{}`$, two significant features are observed. Firstly, there is an overall reduction in the baryonic mass of the objects in the outflows run. This is caused by outflows suppressing accretion onto the host galaxy, a process that is similar to that described in earlier work (Dekel & Silk (1986)) on dwarf galaxies which showed that supernova-driven winds result in a systematic bias toward higher mass-to-light ratios. Secondly, the fiducial run also contains many examples of the “baryonic stripping” mechanism previously identified (Scannapieco, Ferrara & Broadhurst (2000)). In these cases, the momentum in an outflow is sufficient to evacuate gas from nearby overdense regions that would have otherwise later formed into dwarf galaxies, leaving behind empty ‘dark halos’ of cold dark matter. Note also that this overall suppression of objects points to a possible solution to the discrepancy between the number of local group dwarf galaxies and predictions from CDM models that do not include outflows (Klypin et al. (1999)).

Above $`>5\times 10^9M_{}`$, the number-densities in the two models become similar, which was not observed in our previous work (Scannapieco et al. (2001)). Thus this feature may be related to the key difference between these simulations: the inclusion of gas enrichment by supernovae. To study this further, we plot the mean metallicity of the gas particles identified in galaxies as a function of galaxy mass in Figure 2A. Here we see that above $`10^9M_{}`$, the metallicity is roughly constant at $`0.1Z_{}`$. While our models do not include the quiescent mode of star formation that is important in such large objects, this value is suggestive of the $`0.1Z_{}`$ pre-enrichment necessary to explain the lack of low-metallicity G dwarf stars in the solar neighborhood (Ostriker & Thuan (1975)), as well as in other large nearby galaxies (Thomas, Greggio & Bender (1999)), a result that is consistent with semi-analytic studies (Scannapieco & Broadhurst (2001)).

Using a simple estimate for the “cooling time” (see Scannapieco & Broadhurst (2001)) we can get a rough idea of the delay between gas collapse and cooling into a disk for the large objects in these simulations, which is a strong function of metallicity. Referring to Figure 1, we see that above $`2\times 10^{10}M_{}\mathrm{\Omega }_0/\mathrm{\Omega }_b`$ the enrichment run appears to contain more massive galaxies than the non-enrichment run, although there are only two halos in the simulation that are above this threshold and we would prefer better statistics before making definitive conclusions. For the sake of reference we call this mass threshold the ‘crossover’ mass. We note that there was no evidence for a crossover phenomenon in our previous simulation work (Scannapieco et al. 2001).

For galaxies with total masses roughly corresponding to the crossover value the cooling time is $`3.0\times 10^8`$ years for $`0.05Z_{}`$ metallicity gas and $`1.5\times 10^8`$ years for gas at $`0.1Z_{}`$ metallicity. Comparing these values with the age of the universe in this cosmology ($`11.5\times 10^8`$ years at $`z=5`$ and $`14.8\times 10^8`$ years at $`z=4`$) shows that approximately $`3.3\times 10^8`$ years passes between $`z=4`$ and $`z=5`$, which is very slightly larger than the cooling time for the $`0.05Z_{}`$ gas, and about twice the cooling time of the $`0.1Z_{}`$ gas. Both of these cooling times are significantly shorter than the Hubble time at $`z=4`$ which implies that the increased cooling rate can only have an effect on halo gas that is close to the galaxy. Given the small mass fraction of this gas it seems somewhat unlikely that it would contribute to such a noticeable mass increase. An alternative explanation may be that massive halos act as a sink for outflow material from surrounding smaller systems. However, if this were the case we would expect to have seen the effect of this in Scannapieco et al. (2001). Ultimately, an exact determination of the nature of the apparent increase in mass requires a larger simulation box to improve statistics.

Our outflow model also has important implications for galaxies with baryonic masses less than $`10^9M_{}`$, whose metallicities range from $`0.02Z_{}`$ to $`0.1Z_{}`$ in our simulations. This distribution can be compared with metallicities observed in local group dwarf galaxies, which are thought to form at high redshift. To facilitate this comparison, we overplot in Figure 2A the masses and metallicities of these objects as compiled in a review of the local group galaxies (Mateo (1998)), using the iron abundances when available, and assuming a “shifted” solar abundance of $`10^{9.4}`$ in cases in which only the oxygen abundance was measured. Note that this standard shift is also taken in the observational literature (Mateo (1998)) because these objects are rich in $`\alpha `$ elements such as oxygen, which are abundant in the type-II SNe such as those driving our outflows. In each case we assume a standard error in the mass determination of $`50\%`$, due to errors in determining the velocity dispersion and distance to these objects, and apply an overall mass shift of $`\mathrm{\Omega }_b/\mathrm{\Omega }_0.`$ Clearly our model does a good job of reproducing the metallicity-luminosity relation seen in dwarf galaxies, in which metallicity increases as a function of luminosity, but with a scatter of almost an order of magnitude.

Compared to the enrichment values presented in Gnedin & Ostriker (1997) our values are lower, by an approximate factor of 3-4, in the highest density regions. Most of this difference is probably attributable to the difference between the star formation rates, as the inclusion of outflows systematically lowers our SFR (although we lack a photoionization background). We also have significantly less scatter in the observed metal fractions, which is a function of differences between the numerical methods employed. Gnedin & Ostriker adapt the method presented in Cen & Ostriker (1992) which involves converting mass in cells obeying overdensity, cooling and convergence criteria to collisionless stellar particles. This can lead to a widely varying enrichment history depending upon the details of the density of the gas in the protogalaxy. In our model we ignore the inner structure of our halos and instead treat them as uniform systems that undergo global starbursts when sufficient gas has been accumulated, as is done in semi-analytic models. In passing, we mention that our results for protogalaxy metallicity show good agreement with averaged enrichment values presented in Cen & Ostriker (1999), provided that a log-linear extrapolation to z=4 is reasonable.

Despite the increased formation of large objects in the outflows simulation, their comparative rarity at these redshifts causes the total amount of cooled gas to be dominated by the suppression effect. In Figure 2B we plot the overall star formation rate (SFR) as a function of redshift in each of these simulations. The plot shows that the inclusion of outflows results in a decrease in the SFR by approximately a factor of two. We also compare the calculated SFRs with direct observations of star formation (Steidel et al. (1999)) and a broad range of indirect constraints (Hopkins et al. (2001)). While both models are consistent with the indirect constraints, the direct observations are slightly over-predicted by the no-outflows run and slightly under-predicted in the fiducial run. As our starburst model does not include the quiescent mode of star formation, our results should be considered lower limits on the global SFR. Thus, although the data seems to favor the fiducial run, no definitive conclusions can be drawn.

It is worth noting that the star formation rate only becomes significant at redshifts $`<10`$. This redshift is determined by the minimum mass scale in our simulations, chosen to be the minimum mass that can cool effectively without molecular hydrogen, which we assume to be completely dissociated.

We have not included a photoionization background in our simulation primarily because we wanted to observe the effect of the outflow model without having to disentangle it from other physical phenomena. Nonetheless, it is important to assess what effect the inclusion of photoionization field might have. Photoionization is known to reduce the number density of dwarf systems, and is considered to be the primary solution to the ‘substructure problem’ (Klypin et al. (1999); Moore et al. (1999)). Most simulations of photoionization have concentrated on either large scales with a comparatively low mass resolution (eg Katz et al. (1996)) or very small scales at high mass resolution to examine phenomenology of individual halos (eg Navarro & Steinmetz (1997)). On the other hand semi-analytic studies are able to provide both good statistics and high mass resolution (eg Sommerville (2002)), at the cost of lost geometric information and highly phenomenological physics.

Our minimum halo mass resolution is slightly under $`1.5\times 10^8M_{}`$. At this mass scale the $`z=0`$ V-band luminosity function for the semi-analytic model of Sommerville (2002) is about a factor of three times lower when a photoionization background is included. This result is surprisingly insensitive to the actual photoionization epoch, and large differences between different photoionizing models are only found at mass scales that we do not resolve (around $`10^7M_{}`$, i.e. close to the mass of the Draco dwarf spheroidal). Thus, assuming that not all dwarf systems are late forming, we can expect our simulation to over-predict the number of dwarf halos by at most a factor of three. The natural effect of a reduced number density of gaseous halos is a reduced global SFR and reduced enrichment volume factor, although the metallicity values for individual galaxies will be similar. However, in all likelihood, some of the protogalaxies that are destroyed by outflows would be suppressed, or possibly erased, by a photoionization background. Thus an accurate assessment of the effect of a photoionization background on the outflow model is quite difficult. Also if the epoch of reionization is optimistically late ($`z7`$), then we can expect a less significant impact on our model, as a large amount of the evolution will occur without the presence of a photoionizing background. In this case the SFR at $`z=7`$ will undoubtedly plateau, or turn off somewhat, and the enrichment fractions will be reduced. Ultimately, determining the precise impact of a photoionization background requires running a self consistent simulation that includes it.

## 4 Outflows and the Intergalactic Medium

The widespread generation of outflows has an impact over regions much larger than the dense peaks identified as galaxies. Indeed the strongest inference of an epoch of galaxy outflows is drawn from observations of the IGM itself.

There are two methods by which this impact can be quantified. The first of these is by examining the overall mass of gas enriched by stellar material, which we plot as a function of redshift in Figure 3A. While outflows vigorously enrich the densest regions, the IGM transition associated with outflows is incomplete, and only $`20\%`$ of the particles are enriched by redshift 4. In this plot we also see that the majority of enrichment takes place at redshifts below 10, when atomic transitions become an effective coolant.

The onset of this enrichment roughly corresponds with the sharp rise in the number of outflows at $`z>10`$, as shown in Figure 3B. This plot, which displays the outflow rate per unit Mpc3 per Gyr, also hints at the end of the epoch of galaxy outflows as the curve begins to turn over at lower redshift values. The details of this turn-over are quite uncertain, however, and model dependent. The determining factor in this case is the maximum size of the gravitational potential from which SNe can successfully cooperate to drive a large-scale outflow, a quantity that is crudely modeled by our overall upper limit of 1000 particles. In reality this turn-over is dependent on such uncertain quantities as the distribution of high-redshift OB associations, the star formation efficiency, and the stellar initial mass function in high-redshift galaxies (Mori et al. (2001)).

In Figure 3C we examine the overall level of enrichment within the affected regions, and also plot the bounds from analytic models (Pei, Fall & Hauser (1999)) for reference. The slight decrease in these values seen at large redshifts is likely to be due to low-number statistics, as only $`15`$ objects are formed in this simulation by $`z=15`$. Despite the fact that outflows are able to enrich the largest objects to metallicities of $`0.1Z_{}`$, their mean impact on polluted regions of the IGM is only on the order of $`0.02Z_{}`$, and the mean value when averaged over the simulation as a whole is $`0.003Z_{}`$. Thus our model not only results in a patchy distribution of enriched regions, but also a large variance within these regions. This is consistent with observations of C$`\mathrm{IV}`$at $`z3`$ in QSO absorption systems, which can be reproduced by a model in which the overall mean metallicity of the IGM is $`0.003Z_{}`$ with a large scatter of roughly an order of magnitude (Hellsten et al. (1997)).

To get a better handle on this inhomogeneity, we turn our attention to the second method of quantifying this IGM transition and examine the volume impacted by outflows. Figure 4 is a volume-rendered image of the smoothed distribution of metals in our simulation volume. Here we see that the spatial distribution of material from outflows is highly biased, illustrating the necessity of studying this problem through simulations or other techniques that capture this inhomogeneity (Scannapieco & Broadhurst (2001)). While much of the simulation remains metal free, large clouds of $`0.02Z_{}`$ material roughly outline the underlying filaments of dark matter, punctuated by smaller regions of higher metallicity.

In Figure 3A we plot the redshift evolution of the volume enriched to 0.001$`Z_{}`$ or greater. The volume fraction enriched follows the trend of the mass enriched but is reduced because of the large overdensities at which galaxies form. To quantify this further, in Figure 5 we plot the volume filling factor of material with metallicity above a given threshold value, plotted as a function of the threshold value at $`z=4`$. At this redshift the maximum enrichment volume is slightly over 20%, consistent with new semi-analytic calculations (Scannapieco et al. (2002)), and apparently at odds with the claimed (Madau et al. (2001)) 100% value for $`z=10`$.

While Madau, Ferrara, & Rees (2001) only considered outflows from objects at a fixed mass scale of of $`10^8M_{}/h`$ and then multiplied by a Press-Schechter estimate of their number density to compute to overall filling factor, the overall energetics of their wind models are quite similar to ours. Thus it is likely that the primary source of this inconsistency is related instead to their adoption of an sCDM cosmology with parameters $`\mathrm{\Omega }_m=1`$, $`h=0.5`$, $`\sigma _8=0.63`$, $`n=1`$, $`\mathrm{\Omega }_bh^2=0.019`$. Indeed, using the standard Press-Schechter approach to calculate the comoving number density of halos of masses greater than $`10^8M_{}/h`$ at $`z=9`$ we find that the sCDM model of Madau, Ferrara & Rees (2001) contains 5.5 times as many objects per Mpc h-1 as our $`\mathrm{\Lambda }`$CDM model. This enhancement is primarily due to the association of a smaller top-hat filtering scale to sCDM objects as compared with their $`\mathrm{\Omega }_M=0.35`$ counterparts, and is sufficient to explain the factor of 5 difference that we observe between our results. We therefore conclude that the difference between filling factors can be attributed almost entirely to the reduced number counts of dwarf galaxies in the observationally motivated $`\mathrm{\Lambda }`$CDM cosmology compared to the out of favor sCDM model.

Nonetheless, Figure 5 clearly shows that enrichment is far from uniform. Metallicities greater than $`0.01Z_{}`$ are confined to less than 5.5% of the volume, while enrichments in excess of $`0.1Z_{}`$ are confined primarily to larger dwarf galaxies and their halos, which occupy less than 0.1% of the simulation volume. This correlation between the distribution of metals and starbursting galaxies provides a natural explanation for a long standing dichotomy seen in QSO absorption line studies. For while it is clear that the distribution of Ly$`\alpha `$ in the IGM as a whole is well modeled without considering outflows (McDonald et al. (2000)), studies of C$`\mathrm{IV}`$clouds indicate that these enriched objects have been turbulently stirred by events occurring on timescales on the order of $`10^8`$ years (Rauch et al. (2000)).

## 5 Conclusion

In this work, we have studied the role of outflows and metal enrichment in structure formation, comparing our results to a diverse set of observations. We have not attempted to calculate the properties of individual galaxies, but rather have used results from previous studies as input into our model. While based on the approach described in Scannapieco, Thacker, & Davis (2001), this study is a significant advance over our previous work as it includes a more realistic model of metal enrichment as well as a more efficient SPH algorithm that allows us to integrate to directly observable redshifts. This leads to a number of important conclusions, which merit further investigation.

First, while the presence of outflows greatly reduces the number of galaxies formed at lower mass scales as a result of baryonic stripping, this mechanism is effective only for galaxies with total mass below $`5\times 10^9`$ M. At greater mass scales, the galaxy number counts in the outflows case become similar or even exceed that of the no-outflows case as a result of more rapid cooling caused by metal enrichment.

Second, the inhomogeneous distribution of metals by outflows leads naturally to a scatter in metallicity of lower mass objects that is consistent with observations of local group dwarf galaxies. Above $`10^9`$ M, however, the metallicity is roughly constant at 0.1 Z in agreement with pre-enrichment estimates necessary to explain the lack of low-metallicity G-dwarf stars in the solar neighborhood.

Third, the global SFR for our model is in broad agreement with observations, albeit on the low side. This is attributed to the lack of quiescent star formation in our model.

Finally, our global enrichment value of 0.3% Z is consistent with observations of C$`\mathrm{IV}`$at z$``$3 in QSO absorption systems This metal enrichment is far from uniform and concurs with the interpretation that C$`\mathrm{IV}`$clouds are turbulently stirred on small scales while the distribution of Ly$`\alpha `$ clouds remains comparatively unaffected.

Future work needs to address the role of H2 cooling in our simulations as at the moment we make no attempt to calculate production of H2. Progress is difficult on this front as while calculation of rate equations is comparatively straightforward (and has been conducted within in our simulation framework, see Fuller & Couchman (2000)) a fully consistent solution with local photodisassociation is, at least at the moment, an extreme computational challenge (eg Razoumov & Scott (1999)). We also plan to include the fiducial photoionization background associated with the epoch of reionization. We will address both of these issues over the coming months, although we do not anticipate a major revision of our current results.

We thank Andrea Ferrara, Tom Broadhurst and Mordecai Mac Low for helpful discussions, Greg Bryan for making available unpublished results and the anonymous referee for suggesting a number of changes that significantly improved the content of the paper. We also gratefully acknowledge the hospitality of Carlos Frenk and the Physics Department at the University of Durham where this paper was finished. E.S. has been funded in part by an NSF MPS-DRF fellowship. R.J.T. acknowledges funding from the Canadian Computational Cosmology Consortium and NSERC of Canada via the operating grant of Prof H. M. P. Couchman. This project was supported by NSF KDI grant 9872979.


text:
# Nuclear modifications of fragmentation functions and rescaling models Talk presented at nthe “European workshop on the QCD structure of the nucleon” (QCD-N’02), Ferrara (ITA), April 3rd-6th, 2002. This work is (partially) funded by the European Commission IHP program under contract HPRN-CT-2000-00130.

## Abstract

We discuss nuclear modification of fragmentation functions in the context of the so-called “rescaling models”. These models implement partial deconfinement inside nuclei by modifying the fragmentation functions perturbatively. We apply these models to the analysis of nuclear hadron production in deep inelastic scattering processes at the HERMES and EMC experiments.

In deep inelastic scattering a projectile lepton $`\mathrm{}`$ emits a virtual photon $`\gamma ^{}`$, which scatters on a quark $`q`$ from the target $`T`$, in our case a deuteron $`D`$ or a heavier nucleus with atomic number $`A`$. The struck quark fragments into the observed hadron $`h`$, see Fig. 1. We use kinematic variables as summarized in the same figure: $`x`$ is the Bjorken’s scaling variable, $`\nu `$ is the energy of the virtual photon in the target rest frame, and $`z`$ is the fraction of the energy transferred to the produced hadron.

The experimental data on nuclear effects in hadron production are usually presented in terms of the multiplicity ratio

$`R_M^h(z)={\displaystyle \frac{1}{N_A^{\mathrm{}}}}{\displaystyle \frac{dN_A^h}{dz}}/{\displaystyle \frac{1}{N_D^{\mathrm{}}}}{\displaystyle \frac{dN_D^h}{dz}},`$ (1)
where $`N_A^{\mathrm{}}`$ is the number of outgoing leptons in DIS processes on a nuclear target of atomic number $`A`$ and $`dN_A^h/dz`$ is the $`z`$-distribution of produced hadrons in the same processes; the subscript $`D`$ refers to the same quantities when the target is a deuteron. A similar definition for the multiplicity ratio as a function of $`\nu `$ may be written. In absence of nuclear effects the ratio $`R_M^h`$ would be equal to 1. The experimental observation that $`R_M^h1`$  has been explained theoretically in many ways: as an effect of nuclear absorption of the produced hadrons , in a gluon-bremsstrahlung model for leading hadron production , as an higher-twist effect . In this short note we will not discuss these models in detail (see Ref. ) and, instead, will concentrate on the so-called rescaling models .

The starting point of rescaling models is to assume a change in the confinement scale $`\lambda _A`$ in nuclei, compared to the confinement scale $`\lambda _0`$ in free nucleons:

$`\lambda _A>\lambda _0.`$
This assumed partial deconfinement in nuclei affects both the parton distribution functions (PDF) and the fragmentation functions (FF). For what concerns PDF’s , consider a valence quark which carries a momentum $`Q_0`$ when it is confined on a scale $`\lambda _0`$. If the scale changes to $`\lambda _A`$ it carries a corresponding momentum $`Q_A`$. Since there is no other dimensionful scale, the product $`Q\lambda `$ must remain constant, so that

$`Q_0\lambda _0=Q_A\lambda _A.`$ (2)
Therefore, if we take $`Q_0`$ to be the initial scale for the DGLAP evolution of distribution functions, we may set

$`q_f^A\left(x,Q_A={\displaystyle \frac{\lambda _0}{\lambda _A}}Q_0\right)=q_f(x,Q_0),`$ (3)
where $`q_f`$ is the distribution function of a quark of flavour $`f`$ in a free nucleon and $`q_f^A`$ is the same quantity when the nucleon is inside a nucleus. For FF’s  a similar argument holds because constituent quark and effective hadron masses are sensitive to the confinement scale. Therefore if we take $`Q_0`$ to be the physical threshold for hadron production we may set

$`D_f^{h|A}\left(x,Q_A={\displaystyle \frac{\lambda _0}{\lambda _A}}Q_0\right)=D_f^h(x,Q_0),`$ (4)
where $`D_f^h`$ is the fragmentation function of a quark of flavour $`f`$ into a hadron $`h`$ and $`D_f^{h|A}`$ is the nuclear modified fragmentation function. To extend Eqs. (3) and (4) to an arbitrary scale $`Q`$ we apply pQCD evolution. The nuclear structure and fragmentation functions evolve over larger interval in momentum compared to the corresponding functions at the same scale Q, since the starting scale is smaller, see Eq. (2). The final result  is

$`q_f^A(x,Q)`$ $`=q_f(x,\xi _A(Q)Q)`$ (5)
$`D_f^{h|A}(z,Q)`$ $`=D_f^h(x,\xi _A(Q)Q),`$ (6)
where the scale factor $`\xi _A(Q)`$ is greater than one. These two equations are the main tool which we will exploit. There are two models for the computation of the scale factor:

1. the maximal deconfinement model (NP)  assumes the onset of “colour conductivity” in nuclei and takes $`\lambda _A=R_A`$, where $`R_A`$ is the nuclear radius, so that

$$\xi _A(Q)=R_A/R_p;$$
2. the partial deconfinement model (CJRR)  assumes the deconfinement scale $`\lambda _A`$ to be proportional to the degree of overlap of the nucleons inside the given nucleus, and

$$\xi _A(Q)=\left(\lambda _A/\lambda _0\right)^{\frac{1}{2}\frac{\alpha _s(Q_0)}{\alpha _s(Q)}}.$$

Note that the maximal deconfinement model assumes a much larger scale factor, which results in a larger nuclear modification of hadron production at high $`z`$.

We calculate the multiplicity ratio (1) using the rescaled PDF (5) and the rescaled FF (6) in the leading order pQCD computation of $`N_A^{\mathrm{}}`$ and $`dN_A^h/dz`$:

$`{\displaystyle \frac{1}{N_A^{\mathrm{}}}}{\displaystyle \frac{dN_A^h}{dz}}`$ $`={\displaystyle \frac{1}{\sigma ^{\gamma ^{}A}}}{\displaystyle _{\mathrm{exp}.\mathrm{cuts}}}𝑑x𝑑\nu {\displaystyle \underset{f}{}}e_f^2q_f(x,\xi _AQ){\displaystyle \frac{d\sigma ^{\gamma ^{}q}}{dxd\nu }}D_f^h(z,\xi _AQ)`$ (7)
$`\sigma ^{\gamma ^{}A}`$ $`={\displaystyle _{\mathrm{exp}.\mathrm{cuts}}}𝑑x𝑑\nu {\displaystyle \underset{f}{}}e_f^2q_f(x,\xi _AQ){\displaystyle \frac{d\sigma ^{\gamma ^{}q}}{dxd\nu }},`$
where $`e_f`$ is the electric charge of a quark of flavour $`f`$, $`d\sigma ^{\gamma ^{}q}/dxd\nu `$ is the differential cross-section for a $`\gamma ^{}q`$ scattering computed in pQCD at leading order, and $`Q^2=2Mx\nu `$, with $`M`$ the nucleon mass. In the numerical computations we used GRV98 parton distribution funtions  and Kretzer’s parametrization of FF’s at leading order .

The results for the $`\nu `$-dependence of the multiplicity ratio for charged hadrons are shown in Fig. 2a and compared to EMC data . The maximal deconfinement model underestimates the data over nearly the whole range. On the contrary the partial deconfinement model fit the data at high $`\nu `$, but overestimates them at smaller $`\nu `$. Because of Lorentz dilatation, at high $`\nu `$ the hadron is formed mainly outside the nucleus and it is affected by rescaling effects only. On the contrary, at small $`\nu `$ the hadron is formed inside the nucleus and starts interacting with its nucleons, with some probability of being absorbed.

To take into account nuclear absorption, we follow the analysis of Ref. : the quark propagates in the nucleus with a cross-section $`\sigma _q`$ for inelastic scattering on the nucleons. Subsequently, it creates a pre-hadronic states, which has a cross-section $`\sigma _{}`$, and finally the observed hadron is formed, which has a cross-section $`\sigma _h`$. Following the experimental indications of EMC and HERMES, we take $`\sigma _q=0`$ and $`\sigma _{}=\sigma _h=20`$ mbarn for charged hadrons. With the latter assumption we need to consider effectively only the formation time $`\tau _F`$ for the pre-hadronic state :

$`\tau _F=\left({\displaystyle \frac{\mathrm{ln}(z^2)1+z^2}{1z^2}}\right){\displaystyle \frac{z\nu }{\kappa }},`$
where $`\kappa `$ is a parameter which we fix to $`\kappa =0.25`$ GeV/fm. Note that $`\tau _F\frac{1}{\kappa }(1z)\nu `$ as $`z1`$, giving the formation time suggested by the gluon bremsstrahlung model of Ref. . Finally, nuclear absorption effects are included in the computations by multiplying the integrand in (7) by the nuclear absorption factor $`𝒩_A`$, which represent the fraction of hadrons which escape from the nucleus:

$`𝒩_A={\displaystyle d^2b_{\mathrm{}}^{\mathrm{}}𝑑y\rho _A(b,y)\left[S_A(b,y)\right]^{A1}},`$
where $`\rho _A(b,y)`$ is the nuclear density normalized to 1 at transverse and longitudinal coordinates $`(b,y)`$, and $`S_A(b,y)`$ is the survival probability of a hadron produced in $`(b,y)`$:

$`S_A(b,y)=1\sigma _h{\displaystyle _y^{\mathrm{}}}𝑑y^{}\rho _A(b,y^{})\left(1e^{(y^{}y)/\tau _F}\right).`$
For the deuteron we used as a density the sum of the Reid’s soft-core S- and D-wave functions squared . For heavier nuclei we used a Woods-Saxon density with radius $`R_A=1.12A^{1/3}0.86A^{1/3}`$ fm.

The $`\nu `$\- and $`z`$-dependence of the multiplicity ratio at EMC and HERMES after the inclusion of nuclear absorption are shown in Fig. 2b and Fig. 2c. It is clearly seen that while the CJRR model gives a nice description of the data, the NP model is ruled out.

The average values of $`\nu `$ available at HERMES are smaller than at EMC and the absorption effects larger. In Fig. 3 we compare the $`z`$-distributions in the partial deconfinement models at EMC and HERMES before and after inclusion of absorption. As expected, at EMC with a 63Cu target nuclear absorption is marginal and rescaling alone gives a satisfactory description of the data. At HERMES with a 84Kr target, which is comparable to copper, both effects are larger and absorption is dominant, tending to mask rescaling effects. With a 14N target, which is smaller than krypton, both effects are smaller.

In summary, rescaling models (supplemented by nuclear absorption) are shown to be able to describe both HERMES and EMC data on the nuclear modification of hadron production in DIS processes. While the maximal deconfinement model is ruled out by the data - it assumes a too large deconfinement - the partial deconfinement model is shown to be a good one. Further precise data at moderate and high $`\nu `$’s and for light and heavy targets are needed to disentangle rescaling and formation time effects.

Acknowledgments

We are grateful to N. Bianchi, P. di Nezza, B. Kopeliovich and V. Muccifora for stimulating discussions and to R. Fabbri for technical support during the conference.


text:
# Measuring CMB Polarization with ESA PLANCK SubMM-Wave Telescope

## 1 Introduction

The dual-reflector submillimeter-wave telescope on the ESA PLANCK Surveyor is being designed for measuring the temperature anisotropies and polarization characteristics of the cosmic microwave background (CMB).

Our research is concerned with one of its focal plane instruments, Far-IR High Frequency Instrument (HFI) , which will cover six frequency bands centered at $`100`$, $`143`$,$`217`$, $`353`$, $`545`$ and $`857GHz`$ providing the sensitivity of $`\mathrm{\Delta }T/T10^6`$ and the angular resolution down to $`5`$ arcminutes. The HFI consists of an array of $`36`$ horn antenna structures (Fig. 1, a) feeding the bolometric detectors which will be cryogenically cooled to a temperature of $`100mK`$.

One objective of the research is the optimization of the HFI optical design and the computation of the HFI beam patterns. The challenge of the problem is that the telescope is electrically large ($`D/\lambda =4300`$ at $`\lambda =350\mu m`$) and consists of two essentially defocused ellipsoidal reflectors providing a very large field of view at the focal plane.

Another objective is the characterization of the polarization properties of the multi-beam telescope system. The CMB polarization is expected to be at a level of only 10% of the temperature anisotropy quadrupole, and the success of the measurements will depend crucially on the precise knowledge of the polarization properties of the telescope.

## 2 Formulation of the problem

While the performance of the antenna structures can be thoroughly tested in terrestrial conditions, the coupling of the HFI with the telescope is, basically, optimized through the computer simulations.

Among various simulation techniques, physical optics (PO) is the most adequate one for the given purpose. However, conventional implementations of the technique  do not fit the size of the problem. Commercially available packages are also very limited in their capacity to rigorously answer this sort of questions. For example, even the best commercial software requires many hours of the full-scale physical optics computation of the main beam of the telescope at the relatively low frequency of $`143GHz`$, while all conventional physical optics codes collapse at the highest frequencies of $`545857GHz`$.

To solve the problem, we developed a special PO code  that allowed us to overcome the limitations of a generic approach for large multi-reflector systems and perform typical simulations of the telescope in the order of minutes. Generally, it requires only $`1`$ minute for a polarized Gaussian beam at the frequency of $`143GHz`$ and about $`30`$ minutes for the beam of $`30`$ modes at $`857GHz`$ using a PC Pentium III ($`500MHz`$) under the Linux operating system.

We simulated the beams of twelve linearly polarized Gaussian horns operating in the transmitting mode at the nominal frequencies $`143`$, $`217`$, and $`353GHz`$. Each horn is designed for the simultaneous measurement of two orthogonal linear polarizations, ’a’ and ’b’, of the incoming radiation. The polarization is characterized by the tilt $`\psi `$ of the electric field at the beam axis in the sky with respect to the local vertical as defined below. The tilt is $`\psi _a=+45^{}`$ and $`\psi _b=45^{}`$ for the horns 143-1, 143-2, 217-5, 217-6, 353-3, and 353-4 (see Fig. 1, a), while $`\psi _a=0^{}`$ and $`\psi _b=90^{}`$ for the other six polarized horns (the angle $`\psi `$ is measured clockwise from the upward vertical direction as seen from the telescope).

We compared three different models of the horn feed, F1, F2 and F3, when the source field is specified in different manner. The feed model F1 is represented by the far-field power and phase patterns of the actual corrugated horn in $`0^{}/45^{}/90^{}`$ cuts which are used in a somewhat approximate manner for computing the incident complex-vector electromagnetic field on the secondary mirror. The feed F2 is specified by the ’clipped Gaussian’ model distribution of the horn aperture electric field which is accurately propagated to the secondary mirror and further through the system. Finally, the feed F3 is represented by the far-field patterns similar to the model F1 which, however, are computed precisely from the clipped Gaussian aperture field of the model F2 and then used in an approximate manner identical to the feed model F1.

The electric field at the aperture of the Gaussian horn in the feed models F2 and F3 is specified as

$$\stackrel{}{E}(r)=E_0\mathrm{exp}(r^2/w_a^2)\mathrm{exp}(iqr^2)\stackrel{}{e},0ra$$
(1)
where $`r`$ is the radial coordinate, $`w_a`$ is the beam radius at the horn aperture, $`a`$ is the horn aperture radius, $`q=\pi /\lambda L`$, $`L`$ is the curvature radius of the wavefront at the aperture (approximately, $`L`$ is the horn slant length), and $`\stackrel{}{e}=\stackrel{}{e}_{a,b}`$ is the unit polarization vector ($`\stackrel{}{E}(r)=0`$ if $`r>a`$). The aperture model of this kind is quite accurate in representing the field of real Gaussian horns. In particular, both the power and phase patterns of such a horn computed at the frequency $`100GHz`$ coincide perfectly well with the experimental data available for the model Gaussian horn designed for this frequency. However, some discrepancies in the sidelobes at the level of about $`30dB`$ appear for the horns specifically optimized for the higher angular resolution.

The comparison of the feed models F2 and F3 shows that the minor differences can only be observed in the polarization patterns at the periphery of the telescope beam while the power patterns of the main beam are identical for both models. The beams computed with the feed models F1 and F3 are also very similar when the parameters of the model F3 are properly adjusted. In this case, also, it is only the polarization at the periphery of the beam that differs slightly for different models while the power patterns are very similar even at the level below $`30dB`$. Finally, the orientation of the polarization vector on the beam axis in the sky field is precisely the same for all the feed models, and the polarization pattern, in general, is rather independent of both the horn pattern and of the fine features of the propagation model used in the simulations. It proves that the far-field model F1 can safely be used for the simulation of the telescope beams from the real horns despite the approximations used in the model.

## 3 Polarization of the Gaussian beams

Fig. 1, b, shows the power pattern of the telescope beam H-143-1 as projected on the plane normal to the telescope line-of-sight at the $`(0,0)`$ point ($`a_x`$ and $`a_y`$ are the horizontal and vertical axes on the plane, respectively, measured in degrees). The beam axis is at the point $`a_x=1.3645^{}`$, $`a_y=1.1975^{}`$ which is defined as the point of maximum power of the beam. The pattern is computed for the clipped Gaussian feed model F3 with adjusted parameters $`a=3.9mm`$, $`w=2.8mm`$, and $`L=120mm`$, simulating the enhanced horn operating at the frequency $`f=143GHz`$.

The beam is well shaped down to almost $`30dB`$ below the maximum and can be approximated by a Gaussian function (at this level, the pattern does not depend on polarization) with the full beam width of $`W_{min}=7.05`$ arcmin and $`W_{max}=7.50`$ arcmin measured at $`3dB`$. The polarization of the beam is generally elliptical except precisely at the beam axis where it remains linear. In order to achieve the required orientation of the polarization pattern in the sky, we should orient the polarization vector $`\stackrel{}{e}`$ properly on the horn aperture .

For immediate comparison of polarizations measured by different horns when scanning through the sky, we should use easily aligned directions in the sky as equivalent reference polarization axes for different beams. Such directions are the meridians in the spherical frame of the telescope, with the pole being the telescope spin axis tilted by $`\eta =85^{}`$ with respect to the line-of-sight (the meridians define local verticals at various observation points while the parallels are local horizontals that constitute the orthogonal directions).

Now, we should define the reference axis for the polarization vector $`\stackrel{}{e}`$ of a tilted horn. We define the reference axis as the direction of $`\stackrel{}{e}`$ in the horn aperture plane that is projected on the vertical axis of the focal plane. The orientation of $`\stackrel{}{e}`$ is specified by the angle $`\varphi `$ in the horn aperture plane measured from this reference in a clockwise direction when looking from the horn to the secondary mirror.

Using these definitions, we have found that the beam of the H-143-1 horn is polarized at its axis at the required angles of (a) $`\psi _a=+45^{}`$ and (b) $`\psi _b=45^{}`$ if the horn polarization vector $`\stackrel{}{e}`$ is specified by the angle $`\varphi _a=+42.99^{}`$ and $`\varphi _b=47.01^{}`$, respectively.

When the polarization vector $`\stackrel{}{e}`$ is properly oriented, the cross-polarized component of the far field measured along the respective orthogonal direction in the sky is minimized. In such a case, the power pattern of the cross-polarized component is, basically, determined by the power pattern of the semi-minor axis of the polarization ellipse at each observation point of the beam. This pattern is very much the same for any orientation of the polarization vector and can be approximated as follows

$$P_{cr}=P_{cr0}[(\theta /\theta _0)\mathrm{sin}(\phi \phi _0)]^p\mathrm{exp}[(\theta /\theta _0)^q]$$
(2)
where $`P_{cr0},dB`$, is the maximum power of the minor axis component achieved at the points specified by the polar angle $`\theta _0`$ and the azimuthal angles $`\phi _0\pm 90^{}`$ measured from the center of the pattern which is located at the beam axis (the values of the parameters in this approximation depend on the position of the horn considered but the magnitude of $`P_{cr0}`$ is, typically, more than $`30dB`$ below the maximum value of the total power of the beam).

Since the actual frequency bands of the horns are rather wide ($`121166`$, $`182252`$, and $`295411GHz`$), the far-field patterns of real horns differ essentially at the lower and upper frequencies of each band. It results in somewhat different power and polarization patterns of the telescope beam as shown in Figures 2 - 4. Nevertheless, despite these differences, the polarization angles at the beam axis are the same irrelevant of the horn patterns and of the operating frequency.

The values of the aperture polarization angle $`\varphi `$ computed for all the HFI horns as well as the direction cosines (the Cartesian components of the unit vector) of the aperture electric field with respect to the coordinate frame $`M1`$ of the primary mirror are summarized in the Table 1.

## 4 Conclusions

A fast physical optics code has been developed for the analysis of the dual-reflector submillimeter-wave telescope on the ESA PLANCK Surveyor. The code overcomes the limitations of a generic approach for large multi-reflector quasi-optical systems and can perform typical simulations of the telescope in the order of minutes.

Simulations of the telescope beams from the linearly polarized horns have shown that the far-field of the telescope is, generally, elliptically polarized except precisely at the beam axis where the linear polarization is preserved. The magnitude of the minor semi-axis of the polarization ellipse in the telescope beam remains at the level of $`30dB`$ below the maximum total power of the beam even for the most tilted horns located at the edge of the HFI horn array.

When rotating the polarization vector of the horn field about the horn axis, the major axes of the polarization ellipses at the central part of the telescope beam rotate virtually by the same angle about the beam axis while the power patterns of both the co- and cross-polarized components of the far field remain basically unchanged.

Thus, the orthogonal polarization directions at the beam axis in the sky are transformed into the orthogonal directions in the focal plane of the telescope. Some deviations from this basic rule which occur mainly at the periphery of the beam are of minor importance since they virtually do not contribute to the measured total power of the polarized component.

## 5 Acknowledgments

The author is grateful to Bruno Maffei for the power and phase patterns of enhanced horns and to Yuying Longval for providing the updated positions and aiming angles of the HFI horns along with the map of the focal plane. The author would like to acknowledge the support of Enterprise Ireland.