DRAFT January4,2011 PreprinttypesetusingLATEXstyleemulateapjv.5/25/10 BARYON LOADED RELATIVISTIC BLASTWAVES IN SUPERNOVAE Sayan Chakraborti1 & Alak Ray1 DepartmentofAstronomyandAstrophysics,TataInstitute ofFundamental Research, 1HomiBhabhaRoad,Mumbai400005,India DRAFTJanuary 4, 2011 ABSTRACT We provide a new analytic blastwave solution which generalizes the Blandford-McKee solution to 1 arbitraryejectamassesandLorentzfactors. Untilrecentlyrelativisticsupernovaehavebeendiscovered 1 onlythroughtheirassociationwithlongdurationGammaRayBursts(GRB).Theblastwavesofsuch 0 explosions are well described by the Blandford-McKee (in the ultra relativistic regime) and Sedov- 2 Taylor(inthenon-relativisticregime)solutionsduringtheirafterglows,astheejectamassisnegligible n in comparison to the swept up mass. The recent discovery of the relativistic supernova SN 2009bb, a without a detected GRB,opens up the possibility of highly baryonloadedmildly relativistic outflows J which remains in nearly free expansion phase during the radio afterglow. In this work, we consider 3 a massive, relativistic shell, launched by a Central Engine Driven EXplosion (CEDEX), decelerating adiabatically due to its collision with the pre-explosion circumstellar wind profile of the progenitor. ] Wecomputethesynchrotronemissionfromrelativisticelectronsintheshockamplifiedmagneticfield. E This models the radio emission from the circumstellar interaction of a CEDEX. We show that this H model explains the observedradio evolution of the prototypicalSN 2009bband demonstrate that SN . 2009bbhadahighlybaryonloaded,mildlyrelativisticoutflow. Wediscusstheeffectofbaryonloading h on the dynamics and observational manifestations of a CEDEX. In particular, our predicted angular p sizeofSN2009bbisconsistentwithVLBIupperlimitsonday85,butispresentlyresolvableonVLBI - o angular scales, since the relativistic ejecta is still in the nearly free expansion phase. r Subject headings: gamma rays: bursts — supernovae: individual (SN 2009bb) — shock waves — t s radiation mechanisms: non-thermal — techniques: interferometric a [ 2 1. INTRODUCTION pansionofa non-relativisticsupernovablastwave,inter- v Ultrarelativisticbulkmotionofmatterparticlesinas- acting with the surrounding circumstellar medium, was 4 trophysical settings is implied most notably in Gamma found by Chevalier (1982); Nadezhin (1985). Once the 8 RayBursts(GRBs)(seePiran(1999,2004)forreviews). blast wave sweeps up more CSM material than its own 7 Afterglows from GRBs are generated from the emission rest mass, the self-similar solutions of non-relativistic 2 blast waves are described in the Newtonian regime by by relativistic shocks that result from slowing down of a . the Sedov (1946)von Neumann (1963) Taylor(1950) so- 1 relativisticshellbythethemediumsurroundingthepro- lution. 1 genitor star that exploded. Similar interaction of stellar In this paper we provide an analytic solution of the 0 material(ejecta)fromanexplodingstarwiththecircum- standard model of relativistic hydrodynamics (see e.g. 1 stellarmatter (CSM)givesrisetonon-relativisticshocks Piran 1999; Chiang & Dermer 1999) for an adiabatic : in core collapse supernovae. v blastwave. Here, the exploding shell decelerates due to Fluid dynamical treatment of ultra-relativistic spher- Xi ical blast waves mediated by strong shocks2 has been inelastic collision with an external medium. That is, we providethesolutionforanarbitraryLorentzfactorofthe given by Blandford & McKee (1976, 1977). They de- r expanding supernova shell. The need for sucha solution a scribe a similarity solution of an explosion of a fixed whichcanhandleatrans-relativisticoutflowismotivated amount of energy in a uniform medium. This includes by the discovery of SN 2009bb, a type Ibc supernova anadiabaticblastwaveandanimpulsiveinjectionofen- without a detected GRB which shows clear evidence of ergy on a short timescale as well as an explosion where a mildly relativistic outflow poweredby a central engine the totalenergyincreaseswith time, suggestingthat the (Soderberg et al. 2010). SN 2009bb-like objects (Cen- blast wave has a continuous central power supply. An- tral Engine Driven Explosions, hereafter CEDEX) dif- other important model considered by them is that of fer in another significant way from classical GRBs: they a blast wave propagating into a spherically symmetric are highly baryon loaded explosions with non-negligible wind. On the other hand, the initial nearly free ex- ejecta masses. Our new analytic blastwave solution therefore generalizes the Blandford & McKee (1976) re- [email protected],[email protected] sult, in particular their impulsive, adiabatic blast wave 1Also at Institute for Theory and Computation, Harvard- in a wind-like ρ r−2 CSM. SmithsonianCenterforAstrophysics,60GardenSt.,Cambridge, ∝ MA02138, USA The new class of relativistic supernovae without a de- 2 Here the strong shocks are those where the random kinetic tected GRB, e.g. SN 2009bb, relaxes two important energyperparticlebehindtheshockismuchgreaterthantheun- and well-known constraints of the GRBs, namely the shocked medium, i.e., p2/n2 ≫p1/n1, the subscripts 2 and 1 de- compactness problem and baryon contamination. Thus note the shocked gas and unshocked gas respectively and p,n are thepressureandnumberdensities. highly baryon loaded mildly relativistic outflows can re- 2 Chakraborti & Ray main in the nearly free expansion phase during the ra- ation with long duration GRBs. An acceptable model dio afterglow. Therefore we consider the evolution of a for GRBs must find a way to circumvent the compact- massive relativistic shell launched by a CEDEX, experi- ness and baryon contamination problems. The observed encing collisional slowdown due to interaction with the rapid temporal variability in GRBs imply a very com- pre-explosion circumstellar wind of the progenitor. We pact source. The compactness problem pointed out by calculate the time evolution of the radius, bulk Lorentz Ruderman (1975) and Schmidt (1978), indicates that factor and thermal energy of the decelerating blastwave. γ-ray photons of sufficiently high energies will produce Our solution reduces to the Blandford & McKee (1976) electron-positron pairs and will not be able to come out solution, in the ultra-relativistic, negligible ejecta mass due to the resulting high opacity. However, the ob- limit. We also quantify the density of accelerated elec- servation of high energy photons from GRBs was rec- trons,amplifiedmagneticfieldsandhence the radiosyn- onciled with such theoretical constraints by Goodman chrotron emission from the blastwave of a CEDEX. (1986) and Paczynski (1986), allowing the high en- In Section 2 we outline the conditions that require ergy photons to come out from a relativistic explosion. ultra-high bulk Lorentz factors and very low ejecta This requires initial bulk Lorentz factors of the radi- masses in classical Long GRBs. We argue why these ating shell, γ & 102 (Piran 1999). Shemi & Piran 0 constraints are relaxed in the case of mildly relativistic (1990) pointed out that in the presence of even small outflowsdetected inSN2009bblike events,thathaveno amounts of baryons or baryon contamination, essen- detected GRBs. It provides the motivation for the ana- tially the entire energy of the explosion gets locked up lytic solution derived here. In Section 3 we develop the in the kinetic energy of the baryons, leaving little en- analyticsolutionofarelativisticblastwavelaunchedbya ergy for the electromagnetic display. This problem was CEDEX,forthecollisionalslowdownmodeldescribedby solvedbyRees & Meszaros(1992)andindependently by Piran (1999). We use this solution to show that the SN several authors (Narayan et al. 1992; Rees & Meszaros 2009bbblastwaveissubstantiallybaryonloadedandre- 1994; Paczynski & Xu 1994), by considering the recon- mains inthe nearlyfree expansionphase throughoutthe version of the kinetic energy of the fireball into radia- 1 year of observations (Figure 1). In section 4 we dis- tion, due to interaction with an external medium or via ∼ cuss the relativistic blast wave energetics. We quantify internal shocks. However the allowed initial mass is still the amount of shock accelerated electrons and magnetic very small (M 10−6M⊙), as too much baryonic mass ≈ field amplification in a CEDEX blast wave. In Section will slow down the explosion and it will no longer be 5 we use this information to model the radio spectrum relativistic (Piran 1999). Hence, observation of a short andlightcurveofaCEDEXinthe nearlyfreeexpansion bright pulse of γ-ray photons in GRBs require a very phase. In section 6 we provide expressions to deduce small amount of mass to be ejected with a very high the blastwaveparametersfromthe observedradiospec- bulk Lorentz factor. trum. We also compare our blastwave solution to other TheburstintheGRBitselfresultsfromtheconversion knownsolutionsrelevantforrelativisticblastwavesinthe of kinetic energy of ultra-relativistic particles or possi- literature (Section 7). In section 8 we discuss the impli- bly the electromagnetic energy of a Poynting flux to ra- cations of baryon loading in determining the evolution diation in an optically thin region. An inner engine is and observational signatures of a CEDEX. In the Ap- believed to accelerate the outflow to relativistic speeds, pendix we provide the analytic expressions for the tem- although the engine may remain hidden from direct ob- poralevolutionoftheblastwaveparameters(SectionA). servations. The “afterglow” on the other hand results We demonstrate that our solution reduces to the ultra- fromtheslowingdownofarelativisticshellontheexter- relativisticBlandford & McKee(1976)solutionforacon- nal medium surrounding the progenitor star. There can stant velocity wind, in the low mass, ultra-relativistic also be an additional contribution to the afterglow from limit (Section B). In Section C we provide “stir-fry ex- the inner engine that powers the GRB, since the engine pressions”3 for the radio spectrum of a CEDEX with may continue to emit energy for longer duration with a known initial blast wave parameters. In Section D we lower intensity and may produce the earlier part of the inverttheproblemandprovidehandyexpressionsfores- afterglow, say the first day or two in GRB 970228 and timating the initial parametersfor a CEDEXfromradio GRB 970508 (Piran 1999; Katz et al. 1998). observations. In Section E we provide an expression to In supernovae associated with GRBs, after a brief computeangularsizeofaCEDEXfrommulti-bandradio high energy electromagnetic display, the relativistic spectrum. Using the radio spectrum of SN 2009bb, we ejecta continues to power a long lived radio afterglow show that our predicted angular size is consistent with (Kulkarni et al. 1998; Soderberg et al. 2004). Since the the reported upper limits (Bietenholz et al. 2010) from emergence of γ-ray photons, early in its evolution, con- VLBI, but should be resolvable presently. strainstheinitialrelativisticejectamasstobeverysmall 2. RELATIVISTIC OUTFLOWINGRBSANDCEDEX and the initial bulk Lorentz factor to be very large, the relativistic ejecta sweeps up more circumstellar material Relativistic supernovae have, until recently, been dis- than its own rest mass by the time the radio afterglow covered only through their temporal and spatial associ- isdetected. The evolutionofthe radiativeblastwavehas 3Wenotethatasimilarterm“TVDinnerEquations”appeared been described by Cohen et al. (1998). During the ra- intheAstrophysicalliteratureinapaperonGRBs(Rhoads1999) dio afterglow phase if radiative losses do not take away todenoteequationswhere“numericalvaluesforphysicalconstants a significant fraction of the thermal energy, the blast- havebeeninserted,sotheyarereadytousewithoutfurtherprepa- wave may be treated as adiabatic. Under these condi- ration”; in the present work, “stir-fry” expressions would allow tions, the evolution of the blastwave is well described userstopickandchooseusableformulaederivedorapproximated in this work and quickly (re)assemble their own combination of by the Blandford & McKee (1976) solution if the blast- usefulinputparameters. wave remains ultra-relativistic or by the Sedov-Taylor Baryon Loaded Relativistic Blastwaves 3 solution if the blastwave has slowed down into the New- and tonian regime. The interaction of GRBs with their cir- dE =c2(γ 1)dm, (2) cumstellar wind has been discussed by Chevalier & Li − (2000). The spectra and light curves of GRB after- respectively, where dE is the kinetic energy converted glowshavebeencomputedintheultra-relativisticregime into thermal energy, that is the energy in random mo- by Sari et al. (1998) and in the Newtonian regime by tions as opposed to bulk flow, by the infinitesimal colli- Frail et al. (2000). sion. In the adiabatic case this energy is retained in the Soderberg et al.(2010)haverecentlydiscoveredbright shellandwehavetheanalyticrelation(fromPiran1999) radio emission, associatedwith the type Ibc SN 2009bb, requiring a substantial relativistic outflow powered by m(R) = (γ 1)1/2(γ +1)1/2 (3) a central engine. The search for a γ-ray counterpart, 0 0 M − − 0 spatially and temporally coincident with SN 2009bb, γ in data obtained from the all-sky Inter Planetary Net- (γ′ 1)−3/2(γ′+1)−3/2dγ′. work (IPN) of high energy satellites, did not detect ×Zγ0 − a coincident GRB (Soderberg et al. 2010). The ab- At the age of 17 days, the radio luminosity of SN sence of detected γ-raysfrom this event relaxes the con- 2009bb was L 5 1028 ergs s−1 Hz−1 at 8.4 GHz straints of very high Lorentz factors and very low mass ν ∼ × (Soderberg et al. 2010). This implies a radio luminosity in the relativistic ejecta. In fact, the radio emitting out- ofνL 4.2 1038 ergss−1. Comparedwith the energy flow in SN 2009bbis mildly relativistic (Soderberg et al. inrelaνti∼vistic×electronsof1.3 1049ergs(Soderberg et al. 2010) and has a large baryonic mass coupled to it, as 2010), this gives a radiation×timescale of 103 years at evidenced by the nearly free expansion for 1 year ∼ ∼ age 17 days, which is very large compared to the age of (Figure 1). Since the swept up mass is still smaller the object. Hence, it is safe to compute the dynamics of than the rest mass of the relativistic ejecta, the evolu- theobjectintheadiabaticlimit. Wealsonotethatsince tion of the blastwave can neither be described by the the availablethermalenergy initially growswith time as the rapidly decelerating Blandford & McKee (1976) so- E R t (shown later) and the luminosity falls of as lution nor the Sedov-Taylor solution. Dermer & Chiang νL∝ t∝−1 (shown later) the adiabaticity condition only (1998); Chiang & Dermer (1999) have considered the ν ∝ gets better with time. However, the radiation timescale emission from a collisionally decelerating blast wave wouldhavebeensameastheageatt 70seconds,which as it transitions from nearly free expansion to ∼ is comparable to the duration of very Long Soft GRBs the Blandford & McKee (1976) self similar solution. and XRFs. Huang et al. (1999) have discussed the subsequent tran- For a circumstellar medium set up by a steady wind, sition to a Sedov-like phase. whereweexpectaprofilewithρ r−2,wehavem(R)= Hence, the discovery of relativistic supernovae, such ∝ AR where A is the mass swept up by a sphere per unit as SN 2009bb, without detected GRBs, relaxes the con- radialdistance. A M˙ /v can be set up by a steady straints from the compactness problem and baryon con- mass loss rate of M≡˙ withwiandvelocity of v from the tamination. This motivates the study of highly baryon wind pre-explosionCEDEXprogenitor,possiblya WolfRayet loaded mildly relativistic outflows which can remain in star. Note that our definition of A is equivalent to 4πq the nearly free expansion phase during the radio after- in Equation (3 and 5) of Chevalier (1982). Integrating glow. the right hand side then gives us 3. RELATIVISTIC BLASTWAVESOLUTION AR γ γ2 1 = 0 − γ . (4) We use the simple collisional model described by M0 pγ2 1 − 0 Piran (1999); Chiang & Dermer (1999) where the rela- − tivistic ejecta forms a shell which decelerates through Forthephysicallyrelevanptdomainofγ >1thisequation infinitesimal inelastic collisions with the circumstellar has only one analytical root at wind profile. The initial conditions are characterized γ M +AR 0 0 γ = , (5) by the rest frame mass M of the shell launched by a 0 M2+2Aγ RM +A2R2 CEDEX and its initial Lorentz factor γ . In contrast 0 0 0 0 to the self similar solutions, which describe the evolu- which gives the epvolution of γ as a function of R. This tion away from the boundaries of the independent vari- is similar to Equation 8 of Chiang & Dermer (1999) for ables (Barenblatt & Zel’Dovich 1972) and the need for a adiabatic spherical blastwave in a wind like (ρ r−2) proper normalizationto getthe correcttotalenergy,our CSM profile. The amount of kinetic energy con∝verted set of initial conditions directly fixes the total energy as into thermal energy when the shell reaches a particular E0 =γ0M0c2. RcanbeobtainedbyintegratingEquation(2)aftersub- stituting for γ from Equation (5) and dm = AdR, to 3.1. Equations of Motion get The shell slows down by a sequence of infinitesi- E =c2( M AR (6) mal inelastic collisions with the circumstellar matter. 0 − − The swept up circumstellar matter is given by m(R). Conservation of energy and momentum give us (see + M02+2Aγ0RM0+A2R2 . Chiang & Dermer 1999; Piran 1999), q (cid:19) dγ dm 3.2. Evolution in Observer’s Time = , (1) γ2 1 −M − 4 Chakraborti & Ray The evolution of R and γ can be compared with ob- Appendix B) or Sedov like phase, depending upon the servations once we havethe time in the observer’sframe Lorentz factor at that time. thatcorrespondstothecomputedRandγ. Foremission 4. BLASTWAVEENERGETICS alongthelineofsightfromablastwavewithaconstantγ the commonly used expression (Meszaros & Rees 1997) We shall now use the blastwave solution developed in is the previous section to predict the radio evolution of a R t = . (7) CEDEX with a relativistic blastwave slowing down due obs 2γ2c to circumstellar interaction. For the prototypical SN However, Sari (1997) has pointed out that for a deceler- 2009bb,the blastwavewas only mildly relativistic at the ating ultra-relativistic blastwavethe correctt is given time of the observed radio afterglow. In the absence of obs by the differential equation a significant relativistic beaming, the observerwould re- dR ceive emission from the entire shell of apparent lateral dtobs = 2γ2c. (8) extent Rlat at a time tobs given by dR lat dt = , (15) We substitute γ from Equation (5) and integrate both obs βγc sides to get the exact expression which is valid even in the mildly relativistic regime. Be- R(M0+Aγ0R) causeofitsmildlyrelativisticoutflow,allexpressionsde- t = . (9) obs 2cγ (γ M +AR) rivedforSN2009bb,usetheaboveEquationratherthan 0 0 0 Equation 8 which is applicable in the ultra-relativistic Note that, this reduces to the Meszaros & Rees (1997) case. IneithercaseradioobservationsofSN2009bbmea- expression only in the case of nearly free expansion and sure essentially the transverse R not the line of sight lat deviates as the shell decelerates. In the rest of the work R. Integrating term by term gives us the time evolution weusettoindicatethetime t intheobserver’sframe. obs of the lateral radius as Inverting this equation and choosing the physically rele- vant growing branch, gives us the analytical time evolu- R =c γ2 1t (16) tion of the line of sight blastwave radius, as lat 0 − q 1 2 Ac2γ3 γ2 1 t2 R= ( M0+2Acγ0t (10) 0 0 − +O t3 . 2Aγ0 × − − (cid:16) Mp0 (cid:17) The thermal energy available when the sh(cid:0)ell(cid:1)has moved + 8AcM0tγ03+(M0−2Acγ0t)2 , out to a radius R is given exactly by Equation(6), how- q (cid:19) everitisagainconvenienttolookatitsTaylorexpansion, in the ultra-relativistic regime. This can now be substi- tuted into Equations (5 and 6) to get the time evolution E =Ac2(γ 1)R (17) 0 of the Lorentz factor γ and the thermal energy E (see − A2c2 γ2 1 R2 Appendix). This gives us a complete solution for the 0 − +O R3 . blastwavetime evolution,parametrizedbythe values for − 2M (cid:0) (cid:0) 0 (cid:1)(cid:1) γ0, M0 and A. 4.1. Electron Acceleratio(cid:0)n (cid:1) 3.3. Series Expansions Sari et al. (1998) give the minimum Lorentz factor γm of the shock accelerated electrons as Even though an analytical solution is at hand, it is p 2 m p instructive to look at the Taylor expansions in time for γ =ǫ − γ. (18) m e p 1 m the relevant blastwave parameters of radius (cid:18) − (cid:19) e At the mildly relativistic velocities seen in SN 2009bb, 4 Ac2γ3 γ2 1 t2 R=2cγ2t 0 0 − +O t3 , (11) thepeaksynchrotronfrequencyofthelowestenergyelec- 0 − (cid:0) M(cid:0) 0 (cid:1)(cid:1) trons are likely to be below the synchrotron self absorp- Lorentz factor (cid:0) (cid:1) tion frequency. This explains the ν5/2 low frequency be- 2 Acγ2 γ2 1 t γ =γ 0 0 − +O t2 , (12) havior of the spectrum. Hence, considering a electron 0− (cid:0) M(cid:0) 0 (cid:1)(cid:1) distribution with an energy spectrum N0E−pdE, which and thermal energy (cid:0) (cid:1) we assume for simplicity to be extending from γ m c2 m e to infinity, filling a fraction f of the spherical volume of E =2Ac3(γ0−1)γ02t (13) radius R, we need an energy 2A2c4(3γ0−2)γ03 γ02−1 t2 +O t3 . E = 4f(γmmec2)1−pN0πR3. (19) − M0 (cid:0) (cid:1) e 3(p 1) Equation(12) immediately tells us that the b(cid:0)las(cid:1)twaveof − a CEDEX starts out in a nearly free expansion phase, If a fraction ǫ E /E of the available thermal energy e e ≡ but slows down significantly by the time t , when the goes into accelerating these electrons, then for the lead- dec firstnegativetermintheTaylorexpansionforγ becomes ing order expansion of E in R the normalization of the equal to γ , where electron distribution is given by 0 M t = 0 . (14) 3Ac2ǫ (γ 1)(γ m c2)2 dec 2(Acγ (γ2 1)) N e 0− m e , (20) 0 0 − 0 ≃ 2fπR2 This signals the end of the nearly free expansion phase for p=3, as inferred from the optically thin radio spec- and the solution enters a Blandford McKee like (see trum of SN 2009bb (Soderberg et al. 2010). Baryon Loaded Relativistic Blastwaves 5 In the rest of this work we have made the simplifying where c , c and c are constants given by Pacholczyk 1 5 6 assumptionof a constant ǫ . The physics of shock accel- (1970) and D is the distance to the source. Similarly, e eration is unlikely to change during time of interest as the optically thin flux is given by relevant parameters such as the shock velocity does not changemuch. Moreover,we show later in this workthat 4πfR3 ν −(p−1)/2 F = c N B(p+1)/2 , (24) the light curve of a CEDEX is controlled by the com- ν 3D2 5 0 2c petition between a decreasing synchrotron flux at high (cid:18) 1(cid:19) frequencies and a decreasing optical thickness (hence in- In mildly relativistic cases such as SN 2009bb, we may receive radiation from the entire disk projected on the creasing flux) at low frequencies. The major light curve sky, hence R is to be understood as R . Substituting features of a shifting peak frequency and a nearly con- lat forN andB fromEquations(20and22)andtheleading stant peak flux, arise from this effect rather than micro- 0 physical processes that may change ǫ . orderexpansionfortheprojectedlateralradius(Eq. 16), e we have the optically thick flux as 4.2. Magnetic Field Amplification Ifa magnetic fieldofcharacteristicstrengthB fills the νt 5/2 F c2c (γ 1)(γ +1)π (25) same volume, it needs a ma1gnetic energy of ν ≃ 5 0− 0 (cid:18)c1(cid:19) ! EB = 6B2fR3. (21) 4c D2 24AǫB 1/4 . dIfenasfirtya,cttihoennǫtBhe≡chEaBra/cEtergiosteiscimntaogntheteicmfiaegldneitsicgievnenerbgyy , 6 (cid:18)γ0f +f(cid:19) ! Similarly, the optically thin flux can be expressed as c 6Aǫ (γ 1) B 0 B − (22) ≃ Rs f 24A2c3c1c5ǫBǫe(γ0 1)2(γmmec2)2 F − . (26) ν ≃ (νt)D2f γ2 1 In the rest of this work, following the argument in the 0 − previous sub-section, we have made the simplifying as- These equations together prpovide the flux density of a sumption of a constant ǫB. This explains the observed CEDEX as a function of time and frequency. B R−1 behavior (Figure 2) seen in SN 2009bb. These ∝ observations strengthens the case for a nearly constant 5.1. SSA Peak Frequency ǫ . This feature of an explosion within a ρ r−2 wind B ∝ Thetransitionfromtheopticallythicktoopticallythin profile is also seen in several non-relativistic radio su- regime happens at the peak frequency ν . At a fixed pernovae. The normalization is given by the progenitor p observation frequency this is reached in time t . The mass loss parameter, initial Lorentz factor, filling factor p conditionfortheSSApeakmaybeobtainedbyequating of the electrons and the efficiency with which the ther- the optically thin and thick flux, to get mal energy is used in amplifying magnetic fields. Note that this phase would last only as long as the expansion hriseelnnaceteiaoranlyCifnEreDSeNEoXr2t0w0<i9tbhtdbeaca.hriTgguhheleysrebffooarrreyaotnhtdeleocoabd&seedr1voeyudetaflBrow−a.nRd νptp ≃223/1435/14cc6ǫes4 γ0Af9ǫ+5Bfqc71(γ02−1) Note that the highest energy to which a cosmic ray 1 2/7 proton can be accelerated is determined by the B R (γ m c2)2 (27) product (Hillas 1984; Waxman 2005). We have argued × m e πf(γ +1)2 0 (cid:19) elsewhere (Chakraborti et al. 2010) that the mildly rel- All quantities on the right hand side are constants or ativistic CEDEX are ideal for accelerating nuclei to the parameters of the problem which are uniquely fixed for highest energies in sufficient volumetric energy injection a particular relativistic supernova. Hence, during the rates and independent arrival directions to explain the post GZK cosmic rays. The aboveB R−1 dependence nearlyfreeexpansionphaseνp ∝t−p1andthepeakmoves ∝ tolowerandlowerfrequenciesastheplasmaexpandsand seen in Fig 2 implies that a CEDEXlike SN 2009bbwill becomes optically thin with time. This behavior is also continue to accelerate cosmic rays up to a nearly con- seen in the radio spectra of SN 2009bb. stant very high energy for long times, i.e. for the entire The range in the values of ν t is then weakly depen- time that this (nearly free expansion) phase lasts. p p dent on the initial bulk Lorentz factor γ and the mass 0 5. ELECTRONSYNCHROTRONRADIOSPECTRUM loss rate, parametrized by A. The above expression can be used to select the frequency and cadence of radio It is expected that the radio emission from a CEDEX follow-upsoftypeIbcsupernovaefordetectingCEDEXs. will be produced by synchrotron emission of shock ac- celerated electrons in the shock amplified magnetic field 5.2. SSA Peak Flux quantified in the previous section. Once we have the evolution of N and B, we can follow Chevalier (1998) The peak flux density can now be obtained by substi- 0 to express the Synchrotron Self Absorbed (SSA) radio tuting the expression for νptp into that for Fν, to get, spectrum from a CEDEX in the optically thick regime using Rybicki & Lightman (1979) as Fνp ≃ 219/1439/14c19/7c5ǫ9B/14ǫ5e/7(Af(γ0−1))19/14 πR2c ν 5/2 ((cid:16)γ m c2)10/7π2/7 / c2/7D2f2 . (28) F = 5B−1/2 , (23) × m e 6 ν D2 c6 (cid:18)2c1(cid:19) (cid:17) (cid:16) (cid:17) 6 Chakraborti & Ray This peak flux is obtained by keeping only the leading which gives us the initial Lorentz factor γ in terms of 0 orders,foralmostfreeexpansion,isnearlyaconstantfor the early-time(t t )SSApeak frequencyν andthe dec p ≪ for t t . As the blastwave slows down, the decay of peak flux F . In the non-relativistic limit this gives us dec νp ≪ the peak flux can be obtained exactly from the analyti- cal expression for R(t). However for simplicity a Taylor 3c8ǫ (D2F )9 1/19 c v 2 6 B νp 1 , (34) expansion can be written down as ≃ π8c9ǫ f(γ m c2)2 ν t (cid:18) 5 e m e (cid:19) (cid:18) p p(cid:19) FFννpp((0t)) =1− 5(cid:0)AMcγ003(cid:1)t +O(cid:16)t3/2(cid:17) (29) wdǫeeh/niǫccBhe aoissnivnγ∝senαas−int1di/v1ei9s.totTothhbeeeinesqfoeulrvirpeedadrtsγie0tlifohncaospnaasirwsatemeanekttledyrewpαeitn≡h- m The peak flux may also come down if synchrotronlosses Equation (18) giving us γ 1.16 for SN 2009bb. Not 0 fortheelectronsaresignificant. Thiswouldrequireaself ≃ doingso willincur anerrorinestimating the initial bulk consistent modeling of the time dependent acceleration Lorentz factor. Equations (33 or 34) may be used as a andsynchrotronlossesof the relativistic electrons in the direct test of whether a radio supernova has relativistic shock amplified magnetic field of a CEDEX blastwave. ejecta. The peak flux density F is then strongly dependent νp The equipartition parameter α for a given CEDEX onthebulkLorentzfactorγ andtheparametrizedmass 0 may be determined directly from radio observations if lossrate,A. Thissuggestsalargerangeinradiofluxesof a synchrotron cooling break is seen in the broad band CEDEXsandhenceincreasedsensitivityoftheExtended radio spectrum. Chandra et al. (2004a) have used this VeryLargeArray(EVLA)shouldhelpindetectingmore to determine the magnetic field in SN 1993J and its of these enigmatic objects. radial evolution (Chandra et al. 2004b), independent of the equipartition argument. Subsequently A (essentially 5.3. Synchrotron Cooling the scaled mass loss rate) may be obtained by fitting Rybicki & Lightman(1979)givethecharacteristicsyn- Equation (22) to the B-R data (e.g. Figure 2). A direct chrotron frequency of an electron with Lorentz factor methodofestimatingAwouldbetoeliminateγ0between γ 1 in a magnetic field B as equations (27 and 28) to get A 1.2 1012 gram cm−1. e ≫ q B ForatypicalWolfRayetwindve≃locity×of103kms−1 this ν(γ )=γγ2 e . (30) e e2πm c densityprofilemaybesetupbyaconstantmasslossrate Sari et al. (1998) provide the critiecal electron Lorentz of M˙ 1.9 10−6M⊙ yr−1. This is consistent with the ≃ × massmasslossrateofthe SN2009bbprogenitoralready factor γ , above which an electron will loose a signifi- c determined by Soderberg et al. (2010). cant portion of its energy within the age of the object, All CEDEXs at a given distance, in their nearly free as 6πm c e expansionphase,arecharacterizedbyoneoftheconstant γ = . (31) c σ γB2t γ curves (Figure 4) and one of the constant M˙ curves T 0 Substituting this into Equation 30, we get the syn- (again Figure 4). This F vs ν t diagnostic plot for νp p p chrotron cooling frequency as CEDEXs is essentially the relativistic analogue of Fig- 18πm cq ure4ofChevalier(1998)wherenon-relativisticradiosu- e e ν ν(γ )= . (32) c ≡ c σ γB3t2 pernovae of different types are seen to lie on one of the T Substituting thevalueofSN2009bbmagneticfieldat20 constant velocity curves. For example, SN 2009bb, lies days from Table 1 of Chakraborti et al. (2010), we get on a intermediate mass loss parameter A and a mildly the cooling frequency as 2.6 THz. Since the magnetic relativistic γ0. For the same initial γ0 as SN 2009bb, an field decays as B t−1, t∼he cooling frequency will grow objectencounteringlowercircumstellardensity,wouldbe asν t. Hence,∝synchrotroncoolingisunimportantfor significantlyfainterandfasterevolving,suchasthoseoc- c electr∝onsradiatingneartheSSApeakfrequencyν 7.6 curring near the points of intersection of the solid green p GHz, and they are indeed in the slow cooling re∼gime curve and the dashed red curve. These objects, possibly during the timescale of interest. at the faint end of the luminosity function of CEDEXs, The scaling relations for ν and ν show that they maybe discoveredby a dedicatedhighsensitivity search c p should have been comparable at around 1 day and forrelativisticoutflowsfromnearbytype Ibcsupernovae at a frequency of 140 GHz. Very earl∼y time ( 1 with the EVLA. day)millimeterwav∼eobservationsofCEDEXsareth∼ere- The ejecta mass can then be inferred from the decel- fore encouraged. They may revealthe presenceof a syn- eration timescale tdec fitting the time evolution of the chrotron cooling break in the radio spectrum. This may observed radii with the Taylor expansion (Equation 16) be used to determine the magnetic field independent of ofRlatobtainedfromourmodel. InFigure1wecompare the equipartition argument (see Chandra et al. 2004a). predictions from our model with the observed temporal evolutionoftheblastwaveradius. Theobservedradiiare 6. BLASTWAVEMODELINVERSION consistent with a nearly free expansion for an apparent lateralvelocityofγβ =0.527 0.022c. Asthe ejectahas The input parameters of the model are specified by ± not yet slowed down significantly, it is not possible to M , γ , A, ǫ and ǫ . Under the assumption of equipar- 0 0 e B determine M but we can nevertheless put a lower limit tition of the thermal energy between electrons, protons 0 and magnetic fields we have ǫ = ǫ = 1/3. We then on the ejecta mass. Models with M0 below 10−2.5M⊙ e B are ruled out at the 98.3% level by the observations of eliminate A between Equations (27 and 28) to get SN 2009bb. Hence, radio observations of a relativistic 3c8ǫ (D2F )9 2/19 c 2 supernovae can constrain all the physically relevant pa- γ2 1=4 6 B νp 1 , (33) rameters of our model. 0− π8c9ǫ f(γ m c2)2 (ν t )c (cid:18) 5 e m e (cid:19) (cid:18) p p (cid:19) Baryon Loaded Relativistic Blastwaves 7 7. COMPARISONWITHKNOWNSOLUTIONS and remain in equipartition with the electrons through- out the entire remnantvolume, the shock Lorentz factor Blandford & McKee (1976) (hereafter BM) analyze evolvesfor a homogeneousexternalmedium as: Γ r−3 the dynamics of both non-relativistic (NR) and ultra- ∝ (Blandford & McKee 1976). On the other hand, for the relativistic (ER) blast waves in the adiabatic impulsive adiabatic regime where the radiative losses do not tap (AI) approximationas well as steady injection (SI) from the dominant energy reservoir in protons and magnetic a central power supply. A classification of the shock dy- energy but only the electrons are responsible for the ra- namics and the corresponding synchrotron and inverse diation, one has: Γ r−3/2 for a homogeneous medium Compton radiation from the strong relativistic spheri- ∝ (Meszaros et al. 1998). cal shock in an ionized magnetized medium is given by For a NR blast wave of energy E into a medium of Blandford & McKee (1977). The dynamics of both NR density ρ in the Sedov-Taylor phase, it is possible to andERshocksaregovernedbytheenergyE oftheadia- 0 construct a characteristic velocity at a time t from the baticimpulsive(AI)explosion;whileinthecaseofsteady combination(E/ρ t3)1/5. Fortherelativisticproblem,as injectionofenergythedynamicsisgovernedbyLt,where 0 BMpointout,anadditionalvelocitythespeedoflightis L is the luminosity of the central power supply. Other introduced into the problem. The mean energy per par- assumptionsmadebyBMwerethattheshockisastrong ticle in the shocked fluid varies as Γ2, where one factor one (see Footnote 2) and that the magnetic field is not of Γ arise from the increase in energy measured in the dynamically important. co-moving frame and the second arise from a Lorentz The external medium into which the shock wave ex- transformation into the fixed (i.e. the observer) frame. pands can have radial variations of density. It can ei- When most of the energy is resides in recently shocked ther be of uniform density or could be stratified into ρ r−k, with k = 2 for a constant velocity wind. The particles, the total energy is proportionalto Γ2R3. Here ∝ R is the current shock radius. If the total energy con- shock can be either adiabatic or radiative. In adiabatic tained in the shocked fluid is to remain constant, then shocks the radiative mechanisms are slow compared to Γ2 t−3. Onecanconsideramoregeneralcaseinwhich hydrodynamics timescale, while in radiative shocks the ∝ theenergyissuppliedcontinuouslyatarateproportional radiative mechanisms are faster than the hydrodynamic to a power of the time and set: timescale (measured in the observer-frame). A fully ra- diative blast wave may radiate away all the thermal en- Γ2 t−m for m> 1. (35) ergygeneratedbythe shock,if forexample, the external ∝ − medium is composed mainly of electrons and positrons. For the m=3 case, corresponding to an impulsive injec- Cohen et al.(1998)categorizeshocksas“semi-radiative” tion of energy into the blast wave on a time scale short in which the cooling is fast, but only a fraction of the compared with R, the total energy can be shown to be energy is radiated away. This could take place in a col- given by: lisionless shock acceleration,where the shock distributes E =8πw1t3Γ2/17, (36) the internal energy between the electrons and protons. wherew istheenthalpyaheadoftheshock,sothatwith After the material passes behind the shock, there is Γ2 t−13, the total energy is indeed a constant in time no coupling, effectively at the low densities behind the ∝ in the fixed frame. shock, between the electrons and protons. During the BM also provide a solution for an ER adiabatic blast initial stage of the afterglow, the electrons may undergo wave in an external density gradient ρ n r−k but synchrotron cooling or Compton scatter off low energy 1 ∝ 1 ∝ with a cold, pressure-less external medium. This is the photons on timescales which are shorter than the dy- case for a blast wave propagating through a spherically namicaltimescale(Waxman1997b;Meszaros et al.1998; symmetricwind. Forthisexampletheyfindthatfork=2, Sari et al. 1998), but the protons may not radiate and m=1correspondstoanimpulsiveenergyinjectiondueto may remain hot. A narrow cooling layer may form be- ablastwaveinaconstantvelocitywind. Inthiscase,the hindtheshockasinafastcoolingscenarioevenwhenthe analytic solution to the adiabatic impulsive blast waves protonsremainadiabaticandonlytheelectronscool. As canbeobtainedform=3 k > 1,andthetotalenergy the cooling parameter ǫ (the fraction of the energy flux − − in this case is: lost in the radiative layer) increases, the matter concen- trates in a small shell near the shock in the Newtonian E =8πρ Γ2t3/(17 4k). (37) 1 − solution; however in the ER case, even in the adiabatic Againthis energyis constantto the lowestorderin Γ−2. case,thematterisconcentratedinanarrowshellofwidth R/Γ2 (where Γ = √2γ, is the ultra-relativistic shock We compareour solution to other wellknown analytic and numerical schemes in the two final Figures. The wave Lorentz factor in the observer’s frame), but this graphs are not really for SN 2009bb-like systems. As concentration in a narrower denser shell increases with they are for a baryon loaded explosion with a γ of 300 ǫ. 0 going into a wind-like CSM. There may not be any such Classical fireball models of GRBs have a fully radia- spherical explosions with these kind of energies. The tive stage during the (initial) γ-ray event with a radia- comparisoninthesefiguresarenotforanyphysicallyrel- tive efficiency near unity and the implied energy is typ- ically E 1051−52 erg and have a bulk Lorentz factor evantobjects(orscenarios)butratherforthepurposeof of Γ 10∼2 103. When the newly shocked electrons are showing the consistency of our analytical solutions with ∼ − previouslyreportednumericaltechniquesandknownreal radiating near the peak of the initial post-shock energy asymptotic or intermediate asymptotic solutions. γ ξ (m /m )Γ(t) (here ξ = O(1)), they retain high e ∼ e p e e In Figure 5 we compare our solution, in the ultra- radiativeefficiencyevenaftertheGRBoutburstforsome relativistic regime with nearly free expansion or the time (Meszaros et al. 1998). If the protons establish coasting solution at early times and a BM like solution 8 Chakraborti & Ray at intermediate times. The blast wave solution for the hadvery heavy mass loss soon before the explosiontook special case of adiabatic impulsive solution for a con- place,thenaninitiallyfreeexpansionoftheshockquickly stant velocity wind obtained in Section 3 is valid for ar- sweeps up a mass equal to or greater than M and the 0 bitraryvaluesofM andγ andthemasslossparameter explosionenterstheSedov-Taylorlikephase. Alternately 0 0 A. The corresponding solutions in the Blandford Mc- a low ejecta mass initially would also lead to the early Kee analysis can be obtained from our exact solution in onset of the S-T phase. the limit of small initial mass of the high velocity ejecta The longdurationGRB afterglowsare poweredbyex- M m(R)γ where m(R) is the swept-up mass in the plosions that have very small initial ejecta mass, typi- 0 0 wind≪. Thus Equation (B1) shows that our intermedi- cally 10−6 M⊙. Comparatively, the initial mass of the ate asymptotic expressions for R and γ have the same ejecta in SN 2009bb estimated in Section 6 is consider- time dependence as those of the BM solution. Simi- ably largerat M0 >10−2.5 M⊙, but not as largeas non- larly, the opposite limit M m(R)γ gives a nearly relativistic supernovae. These explosions are therefore 0 0 ≫ free expansion of the blast wave. Thus the relativistic baryonloaded. Asbaryonloadingisafactorthataffects CEDEX which had relatively low mass (although the the conversionof the impulsive release of energy into ki- prototypeSN2009bbwasstillsignificantlybaryonloaded neticenergyofthematteraroundthecentralsource,this compared to classical GRBs), underwent nearly free ex- canquenchthe emergenceofthe gamma-raysin aburst. pansion initially. This is the relativistic analogue of The substantial baryon loading in SN 2009bb compared the Chevalier (1982); Nadezhin (1985) phase of the non- to classical GRBs may in fact have played a significant relativistic supernovae. Eventually, an ultra-relativistic roleinthecircumstancethatnogamma-rayswereinfact CEDEX would enter the Blandford-McKee phase (the seen from this CEDEX, despite a thorough search by a relativistic analogue of the Sedov-Taylor phase) when it suite of satellites in the relevant time window. At the has swept up enough mass from the external medium same time, the CEDEX, such as SN 2009bb, are also surrounding the CEDEX. unique in that they span a region of γ which is trans- 0 In Figure 6 we compare our solution, in the mildly- relativistic. So far the LGRBs and the corresponding relativistic regime with a BM like solution at intermedi- outflows were being described in the extreme relativis- ate times and a Snowplough phase at late times. The tic limit as in BM. Here we provide a solution of the BM-like phase becomes unphysical when γ approaches relativistic hydrodynamics equations which is uniquely 1 (hence β becomes 0) and the Snowplough phase is tuned to the CEDEX class of objects like SN 2009bb. feasible only much later when β is less than 1. Both Because of the non-negligible initial ejecta mass of such Figures also compare our solution to a direct numerical a CEDEX, objects like this would persist in the free ex- solution of the equations of motion computed for this pansion phase for quite a long time into their afterglow. work, under under the relevant conditions as prescribed Sweeping up a mass equal to that of the original ejecta byChiang & Dermer(1999). TheseFiguresdemonstrate would take considerable time, unless the mass loss scale that patched solutions are much worse than the solu- in its progenitor was very intense, i.e. it had a large A. tion proposed here or even numerical solutions. While In this paper we also provide an analysis of the peak there can be a patch (though worse than our solution) radiofluxversustheproductofthepeakradiofrequency betweencoastingandBM,thereexistsnopatchbetween and the time to rise to the radio peak. There are loci in BM and Snowplough. While our analysis is for a wind thisplanethatarespannedbylow,intermediateandhigh like (ρ r−2) medium, see Kobayashiet al. (1999) for velocity explosions. Radio SNe, LGRBs and CEDEXs ∝ similar transitions seen in numerical simulations of the are seen to occupy different niches in this plane. Sim- evolutionofanadiabaticrelativisticblastwaveinanuni- ilarly, we provide expressions for the peak fluxes and form medium. peak radio frequencies in terms of mass-loss factors of Ourframeworkmaybe usedinthefuture tostudy the the explosions. We also invert their dependence on γ 0 transition of a relativistic blastwaveinto a non relativis- and M˙ in terms of the peak frequency, peak time and tic Sedov phase. Frail et al. (2000) have used the Sedov peak fluxes to interpret the parameters of the explosion. phaseforverylatetimecalorimetryofGRB970508. The We have also demonstrated that a seed magnetic field relationbetweentheblastwaveradiusandtheobserver’s amplifiedbytheshockasdescribedaboveplaysacrucial frametimealsoleadstotheappropriatenumericalfactor role in the generationof synchrotronradiation and have connectingthe twoforaρ r−2 exterior(See Appendix argued (Chakraborti et. al, in prep) that the accelera- ∝ B). NotethedifferentresultsintheliteratureduetoSari tion of cosmic rays to ultra-high energies is possible in (1997); Waxman (1997a) for an uniform media. SN 2009bb-like objects. 8. DISCUSSIONS Explosion dynamics where adiabatic, ultra-relativistic We thank Alicia Soderberg, Abraham Loeb and or trans-relativistic outflow takes place in a wind-like Poonam Chandra for discussions on mildly relativis- CSM can be classified in a parameter space spanned by tic supernovae. We thank Naveen Yadav for checking M , γ and A, i.e. the initial ejected mass, bulk Lorentz the manuscript and Swastik Bhattacharya for help with 0 0 factoroftheejectaandthemasslossparameteroftheex- mathematics. This work has made use of the Wolfram ternal wind established before the explosion, apart from Mathematica computer algebra system. This research the electron acceleration and magnetic field amplifica- wassupportedbytheTIFR11thFiveYearPlanProject tion efficiencies parametrized by ǫ and ǫ . At low γ , no. 11P-409. We thank the Institute for Theory and e B 0 we have the non relativistic supernovae, which usually Computation,HarvardUniversity forits hospitality. We have considerable ejecta mass (e.g. at least 0.1 M⊙, thankananonymousrefereefordetailedcommentswhich often much larger). If the progenitor of the SN has significantly extended this work, in particular for point- Baryon Loaded Relativistic Blastwaves 9 SN 2009bb MMM0=001==011-002--.125 MMMsssooolll Radio Light Curve 20 18 15 16 R (cm)lat1017 FFνν ((mmJJyy)) 1 005 468111024 2 0 50 100 tobs (days) 150 10 1016 200 1 Frequency (GHz) 101 102 Tobs (days) Fig.1.— Evolution of the blast wave radius Rlat, determined Fig.3.— The flux density Fν as a function of the frequency ν from SSA fit to observed radio spectrum, as a function of the and the time in the observer’s frame tobs, as predicted from the time tobs in the observer’s frame. The evolution is consistent modelproposedinthiswork. Thevaluesusedfortheparameters, with nearly free expansion. Note that the observations require arethoseofSN2009bb. M0&10−2.5M⊙. 1 SN 2009bb Radio SNe uss) Ga d ( Fiel netic 0.1 g Ma 1000 1015 1016 1017 1018 Radius (cm) 100 Fig.2.— Magnetic field as a function of blast wave radius, as determined from SSA fits. Blue dots represent size and magnetic fieldofradiosupernovaefromChevalier(1998)atpeakradiolumi- y) nosity. Redcrosses(with3σerror-bars)givethesizeandmagnetic mJ 10 lfiienledgoifveSsNth2e00b9ebsbtBat∝diffRe−re1n(tEeqpuoacthiso,nfr2o2m)fistp.ectralSSAfits. Red F (νp ing out to us the Chiang and Dermer papers, which we 1 were unaware of. We thank Charles Dermer for sharing their numerical results in these papers. We thank Tsvi SN2009bb 0.1 PiranandCharlesDermerfordiscussionsattheAnnapo- 0.1 1 10 lis GRB2010 meeting. (νp / 10 GHz)(tp / 20 day) Fig. 4.—DiagnosticFνp vsνptp plotofourfiducial CEDEXin the nearly free expansion phase, at a distance of 40 Mpc. Solid linesfromrighttoleftrepresentlociofconstantγ0=1.005(Red), 1.16(Green) and 3.0(Blue), traced by keeping the LHSof Equa- tionD1 fixed at the respective values. Dashed lines from rightto leftrepresent lociof constant M˙ = 0.2(Red), 1.9 (Green) and 20 (Blue)inunitsof10−6 M⊙yr−1,againtracedbykeepingtheLHS of Equation D2 fixed at the respective values. The cross marks thesetofvaluesderivedfortheprototypical SN2009bbspectrum fromSoderbergetal.(2010)at20days,with3σerror-bars. Radio observations canbeusedtodeducetheinitialbulkLorentzfactor γ0 ofaCEDEXandthepre-explosionmasslossrateM˙ fromthis diagram. 10 Chakraborti& Ray 1000 γ(R) Coasting Blandford-McKee Numerical γ 100 10 1011 1012 1013 1014 1015 R (cm) Fig. 5.— EarlytimeradialevolutionofthebulkLorentzfactor γ(R) in a CEDEX, derived in Equation 5 of this work (red solid line), as compared to an initial nearly free expansion or coasting phase (green dotted line) and the intermediate Blandford-McKee like phase (blue dashed line). Magenta squares represent a direct numericalsolutionoftheequationsofmotionfollowingthemethod prescribed by Chiang&Dermer (1999). Note the transition from acoastingphasetotheBlandford-McKeelikephase. β(R) Blandford-McKee 1 Snowplough Numerical β 0.1 0.01 1016 1017 1018 1019 1020 R (cm) Fig. 6.— Late time radial evolution of β(R) = v(R)/c in a CEDEX, derived using Equation 5 of this work (red solid line),ascomparedtoanintermediateBlandford-McKeelikephase (blue dashed line) and a final Snowplough phase (orange dot- ted line). Magenta squares represent a direct numerical solution of the equations of motion following the method prescribed by Chiang&Dermer(1999). NotethetransitionfromtheBlandford- McKeelikephasetoaNewtonianSnowploughphase.