text:
# Magnetars in the Metagalaxy: An Origin for Ultra High Energy Cosmic Rays in the Nearby Universe

## 1 Introduction

I study the possibility that relativistic winds from rapidly rotating magnetars, neutron stars with surface dipole fields on the order of $`10^{15}`$ Gauss (Duncan and Thompson 1992, Paczynski 1992, Kouveliotou et al. 1998, 1999) create the highest energy cosmic rays (ultra-high energy cosmic rays, a.k.a. UHECR). I assume magnetars occur in all normal galaxies which form massive stars, with the UHECR ariving from outside our own Galaxy, except on the occasions (perhaps once per $`10^5`$ years) when a rapidly rotating magnetar is born in our own galaxy. This model has the virtue of having little difficulty in accelerating protons to energies in excess of $`10^{21}`$ eV, with the sources being in all normal (star forming) galaxies and a luminosity density entirely acceptable from the point of view of the (very approximately known) rate of formation of magnetars in our own galaxy, thus offering an explanation of the puzzling air showers with energies above $`10^{19.6}`$ eV without having to introduce major extrapolations of known or suspected energetics.

In §2, I summarize the data on UHECR, the loss processes that affect their transport through intergalactic space, the energetics of UHECR and aspects of the astronomy and physics of magnetized compact objects as known from studies of rotation powered pulsars that are relevant to the present investigation. §3 outlines the calculation of the particle injection spectrum from a newly born magnetar and the effects intergalactic transport have in altering the injected spectrum to the spectrum received at the Earth. In §4.1, I discuss the escape of the relativistic wind from a magnetar’s natal supernova, and in §4.2 I outline how the wind drives a relativistic blast wave containing a Magnetar Wind Nebula (MWN) into the surrounding interstellar medium. The accceleration mechanism of the UHE ions in the magnetar’s wind is addressed in §4.3, and their escape from the expanding MWN nebula is outlined in §4.4. The effect the electromagnetic energy lost from a magnetar has in creating HI supershells in the interstellar medium is discussed in §5, with results that are used to suggest that more than 90% of a newly born, rapidly rotating magnetar’s $`5\times 10^{52}`$ ergs of rotational energy gets lost as gravitational radiation. The typical strains of such gravity waves, from relatively nearby events that could contribute particles in the UHECR spectrum with energies above $`10^{20}`$ eV, are estimated in §6. The possibility of observing bursts of UHECR associated with the birth of individual magnetars gets attention in §7. I discuss the relation of this study to other work in §8. I draw my conclusions in §9.

The most prominent difference between this model for UHECR and most others which appeal to sources more or less uniformly distributed throughout the Universe is that the underlying objects create an extremely flat injection spectrum, with particles per unit energy range injected at a rate $`E^1[1+(E/E_g)^s]^1`$, with $`E_g`$ reflecting the strength of gravity wave emission on the magnetars’ early spin down; for the model developed in detail here, $`s=1`$. As a result, the usual Greisen-Zatsepin-Kuzmin (GZK) cutoff disappears. This cutoff, apparent in models that assume power injection spectra $`E^s,s2`$ (see Berezinsky et al. 2002 for a recent example of such “conservative” models), becomes a flattening of the observed $`E^3J(E)`$ spectrum between $`10^{19}`$ and $`10^{20}`$ eV; the flatttening becomes a moderate replica of a GZK cutoff if the gravitational wave losses are strong enough.

## 2 Background

Recent observations of ultra-high energy cosmic rays (UHECR), manifested as atmospheric air showers with energy above $`10^{18.8}`$ eV, suggest an origin in the metagalaxy (Bird et al. 1995, Yoshida et al. 1995, Takeda et al. 1999, Abu-Zayyad et al. 2002a,b); see Nagano and Watson (2000), Watson (2002) and Halzen and Hooper (2002) for recent reviews. This interpretation relies on three basic facts.

1) The total cosmic ray energy spectrum hardens above $`E_a10^{18.8}`$ eV, with isotropic intensity $`J(E)E^{3.2\pm 0.05}`$ below $`E_a`$ to $`JE^{2.8\pm 0.2}`$ for $`E>E_a`$, a hardening of the spectrum known as the cosmic ray “ankle.” See Nagano and Watson (2000) for a thorough discussion of the strong evidence for the spectral hardening at the highest energies, including the uncertainties in the value of $`E_a`$. Such hardening of the spectrum indicates an origin for the high energy particles causing the showers at energies above $`E_a`$ differing from the lower energy particle sources. Numerically, for $`E>E_a`$

$$J(E)10^{36}\left(\frac{E_a}{E}\right)^{2.8\pm 0.2}\mathrm{eV}^1\mathrm{cm}^2\mathrm{s}^1\mathrm{ster}^1.$$
(1)
This flux yields the integral intensity above $`E_a`$, $`J(>E_a)3\times 10^{18}`$ cm-2-sec-1-ster-1, corrsponding to particle number density $`n(>E_a)1.5\times 10^{27}`$ cm-3 and energy density $`U(>E_a)2\times 10^8`$ eV cm-3.

2) Correlations of event directions on the sky with the supergalactic equator may exist at energies above $`10^{19.5}`$ eV (Uchiori et al. 2000), although the results are also marginally consistent with isotropy, as are other analyses (Takeda et al. 1999 and references therein.) An extragalactic origin remains consistent with either interpretation of the isotropy investigations. The few events with energy above $`10^{20}`$ eV say nothing useful about anisotropy.

3) The shower data are consistent with light nuclei or protons as the initiating events above $`E_a`$. At lower energies, in the cosmic ray “shin”, lying between the “knee” at $`10^{15}`$ eV and the ankle, the data suggest a heavy nuclei composition (Bird et al. 1994).

If light charged particles such as protons initiate the UHE showers, the Larmor radii of such particles in the interstellar magnetic field, $`r_L=E/ZeB=100(E/10^{20}\mathrm{e}V)(10^6\mathrm{Gauss}/ZB)`$ kpc, exceed the radial extent of the gaseous disk of our galaxy (for $`Z=1`$), and greatly exceed the vertical thickness of the disk - indeed, $`r_L(100\mathrm{EeV})`$ greatly exceeds the thickness of the galactic synchrotron halo. Even if $`Z=26`$, as in the galactic pulsar source model of Blasi et al. (2000), the Larmor radius of 100 EeV particles still exceeds the thickness of the known magnetized disk (the synchrotron halo, with scale height $`1`$ kpc), implying a substantial reduction of the flux estimated by Blasi et al.111From what is known of the top layers of a pulsar’s atmosphere from X-ray observations, the composition almost certainly is not iron. Helium may be a more probable candidate, based on interpretaion of neutron stars’ X-ray spectra \[including line features in one case, Sanwal et al. 2002 and references therein, or perhaps O and Ne (Hailey and Mori 2002).\] Lighter elements would have floated up to the X-ray photosphere (Zavlin et al. 1998, Pavlov and Zavlin 2000). Such high altitude ions (perhaps a few centimeters above the crust) are the reservoir for extraction from the star by the electric fields parallel to the almost vertical $`B`$ fields, to be injected into the escaping wind, not the iron thought to dominate the crustal composition. Therefore, a charged particle model of the UHECR based on acceleration driven by neutron stars more plausibly should accelerate protons or $`\alpha `$ particles, as indicated by the shower structures at energies above $`E_a`$, which in turn requires consideration of an extragalactic model, to be consistent with the isotropy data. The observed isotropy also argues against the UHECR being galactic or local group neutrons, which can survive for a distance $`9(E/10^{18}\mathrm{eV})`$ kpc from their sources; neutrons at 100 EeV certainly cannot arrive from the whole local supercluster.

Extragalactic ultra-high energy charged particles cannot arrive from truly cosmological distances. As first noted by Greisen (1966) and by Zatsepin and Kuz’min (1966), interaction of relativistic charged particles with the cosmic microwave background (CMB) degrades their energy at a rapid rate through production of $`\pi `$ mesons when the CMB photons have sufficient energy in the rest frame of the cosmic ray. The attenuation length for protons due to pion losses (Protheroe and Johnson 1996, Stanev et al. 2000, Fodor and Katz 2000, Berezinsky et al. 2002) can be adequately represented as

$`l_{GZK}`$ $`=`$ $`c{\displaystyle \frac{E}{\dot{E}}}=cT_{GZK}\left[1+\left({\displaystyle \frac{E_G}{E}}\right)^2\right],E_G=2.8\times 10^{20}\mathrm{eV},`$ (2)
$`D_{GZK}`$ $``$ $`cT_{GZK}=20\mathrm{Mpc},T_{GZK}=6.5\times 10^7\mathrm{years}.`$ (3)
Expression (2) gives the largest loss rate above $`4\times 10^{19}`$ eV. Below this energy, photopion losses are weaker than those due to $`e^\pm `$ pair creation and to cosmological expansion. For the purposes of the present study, pair creation and cosmic sxpansion can be adequately represented as having constant attenuation lengths222More accurately, $`l_\pm =500[1+(E_\pm /E)^{0.4}]`$ Mpc for $`EE_\pm =4\times 10^{19}`$ eV, but the slow increase in path length with decreasing energy has no effect on the observable aspects of the model outlined here, at its present level of development., $`l_\pm =1000`$ Mpc, $`T_\pm =3.3`$ Gyr, and $`l_H=cT_H=4300`$ Mpc, $`T_H=14`$ Gyr. Therefore, the average loss rate is

$$\frac{\dot{E}}{E}=\frac{1}{T_{loss}(E)}=\frac{1}{T_H}\frac{1}{T_\pm }\frac{1}{T_{GZK}\left[1+\left(\frac{E_G}{E}\right)^2\right]}.$$
(4)
Many hypotheses have appeared concerning the origin of these air showers. The most puzzling aspect of the problem has been the observation of air showers at energies so large (in excess of $`10^{20}`$ eV) that charged particles coming directly from the accelerator must have started no more than 50 Mpc away, due to the radiation losses in the CMB. Within this local patch of the Universe, there is a lack of the most commonly hypothesized accelerators, disks and jets associated with giant black holes in the most active galactic nuclei, an idea with a long history going back to Lynden-Bell (1969).

This scarcity has led to a plethora of explanations. Some are of the “bottom up” variety as in this paper, based on acceleration of normal charged particles to extreme energies, usually in objects not otherwise known to exhibit the ability to accelerate extremely high energy particles (e.g., Torres et al. 2002 and references therein). Others are of the “top-down” variety, based on the behavior of particles left over from the early Universe - see Sarkhar (2002) for a recent review of such ideas. The Z-burst model (Singh and Ma 2003 and references therein) offers a mixed alternative, with extremely high energy charged particle acceleration and ultra high energy neutrino production moved to high redshift (where hyper-active AGN are more common), while the high energy protons we observe are created locally through interactions between the UHE neutrinos and low energy neutrinos left over from the Big Bang (clustered in dark matter halos, if the background neutrinos have sufficient mass). The main point of the model proposed here is to argue that adequate accelerators may be common in the local Universe, since they may be in all normal galaxies with active star formation. In particular, the accelerators may be associated with magnetars, a relatively rare form of core-collapse supernova debris.

The required galactic luminosity in UHECR gives a strong hint that the UHECR sources have something to do with massive stars and supernovae in ordinary galaxies. This luminosity is, using (1) and (4),

$$L_g^{UHE}=\frac{4\pi }{cn_g}_{E_a}^{3\times 10^{20}\mathrm{eV}}𝑑EE\frac{J(E)}{T_{loss}(E)}=0.8\times 10^{39}n_{g2}^1\mathrm{ergs}/\mathrm{sec}.$$
(5)
Here $`n_g=0.02n_{g2}`$/Mpc3 is the space density of galaxies (Blanton et al. 2001), the majority of which actively form stars and are not in large clusters. The particles’ intergalactic residence time, weighted by the energy in the observed spectrum in (5), is $`T_{loss}T_H/6.92.1`$ Gyr; when weighted by particle number, the mean residence time is $`T_H/3.8=3.7`$ Gyr. This luminosity per galaxy coresponds to the metagalactic luminosity density $`5\times 10^{44}`$ ergs/Mpc3-year. These energetics are consistent with a supernova related source, as is the case with normal cosmic rays. The high energy per particle forms the distinctive peculiarity of the UHECR.

The basic model for the accelerator proposed here is a variant of that advanced by Blasi et al. (2000), who assumed the relativistic winds from pulsars in our own galaxy could be the acceleration site, granted that these stars begin their lives as almost maximally rapid rotators. Relativistic magnetic rotators have magnetospheric voltage drops across the magnetic field with magnitude

$$\mathrm{\Phi }_{\mathrm{mag}}=\frac{\mathrm{\Omega }^2\mu }{c^2}=\mathrm{\Psi }_{\mathrm{mag}}\frac{\mathrm{\Omega }}{c}=3\times 10^{22}\frac{\mu }{10^{33}\mathrm{cgs}}\left(\frac{\mathrm{\Omega }}{10^4\mathrm{s}^1}\right)^2\mathrm{Volts},$$
(6)
which is in excess of the energy/charge exhibited in the UHECR, if the underlying star has large angular velocity ($`\mathrm{\Omega }>10^3`$ s-1) and large dipole moment ($`\mu >10^{32}`$ cgs, corresponding to a surface dipole field above $`10^{14}`$ Gauss)333Indeed, a neutron star which had magnetic energy in the poloidal field comparable to the gravitational binding energy of the star ($`B_{\mathrm{dipole}}10^{18}`$ Gauss) would be able to accelerate charged particles to $`10^{25}`$ eV, if it rotated close to breakup ($`\mathrm{\Omega }>10^4\mathrm{sec}^1`$). Of course, such an object would look nothing like a sphere and would lose most of its rotational energy in gravitational waves (see below), should it be so fortunate to survive the host of instabilities plaguing its existence. Magnetized fallback disks (e.g., Macfadyen and Woosley 1999) might be functional equivalents of magnetically supported, rapidly rotating neutron stars; see §8.. Here $`\mathrm{\Psi }_{\mathrm{mag}}=R_L^2B(R_L)`$ is the magnetic flux in the “open” field that extends beyond the light cylinder, located at $`r=R_Lc/\mathrm{\Omega }=30/\mathrm{\Omega }_4`$ km that connects the magnetosphere to the outside world, and $`\mathrm{\Omega }=10^4\mathrm{\Omega }_4\mathrm{s}^1`$.

Acceleration in the wind avoids the catastrophic radiation losses associated with acceleration inside a rapidly rotating neutron star’s magnetosphere. Charged particles moving on paths with radius of curvature $`\rho _c`$ have radiative power $`P=(2/3)(q^2c/\rho _c^2)(E/mc^2)^4`$. Radiation losses would therefore prevent protons from accelerating to energies above
$`E_{rad}10^{16.5}(A^2/Z)^{1/3}\mathrm{\Omega }_4^{1/3}(r/R_L)^{1/3}(\rho _c/r)^{2/3},`$ eV, where a particle’s mass is $`Am_p`$ and its charge is $`Ze`$. Within and near the magnetosphere, the strong magnetic field guides particles along $`B`$. Therefore, even if strong parallel electric fields can form, the radius of curvature of a particle’s orbit is $`\rho _crR_L`$ and all the energy expended in accelerating the putative cosmic ray gets lost as gamma rays (and $`e^\pm `$ pairs.)

At radii $`rR_L`$, the physical situation differs. I shall argue that at large radii, acceleration can and does occur in an angular sector around a magnetar’s rotational equator, where the radii of curvature of highly accelerated ions’ orbits do greatly exceed the local radius $`r`$ at each point on their orbits (see §§4.3 and 4.4), so that radiative losses are negligible.

Relativistic magnetic rotators, appearing as rotation powered neutron stars and as black holes with magnetized disks losing rotational energy in jets, generally create non-radiative relativistic outflows. These outflows - relativistic “winds” - include wound up magnetic fields and outflowing plasma, which support the “MHD” electric field $`𝐄_{}=(1/c)𝐯\times 𝐁`$. The magnetic field has magnitude $`BB(R_L)(R_L/r)`$, an estimate appropriate both to quasi-spherically expanding winds and to (current-carrying) jets. At radii large compared to $`R_L`$ , the cross field voltage in the wind is ($`v/c1`$)

$$\mathrm{\Phi }_{wind}rE_{}(r)=rB(r)=R_LB(R_L)=\frac{\mathrm{\Psi }_{\mathrm{mag}}}{R_L}=\mathrm{\Phi }_{\mathrm{mag}}=\frac{\mu }{R_L^2}=\frac{\mu \mathrm{\Omega }^2}{c^2}$$
(7)
The last two forms of (7) are appropriate for underlying rotators which have magnetic dipole moments, such as neutron stars, while expressing $`\mathrm{\Phi }_{wind}`$ in terms of the open magnetic flux $`\mathrm{\Psi }`$ allows application of the same physical considerations to outflows from the magnetized disks thought to orbit black holes.

The estimate given by equation (7) assumes the relativistic wind does not suffer strong radiative dissipation, which is true in the outflows from rotation powered pulsars, and is also the case in electromagnetic models of jets from underluminous galactic nuclei in radio galaxies, galactic microquasars and gamma ray bursts (GRBs). Then the voltage generated by the magnetized rotator is available to accelerate particles at radii $`rc/\mathrm{\Omega }`$ \- for protons, particles of energy greater than $`10^{20.3}`$ eV could be accelerated at distances $`r>0.3(r/\rho _c)^2`$ pc, a dimension easily achievable in the winds from magnetized rotators. If the acceleration process approximates that of a linear accelerator, rather than of a synchrotron, as may be the case if “surf-riding” acceleration (e.g., Tajima and Dawson 1979, Chen et al. 2002, §4.3) applies, then $`r/\rho _c1`$ and radiation losses are even less of a constraint on the model (see §4.4).

