Particle acceleration and magnetic field amplification in massive young stellar object jets
MMNRAS , 1–16 (2015) Preprint 24 February 2021 Compiled using MNRAS L A TEX style file v3.0
Particle acceleration and magnetic field amplification in massiveyoung stellar object jets
Anabella T. Araudo, , (cid:63) Marco Padovani, and Alexandre Marcowith ELI Beamlines, Institute of Physics, Czech Academy of Sciences, 25241 Doln´ı Bˇreˇzany, Czech Republic Astronomical Institute, Czech Academy of Sciences, Boˇcn´ı II 1401, CZ-141 00 Prague, Czech Republic INAF-Osservatorio Astrofisico di Arcetri, Largo E. Fermi, 5 - 50125 Firenze, Italy Laboratoire Univers et Particules de Montpellier (LUPM) Universit´e Montpellier, CNRS/IN2P3, CC72, place Eug`ene Bataillon,34095, Montpellier Cedex 5, France
Accepted XXX. Received YYY; in original form ZZZ
ABSTRACT
Synchrotron radio emission from non-relativistic jets powered by massive protostarshas been reported, indicating the presence of relativistic electrons and magnetic fieldsof strength ∼ − − . We show thatthe magnetic field in the synchrotron emitter can be amplified by the non-resonanthybrid (Bell) instability excited by the cosmic-ray streaming. By combining the syn-chrotron data with basic theory of Bell instability we estimate the magnetic field inthe synchrotron emitter and the maximum energy of protons. Protons can achievemaximum energies in the range . − . TeV and emit γ rays in their interactionwith matter fields. We predict detectable levels of γ rays in IRAS 16547-5247 andIRAS 16848-4603. The γ ray flux can be significantly enhanced by the gas mixing dueto Rayleigh-Taylor instability. The detection of this radiation by the Fermi satellite inthe GeV domain and the forthcoming Cherenkov Telescope Array at higher energiesmay open a new window to study the formation of massive stars, as well as diffu-sive acceleration and magnetic field amplification in shocks with velocities of about1000 km s − . Key words: acceleration of particles – radiation mechanisms: non-thermal – shockwaves – stars: jets – gamma-rays: general
Stars are formed within dense molecular clouds, accretingmatter onto the central protostar with the formation of acircumstellar disc and bipolar jets. These ejections are col-limated flows of disc/stellar matter accelerated by magneticfield lines (Blandford & Payne 1982; Shu et al. 1994), andmoving with speeds v j ∼ − km s − into the ambi-ent molecular cloud. Molecular matter from the cloud is en-trained by the jet, forming molecular outflows. Protostellarjets are thermal radio emitters in most of the cases (Angladaet al. 2018). However, radio emission with negative spectralindices has been detected in several protostellar jets suggest-ing a non-thermal nature of the emission (e.g. Garay et al.2003; Rodr´ıguez-Kamenetzky et al. 2017; Purser et al. 2016;Obonyo et al. 2019). In the particular case of the jet as-sociated with the massive protostar IRAS 18162-2048 andits famous Herbig-Haro (HH) objects HH 80 and HH 81, (cid:63) E-mail: [email protected]
Carrasco-Gonz´alez et al. (2010) reported on polarized radioemission and confirmed that the radiation is produced by thesynchrotron process. The detection of synchrotron radiationis an evidence that there is a population of mildly-relativisticor relativistic electrons in the jet.The cosmic-ray (CR) energy density in molecular cloudsis poorly known, but it should be at least of the same or-der of the CR energy density in the interstellar medium,i.e. ∼ − erg cm − (e.g. Ferri`ere 2001) . Therefore, CRswith energy densities larger than − erg cm − (see Sect. 3)have to be locally accelerated at the jet termination shocksor/and in the internal jet shocks. This result is in agreementwith Padovani et al. (2015, 2016), who have shown that thelarge ionization rate estimated from molecular line ratios inL1157-B1 (Podio et al. 2014) and OMC2 FIR4 (Ceccarelliet al. 2014; Fontani et al. 2017; Favre et al. 2018) can be Unless a strong source of CRs contributes to a local enhance-ment, this number is actually an upper limit as CRs are subjectto strong ionization losses while propagating in molecular clouds. © a r X i v : . [ a s t r o - ph . H E ] F e b A. T. Araudo et al. due to local acceleration of thermal particles to a relativis-tic regime in protostellar jets.Diffusive Shock Acceleration (DSA) is the most estab-lished mechanism to accelerate particles in sources whereshock waves are present, from solar flares (e.g. Chen et al.2015) to the outskirts of galaxy clusters (e.g. Petrosian2001). Magnetic turbulence near the shock allows particlesto diffuse back and forth the shock front gaining energy ateach crossing cycle (Axford et al. 1977; Krymskii 1977; Bell1978a; Blandford & Ostriker 1978). The maximum energythat particles can achieve in non-relativistic parallel shocksis usually determined by radiative losses or advection escapefrom the acceleration region. Araudo et al. (2007) showedthat completely ionized and fast protostellar jets can accel-erate particles up to TeV energies if Bohm diffusion applies.However, in a non-completely ionized medium, damping ofAlfv´en (or magnetohydrodynamic) waves by ion-neutral col-lisions can reduce the particle maximum energy (Drury et al.1996; Padovani et al. 2015, 2016).The detection of molecular emission lines such as H α ,[NII], [SII], and [OIII] in HH objects indicate that pro-tostellar jets are not completely ionized. Emission linesare usually associated with internal (low-velocity) shocksin jets emanating from low-mass protostars. However, insources like HH 80 and HH 81 molecular lines are co-spatialwith synchrotron radiation of locally accelerated electrons(Rodr´ıguez-Kamenetzky et al. 2019a). Therefore, strongshocks and scattering waves are present in non-completelyionized plasmas. As a result, short wavelength Alfv´en wavesare unlikely to scatter particles back and forth the shocksdue to their damping by ion-neutral collisions. Bell (2004) re-alized that under certain conditions non-resonant (hereafterNR) hybrid waves can grow faster than Alfv´en (resonant)waves. The NR (Bell) instability has two main advantages: ( i ) the magnetic field is amplified by orders of magnitude and ( ii ) NR waves are not strongly damped by ion-neutral colli-sions, even if the maximum growth rate is reduced (Revilleet al. 2007).In this study we show that the magnetic field in thesynchrotron emitter of high-mass protostellar jets can beamplified by the streaming of CRs (Bell 2004), as in super-nova remnants (e.g. Vink & Laming 2003) and jets in ActiveGalactic Nuclei (Araudo et al. 2015). By assuming that thenumber of non-thermal electrons and protons is the same, wefind that the energy density in non-thermal protons is largeenough to drive the Bell instability in the termination regionof protostellar jets. In addition, we show that detectable lev-els of gamma rays in the GeV domain are expected from pro-tostellar jets. We consider the sample of non-thermal lobes inthe termination region of high-mass protostellar jets stud-ied by Purser et al. (2016) with the Australian TelescopeCompact Array (ATCA) facility. We select the sources withradio spectral index . < α < . , where the radio flux den-sity at frequency ν is S ν ∝ ν − α . In addition, Obonyo et al. The Bohm diffusion regime is obtained when the mean-free pathof the particle imposed by angular scattering is equal to its Lar-mor radius. The term hybrid comes from the nature of the instability whichcan be described using both magnetohydrodynamic and kinetictheory.
Table 1.
Non-thermal MYSOs selected from the observationalstudies carried out by Purser et al. (2016) and Obonyo et al.(2019). From left to right we list the name of the source in theGreen and IRAS catalogues, the distance d , and the ionized massloss rate ˙ M i = X ˙ M j , where X is the ionization fraction and ˙ M j isthe total jet mass loss rate provided in the literature mentionedbefore. The values of X and jet velocity fixed by Purser et al.(2016) and Obonyo et al. (2019) to compute ˙ M j are X = . and0.4, and v j = km s − . The jet kinetic luminosity computed as L kin = ˙ M i v / is also listed.Source IRAS d ˙ M i X L kin name [kpc] [M (cid:12) yr − ] [erg s − ]G263.7434 08470-4321 0.7 6.4 × − . × G263.7759 08448-4343 0.7 1.3 × − . × G310.1420 13484-6100 5.4 1.1 × − . × G313.7654 14212-6131 7.8 4.1 × − . × G339.8838 16484-4603 2.7 5.3 × − . × G343.1261 16547-4247 2.8 3.5 × − . × G114.0835 23262+5834 4.2 4.4 × − . × (2019) found that 6 of the 15 massive young stellar ob-jects (MYSOs) observed with the Jansky Very Large Ar-ray (JVLA) show clear evidence of non-thermal emission,with spectral indices α > . . From this last sample we se-lect the source IRAS 23262+5834. In Table 1 we list theMYSOs selected for our study. They are located at a dis-tance . ≤ d / kpc ≤ . from Earth. The jet ionized mass-lossrate is × − ≤ ˙ M i / ( M (cid:12) yr − ) ≤ × − and the kineticluminosity is . × ≤ L kin / ( erg s − ) ≤ . × .The paper is organized as follows: in Sect. 2 we describethe physical properties of protostellar jets and the conditionfor them to be non-thermal emitters. In Sect. 3 the magneticfield and the content of relativistic particles is inferred fromsynchrotron emission. In Sect. 4 we discuss both the resonantand non-resonant instabilities in the context of young stellarjets . In Sect. 5 we briefly describe the non-linear saturationof the amplified magnetic field and determine its strength.In Sect. 6 we compute the maximum energy of protons andin Sect. 7 we estimate the γ -ray flux from our source sample.In Sect. 8 we present our conclusions. Gaussian-cgs units areused throughout the paper. The detection of thermal (free-free) and synchrotron emis-sion from protostellar jets indicates that ( i ) there are shocksthat heat up the plasma and accelerate particles, ( ii ) theregions in the jet where thermal and non-thermal emis-sion is radiated are at least partially ionized, and ( iii ) mag-netic fields of about 1-10 mG are needed to explain the de-tected synchrotron flux if equipartition arguments apply (seeSect. 3.1). In this section we provide some constraints on thetemperature T u , the ionization fraction X u , the ion density n i , and the magnetic field B u upstream of the shock of speed v sh by using the most updated data available in the liter-ature (see Table 2) and conservation equations. We define X u = ( n e + n i ) / n j , where n j = n i + n e + n n , n e and n n are the free We note that this section is plasma-physics oriented. It is writ-ten is a self-contained way. Unfamiliar reader can skip it.MNRAS , 1–16 (2015)
FA in jets from MYSOs electrons and neutral densities, respectively, and we assume n e = n i hereafter. Bacciotti & Eisl¨offel (1999) estimated the ionization frac-tion, X , in a sample of low-mass protostellar jets. They foundthat X generally decreases almost monotonically along thejet, hence the effects of the shocked materials over X seems tobe minor and we may identify X with X u . Typical values of X u between 0.02 and 0.4 have been derived by the authors. Mau-rri et al. (2014) analyzed Hubble Space Telescope arc-secondimages of the DG Tau jet, hence at a much closer locationto the central source than Bacciotti & Eisl¨offel (1999). Theyfound a low ionization fraction ( ≈ . ) at the base of thejet, then increasing to a plateau value of ≈ . at a distanceof 3 (cid:48)(cid:48) from the central source. The ionization fraction seemsthen to decrease again. By modeling the optical [OII]/[OI]line ratio in shocks at the base of the HN Tau jet, Har-tigan et al. (2004) find X u ≈ . . Te¸sileanu et al. (2012) de-veloped a magnetohydrodynamic (MHD) model of low-massprotostellar jets and adapted it to the source RW Aur. Theyfound that X u ≈ . − . in the central spine correctly fitsthe surface brightness fluxes in several forbidden lines. Theyattributed this ionization fraction to irradiation by X-raysfrom the protostellar magnetic corona.The temperature of the unshocked gas is difficult toconstrain as line brightness or line ratios are only weaklydependent on it. Garcia et al. (2001) investigated low-massMHD jets launched from accretion discs including ambipolardiffusion heating and several sources of ionization (X-raysand UV radiation as well as collisions). They found a typicaljet temperature T u ≈ K.In high-mass protostars both X u and T u are barelyknown or even constrained. Following Garcia et al. (2001) wecan reasonably expect the jets from high-mass protostars tohave higher mass loss rates and to be denser and cooler thanjets from low-mass protostars. Values of X u as small as − are found for jets with T u ≈ K. However, other sources ofheating like turbulence, MHD wave dissipation, magnetic re-connection, or weak shock heating have not been consideredin this work. In a very recent study, Fedriani et al. (2019)found that the ionization fraction in the jet of the high massprotostar G35.20-0.74N is ∼ . , similar to the values forlow-mass protostars. Purser et al. (2016) assumed X u = . for jets powered by high-mass protostars. For our study weselect . ≤ X u ≤ and T u = K in MYSO jets. Hotter(colder) jets likely have higher (lower) ionization fraction.
Internal shocks due to a variable ejection velocity is the mostpopular model to explain the chain of knots in protostellarjets such as HH 111 (Raga et al. 1990). In most of the cases,knots are bright H α and [SII] emitters in the optical bandindicating shock speeds of few tens of km s − , much smallerthan typical proper motions with speeds of hundreds km s − .Since velocity perturbations are of the order of 10% of thejet speed, v j (Reipurth & Bally 2001), internal shocks arenot fast enough to accelerate particles up to very high en-ergies, which is consistent with the fact that most of the knots (and HH objects) are thermal emitters (however, seeRodr´ıguez-Kamenetzky et al. 2017; Osorio et al. 2017, forsome exceptions).Non-thermal emission has been detected in the outerlobes of a handful of protostellar jets, which is more con-sistent with the presence of a strong adiabatic shock at thelocation of the non-thermal emitter where synchrotron emit-ting electrons are accelerated. In the jet termination region,where the jet impacts against the external medium (seeFig. 1), the bow shock moves into the molecular cloud ata speed v bs ∼ v j / ( + / β √ χ ) , where χ ≡ n j / n mc is the jet ( n j )to molecular cloud ( n mc ) density contrast, and β ≡ R h / R j ,being R j and R h the radius of the jet beam and head, respec-tively (e.g. Raga et al. 1998; Hartigan 1989; Chernin et al.1994). The reverse shock (or Mach disc) in the jet moves at v rs = v j − v bs / . In “heavy” jets ( χ > ), the bow shock isfaster than the Mach disc , whereas in “light” jets ( χ < ),the reverse shock is faster than the bow shock. In particular, v rs ∼ v j and v bs ∼ v j √ χ when χ (cid:28) , whereas v bs ∼ v j when χ (cid:29) . As an example, v rs and v bs are ∼
750 and 250 km s − ,respectively, when v j = km s − and χ = . . Propermotions of about − km s − have been observed inprotostellar jets (e.g. Marti et al. 1993; Masqu´e et al. 2012),and therefore we can expect jet velocities comparable to orlarger than these values.The detection of GeV photons from classical novae(Ackermann et al. 2014) indicates that radiative shockswith velocities ≤ km s − can accelerate particlesup to ∼
10 GeV (see e.g. Vurm & Metzger 2018). Thepresent study concentrates, however, on adiabatic (i.e. non-radiative) shocks. The condition for a shock with velocity v sh to be adiabatic is q ≡ l th / R j > (Blondin et al. 1990),where we define R j as the radius of the section of the jet atthe position of the hotspot. The thermal cooling length is l th cm (cid:39) . × (cid:16) n u cm − (cid:17) − (cid:16) v sh − (cid:17) (1)(e.g. Raga et al. 2002), where n u is the density upstream ofthe shock. The condition q > can be rewritten as v sh > v sh , ad ,where v sh , ad km s − (cid:39) (cid:16) n u cm − (cid:17) (cid:18) R j cm (cid:19) . (2)The reverse- to bow-shock cooling parameter ratio is q rs / q bs = l rs / l bs ∼ χ − . indicating that the termination re-gion of light jets ( χ < ) is composed by an adiabatic Machdisc and a radiative bow shock.The jet ion density can be estimated as n ˙ M cm − ≈ (cid:18) ˙ M i − M (cid:12) yr − (cid:19)(cid:16) v j − (cid:17) − (cid:18) R j cm (cid:19) − (3)(e.g. Rodr´ıguez-Kamenetzky et al. 2017). Owing to the largeuncertainties in the different parameters we neglect the con-tribution of helium in Eq. (3). In Fig. 2 we plot n ˙ M for R j = R / and R / , where R is the average linear size of the This situation is very similar to the cloudlet model (Norman &Silk 1979) where a fast cloud moves into a molecular cloud (seee.g. Hartigan 1989) or into a slow jet (Yirak et al. 2012).MNRAS000
10 GeV (see e.g. Vurm & Metzger 2018). Thepresent study concentrates, however, on adiabatic (i.e. non-radiative) shocks. The condition for a shock with velocity v sh to be adiabatic is q ≡ l th / R j > (Blondin et al. 1990),where we define R j as the radius of the section of the jet atthe position of the hotspot. The thermal cooling length is l th cm (cid:39) . × (cid:16) n u cm − (cid:17) − (cid:16) v sh − (cid:17) (1)(e.g. Raga et al. 2002), where n u is the density upstream ofthe shock. The condition q > can be rewritten as v sh > v sh , ad ,where v sh , ad km s − (cid:39) (cid:16) n u cm − (cid:17) (cid:18) R j cm (cid:19) . (2)The reverse- to bow-shock cooling parameter ratio is q rs / q bs = l rs / l bs ∼ χ − . indicating that the termination re-gion of light jets ( χ < ) is composed by an adiabatic Machdisc and a radiative bow shock.The jet ion density can be estimated as n ˙ M cm − ≈ (cid:18) ˙ M i − M (cid:12) yr − (cid:19)(cid:16) v j − (cid:17) − (cid:18) R j cm (cid:19) − (3)(e.g. Rodr´ıguez-Kamenetzky et al. 2017). Owing to the largeuncertainties in the different parameters we neglect the con-tribution of helium in Eq. (3). In Fig. 2 we plot n ˙ M for R j = R / and R / , where R is the average linear size of the This situation is very similar to the cloudlet model (Norman &Silk 1979) where a fast cloud moves into a molecular cloud (seee.g. Hartigan 1989) or into a slow jet (Yirak et al. 2012).MNRAS000 , 1–16 (2015)
A. T. Araudo et al. D E C L I NA T I O N ( J2000 ) RIGHT ASCENSION (J2000)16 58 18.0 17.8 17.6 17.4 17.2 17.0 16.8 16.6 16.4-42 51 505552 0005101520
CentralSource N-1S-1 jetjet jet reverse shock V rs ≈ V j bow shock V bs ≪ V j Molecular cloud contact discontinuity
Figure 1.
Top: The radio image of the source IRAS 16547-4247(G343.1261) at 14.9 GHz (adapted from Rodr´ıguez et al. 2005).Bottom: Sketch of the jet termination region in light protostellarjets. Here, V bs , V rs , V j refer to the bow shock, the reverse shock,and the jet velocity, respectively. Ionized mass loss rate [M sun yr -1 ] D e n s it y [ c m - ] n M (R j = R/2)n M (R j = R/6)n ff (R j = R/2)n ff (R j = R/6)(1) (2) (11) (3,4) (9,10)(5,6)(7,8) Figure 2.
Jet ion density, n ˙ M , for the sources listed in Table 2 asa function of the ionized mass loss rate assuming R j = R / (greentriangles up) and R / (green triangles down). Red squares andcircles indicate the value of n ff when R j = R / and R / , respectively. hotspot also listed in Table 2. The jet mass loss rate of ion-ized matter ˙ M i ∝ v j is assumed to be constant along the jetand then n ˙ M turns out to be independent of the jet velocity.By inserting Eq. (3) in Eq. (2) we find that v sh , ad ∝ R − / .In Table 2 we list n ˙ M and v sh , ad assuming R j = R / . We notethat n ˙ M is a rough estimation of the jet ion density given theuncertainties in the values of ˙ M i and R j . In particular, lightadiabatic jets form a cocoon and therefore the size of thenon-thermal lobe at the jet head is expected to be R h > R j .Krause (2003) found that ( R j / R h ) ∼ . − when χ > . ;this gives an increment in the jet density by a factor of ∼ when R j = R / instead of R / (see Fig. 2). Nevertheless, thederived jet density is likely lower than the typical densityvalues of ∼ cm − in the molecular clouds where mas-sive stars form (Hennebelle & Falgarone 2012) giving χ < .Therefore q rs > q bs and we should rather expect non-thermalsources in light jets ( χ < ), where the reverse shock is fasterthan the bow shock.In order to check if synchrotron emission comes fromthe reverse-shock rather than from the bow-shock down-stream region, we compare the synchrotron ( ε synchr , ν ) andfree-free ( ε ff , ν ) emissivities. The temperature of the plasmaimmediately downstream of the shock with compression ra-tio r = is T d ∼ . × ( v sh / − ) K. The free-free emissivity of the shocked plasma emitted at frequency ν (cid:28) k B T d / h ∼ ( v sh / − ) GHz is ε ff , ν ergcm − s − Hz − ≈ . × − g ( ν , T d ) × (cid:16) n e cm − (cid:17) (cid:16) v sh − (cid:17) − (4)(Lang 1974), where k B and h are the Boltzmann andPlanck constants, respectively. The Gaunt factor is g ( ν , T d ) ∼ . Λ and Λ ≈ . × (cid:16) v sh − (cid:17) (cid:16) ν GHz (cid:17) − (5)when πν (cid:29) ω p and T d > . × K. Here, ω p ∼ . × MNRAS , 1–16 (2015)
FA in jets from MYSOs (cid:112) n e / cm − rad s − is the electron plasma frequency.We note that g ( ν , T ) ∼ . when v sh = km s − and ν ∼ GHz.On the other hand, the synchrotron emissivity of asource located at distance d is ε synchr , ν = π d S ν / V e , where S ν is the synchrotron flux measured at frequency ν and V e ∼ R is the volume of the synchrotron emitter (see Eq. A1). Byimposing ε synchr , ν > ε ff , ν , we find that the density of thermalelectrons downstream of the shock has to be smaller than n ff , where n ff cm − ≈ . × (cid:18) d kpc (cid:19)(cid:18) S ν mJy (cid:19) (cid:18) R j cm (cid:19) − (cid:16) v sh − (cid:17) (6)for the non-thermal emission to dominate at frequency ν (seeHenriksen et al. 1991). In Fig. 2 we plot n ff and in Table 2we list the values of n ff for the case v sh = km s − and R j = R / . We can see that molecular clouds with densities ∼ cm − will not (or marginally) satisfy the condition n mc < n ff for hosting the synchrotron emitter in the shockedregion downstream of the bow shock. On the other hand, n ˙ M (cid:28) n ff for all the sources in the sample.Hereafter, we will consider light jets ( χ < ) for whichthe synchrotron radiation dominates over the free-free emis-sion and identify the reverse shock with jet speed ( v j = v rs ).In order to simplify the notation we will name it v sh and use km s − as a characteristic value . The jet ion densitywill be n ˙ M ≤ n i ≤ n ff with a mean value (cid:104) n i (cid:105) = √ n ˙ M n ff (seeTable 2). We note that (cid:104) n i (cid:105) is similar to n ˙ M when R j = R / .In Fig. 1 we show a sketch of the scenario considered in thisstudy.Setting constraints on n i and v sh is highly important inour study given that the jet kinetic energy density U kin ergcm − (cid:39) . × − X u (cid:16) n i cm − (cid:17)(cid:16) v sh − (cid:17) , (7)and the jet kinetic luminosity L kin ergs − (cid:39) . × X u (cid:16) n i cm − (cid:17)(cid:16) v sh − (cid:17) (cid:18) R j cm (cid:19) , (8)represent the energy budget to accelerate particles at theshock. Another very important parameter in our study isthe magnetic field. Maser emission provides information about the orientationand strength of magnetic fields in the molecular outflowsassociated to the central jet. Goddi et al. (2017) derivedmagnetic field strengths between 100 and 300 mG in theprotostellar jet W3(H O), which likely results from stronggas compression behind shocks associated with the outflowexpansion in the ambient molecular cloud. However, it is This value probably overestimates the effective v rs but, on theother hand, v j deduced from proper motions is underestimatedgiven that it corresponds to the velocity projected on the planeof the sky. difficult to infer the magnetic field orientation and strengthin the central jet itself where the origin of the magnetic fieldis more likely connected to the ejection process from theaccretion disk. Lee et al. (2018) reported on the detection ofSiO line polarization in the HH 211 protostellar jet, with anestimated magnetic field of about 15 mG at ∼
300 AU fromthe central protostar.Theoretical studies of magnetically-accelerated jets re-quire the action of a poloidal component to accelerate theplasma and a toroidal component to confine it (Casse &Keppens 2002). The toroidal component evolves with thedistance z from the star along the jet as B φ ∝ / z . At thejet base, on the edge of the disk, typical values of themagnetic field strength are obtained from the local betaplasma parameter β p = P g / P m , where P g and P m are thegas and magnetic pressures, respectively. Jet launching re-quires β p ≥ (Casse & Keppens 2002) giving an upper limitfor the magnetic field strength at the base of the jet of . ( n disc / cm − )( T disc / K ) G where n disc and T disc arethe gas density and temperature in the disc, respectively.Using the model of Combet & Ferreira (2008) we caninfer typical values of the vertical magnetic field strengthin jet emitting disks B z ∼ . – G at 1 AU from the stardepending on the central mass of the object and the accre-tion rate. We will henceforth consider typical magnetic fieldstrengths at the base of the jet in the range 10 mG–10 Gto be conservative. By using a / z dilution factor, we ob-tain typical magnetic field strengths of B j ∼ µ G – mGat distances of the jet termination shock of ∼ AU or tentimes smaller at a distance of AU. By considering theflux-freezing condition, Hartigan et al. (2007) showed thatthe magnetic field strength of variable MHD jets can be de-scribed by the relation ( B j / µ G ) ∼ ( n i / − ) . givingvalues in agreement with the above estimate. However, atdistances larger than ∼ AU from the protostar, densi-ties and magnetic fields are mostly influenced by shocks andrarefaction waves.The magnetic field topology at the jet termination shockin high-mass protostars is difficult to assess, but it is likelyhelical due to a combination of a toroidal and a poloidalcomponent (Cerqueira et al. 1997; C´ecere et al. 2016), hencethe orientation of the magnetic field with respect to the nor-mal of the termination shock front is likely oblique. On theother hand, numerical studies of non-relativistic and mag-netically driven jets show that kink instabilities destroy theordered helical structure of the magnetic field. Moll (2009)showed that at z (cid:38) – R A , the toroidal field is dissipated andthe poloidal component dominates. Here, R A is the Alfv´enradius, where the poloidal jet and the Alfv´en speeds areequal. For the objects of interest, this length scale is typi-cally less than a few tens of AU (Pelletier & Pudritz 1992).A dominant poloidal field is also supported by laboratoryexperiments (Albertazzi et al. 2014).In this work, acknowledging for the above uncertainties,we adopt B j = µ G as a characteristic value of the jetmagnetic field in the termination region. The Alfv´en speedin the jet is given by v A kms − (cid:39) . (cid:18) B j µ G (cid:19)(cid:16) n i cm − (cid:17) − , (9) MNRAS000
300 AU fromthe central protostar.Theoretical studies of magnetically-accelerated jets re-quire the action of a poloidal component to accelerate theplasma and a toroidal component to confine it (Casse &Keppens 2002). The toroidal component evolves with thedistance z from the star along the jet as B φ ∝ / z . At thejet base, on the edge of the disk, typical values of themagnetic field strength are obtained from the local betaplasma parameter β p = P g / P m , where P g and P m are thegas and magnetic pressures, respectively. Jet launching re-quires β p ≥ (Casse & Keppens 2002) giving an upper limitfor the magnetic field strength at the base of the jet of . ( n disc / cm − )( T disc / K ) G where n disc and T disc arethe gas density and temperature in the disc, respectively.Using the model of Combet & Ferreira (2008) we caninfer typical values of the vertical magnetic field strengthin jet emitting disks B z ∼ . – G at 1 AU from the stardepending on the central mass of the object and the accre-tion rate. We will henceforth consider typical magnetic fieldstrengths at the base of the jet in the range 10 mG–10 Gto be conservative. By using a / z dilution factor, we ob-tain typical magnetic field strengths of B j ∼ µ G – mGat distances of the jet termination shock of ∼ AU or tentimes smaller at a distance of AU. By considering theflux-freezing condition, Hartigan et al. (2007) showed thatthe magnetic field strength of variable MHD jets can be de-scribed by the relation ( B j / µ G ) ∼ ( n i / − ) . givingvalues in agreement with the above estimate. However, atdistances larger than ∼ AU from the protostar, densi-ties and magnetic fields are mostly influenced by shocks andrarefaction waves.The magnetic field topology at the jet termination shockin high-mass protostars is difficult to assess, but it is likelyhelical due to a combination of a toroidal and a poloidalcomponent (Cerqueira et al. 1997; C´ecere et al. 2016), hencethe orientation of the magnetic field with respect to the nor-mal of the termination shock front is likely oblique. On theother hand, numerical studies of non-relativistic and mag-netically driven jets show that kink instabilities destroy theordered helical structure of the magnetic field. Moll (2009)showed that at z (cid:38) – R A , the toroidal field is dissipated andthe poloidal component dominates. Here, R A is the Alfv´enradius, where the poloidal jet and the Alfv´en speeds areequal. For the objects of interest, this length scale is typi-cally less than a few tens of AU (Pelletier & Pudritz 1992).A dominant poloidal field is also supported by laboratoryexperiments (Albertazzi et al. 2014).In this work, acknowledging for the above uncertainties,we adopt B j = µ G as a characteristic value of the jetmagnetic field in the termination region. The Alfv´en speedin the jet is given by v A kms − (cid:39) . (cid:18) B j µ G (cid:19)(cid:16) n i cm − (cid:17) − , (9) MNRAS000 , 1–16 (2015)
A. T. Araudo et al.
Table 2.
Observed and derived parameters of non-thermal lobes associated with MYSOs. From left to right we list the name of the sourceand the non-thermal component in the jet, the observed frequency ν , the measured flux density S ν , the radio spectral index α , the major( θ maj ) and minor axis ( θ min ) deconvolved dimensions of the component and its average linear size R , the synchrotron emissivity ε synchr ,the lower limit on the shock velocity v sh , ad , the density n ˙ M , the upper limit on the lobe density given by Eq. (6), and the average density (cid:104) n i (cid:105) . We define R = ( R maj + R min ) / as the characteristic size of the hotspot, where R min , maj = . × − d θ min , maj and d is the source distancegiven in Table 1.Source ν S ν α θ maj × θ min R ε synchr v sh , ad n ˙ M n ff (cid:104) n i (cid:105) [GHz] [mJy] [arcsec ] [cm] [erg cm − s − Hz − ] [km s − ] [cm − ] [cm − ] [cm − ](1) G263.7434 N 9 0.56 .
45 1 . × .
36 4 . × . × − . × . × . × (2) G263.7759 NW 17 2.16 .
78 1 . × .
34 5 . × . × − . × . × . × (3) G310.1420 A4 9 0.93 .
70 1 . × .
89 3 . × . × − . × . × . × (4) D 9 1.93 .
46 1 . × .
05 4 . × . × − . × . × . × (5) G313.7654 A2 9 0.14 .
68 0 . × .
38 6 . × . × − . × . × . × (6) D 5.5 0.15 .
32 0 . × .
41 3 . × . × − . × . × . × (7) G339.8838 NE 9 2.36 .
39 0 . × .
38 8 . × . × − . × . × . × (8) SW 9 1.50 .
72 5 . × .
85 6 . × . × − . × . × . × (9) G343.1261 N4 17 1.80 .
67 0 . × .
35 1 . × . × − . × . × . × (10) S1 17 4.72 .
45 0 . × .
35 1 . × . × − . × . × . × (11) G114.0835 B 1.5 0.21 .
42 2 . × .
20 1 . × . × − . × . × . × Average values 0.55 . × .
57 3 . × . × − . × . × . × Observational data Section 2.2 and the Alfv´en Mach number is M A = v sh v A (cid:39) (cid:16) v sh − (cid:17)(cid:18) B j µ G (cid:19) − (cid:16) n i cm − (cid:17) . (10)To keep the calculations simple, we derive most of estimatesin the parallel shock configuration (i.e. with a purely poloidalmagnetic field), but we discuss in Sect. 4.1 the impact of themagnetic field obliquity on our results. The energy distribution of particles downstream of a shockis a Maxwellian (thermal) distribution with a power-law(non-thermal) tail at a high energies. The transition en-ergy between both distributions is unknown. By using hy-brid simulations of parallel non-relativistic shocks, Caprioli& Spitkovsky (2014a) have shown that immediately behindthe shock there is a bridge of supra-thermal ions smoothlyconnecting the thermal peak with the power-law, while fardownstream there is a sharp boundary between thermal andnon-thermal protons at E inj , p ∼ E th , where E th = m p v / ∼ . ( v sh / − ) keV, in agreement with Bell (1978b).Park et al. (2015) have shown that the typical energy neededfor electron injection into the shock acceleration is compara-ble with E inj , p . Achieving such injection energy E inj is muchmore difficult for electrons, as their initial momentum is afactor m p / m e smaller. In dense and non-completely ionizedmedium, ionization and Coulomb losses can suppress the ac-celeration of low-energy particles. Ionization and Coulomblosses are almost constant at energies below E = − m p c and E th = × − ( T / K ) m p c , respectively, and thereforethey cannot quench the acceleration at the injection energyif n i cm − < . × κ (cid:18) B j µ G (cid:19)(cid:16) v sh − (cid:17) F , (11) where κ ≥ indicates the departure from the Bohm diffu-sion regime, F = max [ X u / (cid:112) T u , , ( − X u )] , and T u , = T u / K(Drury et al. 1996). For typical parameters in jets fromhigh mass protostars, and X u > . neither ionization norCoulomb losses quench the acceleration at any energy. Weconsider that the spectrum of accelerated particles is givenby Eq. (B1). Synchrotron emission at frequency ν in a magnetic field B s is mostly produced by electrons with energy E e MeV (cid:39) (cid:16) ν GHz (cid:17) / (cid:18) B s mG (cid:19) − / , (12)indicating that synchrotron photons in the range 100 MHz-10 GHz correspond to relativistic electrons with energies inthe range 150 MeV-1.5 GeV, when B s = mG. The energydensity in non-thermal electrons following a power-law dis-tribution N e = K e E − s e is U e , tot = K e f e , where f e is defined inEq. (B2) and s = α + . For the present study we selectnon-thermal lobes with . < α < . and then . < s < . .The normalization factor K e can be determined from themeasured synchrotron flux S ν at a particular frequency ν ,as described in Appendix A. By combining Eqs. (A5) and(B4) we find that U e , tot ergcm − ≈ . × − ξ K ( s ) (cid:18) f e (cid:19)(cid:16) ν GHz (cid:17) s − × (13) (cid:18) ε syn , ν − erg s − cm − Hz − (cid:19)(cid:18) B s mG (cid:19) − ( s + ) , where ξ K ( s ) is given in Eq. (A6) and B s is the magneticfield in the synchrotron emitter, i.e. mostly coming from theshock downstream medium . The magnitude of B s is un-known. It is commonly assumed that B s is in equipartition Depending of the properties of the turbulence around the shock,the magnetic field strength in the synchrotron emitting region hasMNRAS , 1–16 (2015)
FA in jets from MYSOs Magnetic field in the shock downstream region B s [mG] U e , t o t [ e r g c m - ] n i = n Mdot n i =
Non-thermal electron energy density ( U e , tot ) needed toexplain the synchrotron flux at the observed frequency ν in thesources listed in Table 3. The magnetic field strength in the syn-chrotron emitter is B s ≤ B eq . Orange and violet circles indicate thevalues of U e , tot and B s that satisfy the condition U nt = U kin . with non-thermal particles, satisfying the minimum energyrequirement to explain the synchrotron emission. By setting U e , tot = B / ( π ) we find that the magnetic field in equiparti-tion with relativistic electrons is B eq , e mG ≈ ξ eq ( s ) (cid:18) f e (cid:19) s + × (cid:18) ε syn , ν − erg s − cm − Hz − (cid:19) s + (cid:16) ν GHz (cid:17) s − s + , (14)where ξ eq ( s ) (cid:39) ( . ξ K ( s )) / ( s + ) is plotted in Fig. A1. The energy density in protons with energy E p ≥ m i c is U p , tot = aU e , tot , where a = ( m p / m e ) ( s − ) / f p / f e (see Ap-pendix B). The non-thermal energy density is U nt = U p , tot + U e , tot = ( + a ) U e , tot , and then the magnetic field in equipar-tition with non-thermal electrons and protons is B eq = ( + a ) B eq , e . (15)In Table 3 we list the values of s , a , and B eq for all thesources in Table 2 and considering E e , max = . TeV in f e .However, we note that B eq ∝ f / ( s + ) e and therefore the de-pendence on E e , max is negligible. We also point out that B eq is an upper limit for the magnetic field in the synchrotronemitter. In Fig. 3 we plot U e , tot for B s ≤ B eq . Orange andviolet circles corresponds to the case where the total accel-eration efficiency in non-thermal electrons and protons is U nt / U kin = and we fixed v sh = km s − and n i = (cid:104) n i (cid:105) and n ˙ M , respectively, to compute U kin .We define the total proton acceleration efficiency η p , tot = to account for the residence times of the particles downstream andupstream. The exact expression is given by Eq. (18) in Parizotet al. (2006). U p , tot / U kin . These protons can excite various types of insta-bilities at different scales (see e.g. Marcowith et al. 2016,for a review). In particular they drive a current that com-bined with small perturbations present in the magnetic fieldcan excite NR waves (Bell 2004). Previous studies of DSAin protostellar jets consider that particles diffuse back andforth the shock due to (resonant) Alfv´en waves (e.g. Crusius-Watzel 1990; Henriksen et al. 1991; Araudo et al. 2007;Bosch-Ramon et al. 2010; Padovani et al. 2015). However,we will see that, if the jet (unperturbed) magnetic field B j is smaller than a certain value and the driving param-eter (related to the CR current in the upstream medium) ζ = η p v sh / c is high enough, the dominant instability is non-resonant. The acceleration efficiency η p of protons with par-ticular energy E p and energy density U p in Eq. (B4) is definedas η p = U p / U kin = η p , tot / f p , where f p is defined in Eq. (B2)and plotted in Fig. B1. The current j p = n p ev sh of relativistic protons with energy E p and a number density n p = η p U kin / E p in Eq. (B3) candrive MHD turbulence due to the force j p × B added in themomentum equation (Bell 2004, 2005). The parameter η p represents the fraction of the kinetic energy imparted intoprotons with energy E p driving the instabilities. We considera parallel shock in the z -direction with a small perturbationin the plasma. Wave solutions ∝ exp ( ikz − ω t ) of the first or-der MHD perturbed equations lead to the dispersion relation ω = k v + ε k ζ v r g ( σ ( x , s ) − ) , (16)where ε = ± , r g cm (cid:39) . × (cid:18) E p GeV (cid:19)(cid:18) B j µ G (cid:19) − (17)is the proton Larmor radius and σ ( x , s ) = when x ≡ kr g >> (Bell 2014). The solutions of the dispersion relation inEq. (16) are plotted in Fig. 4 for s = , v sh = and km s − and n i = and cm − . We plot Im( ω )and Re( ω ) for ε = (Bell or NR branch, bottom panel)and ε = − (Alfv´en or resonant branch, top panel). Growingmodes correspond to Im ( ω ) > . Im ( ω ) and Re ( ω ) are nor-malized to / ( t (cid:107) ) , where t (cid:107) = r g c / v is the time it takes aparallel shock to travel through the layer within which pro-tons of speed ∼ c are confined and the upstream diffusioncoefficient is κ Bohm = r g c / (Bohm diffusion). More generally,this layer is controlled by the upstream diffusion coefficient κ u which depends on the shock obliquity and we define theadvection time as t adv = κ u / v (see Sect. 4.1).We can see in Fig. 4 (top panel) that resonant wavesare unstable when kr g ≈ and therefore this instability isresonant. On the other hand, the bottom panel shows thatin the NR regime ( kr g > ) Im ( ω ) (cid:29) Re ( ω ) and then the NRmode is almost purely growing. The maximum growth rate Γ max ≡ max ( Im ( ω )) moves to the resonant regime ( kr g ≈ )when v sh decreases. However, even when v sh is as small as300 km s − , the NR mode is still dominant due to the largeion density of protostellar jets. Both modes are present in the MNRAS000
Non-thermal electron energy density ( U e , tot ) needed toexplain the synchrotron flux at the observed frequency ν in thesources listed in Table 3. The magnetic field strength in the syn-chrotron emitter is B s ≤ B eq . Orange and violet circles indicate thevalues of U e , tot and B s that satisfy the condition U nt = U kin . with non-thermal particles, satisfying the minimum energyrequirement to explain the synchrotron emission. By setting U e , tot = B / ( π ) we find that the magnetic field in equiparti-tion with relativistic electrons is B eq , e mG ≈ ξ eq ( s ) (cid:18) f e (cid:19) s + × (cid:18) ε syn , ν − erg s − cm − Hz − (cid:19) s + (cid:16) ν GHz (cid:17) s − s + , (14)where ξ eq ( s ) (cid:39) ( . ξ K ( s )) / ( s + ) is plotted in Fig. A1. The energy density in protons with energy E p ≥ m i c is U p , tot = aU e , tot , where a = ( m p / m e ) ( s − ) / f p / f e (see Ap-pendix B). The non-thermal energy density is U nt = U p , tot + U e , tot = ( + a ) U e , tot , and then the magnetic field in equipar-tition with non-thermal electrons and protons is B eq = ( + a ) B eq , e . (15)In Table 3 we list the values of s , a , and B eq for all thesources in Table 2 and considering E e , max = . TeV in f e .However, we note that B eq ∝ f / ( s + ) e and therefore the de-pendence on E e , max is negligible. We also point out that B eq is an upper limit for the magnetic field in the synchrotronemitter. In Fig. 3 we plot U e , tot for B s ≤ B eq . Orange andviolet circles corresponds to the case where the total accel-eration efficiency in non-thermal electrons and protons is U nt / U kin = and we fixed v sh = km s − and n i = (cid:104) n i (cid:105) and n ˙ M , respectively, to compute U kin .We define the total proton acceleration efficiency η p , tot = to account for the residence times of the particles downstream andupstream. The exact expression is given by Eq. (18) in Parizotet al. (2006). U p , tot / U kin . These protons can excite various types of insta-bilities at different scales (see e.g. Marcowith et al. 2016,for a review). In particular they drive a current that com-bined with small perturbations present in the magnetic fieldcan excite NR waves (Bell 2004). Previous studies of DSAin protostellar jets consider that particles diffuse back andforth the shock due to (resonant) Alfv´en waves (e.g. Crusius-Watzel 1990; Henriksen et al. 1991; Araudo et al. 2007;Bosch-Ramon et al. 2010; Padovani et al. 2015). However,we will see that, if the jet (unperturbed) magnetic field B j is smaller than a certain value and the driving param-eter (related to the CR current in the upstream medium) ζ = η p v sh / c is high enough, the dominant instability is non-resonant. The acceleration efficiency η p of protons with par-ticular energy E p and energy density U p in Eq. (B4) is definedas η p = U p / U kin = η p , tot / f p , where f p is defined in Eq. (B2)and plotted in Fig. B1. The current j p = n p ev sh of relativistic protons with energy E p and a number density n p = η p U kin / E p in Eq. (B3) candrive MHD turbulence due to the force j p × B added in themomentum equation (Bell 2004, 2005). The parameter η p represents the fraction of the kinetic energy imparted intoprotons with energy E p driving the instabilities. We considera parallel shock in the z -direction with a small perturbationin the plasma. Wave solutions ∝ exp ( ikz − ω t ) of the first or-der MHD perturbed equations lead to the dispersion relation ω = k v + ε k ζ v r g ( σ ( x , s ) − ) , (16)where ε = ± , r g cm (cid:39) . × (cid:18) E p GeV (cid:19)(cid:18) B j µ G (cid:19) − (17)is the proton Larmor radius and σ ( x , s ) = when x ≡ kr g >> (Bell 2014). The solutions of the dispersion relation inEq. (16) are plotted in Fig. 4 for s = , v sh = and km s − and n i = and cm − . We plot Im( ω )and Re( ω ) for ε = (Bell or NR branch, bottom panel)and ε = − (Alfv´en or resonant branch, top panel). Growingmodes correspond to Im ( ω ) > . Im ( ω ) and Re ( ω ) are nor-malized to / ( t (cid:107) ) , where t (cid:107) = r g c / v is the time it takes aparallel shock to travel through the layer within which pro-tons of speed ∼ c are confined and the upstream diffusioncoefficient is κ Bohm = r g c / (Bohm diffusion). More generally,this layer is controlled by the upstream diffusion coefficient κ u which depends on the shock obliquity and we define theadvection time as t adv = κ u / v (see Sect. 4.1).We can see in Fig. 4 (top panel) that resonant wavesare unstable when kr g ≈ and therefore this instability isresonant. On the other hand, the bottom panel shows thatin the NR regime ( kr g > ) Im ( ω ) (cid:29) Re ( ω ) and then the NRmode is almost purely growing. The maximum growth rate Γ max ≡ max ( Im ( ω )) moves to the resonant regime ( kr g ≈ )when v sh decreases. However, even when v sh is as small as300 km s − , the NR mode is still dominant due to the largeion density of protostellar jets. Both modes are present in the MNRAS000 , 1–16 (2015)
A. T. Araudo et al.
Table 3.
From left to right we list the source names, the relativistic electrons and protons energy power-law slope s = α + , the totalprotons-to-electrons energy density ratio ( a , for the case E e , max = E p , max = TeV), the equipartition magnetic field ( B eq ), the amplified fieldin the shock downstream region ( B sat , NR ), the total energy density in protons ( U p , tot ) and their acceleration efficiency ( η p , tot ), the protonsmaximum energy ( E p , max ), and the electrons ( K e ) and protons ( K p ) normalization constants. We note that η p , tot and E p , max are computedfor n i = n ˙ M (Eq. 3) and n i = (cid:104) n i (cid:105) = √ n ˙ M n ff (Eqs.3,6).Source s a B eq B sat , NR U p , tot η p , tot E p , max K e K p [mG] [mG] [erg cm − ] [TeV] [erg s − cm − ] [erg s − cm − ] n i = n ˙ M n i = (cid:104) n i (cid:105) n i = n ˙ M i n i = (cid:104) n i (cid:105) (1) G263.74 N 1.90 18.81 4.79 0.55 . × − . × − . × − (2) G263.7759 NW 2.56 5.11 7.70 0.92 . × − . × − . × − (3) G310.1420 A4 2.40 8.90 3.27 0.40 . × − . × − . × − (4) D 1.92 19.24 3.65 0.42 . × − . × − . × − (5) G313.7654 A2 2.36 10.14 2.25 0.28 . × − . × − . × − (6) D 1.64 10.14 2.01 0.22 . × − . × − . × − (7) G339.8838 NE 1.78 15.06 14.6 1.63 . × − . × − . × − (8) SW 2.43 7.80 1.60 0.19 . × − . × − . × − (9) G343.1261 N4 2.34 10.80 10.2 1.22 . × − . × − . × − (10) S1 1.89 63.07 9.93 1.13 . × − . × − . × − (11) G114.0853 B 1.84 17.13 0.82 0.09 . × − . × − . × − Average values 2.13 5.52 0.19 . × − . × − . × − Section 3 Section 5 Section 6 Section 7 plasma, although the dominance of one over the other de-pends on the competition between the Alfv´en (that contains v A ) and the driving (that contains ζ ) terms in Eq. (16). Wenote that the NR mode exists when / r g ≤ k ≤ k max , where k max = Γ max / v A . Therefore, the condition for the growth ofthe NR mode is k max r g > , which can be written as ζ M > ,where ζ M (cid:39) (cid:16) η p . (cid:17)(cid:16) n i cm − (cid:17)(cid:18) B j µ G (cid:19) − (cid:16) v sh − (cid:17) . (18)It is clear from Eq. (18) that the NR mode is more easilytriggered by the fastest shocks. The condition ζ M > canbe written also as U p ( v sh / c ) > B / ( π ) , or equivalently U p ergcm − > . × − (cid:18) B j µ G (cid:19) − (cid:16) v sh − (cid:17) − . (19)In the limit ζ M > , the NR instability dominates. Inthe regime k max , NR ≥ k ≥ / r g , σ ≈ and the dispersion rela-tion in Eq. (16) reduces to ω = k v A − k ζ v r g . (20)The maximum growth rate of the NR instability excited byprotons with energy E p and carrying a current of strength j p is Γ max = . ( j p / c ) (cid:112) π / n i m i , giving Γ max , NR s − (cid:39) . × − (cid:16) η p . (cid:17)(cid:16) v sh − (cid:17) × (cid:16) n i cm − (cid:17) (cid:18) E p GeV (cid:19) − . (21)We note that Γ max , NR is independent of B j , however thewavenumber k max , NR = Γ max , NR / v A = ζ M / r g of the fastestgrowing modes does depend on B j (see Eq. 18). We alsonote that as Γ max , NR and k max , NR decrease with v sh . Reso-nant modes should take over as it is the case in low massYSOs (Padovani et al. 2016). In the limit ζ M A < , the Alfv´en instability dominatesand the dispersion relation in Eq. (16) reduces to ω = k v A .The maximum growth rate of resonant modes occurs for par-allel propagating modes (Achterberg 1981) at k max , R ≈ √ / r g and therefore k max , R r g ≈ , as shown in Fig. 4 (top panel).The growth rate reaches its maximum value at Γ max , R ∼ Γ max , NR (e.g. Pelletier et al. 2006; Amato & Blasi 2009).Alfven instabilities dominates over Bell instabilities when v sh kms − < (cid:16) η p . (cid:17) − (cid:18) B j µ G (cid:19)(cid:16) n i cm − (cid:17) − . (22)This regime was studied by Padovani et al. (2015, 2016) inthe context of low mass protostars, where the jet typical ve-locities are ∼ km s − . The important aspect of destabi-lizing perturbations in resonance with CRs, i.e. when kr g ≈ ,is that the scattering process is more efficient. However, thisdoes not guarantee to have a strong amplification. In thepresent paper we mostly focus on the NR instability. In the perpendicular shock case the wave speed in Eq. (20)now includes a contribution of the compression modes char-acterized by the local sound speed c s (Bell 2005; Marcowithet al. 2018). If c s < v A the non-resonant growth rate isunmodified (Bell 2005; Matthews et al. 2017). If c s > v A the growth rate drops because of thermal effects (see Sec-tion 4.2.1).The necessary condition for the streaming instabilitiesto grow is also modified. In perpendicular shocks, the precur-sor size is shorten because the particle transport is controlledby the perpendicular diffusion κ u = κ ⊥ = κ Bohm κ / ( + κ ) (Forman & Gleeson 1975) where κ = κ (cid:107) / κ Bohm . If angu-lar diffusion proceeds at a smaller rate than Bohm then κ (cid:29) and κ ⊥ (cid:39) κ Bohm / κ . We derive the general expres-sion Γ max , NR t adv ∼ η p ( v sh / v A ) κ b / , where b = and − (with κ > ) for parallel and perpendicular shocks, respectively. Inperpendicular shocks, Γ max , NR t adv > gives η p > κ / M A . High MNRAS , 1–16 (2015)
FA in jets from MYSOs k [1/r g ] I m ( ω ) , R e ( ω ) [ v s h2 / r g c ] v sh =1000 km s -1 , n i =10 cm -3 v sh =1000 km s -1 , n i =10 cm -3 v sh =300 km s -1 , n i =10 I m ( ω ) R e ( ω ) ε = -1 k [1/r g ] I m ( ω ) , R e ( ω ) [ v s h2 / r g c ] v sh =1000 km s -1 , n i =10 cm -3 v sh =1000 km s -1 , n i =10 cm -3 v sh =300 km s -1 , n i =10 I m ( ω ) R e ( ω ) ε = 1 Figure 4.
Dispersion relation ω ( k ) for the resonant (upper andlower panels for kr gm ≤ ) and non-resonant (lower panel only for kr gm ≥ ) modes destabilized by CRs following a power-law energydistribution with s = . Re ( ω ) – dashed lines – and Im ( ω ) – solidlines – are plotted for v sh = and km s − and n i = and cm − . We fix B j = µ G and η p = . . κ values probably quench the development of fast streaminginstabilities in the configuration of a weakly CR modifiedshock. In the mean time, if κ > then η p > . and thedestabilization of streaming instabilities occurs in the regimeof strongly CR modified shocks, a case beyond the scope ofthe present linear analysis. Different effects may reduce the growth rate of the NR in-stability. In this paper we discuss thermal effects and ion-neutral collisions in non-completely ionized jets.
Thermal effects are important when ( v A / v i ) < ( n p / n i )( v sh / v i ) , where the ion speed is v i = k B T u / m p (Zweibel & Everett 2010). This leads to the condition T u K > . × (cid:16) η p . (cid:17) − (cid:18) E p GeV (cid:19) × (cid:16) v sh − (cid:17) − (cid:18) B j µ G (cid:19) (cid:16) n i cm − (cid:17) − (23)in which case an extra term ∝ k T r g ω needs to be addedin the dispersion relation in Eq. (16). The maximum growthrate of the thermally modified Bell instability is Γ max , th (cid:39) ω ci (cid:18) n p n i (cid:19) (cid:18) v sh v i (cid:19) , (24)where ω ci = eB j / ( m p c ) is the ion cyclotron frequency (Revilleet al. 2007). For typical values in protostellar jets, Γ max , th Γ max , NR (cid:39) . (cid:16) η p . (cid:17) − (cid:18) E p GeV (cid:19) (cid:18) T u K (cid:19) − × (cid:16) v sh − (cid:17) − (cid:18) B j µ G (cid:19)(cid:16) n i cm − (cid:17) − . (25)We stress that Γ max , th < Γ max , NR when the temperature iswithin the range where the condition in Eq. (23) is satis-fied. The thermally modified Bell instability is damped when Γ max , th t adv < , i.e. when T u K > . × (cid:16) η p . (cid:17) (cid:18) E p GeV (cid:19) , (26)where we have assumed t adv = t (cid:107) . When thermal effects aretaken into account, they can reduce the growth rate of theNR instability but we do not expect strong damping bythermal effects over the development of the NR streamingmodes, unless the proton acceleration efficiency η p is unrea-sonably small in which case the NR modes are not destabi-lized. In a partially ionized protostellar jet (i.e. X u < ) the frictionarising between charged and neutral particles can quenchthe growth of CR driven instabilities at shock precursors,and therefore DSA is less efficient (e.g. Drury et al. 1996;Reville et al. 2007). When ion-neutral collisions are takeninto account in the MHD equations, the dispersion relationof the CR-driven instability becomes ω + i ω ν in (cid:18) − X u (cid:19) + ω (cid:32) εζ v r gm k ( − σ ) − k v (cid:33) + i ν in (cid:18) X u − X u (cid:19)(cid:32) εζ v r gm k ( − σ ) − k v (cid:33) = , (27)where the ion-neutral collision frequency is given by ν in s − (cid:39) . × − (cid:18) T u K (cid:19) . (cid:16) n n cm − (cid:17) , (28)in the shock upstream region with T u > K and neutraldensity n n = ( / X u − ) n i (Jean et al. 2009). We point outthat in a completely ionized plasma X u = , n n = and then ν in = . In such a case the dispersion relation in Eq. (27) isidentical to Eq. (16). We solve Eq. (27) for ε = and dif-ferent values of X u . As pointed out by Reville et al. (2007),ion-neutral collisions are unable to stabilize the Bell modes MNRAS000
Thermal effects are important when ( v A / v i ) < ( n p / n i )( v sh / v i ) , where the ion speed is v i = k B T u / m p (Zweibel & Everett 2010). This leads to the condition T u K > . × (cid:16) η p . (cid:17) − (cid:18) E p GeV (cid:19) × (cid:16) v sh − (cid:17) − (cid:18) B j µ G (cid:19) (cid:16) n i cm − (cid:17) − (23)in which case an extra term ∝ k T r g ω needs to be addedin the dispersion relation in Eq. (16). The maximum growthrate of the thermally modified Bell instability is Γ max , th (cid:39) ω ci (cid:18) n p n i (cid:19) (cid:18) v sh v i (cid:19) , (24)where ω ci = eB j / ( m p c ) is the ion cyclotron frequency (Revilleet al. 2007). For typical values in protostellar jets, Γ max , th Γ max , NR (cid:39) . (cid:16) η p . (cid:17) − (cid:18) E p GeV (cid:19) (cid:18) T u K (cid:19) − × (cid:16) v sh − (cid:17) − (cid:18) B j µ G (cid:19)(cid:16) n i cm − (cid:17) − . (25)We stress that Γ max , th < Γ max , NR when the temperature iswithin the range where the condition in Eq. (23) is satis-fied. The thermally modified Bell instability is damped when Γ max , th t adv < , i.e. when T u K > . × (cid:16) η p . (cid:17) (cid:18) E p GeV (cid:19) , (26)where we have assumed t adv = t (cid:107) . When thermal effects aretaken into account, they can reduce the growth rate of theNR instability but we do not expect strong damping bythermal effects over the development of the NR streamingmodes, unless the proton acceleration efficiency η p is unrea-sonably small in which case the NR modes are not destabi-lized. In a partially ionized protostellar jet (i.e. X u < ) the frictionarising between charged and neutral particles can quenchthe growth of CR driven instabilities at shock precursors,and therefore DSA is less efficient (e.g. Drury et al. 1996;Reville et al. 2007). When ion-neutral collisions are takeninto account in the MHD equations, the dispersion relationof the CR-driven instability becomes ω + i ω ν in (cid:18) − X u (cid:19) + ω (cid:32) εζ v r gm k ( − σ ) − k v (cid:33) + i ν in (cid:18) X u − X u (cid:19)(cid:32) εζ v r gm k ( − σ ) − k v (cid:33) = , (27)where the ion-neutral collision frequency is given by ν in s − (cid:39) . × − (cid:18) T u K (cid:19) . (cid:16) n n cm − (cid:17) , (28)in the shock upstream region with T u > K and neutraldensity n n = ( / X u − ) n i (Jean et al. 2009). We point outthat in a completely ionized plasma X u = , n n = and then ν in = . In such a case the dispersion relation in Eq. (27) isidentical to Eq. (16). We solve Eq. (27) for ε = and dif-ferent values of X u . As pointed out by Reville et al. (2007),ion-neutral collisions are unable to stabilize the Bell modes MNRAS000 , 1–16 (2015) A. T. Araudo et al.
100 1000 k [1/r g ] I m ( ω ) , R e ( ω ) [ v s h2 / r g c ] x u = 1x u = 0.8x u = 0.6x u = 0.4x u = 0.2 I m ( ω ) v sh = 1000 km s -1 n i = 10 cm -3 T u = 10 Ks = 2
Figure 5.
Growth rate for different values of the ionization fractionupstream of the shock. although the maximum growth rate decreases with X u , aswe illustrate in Fig. 5. This behaviour for Γ max , inc (the max-imum growth rate in the incompletely ionized jet) indi-cates that longer times are required to satisfy the condition Γ max , inc t adv > to excite the NR modes. Magnetic fields in the synchrotron emitter in protostellarjets have strengths B s ≈ . – mG (see Sect. 3). These val-ues are larger than the expected magnetic field in the jettermination region, B j ≈ – µ G, requiring amplificationof the magnetic field, given that compression by the shockdoes not produce strong enough fields. In this study we con-sider the magnetic field amplification by the Bell instabilitydiscussed in the previous section, but in the non-linear case(Bell 2004, 2005).In the linear regime analyzed in Sect. 4, the magneticfield increases exponentially with time until it reaches avalue B ∼ B j , after that the amplification enters in a non-linear regime and the magnetic field growth becomes linearwith time. At a particle energy E p , the amplified magneticfield saturates at B π ∼ η p n i m i v c = U p (cid:16) v sh c (cid:17) (29)(Bell 2004; Pelletier et al. 2006; Zirakashvili & Ptuskin2008). The above expression refers to magnetic field am-plification by protons with energy E p . However, the distri-bution of relativistic protons driving the current spans from ≈ GeV to E p , max ≈ TeV (see Sect. 6) and therefore the (to-tal) saturated magnetic field B sat , NR immediately upstreamof the shock is better estimated by considering the total pro-ton population with energy density U p , tot and accelerationefficiency η p , tot , which gives B sat , NR mG ≈ . (cid:18) U p , tot − ergcm − (cid:19) (cid:16) v sh − (cid:17) ≈ . (cid:16) η p , tot . (cid:17) (cid:16) n i cm − (cid:17) (cid:16) v sh − (cid:17) . (30) We note that B sat , NR does not depends on B j .Once the magnetic field strength reaches the saturationvalue B sat , NR , the Alfv´en velocity becomes larger and η p M decreases. At the same time the Larmor radius of particlesdecreases down to values r g ∼ / k and then the NR Bellinstability becomes subdominant. Caprioli & Spitkovsky(2014b) proposed a refined model for the growth of thenon-resonant modes including a non-linear stage (when en-ergy densities in the turbulent and background magneticfields become similar). They find a maximum magnetic fieldstrength B max / B j ∼ M / √ , where M = v sh / c s is the sonicMach number and c s is the sound speed. In our case wehave moderate M with a typical range 5-50, which makesthe result obtained in Eq. (30) likely optimistic. It shouldbe also noted that Caprioli & Spitkovsky (2014b) con-duct their study in the regime of beta plasma parameter β p = ( c s / v A ) ∼ whereas MYSO jets have rather β p > ifnot (cid:29) because they are relatively hot and dense with aweak magnetic field. The perturbations produced by the NR instability can leadto various effects in the shock downstream region (Giacalone& Jokipii 2007). As NR modes are not normal modes of theplasma, they are expected to be damped rapidly once thesource of excitation disappears. Weibel instabilities wouldalso contribute to produce magnetic field fluctuations up-stream of the shocks. However, Weibel modes are expectedto decay over a few plasma skin depth downstream, so overmuch smaller scales than R j . Conversely, dynamo processescan further amplify the magnetic field produced at the shockfront (Bell 2004; Pelletier et al. 2006; Inoue et al. 2009; Mar-cowith & Casse 2010; van Marle et al. 2018; Tzeferacos et al.2018), even if the upstream region is partially ionized (Xu &Lazarian 2017). This aspect involves a model of non-linearevolution of the plasma and it is beyond the scope of thispaper. For simplicity, we consider that the amplified mag-netic field in the shock upstream region is compressed bythe shock and maintained on a spatial scale ∼ R (the size ofthe synchrotron emitter).In the shock downstream region, the (turbulent)isotropic upstream random magnetic field is compressed bythe shock by a factor r B = (cid:112) ( r + ) / ≈ . when theshock compression ratio is r = (e.g. Parizot et al. 2006;Zirakashvili & Ptuskin 2008) . Therefore, the amplified fielddownstream of the shock is B sat , d = r B B sat , NR . On the otherhand, the magnetic field in the synchrotron emitter B s , i.e.mostly coming from the shock downstream region, is a func-tion of the energy density in non-thermal electrons. In Fig. 6we plot U p , tot = aU e , tot (green-solid lines), where U e , tot is givenby Eq. (13), and U p , tot = ( c / v sh ) B , NR / ( π ) by Eq. (30) (red-dashed line). By assuming that B sat , d = B s we find that B sat , d B eq = (cid:20) r (cid:18) a + a (cid:19)(cid:16) v sh c (cid:17)(cid:21) s + ∼ . ξ sat ( s ) (cid:16) v sh − (cid:17) s + (31) Hereafter we use r = and r B = . . MNRAS , 1–16 (2015) FA in jets from MYSOs Magnetic field in the shock downstream region B s [mG] U p , t o t [ e r g c m - ] n i = n Mdot n i =
Non-thermal protons energy density ( U p , tot = aU e , tot )needed to explain the synchrotron flux at the observed frequency ν in the sources listed in Table 3 (green-solid lines). The magneticfield strength in the synchrotron emitter is B s ≤ B eq , as in Fig. 3.The red dashed line represents the energy density in non-thermalprotons needed to amplify the magnetic field up to the saturationvalue B sat , NR as indicated in Eq. (30). Black squares indicate thevalue of the magnetic field when v sh = km s − . Orange andviolet circles indicate the values of U p , tot and B s that satisfy thecondition U nt = U kin . and U p , tot ergcm − = . × − ξ (cid:16) v sh − (cid:17) − s + s + (cid:18) B eq mG (cid:19) , (32)where ξ sat = . ( s − ) / ( s + ) (see Fig. A1) and we have fixed r B = . and a / ( + a ) ∼ . In Fig. 6 we indicate B s and U p , tot for the cases with v sh = km s − (black squares). In Ta-ble 3 we list B sat , NR = B s / . and U p , tot for v sh = km s − .We list also η p , tot = U p , tot / U kin for n i = n ˙ M and (cid:104) n i (cid:105) . We re-mark that sources G263.7759 NW (2) and G310.1420 A4,D(3,4) have η p , tot ∼ when B s = . B sat , NR , n i = n ˙ M and v sh = km s − . Therefore, we conclude that in these sourcesthe magnetic field either cannot be amplified by Bell insta-bilities or n i > n ˙ M .Contrary to the case of supernova remnants where themagnetic field in the shock downstream region is calculatedfrom the synchrotron cooling length of X-ray emitting elec-trons, in the non-thermal lobes considered in the presentstudy the radio emitting electrons are not in the fast cool-ing regime. We have therefore estimated the magnetic fieldin the non-thermal hotspots by assuming that upstream ofthe shock the magnetic field is amplified by the Bell’s insta-bility up the the saturation regime given by Eq. (30) andby equating B s = . B sat , NR . By assuming that the down-stream magnetic field is B s = B sat , d = . B sat , NR we determine B s , U p , tot , and η p , tot . Once we know the protons accelerationefficiency and the magnetic field, we compute in the follow-ing sections the maximum energy of particles (Sect. 6) andthe γ -ray emission (Sect. 7), respectively. If the non-resonant streaming instability is important andthe magnetic field is amplified, the maximum energy of pro-tons is determined by the amount of protons that escapefrom the shock upstream region (Zirakashvili & Ptuskin2008; Bell et al. 2013). Given that only the most energeticprotons can penetrate far upstream from the shock and am-plify the magnetic field in the shock precursor, the avail-able time to accelerate these particles is ∼ / Γ max , NR ( E p , max ) ,where Γ max , NR ( E p , max ) is the maximum growth rate of NRmodes driven by protons with an energy E p , max . By equat-ing / Γ max , NR ( E p , max ) = R j / v sh and using Eq. (B3) with E p = E p , max we find E p , max m p c = ( − s ) F s < (cid:16) E p , max GeV (cid:17) − F s = (cid:104) ( s − ) m p c F (cid:105) s − s > (33)(Bell et al. 2013; Schure & Bell 2014), where F (cid:39) . (cid:18) U p , tot − ergcm − (cid:19)(cid:18) R j cm (cid:19)(cid:16) n i cm − (cid:17) − . (34)Values of E p , max for the non-thermal lobes in our study arelisted in Table 3 for the case v sh = − and R j = R / .We also considered n i = n ˙ M i and (cid:104) n i (cid:105) . Notice that around theshock the proton mean free path is κ r g . Unless κ (cid:29) , thismean free path at E p , max is always much smaller than R j . Thecondition Γ max , NR ( E p , max ) R j / v sh > was introduced by Bellet al. (2013) for the case of SNR and where they consideredthe size of the spherical shock (instead of R j ). However, Zi-rakashvili & Ptuskin (2008) found a similar expression forplanar non-relativistic shocks. We also note that the samelimit applies when Alfv´en waves dominate given that themaximum growth rate is almost the same (Zirakashvili &Ptuskin 2008).For large temperatures or low magnetic fields (seeEq. (23)), thermal effects are important and therefore E p , max is obtained by equating / Γ max , th = R j / v sh . If thejet is not completely ionized, the maximum energy ofprotons due to escape upstream is reduced by a factor ( Γ max , inc / Γ max , NR ) / ( s − ) . Previous models consider Alfv´en turbulence and Bohm dif-fusion (see e.g. Araudo et al. 2007; Bosch-Ramon et al. 2010;Padovani et al. 2015, 2016). In the present study, electronsdiffuse in the turbulence self generated by protons. Theyare not the main drivers of the turbulence and can be con-sidered as test particles. If the non-resonant streaming in-stability gets into the non-linear phase, and possibly otherinstability can contribute to generate longer wavelength per-turbations, we can expect to have the turbulence coherencelength (cid:96) c limited by the Larmor radius of the protons at themaximum energy E p , max in the amplified field B s (see thediscussion in e.g. Reville et al. 2009; Caprioli & Spitkovsky2014b).If at first we assume that the maximum electrons MNRAS000
Non-thermal protons energy density ( U p , tot = aU e , tot )needed to explain the synchrotron flux at the observed frequency ν in the sources listed in Table 3 (green-solid lines). The magneticfield strength in the synchrotron emitter is B s ≤ B eq , as in Fig. 3.The red dashed line represents the energy density in non-thermalprotons needed to amplify the magnetic field up to the saturationvalue B sat , NR as indicated in Eq. (30). Black squares indicate thevalue of the magnetic field when v sh = km s − . Orange andviolet circles indicate the values of U p , tot and B s that satisfy thecondition U nt = U kin . and U p , tot ergcm − = . × − ξ (cid:16) v sh − (cid:17) − s + s + (cid:18) B eq mG (cid:19) , (32)where ξ sat = . ( s − ) / ( s + ) (see Fig. A1) and we have fixed r B = . and a / ( + a ) ∼ . In Fig. 6 we indicate B s and U p , tot for the cases with v sh = km s − (black squares). In Ta-ble 3 we list B sat , NR = B s / . and U p , tot for v sh = km s − .We list also η p , tot = U p , tot / U kin for n i = n ˙ M and (cid:104) n i (cid:105) . We re-mark that sources G263.7759 NW (2) and G310.1420 A4,D(3,4) have η p , tot ∼ when B s = . B sat , NR , n i = n ˙ M and v sh = km s − . Therefore, we conclude that in these sourcesthe magnetic field either cannot be amplified by Bell insta-bilities or n i > n ˙ M .Contrary to the case of supernova remnants where themagnetic field in the shock downstream region is calculatedfrom the synchrotron cooling length of X-ray emitting elec-trons, in the non-thermal lobes considered in the presentstudy the radio emitting electrons are not in the fast cool-ing regime. We have therefore estimated the magnetic fieldin the non-thermal hotspots by assuming that upstream ofthe shock the magnetic field is amplified by the Bell’s insta-bility up the the saturation regime given by Eq. (30) andby equating B s = . B sat , NR . By assuming that the down-stream magnetic field is B s = B sat , d = . B sat , NR we determine B s , U p , tot , and η p , tot . Once we know the protons accelerationefficiency and the magnetic field, we compute in the follow-ing sections the maximum energy of particles (Sect. 6) andthe γ -ray emission (Sect. 7), respectively. If the non-resonant streaming instability is important andthe magnetic field is amplified, the maximum energy of pro-tons is determined by the amount of protons that escapefrom the shock upstream region (Zirakashvili & Ptuskin2008; Bell et al. 2013). Given that only the most energeticprotons can penetrate far upstream from the shock and am-plify the magnetic field in the shock precursor, the avail-able time to accelerate these particles is ∼ / Γ max , NR ( E p , max ) ,where Γ max , NR ( E p , max ) is the maximum growth rate of NRmodes driven by protons with an energy E p , max . By equat-ing / Γ max , NR ( E p , max ) = R j / v sh and using Eq. (B3) with E p = E p , max we find E p , max m p c = ( − s ) F s < (cid:16) E p , max GeV (cid:17) − F s = (cid:104) ( s − ) m p c F (cid:105) s − s > (33)(Bell et al. 2013; Schure & Bell 2014), where F (cid:39) . (cid:18) U p , tot − ergcm − (cid:19)(cid:18) R j cm (cid:19)(cid:16) n i cm − (cid:17) − . (34)Values of E p , max for the non-thermal lobes in our study arelisted in Table 3 for the case v sh = − and R j = R / .We also considered n i = n ˙ M i and (cid:104) n i (cid:105) . Notice that around theshock the proton mean free path is κ r g . Unless κ (cid:29) , thismean free path at E p , max is always much smaller than R j . Thecondition Γ max , NR ( E p , max ) R j / v sh > was introduced by Bellet al. (2013) for the case of SNR and where they consideredthe size of the spherical shock (instead of R j ). However, Zi-rakashvili & Ptuskin (2008) found a similar expression forplanar non-relativistic shocks. We also note that the samelimit applies when Alfv´en waves dominate given that themaximum growth rate is almost the same (Zirakashvili &Ptuskin 2008).For large temperatures or low magnetic fields (seeEq. (23)), thermal effects are important and therefore E p , max is obtained by equating / Γ max , th = R j / v sh . If thejet is not completely ionized, the maximum energy ofprotons due to escape upstream is reduced by a factor ( Γ max , inc / Γ max , NR ) / ( s − ) . Previous models consider Alfv´en turbulence and Bohm dif-fusion (see e.g. Araudo et al. 2007; Bosch-Ramon et al. 2010;Padovani et al. 2015, 2016). In the present study, electronsdiffuse in the turbulence self generated by protons. Theyare not the main drivers of the turbulence and can be con-sidered as test particles. If the non-resonant streaming in-stability gets into the non-linear phase, and possibly otherinstability can contribute to generate longer wavelength per-turbations, we can expect to have the turbulence coherencelength (cid:96) c limited by the Larmor radius of the protons at themaximum energy E p , max in the amplified field B s (see thediscussion in e.g. Reville et al. 2009; Caprioli & Spitkovsky2014b).If at first we assume that the maximum electrons MNRAS000 , 1–16 (2015) A. T. Araudo et al. energy E e , max exceeds E p , max , then the turbulence expe-rienced by electrons will be in the small-scale turbu-lence regime and the diffusion coefficient would be κ u (cid:39) κ D Bohm ( E p , max , B s )( E e / E p , max ) . This rapid increase with theenergy would limit E e , max to E p , max unless radiative lossesdominate and limit it to smaller values. Then, one shouldexpect to have E e , max ≤ E p , max .By inserting κ u into the acceleration timescale ( t acc ∝ κ u / v ) and the diffusive loss timescale ( R / κ u ) we finda maximum electron energy E e , max , diff to be within a fac-tor of a few of the maximum electron energy fixed by syn-chrotron losses E e , max , rad for the different sources in our sam-ple. Then, one can argue that the statement E e , max ≤ E p , max is reasonably verified accounting for the uncertainties onthe parameters controlled by the microphysics of turbu-lence generation at fast shocks, namely (cid:96) c , κ , B s , and onthe macroscopic jet parameters v sh and R j . Hereafter we as-sume that E e , max = E p , max . ( E e , max , rad < E p , max never happensin our source sample with the derived values of the magneticfield strength B s .) We left to a future work more precise cal-culation of E e , max . TeV electrons and protons can emit gamma-rays by theirinteraction with ambient cold protons through relativis-tic Bremsstrahlung and proton-proton (pp) collisions. In-verse Compton scattering is a mechanism contributing togamma-ray emission as well. In the jet termination region,at ∼ AU from the central protostar, the stellar pho-ton field is very diluted and therefore the Inverse Comptonscattering is not expected to be very relevant. However, westress that thermal photons from the bow-shock downstreamregion can be boosted to the γ -ray domain. The γ -ray fluxby up-scattering of photons from the bow shock will be com-puted in a following paper, together with the diffuse emissionin the molecular cloud.Proton-proton and relativistic Bremsstrahlung coolingtimescales are comparable. However, the number of rela-tivistic protons per unit energy is larger than that of elec-trons, i.e. N p > N e , given that K p > K e (see Appendix B andFig. A1) and we assume that the slope is the same in bothelectrons and protons energy distributions. In particular,the distribution of relativistic protons in the jet hotspot is N p = K p E − s p , where K p = K e ( m p / m e ) ( s − ) / and K e is computedfrom Eq. (A5) by fixing B s = . B sat , NR . In Table 3 we listthe values of K e and K p . As a consequence, the emission inthe γ -ray domain will be dominated by hadrons. In the one-zone model approximation, the specific luminos-ity of gamma rays due to π -decay can be written as E γ L γ , pp = E γ q γ ( E γ ) n lobe V e , (35)where q γ is the emissivity (see Eqs. (16) and (17) in Araudoet al. (2007)). We assume that the volume of the γ -ray emit-ter is the same as the synchrotron emission volume ( V e ) andthat the thermal ion density in the lobe is n lobe = n ˙ M . InFig. 7 we plot the flux F γ , pp = E γ L γ , pp / ( π d ) for all thesources from our sample and the Fermi -LAT sensitivity. We E γ [erg] F γ , pp [ e r g c m - s - ] Fermi (10 yr) (3) (5)(6)(2) (8) (9) (10) (7)(11) (4)(1) n lobe = 4n M dot Figure 7. π -decay for all the sources in our sample for the case of n i = n ˙ M . Black dashed lines indicates the Fermi sensitivity for 10 yrof observation. Numbers indicate the source name (see Table 3). E γ [erg] F γ [ e r g c m - s - ] IRAS 16547 N4 (9)IRAS 16547 S1 (10)IRAS 16484 NE (7)
Fermi (10 yrs)Synchrotron Rel. Bremsstrahlungproton-proton
Figure 8.
Spectral energy distribution from radio to gamma-raysfor the three hot spots: IRAS 16547 N4, IRAS 16547 S1, and IRAS16848 NE. Black squares indicate synchrotron data (Purser et al.2016) and the black dashed line represent the Fermi sensitivityfor 10 years of observation. can see that the predicted emission in G339-8838 NE (7),G34301261 N4 (9), and S1 (10) is detectable by
Fermi -LAT after 10 years of observation (with σ -confidence)when n lobe = n ˙ M . For these three sources we plot in Fig. 8the spectral energy distribution (SED) including also syn-chrotron emission and relativistic Bremsstrahlung and usingthe formulation in Blumenthal & Gould (1970). We assume E e , max = E p , max as it is calculated in Eq. (33).We note that F γ , pp ∝ n lobe and therefore the interactionof relativistic protons with clumps denser than the jet willincrease the γ -ray flux. The termination region of non-relativistic light jets ( n j < n mc ) is expected to be a combination of an adiabatic reverse MNRAS , 1–16 (2015)
FA in jets from MYSOs shock and a radiative bow shock, leading to a large densityratio at the contact discontinuity between both shocks (seee.g. Rodr´ıguez-Kamenetzky et al. 2019b). The density of theplasma downstream of the radiative bow shock is n (cid:48) mc n mc = (cid:16) v bs − (cid:17) (cid:18) T K (cid:19) − (36)(e.g. Blondin et al. 1990) when the plasma is cooled downto a temperature T , making the density contrast at the con-tact discontinuity χ (cid:48) = n j / n (cid:48) mc (cid:28) χ . As a consequence, thecontact discontinuity is unstable to dynamical and thermalinstabilities. A dense layer of density n (cid:48) mc located at distance l th downstream of the bow shock (see Eq. 1) fragments intoseveral clumps (e.g. Calder´on et al. 2020). In this case, theeffective density downstream the reverse shock will increaseup to a value ∼ f v , clump n (cid:48) mc , where f v , clump ≤ is the volumefilling factor of clumps in the lobe with volume ∼ R . How-ever, we draw attention to the fact that the component ofthe magnetic field in the ambient molecular cloud parallelto the bow shock front ( B mc , ⊥ ) limits the compression factorto a maximum value (Blondin et al. 1990) n (cid:48) max n mc ∼ (cid:16) n mc cm − (cid:17) (cid:16) v bs − (cid:17)(cid:18) B mc , ⊥ . (cid:19) − . (37)The magnetic field in molecular clouds is B mc ∼ mG(Crutcher 2012) and we have assumed that B mc , ⊥ ∼ . B mc .In the case of light jets, v bs ∼ v j √ χ = v j (cid:112) n j / n mc (see Sect. 2.2)and therefore Eq. (37) can be written as n (cid:48) max n mc ∼ (cid:16) n j cm − (cid:17) (cid:16) v j − (cid:17)(cid:18) B mc , ⊥ . (cid:19) − . (38)This indicates that significant enhancement in the plasmadensity downstream of the reverse shock is feasible if insta-bilities grow fast enough to fragment the dense shell andform the clumps. However, we note that even when clumpsare not formed, protons accelerated at the jet reverse shockcan diffuse down to the dense shell and radiate there. Kelvin-Helmoltz (KH) and Rayleigh-Taylor (RT) instabili-ties can grow in the contact discontinuity due to the velocityshear and the force exerted by the downstream material ofthe reverse shock on that of the bow shock. If the bow shockis radiative, the formation of a shell much denser than thejet ( n j ) and the molecular cloud ( n mc ) at the contact discon-tinuity makes the working surface unstable even in the caseof light jets. Following the analysis in Blondin et al. (1990),the acceleration of the dense shell with a width w shell can bewritten as a ∼ n j v / n (cid:48) mc w shell . For a characteristic dynamicaltimescale t dyn = R j / v j we find t RT t dyn ∼ . (cid:18) T K (cid:19) − (cid:16) v j − (cid:17) , (39)where t RT ∼ / √ ak is the growing time of RT instabilitiesand we have assumed that k ∼ / w shell and w shell ∼ R j / . The detection of synchrotron emission from jets poweredby high-mass protostars indicates that electrons are in-
Table 4.
Average values of the most relevant parameters in ourstudy (see Tables 2 and 3). From observationsRadio spectral index α = . Synchrotron emissivity ε syn , ν = × − erg cm − s − Hz − Size of the non-thermal lobe R = . × cmAssumed valuesShock speed v sh = km s − Jet magnetic field B j = µ GJet temperature T u = KComputed valuesJet ion density n i = n ˙ M = . × cm − n i = (cid:104) n i (cid:105) = . × cm − Minimum shock speed v sh , ad = km s − Equipartition magnetic field B eq = . mGDownstream magnetic field B s = . mGProtons total energy density U p , tot = . × − erg cm − Protons total acceleration efficiency η p , tot = . Protons maximum energy E p , max = . TeVElectrons normalization constant K e = . × − erg s − cm − Protons normalization constant K p = . × − erg s − cm − situ accelerated. We consider a sample of 11 non-thermallobes in MYSOs as referenced by Purser et al. (2016) andObonyo et al. (2019). From the observational data in theabove mentioned papers, the lobes selected for our studyhave an average radio spectral index α = . , emissivity ε synchr = . × − erg cm − s − Hz − , and linear size R = . × cm (see Table 4). Magnetic fields B s ∼ . − mGare needed to explain the synchrotron radio flux. These largevalues of B s are difficult to explain by a simple compressionin an adiabatic strong shock at the jet termination region,at ∼ . pc from the central protostar. We note that strongcompression of the magnetic field in the molecular cloud bythe radiative bow shock could explain values of B s in therange . − mG , and this scenario will be analysed in afuture study. In the present paper we focus on magnetic fieldamplification by streaming instabilities.CR streaming instabilities can generate magnetic fieldperturbations necessary for the particles to diffuse back andforth the shock. The CR streaming can amplify resonantAlfv´en waves with wavelength of the order of the CR Lar-mor radius r g but also excite the NR (Bell) instability whichproduces perturbations at scales much smaller than r g (Bell2004). Bell’s instability can amplify small perturbations inthe magnetic field up to values much larger than the unper-turbed jet magnetic field B j if ζ M A > where ζ M = (cid:18) U p − ergcm − (cid:19)(cid:18) B j µ G (cid:19) − (cid:16) v sh − (cid:17) (40)is an equivalent expression for Eq. (18). We remark that evenin the case where ζ is small due to the slow speed of theshock (compared e.g. with supernova remnants), the largevalues of jet densities makes M A significantly large to satisfythe condition ζ M (cid:29) over a large range of parameters.However, large ion densities in a non-completely ionized jetwould increase the ion-neutral collisions. Alfv´en waves canbe damped by ion-neutral collisions (e.g. Drury et al. 1996).Conversely, Bell modes are not as heavily damped, althoughthe maximum growth rate decreases (Reville et al. 2007). MNRAS000
Average values of the most relevant parameters in ourstudy (see Tables 2 and 3). From observationsRadio spectral index α = . Synchrotron emissivity ε syn , ν = × − erg cm − s − Hz − Size of the non-thermal lobe R = . × cmAssumed valuesShock speed v sh = km s − Jet magnetic field B j = µ GJet temperature T u = KComputed valuesJet ion density n i = n ˙ M = . × cm − n i = (cid:104) n i (cid:105) = . × cm − Minimum shock speed v sh , ad = km s − Equipartition magnetic field B eq = . mGDownstream magnetic field B s = . mGProtons total energy density U p , tot = . × − erg cm − Protons total acceleration efficiency η p , tot = . Protons maximum energy E p , max = . TeVElectrons normalization constant K e = . × − erg s − cm − Protons normalization constant K p = . × − erg s − cm − situ accelerated. We consider a sample of 11 non-thermallobes in MYSOs as referenced by Purser et al. (2016) andObonyo et al. (2019). From the observational data in theabove mentioned papers, the lobes selected for our studyhave an average radio spectral index α = . , emissivity ε synchr = . × − erg cm − s − Hz − , and linear size R = . × cm (see Table 4). Magnetic fields B s ∼ . − mGare needed to explain the synchrotron radio flux. These largevalues of B s are difficult to explain by a simple compressionin an adiabatic strong shock at the jet termination region,at ∼ . pc from the central protostar. We note that strongcompression of the magnetic field in the molecular cloud bythe radiative bow shock could explain values of B s in therange . − mG , and this scenario will be analysed in afuture study. In the present paper we focus on magnetic fieldamplification by streaming instabilities.CR streaming instabilities can generate magnetic fieldperturbations necessary for the particles to diffuse back andforth the shock. The CR streaming can amplify resonantAlfv´en waves with wavelength of the order of the CR Lar-mor radius r g but also excite the NR (Bell) instability whichproduces perturbations at scales much smaller than r g (Bell2004). Bell’s instability can amplify small perturbations inthe magnetic field up to values much larger than the unper-turbed jet magnetic field B j if ζ M A > where ζ M = (cid:18) U p − ergcm − (cid:19)(cid:18) B j µ G (cid:19) − (cid:16) v sh − (cid:17) (40)is an equivalent expression for Eq. (18). We remark that evenin the case where ζ is small due to the slow speed of theshock (compared e.g. with supernova remnants), the largevalues of jet densities makes M A significantly large to satisfythe condition ζ M (cid:29) over a large range of parameters.However, large ion densities in a non-completely ionized jetwould increase the ion-neutral collisions. Alfv´en waves canbe damped by ion-neutral collisions (e.g. Drury et al. 1996).Conversely, Bell modes are not as heavily damped, althoughthe maximum growth rate decreases (Reville et al. 2007). MNRAS000 , 1–16 (2015) A. T. Araudo et al.
Thermal effects only weakly affect the NR instability growthrate (see Fig. 5). High magnetic field obliquity can also leadto a decrease of the streaming instability growth rates oreven a complete quenching if the precursor length becomestoo short (see Sect. 4.1).By assuming that the large magnetic field in the syn-chrotron emitter ( B s ) is due to amplification through theBell’s instability and fixing v sh = km s − we estimate B s ∼ . B eq and the energy density in non-thermal electrons U e , tot . Then, the energy density in non-thermal protons is U p , tot = aU e , tot . Under the assumption that the jet ionizeddensity is n i = (cid:104) n i (cid:105) , we estimate the proton acceleration ef-ficiency η p , tot to be ∼ . . We stress that this method isdifferent with respect to that used in supernova remnants,where the magnetic field is usually estimated by comparingthe width of X-ray filament profiles with the synchrotroncooling length.By knowing η p , tot and B s we estimate the maximum en-ergy of protons accelerated in the jet reverse shock and the γ -ray emission that they produce. By considering the ampli-fication timescale of the magnetic field in the shock upstreamregion we find E p , max ∼ . TeV. These protons can emit γ rays through their collisions with thermal ions. We computethe γ -ray flux F γ , pp for all the sources in our sample. We notethat F γ , pp ∝ n lobe and the jet ionized density at the termina-tion region is n ˙ M ≤ n i ≤ n ff , where n ff / n ˙ M ∼ (see Fig. 2).We find that having a density n lobe = n ˙ M in the non-thermalemitter is high enough to reach detectable levels of γ rayswith Fermi in IRAS 16547 N4, IRAS 16547 S1, and IRAS16848 NE, as we show in Fig. 7. Although there is no claimof detection of these sources by Fermi, we expect that ourresult will motivate a future study on the Fermi data at thelocation of these sources. In addition, these sources will beperfect targets for a point-source mode. We also note thatmixing due to dynamical instabilities can significantly en-hance the density of targets in the lobes. The very denseshell formed downstream of the radiative bow shock is un-stable and fragmented in clumps with density ∼ − times the density of the ambient medium (i.e. the molecu-lar cloud). By achieving F γ , pp ∼ × − erg cm − s − at0.1 TeV, CTA will be of great importance to measure thecut-off of the spectrum. Also, the spatial resolution of CTAis expected to be better than Fermi. The detection of γ raysfrom protostellar jets will open a new window to study stellarformation, as well as the efficiency of DSA in the high density( n i ∼ − cm − ) and low velocity ( v sh ∼ km s − )regime. In particular, the detection of diffuse gamma-rayemission in molecular clouds where MYSOs are embeddedwill be a piece of evidence of proton acceleration in proto-stellar jets.In a following paper we will perform more detailedcalculations of the multi-wavelength lepto-hadronic spec-tral energy distribution, from radio to γ rays. The diffuse γ -ray emission of particles accelerated in the jet termi-nation shocks and interacting with ions in the molecularcloud where the protostar is embeded will be also mod-eled, as well as the emission of secondary particles. Wewill select the most promising candidates from the sam-ple of sources in the present study, as well as other sourcessuch as IRAS 18162-2048, the powering massive protostar ofthe Herbig Haro objects HH80, HH81 and HH80N (see e.g. Rodr´ıguez-Kamenetzky et al. 2017; Rodr´ıguez-Kamenetzkyet al. 2019b).High sensitivity radio observations in the GHz do-main are very important to model the synchrotron emissionand constrain the energy density in non-thermal particles(Sect. 3). To this purpose, the detection of polarized emis-sion will be crucial to disentangle thermal contamination inthe GHz domain, as well as to investigate the morphologyof the magnetic field near the shock. The Next GenerationVery Large Array (ngVLA) will play a fundamental role onthat (Galv´an-Madrid et al. 2018; Hull et al. 2018). Movingto lower frequencies, the low-energy cutoff of the electronsdistribution is an important parameter related with the effi-ciency of shocks to inject electrons from the thermal pool tothe high energy tail. In this sense, the forthcoming SquareKilometre Array (SKA) will be extremely important to ob-serve the low-energy cutoff (see e.g. Feeney-Johansson et al.2019).Finally we note that protostellar jets have been wellstudied through several laboratory experiments (e.g. Nico-la¨ı et al. 2008; Liang et al. 2018; Suzuki-Vidal et al. 2012).In particular, Suzuki-Vidal et al. (2015) have shown thatlaboratory bow shocks formed by the collision of two coun-terstreaming and supersonic plasma jets is fragmented dueto the rapid growth of thermal instabilities. The formationof collisionless shocks (Li et al. 2019), the development ofplasma instabilities, and the acceleration of particles in laserplasma (Reville et al. 2013) is going to open a fascinatingera of laboratory astrophysics in synergy with high energy-astrophysics. ACKNOWLEDGEMENTS
The authors thank the anonymous referee for his/her com-ments that helped us to improve the quality of our pa-per. The authors thank C. Carrasco-Gonz´alez and S. Cabritfor insightful discussions on protostellar jets. A.T.A. thanksTony Bell and K. Blundell for encouragement and motivat-ing discussions on jets, particle acceleration, and streaminginstabilities. A.T.A. thanks M.V. del Valle and R. Santos-Lima for their help with numerical calculations. A.T.A.thanks the Czech Science Foundation under the grant GAˇCR20-19854S titled “Particle Acceleration Studies in Astro-physical Jets”. M.P. acknowledges funding from the INAFPRIN-SKA 2017 program 1.05.01.88.04, by the Italian Min-istero dell’Istruzione, Universit`a e Ricerca through the grantProgetti Premiali 2012–iALMA (CUP C52I13000140001).and by the project PRIN-INAF-MAIN-STREAM 2017“Pro-toplanetary disks seen through the eyes of new-generationinstruments”. This work has been carried out thanks to thesupport of the OCEVU Labex (ANR-11-LABX-0060) andthe A*MIDEX project (ANR-11- IDEX-0001-02) funded bythe ”Investissements d’Avenir” French government programmanaged by the ANR.
REFERENCES
Achterberg A., 1981, A&A, 98, 161Ackermann M., et al., 2014, Science, 345, 554Albertazzi B., Ciardi A., Nakatsutsumi M., et al. 2014, Science,346, 325 MNRAS , 1–16 (2015)
FA in jets from MYSOs Amato E., Blasi P., 2009, MNRAS, 392, 1591Anglada G., Rodr´ıguez L. F., Carrasco-Gonz´alez C., 2018,A&ARv, 26, 3Araudo A. T., Romero G. E., Bosch-Ramon V., Paredes J. M.,2007, A&A, 476, 1289Araudo A. T., Bell A. R., Blundell K. M., 2015, ApJ, 806, 243Arbutina B., Uroˇsevi´c D., Andjeli´c M. M., Pavlovi´c M. Z., Vukoti´cB., 2012, ApJ, 746, 79Axford W. I., Leer E., Skadron G., 1977, International CosmicRay Conference, 11, 132Bacciotti F., Eisl¨offel J., 1999, A&A, 342, 717Beck R., Krause M., 2005, Astronomische Nachrichten, 326, 414Bell A. R., 1978a, MNRAS, 182, 147Bell A. R., 1978b, MNRAS, 182, 443Bell A. R., 2004, MNRAS, 353, 550Bell A. R., 2005, MNRAS, 358, 181Bell A. R., 2014, Brazilian Journal of Physics, 44, 415Bell A. R., Schure K. M., Reville B., Giacinti G., 2013, MNRAS,431, 415Blandford R. D., Ostriker J. P., 1978, ApJ, 221, L29Blandford R. D., Payne D. G., 1982, MNRAS, 199, 883Blondin J. M., Fryxell B. A., Konigl A., 1990, ApJ, 360, 370Blumenthal G. R., Gould R. J., 1970, Reviews of Modern Physics,42, 237Bosch-Ramon V., Romero G. E., Araudo A. T., Paredes J. M.,2010, A&A, 511, A8Calder´on D., Cuadra J., Schartmann M., Burkert A., Prieto J.,Russell C. M. P., 2020, MNRAS, p. 91Caprioli D., Spitkovsky A., 2014a, ApJ, 783, 91Caprioli D., Spitkovsky A., 2014b, ApJ, 794, 46Carrasco-Gonz´alez C., Rodr´ıguez L. F., Anglada G., Mart´ı J.,Torrelles J. M., Osorio M., 2010, Science, 330, 1209Casse F., Keppens R., 2002, ApJ, 581, 988Ceccarelli C., Dominik C., L´opez-Sepulcre A., Kama M.,Padovani M., Caux E., Caselli P., 2014, ApJ, 790, L1C´ecere M., Vel´azquez P. F., Araudo A. T., De Colle F., EsquivelA., Carrasco-Gonz´alez C., Rodr´ıguez L. F., 2016, ApJ, 816,64Cerqueira A. H., de Gouveia dal Pino E. M., Herant M., 1997,ApJ, 489, L185Chen B., Bastian T. S., Shen C., Gary D. E., Krucker S., GlesenerL., 2015, Science, 350, 1238Chernin L., Masson C., Gouveia dal Pino E. M., Benz W., 1994,ApJ, 426, 204Combet C., Ferreira J., 2008, A&A, 479, 481Crusius-Watzel A. R., 1990, ApJ, 361, L49Crutcher R. M., 2012, ARA&A, 50, 29Drury L., Duffy P., Kirk J. G., 1996, A&A, 309, 1002Favre C., et al., 2018, ApJ, 859, 136Fedriani R., et al., 2019, arXiv e-prints, p. arXiv:1908.05346Feeney-Johansson A., Purser S. J. D., Ray T. P., Eisl¨offel J., HoeftM., Drabent A., Ainsworth R. E., 2019, ApJ, 885, L7Ferri`ere K. M., 2001, Reviews of Modern Physics, 73, 1031Fontani F., Ceccarelli C., Favre C., Caselli P., Neri R., Sims I. R.,et al. 2017, A&A, 605, A57Forman M. A., Gleeson L. J., 1975, Ap&SS, 32, 77Galv´an-Madrid R., Beltr´an M., Ginsburg A., Carrasco-Gonz´alezC., Liu H. B., Rodr´ıguez L. F., Kurtz S., 2018, preprint,( arXiv:1806.10225 )Garay G., Brooks K. J., Mardones D., Norris R. P., 2003, ApJ,587, 739Garcia P. J. V., Ferreira J., Cabrit S., Binette L., 2001, A&A,377, 589Giacalone J., Jokipii J. R., 2007, ApJ, 663, L41Goddi C., Surcis G., Moscadelli L., Imai H., Vlemmings W. H. T.,van Langevelde H. J., Sanna A., 2017, A&A, 597, A43Hartigan P., 1989, ApJ, 339, 987Hartigan P., Edwards S., Pierson R., 2004, ApJ, 609, 261 Hartigan P., Frank A., Varni´ere P., Blackman E. G., 2007, ApJ,661, 910Hennebelle P., Falgarone E., 2012, A&ARv, 20, 55Henriksen R. N., Mirabel I. F., Ptuskin V. S., 1991, A&A, 248,221Hull C. L. H., Carrasco-Gonz´alez C., Williams P. K. G., Gi-rart J. M., Robishaw T., Galv´an-Madrid R., Bourke T., 2018,preprint, ( arXiv:1806.06313 )Inoue T., Yamazaki R., Inutsuka S.-i., 2009, ApJ, 695, 825Jean P., Gillard W., Marcowith A., Ferri`ere K., 2009, A&A, 508,1099Krause M., 2003, A&A, 398, 113Krymskii G. F., 1977, Soviet Physics Doklady, 22, 327Lang K. R., 1974, Astrophysical formulae: A compendium for thephysicist and astrophysicistLee C.-F., Hwang H.-C., Ching T.-C., Hirano N., Lai S.-P., RaoR., Ho P. T. P., 2018, Nature Communications, 9, 4636Li C. K., Tikhonchuk V. T., Moreno Q. e. a., 2019, Phys. Rev.Lett., 123, 055002Liang G. Y., et al., 2018, ApJ, 868, 56Marcowith A., Casse F., 2010, A&A, 515, A90Marcowith A., et al., 2016, Reports on Progress in Physics, 79,046901Marcowith A., Dwarkadas V. V., Renaud M., Tatischeff V., Gi-acinti G., 2018, MNRAS, 479, 4470Marti J., Rodriguez L. F., Reipurth B., 1993, ApJ, 416, 208Masqu´e J. M., Girart J. M., Estalella R., Rodr´ıguez L. F., Beltr´anM. T., 2012, ApJ, 758, L10Matthews J. H., Bell A. R., Blundell K. M., Araudo A. T., 2017,MNRAS, 469, 1849Maurri L., Bacciotti F., Podio L., Eisl¨offel J., Ray T. P., MundtR., Locatelli U., Coffey D., 2014, A&A, 565, A110Moll R., 2009, A&A, 507, 1203Nicola¨ı P., Stenz C., Kasperczuk A. e. a., 2008, Physics of Plas-mas, 15, 082701Norman C., Silk J., 1979, ApJ, 228, 197Obonyo W. O., Lumsden S. L., Hoare M. G., Purser S. J. D.,Kurtz S. E., Johnston K. G., 2019, MNRAS, 486, 3664Osorio M., D´ıaz-Rodr´ıguez A. K., Anglada G., et al. 2017, ApJ,840, 36Padovani M., Hennebelle P., Marcowith A., Ferri`ere K., 2015,A&A, 582, L13Padovani M., Marcowith A., Hennebelle P., Ferri`ere K., 2016,A&A, 590, A8Parizot E., Marcowith A., Ballet J., Gallant Y. A., 2006, A&A,453, 387Park J., Caprioli D., Spitkovsky A., 2015, Physical Review Let-ters, 114, 085003Pelletier G., Pudritz R. E., 1992, ApJ, 394, 117Pelletier G., Lemoine M., Marcowith A., 2006, A&A, 453, 181Petrosian V., 2001, ApJ, 557, 560Podio L., Lefloch B., Ceccarelli C., Codella C., Bachiller R., 2014,A&A, 565, A64Purser S. J. D., et al., 2016, MNRAS, 460, 1039Raga A. C., Binette L., Canto J., Calvet N., 1990, ApJ, 364, 601Raga A. C., Canto J., Cabrit S., 1998, A&A, 332, 714Raga A. C., Noriega-Crespo A., Vel´azquez P. F., 2002, ApJ, 576,L149Reipurth B., Bally J., 2001, ARA&A, 39, 403Reville B., Kirk J. G., Duffy P., O’Sullivan S., 2007, A&A, 475,435Reville B., Kirk J. G., Duffy P., 2009, ApJ, 694, 951Reville B., Bell A. R., Gregori G., 2013, New Journal of Physics,15, 015015Rodr´ıguez-Kamenetzky A., et al., 2017, ApJ, 851, 16Rodr´ıguez-Kamenetzky A., Carrasco-Gonz´alez C., Gonz´alez-Mart´ın O., Araudo A. T., Rodr´ıguez L. F., Vig S., HofnerP., 2019a, MNRAS, 482, 4687MNRAS000
FA in jets from MYSOs Amato E., Blasi P., 2009, MNRAS, 392, 1591Anglada G., Rodr´ıguez L. F., Carrasco-Gonz´alez C., 2018,A&ARv, 26, 3Araudo A. T., Romero G. E., Bosch-Ramon V., Paredes J. M.,2007, A&A, 476, 1289Araudo A. T., Bell A. R., Blundell K. M., 2015, ApJ, 806, 243Arbutina B., Uroˇsevi´c D., Andjeli´c M. M., Pavlovi´c M. Z., Vukoti´cB., 2012, ApJ, 746, 79Axford W. I., Leer E., Skadron G., 1977, International CosmicRay Conference, 11, 132Bacciotti F., Eisl¨offel J., 1999, A&A, 342, 717Beck R., Krause M., 2005, Astronomische Nachrichten, 326, 414Bell A. R., 1978a, MNRAS, 182, 147Bell A. R., 1978b, MNRAS, 182, 443Bell A. R., 2004, MNRAS, 353, 550Bell A. R., 2005, MNRAS, 358, 181Bell A. R., 2014, Brazilian Journal of Physics, 44, 415Bell A. R., Schure K. M., Reville B., Giacinti G., 2013, MNRAS,431, 415Blandford R. D., Ostriker J. P., 1978, ApJ, 221, L29Blandford R. D., Payne D. G., 1982, MNRAS, 199, 883Blondin J. M., Fryxell B. A., Konigl A., 1990, ApJ, 360, 370Blumenthal G. R., Gould R. J., 1970, Reviews of Modern Physics,42, 237Bosch-Ramon V., Romero G. E., Araudo A. T., Paredes J. M.,2010, A&A, 511, A8Calder´on D., Cuadra J., Schartmann M., Burkert A., Prieto J.,Russell C. M. P., 2020, MNRAS, p. 91Caprioli D., Spitkovsky A., 2014a, ApJ, 783, 91Caprioli D., Spitkovsky A., 2014b, ApJ, 794, 46Carrasco-Gonz´alez C., Rodr´ıguez L. F., Anglada G., Mart´ı J.,Torrelles J. M., Osorio M., 2010, Science, 330, 1209Casse F., Keppens R., 2002, ApJ, 581, 988Ceccarelli C., Dominik C., L´opez-Sepulcre A., Kama M.,Padovani M., Caux E., Caselli P., 2014, ApJ, 790, L1C´ecere M., Vel´azquez P. F., Araudo A. T., De Colle F., EsquivelA., Carrasco-Gonz´alez C., Rodr´ıguez L. F., 2016, ApJ, 816,64Cerqueira A. H., de Gouveia dal Pino E. M., Herant M., 1997,ApJ, 489, L185Chen B., Bastian T. S., Shen C., Gary D. E., Krucker S., GlesenerL., 2015, Science, 350, 1238Chernin L., Masson C., Gouveia dal Pino E. M., Benz W., 1994,ApJ, 426, 204Combet C., Ferreira J., 2008, A&A, 479, 481Crusius-Watzel A. R., 1990, ApJ, 361, L49Crutcher R. M., 2012, ARA&A, 50, 29Drury L., Duffy P., Kirk J. G., 1996, A&A, 309, 1002Favre C., et al., 2018, ApJ, 859, 136Fedriani R., et al., 2019, arXiv e-prints, p. arXiv:1908.05346Feeney-Johansson A., Purser S. J. D., Ray T. P., Eisl¨offel J., HoeftM., Drabent A., Ainsworth R. E., 2019, ApJ, 885, L7Ferri`ere K. M., 2001, Reviews of Modern Physics, 73, 1031Fontani F., Ceccarelli C., Favre C., Caselli P., Neri R., Sims I. R.,et al. 2017, A&A, 605, A57Forman M. A., Gleeson L. J., 1975, Ap&SS, 32, 77Galv´an-Madrid R., Beltr´an M., Ginsburg A., Carrasco-Gonz´alezC., Liu H. B., Rodr´ıguez L. F., Kurtz S., 2018, preprint,( arXiv:1806.10225 )Garay G., Brooks K. J., Mardones D., Norris R. P., 2003, ApJ,587, 739Garcia P. J. V., Ferreira J., Cabrit S., Binette L., 2001, A&A,377, 589Giacalone J., Jokipii J. R., 2007, ApJ, 663, L41Goddi C., Surcis G., Moscadelli L., Imai H., Vlemmings W. H. T.,van Langevelde H. J., Sanna A., 2017, A&A, 597, A43Hartigan P., 1989, ApJ, 339, 987Hartigan P., Edwards S., Pierson R., 2004, ApJ, 609, 261 Hartigan P., Frank A., Varni´ere P., Blackman E. G., 2007, ApJ,661, 910Hennebelle P., Falgarone E., 2012, A&ARv, 20, 55Henriksen R. N., Mirabel I. F., Ptuskin V. S., 1991, A&A, 248,221Hull C. L. H., Carrasco-Gonz´alez C., Williams P. K. G., Gi-rart J. M., Robishaw T., Galv´an-Madrid R., Bourke T., 2018,preprint, ( arXiv:1806.06313 )Inoue T., Yamazaki R., Inutsuka S.-i., 2009, ApJ, 695, 825Jean P., Gillard W., Marcowith A., Ferri`ere K., 2009, A&A, 508,1099Krause M., 2003, A&A, 398, 113Krymskii G. F., 1977, Soviet Physics Doklady, 22, 327Lang K. R., 1974, Astrophysical formulae: A compendium for thephysicist and astrophysicistLee C.-F., Hwang H.-C., Ching T.-C., Hirano N., Lai S.-P., RaoR., Ho P. T. P., 2018, Nature Communications, 9, 4636Li C. K., Tikhonchuk V. T., Moreno Q. e. a., 2019, Phys. Rev.Lett., 123, 055002Liang G. Y., et al., 2018, ApJ, 868, 56Marcowith A., Casse F., 2010, A&A, 515, A90Marcowith A., et al., 2016, Reports on Progress in Physics, 79,046901Marcowith A., Dwarkadas V. V., Renaud M., Tatischeff V., Gi-acinti G., 2018, MNRAS, 479, 4470Marti J., Rodriguez L. F., Reipurth B., 1993, ApJ, 416, 208Masqu´e J. M., Girart J. M., Estalella R., Rodr´ıguez L. F., Beltr´anM. T., 2012, ApJ, 758, L10Matthews J. H., Bell A. R., Blundell K. M., Araudo A. T., 2017,MNRAS, 469, 1849Maurri L., Bacciotti F., Podio L., Eisl¨offel J., Ray T. P., MundtR., Locatelli U., Coffey D., 2014, A&A, 565, A110Moll R., 2009, A&A, 507, 1203Nicola¨ı P., Stenz C., Kasperczuk A. e. a., 2008, Physics of Plas-mas, 15, 082701Norman C., Silk J., 1979, ApJ, 228, 197Obonyo W. O., Lumsden S. L., Hoare M. G., Purser S. J. D.,Kurtz S. E., Johnston K. G., 2019, MNRAS, 486, 3664Osorio M., D´ıaz-Rodr´ıguez A. K., Anglada G., et al. 2017, ApJ,840, 36Padovani M., Hennebelle P., Marcowith A., Ferri`ere K., 2015,A&A, 582, L13Padovani M., Marcowith A., Hennebelle P., Ferri`ere K., 2016,A&A, 590, A8Parizot E., Marcowith A., Ballet J., Gallant Y. A., 2006, A&A,453, 387Park J., Caprioli D., Spitkovsky A., 2015, Physical Review Let-ters, 114, 085003Pelletier G., Pudritz R. E., 1992, ApJ, 394, 117Pelletier G., Lemoine M., Marcowith A., 2006, A&A, 453, 181Petrosian V., 2001, ApJ, 557, 560Podio L., Lefloch B., Ceccarelli C., Codella C., Bachiller R., 2014,A&A, 565, A64Purser S. J. D., et al., 2016, MNRAS, 460, 1039Raga A. C., Binette L., Canto J., Calvet N., 1990, ApJ, 364, 601Raga A. C., Canto J., Cabrit S., 1998, A&A, 332, 714Raga A. C., Noriega-Crespo A., Vel´azquez P. F., 2002, ApJ, 576,L149Reipurth B., Bally J., 2001, ARA&A, 39, 403Reville B., Kirk J. G., Duffy P., O’Sullivan S., 2007, A&A, 475,435Reville B., Kirk J. G., Duffy P., 2009, ApJ, 694, 951Reville B., Bell A. R., Gregori G., 2013, New Journal of Physics,15, 015015Rodr´ıguez-Kamenetzky A., et al., 2017, ApJ, 851, 16Rodr´ıguez-Kamenetzky A., Carrasco-Gonz´alez C., Gonz´alez-Mart´ın O., Araudo A. T., Rodr´ıguez L. F., Vig S., HofnerP., 2019a, MNRAS, 482, 4687MNRAS000 , 1–16 (2015) A. T. Araudo et al.
Rodr´ıguez-Kamenetzky A., Carrasco-Gonz´alez C., Gonz´alez-Mart´ın O., Araudo A. T., Rodr´ıguez L. F., Vig S., HofnerP., 2019b, MNRAS, 482, 4687Rodr´ıguez L. F., Garay G., Brooks K. J., Mardones D., 2005,ApJ, 626, 953Schure K. M., Bell A. R., 2014, MNRAS, 437, 2802Shu F., Najita J., Ostriker E., Wilkin F., Ruden S., Lizano S.,1994, ApJ, 429, 781Suzuki-Vidal F., Bocchi M., Lebedev e. a., 2012, Physics of Plas-mas, 19, 022708Suzuki-Vidal F., Lebedev S. V., Ciardi A., Pickworth L. A., Ro-driguez R., Gil J. M., Espinosa G., Hartigan P., 2015, ApJ,815, 96Te¸sileanu O., Mignone A., Massaglia S., Bacciotti F., 2012, ApJ,746, 96Tzeferacos P., Rigby A., Bott A. F. A., Bell A. R., Bingham R.,Casner A., Cattaneo F., et al 2018, Nature Communications,9, 591Vink J., Laming J. M., 2003, ApJ, 584, 758Vurm I., Metzger B. D., 2018, ApJ, 852, 62Xu S., Lazarian A., 2017, ApJ, 850, 126Yirak K., Schroeder E., Frank A., Cunningham A. J., 2012, ApJ,746, 133Zirakashvili V. N., Ptuskin V. S., 2008, ApJ, 678, 939Zweibel E. G., Everett J. E., 2010, ApJ, 709, 1412van Marle A. J., Casse F., Marcowith A., 2018, MNRAS, 473,3394
APPENDIX A: SYNCHROTRON EMISSION
The synchrotron emissivity per unit frequency ν emitted bya source with volume V e located at distance d from Earth is ε syn , ν = π d S ν / V e , giving ε syn , ν ergHz − s − cm − = . × − f Ve (cid:18) d kpc (cid:19) (cid:18) S ν mJy (cid:19)(cid:18) R cm (cid:19) − , (A1)where S ν is the flux at frequency ν . We have defined V e = f Ve R , where V e is the volume filled in with non-thermal elec-trons and f Ve is the volume filling factor. We consider that ε syn , ν = πε ν , where the emission coefficient for synchrotronradiation is ε ν = c K e (cid:104) sin Θ s + (cid:105) B s + s (cid:18) ν c (cid:19) − s , (A2) Θ is the electron pitch angle, c = . × Hz and c = . × − ( s + ) Γ (cid:18) s − (cid:19) Γ (cid:18) s + (cid:19) ergG − sterad − (A3)(Arbutina et al. 2012; Beck & Krause 2005). We assume anisotropic distribution for the orientation of Θ and therefore (cid:104) sin Θ s + (cid:105) = √ π Γ (cid:16) s + (cid:17) Γ (cid:0) s + (cid:1) . (A4)By combining Eqs. (A1)-(A4) we find that K e erg s − cm − ≈ . × − ξ K ( s ) f Ve (cid:18) d kpc (cid:19) (cid:18) S ν mJy (cid:19)(cid:18) R cm (cid:19) − × (cid:16) ν GHz (cid:17) s − (cid:18) B s mG (cid:19) − s + , (A5) Non-thermal particles spectral index s ξ ξ K ξ eq ξ sat K p / K e K p /K e Figure A1.
Functions ξ K (Eq. A6), ξ eq , and ξ sat (left y-axis) and K p / K e (right y-axis) as a function of the spectral index of non-thermal particles. where ξ K ( s ) (cid:39) − . ( s − ) Γ (cid:16) s + (cid:17) Γ (cid:0) s − (cid:1) Γ (cid:0) s + (cid:1) Γ (cid:0) s + (cid:1) ( s + ) . (A6)The function ξ K ( s ) is plotted in Figure A1. APPENDIX B: PARTICLE ENERGY DISTRIBUTION INTHE SHOCK DOWNSTREAM REGION
Non-thermal particles accelerated in a non-relativistic shockfollow a power-law energy distribution N e , p = (cid:40) K e , p m e , p c ( − s ) E − s + e , p E inj < E e , p < m e , p c , K e , p E − s e , p m e , p c < E e , p < E e , p , max , (B1)where e and p stands for electrons and protons, respectively,and E inj ∼ m p v . We note that the condition N e = N p at E < m e c gives K p / K e = ( m p / m e ) ( s − ) / (see Fig. A1). Thespectrum is flatter in the non-relativistic regime, and there-fore non-relativistic particles do not contribute significantlyto the total energy density U e , p and hence most of thenon-thermal number and energy density is due to particleswith E e , p (cid:38) m e , p c . In the case of negligible energy losses,the total energy density stored in non-thermal particles is U e , p , tot = K e , p f e , p ( s ) , where f e , p ∼ (cid:0) − s (cid:1) E − se , p , max if s < (cid:16) E e , p , max m e , p c (cid:17) if s = (cid:0) s − (cid:1) ( m e , p c ) − s if s > . (B2)We note that the ratio a = U p , tot / U e tot can be written as a =( m p / m e ) ( s − ) / f p / f e , and for s > it can be approximated as a ∼ ( m p / m e ) ( − s ) / . In Figure B1 we plot f p , f e , and a for . ≤ s ≤ . and E e , max = E p , max = . , 1, and 10 TeV.The energy and number density of particles with a cer-tain energy E e , p are n e , p = K e , p E − se , p and U e , p = K e , p E − se , p , re-spectively. We are particularly interested in protons, giventhat they are responsible for the CR streaming instability. MNRAS , 1–16 (2015)
FA in jets from MYSOs Non-thermal particles spectral index s f e , p E max = 1 TeVE max = 0.1 TeVE max = 10 TeV 0510152025 a E l ec t r on s ( E e , m i n = m e c ) P r o t o n s ( E p , m i n = m p c ) Figure B1.
Functions f e and f p (left axis) and a (right axis) as afunction of the spectral index of non-thermal particles for differentvalues of E e , max and E p , max . The mumber and energy density of protons with energy E p are n p = U p , tot f p E − sp ∼ (cid:16) U p , tot E p , max (cid:17) ( − s ) (cid:16) E p E p , max (cid:17) − s if s < U p , tot m p c log (cid:16) Ep , max mc (cid:17) (cid:16) E p m p c (cid:17) − if s = (cid:16) U p , tot m p c (cid:17) ( s − ) (cid:16) E p m p c (cid:17) − s if s > (B3)and U p = U p , tot f p E − sp ∼ U p , tot ( − s ) (cid:16) E p E p , max (cid:17) − s if s < U p , tot log (cid:16) Ep , max mpc (cid:17) if s = U p , tot ( s − ) (cid:16) E p m p c (cid:17) − s if s > , (B4)respectively. This paper has been typeset from a TEX/L A TEX file prepared bythe author.MNRAS000