Inhomogeneous extragalactic magnetic fields and the second knee in the cosmic ray spectrum

0.47 MB·English

by  Kumiko Kotera

APS/123-QED
Inhomogeneous extragalactic magnetic fields
and the second knee in the cosmic ray spectrum
Kumiko Kotera∗ and Martin Lemoine†
Institut d’Astrophysique de Paris
UMR7095 - CNRS, Universit´e Pierre & Marie Curie,
98 bis boulevard Arago
F-75014 Paris, France
(Dated: February 5, 2008)
Various experiments indicate the existence of a second knee around energy E = 3 1017 eV in
×
thecosmic ray spectrum. This feature could bethe signature of theend of thegalactic component
and of the emergence of the extragalactic one, provided that the latter cuts off at low energies.
Recent analytical calculations haveshown that thiscut-off could be a consequenceof theexistence
8 of extragalactic magnetic fields (Refs. [1, 2]): low energy protons diffuse on extragalactic magnetic
0 fieldsandcannot reach theobserverwithinagiven time. Westudytheinfluenceofinhomogeneous
0 magneticfieldsonthemagnetichorizon,usinganewsemi-analyticalpropagation code. Ourresults
2 indicate that, at a fixed value of the volume averaged magnetic field B , the amplitude of thelow
n energy cut-off is mainly controled by the strength of magnetic fieldshinithe voids of the large scale
a structuredistribution.
J
PACSnumbers: 98.70.Sa,98.62.En,98.65.Dx
4

]
h I. INTRODUCTION extragalacticcomponentappearsatthe ankle,oneneeds
p to postulate the existence of a third cosmic ray compo-
- nentbetween 1017eVand 1019eV(seeforinstance
o Recent developments in our undertanding of the ori- ∼ ∼
Refs. [6] for a recent proposal), or to assume that super-
r
gin of cosmic rays have led to the suggestion that the
t novaeareabletoaccelerateparticlesuptotheankle(see
s cosmic ray spectrum might comprise only two compo-
a for instance [7]).
[ nents, one of galactic origin, dominant in the energy
range E . 1017eV, and another of extragalactic origin Themoderninterpretationofatransitionbetweenthe
2 in the energy range E & 1018eV [3]. In this interpreta- galactic and the extragalactic component at the second
v tion, the cross-over between the two components marks knee does not come without flaws, however. In particu-
1
the so-called“secondknee” in the allparticle cosmic ray lar, the smooth matching of the galactic and extragalac-
9
8 spectrum (see Ref. [4]for a review oncosmic ray data at tic components at E ∼ 1017eV bears the unaesthetic
1 highenergies),whilethe “firstknee”atE 2 1015eV look of fine-tuning. Of course, this problem is generic
. wouldbeassociatedwithachangeofpropa∼gatio×nregime to the matching of two distinct components at a point
6
0 orthemaximalenergyofprotonsatthesource. Thethird wheretheslopesteepens;theintroductionofathirdcos-
7 feature in the cosmic ray spectrum, i.e. the “ankle” at mic ray component would not help in this respect. As
0 E 1019eV would be a consequence of pair production farasthe galacticcomponentisconcerned,the fall-offat
: loss∼esofultra-highenergyprotonsacceleratedinsources E & 1017eV arises as a direct consequence of the ob-
v
i located at cosmologicaldistances [3]. servationof the firstknee, as mentionedabove. However
X one must explain why the extragalactic component van-
On purely phenomenological grounds, this modern
r ishes at energies below the second knee. In the original
a view is appealing when compared to the more tradi- scenario of Berezinsky et al. [3], this low energy cut-off
tionalinterpretationinwhichtheankleisassociatedwith
wasattributedtophysicsatthesource. InRef.[8],itwas
the emergence of an extragalactic cosmic ray compo-
suggested to interpret it as the modulation of the extra-
nent out of a more steeply falling spectrum at energies
galacticfluxduetoagalacticmagnetizedwind,although
E . 1019eV. One of these major advantages certainly
the calculations of Ref. [9] bring down this cut-off to a
is the “economy” of sources. It is indeed notoriously too low energy, E 1015 1016eV.
difficult to accelerate particles beyond E 1015eV in ∼ −
supernovaeremnants [5]. Therefore,if the fi≈rstknee cor- Anotherpossibility,advocatedinRefs.[1,2]istorelate
respondstothemaximalenergyofprotonsatthesource, this cut-off with the influence of extragalactic magnetic
thefall-offofthegalacticcomponentat 1017eVwould fields. As demonstrated in these studies, if the inten-
naturally be associated with the maxim∼al energy of the sityofextragalacticmagneticfieldsisrathermodest,say
iron component at the source. On the contrary, if the B 10−9G, the diffusion time of particles with en-
hergiy ∼E . 1017eV from the closest sources (located at,
say 50 100Mpc)becomeslongerthanthe ageofthe
∼ −
Universe. Thisproducesalowenergycut-offintheprop-
∗Electronicaddress: [email protected] agated spectrum at the required location, which allows
†Electronicaddress: [email protected] to reproduce a smooth transition at the second knee in
2
agreement with observational data (see [1]). on the results. As a by-product of the present study,
It has also been argued that if ultra-high energy cos- we thus propose a simple and new recipe to build semi-
mic rays comprised a significant fraction of heavy nu- realistic magnetic field distributions out of dark matter
clei (& 20%), the scenario of a transition at the second simulations (which can be obtained at a lesser cost than
knee would loose its merits as the energy losses would MHD numerical simulations) as well as a new transport
no longerbe able to reproducethe ankle feature[10, 11]. schemewhichismoreefficientthanexistingcodesinsev-
However even a solar type (or Galactic cosmic ray type) eral respects. In particular, it allows to enlarge artifi-
chemical composition, with 10% helium and only cally the range of scales on which the magnetic field is
traces of heavier elements, allo∼ws to fit the existing data distributed, hence to model the influence of intergalac-
at the ankle with a single powerlaw spectrum at injec- tic magnetized turbulence on particle transport. These
tion [3]. Furthermore, there is no particular reason to techniques,which aredevelopedinSection II andin Ap-
expect the source composition to be enriched in metals. pendix A, allow us to bracket the possible distributions
In any case, future measurements of the chemical com- of extragalactic magnetic fields at the present time and
position will tell [10, 12]. their impact on the ultra-high energy cosmic ray spec-
trum. Our results indicate that, at a fixed value of the
Quantizing the influence of extragalactic magnetic
volume averaged B , the amplitude of the low energy
fields on the spectrum of cosmic rays with energy E
h i
1017eV is not an easy task as the propagationtimes be∼- cut-off is controlled by the strength Bvoid of magnetic
fields in the voids of the large structure distribution and
come of the order of a Hubble time, hence one must ac-
count for the effects of expansion. For the sake of sim- thesourcedistancescalens−1/3 = 50Mpc. Thefactthat
plicity, Refs. [1, 2] have thus assumed the magnetic field ourconclusionsdependmoreweaklyonotherparameters
power to be distributed homogeneously in space. How- characterizingthemagneticfielddistributionprovidesan
ever, this approximation deserves to be refined since the adequate a posteriori justification for our semi-analytic
magnetic field is most likely distributed as the charged construction. We also argue that this simulation tech-
baryonic plasma. Since the scale of inhomogeneity of nique offers various advantages over existing full-blown
largescalestructureintheUniverseiscomparabletothe MHD simulations of large scale structure formation, at
distance to the closest sources, 50 100Mpc, the in- least as far as cosmic ray propagation is concerned.
∼ −
homogeneityofthe magneticfield mayaffectthe conclu- The paper is organized as follows. In Section II, we
sions of Refs. [1, 2]. The objective of the present paper presentourschemeofinhomogeneousmagneticfieldsim-
is precisely to address this issue and to study the sce- ulationandthe numericaltechnique ofcosmic raytrans-
nario put forward in these references in a more realistic port. We compare these techniques to existing simula-
extragalactic magnetic field configuration. tions and discuss the advantages and drawbacksof each.
Thisimmediatelybringsforwardthedifficultyofdefin- InSectionIII weaddressthe issue ofthe low energycut-
ing a realistic distribution of large scale extragalactic offinvariousmodels of extragalacticmagneticfields dis-
magnetic fields, including the shape and amplitude of tributions, compute the spectra and compare them to
a turbulent magnetic cascade. From an observational experimental data. Section IV discusses the limitations
pointof view,one has beenable to measurethe strength of our approach and possible future avenues of research.
ofextragalacticmagneticfields“only”inthecoreofclus- Finally, Section V summarizes our findings.
ters of galaxies [13]; a bridge of synchrotronemission on
Mpc scales has been observed in the Coma cluster [14].
Hopefully the SKA project will enlarge considerably the II. PROPAGATION OF HIGH ENERGY
datasetonextragalacticmagnetic fields [15] but it is not COSMIC RAYS IN EXTRAGALACTIC
expected to enter operation before 2017. In the mean- MAGNETIC FIELDS
time, one thus has to rely on theory. Unfortunately, the
very origin of extragalactic magnetic fields is unknown, The straightforward way to study the influence of
see Ref. [16]for a review. Furthermore,evenif one knew extragalactic magnetic fields boils down to performing
exactly the initial conditions that set the configuration MonteCarlosimulationsofparticlepropagationinasim-
of magnetic fields at a high redshift, the simulation of ulatedmagnetizedUniverse. This,however,bringsintwo
theirevolutionthroughoutcosmichistorytothepresent, majordifficulties,whichwerealludedtoearlierbutwhich
carrying sufficient accuracy on a large dynamic range of arerarelydiscussedin the literature: (i) anaccurate nu-
spatial scales, remains a formidable task for numerical merical modeling of the transport of charged particles
computing. in magnetic fields; (ii) an accurate numerical modeling
In regards of all these uncertainties on the origin of of the magnetized volume, including magnetized turbu-
extragalacticmagneticfields, ontheir distributioninthe lence.
present Universe, on the nature and shape of magnetic Point (i) deals with the theory of cosmic ray diffu-
turbulence as well as on the transport properties of par- sion, which in spite of a long history and recent major
ticles in chaotic magnetic fields, we adopt a simplified progress, has not yet reached a consensus on the trans-
and parametrized description which allows us to eval- port of cosmic rays in MHD turbulence (see Ref. [17]
uate the effects of the various sources of uncertainties for a recent review). Actually, the simulation of parti-
3
cle transport in a well-defined MHD environment is not
trivial even from a purely numerical point of view. For
example, Ref. [18] has demonstrated that the interpola-
tion of the magnetic field from a numerical grid gives an
erroneousdescription of particle transportif the Larmor
radius r . l , where l represents the grid size,
L min min
i.e. the minimum scale of the turbulence inertial range.
Point (ii) deals with the problem of simulating realistic
MHD flows on a large range of spatial scales, which also
constitutes a field of research in its own right.
A. Magnetic field modeling
Severalpioneeringworkshavestudiedthe propagation
of cosmic rays in so-called “realistic” magnetized envi-
ronments [19, 20, 21]. These studies have constructed
the magnetized cube out of a hydrodynamical simula-
tion of largescale structure formation, which follows the
magnetic field in a passive way for Ref. [20], and with
feedback effects on the matter evolution for Ref. [21].
The initial conditions for this magnetic field have been
set at a highredshift (althoughRef. [20] also models the FIG.1: Volumefillingfactorofthemagneticfieldindifferent
production of magnetic fields at accretion shock waves) scenarios. Indot-dashedline,themagneticfieldsimulatedby
and the overall amplitude of this field has been rescaled Dolag and coauthors in Ref. [19]; in long-dashed line, that
simulated by Sigl and coauthors [20]. In solid line, the semi-
at the end of the simulation so as to reproduce the ob-
analytic model with B ρ2/3, in orange dotted line, B
served strength of magnetic fields in the core of clusters ρ. In red dashed line, th∝e model B ρˆ1+(ρ/ρ)−2˜; th∝is
of galaxies. This ingenious procedure allows to fix the ∝
model simulates a volume with unmagnetized voids. In all
volumeaveragedmagneticfieldindependently ofthe ori-
cases theproportionality factor B0 =2 nG.
gin of the magnetic field, although the volume averaged
magnetic field now depends on the details of the ampli-
fication process during cluster formation.
InRef.[20],theauthorsfollowthetrajectoryofcosmic phistication of the numerical codes used, cannot truly be
rays using Monte Carlo methods while the authors of deemed as realistic. Italsoindicatestheneedforalterna-
Ref.[19]deriveanupperboundonthetypicalcosmicray tivemethodstostudythetransportofhighenergycosmic
deflection using a semi-analytic transport scheme. Their rays in extragalactic magnetic fields, in order to provide
conclusions are radically different: the former authors new angles of attack on this difficult problem. This con-
derive a typical deflection of 10 20o above 1020eV stitutes one major motivation of the present work, in
while the latter find a deflect∼ion le−ss than a degree at which we develop one such method and apply it to the
these energies. This discrepancy illustrates the inherent study of the low energy cut-off at energies close to the
complexity of such simulations. The complexity and the second knee.
cost of such numerical simulations are such that it has Ourmagnetizedvolumeisconstructedinasimpleway
not been possible to elucidate the precise origin of this ascomparedtoRefs.[19,20],thissimplicityofferingvar-
discrepancy yet. It is likely that most of this difference ious advantages (and admittedly, several drawbacks) as
is to be attributed to the modeling of the extragalactic discussed further below. The core of our method is to
magnetic field, and to a lesser degree, to the transport map the magnetic field strength over the gas density us-
scheme. ingananalyticalrelationB(ρ)(tobespecifiedlater)and
Figure 1, which presents the volume filling factors of to distribute randomly the magnetic field orientation in
the magnetic field strength obtained in these numeri- cells of coherence length l . The gas density itself is
c
cal simulations is particularly instructive (the model of obtained from a high resolution dark matter numerical
Ref.[19]isshownasthedot-dashedline,whilethemodel simulation of large scale structure formation (with stan-
of Ref. [20] is given by the long-dashed line). It reveals dard cosmological parameters Ω = 0.7, Ω = 0.3 and
Λ m
largedifferencesinthevolumeaveragedmagneticfieldas H = 70km/s/Mpc). Once the volume has been set up,
0
well as in the spatial distribution of these fields (which cosmic ray trajectories are simulated as follows. At each
translates in this figure as a difference in the slopes of step, the cosmic ray is supposed to enter a spherical cell
thevolumefillingfactor). Again,theoriginofthisdiffer- ofcoherenceofthemagneticfielddefinedbyitsdiameter
ence is not understood. This figure clearly demonstrates l ,inwhichthemagneticfieldorientationisrandom. The
c
thatthesimulatedmagnetizedvolumes,despitealltheso- timespentinthe cellandthe directionofexitofthecos-
4
mic ray are then drawn from semi-analytic distributions along one (i.e. collapse on a wall) or two (i.e. collapse
whichsimulatethetransportoftheparticleinMHDtur- on a filament) spatial directions. The law B ρ0.9 has
∝
bulence, accordingto studies carriedoutinRefs.[18,22] indeed been observed in the simulations of Dolag and
(see Appendix A for a detailed discussion). The particle co-authors [19].
is then movedto another coherentcell andthe next step Viscosityandshearflowsduringcollapsemayalsoam-
is simulated. As we explain in Section IIIB, we finally plify further the magnetic field, leading to departures
compute the propagated spectrum in the low energy re- from the law B ρ between regions of very different
gion E . 1017eV in a semi-analytic way which allows density, in particu∝lar the voids and the structures (fila-
usto modelthe effectofcosmologicalexpansionoverthe ments and pancakes). For instance both simulations of
course of propagation from the source to the detector. Refs.[19,20]leadto weakervaluesofthe magneticfields
In detail, we use the Monte Carlo simulations of particle in the voids than would be expected from an extrapola-
propagation in the extragalactic magnetic fields at zero tion of the law B ρ to regions of low density. In order
∝
redshiftinordertomeasurethediffusioncoefficients,and to bracket these different effects, we consider several re-
use existing analyticalformulae in order to calculate the lations B(ρ):
propagatedspectra from the diffusion equationin an ex-
B ρ2/3 , (3)
panding space-time. More details on this latter step are
∝
provided in Section IIIB. B ρ0.9 , (4)
∝
Note that the density field of dark matter provides a ρ −2
good approximation to that of the gas density in the in- B ρ 1+ . (5)
tergalactic medium (IGM) on scales larger than a few ∝ " (cid:18)hρi(cid:19) #
hundred kpc, corresponding to the baryon Jeans length.
The last model is an ad-hoc modeling of the suppression
Therefore the overall baryonic gas density field can be
of magnetic fields in the voids of large structure which
obtainedbysmoothingthedarkmatterdistributionbya
leavesunchangedthedistributioninthedenseintergalac-
windowfunctionofthissize. Ourdarkmattersimulation,
tic medium (meaning ρ> ρ ).
runwiththehydrodynamicalcodeRAMSES[23],is5123, h i
with extent 280Mpc, hence with a grid size 560kpc:
≃
in this case the minimum scale of the simulation then
B. Why a semi-analytical propagation method?
plays the role of the window function and no smoothing
is required. Of course, this treatment does not provide
At this stage, one should compare the respective mer-
a perfect description of the gas distribution but we have
itsanddrawbacksofthisnewmethodwithotherexisting
checkedthat decreasingthe resolutionby a factor2 does
techniques. Concerning the magnetic field distribution,
not affect our results. Moreover the gas density field
our method obviously neglects subtle effects such as the
serves only as a marker of the magnetic field distribu-
amplification of the magnetic field in the vicinity of ac-
tion, so that the above error is negligible in comparison
cretion shock waves of large scale structure. However, it
to the uncertainty surrounding the strength and config-
shouldbeclearthatnonumericalsimulationcanclaimto
uration of the magnetic field. The essential is rather the
simulate with accuracy the magnetic field in the vinicity
law B(ρ) which provides the mapping between the mag-
ofcosmologicalshockwavesduetotheintricacyofMHD
netized volume and the density field.
physics at play. The amount of amplification, the coher-
In the case of isotropic collapse, it is well known that
ence length and the shape of the turbulence spectrum
B ρ2/3inaplasmaofinfiniteconductivity. Thislawis
∝ remain open questions (see however [25] for a detailed
slightlyoversimplisticbecauseitignorestheanisotropyof
discussion of the Weibel instability operating at inter-
collapse in large scale structures, which results in an en-
galactic shocks).
hanced amplification of the magnetic field by shear and
ByconsideringaonedimensionallawB(ρ),oursimula-
anisotropic compressive flows [24]. In detail, the equa-
tionapparentlyneglectstheinfluenceofthevelocityfield
tions of ideal MHD lead to the conservation law:
onthe magnetic fieldamplification. Indeed, B shouldbe
d B B a multidimensional function that depends on ρ as well
= v , (1)
dt ρ ρ ·∇ as on the velocity field in order to take into account dy-
(cid:18) (cid:19) (cid:18) (cid:19)
namo and shear effects. The three models alluded to
which allows to derive the magnification of B from the
earlier do actually account for these effects up to some
deformationofthedensityfield. Iftheseparationδq<be-
extent, as they reproduce the characteristic features ob-
tweentwopointsismappedintoδq> = δq< through
D· tainedinthe numericalsimulationsthatinclude dynamo
deformation, where indicates the deformation tensor,
then B</ρ< is mappDed into B</ρ<, so that: and shear effects. The choice of a random orientation of
D· themagneticfieldineachcoherencecellalsoneglectsthe
B> ρ> B< influence of large scale motions. Simulations which fol-
= |D· | . (2) lowexplicitelythemagneticfieldindicatethatthislatter
B< ρ< B<
tends to be aligned with the principal directions of large
From Eq. (2), it is then easy to derive the law B ρ2/3 scale structure, i.e. the axis of the filament for exam-
∝
for isotropic collapse, or B ρ for anisotropic collapse ple[26]. However,oneshouldnotethatthesimulationsof
∝
5
Refs.[19,20,26]assumeaninitialmagneticfieldwithin- (which sets the initial l ) as well as to the velocity fields
c
finitecoherencelength,sothatthefinalcoherencelength whichdistortthe fieldduringtheevolution,inparticular
of the final magnetic field along the axis of the filament withthe upbringingofMHD turbulence. Hereaswell, it
equalsthelengthofthefilament(seeFig.4ofRef.[26]for should be clear that no simulation can claim to simulate
anillustrationofthiseffect). Thisofcourseisunrealistic these various effects with accuracy.
(unless one assumes an acausal origin for B) since the The issue of turbulence in the IGM is delicate, be-
coherence length of the magnetic field should not exceed cause the Reynolds number in the intergalactic medium
1Mpc, the typical turn-around time of an intergalac- may take large or moderate values depending on the en-
∼
tic eddy of this size being comparable to the age of the vironment. For turbulent excitation on a length scale L
Universe [27, 28]. at velocity v, this Reynolds number reads [33]:
In contrast, our simulation presents the advantage of
−5/2
simulating this multiple field reversalalongthe filament. L v T ρ
e 105 .
The alignment of the magnetic field direction along the R ≃ 1Mpc 300km/s 105K ρ
filamentshouldtakeplaceifthecoherencelengthislarger (cid:18) (cid:19)(cid:18) (cid:19)(cid:18) (cid:19) (cid:18)h i((cid:19)6)
thanthetransversesizeofthefilament(moreexactlythe Turbulenceisthusprobablyfullydevelopedinmostof
typicalscaleheightofthedensitygradient). If,asismore the IGM, except in the high temperature regions repre-
likely,the coherencelength is smaller,then the compres- sentative of clusters of galaxies. In these regions of high
sion is similar to planar collapse as far as the magnetic kinematic viscosity, the shape and extent of the inertial
fieldinacellisconcerned,hencethefieldbecomesaligned range of turbulence is rather complex, and most likely
transversely to the density gradient, and not necessarily influenced by the strong magnetic field, see Ref. [34] for
along the filament axis. Our method offers the means detailed discussions. Turbulence plays a fundamental
to include this effect but we leave this investigation to role in the transport of charged particles as well as
future work for simplicity. in the reshaping of the distribution of the magnetic
A last point concerning the distribution of the mag- field, but its incorporation in numerical simulations of
netic field in the intergalactic medium is related to its the gas density is extremely complex. As the largest
origin. The simulations of Refs. [19, 20, 26] set the ini- scale of the turbulence cannot exceed a few hundred
tial conditions at high redshift and ignore other sources kiloparsecs or a megaparsecs, taking into account just
of magnetic pollution of the intergalactic medium (see one or two decades of inertial range necessitates an
Ref. [16] for a review on the origin of cosmic mag- unrealistically high resolution since the cube size must
netic fields), such as galactic outflows, AGN pollu- remain larger than the inhomogeneity scale 100Mpc.
∼
tion [29, 30, 31], the amplification of magnetic fields in In this respect, our simulations provide more flexibility
accretionshocksoflargescalestructure(exceptRef.[20] because our transport scheme in the simulated magnetic
which uses a model of Biermann battery effects), tur- field allows to simulate the influence of a turbulence
bulent amplification in the IGM [32], etc. These sources spectrum down to scales well below rL, see Appendix A.
shouldinfluencethetransportofhighenergycosmicrays
in two ways: by modifying the relationship between the The present simulation technique thus combines sim-
volume averaged B and the value observed in cluster plicity,efficiencywithflexibility,andtheapproximations
h i
cores, and by adding additional scattering centers which on which it rests appear reasonablein regardsof the un-
have been omitted in these simulations. Our simulation certainties surrounding the origin of extragalactic mag-
techniqueoffersthefreedomtoincludesuchlocalizedpol- neticfields,thenatureofMHDturbulenceandtheprop-
lution effects in a simple and efficient way: one could in- erties of cosmic ray transport in such turbulence. Most
clude these highly magnetized regions in our simulation importantly, its parametrized description allows us to
cubebysamplingthemaccordingtothelocalmatterden- test the influence of the various parameters on the re-
sity. For the sake of simplicity, in a first step we ignore sults in contrast with most other works on this topic.
theseadditionalsourcesandpostponetheirstudy tofur-
ther work. This choice is conservative in so far as the
inclusion of localized regions of enhanced magnetic field III. RESULTS
would tend to amplify the magnitude of the low energy
cut-off, all things being equal. The following results were obtained by computing the
Finally, our simulation assumes for simplicity that the trajectories of 103 protons in inhomogeneous magnetic
coherence length l is uniform in space, whereas it is fields mapped accordingto four models,for many sets of
c
likely to evolve as a function of the density and veloc- energies E, magnetic field characteristic values B and
0
ity fields. However, this brings in additional parameters coherencelengthsl . Wewilllabelinwhatfollows“mod-
c
which enlarge the parameter space. We believe that at els 1 3” our modeling of B(ρ) presented in Eqs. (4 5).
− −
this stage,itis morereasonableto study the influence of We add to these models a last one (model 4) for which
l by performing different runs with different values of l B ρ2/3 and the level of turbulence η = δB2 / B2
c c
∝ h i h i≪
andcomparingtheresults. Furthermore,theactualvalue 1, where δB is the inhomogeneous perturbation compo-
l isintimately relatedto the originofthe magneticfield nent of B (defined such as: B = B +δB).
c
h i
6
FIG. 3: Root mean square of the distance of 103 particles
to their source after one Hubble time (tH 13.9 Gyr) as a
∼
FIG. 2: Trajectories of protons of different energies (solid function of their energy, for B0 = 2 nG and lc = 300 kpc.
line: E = 1017 eV, dotted line: E = 1019 eV) in a slice of Thesolidlinerepresentstherootmeansquareofthedistance
andthesurroundingcolor band itsvariance. Thedotted line
simulateduniverse. Thecharacteristicmagneticfieldistaken
showsthevaluesobtainedfromanalyticalcalculationsinaho-
as B0 =1nGand thecoherencelength as lc =100 kpc. The
mogeneousmagneticfield(Eq.9)andthedot-dashedlinethe
colorbar on the side indicates the intensity of the magnetic
threshold energy above which the energy loss time becomes
field (in log).
.tH/2.
Though B and B have quite similar numerical val-
0
h i
ues,theyarenotstrictlyequal(theydifferapproximately diffusive trajectory, one notices the expected intuitive
by a factor 1.5). B represents the volume averaged correlation between the fluffiness of the trajectory and
h i
magnetic field and B0 is the proportionality factor in the clustered regions.
models (1 4), so that B = B f(ρ), where f(ρ) is
0
− ×
dimensionless.
Figure3illustratesthesecommentsinamorequantita-
Theparticlesareemittedfrom10differentsourcescho-
tive way. It shows the root mean square of the distance
sen randomly amongregions of high baryonicdensity. A
of 103 particles to their source after one Hubble time
detailed description of our code is given in appendix A.
(13.9 Gyr) as a function of their energy, for a character-
Wefirstspottheexistenceofamagnetichorizonusing
istic magnetic field B = 2 nG and a coherence length
the isotropic collapse magnetic field model (Eq. 4). We 0
l = 300 kpc. Our results are no longer valid beyond
then move on to other models, calculate their resulting c
the dot-dashed line which represents the threshold en-
propagation spectra and study the influence of our two
ergy above which the energy loss time becomes . t /2,
main parameters: B and l . H
0 c as our simulations do not compute energy losses.
The first striking remark is that particles of energy
below E 3 1017 eV cannot travel farther than a dis-
A. Existence of a magnetic horizon ∼ ×
tanceofa hundredmegaparsecsfromtheir sources. This
corroborates the scenario of Refs. [1, 2] on the existence
Figure 2 shows two examples of proton trajectories
of a magnetic horizon, and extends it to the case of a
(solid line: E =1017 eV, dotted line: E =1019 eV) in a
inhomogeneous magnetic field.
slice of simulated universe. The characteristic magnetic
The curve presented in Fig. 3 comprises three distinct
field is taken as B = 1 nG and the coherence length as
0
parts: a diffusive part with a slope of 1/6, a semi-
l =100 kpc.
c ∼
Obviously the particles at E = 1019 eV and diffusive part (slope ∼ 1) and a quasi-rectilinear part
E = 1017 eV evolve completely differently: the former with a slope tending towards zero. These trends can be
naturallyexplainedbyanalysingthepropagationregimes
travel in a rectilinear regime without being affected
at different energies, in homogeneous magnetic fields.
by changes in density, whereas the latter experience
a diffusive propagation. Taking a closer look at the The quantity r2 1/2 plotted in Fig. 3 can be easily
h i
7
related to the diffusion coefficient through the equation: r and then D using (A2). It is then easy to obtain Rˆ
L
and calculate T using (B3).
r2 =2Dt , (7)
H Fortheisotropiccollapsemodel(model1)andtherep-
h i
resented parameters (B =2 nG, l = 100 kpc), there is
where t is the Hubble time and D the diffusion coeffi- 0 c
H anoticeabledifferencebetweenthehomogeneousandthe
cient. Thisequationfollowsstraightlyfromthedefinition
inhomogeneouscases. The cut-offoccursatlowerenergy
of D ∆x2 /2∆t, where ∆x represents the displace-
≡ h i for the inhomogeneous case, probably due to voids that
ment during the time interval ∆t. Thus our computa-
enable particles to travel farther.
tionofD inoursimulationshasadirectinfluence onthe
The previous remark does not stand for a travelled
shape of the curve observed in Fig. 3.
distance of10Mpc (blue dashedlines). Onthe contrary,
As explained in appendix A, our diffusion coefficient
the transmission factor in an inhomogeneous magnetic
is calculated following the results of Casse et al. [18]. It
field is lower than in the homogeneous case. This is due
accounts for both diffusive (r l ) and semi-diffusive
L ≪ c to the influence of the dense environment of the source
(r >l )regimes. Inthecaseofadiffusiveregime,equa-
L c whereparticleswereemitted. Onasmallscaleof10Mpc,
tion (A2) becomes
low energy particles have just escaped the high density
D r1/3l2/3, (8) regionsurroundingtheir sourceandcannotpropagateas
∝ L c far as in the homogeneous case. The influence of the
whichcorrespondsto the standardKolmogorovdiffusion environment will be discussed in section IVB.
regime. Besides,whentheLarmorradiusr issomewhat
L
greater than the coherence length of the magnetic field,
we have the well known dependence B. Calculated spectra
D r2l−1. (9)
∝ L c In order to compare our results with observational
data, we derive spectra from our simulations in the fol-
Knowing that r EB−1, we get from (8) and (9):
r2 1/2 E1/6 aLtlo∝wenergiesand r2 1/2 E forhigher lowing way. Ref. [1] has shown that the solution to the
h i ∝ h i ∝ diffusion equationin an expanding Universe, assuming a
energies.
constant comoving distance between scattering centers,
It is quite surprising that these slopes, expected for
and limiting itself to energy losses by expansion (which
homogeneous magnetic fields, are also observed in our
is correct at energies below 1018.3eV) takes the form:
inhomogeneous simulations for the magnetic field model ∼
and the set of parameters presented in Fig. 3. We will c e−ri/(4λ2) dE (t,E)
show in the following section that this is not true for J = dt g Q(E (t,E)).
other models and parameters. diff 4π (4πλ2)3/2 dE g
Z i
When we get to very high energy (E 1018.5), the X (10)
slope of r2 1/2 versus E gets weaker, as p∼articles enter
h i
the quasi-rectilinearregime. Eq.(7) is no longervalid as This solution agrees with Ref. [35], which derives the
particles never reach the diffusion regime. general diffusion equation in an expanding universe and
Another illustration of the existence of the magnetic
horizon is presented in figure 4. The transmission factor
is plotted as a function of particle energy for three dis-
tances to the source (dashed lines: 10 Mpc, solid lines:
100Mpc,dottedlines: 300Mpc). Givenaninitialsource
position, we propagate protons over one Hubble time.
At a distance R from the source, we calculate the trans-
mission factor by taking the ratio between the number
of particles situated beyond R and the total number of
particles that were emitted.
Figure 4 clearly indicates the presence of a magnetic
horizon: for energies below 2 1017 eV, only half of
∼ ×
the emitted particles reach a distance of 100 Mpc in a
Hubble time. The cut-off energy is lower than for the
case representedin Fig.3 due to the lowervalue ofl , as
c
will be explained in section IIIB.
Thinlinesrepresenttheanalyticaltransmissionfactors FIG.4: Particletransmissionfactoratvariousdistancesfrom
calculated in appendix B using the diffusion coefficient the source, as a function of particle energy. Thick lines
implemented in our code (Eq. A2), for the homogeneous are results from the simulation run with B0 = 2 nG and
case. Foragivenenergywithaparticularsetofparame- lc = 100 kpc. Thin lines represent the analytical transmis-
sion factor for a homogeneous magnetic field (Eq. B3).
ters, one can calculate the corresponding Larmor radius
8
presents the solutions for various energy losses.
InEq.(10)above,r representsthe comovingdistance
i
to source i, E (t,E) the required energy at time t in
g
order to have an energy E at t given the energy losses,
0
Q(E ) the emission rate per source at energy E and λ
g g
the comoving “path length”. λ is defined as:
t0 dt a E
λ2 = D e e , (11)
a(t) a(t)
Zte (cid:20) (cid:21)
where a is the scale factor at emission and D the diffu-
e
sion coefficient. Physically, λ represents the typical dis-
tancetravelledbydiffusion,accountingforenergylosses.
In order to calculate λ, we first study the dependence
of D on E, B and l using our simulations. We find
0 c
that each set of parameters corresponds to a different
functionD(E,B ,l ). Wethenparametrizetheevolution
0 c
of the magnetic configuration as done by Berezinsky &
Gazizov [36], as:
l (z)=l (1+z) and B (z)=B (1+z)2−m, (12)
c c 0 0
where m characterizes the MHD amplification of the
field. For simplicity, we set m to 0 in our calculations.
This toy model corresponds to a constant comoving dis-
tance between scattering centers and ignores magnetic
field amplification during structure formation. In this
way, we obtain the required dependence of D over t and
thus the function D[a E /a(t)].
e e
The function dE (t,E)/dE is calculated by inte-
g
FIG. 5: Upper panel: influence of B0 on the spectra for lc = grating the energy losses, following the calculations
100 kpc and model 1. Middle panel: influence of lc on the of Berezinsky et al. [3]. The injection spectrum ex-
spectraforafixedvalueofB0 =2nGandformodel1. Lower tends from 1016 eV to E = 1020 eV. The function
panel: influenceofdependenceofBoverρonspectra(models max
Q(E ) = K(E /E )−γ gives the emission rate per
1 4), for fixedvalues of B0 =2 nGand lc=100 kpc. g g max
− source at energy E , K being a normalisation factor
g
such that dEEQ(E)=L, with L the total luminosity,
which is assumed to scale as the cosmic star formation
R
rate from Ref. [37]. We will assume in our calculation
a spectral index of γ = 2.6. In any case, it should
be pointed out that the choice of the star formation
rate has little influence on our spectra, since the effects
of the magnetic horizon dominates those of the star
formationhistoryonthelowenergypartofthespectrum.
At higher energies, when the comoving light cone dis-
tance r(t) = t0dt′/a(t′) becomes smaller than λ(t,E),
t
the propagationis no longer diffusive andenters the rec-
R
tilinear regime. In this case, the propagatedspectrum is
given by:
c 1 1 dE (t ,E)
g i
J (E)= Q(E (t ,E)),
rect 4π 4πr2 1+z dE g i
i i i
X
(13)
where t is related to r by r = t0dt′/a(t′), r and z
i i i t i i
denoting the comoving distance and redshift of the ith
FIG. 6: Same as Fig. 3, for B0 =2 nG and lc =100 kpc, for R
source. The factor 1/(1+ z ) was omitted in Ref. [1],
models 1 4. Variances are not represented. i
− but it has no influence whatsoever as z 1 when the
≪
9
FIG. 7: Total spectra (galactic + extragalactic) compared to data. Each row corresponds toa model and a set of parameters.
Caution: n = 10−6 Mpc−3 for the third row. The left panels show KASKADE, Akeno and AGASA data. The right panels
s
showKASKADE,HiRes-1,HiRes-2andFly’sEyesdata. Solid linesrepresentsthemedian valuesofthetotal flux,dot-dashed
linestheseparate galactic and extragalactic componentsandthedottedlines theupper75th andlower 25th percentilesforthe
magnetic cut-off of theextragalactic flux.
rectilinear regime is reached. the first hundred sources were uniformly sampled. For
farther sources, the continuous source approximation is
Figure5presentstheinfluenceoftheB(ρ)-modelsand valid and it was used numerically.
of parameters B and l on the magnetic cut-off. Only The upper panel shows the intuitive result that the
0 c
the diffusive part of the spectra is represented there and greater the mean magnetic field, the steeper the cut-off.
the fall-off of the curve around 1018 eV corresponds Ofcoursethis lawis notrestrictedtomodel1but isalso
∼
to the transition between the diffusive and rectilinear valid for models 2 4.
−
propagation regimes. We assume continuously emitting The middle panel shows interesting features that are
sources with density n = 10−5 Mpc−3 and plot the in agreement with Eqs. (8) and (9). For a fixed value of
s
median spectrum obtained over 100 realisations of the B , for low energies, particles are in the diffusive regime
0
source locations. For each realisation, the location of [see Eq. (8)] and r2 scales with the coherence length
h i
10
as l2/3. For higher energies particles are in the semi- We use the data of six major experiments that
c
diffusive regime [see Eq. (9)] and r2 l−1. In other measuredthe cosmicrayfluxesinourregionsofinterest:
h i ∝ c
words,thespectrumcutsoffmoresteeplyforlowervalues KASKADE (2004 data), with an energy range going
of l for low energies and for greatervalues of l for high from 1015 to 1017 eV [38], Akeno from 1015 to 1018.6 eV
c c
energies. This can be seen on figure 5: for l =300 kpc, [39],AGASAfrom1018.5to1020.5 eV[40],HiResIandII
c
the spectrum cuts off at high energy but the slope is from1017.3 to 1020 eV [41] andFly’s Eyes from1017.3 to
shallow for low energies, whereas for l = 30 kpc, the 1020 eV [42]. We split these data in two sets in order to
c
slope is steep at low energies but the cut-off starts at accountforthediscrepancybetweenHiResandAGASA.
lower energies. This enables us to have two different normalizations for
In view of these trends, one will have to find a good the extragalactic flux on the left and right panels. The
compromiseinordertoobtainsatisfactoryfits totheob- normalization of KASKADE data remains the same for
servational spectra. both sets.
The lower panel illustrates the shape of the cut-off
for the four models previously described. Models 3 and
4 present a much shallower slope compared to models
1 and 2. The almost total absence of magnetic field
in the large scale structure voids for model 3 and the Four main points emerge from Fig. 7. (i) The sec-
cancellation of turbulence in model 4 can explain this. ond knee feature appears more or less clearly in the four
Wealsonoticethatmodel1whichhasahighermagnetic models, but ultimately remains quite robust to model
field intensity in the voids cuts off in a steeper way than changes. (ii) However, again, the influence of the mag-
for model 2. neticfieldintensityinvoidsisobvious: evenwithasource
densityofn =10−6 Mpc−3,thegoodnessoffitofmodel
s
For a better understanding of the trends seen in the 3withtheobservedspectraisonlymarginal. This situa-
lower panel of Fig. 5, we plot in figure 6 the root mean tionis clearlyimprovedinthe othermodels, especially if
square of the distance of 103 particles to their source weconsidertheuncertaintyonthepositionoftheclosest
after one Hubble time as a function oftheir energy,asin sources. (iii)Onemightalsonoticethatthislastelement
Fig. 3, for models 1 4. Variances have the same width has a considerable impact on the cut-off energy, much
−
for all models; we did not represent them for clarity. morethaninthe caseofthe homogeneousmagneticfield
As already mentioned, the functions represented in this of Ref. [1]. This is due to the presence of the diffusive
figure are closely related to the diffusion coefficient D regime at the low energy tail. One can indeed observe
[see Eq. (7)], which is required to calculate the spectra in Fig. 6 the flat diffusive locus at low energies for mod-
in the diffusive regime. els 1 and 4. Phenomenologically, one understands that
for these models, a slight change in the closest source
In figure 7 we present the total spectra (galac- distance caninfluence greatlythe flux oflowenergypar-
tic+extragalactic)comparedto the data, for our param- ticles. (iv) Finally, comparing our plots for AGASA and
eter fit for each model. As for figure 5, we draw the HiRes data, we conclude that the fits are better for the
medianspectrum(dot-dashedline)obtainedover100re- latter. The higher slope above the second knee break
alisations of the source locations. The upper and lower point in the HiRes data as well as the gap of data be-
dotted curves show the 75th and 25th percentiles around tweenthe KASCADE andHiRes rangesmakethe fitting
this prediction, meaning that only 25% of spectra are easier.
higherorlowerrespectivelythanindicatedbythecurves.
This uncertainty is related to the location of the closest
sources. As explained in Ref. [1], we draw a straight
dashed line in the region slightly above 1018 eV, where One should emphasize, however, that the above fits
the propagation is neither rectilinear nor diffusive (see were obtained by hand, not by any optimization proce-
Ref. [1] for more details on this transition zone). dureduetothecomputingtimerequiredtocomputeone
The galactic cosmic ray component is modeled as fol- spectrum. Therefore, the spectra shown above do not
lows. Supernovae are accepted as standard acceleration strictly speaking represent the best fit to the data. Fur-
sites, yet it is notoriously difficult to explain acceler- thermore,oneshouldalsoexertsomecautionwhencom-
ation up to maximal energy 1018 eV. Thus it is as- paring datasets from different experiments. In Fig. 7,
sumed that the knee sets the maximal acceleration en- we chose to plot separately the AGASA and HiRes data
ergyfor galacticcosmic rays: inthis conservativemodel, because of the well-known discrepancy, but one cannot
the spectrum of species i with charge Z takes the form excludeadiscrepancybetweentheenergyscalesofKAS-
j (E) (E/E )−γiexp( E/E ), with γ 2.4 2.7, a CADE and HiRes for instance, which would shift one
Z Z Z i
species≃dependent spectra−l index, E Z∼ 2 1−015 eV dataset with respect to the other. Given all these un-
Z
≃ × ·
[38]. Thetotalgalacticcomponentisobtainedasthesum certainties, the fits shown in Fig. 7 appear satisfactory,
of elemental spectra, each adjusted to KASKADE data exceptmaybethatofmodel3inwhichthecut-offalways
as described in Ref [1]. appears too mild.