The new aspect of the idea pursued here is to apply the same physics to magnetars in all galaxies. These objects have no difficulty (in principle) with accelerating light nuclei to the required energies, and by hypothesis, being located in all galaxies, have no difficulty in supplying particles with an isotropic distribution of source directions. The initial spin rates of magnetars are unconstrained, except from theoretical considerations of the stars’ dynamical stability, which limits the rotational energy to be less than 25% of their gravitational energy ($`\mathrm{\Omega }<1.2\times 10^4`$ s-1 for a uniform sphere).

Therefore, one can assume that the stars are born with sub-millisecond initial rotation periods. The limited evidence on the initial spins of rotation powered galactic pulsars suggests some begin their lives rotating at less than their maximal rate, if electromagnetic torques dominate their spindown from birth (e.g., Murray et al. 2002, Kaspi and Helfand 2002). As I will show, other constraints suggest that at least for magnetars, gravitational wave torques may govern the initial spindown, which alleviates concerns over whether all such objects might be born as slow rotators - in any case, only 5-10% of the magnetars must be born as rapid rotators, for the UHECR source model proposed here to be viable.

Magnetars are of particular interest as an acceleration site since a) they exist (Kouveliotou et al. 1998, 1999; see the summary in Baring and Harding’s 2001 Table 1), and b) if they rotate fast enough sometime in their lives, they easily can have voltage drops well in excess of the energy per unit charge observed so far in UHECR (expression 6). $`\mathrm{\Omega }_4=1`$ corresponds to a rotation period of 0.63 msec, about half the centrifugal breakup angular velocity for a $`1.4M_{}`$ neutron star with 10 km radius, while $`\mu _{33}=1`$ corresponds to a polar surface dipole field of $`2\times 10^{15}\mu _{33}R_{10}^3`$ Gauss, with $`R_{10}=R_{}/10`$ km; these values corrspond to a charateristic voltage of $`3\times 10^{22}`$ V. Clearly, rapidly rotating, magnetized neutron stars can provide the energy/particle seen in UHECR, in principle, so long as the wind carrying the voltage survives to distances far from the star.

The energetics of a magnetar acceleration model are simple. Suppose the particles each gain the energy

$$E(\mathrm{\Omega })=q\eta \mathrm{\Phi }_{\mathrm{mag}}=q\eta \frac{\mathrm{\Omega }^2\mu }{c^2}=3\times 10^{21}Z\eta _1\mathrm{\Omega }_4^2\mu _{33}\mathrm{eV},\eta _1\eta /0.1.$$
(8)
Here $`\eta `$ is the fraction of the open field line voltage experienced by each particle on its way from the star to the outside world. If the magnetosphere is relativistic, the usual electromagnetic rate of rotational energy loss

$$\dot{E}_{EM}=\frac{4}{9}\frac{\mathrm{\Omega }_{}^4\mu ^2}{c^3}$$
(9)
applies444Expression (9) is the standard result for vacuum dipole radiation, after averaging the geometric factor $`(2/3)\mathrm{sin}^2i`$ over the sphere; $`i`$ is the angle between the dipole moment and the rotation axis. All electromagnetic spindown theories are thought to lead to comparable rates of energy loss, with $`\dot{E}_{EM}=c\mathrm{\Phi }_{\mathrm{mag}}^2`$..

I adopt the view that electric conduction currents are a major part of the support for the electromagnetic fields which exert torques on the rotating objects. I also adopt the assumption that charged particle beams extracted from the rotating object by electric fields compose these currents; see Arons (2002) for a recent discussion of these currents’ significance, and a discussion (briefly recapitulated below) for the evidence that such beams actually exist.

If conduction currents are the sole source of the electromagnetic fields, and if a charged particle beam carries all of the current, $`\dot{E}_R=I\mathrm{\Phi }_{\mathrm{mag}}=qc\dot{N}\mathrm{\Phi }_{\mathrm{mag}}`$, and

$`\dot{N}`$ $`=`$ $`c\left[{\displaystyle \frac{\mathrm{\Omega }B_{,\mathrm{dipole}}}{2\pi c}}\right](2\pi A_{\mathrm{cap}}),`$ (10)
$`=`$ $`{\displaystyle \frac{\mathrm{\Omega }^2\mu }{|q|c}}={\displaystyle \frac{\dot{E}_R}{|q|\mathrm{\Phi }_{\mathrm{mag}}}}`$ (11)
One recognizes $`A_{\mathrm{cap}}\pi (\mathrm{\Omega }R_{}/c)R_{}^2`$ to be the area of a magnetic polar cap, if the closed field lines extend to the light cylinder, and the expression in square brackets in (10) to be the Goldreich-Julian charge density $`\rho _{GJ}`$ at the stellar surface, the number of elementary charges per unit volume such that parallel electric fields are shorted out in the magnetosphere.

The electric current

$$I_{GJ}=2A_{cap}c\rho _{GJ}=\frac{\mathrm{\Omega }_{}^2\mu }{c}=c\mathrm{\Phi }_{\mathrm{mag}}$$
(12)
is a fundamental item in the theory of the spindown of relativistic magnetized rotators. Since $`\dot{E}_{EM}=I_{GJ}\mathrm{\Phi }_{\mathrm{mag}}`$, rotators which have such currents linking the star to the outside world can spin down at a rate comparable to that expected from vacuum (magnetic dipole) radiation losses. Both contributions to spindown are expected to be comparable, in oblique rotators - there is substantial evidence, from the lack of dependence of radio pulsars’ spindown torques on the angle between the rotation and dipole axes (Lyne and Manchester 1988) that this assumption is true.

Such currents, in stars with dipolar magnetic fields, include a primary charged particle beam emanating from the polar cap - electrons, when $`i\mathrm{}(\mu ,𝛀)<90^{}`$. That current must be balanced, at least on average, by an electric return current. The most traditional ideas assume this return current is an outflowing beam of charges of the opposite sign - ions, if $`i<90^{}`$ \- usually assumed to occur on the ring of field lines around the polar cap which form the boundary layer separating the closed and open field line zones of the magnetosphere (Goldreich and Julian 1969); in force free models (Michel 1975, Contopoulos et al. 1999), the return current appears as a current sheet. Such “auroral” field lines map to the rotational equator, in oblique as well as aligned rotators, where the return current forms a (corrugated) current sheet, with spatial oscillations imposed by the rotation of the underlying non-aligned dipole magnetic field. The model developed here for UHECR assumes the ion injection rate per magnetar to be that of this electrodynamically fundamental return current, based on an analogy to what has been learned from Pulsar Wind Nebulae.

$`\dot{E}_{EM}=I_{GJ}\mathrm{\Phi }_{open}`$ does not mean that such an object automatically is a maximal particle accelerator, with all the energy going into the energy of the particles in the electric currents. For that to be true, the particles in the currents would have to transport all of the rotational energy lost, and each would have have to reach the maximum energy per particle $`q\mathrm{\Phi }_{\mathrm{mag}}`$ \- other particle constituents would have to be largely absent, and Poynting fluxes would have to be unimportant. Pulsar theory, and related theories of magnetized disks, begins with the opposite assumption - that a large scale Poynting flux carries almost all of the outlowing energy, with the current carriers gaining very little of the potential energy available until dissipation occurs in the region where the outflow from the star terminates (Goldreich and Julian 1969). Even then, the energy dissipated might go into bulk expansion of the surroundings, rather than into acceleration of a few high energy particles. This approximation - represented in its most extreme form by force free models of the outflow (e.g., Contopoulos et al. 1999) - gets its emprirical support from the lack of intense photon emission from many known relativistic outflows. Therefore, one cannot blithely assume the existence of an energy/particle given by (8) with $`\eta 1`$, even when the ion loss rate is given by (11).

In addition, the study of pair creation in the magnetosphere has led to the belief that a dense electron-positron plasma provides the dominant constituent by number in the outflow. The feeding of the X-ray source in the Crab Nebula, in the X-ray nebula around the Vela pulsar (Helfand et al. 2001) and in the nebula G320.4 around PSR 1509-58 (Gaensler et al. 2002) strongly supports the existence of such a dense plasma outflow - if the flow leaving the neutron star has the same composition on all field lines, the particle input to the Crab ($`\dot{N}_{\mathrm{total}}10^{38}`$ pairs/second) suggests $`\eta _{pairs}<10^3`$ \- the pairs carry much of the energy flow at large distances from the star, but being a dense plasma, they have energy per particle well less than (8).

However, study of the structure of the shock wave terminating the equatorial relativistic wind in the Crab (e.g. Hester et al. 1995), in G320.4 (Gaensler et al. 2002) and perhaps in the Vela nebula (Helfand et al. 2001), where the energy in the pulsar wind gets transfered to the plasma forming the observed synchrotron nebulae, suggests that an interestingly large fraction of the total voltage does get applied to a Goldreich-Julian flux of heavy ions on the equatorial flow lines from rotation powered pulsars (Arons 2002). These are probably light nuclei such as protons or alpha particles (Gallant and Arons 1994). Detailed modeling of this interaction zone (Gallant and Arons 1994, Spitkovsky and Arons 2000 and submitted to ApJ) suggests that the total ion injection rate is indeed the Goldreich-Julian rate (11) and that these these particles, flowing out in a latitudinal sector around the stars’ rotational equators, have experienced $``$ 10% of the full potential drop (i.e., $`\eta 0.1`$), and thus carry an energy loss from the rotator $`\dot{E}_{ions}=I_{GJ}\eta \mathrm{\Phi }_{mag}=\eta \dot{E}_R0.1\dot{E}_R`$ of each neutron star’s rotational energy loss in an equatorial sector filling $`1/5`$ of the sky around each of these pulsars, exactly the region where the models hypothesize the flow of the return current in an equatorial current sheet. Outside this equatorial sector, the wind is more likely to be Poynting flux (AC and DC) dominated.

The specifics of the ion accelerator are not known - most recent suggestions for the acceleration site have focused on non-ideal processes in the equatorial wind, well outside the light cylinder. The ions enter the wind, having been drawn up from the star along the boundary layer separating the closed and open field lines, leave the star in the plane of the star’s rotational equator, and gain their high energy per particle at some large distance from the star, perhaps at the interesting radius where the shortage of charge carriers causes dissipation of the corrugations in the equatorial current sheet (Michel 1971, Coroniti 1990, Michel 1994), or transform to electromagnetic waves, either Alfvenic (Bellan 1999, 2001) or vacuum-like (Melatos and Melrose 1996; Melatos 1998), whose pomderomotive force can accelerate the particles, through surf riding in the EM fields, to energies comparable to $`q\mathrm{\Phi }_{mag}`$. I discuss surf riding acceleration mechanism further in §4.3, but a detailed study of that physics is outside the scope of the present investigation.

In what follows, I assume that such ion acceleration also occurs in magnetars that were rapidly rotating at birth, and explore the consequences of this assumption for an extragalactic origin of UHECRs from normal galaxies. I assume all galaxies which have core collapse supernovae form magnetars with a birthrate $`\nu _m`$ comparable to that in our own galaxy, where $`10^3>\nu _m>10^5`$ years-1 (Gaensler et al. 2001), and that some are born with high angular velocities, $`\mathrm{\Omega }_i10^4`$ s-1. During their spindown, they emit heavy ions (taken to be protons in what follows) with an ion loss rate per magnetar given by (11), with the energy/particle in the emitted equatorial “beam” declining as the rotator ages, as specified by (8), with $`\eta 0.1`$ in analogy to what has been inferred for pulsars. Therefore, the highest energy particles come from new born objects. The particles cannot be contained within normal galaxies, therefore the flux at the earth reflects the contributions from all galaxies within the energy dependent GZK volume. That volume expands to encompass all galaxies in the Hubble volume, as the particle energy decreases, but below the ankle energy $`E_a`$, this extragalactic source becomes swamped by the softer particle spectrum thought to be due to acceleration in our own galaxy.

A simple estimate shows that the idea is a viable candidate for the origin of the UHE spectrum. The star emits ions at the rate (11). Suppose electromagnetic torques are the sole means of spinning down the neutron star. The strong $`\mathrm{\Omega }`$ dependence says that the star emits most of the particles during the initial loss of rotational energy, when $`\mathrm{\Omega }=\mathrm{\Omega }_i=10^4\mathrm{\Omega }_4`$. The initial electromagnetic spindown time is

$$\tau _{EM}=\frac{\frac{1}{2}I\mathrm{\Omega }_i^2}{\dot{E}_{EM}}=\frac{9}{8}\frac{Ic^3}{\mu ^2\mathrm{\Omega }_i^2}=5\frac{I_{45}}{\mu _{33}^2\mathrm{\Omega }_4^2}\mathrm{minutes}.$$
(13)
Therefore, each magnetar injects

$$N_i\dot{N}_{GJ}(\mathrm{\Omega }_i)\tau _{EM}(\mathrm{\Omega }_i)=\frac{Ic^2}{Ze\mu }=4\times 10^{42}\frac{I_{45}}{Z\mu _{33}}$$
(14)
ions into the metagalaxy, each with the energy (8). The existing observations of a few events at energies $`>10^{20}`$ eV require the birth of some magnetars with $`\eta _1\mathrm{\Omega }_4^2\mu _{33}>0.1/Z`$.

With a magnetar birth rate per galaxy $`\nu _m=10^4\nu _{m4}\mathrm{yr}^1`$ and a galaxy density $`n_g=0.02n_{g2}`$ Mpc-3, the number of particles per unit volume from magnetars in all the galaxies is $`n_{UHE,magnetars}^{(est)}=N_in_g\nu _mT_{loss}(E_a),`$ where $`T_{loss}(E_a)`$ is the residence time of a particle at the ankle initially injected at higher energy; very roughly, this is the lifetime $`T_\pm (E_a)10`$ Gyr. Then

$$J^{(est)}(>E_a)\frac{c}{4\pi }N_i\nu _mn_gT_\pm (E_a)7\times 10^{18}\left(\frac{I_{45}}{Z\mu _{33}}\right)\nu _{m4}\mathrm{cm}^2\mathrm{sec}^1\mathrm{ster}^1,$$
(15)
comfortably above the observed $`J(>E_a)=3\times 10^{18}\mathrm{cm}^2\mathrm{sec}^1\mathrm{ster}^1`$.

The required luminosity in UHE particles does not unduly tax the energetics of rapidly rotating magnetars in galaxies, as can be seen by writing $`L_g^{UHE}=ϵ_m\nu _mE_R^{(i)}`$. With $`L_g^{UHE}`$ from (5) and $`E_R^{(i)}=(1/2)I\mathrm{\Omega }_i^2=5\times 10^{52}I_{45}\mathrm{\Omega }_4^2`$ ergs, the data require a fraction $`ϵ_m=0.005/\nu _{m4}n_{g2}I_{45}\mathrm{\Omega }_4^2`$ of the initial rotational energy to go into UHE particle acceleration. I show in §5, using the results of the more detailed theory in §3, that the birth rate of fast magnetars, those with initial voltages large enough to accelerate particles with energy greater than $`E_a`$, is $`\nu _m^{(fast)}=(0.050.1)\nu _m`$. Then the coefficient in $`ϵ_m`$ changes from 0.005 to .05-0.1, still not a taxing demand on the model’s energetics. Put slightly differently, the luminosity density in rapidly rotating magnetars deduced from the fit of the model to the observations described in §3 is $`\nu _m^{(fast)}n_gE_R^{(i)}=10^{46}(\nu _m^{(fast)}/10^5\mathrm{yr}^1)I_{45}\mathrm{\Omega }_4^2\mathrm{ergs}/\mathrm{Mpc}^3\mathrm{year}`$. Given that 10% of the accelerating voltage gets applied to each accelerated particle, and that the injection rate per object is given by (11), the model has the required luminosity density $`10^{45}`$ ergs/Mpc3-year in UHECRs.

As I shall show, more careful treatment of the spectrum, plus approximate accounting for the losses sustained as the magnetized wind bursts out of the supernova that gives rise to the neutron star, improves the correspondence between the model and the data, and makes some interesting predictions for future UHECR observations.

## 3 The Spectrum of UHE Cosmic Rays from Bare Metagalactic Magnetars

### 3.1 The Injection Spectrum from an Isolated Magnetar

During a time $`dt`$, the magnetar spins down from angular velocity $`\mathrm{\Omega }(t)`$ to $`\mathrm{\Omega }(t+dt)=\mathrm{\Omega }(t)+\dot{\mathrm{\Omega }}dt`$. According to (8) and (11), the energy/particle $`E(t)=m_ic^2\gamma (t)`$ and the instantaneous injection rate $`\dot{N_i}`$ each decay $`\mathrm{\Omega }^2`$. Therefore, the number of particles accelerated as the star spins down in time $`dt`$ with energy between $`mc^2\gamma `$ and $`mc^2(\gamma d\gamma )`$ is

$$dN_i=\dot{N}_i[\mathrm{\Omega }(t)]dt=\dot{N}_i\frac{dt}{d\mathrm{\Omega }}\frac{d\mathrm{\Omega }}{d\gamma }d\gamma ,$$
(16)
with

$$\frac{d\mathrm{\Omega }}{dt}\dot{\mathrm{\Omega }}=\frac{\dot{E}_R}{I\mathrm{\Omega }}=\frac{\dot{E}_{EM}+\dot{E}_{grav}}{I\mathrm{\Omega }}.$$
(17)
Here $`I`$ is the moment of inertia, and (9) expresses the electromagnetic energy losses.

The gravitational wave losses may be substantial. If the star has a static equatorial ellipticity $`ϵ`$, gravitational waves extract energy at the rate (Ostriker and Gunn 1969)

$$\dot{E}_{grav}=\frac{32}{5}\frac{GI^2ϵ^2\mathrm{\Omega }^6}{c^5}=1.8\times 10^{51}I_{45}^2\mathrm{\Omega }_4^6ϵ_2^2\mathrm{ergs}/\mathrm{sec},$$
(18)
where $`I_{45}=I/10^{45}`$ g-cm2 and $`ϵ_2=ϵ/10^2`$. Note that $`ϵ_21`$ is a substantial but not ridiculous ellipticity; interior magnetic fields can create ellipticity (Usov 1992, Bonnazola and Gourgoulhon 1996) with magnitude (Ostriker and Gunn 1969)

$$ϵ=\frac{25}{24}\frac{R^4}{GM^2}(3B_{poloidal}^2B_{toroidal}^2)10^2\frac{3B_{poloidal}^2B_{toroidal}^2}{(4\times 10^{16}\mathrm{Gauss})^2};$$
(19)
the numerical value assumes $`M=1.4M_{}`$. The initial spindown time due to gravitational radiation is

$$\tau _{GW}=\frac{\frac{1}{2}I\mathrm{\Omega }_i^2}{\dot{E}_{grav}}=\frac{28.5}{I_{45}ϵ_2^2\mathrm{\Omega }_4^4}\mathrm{seconds}.$$
(20)
Interior fields an order of magnitude stronger than the surface field are by no means incredible. Formation of interior fields as large as $`10^{17.5}`$ Gauss in magnetars has been proposed, as a consequence of dynamo activity in the first few seconds of the proto-neutron star’s life (Duncan and Thompson 1992, Wheeler et al. 2000.) - these stars are thought to have non-dipolar surface fields much stronger than the already large dipole field. Once formed, such an ellipticity would most likely persist throughout the neutron star’s spindown, certainly through the first day or so of the spindown germane to the model discussed here.

Other sources of ellipticity or its equivalent may be dynamical instability of the star, if the rotational energy exceeds 25% of the newly formed star’s gravitational energy, and secular instabilites of F- and R-modes, driven by the gravitational radiation itself (e.g., Lai 2001, Ushomirsky 2001). In the interest of simplicity, I use (19) as motivation to choose $`ϵ_21`$ as a parameter throughout the subsequent discussion, and neglect other versions of gravitational wave spindown. In the appendix, I briefly summarize the possible effects of gravity waves emitted because of the R mode instability.

I rewrite the spin down rate as

$$\dot{\mathrm{\Omega }}=\frac{4}{9}\frac{\mu ^2\mathrm{\Omega }^3}{Ic^3}\left[1+\left(\frac{\mathrm{\Omega }}{\mathrm{\Omega }_g}\right)^2\right]$$
(21)
with $`\mathrm{\Omega }_g`$ the angular velocity at which gravity wave and electromagnetic losses are equal,

$$\mathrm{\Omega }_g\left(\frac{5}{72}\frac{c^2\mu ^2}{GI^2ϵ^2}\right)^{1/2}=3\times 10^3\frac{\mu _{33}}{I_{45}ϵ_2}\mathrm{s}^1.$$
(22)
Expressions (8), (11), (16) and (21) then yield the particle spectrum accelerated by a magnetar during its spindown to be

$$\frac{dN_i}{d\gamma }=\frac{dN_i}{dt}\left(\frac{dt}{d\mathrm{\Omega }}\right)mc^2\frac{d\mathrm{\Omega }}{dE}=\frac{9}{4}\frac{c^2I}{Ze\mu \gamma }\left(1+\frac{\gamma }{\gamma _g}\right)^1,$$
(23)
where

$$E_g=m_ic^2\gamma _g=\frac{Z\eta e\mu }{c^2}\mathrm{\Omega }_g^2=\frac{5}{72}\frac{Z\eta e\mu ^3}{GI^2ϵ^2}=3\times 10^{20}\frac{Z\eta _1\mu _{33}^3}{I_{45}^2ϵ_2^2}\mathrm{eV}.$$
(24)
The star promptly emits the energy that goes into the wind. The time to spin down to a voltage such that a particle of energy $`E`$ can be accelerated is

$`{\displaystyle \frac{t_{spin}(E)}{\tau _g}}`$ $`=`$ $`\left({\displaystyle \frac{\mathrm{\Omega }_g}{\mathrm{\Omega }}}\right)^2\left({\displaystyle \frac{\mathrm{\Omega }_g}{\mathrm{\Omega }_i}}\right)^2\mathrm{ln}\left[{\displaystyle \frac{\mathrm{\Omega }_i^2}{\mathrm{\Omega }^2}}{\displaystyle \frac{1+\frac{\mathrm{\Omega }^2}{\mathrm{\Omega }_g^2}}{1+\frac{\mathrm{\Omega }_i^2}{\mathrm{\Omega }_g^2}}}\right]`$ (25)
$`=`$ $`{\displaystyle \frac{E_g}{E}}{\displaystyle \frac{E_g}{E_i}}\mathrm{ln}\left({\displaystyle \frac{E_i}{E}}{\displaystyle \frac{1+\frac{E}{E_g}}{1+\frac{E_i}{E_g}}}\right),`$ (26)
where

$$\tau _g=\frac{\frac{1}{2}I\mathrm{\Omega }_g^2}{\frac{4}{9}\frac{\mu ^2\mathrm{\Omega }_g^4}{c^3}}=\frac{81}{5}\frac{GI^3cϵ^2}{\mu ^4}=0.9\frac{I_{45}^3ϵ_2^2}{\mu _{33}^4}\mathrm{hours}$$
(27)
is the electromagnetic spin-down time, for a star spinning at angular frequency $`\mathrm{\Omega }_g`$. (13) and (20) give the initial spin-down times for pure electromagnetic and pure gravity wave losses. Expression (25) yields the spindown history when gravity wave losses dominate, $`\mathrm{\Omega }(t)=\mathrm{\Omega }_i[1+(2t/\tau _{GW})]^{1/4}`$, if the limit $`\mathrm{\Omega }_g\mathrm{\Omega }<\mathrm{\Omega }_i`$ exists, and yields pure electromagnetic spindown, $`\mathrm{\Omega }(t)=\mathrm{\Omega }_i[1+(t/\tau _{EM})]^{1/2}`$, in the limit $`\mathrm{\Omega }\mathrm{\Omega }_g`$. For $`t<2\tau _g`$, (26) simplifies to

$$t_{spin}(E)=2\tau _{GW}\left[\left(\frac{\mathrm{\Omega }_i}{\mathrm{\Omega }}\right)^41\right]1.8\frac{Z^2\eta _1^2\mu _{33}^2}{I_{45}ϵ_2^2}\left(\frac{10^{20.5}\mathrm{eV}}{E}\right)^2\mathrm{hours},E_g<E<E_{max}.$$
(28)
Here $`E_{max}=Ze\eta \mathrm{\Phi }_{\mathrm{mag}}(\mathrm{\Omega }_i)=3.3\times 10^{21}Z\eta _1\mu _{33}\mathrm{\Omega }_4^2`$ eV. For $`t>2\tau _g,E<E_g`$, electromagnetic losses alone accurately describe the spindown. In this limit,

$$t_{spin}(E)=2\tau _{EM}\left(\frac{\mathrm{\Omega }_i^2}{\mathrm{\Omega }^2}1\right)1.8\frac{Z\eta _1I_{45}}{\mu _{33}^2}\frac{10^{20.5}\mathrm{eV}}{E}\mathrm{hours},E<E_g.$$
(29)
The general spin history can be written in the simplified form, valid for $`t\tau _{GW}`$,

$$\mathrm{\Omega }=\mathrm{\Omega }_g\left(\frac{2\tau _g}{t}\right)^p,$$
(30)
with $`p=1/4`$ for $`\tau _{GW}t<2\tau _g`$ and $`p=1/2`$ for $`t>2\tau _g`$.

Expression (23) is the spectrum that would be observed from a single nearby event by an experiment which accumulates events of energy $`E`$ and higher, over the time (26). An experiment which resolves a single nearby magnetar birth and spindown would see an instantaneous, approximately monoenergetic spectrum, with the energy/particle declining in proportion to $`\mathrm{\Omega }^2(t)`$.

To the magnetar’s host galaxy and to the metagalaxy, the energy input to the relativistic wind appears as an impulsive burst; for the energies of all the observed particles ($`E_{max}E>E_a=6\times 10^{18}`$ eV), $`t_{spin}(E)<t_{spin}(E_a)=42Z\eta _1I_{45}/\mu _{33}`$ hours. Note that this is the time taken to put energy into the wind, corresponding to voltages high enough to accelerate the observed ultra-high energy cosmic rays; the actual time of acceleration is later. How much later depends on the time it takes for a region of the wind to flow to the radius where acceleration occurs. In §4.3, I suggest this flow time probably is a few hours, for typical parameters.

### 3.2 Intergalactic Particle Transport and the Observed Spectrum from Bare Magentars

#### 3.2.1 Scatter Free Transport

The magnetars forming in galaxies emit UHE particles which cannot be contained within the galaxies. Therefore, normal galaxies inject particles into intergalactic space at the average rate per unit volume $`n_g\nu _mdN_i/d\gamma `$. Assuming scattering in intergalactic magnetic fields has negligible effect on particle transport, the intergalactic spectrum $`N(\gamma )`$ may be determined from

$$\frac{}{\gamma }(\dot{\gamma }N)=W_{geom}\nu _mn_g\frac{dN_i}{d\gamma }.$$
(31)
Here $`\dot{\gamma }=\gamma /T_{loss}(\gamma )`$ and $`W_{geom}=0.5`$ represents the geometric fact that only half the magnetars are expected to have the magnetic geometry appropriate to injecting ions from the stars atmospheres into the wind in the rotational equator (§4.4).

Integrating equation (31) with (4) and converting from $`N(\gamma )`$ to $`J(E)`$ through $`m_ic^2J=cN/4\pi `$ yields

$$J(E)=\frac{c}{4\pi m_ic^2}N(\gamma )=W_{geom}\frac{9}{16\pi }\frac{Ic^3}{Ze\mu }\frac{n_g\nu _mT_{loss}(E)}{E}\mathrm{ln}\left(\frac{E_{max}}{E}\frac{1+\frac{E}{E_g}}{1+\frac{E_{max}}{E_g}}\right).$$
(32)
The spectral form has the $`1/E`$ shape of the injection spectrum below 10 EeV, steepens to become proportional to $`E^3`$ in the observed UHE range $`10\mathrm{EeV}<E<100\mathrm{EeV}`$, then flattens to return to a shape approaching the $`E^1`$ spectrum above 300 EeV ($`E^2`$ if gravitational wave losses during spin down truncate the number of particles produced at energies well below $`E_{max}`$, i.e., when $`E_gEE_{max}`$), where GZK losses saturate. In the strong GR case, gravity wave emission shortens the time at which each magnetar can contribute particles at a given energy, which reduces the flux even in the better studied region $`E<5\times 10^{19}`$ eV, for parameters other than ellipticity fixed.

These ultra high energy particles have little trouble entering our galaxy555The ion Larmor radius $`r_L=E/ZeB=2.3(E/10^{18.8}\mathrm{eV})(3\mu \mathrm{Gauss}/ZB_0\mathrm{exp}(z/2h)`$ kpc, where $`z`$ is the vertical distance from the galactic midplane and $`B_0`$ is the midplane magnetic field, is larger than the vertical extent of the synchrotron halo of our galaxy, whose scale height is $`h`$ 1.8 kpc (Beuermann et al. 1985) implying easy entry even at the low end of the UHE spectrum. An extended, quasi-spherical magnetized halo, if it exists, probably has a magnetic field at least an order of magnitude weaker (with Larmor radius $``$ 20 kpc), an assumption supported by radio synchrotron emission observed in other disk galaxies (e.g, Hummel et al. 1991, Lisenfeld and Völk 2000). Thus even at the lowest observed UHECR energies, individual particles can reach our galaxy’s midplane following almost rectilinear orbits, in times not much exceding the light transit time - at the lowest energies, the guiding center drift velocity might drop a bit below c, to $`c(r_L/R_{halo})`$, if an extended halo with radius larger than 20 kpc exists. Since energy losses are negligible on this scale (and they would be negligible, even if the particles had to diffuse in the galaxy’s magnetic field), the particle density and spectrum inside our galaxy comes to equilibrium with the intergalactic spectrum at the same density and spectrum as is in the metagalaxy (Ginzburg and Syrovatskii 1964, Longair 1994). However, magnetic deflection in halo magnetic fields probably does corrupt the accuracy with which one can infer a possible source’s location from the look back direction along an air shower (itself known to $`1^{}`$ accuracy) for energies below $`5\times 10^{19}`$ eV.. Therefore, expression (32) can be compared directly to the observations. Multiplying (32) by $`E^3`$ and overlaying the results on the data yields the result shown in Figure 1. Given the substantial errors, fitting by eye suffices, with the conclusion that this elementary theory provides a satisfactory fit to the data for the values of $`K_0=W_{geom}\nu _mIn_gT_H/\mu `$ shown in Table 1.

At the relatively low energies now observed, the presence or absence of gravitational wave losses makes little difference unless the gravity wave losses are large. Below 4 EeV, this extragalactic particle spectrum hides below the flux of (probably) galactic particles at energies below the ankle.

The spectrum in Figure 1 represents an average over many particle bursts, each lasting a couple of days (for $`E>E_a`$). These bursts should be distinguishable as such, with sufficient collecting area. One can readily show, using the time dependent form of (31) with its right hand side set to zero and with (23) as the initial condition, that the fluence $``$ of a single burst has the same spectral form as is shown in Figure 1, but with the upper cutoff energy now determined by the losses in the CMB degrading the maximum injection energy to the value determined by the source’s distance.

Because the injection spectrum is extremely flat ($`E^1`$, if gravitational wave losses are negligible), the losses in the CMB flatten $`E^3J`$ but do not produce the GZK pileup characteristic of steeper injection spectra (e.g., Berezinsky et al. 2002.) If the claimed observation of such pile up by the Hi-Res experiment (Abu-Zayyad et al. 2002 a,b) is confirmed and the discrepancy with the AGASA results (Hayashida et al. 1994, Takeda et al. 1999; see http://www-akeno.icrr.u-tokyo.ac.jp/AGASA/ for results through Fall 2001) is resolved, the existence of such a feature would, within the context of the model desrcibed here, point to very rapid early spindown due to gravitational wave losses, which has the effect of steepening the highest energy part of the injection spectrum - see the “strong GR” curve in Figure 1.

### 3.3 Predicted Spectrum at Higher Energies

Unless current experiments are so fortunate as to have observed the upper cutoff, the spectrum of high energy particles should continue into the higher energy region where current detectors lack sensitivity. The specific model used here suggests the spectrum $`E^3J(E)`$ should resume rising with $`E`$ above the energy $`E_G=2.8\times 10^{20}`$ eV where the GZK loss rate becomes approximately energy independent (unless the gravitational wave losses are large, as shown in Figure 1.)

This prediction is simplified, in its assumption that all magnetars have exactly the same starting voltage ($`10^{22.5}`$ Volts) and single particle efficiency $`\eta 0.1`$. The inferred value of $`K_0`$ in Table 1 being small compared to unity in the weak and moderate GR cases has its most natural interpretation in the initial voltages (initial rotation periods, if all objects have the same dipole moment) having a distribution of values, with only a small fraction having voltages such that particles of energy greater than $`E_a`$ can be accelerated. Since the results of §4.4 yield $`W_{geom}=0.5`$, Table 1 suggests that the birth rate of magnetars with voltages high enough to contribute to the spectrum in Figure 1 is only 6% of the overall magnetar birth rate, in the “no GR” case (a “fast” magnetar formation rate of once every 170,000 years/galaxy); for “moderate GR”, the required rate is 10% of the overall magnetar birth rate (once every $`10^5`$ years/galaxy).

A distribution of initial voltages can have potentially observable effects on the predicted UHECR spectrum. To illustrate, suppose the distribution of birth rates as a function of starting voltage is a power law,

$$\frac{d\nu _m}{d\mathrm{\Phi }_i}=\frac{\nu _m}{\mathrm{\Phi }_{i,max}}\frac{s1}{\left(\frac{\mathrm{\Phi }_{i,max}}{\mathrm{\Phi }_{min}}\right)^{s1}1}\left(\frac{\mathrm{\Phi }_i}{\mathrm{\Phi }_{i,max}}\right)^s,$$
(33)
with $`\mathrm{\Phi }_{min}\mathrm{\Phi }_i\mathrm{\Phi }_{i,max}`$, $`\mathrm{\Phi }_{i,max}=3\times 10^{22}\mu _{33}(\mathrm{\Omega }_i/10^4\mathrm{sec}^1)^2`$ and $`\mathrm{\Phi }_{min}=10^{14}`$ Volts, the potential drop characteristic of the known magnetars - as with pulsars, initial rotation periods comparable to the observed periods are possible, so that initial rotation periods as long as 10 seconds must be considered. Expression (32), with $`d\nu _m/d\mathrm{\Phi }_i`$ from (33) replacing $`\nu _m`$, now gives $`E^3J(E;\mathrm{\Phi }_i)/\mathrm{\Phi }_i`$, with $`\gamma _{max}`$ replaced by $`\gamma _iZe\eta \mathrm{\Phi }_i/mc^2`$. Integrating over all initial voltages has the effect, when $`s>1`$, of replacing the logarithmic cutoff in (32) by $`(E/E_{i,max})^{(s1)}S(E/E_{i,max},s)`$, where $`S`$ is a slowly varying function which produces a more complicated logarithmic cutoff of the spectrum at $`E_{i,max}Ze\eta _1\mathrm{\Phi }_{i,max}`$. If $`s<1`$, there is no substantial change from the single initial voltage model - when the spectrum of birth potentials is flat ($`s<1`$), essentially all the stars would have starting voltage equal to $`\mathrm{\Phi }_{i,max}`$.

Fitting the resulting modified spectrum to the data (again by eye) yields the required reduction of the normalization if $`s1.2`$, in the “no GR” case. The extra $`E^{0.2}`$ now introduced in $`E^3J(E)`$ at energies below $`10^{20}`$ eV is undetectable in the current experiments, but with increased data above $`5\times 10^{19}`$ eV the shape of the spectrum will provide significant constraints on such refinements of the model. Above $`10^{20}`$ eV, the decreased number of contributing sources can substantially reduce the high peaks at high energy shown in Figure 1, by as much as a factor of 3. However, the basic conclusion, that the model predicts something to see, roughly at the level one would guess by extrapolating the observed $`E^3J(E)`$ above the formal GZK cutoff, remains robust. With $`s=1.2`$, 5% of all the magnetars have $`Ze\eta \mathrm{\Phi }_i`$ above $`E_a`$, consistent with the simple inference of the fast magnetar birthrate described above.

Also, the highest energy part of the spectrum probes acceleration in the wind at the earliest times. As discussed in §4.1, breakout of the wind from a magnetar’s natal supernova dissipates a fraction $`1W_{blowout}0.3`$ of the initial energy injected, with voltages in the freely expanding wind available for particle acceleration only after spindown from $`\mathrm{\Omega }_i`$ to $`\mathrm{\Omega }_iW_{blowout}^{1/4}`$, suggesting a more realistic upper cutoff for the spectrum to be at $`E_{max}^{(eff)}=E_{max}W_{blowout}^{1/2}=8\times 10^{21}\eta _1\mu _{33}\mathrm{\Omega }_4^2(W_{blowout}/0.7)^{1/2}`$ eV.

Figure 1 illustrates another robust conclusion, that the highest energy part of the spectrum is sensitive to possible spindown losses due to gravity waves. The shortening of the time in which particles can be accelerated to energies above $`E_g`$ suggests the interesting possibility that observations of UHECR at energies above $`10^{20}`$ eV can point the way to sources of gravitational radiation, and vice versa (see §6).

While the agreement of this simple model with the observations is encouraging, one must consider the possible modifications due to losses. Blowout from the source (in this case, the supernova initially containing the relativistic wind), radiative and adiabatic expansion losses are the most prominent pitfalls for any compact object model for high energy particles observed in any part of the cosmic ray spectrum. Likewise, modification of the spectrum due to possible diffusive transport in the intergalactic medium must be considered. I consider scattering in intergalactic magnetic fields first.

### 3.4 Effect of Particle Scattering in an Integalactic Magnetic Field

There is no positive information pointing to the existence of a general intergalactic magnetic field. Various limits suggest the magitude of such a field does not exceed $`10^9`$ Gauss (e.g., Barrow et al. 1997), although a number of authors have argued for larger fields injected from galaxies and spread around by ill-understood processes (e.g., Farrar and Piran 2000, Kronberg et al. 2001); see also Eilek and Owen (2002). If such a field exists, and if it causes high energy particles to flow diffusively through intergalactic space (Adams et al. 1997), significant alterations of the UHE spectrum received from sources in galaxies would occur.

Consider what happens in a magnetic field strong enough, and with sufficient magnetic fluctuations, to cause particles to have a scattering mean free path small compared to the GZK length at energy $`E`$. Assume the intergalactic field and its fluctuations fill all space (no intermittency), The mean free path for a particle of energy $`E`$ in the fluctuating field is (Kulsrud and Pearce 1969) $`\lambda _Br_L(E)(B/\delta B)^2,`$ where $`r_L(E)=E/ZeB=341E_{20.5}/ZB_9`$ Mpc, and $`\delta B^2k\delta B_k^2_{kr_L=1}.`$ Here $`E_{20.5}=E/10^{20.5}`$ eV and $`B_9=B/10^9`$ Gauss. If $`\lambda _B(E)l_{loss}`$ (see expression 2), particles come only from a distance $`D_{scatter}(E)\sqrt{\lambda _Bl_{loss}}<l_{loss}`$. This diffusive reduction in the observable volume can occur at the highest energies, including those above the saturation energy of the pion losses $`E_G=2.8\times 10^{20}`$ eV, only if the general intergalactic magnetic field is strong, $`B>2\times 10^8`$ Gauss, and the magnetic fluctuations are intense, $`\delta B/B1`$, conditions perhaps unlikely for a general intergalactic field.

But at lower energies, the volume contributing to the observed particle flux shrinks when scattering is significant, compared to that observed if scattering is negligible. This shrinkage depends on particle energy, which means that diffusive transport can have a significant effect on the observed particle spectrum between the ankle at $`E_a=6`$ EeV and $`E_B=1.1\times 10^{20}(ZB_9)^{1/3}(\delta B/B)^{2/3}`$ eV. The shape of the spectrum can be readily estimated in such a diffusive transport regime. The distance $`D(E)`$ from sources to the observer is $`D_{eff}=(\lambda _B(E)cT_{loss})^{1/2}=cT_{eff}<cT_{loss}`$. This contraction of the observable volume leads to the spectrum in (32) becoming $`J(E)(c/4\pi )T_{eff}n_g\nu _mdN_i/d(m_ic^2\gamma )`$ $`E^{3/2}`$ instead of $`E^3E^3JE^{3/2},`$ in the energy range where pion losses dominate, $`E>4\times 10^{19}`$ eV, instead of $`E^0`$, contrary to the observations. Solving the diffusion equation with systematic energy losses in the short mean free path regime confirms this estimate.

Therefore, if there is a general intergalactic magnetic field, this and all other metagalactic source models which explain the absence of the traditional GZK cutoff by the existence of relatively nearby hard spectrum particle sources can work only if $`E_B<E_a`$, which implies $`ZB_9(\delta B/B)^2<1.5\times 10^4`$ \[equivalently, $`\delta B/B<0.01/(ZB_9)^{1/2}`$\] - either the general intergalactic field is weak, or the amplitude of the magnetic turbulence is small, or both.

Compound diffusion, in which particles scatter strongly in localized regions, with scatter free transport between the scattering sites, yields a similar contradiction with the data. If the average distance between scattering centers is an energy independent distance $`l`$ and the particles lose negligible energy while trapped in the scattering sites, then $`D_{eff}=(lcT_{loss})^{1/2}`$ and and one finds $`E^3JE`$ between $`E_a`$ and the highest energy that can be trapped in a scattering site, also contrary to observations of the well-determined spectrum at 10 to 50 EeV.

The most likely candidates for such localized scattering sites are galaxy clusters, with magnetic fields perhaps a few $`\times 10^7`$ Gauss if smoothly distributed in the clusters (see Eilek and Owen 2002 and references therein). The particles which can be diffusively trapped in a (large) cluster have energy less than $`10^{20}ZB_7(\delta B/B)^2`$. Since $`E^3J`$ is flat above $`E_a`$, one can conclude that if UHE cosmic rays arise in normal galaxies with active star formation, as is the case for the model proposed here, most of the particles get from the sources to us without becoming diffusively trapped in large galaxy clusters - which is not surprising, since galaxies with active star formation are not in rich clusters themselves, and rich clusters occupy only a small fraction of the metagalaxy.

## 4 Acceleration and Escape of High Energy Particles from the Magnetar Wind and the Magnetar Wind Nebula

If magnetars were born naked, the theory would be complete, in terms of its energetics. The model’s basic results appear in Figure 1 and Table 1. For fiducial values of the input parameters and weak or moderate gravitational wave losses, the theory does too well - the source strength $`K_0`$ is too large, leading to the inference that only 5-10% of the magnetars are born with voltages sufficient to contribute to the UHECR. The aim of this section is to estimate the losses imposed on the winds as they escape their natal supernovae, to address the acceleration mechanism, and to argue that only half the magnetars have magnetic geometry appropriate to their being UHECR sources. The losses associated with breakout from the supernova (estimated below to be $``$ 30% of the initial rotation energy) affect the energy of the upper cutoff of the spectrum; the geometric factors affect the overall normalization. I also argue that in the equatorial outflow geometry appropriate to the model, adiabatic and radiation losses are negligible.

### 4.1 Escape of the Wind from the Supernova

Magnetars are neutron stars. Neutron stars form in core collapse supernovae (Filippenko 2001 and references therein). The injection spectrum (23) can form only if the relativistic wind survives with negligible energy loss to large radii, and only if the high energy particles accelerated in the wind can escape into interstellar and intergalactic space with negligible adiabatic and radiation losses. In this and the following two sections, I outline the physics of the wind and particle escape which justifies adoption of the bare magnetar spectrum (32) as a reasonable “bottoms-up” model for the origin of UHE cosmic rays.

Particles with energy $`EE_{max}`$, an energy well above existing observations, must come from the outer limits of the relativstic wind, which at time $`t`$ lie between $`r=ct`$ and $`r=c(tt_{spin,i})`$, with $`t_{spin,i}=\mathrm{min}(\tau _{EM},\tau _{GW})`$ \- the outermost 5 light minutes or less. If the energy in this region instead is dissipated through doing work on the plasma and radiation in the supernova and the presupernova environment, the spectrum would be truncated at the highest energies, and possibly suppressed at lower energies if the environment dissipates wind energy at later times. The observed UHECR, with energies $`E<10^{20.5}E_{20.5}`$ eV, require escape of the nonradiative wind to large radii only for $`\mathrm{\Omega }<3\times 10^3(E_{20.5}/Z\eta _1\mu _{33})^{1/2}\mathrm{s}^1`$; from (30), such rotation rates occur 45 minutes and more after the magnetar’s formation.

The essential issue arises in the initial spindown time being short or comparable to the time for the explosive ejection of the supernova envelope666Many of the topics discussed in this section have also been considered by Inoue et al. (2002) using a rapidly rotating pulsar with a $`10^{13}`$ Gauss magnetic field, in the context of the supranova model for gamma ray bursts, assumptions which yield rotational generation of electromagnetic energy on a time scale long compared to the supernova explosion itself.. Magnetars, if born rapidly rotating, necessarily get rid of much of their rotational energy before the extended envelope of an isolated Type II SNe can be fully ejected, which opens the possibility of the energy of the relativistic wind being dissipated in the dense matter around the newly formed magnetar, rather than being made available to relativistic expansion and particle acceleration at large radii.

If core collapse SNe form in binaries, or shed their extended envelopes before the core collapses through large mass loss, the presupernova star may be the compact ($`R_{}10^{5.3}R_{5.3}`$ km) stripped helium core of the pre-supernova massive star, with envelope mass after core collapse $`M4M_4`$ solar masses (Wheeler et al. 2000, Mazzali et al. 2002). Compact progenitors for core collapse SNe occur in Type Ib and Ic supernovae, and it is interesting to note that the poorly known Ib/c supernova rates (Caparello et al. 1999) are comparable to the even more poorly known magnetar birth rate (Gaensler et al. 2001). Then the dynamical time of the envelope surrounding the newly formed neutron star is fairly short, $`t_d2.3R_{5.3}^{3/2}/M_4^{1/2}`$ minutes. The gravitational binding energy of the envelope is small, $`W_g2\times 10^{50}M_4^2/R_{5.3}`$ ergs, a figure which includes the contribution from the neutron star’s gravity. If there were no macroscopic electromagnetic energy input from the newly formed magnetar, the supernova shock, generated either from core bounce or neutrino driven convection and expanding through the envelope at speed $`v_s`$, would begin ejection of the envelope in $`t_{SN}=R_{}/v_s18R_{5.3}/v_{s30}`$ seconds $`t_d`$, where $`v_{s30}=v_s/30,000`$ km/s. The envelope thus forms a moving cavity containing the magnetar, expanding with energy $`4\times 10^{52}M_4v_{s30}^2`$ ergs (Wheeler et al. 2000).

The newly formed magnetar initially emits electromagnetic energy at a prodigious rate - $`\dot{E}_{EM}2\times 10^{50}\mu _{33}^2\mathrm{\Omega }_4^4`$ ergs/s, declining to $`10^{47}`$ ergs/s after a few days. Thus, within the time used by the supernova shock to set the envelope in motion, the spinning magnet in the middle fills the cavity with “lightweight” relativistic energy, mostly in the form of electromagnetic fields, since at this early phase radiation losses bleed off the energy of any accelerated particles777The photon energy created contributes to the lightweight energy inside the cavity, until the Rayleigh-Taylor instabilities which shred the envelope allow all the photons, electromagnetic fields and relativistic particles to escape into the surrounding circumstellar/interstellar (ISM/CSM) medium.. Filling the cavity with buoyant “lightweight” energy leads to prompt “shredding” of the envelope through Rayleigh-Taylor instabilities.

Such shredding happens quickly. Consider the “radiation pressure” exerted on the envelope by the electromagnetic fields spun off from the newly formed magnetar. This pressure creates an effective gravity at the expanding envelope’s inner surface $`gU_{EM}A_{envelope}/M_{envelope}`$. Since this acceleration of the envelope’s inner boundary is Rayleigh-Taylor unstable, buoyant “bubbles” of light weight energy rise through the shell at the generalized Alfven speed $`v_A=(U_{EM}/\rho _{envelope})^{1/2}[\mathrm{\Delta }E_{EM}(t)/M]^{1/2}`$ (adapting the results of Davies and Taylor 1950 and Layzer 1954) and burst out, shredding the envelope in the time $`t_{shred}`$. The total electromagnetic energy contained after spin down to angular velocity $`\mathrm{\Omega }`$ is

$$\mathrm{\Delta }E_{EM}=\frac{1}{2}I\mathrm{\Omega }_g^2\mathrm{ln}\left(\frac{1+\frac{\mathrm{\Omega }_i^2}{\mathrm{\Omega }_g^2}}{1+\frac{\mathrm{\Omega }^2}{\mathrm{\Omega }_g^2}}\right)\frac{1}{2}I\mathrm{\Omega }_g^2\mathrm{ln}\mathrm{\Lambda }(\mathrm{\Omega })1\times 10^{52}\frac{\mu _{33}^2}{I_{45}ϵ_2^2}\frac{\mathrm{ln}\mathrm{\Lambda }}{2.4}\mathrm{ergs}.$$
(34)
For $`\mathrm{\Omega }\mathrm{\Omega }_g(ϵ_2=1),E10^{20.5}`$ eV, $`\mathrm{ln}\mathrm{\Lambda }1.7`$, while at late times ($`\mathrm{\Omega }\mathrm{\Omega }_g,E10^{20.5}`$ eV) $`\mathrm{ln}\mathrm{\Lambda }\mathrm{ln}\mathrm{\Lambda }(\mathrm{\Omega }=0)2.4`$. Expression (30) gives an adequate approximation for the time to spin down to angular velocity $`\mathrm{\Omega }`$; the general expression is in (26). Interstellar medium constraints (see §5) suggest gravitational radiation losses are substantial, thus requiring $`\mathrm{\Omega }_g<\mathrm{\Omega }_i`$ and the use of expression (34) in evaluating the electromagnetic input to the supernova.

Rayleigh-Taylor instability, an adiabatic process, typically focuses the heavy fluid (the material in the envelope) into elongated “spikes” presenting a cross sectional size $`R/\kappa `$ with respect to the outflow, covering a fraction $`1/\kappa ^2`$ of the initially unstable surface. The dominant scale of the magnetic Rayleigh-Taylor instability (Kruskal and Schwarschild 1954, Hester et al. 1996) as it becomes nonlinear is the largest unstable mode that will fit in the thick envelope, $`\lambda _iR/2`$. Therefore, the number of envelope fragments is $`N_{fragment}A_{envelope}/\pi \lambda _i^216`$. Electromagnetic compression rapidly collapses the fragments into a small percentage of the volume \- simulations of magnetic Rayleigh-Taylor instabilities (e.g., Jun et al. 1995) suggest the fragments would have typical dimensions $`0.5\lambda _i(\kappa 4)`$ .

Envelope shredding happens when the bubbles break through the expanding stellar envelope. The bubbles rise faster in the expanding envelope as the electromagnetic energy accumulates. Most of the envelope’s material is in a shell of thickness $`\mathrm{\Delta }R`$ not greater than the presupernova star’s radius $`R_{}`$. Therefore $`t_{shred}\mathrm{\Delta }R/v_A(t_{shred})R_{}/v_A(t_{shred})`$. Assume $`\tau _gt_{shred}\tau _{GW}`$, an assumption found to be self-consistent a posteriori. Then $`\mathrm{ln}\mathrm{\Lambda }1`$ in (34), $`\mathrm{\Delta }E_{EM}(t_{shred})(1/2)I\mathrm{\Omega }_g^2\mathrm{ln}\mathrm{\Lambda }=4.5\times 10^{51}(\mu _{33}^2/I_{45}ϵ_2^2)\mathrm{ln}\mathrm{\Lambda }`$, $`v_A7500(\mu _{33}/ϵ_2)\sqrt{\mathrm{ln}\mathrm{\Lambda }/I_{45}M_4}`$ km/sec and

$$t_{shred}27R_{5.3}\frac{ϵ_2}{\mu _{33}}\left(\frac{M_4I_{45}}{\mathrm{ln}\mathrm{\Lambda }}\right)^{1/2}\mathrm{seconds};$$
(35)
this estimate assumes the envelope thickness is comparble to its radius.

The energy $`\mathrm{\Delta }E_{shred}`$ expended in shredding the shell into fragments which cover a fraction of the area around the newly forming relativistic bubble is not large. The mechanical work that has to be done to shred the shell is

$`\mathrm{\Delta }E_{shred}`$ $``$ $`M_{envelope}gv_A(t_{shred})t_{shred}\mathrm{\Delta }E_{EM}{\displaystyle \frac{v_A}{v_s}}`$ (36)
$`=`$ $`\mathrm{\Delta }E_{EM}\left(\mathrm{\Delta }E_{EM}/Mv_s^2\right)^{1/2}|_{t=t_{shred}}(1W_{shred})\mathrm{\Delta }E_{EM}|_{t=\mathrm{}}.`$
I used $`g=U_{EM}A/M\mathrm{\Delta }E_{EM}/MR(t)\mathrm{\Delta }E_{EM}/Mv_st`$ to obtain this expression.

Shredding of the envelope leaves most of the electromagnetic energy free to escape into the surrounding interstellar medium. From expression (36) the fraction

$$W_{shred}=1\frac{\mathrm{\Delta }E_{shred}}{\mathrm{\Delta }E_{EM}|_{t=\mathrm{}}}=10.1\frac{2.4}{\mathrm{ln}\mathrm{\Lambda }(\mathrm{\Omega }=0)}\left(\frac{\mu _{33}^2}{I_{45}M_4ϵ_2^2v_{30}^2}\right)^{1/2}$$
(37)
of the electromagnetic energy can escape to form a relativstic wind and Magnetar Wind Nebula, assuming mixing of the nonrelativistic material from the shredded supernova envelope does not create too much mass loading.

The relativistic fields and particles emitted by the magnetar behave as a relativstic “fluid” only so long as the non-relativistic material in the envelope remains unmixed with the relativistic gas. The Rayleigh-Taylor “fingers” are themselves subject to disruption from subsequent Kelvin-Helmholtz instabiltities, whose final consequence is to mix the non-relativistic material of the presupernova star’s envelope with the relativistic outflow immediately in contact with the nonrelativistic material, creating local “mass loading”. The Kelvin-Helmholtz mixing speed when one of the fluids is relativistic is (Arons and Lea 1980) $`v_{mix}0.1\gamma _{KH}/k0.1(U_{EM}/\rho _{finger})^{1/2}=.1v_A/\kappa ^{3/2}`$, where $`v_A=(\mathrm{\Delta }E_{EM}/M)^{1/2}`$ as in the discussion of interchange motions. For the environment considered here, the mixing times are on the order of an hour, during which time the envelope fragments move to several hundred times the size of the presupernova star, while the relativistic gas not in contact with the fragments escapes at the speed of light. The nonrelativistic matter mixes into only a small fraction, $`\kappa ^2`$ 10-20%, of the relativistic outflow streamlines, thus leaving almost all of the relativistic energy that survives envelope shredding and fed in by the magnetar at later times to be free of mass loading and able to accelerate high energy particles at radii large compared to the size of the ejecta. The fraction of the energy surviving mass loading thus is given by the covering factor of the shredded ejecta, $`W_{load}1N_{fragment}(A_{fragment}/A_{envelope})=1(N_{fragment}/4\kappa ^2)0.75.`$ Thus, $`W_{blowout}=W_{shred}W_{load}0.60.7`$; 60-70% of the lightweight energy survives to break out of the confining supernova to blow a cavity in the surrounding interstellar/circumstellar medium (ISM/CSM), before the envelope ejecta get a chance to drive the standard supervova blast wave. Instead, the relativistic bubble (large scale, low frequency EM fields and particles) drives a relativistic blast wave into the ISM/CSM.

The losses involved in the relativstic wind pushing out of the supernova come from the initial rotational energy loss put into the wind. Therefore, the highest voltage parts of the wind would be lost, which would reduce the maximum particle energy that can be achieved by a factor $`W_{blowout}^{1/2}0.8`$.

### 4.2 Relativistic Blast Wave and Magnetar Wind Nebula in the Interstellar Medium

The particles accelerated in the wind still must escape the Magnetar Wind Nebula without appreciable loss, even after the wind and the MWN escape the supernova, if the bare magnetar spectrum (23) is to represent the rate ultra-high energy ions are injected into the metagalaxy. For times greater than $`t_{shred}`$, the relativistic energy escaping from the supernova drives a shock wave directly into the surrounding interstellar or circumstellar medium, superseding the usual shock driven by the supernova ejecta - the ejecta are left behind. For simplicity, I consider only the case of an explosion into a uniform external medium of mass density $`\rho _1`$; the results for an explosion into a cavity excavated by mass loss from the supernova’s progenitor, with mass density $`r^a`$, are similar to the uniform case.

Since $`\mathrm{\Delta }E_{EM}/(ct)^3\rho _1c^2`$ for many months after the explosion, for all reasonable densities of the interstellar medium, the relativistic energy drives a relativistic blast wave, whose Lorentz factor and radius are

$`\mathrm{\Gamma }_b`$ $`=`$ $`\left[{\displaystyle \frac{17}{8\pi }}{\displaystyle \frac{\mathrm{\Delta }E_{EM}}{\rho _1c^2(ct)^3}}\right]^{1/2}`$ (38)
$`R_b`$ $`=`$ $`c\left[t{\displaystyle _0^t}𝑑t^{}{\displaystyle \frac{1}{2\mathrm{\Gamma }_b^2(t^{})}}\right]=ct\left[1{\displaystyle \frac{1}{8\mathrm{\Gamma }_b^2(t)}}\right].`$ (39)
The dependence of $`\mathrm{\Gamma }_b`$ on $`\mathrm{\Delta }E_{EM}`$ and $`\rho _1`$ comes from simple consideration of momentum conservation; the factor $`17/8\pi `$ comes from the similarity solution (Eltgroth 1972, Blandford and McKee 1976). The shocked interstellar medium occupies a thin shell, whose thickness is $`\mathrm{\Delta }R_{ism}=R_b/6\mathrm{\Gamma }_b(t)`$. Within this thin shell, the interstellar magnetic field is greatly amplified, $`B_{2ism}=3B_1\mathrm{\Gamma }(t).`$ Flux conservation since the start of the blast yields the coeffcient $`3`$ in this expression; the instantaneous jump conditions replace $`3`$ with $`2`$, which applies to the $`B`$ field just behind the blast wave. The relativstic expansion continues until the time $`T_{rel}`$ such that $`\mathrm{\Gamma }_b(T_{rel})1`$, which yields

$$T_{rel}=\left(\frac{17}{8\pi }\frac{\mathrm{\Delta }E_{EM}}{\rho _1c^5}\right)^{1/3}=1.75\left(\frac{\mathrm{\Delta }E_{52}}{n_1}\right)^{1/3}\mathrm{years};$$
(40)
$`n_1`$ is the baryon number density upstream of the blast wave in units of 1 cm-3. For $`t>T_{rel}`$, the blast enters the nonrelativistically expanding, energy conserving Sedov phase.

The relativistic fields and particles from the magnetar fill the volume out to a contact surface located at $`R_c=R_b\mathrm{\Delta }R_{ism}`$, which is close to $`R_b`$ when $`t<T_{rel}`$; it suffices to take $`R_c`$ as being equal to $`R_b`$. This volume forms the Magnetar Wind Nebula (MWN) in its relativstic expansion phase. The freely expanding magnetar’s wind terminates in a relativistic reverse shock wave located at radius $`R_s`$ where the momentum flux in the wind matches the pressure in the magnetar wind nebula (almost all of which is in the initial energy pulse $`W_{blowout}\mathrm{\Delta }E_{EM}`$). I consider only times greater than $`t_{shred}>\tau _{GW}`$, since at earlier times the relativistic blast wave has not yet formed.

The pressure balance condition is

$$\frac{\dot{E}_{EM}(tR_s/c)}{4\pi R_s^2c}=U_{MWN}\frac{3}{4\pi }\frac{\frac{1}{2}I\mathrm{\Omega }_g^2W_{blowout}\mathrm{ln}\mathrm{\Lambda }(\mathrm{\Omega })}{R_b^3}.$$
(41)
One readily finds, upon using (30) and assuming $`R_s`$ is not almost equal to $`ct`$, that the reverse shock forms where

$$\left(\frac{R_s}{R_b}\right)^2=\frac{2}{3W_{blowout}\mathrm{ln}\mathrm{\Lambda }}\left(\frac{2\tau _g}{t}\right)^{4p1}$$
(42)
with $`p=1/4`$ for $`\tau _{GW}t<2\tau _g`$, and $`p=1/2`$ for $`t>2\tau _g`$. At early times, when gravity waves control the spindown, the reverse shock expands with the MWN ($`4p1=0,R_s=(2\sqrt{W_{blowout}\mathrm{ln}\mathrm{\Lambda }}/3)ct0.63ct`$), but for times greater than $`2\tau _g,p=1/2`$ and the blast wave leaves the reverse shock behind888If gravity wave losses are negligible, pure electromagnetic spindown yields $`(R_s/R_b)^2=(2/3W_{blowout})(2\tau _{EM}/t),R_s=\sqrt{2/3W_{blowout}}c(2\tau _{EM}t)^{1/2}=4\times 10^{13}(I_{45}/\mu _{33}^2\mathrm{\Omega }_4^2)^{1/2}t_{hours}^{1/2}\mathrm{cm}`$ when $`t2\tau _{EM}t_{shred}`$. Since $`t_{shred}2\tau _{EM}`$, at very early times $`R_s(2/3W_{blowout})^{1/2}ct=3\times 10^{13}(0.7/W_{blowout})^{1/2}(t/10\mathrm{min})\mathrm{cm},t<2\tau _{EM}`$., with

$$R_s=\left(\frac{2}{3W_{blowout}\mathrm{ln}\mathrm{\Lambda }}\right)^{1/2}c(2t\tau _g)^{1/2}=6.5\times 10^{13}\left(\frac{0.7}{W_{blowout}}\frac{2.4}{\mathrm{ln}\mathrm{\Lambda }}\frac{I_{45}^3ϵ_2^2}{\mu _{33}^4}\right)^{1/2}t_{hours}^{1/2}\mathrm{cm}.$$
(43)
Because the wind, at $`r<R_s`$, and the MWN, between the $`R_s`$ and $`R_c`$, both expand relativistically, the electromagnetic fields everywhere ($`0<r<R_c`$) have amplitude proportional to $`1/r`$,

$$B(r,t)=\sqrt{\frac{\sigma }{1+\sigma }}\frac{\mu \mathrm{\Omega }^2(tr/c)}{c^2r}f(tr/c)=\sqrt{\frac{\sigma }{1+\sigma }}\frac{\mathrm{\Phi }(tr/c)}{r}f(tr/c),$$
(44)
where $`\sigma `$ is the ratio of the electromagnetic energy flux to the plasma kinetic energy flux in the wind interior to $`R_s`$. $`f(r,t)`$ describes short wavelength \[$`\lambda 2\pi c/\mathrm{\Omega }(tr/c)`$\] structure in the wind, as is the case if the wind has wave structure similar to vacuum strong waves or, if an MHD wind, is “striped” (Michel 1971, Coroniti 1990, Michel 1994, Bogovalov 1999), with the B field direction reversing with wavelength $`\lambda `$ while the overall amplitude declines inversely with $`r`$.

Pulsar Wind Nebulae (PWN) observed in our own galaxy expand nonrelativistically, due to the fact that the initially injected relativistic energy $`\mathrm{\Delta }E_{EM}`$ is small compared to the kinetic energy of the conventional supernova. The resulting “bags” of electromagnetic energy and particles injected by the pulsar are confined by the inertia of the supernova ejecta until times large compared to their ages when we see them, as is the case in the Crab Nebula (Hester et al. 1996), or by the surrounding CSM or ISM, as might be the case for the nebula around PSR 1509-58 (e.g., Gaensler et al. 2002). In contrast, the hypothesized MWN described here expands relativistically at all radii $`R_srR_c=R_b\mathrm{\Delta }R_{ism}`$, for $`t<T_{rel}`$. Therefore, the post-shock velocity at $`R_s`$ must be relativistic, and $`\sigma 1`$ \- the magnetar’s wind and the MWN must be everywhere electromagnetically dominated and “clean” - throughout the relativstic expansion phase, the field structure will not be disturbed by interaction with clumps of nonrelativstic matter, since all such material lies in the shreds of supernova material through which the wind escapes, and in the shocked shell of ISM/CSM right behind the blast wave.

If magnetars behave like pulsars and have wind outflows with a dense $`e^\pm `$ plasma $`(\dot{N}_\pm \dot{N}_{GJ}`$) in addition to the Goldreich-Julian fluxes associated with the electric currents (the high energy ions form one of these electric current flows), the wind and MWN would have a structure perhaps described by MHD. $`\dot{N}_\pm `$ for a magnetar is not well understood, either theoretically or observationally. Observed magnetars are slow rotators with voltages below the level at which pulsars stimulate observable PWN, so prominent nonrelativistic MWN around observed galactic magnetars would not be expected. Theoretically, pair production in the ultrastrong magnetic fields of magnetars has been controversial, with the inhibitions of pair creation by photon splitting being uncertain (see Baring and Harding 2001). However, since $`\sigma `$ must be large, the MWN would not be an efficient radiator, if all the acceleration of pairs is associated with the reverse shock at $`R_s`$, until $`t>T_{rel}`$ 1 year.

### 4.3 Acceleration of UHE Ions in the Magnetar Wind

The metagalactic magnetar model for UHE cosmic rays has been developed using the phenomenology and theory of relativistic winds in galactic Pulsar Wind Nebulae (PWN) as the basic physics input. So far I have not described a specific acceleration mechanism, through which the ions in the equatorial Goldreich-Julian current actually tap the voltage available in the wind; the same lacuna exists in the theory of ions in PWN. A detailed investigation of this subject is outside the scope of this paper. However, a few remarks are relevant, addressed to the most likely mechanism in this context, the “surf-riding” of charged particles in the relativistically outflowing electromagnetic fields.

As is clear from expression (7), the ions would gain energy comparable to the total available voltage if they can cross magnetic field lines for a distance comparable to one decade in radius. Such field line crossing can occur in either radius or in latitude, or both. To see how the simplest form of surf-riding works, consider particle motion in the asymptotic wind of the aligned rotator (Buckley 1977, Contopoulos and Kazanas 2002).

Substantial theoretical evidence has accumulated that when the electromagnetic fields dominate the outflow’s energy density, the structure of the dipolar magnetosphere and its outflow can be approximated as force-free \[$`\rho _q𝐄+(1/c)𝐉\times 𝐁=\mathrm{𝟎}`$, with $`\rho _q`$ the charge density and $`𝐉`$ the current density\], and that under the magnetohydrodynamic conditions appropriate to PWN ($`n\rho _q/q`$, where n is the total quasi-neutral plasma density), the structure of the aligned rotator’s wind at radii larger than a few $`R_L`$ closely approximates that of the split monopole (Michel 1974 , Contopoulos, Kazanas and Fendt 1999), except perhaps in narrow cones around the rotation axis. These fields have the form (Michel 1973, Bogovalov 1997)

$`B_r`$ $`=`$ $`\pm {\displaystyle \frac{M}{r^2}},B_\varphi ={\displaystyle \frac{M\mathrm{sin}\theta }{rR_L}}`$
$`E_\theta `$ $`=`$ $`B_\varphi `$ (45)
$`B_\theta `$ $`=`$ $`E_r=E_\varphi =0.`$
The magnitude of the monopole moment is $`M=k\mu /R_L`$, where $`\mu `$ is the dipole moment and $`k`$ is a constant on the order of unity. For the self-consistent numerical solution obtained by Contopoulos et al. (1999), $`k=1.36`$. The coordinate system is spherical, with $`\theta =0`$ being the colatitude of the angular momentum vector, and with the azimuth $`\varphi `$ being measured with respect to an arbitrary $`x`$ axis. The plus and minus signs apply to the opposite (northern/southern or southern/northern, with respect to the rotation axis) hemispheres. The fields reverse direction abruptly across the rotational equator ($`\theta =\pi /2`$), which implies the existence of an equatorial current sheet where the return current flows. In the body of the wind, the volume current is radial.

The electromagnetic energy propagates with the field line velocity

$`𝒗_E`$ $`=`$ $`c{\displaystyle \frac{𝑬\times B}{B^2}}=c\left(\widehat{𝒓}{\displaystyle \frac{x^2}{1+x^2}}+\widehat{\mathit{\varphi }}{\displaystyle \frac{x}{1+x^2}}\right)`$ (46)
$`v_E`$ $`=`$ $`|𝒗_E|=c{\displaystyle \frac{x}{\sqrt{1+x^2}}}`$
$`\gamma _E`$ $`=`$ $`\left[1\left({\displaystyle \frac{v_E}{c}}\right)^2\right]^{1/2}=(1+x^2)^{1/2},x{\displaystyle \frac{r\mathrm{sin}\theta }{R_L}}.`$ (47)
Charged particles flow out with with velocity $`𝒗=𝒗_E+v_{}𝒃`$, where $`𝒃`$ is the unit vector along the magnetic field. When particles enter the wind with initial energy $`mc^2\gamma _L`$ (corresponding to initial flow 4-velocity $`c\beta _L\gamma _L,\gamma _L=(1\beta _L^2)^{1/2}`$), their parallel speed is, when $`\gamma _L1`$ (Arons, in preparation)

$$\beta _{}=\frac{\gamma _L}{\gamma _E\sqrt{\gamma _L^2+\gamma _E^2}}.$$
(48)
Therefore, the wind’s Lorentz factor is

$$\gamma _w=(1\beta _E^2\beta _{}^2)^{1/2}=(\gamma _L^2+\gamma _E^2)^{1/2}=\left[1+\gamma _L^2+\left(\frac{r}{R_L}\right)^2\mathrm{sin}^2\theta \right]^{1/2}$$
(49)
(Contopoulos and Kazanas 2002). When $`r\mathrm{sin}\theta \gamma _LR_L`$, $`\gamma r\mathrm{sin}\theta /R_L`$ (Buckley 1977); the wind is a linear accelerator, with each particle’s energy increasing by one factor of 10 for each decade in radius. At large radii, the particles move with the field lines, with negligible motion parallel to $`𝑩`$, that is, they “surf-ride” on the electromagnetic field.

This result applies to quasi-neutral winds, whose total density is large compared to the Goldreich-Julian density; it does not apply to the charge separated force free split monopole (Arons, in preparation). Pair creation within the magnetosphere supplies densities large compared to $`n_{GJ}`$, with injection energies $`\gamma _L1010^3`$; the precise values of $`n`$ and $`\gamma _L`$ depend on the specifics of the magnetospheric parameters (surface magnetic field strength, voltage and whether curvature or inverse Compton gamma rays radiated by electrons or positrons accelerated within the magnetosphere initiate the pair cascades.)

Such surf-riding acceleration may seem surprising, since $`𝒗E=0`$. The physics of the acceleration appears when one considers the (small) effect of inertia on the velocity transverse to $`B`$, by solving the equation of motion to first order in the rest mass $`m`$. The solution for the velocity perpendicular to $`𝑩`$ at large radius, where the velocity becomes equal to the radial component of the $`𝑬\mathbf{\times }𝑩`$ drift, is

$$𝒗=𝒗_{Er}+\frac{mc}{q}\frac{𝑩\times D(\gamma 𝒗)/Dt}{B^2}𝒗_{Er}+\frac{mc}{q}\frac{𝑩\times D(\gamma _E𝒗_{Er})/Dt}{B^2}=𝒗_{Er}+𝒗_{pol}.$$
(50)
Here $`D/Dt=/t+𝒗/t+v_{Er}(/r)`$. The polarization drift velocity $`𝒗_{pol}`$ is parallel to the electric field:

$$𝒗_{pol}=\frac{mc}{qB^2}𝑩\times \frac{D(\gamma _w𝒗_{Er})}{Dt}=\frac{mc^3}{q}\frac{r}{M}\widehat{𝜽},$$
(51)
where $`\widehat{𝜽}`$ is the unit vector in the $`\theta `$ direction, and I have assumed $`r\mathrm{sin}\theta R_L`$. $`𝒗_{pol}`$ always has a direction such that $`q𝒗_{pol}𝐄`$ is positive - particles gain energy from the electric field in an accelerating flow. Then energy conservation takes the form

$$mc^2\frac{D\gamma }{Dt}=q𝒗_{pol}𝑬_\theta =+\frac{mc^3\mathrm{sin}\theta }{R_L}.$$
(52)
Integrating with $`D/Dt=c(/r)`$ yields expression (49) in the limit $`r\mathrm{sin}\theta \gamma _LR_L`$. This correspondence can be extended to all radii, at the expense of much more complex algebra.

Thus relativistic force free winds can be linear accelerators, with the particles “surf-riding” on the electromagnetic fields. For the aligned rotator, these fields form a (spherical) DC transmission line, with the fluid accelerating following the increase of the energy transmission speed toward $`c`$ as the magnetic field approaches a purely transverse structure with increasing radius - the wind 4-velocity $`\gamma _w`$ is not constant but is asymptotially $`r`$. For $`c/\mathrm{\Omega }=30/\mathrm{\Omega }_4`$ km and $`\gamma _L10^3`$, such linear acceleration sets in at radii not larger than $`(30,000/\mathrm{\Omega }_4)(\gamma _L/10^3)`$ km.

The result (49) applies to force free conditions, which require $`\sigma B^2/4\pi \rho \gamma _wc^21`$, where $`\rho `$ is the rest mass density measured in the frame where the neutron star’s center of mass is at rest. Assume plasma at nonzero latitude $`|\lambda |=|\theta \pi /2|>0`$ crosses the light cylinder with initial energy/particle $`mc^2\gamma _L=E_L`$ and rest mass density $`\rho _L=(m_i+2\kappa _\pm m_\pm )n_{GJ}(R_L)m_{eff}\mathrm{\Omega }B(R_L)/2\pi cqm_{eff}\mathrm{\Omega }\mu /2\pi cqR_L^3`$, where I have assumed the heavy ions have rest mass $`m_i`$ and density $`n_i=n_{GJ}`$ and are mixed with a pair plasma of density $`n_++n_{}=2\kappa _\pm n_{GJ}n_i`$. Then

$`\sigma _L`$ $`=`$ $`\sigma (R_L)={\displaystyle \frac{B_L^2}{4\pi m_{eff}(\mathrm{\Omega }B_L/2\pi cq)\gamma _Lc^2}}={\displaystyle \frac{1}{2}}{\displaystyle \frac{qB_L}{m_{eff}c\gamma _L\mathrm{\Omega }}}={\displaystyle \frac{1}{2}}{\displaystyle \frac{q\mathrm{\Phi }_{mag}}{m_{eff}c^2\gamma _L}},`$ (53)
$`m_{eff}`$ $``$ $`m_i+2\kappa _\pm m_\pm .`$ (54)
In the absence of magnetic disspation and for $`r\mathrm{sin}\theta \gamma _Lc/\mathrm{\Omega }`$, $`\sigma =\sigma _L\gamma _L(R_L/r\mathrm{sin}\theta )`$, with the force free region ending at $`r=R_1=R_L\gamma _L\sigma _L/\mathrm{sin}\theta `$, where $`\sigma `$ drops to unity. Beyond $`R_1`$, the plasma probably coasts, with $`\gamma _w=\gamma _1\gamma _L\sigma _L`$, in regions where $`\mathrm{sin}\theta 1`$. If the ions are absent, $`m_{eff}2m_\pm \kappa _\pm `$.

The polarization drift velocity (51) causes the ions (and other particles) to drift across B (that drift is the origin of the acceleration), moving off their radial orbits toward the rotation pole a distance $`\mathrm{\Delta }l(r)=(m_{rest}c^2/q\mathrm{\Phi }_{mag})(r^2/R_L)`$. Since the acceleration ends at radii no larger than $`r=R_1`$, the maximum fractional energy gain for an ion from polarization drift alone is $`E/q\mathrm{\Phi }_{mag}=m_i/m_{eff}`$. If $`m_{eff}210m_i`$ and $`\gamma _L`$ a few (as is the case when pair creation cascades are strong, as in the young pulsars) in the return current region, 10% - 50% of the magnetsopheric potential on open field lines could go into accelerating the ions due to polarization drift along $`E`$. This fraction is lower if the pressure of the MWN terminates the wind interior to $`R_1`$.

The accelerator closely approximates a purely radial motion - curvature (a.k.a. synchrotron) radiation losses are negligible, so long as the particles surf-ride ($`𝒗𝒗_{Er}`$).

The energy that can be achieved within $`r=R_1`$ thus depends on the plasma composition. For an aligned rotator with magnetic moment parallel to the angular momentum, the ions in the wind would flow in the equatorial current sheet, while at higher latitudes, the plasma would be composed solely of pairs. The theory of pair creation in ultra-magnetized neutron stars is not as well understood as in standard pulsars ($`\mu 10^{30}`$ cgs), since a variety of strong field QED effects can intervene to inhibit pair creation. A traditional estimate of $`\kappa _\pm `$ which ignores these effects (e.g., Ruderman and Sutherland 1975, Arons and Scharlemann 1979) and that takes account of the limitations set by curvature gamma ray emission on the energies of electrons in the primary polar particle beam yields $`\kappa _\pm 10^4(\mu _{33}R_{10}^2)^{1/4}`$. Taking account of photon splitting, and pair production directly into the lowest Landau levels, effects that reduce the contributions of the synchrotron cascade to pair creation (Baring and Harding 2001), suggests $`\kappa _\pm 10100`$. This polar flow would reach the light cylinder with Lorentz factor $`10100`$, yielding $`E_1=m_\pm c^2\gamma _1=m_\pm c^2\gamma _L\sigma _L=3\times 10^{21}\mu _{33}\mathrm{\Omega }_4^2(10/\kappa _\pm `$ for the asymptotic 4 velocity of an ideal MHD, pair dominated wind outside the equatorial current sheet, assuming no magnetic dissipation interior to $`R_1`$.

In the aligned rotator models contructed to date, the dynamics of charges in the current sheet, where the ions flow as the return current, are not specified. Simple models of the current sheet’s structure (Arons, in preparation) suggest it to be a force-free rotational discontinuity separating the hemispheres. This structure causes the ions to be frozen to the fields which partcipate in the general linear acceleration, and therefore also obey the linear acceleration law (49). If the geometry has the dipole moment antiparallel to the rotation axis, then the ions are themselves the charge carriers in the $`|\lambda |>0`$ hemispheres and participate in the generaL acceleration, while the current sheet carries a negatively charged outflowing return current.

In general, magnetized rotators are not expected to have their dipole moments exactly aligned or counter-aligned to the angular velocity. While a solution of the force-free oblique rotator with dipolar field structure inside the light cylinder has not yet appeared, the experience with the aligned rotator suggests the wind emerging from an oblique rotator will have the structure of the oblique split monopole at radii larger than a few times $`R_L`$. An exact solution for the oblique split monopole has been obtained by Bogovalov (1999), whose results for the fields are shown in Figure 2.

The magnetic and electric field strengths in the oblique split monopole are the same as in the aligned split monopole; the only difference is that in the equatorial region (latitudes $`i<\lambda <i`$, where $`i`$ is the angle between the magnetic and rotation axes), the directions of $`E`$ and $`B`$ reverse every half-wavelength, where the wavelength is $`2\pi R_L`$, forming a “striped” wind (Michel 1971, Coroniti 1990). The current sheet is now wrinkled and is frozen into the flow, advecting outwards with the same 4-velocity as the general wind. That 4-velocity is exactly the same as in the aligned split monopole, given by (49). Thus, under the plausible assumption that the oblique dipole’s wind takes the asymptotic form of the oblique split monopole, the maximum ion energy/particle achieved in all magnetic geometries of rapidly rotating magnetars through the action of polarization drift of the ions along the MHD electric field is

$$E_{max}=m_ic^2\gamma _1=m_ic^2\gamma _L\sigma _L=q\mathrm{\Phi }_{mag}\frac{m_i}{m_{eff}},$$
(55)
where $`m_{eff}=m_i+2\kappa _\pm m_\pm `$ and $`\kappa _\pm `$ is the pair multiplicity in the current sheet.

If $`\kappa _\pm 10^3`$ in the current sheet, $`mm_i`$ and this elementary form of surf-riding yields $`\sigma _L=0.5Ze\mathrm{\Phi }_{mag}/Am_pc^2\gamma _L=6\times 10^{12}(Z/A)\mu _{33}\mathrm{\Omega }_4^2(10/\gamma _L)`$. Then $`R_1=R_L\gamma _L\sigma _L=3.5(Z/A)\mu _{33}\mathrm{\Omega }_4^2`$ pc, and $`E_{max}=3.3\times 10^{22}Z\mu _{33}\mathrm{\Omega }_4^2`$ eV.

Linear acceleration due to polarization drift cannot not go to completion in the freely expanding ($`r<R_s`$) magnetar wind since $`R_sR_1`$. The MWN exterior to the large $`\sigma `$ shock ending the wind also expands relativistically but with a decelerating field line velocity. Therefore, if the magnetar wind model is to explain UHECR, there must be a still stronger acceleration gradient in the wind; laminar surf-riding in the decelerating MWN saps particles of their energy999In pulsar wind nebulae, large-scale surf-riding in the wind is an adequate explanation of the infered wind four velocities (Contopoulos and Kazanas 2002, Spitkovsky and Arons 2003).

The polarization drift based linear accelerator is only a lower limit to the rate at which surf-riding in the magnetar’s wind can can accelerate the charges. If the fields become dissipative within the wind, more rapid acceleration can occur. The low magnetization inferred for the equatorial sectors of pulsars’ winds suggests that such dissipation must exist.

Figure 2 shows that for an obique rotator, much of the Poynting flux is tied up in waves. In the ideal MHD outflow, these waves are frozen in current sheets. These structures become thinner as the wind expands, and are subject to inductive dissipation (Coroniti 1990, Michel 1994, Lyubarsky and Kirk 2001, Lyutikov 2002) at and beyond the disipation radius of the crinkled current sheet $`R_{diss}2\kappa _\pm R_LR_s`$. They are also subject to instability with respect to emission of kinetic Alfven waves (Bellan 1999, 2001, and Arons, in preparation) which propagate with respect to the plasma in the wind frame. Current driven instability with respect to wave emission is a faster process than the inductive dissipation mechanisms discussed to date.

Assume wave emission is the dominant process dissipating the wrinkles in the current sheet, and that such emission leads to the formation of large amplitude electromagnetic waves with dimensionless amplitude $`a=q\delta B/m_ic\mathrm{\Omega }1`$. These waves can accelerate ions faster than does polarization drift in $`E_\theta `$. Even when $`r<R_{diss}`$, the Alfven waves can act as powerful ponderomotive accelertors (Chen et al. 2002), while for $`rR_{diss}`$ these waves become superluminal “strong waves” (Melatos 1998), whose ponderomotive force also yields surf-riding acceleration. When the waves are kinetic Alfven modes, the ponderomotive force applies preferentially to the ions.

An estimate of the ponderomotive work done suffices for the present discussion.

This single particle force arises since the transverse motion of a charged particle in the wave electric field gives rise to $`𝒗\times B`$ force in the wave magnetic field which points radially outwards, does not average to zero over a wave oscillation, and has magnitude

$$𝐟_{pond}=\frac{q}{c}\delta 𝐯\times \delta 𝐁=\frac{q}{mc}\frac{\delta 𝒑_{}}{\gamma _{osc}}\times \delta 𝑩=\frac{q^2}{mc\mathrm{\Omega }}\frac{\delta 𝑬}{\gamma _{osc}}\times \delta 𝑩=mc\mathrm{\Omega }\frac{a^2}{\sqrt{1+a^2}}\widehat{𝐫},$$
(56)
where a particle’s transverse oscillation momentum in the strong wave is $`\delta 𝒑_{}q\delta 𝑬/\mathrm{\Omega }`$, the oscillation energy per particle is $`m_ic^2\gamma _{osc}=m_ic^2[1+(\delta p_{}/mc)^2]^{1/2}=m_ic^2\sqrt{1+a^2}`$, $`\widehat{𝒓}`$ is the unit vector in the radial direction and $``$ indicates an average over wave phase. The work done on a particle as it moves a distance $`\mathrm{\Delta }r`$ through the wave is then, when $`a1`$,

$$W_{pond}=f_{pond}\mathrm{\Delta }r=(mc\mathrm{\Omega }a)\mathrm{\Delta }r=qBr\frac{\delta B}{B}\frac{\mathrm{\Delta }r}{r}q\mathrm{\Phi }_{mag}\eta ,\eta \frac{\mathrm{\Delta }r}{r},$$
(57)
since $`\delta E\delta B`$ and $`\delta B/B1`$ in the rotational equator - most of the energy flow is in the wave Poynting flux. The evidence from Pulsar Wind Nebulae is that $`\eta 0.1`$ and therefore $`\mathrm{\Delta }r0.1r`$. The fact that $`\mathrm{\Delta }rr`$ and the particles gain only a relatively small part of the total voltage available most likely is due to the demise of the strong waves coherence, as a result of parametric instabilities in these electromagnetic structures - see, for example, Leboeuf et al.. (1982).

This argument suggests that such surf riding may be the source of the acceleration hypotheiszed in the model for UHECR studied here. A detailed demonstration that this is so is a task for another study. If so, the ion acceleration is prompt, beginning a few seconds and continuing for a few days after the escape of the relativistic wind from the natal supernova - that escape occurs $``$ a minute after the neutron star’s formation - see expression (35).

The standard alternate for the acceleration mechanism of UHECR is some form of shock acceleration, usually invoked as diffusive Fermi acceleration in relativistic shocks behaving in a manner similar to nonrelativistic shocks (e.g., Waxman 1995, Bahcall and Waxman 2003); see Ostrowski and Bednarz (2002) and references therein for a discussion and critique of this mechanism’s plausibility in relativistic shocks.) In the current context, the fact that the electromagnetic fields dominate the energetics of the outflow militates against the relevance of that mechanism. Alternate, “flare” mechanisms, based on direct dissipation of the currents, are even more poorly understood. In this regard, identifying bursts of UHE cosmic rays as coming from the same direction as bursts of quasi-coherent gravitational radiation, with correlated arrival times, would give substantial empirical input into constraining the specific physical processes behind the acceleration of the highest energy cosmic rays.

### 4.4 Escape of UHE ions from the Magnetar Wind Nebula

In accord with the discussion in §4.3, I assume the high energy ions achieve their energy in the relativistic wind through surf-riding acceleration, at a radius $`R_aR_{diss}\kappa _\pm R_L=(10^310^4)R_L=(30,000300,000)/\mathrm{\Omega }_4`$ km. Beyond $`R_a`$, dissipation of the strong waves reduces the structure of the magnetic field to that of the aligned rotator (Blandford 2002), with fields given by (45) but now with $`M\mu /R_LM_{}\mu \mathrm{cos}i/R_L`$. The acceleration injects the ions into the angular sector $`i<\lambda <i`$, with a 4-velocity large compared to the wind’s 4-velocity (49). One can readily show, either from the equation of motion of a UHE ion in the crossed $`E`$ and $`B`$ fields of the wind, or from the motion of an ion in the pure B field in the $`E=0`$ frame of the wind, that the radius of curvature of an accelerated ion’s orbit is

$$\rho _c=2r_L\gamma _w^2=(2\eta \gamma _w^2/\mathrm{cos}i)r,$$
(58)
where $`r_L=m_ic^2\gamma /qB=\eta r`$ is an ion’s formal Larmor radius. When $`\gamma _wr/R_l1`$, as it is at radii comparable to and larger than $`R_a`$, $`\rho _cr`$.

This large radius of curvature of ions’ orbits, applicable wherever the flow is relatistically expanding (both within the wind and in the MWN) means that both synchrotron and adiabatic losses of the ions are negligible. From the power radiated by a charged particle with orbital radius of curvature $`\rho _c,P=(2q^2c/\rho ^2)\gamma ^4`$, $`\gamma _wr/R_L`$, $`\gamma =\eta q\mathrm{\Phi }/m_ic^2`$, one readily finds

$$T_{rad}=\frac{6}{\eta \mathrm{\Omega }}\frac{A^4}{Z^5}\left(\frac{m_pc^2}{e\mathrm{\Phi }}\right)^3\frac{m_pc^2}{e^2}R_L\mathrm{sec}^2i\left(\frac{r}{R_L}\right)^6=\frac{3.6}{\eta _1\mathrm{\Omega }_4^2}\frac{A^4}{Z^5}\frac{\mathrm{sec}^2i}{\mu _{33}^3}r_{10}^6\mathrm{seconds},$$
(59)
where $`\eta _1=\eta /0.1`$ and $`r_{10}=r/10^{10}`$ cm; one also obtains this same radiation time by applying ordinary synchrotron losses in the fluid frame of the wind where $`𝐄=\mathrm{𝟎}`$. Thus, for $`rR_a>10^{10}`$ cm, radiation losses are negligible ($`T_{rad}r/c`$). Also, since $`\rho r`$, the high energy ions are not coupled to the field lines and do not suffer adiabatic expansion losses, so long as they escape the wind and the MWN while the MWN’s expansion is relativistic.

The particles leave the wind and move through the MWN with speed $`v(t)=c(11/2\gamma ^2)c`$. For the observed UHE cosmic rays, with $`E_a<E<2\times 10^{20}\mathrm{eV}E_{max}`$, the particles at radii large compared to $`R_s`$ were accelerated at radius $`R_a(R_a<R_s`$ by assumption) and at the time $`T_a=R_a/c+(E_a/E)t_a`$, where $`t_a=(9/8)(Ze\eta /\mu E_a)2Z\eta _1I_{45}/\mu _{33}E_{a,18.8}`$ days is the time at which the neutron star’s voltage was at the value corresponding to acceleration of ions at the ankle in the spectrum. Here $`E_{a,18.8}E_a/6\times 10^{18}`$ eV. For particles in the UHECR spectrum measured to date, $`T_at_a(E_a/E)`$, so long as $`R_aR_s`$. At time $`t`$, particles of energy $`E`$ have reached $`R_{particle}=c(tT_a)+R_a`$; for all particles of energy $`E>E_a`$, $`\gamma =E/mc^2\mathrm{\Gamma }_b(t)`$.

The outer boundary of the MWN lies at $`R_b(t)`$ given by (39). All particles with energy $`EE_a`$ eV catch up with the blast wave at the time $`T_{catchup}`$ when $`R_{particle}`$ first exceeds $`R_b(t)`$. Equating the radii yields

$$T_{catchup}=\left(\frac{17}{\pi }\frac{\mathrm{\Delta }E_{EM}}{\rho _1c^5}T_a\right)^{1/4}=8.3\left(\frac{\mathrm{\Delta }E_{EM,52}}{n_1}\frac{I_{45}}{\mu _{33}E_{a,18.8}}\frac{E_a}{E}\right)^{1/4}\mathrm{months}.$$
(60)
Here $`\mathrm{\Delta }E_{EM,52}=\mathrm{\Delta }E_{EM}/10^{52}`$ ergs.

Comparison of (60) to (40) shows that all the UHE cosmic rays reach the blast wave before the MWN becomes nonrelativistic. The 4 velocity of the blast wave at $`T_c`$ is $`\beta _b\mathrm{\Gamma }_b(T_c)=4(E_{52}/n_1)^{1/8}(\mu _{33}E_{a,18.8}/Z\eta _1I_{45})^{3/8}(E/E_a)^{3/8}`$. The cosmic rays have no difficulty getting through the plowed up interstellar medium and its compressed magnetic field just behind the blast wave, once they arrive at the MWN’s outer boundary - the cosmic rays’ Larmor radii in the shocked interstellar magnetic field are $`r_{L2}=540(E_{a,18.8}/B_{1,\mu })(4/\mathrm{\Gamma }_b[T_c(E_a)])(E/E_a)^{5/8}`$ parsecs, which is far larger than the thickness of the shocked shell of interstellar matter. Thus the UHE ions avoid substantial adiabatic energy loss in escaping the MWN. Here $`B_{1,\mu }=B_{ISM}/1\mu \mathrm{Gauss}`$.

Other radiation losses in the MWN might occur if the high energy ions encounter a sufficiently dense external radiation field. The supernova’s own optical and ultraviolet luminosity is the most prominent candidate for such radiation. A reasonable upper limit for that luminosity is $`10^{44}`$ ergs/sec, emitted over the first few days, i.e., over the same time as the UHECR acceleration occurs in the wind. At the radii where the acceleration is likely to occur (see §4.3), photopion losses in the stellar radiation field might limit the maximum ion energy that can be achieved to no more than $`10^{22}`$ eV. Losses due to Compton scattering and pair creation are negligible, as is energy loss in the dense matter of the supernova, since the ejecta have only a small filling factor, and move out slowly behind the high energy particles. Losses in the radiation field of the MWN are unlikely to be of significance, since this highly magnetized flow is unlikely to be radiatively dissipative until after the outer blast wave becomes nonrelativistic, which occurs after the UHE ions have left the system. Radiation losses in the interstellar environment surrounding the MWN (pion and pair production, synchrotron radiation, inverse Compton scattering) are all negligible, even if the MWN forms within a cluster of hundreds of O and B stars, and the magnetic field strength is that of molecular clouds (milligauss). Interactions with the circumstellar matter that might surround the presupernova star also can be neglected.

Radiationless escape of the UHE ions from the MWN works only when the ions form the equatorial return current. This geometry applies to rotators with dipole axes making an acute angle with respect to the rotation axis. Thus, only 50% of the magnetars can supply UHE cosmic rays, contributing the factor $`W_{geom}=0.5`$ to the flux of UHECR (expression 32).

If the relativistic wind has a dense pair plasma, as well as the ion return current in the equator and the polar electron current, Compton scattering from the pairs might create a high energy photon flash, which may be an interesting candidate for a prompt Gamma Ray Burst. This is a topic for a separate investigation. For the present purposes, it suffices to observe that the supernova photon luminosity is small compared to the magnetar’s spin down luminosity, therefore Compton drag can have little effect in dissipating the electromagnetically dominated outflow.

## 5 Interstellar Supershells and Supernova Remnants

Formation of interstellar HI shells and supershells provide an interesting constraint pointing toward only a small fraction of the newly born magnetars having voltages large enough to accelerate UHECR, and perhaps to large gravitational wave emission by such objects.

If the magnetars lost energy purely electromagnetically, they would deliver energy on the order of $`(1/2)I\mathrm{\Omega }_i^2W_{blowout}=3.5\times 10^{52}I_{45}\mathrm{\Omega }_4^2(W_{blowout}/0.7)`$ ergs with each event. These explosions would generate supershells, large expanding structures in the interstellar medium with energy in excess of $`10^{52}`$ ergs (Heiles 1979, 1984) and dimensions in excess of 100 pc. From Chevalier’s (1974) expression (26), with a random velocity of clouds in the ISM of $`20`$ km/s, electromagnetic energy deposition with $`\mathrm{\Omega }_41`$ would create shells with radii $`230`$ pc. Further investigation over the last 20 years has led to the conclusion that all such events in our galaxy occur at a rate $`\nu _{supershell}0.8\times 10^5\mathrm{yr}^1`$ (Ehlerova and Palous 1996), about a factor of 12 slower than the current estimate of the magnetar birth rate; other relations between explosion energy and remnant size (e.g. McKee and Cowie 1975) lead to similar conclusions. Shells as large as 200 pc form only a small, but ill determined, fraction of the total supershell population. Therefore, if the theory outlined here is to be a viable explanation of UHECR, either the magnetar rate is overestimated, or the energy per magnetar delivered to the interstellar medium is smaller (or both).

The approximate correspondence shown in Figure 1 between the theory outlined in this paper and the observed UHECR spectrum requires $`K_0=W_{geom}\nu _{m4}(I_{45}/\mu _{33})n_{g2}T_{H14}0.03`$ if the spindown is purely electromagnetic. This result comes from requiring the model to explain the well determined spectrum between $`10^{19}`$ and $`6\times 10^{19}`$ eV. If the particle spectrum measured at energies below the ankle extrapolates into the UHE regime, $`K_0`$ drops to 0.02 for the pure electromagnetic case. If one restricts oneself to magnetars with parameters defined by current observations, the main freedom is in the fact that we don’t know whether all such objects are born with the same initial rotation rate. Magnetars can succeed as a model for UHECR if their voltages are sufficiently large, requiring initial periods $`P_i<4(\eta _1\mu _{33})^{1/2}`$ milliseconds for at least some of the magnetar population.

With $`W_{geom}=0.5`$, and $`n_g/0.02\mathrm{Mpc}^31`$ , the inferred value of $`K_0`$ tells us that the birth rate of such rapidly rotating objects, which by hypothesis are subject only to electromagnetic torques, is $`\nu _m^{(fast)}0.1\nu _m=10^5\nu _{m4}\mathrm{yr}^1`$ (see §3.3). Then most of the supershells in our galaxy would be the result of these “electromagnetic bombs” in the interstellar medium, and most would be of the largest variety. Loeb and Perna (1998) made the related suggestion that Gamma Ray Bursts drive supershells.

If one insists that most supershells do not come from exotic sources (probably true), and that their numbers not be dominated by the largest variety (certainly true), then the energy delivered to the interstellar medium in each event is less than the $`(1/2)I\mathrm{\Omega }_i^2=5\times 10^{52}I_{45}\mathrm{\Omega }_4^2`$ ergs. For pure electromagnetic spindown, initial rotation rates less than $`10^4`$ sec-1 are the only means of achieving that reduction. Creating shells smaller than 100 pc, which would cause them to be missed in Heiles’s (1984) census, requires the energy deposited per event to be less than $`4\times 10^{51}`$ ergs, in turn requiring $`\mathrm{\Omega }_i<3\times 10^3\mathrm{sec}^1.`$ Then the maximum cosmic ray energy would be $`E_{max}<2\times 10^{20}Z\eta _1\mu _{33}(W_{blowout}/0.7)^{1/2}`$ eV, which is too low, for $`Z1`$.

Furthermore, ordinary supernova remnants in galaxies appear to have sizes typically less than $`50`$ pc (see Matonick and Fesen 1997 and references therein). Therefore, from Chevalier’s energy-size relation, the energy per event would have to be $`10^{50.6}`$ ergs - indeed, the standard result, that supernovae create explosions with energy $`10^{51}`$ ergs, comes from precisely these considerations. Requiring the rotational energy per magnetar to be $`10^{51}`$ ergs or less requires initial angular velocities no greater than $`1.4\times 10^3I_{45}^{1/2}\mathrm{sec}^1`$, corresponding to a maximum cosmic ray energy not greater than $`7\times 10^{19}Z\eta _1/I_{45}`$ eV, which is insufficient to explain the highest energy cosmic rays (if $`Z=12`$).

Therefore, the most plausible conclusion is that gravitational radiation must be the recipient of most of the initial rotational energy. The ultra-strong magnetic fields in these stars make this a viable possibility. From (34), the energy delivered would then be $`E=(1/2)I\mathrm{\Omega }_g^2\mathrm{ln}\mathrm{\Lambda }(\mathrm{\Omega }0)W_{blowout,GR}=0.75\times 10^{52}(\mu _{33}^2/I_{45}ϵ_2^2)(\mathrm{ln}\mathrm{\Lambda }/2.4)(W_{blowout}/0.7)`$ ergs. Requiring the interstellar shells created by these electromagnetic explosions to all have radii less than 100 pc demands that the neutron stars have equatorial ellipticities $`ϵ>2\times 10^2`$; requiring the interstellar remnants of the explosions to have sizes less than 50 pc demands $`ϵ>4\times 10^2`$. From expression (19), such ellipticities would be expected for interior magnetic fields $`>(68)\times 10^{16}`$ Gauss, large but within the realm of recent theoretical suggestions. Then more than 90% of the initial rotational energy of the fast rotators radiates away as gravitational waves, which still leaves more than enough energy for the UHECR.

In passing, it is interesting to note that Matonick and Fesen (1997) found a small number of unusually large (diameters exceeding 100 pc) supernova remnants in their sample, with velocities too large to be readily explained as the consequence of several ordinary but coeval supernova explosions. Such a species of supernova remnant clearly is a candidate for identification with the electromagnetically driven explosions discussed here. However, given the selection effects, it is unclear whether the number of these unusual events is consistent with electromagnetically driven explosions from magnetars (or any other exotic explanation) as their origin.

## 6 Gravitational Radiation

Since this model requires significant influence on the magnetars’ spindown by gravitational radiation losses, with release of most of the initial rotational energy going into this so far undetected carrier, the model makes reasonably definite predictions for the detection of almost coherent signals by LIGO, VIRGO and other gravity wave detectors optimized for khz signals. The dimensionless strain of the narrow band oscillatory signal, measured at distance $`D`$ from the source and observed for a time longer than the coherence time of the oscillator $`\tau _{coh}`$, so that $`n\mathrm{\Omega }\tau _{coh}`$ oscillations can be coherently counted, is (Thorne 1997, Brady et al. 1998)

$$h_n=2.9\frac{GIϵ\mathrm{\Omega }^2}{Dc^4}\sqrt{\mathrm{\Omega }\tau _{coh}};$$
(61)
for simplicity, I have averaged over all orientations of the rotation axis with respect to the observer, and I have included only the wave modes at the star’s rotation frequency, which exist for magnetically distorted stars (Bonazzola and Gourgoulhon 1996). The coherence time is $`\tau _{coh}=\mathrm{\Omega }/|\dot{\mathrm{\Omega }}|`$.

Since gravitational radiation controls the early spindown ($`\tau _{EM}\tau _{GW}`$), $`\tau _{coh}=4\tau _{GW}[1+(t/\tau _{GW})]`$ and $`\mathrm{\Omega }=\mathrm{\Omega }_i/(1+t/\tau _{GW})^{1/4}`$ for $`t<2\tau _g=1.8I_{45}^3ϵ_2^2/\mu _{33}^4`$ hours. Then

$`h_n^{(GR)}`$ $``$ $`4\times 10^{24}{\displaystyle \frac{I_{45}ϵ_2\mathrm{\Omega }_4^2}{D_{20}[1+(t/\tau _{GW})]^{1/2}}}\sqrt{\mathrm{\Omega }\tau _{coh}}={\displaystyle \frac{(\mathrm{\Omega }_4I_{45})^{1/2}}{D_{20}\left(1+\frac{t}{\tau _{GW}}\right)^{1/8}}}`$ (62)
$`=`$ $`2.5\times 10^{21}{\displaystyle \frac{I_{45}^{3/8}}{D_{20}ϵ_2^{1/4}t_{hours}^{1/8}}};`$
here $`D_{20}=D/20`$ Mpc. For $`t>2\tau _g`$ (rotation periods exceeding 4 msec), the electromagnetic torques take over and

$$h_n^{(EM)}=6.5\times 10^{25}\frac{I_{45}^2}{D_{20}}\frac{ϵ_2}{\mu _{33}^2t_{hours}}\sqrt{\mathrm{\Omega }\tau _{coh}}=2.9\times 10^{21}\frac{I_{45}^{9/4}ϵ_2}{D_{20}t_{hours}^{3/4}}.$$
(63)
These strain levels are potentially observable, although the problems of finding such quasi-coherent signals in blind searches, with no advance notice of the period, are quite formidable. Nevertheless, since the site of the cosmic ray acceleration is at small enough radii ($`r<R_s`$ 1 AU), the UHE cosmic rays from an individual outburst should arrive before the gravity waves fade from view - conceivably, a UHE cosmic ray event, if observed with sufficient collecting area to identify a burst of particles, might serve as a marker for searches of gravity wave data for quasi-coherent signals with periods in the one to few millisecond regime.

## 7 Individual Events

Charged particles of energy in excess of $`10^{20}`$ eV must come from nearby. For such ultra-high energy particles, the GZK losses limit the observable volume to size

$`D_E`$ $`=`$ $`\mathrm{min}\{D_{GZK}\left[\mathrm{max}(\mathrm{ln}{\displaystyle \frac{E_{max}}{E}},{\displaystyle \frac{E_G^2}{E^2}})\right],D_{H\pm }\},`$ (64)
$`D_{H\pm }`$ $`=`$ $`\mathrm{min}[cT_H,l_\pm \left({\displaystyle \frac{E_\pm }{E}}\right)^{0.4}].`$
$`D_E50`$ Mpc, for maximum energies appropriate to the model discussed here. See expression (4) and the immediately preceeding discussion for the values of the characteristic lengths and energies. The number of source events per unit time which we can detect, in principle, in the ultra-high energy particle spectrum then is $`\dot{N}_{UHE}=(4\pi /3)D_E^3n_g\nu _m^{(fast)}`$. This rate is $``$ eight thousand per year at $`E_a=6\times 10^{18}`$ eV ($`D_E2.1`$ Gpc), drops to 800/year at $`4\times 10^{19}`$ eV ($`D_E1`$ Gpc), but then declines rapidly, at a rate $`E^6`$ to a few times per decade above $`3\times 10^{20}`$ eV, assuming the fiducial event rate per unit volume $`n_g\nu _m^{(fast)}=2\times 10^7`$ Mpc-3 yr-1. The particle fluence per event is

$$(>E)=\frac{9}{32\pi D^2}\frac{c^2I}{Ze\mu }f(E)4\times 10^2\frac{I_{45}}{Z\mu _{33}D_{100}^2}\frac{f}{2.5}\mathrm{km}^2;$$
(65)
here $`D_{100}=D/100`$ Mpc. $`f23`$ is a logarithmic function of $`E`$ and $`D`$.

The model proposed here requires the particles to be accelerated and emitted in an equatorial slice with opening angle $`2i`$, with $`i`$ the angle between the rotation axis and the magnetar s magnetic moment, into the surrounding Universe - the particles are physically beamed. The fraction of the sky into which the particles are beamed is $`b=\mathrm{sin}i`$. If $`i`$ is uniformly distrubuted between $`0`$ and $`\pi /2`$, the average beaming fraction is $`b=2/\pi =0.64`$. The individual fluence/event would be higher by a factor $`1/b`$. Of course, the rate at which events are observed is lower by a factor $`b`$, which leads to the average spectrum shown in Figure 1 being independent of $`b`$.

Occasionally (once in $`1/\nu _m^{(fast)}10^5`$ years), such events must occur in our own galaxy, which would create a primary particle fluence at the earth of order $`10^7`$ km-2, for an average distance of 5 kpc. Such events might leave tracers in the geochemical record, a topic beyond the scope of this study.

## 8 Discussion: Related Models

Blasi et al. (2000) presented the model most closely related to the theory outlined here. They were the first to identify the relativistic wind of a magnetized rotator as a potential site of UHE cosmic ray acceleration. They assumed the accelerators to be initially rapidly rotating pulsars in our own galaxy ($`\mathrm{\Omega }_i3000\mathrm{s}^1`$), and made the same assumption of the charge lost per second, $`I=q\dot{N}_i=\mu \mathrm{\Omega }^2/c=c\mathrm{\Phi }`$. Since they assumed normal (if rapidly rotating) pulsars, they produced particles with energy $`100`$ EeV only by assuming the high energy ions to be fully stripped iron. They made the same construction of the $`E^1`$ particle spectrum as in §3.1 above, for the pure electromagnetic spindown case. They did not include possible gravitational wave losses. Since the pulsar spin down time for their assumed parameters is long compared to the expansion time of the supernova envelope, the pulsar’s rotational energy loss does not substantially affect the dynamics of the supernova, and most of the UHE particles form at times late enough to avoid substantial radiation and stopping losses in the supernova envelope. They gave no consideration to the disruption of the envelope due to the pulsars’ magnetic pressure on the expelled envelopes. They did not use supershells and the sizes of supernova remnants, to constrain the event rate. The large charge/particle assumed allowed them to argue that scattering in the galactic magnetic field might allow this hypothesized source to be in reasonable accord with the observed UHECR isotropy.

However, UHECR observations support light nuclei, rather than Fe, as the particles responsible for the UHE showers (e.g., Bird et al. 1994, Watson 2002). X-ray observations of pulsars suggest the atmospheres from which the ions would be extracted are primarily H or He, as discussed in §2, not Fe. The energies of protons or $`\alpha `$-particles would be too low to explain the super-GZK particles in their model. Also, the Larmor radii would be too large, to allow isotropization in the galactic magnetic field. Finally, within the scheme actually outlined by Blasi et al., with each pulsar’s magnetic field contained within a closed shell of supernova ejecta, the supposed UHE cosmic rays suffer catastropic adiabatic losses. Because of the relatively small rotational energy per pulsar assumed in their model, and the slow release of that energy ($`\tau _{EM}1`$ year for their parameters), disruption of the shell and escape of the wind and its high energy particles takes too long for the UHE particles to survive, with ordinary (if rapidly rotating) pulsars as the drivers.

The metagalactic magnetar model outlined in this paper avoids all of these difficulties, as well as adding consideration of the specific mechanism of acceleration in the wind and of gravitational wave losses on the rotation history of the pulsar, which can have observable effects on the UHE particle spectrum - see Figure 1.

The theory of supernova envelope disruption and relativistic expansion of a Magnetar Wind outlined in §§4.1 \- 4.2 have close kinship, and owe much to, the Gamma Ray Burst models of Wheeler et al. (2000) and of Inoue et al. (2002). The principal difference from the Wheeler et al. scenario is in the disruption of the supernova envelope by Rayleigh-Taylor instabilities, rather than assuming a jet punches through the envelope. Inoue et al. asummed a rather lower energy rotator, similar to Blasi et al.’s fiducial objects, so that Rayleigh-Taylor disruption of the ejecta happens at rather later times. They also assumed the Pulsar Wind Nebula in their model would be radiatively efficient so as to construct a Gamma Ray Burst model based on nonthermal emission from the Pulsar Wind Nebula, even though they also assumed the outflow to be in the form of an electromagnetically driven relativistic blast wave, an outlow which may be an inefficient radiator, at least until it decelerates to nonrelativistic velocities. In the present model, while there may be an initial flash of Compton scattering as the electromagnetic fields disrupt the supernova shell, the subsequent flow is highly magnetized. Such flows have weak shock dissipation (Kennel and Coroniti 1984), and may have weak magnetic dissipation (Lyubarsky and Kirk 2001), so the conversion of outflow energy to efficiently radiating accelerated particles is likely to be weak. However, other processes in a strongly magnetized medium associated with intense electric current flow may allow conversion of magnetic energy to emergent photons from the wind with unknown efficiency \- see Blandford (2002) for comments on some of the possibilities.

The considerations of laminar surf-riding acceleration outlined in §4.3 have their roots in the work of Buckley (1977) and of Contopoulos and Kazanas (2002). The relation of these accelerating, force-free wind results to MHD analyses of relativistic winds will be discussed elsewhere.

Magnetars behaving as outlined in this paper probably cannot co-exist with long lived ($`t_{disk}\tau _{GW}`$) fall back disks, as are invoked in collapsar and supranova models of GRBs - such structures would probably suppress or greatly modify the flow and acceleration of the ion return currents in the rotational equator. However, if that disk is itself strongly magnetized, in principle it could play the same role as a UHECR accelerator. The open magnetic flux in the magnetar model, $`\mathrm{\Psi }=R_L\mathrm{\Phi }10^{29}`$ Maxwells, is similar to that expected in many models of Gamma Ray Burst sources’ disks, either due to binary mergers or to fallback in a supernova (Narayan et al. 1992, Usov 1994, Paczynski 1998, Vietri and Stella 1999, Macfadyen and Woosley 1999). Many of these models involve magnetic fields with magnitude $`>10^{16}`$ Gauss rotating at Keplerian angular velocity at $`r6R_{Schwarschild}=179M/10M_{}`$ km. These objects would have $`\mathrm{\Psi }10^{30}`$ Mx and $`\mathrm{\Phi }10^{24.7}`$ Volts, clearly making them available as candidates for UHECR acceleration. Time dependent decay of the disk might play a role similar to a magnetar’s spindown, in creating a broad energy spectrum of high energy particles. Such a model has the advantage that the $`>`$ petagauss magnetic fields in the disks, if both the disks and the strong fields exist, are surely rapidly rotating. The magnetar model has the advantage that neutron stars with petagauss surface fields do exist, but we do not have independent evidence on how many are born rapidly rotating.

## 9 Conclusions

The model’s basic results and predictions are shown in Figure 1. The upper cutoff of the spectrum, at $`E=10^{21.5}Z\eta _1\mu _{33}\mathrm{\Omega }_4^2`$ eV, is likely to be reduced by dissipation of the wind as it first breaks free of the magnetar’s natal supernova. The estimates of §4.1 suggest that escape from the supernova consumes a fraction $`1W_{blowout}0.30.4`$ of the magnetar’s intial energy loss, in turn suggesting that $`80`$% of the highest voltages are available for cosmic ray acceleration. If so, then the highest energies in the spectra shown in Figure 1 should be be taken seriously.

Agreement with the well determined UHE cosmic ray flux below $`5\times 10^{19}`$ eV constrains the combination $`W_{geom}n_{g2}\nu _{m4}/\mu _{33}`$ to be approximately 0.02-0.06, where $`\mu _{33}`$ is the magnetic moment in units of $`10^{33}`$ cgs, $`\nu _{m4}`$ is the magnetar birth rate in units of $`10^4`$ yr-1, $`n_{g2}=n_{galaxy}/0.02\mathrm{Mpc}^3`$ and $`W_{geom}=0.5`$ is the inefficiency factor due to the fact that only 50% of the magnetic geometries are appropriate to the ion current being emitted into the rotational equator, where radiation and adiabatic losses on the UHECR would be weak. The small value of this normalization required to get a fit, appropriate to the “no GR” and “moderate” GR cases considered in Figure 1, implies that the birth rate of magnetars with voltages high enough to create UHECR is 5-10% of the overall magnetar birth rate, a conclusion consistent with not overpopulating our and other galaxies with too many supershells and overly large supernova remnants.

The high energy ions suffer negligible adiabatic and radiation losses in escaping the magnetars’ winds and their surrounding, relativistically expanding nebulae. The existing observations require voltages $`\mathrm{\Phi }>2\times 10^{21}/Z\eta _1`$ Volts, or $`\mu _{33}\mathrm{\Omega }_4^2>0.06/\eta _1Z(P_{initial}<2.6(Z\eta _1\mu _{33})^{1/2}`$ ms), with $`\eta _1=\eta /0.1`$ the fraction of the voltage actually sampled by each ion. I suggested in §4.3 that electromagnetic surf-riding in the relativistically strong electromagnetic waves in the wind, generated by the oblique rotator, is responsible for the ion acceleration, occuring at or within the radius ($`0.1`$ AU) where the wave structures no longer are frozen into the expanding plasma. That acceleration site is consistent with small radiative losses.

This model works only if the ions travel from their source galaxies on more or less straight lines. That requires intergalactic magnetic turbulence to be small amplitude, or the intergalactic magnetic field to be weak, or both - $`\delta B/B<10^2(10^9\mathrm{Gauss}/ZB_{IGM})^{1/2}`$3.4). Such weak scattering introduces negligible time delays into the particles’ transport from source to observer.

If the current data are interpreted as providing evidence for the conventional GZK cutoff, i.e. the AGASA results are disregarded (Bahcall and Waxman 2003) or are renormalized to bring them into accord with the Hi-Res results (Abu-Zayyad et al. 2002b), agreement of this model with existing observations requires rather strong losses of rotational energy due to gravitational radiation. The “strong GR” curve in Figure (1) shows the effect of the rapid spindown depopulating the highest energy end of the injection spectrum. The model can replicate the high energy end of the spectrum observed by AGASA only if gravitational wave losses are small. The metagalactic magnetar model has no trouble with a spectrum that compromises between the highest energy AGASA results and the Hi-Res (without the 1995 Bird et al. event) results - including that data point favors such a compromise.

Even when one assumes a birth rate of fast magnetars $``$ 5-10% of the overall magnetar birth rate, the requirement that the model not overproduce supershells and large supernova remnants in our and other galaxies’ interstellar media suggests gravitational wave losses are an important drain of rapidly rotating magnetars’ initial rotational energy. Therefore, one expects bursts of almost coherent gravitational waves with millisecond periods and strains with magnitude $`10^{21}`$ if observed in an experiment designed to find almost coherent signals, lasting for several hours and overlapping the arrival of the higher energy particles.

The model suggests the UHE cosmic rays come from sources whose distribution should mimic that of luminous baryons. More specifically, the sources should follow luminous matter with stars that are progenitors of core collapse SNe, thus should anticorrelate with large galaxy clusters. Events with energy above $`10^{20}`$ eV come from small distances ($`D<50`$ Mpc), which may allow some imprint of galaxies’ large scale structure to appear on the isotropy of such events. Starburst galaxies are obvious suspects for especially luminous UHECR sources.

In conclusion, I have shown that bare magnetars in normal galaxies provide a possible source for the origin of UHECR, both in total flux and in form of the spectrum. I have also given plausible arguments for how the particles accelerated by a magnetar embedded in its supernova and its magnetar wind nebula can escape without catastrophic losses, to contribute to the particle spectrum observed at the Earth.

This subject will remain observationally driven, most prominently by the AGASA (e.g. Takeda et al. 1998), Hi-Res (e.g. Sokolosky 1998, Abu-Zayyad et al. 2002b) and Auger (e.g. Boratov 1997, Cronin 2001a,b) experiments, which provide a window into the spectrum at energies above $`10^{20.5}`$ eV. The theory described here suggests that there will be something to see at these super-GZK energies. Also, Auger, with several thousand km2 collecting area and with fluorescence and particle detectors recording events simultaneously at the same site, will have the opportunity to resolve individual UHE particle bursts, with enough statistics to measure a spectrum for one event, at least up to $`10^{20}`$ eV - see §7 for bookeeping on this model’s predictions. On the theoretical front, the theory of how the acceleration actually occurs, in the model presented here and in other schemes, urgently needs substantial improvements.

## 10 Acknowledgments

I am indebted to L. Bildsten, E. Brown, G. Farrar, A. Filippenko, C. Heiles, Y. Levin, C.-P. Ma, C. Max, C. McKee, A. Olinto, P.B. Price, S. Singh, A. Spitkovsky, and S. Woosley for enlightening discussions. I acknowledge partial support from NASA grants NAG5-12031 and HST-AR-09548.01-A, and from the Miller Institute for Basic Research in Science. As always, I thank California’s taxpayers for their indulgence.

## Appendix A Appendix - R-mode Gravitational Wave Losses

Gravitational wave losses modify the accelerated spectrum by shortening the time each star spends at a given rotation frequency, thus reducing the flux at energy $`E(\mathrm{\Omega })`$. I used gravitational wave spindown rates owing to permanent nonaxisymmetric quadupolar distortions of the stars in order to arrive at the spectra shown in Figure 1, and argued that these could be generated by strong internal magnetic fields. The R-mode instability (Andersson 1998) is another possible source of gravitational wave emission from a rapidly rotating neutron star. With a dimensionless mode amplitude $`\alpha `$, the spindown rate due to gravitational radiation is (Owen et al. 1998)

$$\dot{\mathrm{\Omega }}_{Rmode}=\frac{512}{315}\frac{GI\mathrm{\Omega }^5}{c^5}\left(\frac{\mathrm{\Omega }R}{c}\right)^2\alpha ^2,$$
(A1)
where $`R`$ is the stellar radius. The same calculation as in §3, with R mode losses replacing those due to a static quadrupole moment, leads to

$$\frac{dN_i}{d\gamma }=\frac{dN_i}{dt}\left(\frac{dt}{d\mathrm{\Omega }}\right)mc^2\frac{d\mathrm{\Omega }}{dE}=\frac{9}{4}\frac{c^2I}{Ze\mu \gamma }\left[1+\left(\frac{\gamma }{\gamma _g}\right)^2\right]^1,$$
(A2)
replacing expression (23), assuming $`\alpha =`$ constant during the spindown, as in the model of Owen et al.. In (A2), the characteristic energy above which gravitational wave losses affect the spectrum now is

$$E_g=m_ic^2\gamma _g=\frac{Ze\eta \mu }{c^2}\left(\frac{104}{384}\frac{\mu ^2c^4}{GI^2R^2\alpha ^2}\right)^{1/2}=6\times 10^{20}\frac{Z\eta _1\mu _{33}^2}{I_{45}R_{10}\alpha }\mathrm{eV},$$

with $`R_{10}=R/10`$ km. If $`\alpha `$ were large, $`E_g(\alpha )`$ might be low enough to yield a steepened high energy injection spectrum and a GZK pileup.

However, recently Arras et al. (2003) provided the theory for $`\alpha `$ as a function of $`\mathrm{\Omega }`$, showing that damping in the neutron star’s interior, through coupling between R modes and internal g modes, keeps $`\alpha `$ small. In the notation used here, their theory provides $`\alpha =0.025\mathrm{\Omega }_4^{5/2}`$. The gravitational wave spindown rate is now proportional to $`\mathrm{\Omega }^{12}`$, instead of $`\mathrm{\Omega }^5`$ (static quadupole) or $`\mathrm{\Omega }^7`$ (Rossby modes with constant $`\alpha `$.) The injection spectrum still is given by (A2), but now with $`(\gamma /\gamma _g)^2`$ replaced by $`(\gamma /\gamma _g)^{9/2}`$, and with $`E_g=m_ic^2\gamma _g=10^{21}Z\eta _1\mu _{33}^{13/9}/(I_{45}R_{10})^{4/9}`$ eV. Gravitational wave losses due to the R mode instability therefore have negligible effect on the ability of the magnetar model to account for existing UHECR observations. They do suggest a very hard cutoff of the spectrum, above $`10^{21}`$ eV.