Blazar jets launched with similar energy per baryon, independently of their power
Jesus M. Rueda-Becerril, Amanda O. Harrison, Dimitrios Giannios
MMNRAS , 1–11 (2015) Preprint 7 September 2020 Compiled using MNRAS L A TEX style file v3.0
Baryon loading of blazar jets independent of accretionrate, not so their luminosity
Jes´us M. Rueda-Becerril (cid:63) , Amanda O. Harrison , Dimitrios Giannios Department of Physics, Purdue University, 525 Northwestern Avenue, West Lafayette, IN, 47907, USA
Accepted XXX. Received YYY; in original form ZZZ
ABSTRACT
The most extreme active galactic nuclei (AGN) are the radio active ones whose rel-ativistic jet propagates close to our line of sight. These objects were first classifiedaccording to their emission line features into flat-spectrum radio quasars (FSRQs) andBL Lacertae objects (BL Lacs). More recently, observations revealed a trend betweenthese objects known as the blazar sequence , along with an anti-correlation betweenthe observed power and the frequency of the synchrotron peak. In the present work,we propose a fairly simple idea that could account for the whole blazar population:all jets are launched with similar energy per baryon, independently of their power. Inthe case of FSRQs, the most powerful jets, manage to accelerate to high bulk Lorentzfactors, as observed in the radio. As a result, they have a rather modest magnetizationin the emission region, resulting in magnetic reconnection injecting a steep particleenergy distribution and, consequently, steep emission spectra in the γ -rays. For theweaker jets, namely BL Lacs, the opposite holds true; i.e., the jet does not achievea very high bulk Lorentz factor, leading to more magnetic energy available for non-thermal particle acceleration, and harder emission spectra at frequencies (cid:38) GeV. Inthis scenario, we recover all observable properties of blazars with our simulations, in-cluding the blazar sequence for models with mild baryon loading ( (cid:46) µ (cid:46) ). Thisinterpretation of the blazar population, therefore, tightly constrains the energy perbaryon of blazar jets regardless of their luminosity. Key words: galaxies: BL Lacertae objects: general – magnetic reconnection – accel-eration of particles – accretion, accretion discs – radiation mechanisms: non-thermal– methods: numerical
Blazars are a subclass of radio-loud AGNs with a relativisticjet propagating close to the line of sight of the observer. Theemission from these objects covers all frequencies of the elec-tromagnetic spectrum, producing a double bump structure.The peak of the low frequency bump ranges from infraredto X-ray, whereas the high frequency one peaks in the γ -ray.Blazars have been classified into two subclasses based on theproperties of their emission lines: FSRQs and BL Lacs (Urry& Padovani 1995). Blazar science has greatly advanced, dur-ing the last decade, thanks to dedicated monitoring pro-grams at different wavelengths (e.g., Ghisellini et al. 2010;Blinov et al. 2015; Lister 2016; Jorstad & Marscher 2016;Ackermann et al. 2011; Rani et al. 2017). In part because Fermi -LAT allowed, for the first time, for the systematicstudy of the populations as a whole by following an un- (cid:63)
E-mail: [email protected] precedented number of sources in γ -rays (Ackermann et al.2011; Ajello et al. 2014). Therefore, we can now move beyondthe case-by-case studies and attempt a holistic approach inunderstanding the physical processes involved. One of theclear trends identified by Fermi -LAT is that BL Lac objectsare characterized, on average, by harder spectra than FS-RQs (Ghisellini et al. 2009). As a result, BL Lac objectsare the most extreme TeV emitters (Ajello et al. 2014).BL Lacs are also typically characterized by a synchrotronpeak at higher energies (as high as X-rays). Not surpris-ingly, modeling of the spectrum of blazars requires electronsinjected with much higher energies in BL Lacs than in FS-RQs (Celotti & Ghisellini 2008).The systematic differences of the two blazar classes arenot limited to their γ -ray properties. Radio programs likeMOJAVE have shown that FSRQs are characterized by ex-treme apparent speeds ( β app ∼ tens) in contrast to thoseof BL Lacs ( β app ∼
1) (Kovalev et al. 2009; Lister et al.2011; Homan et al. 2009; Lister et al. 2009). Also, BL Lacs © a r X i v : . [ a s t r o - ph . H E ] S e p J.M. Rueda-Becerril, A.O. Harrison and D. Giannios are likely associated with less powerful jets (FR-I equiva-lent) in contrast to FSRQs (FR II equivalent) (Ghisellini& Celotti 2001; Giommi et al. 2012, 2013; Giustini & Proga2019). It has also been pointed out that the luminosity of thebroad-line region (BLR) may be a distinctive between thetwo kinds of blazars (e.g., Ghisellini & Celotti 2001; Ghis-ellini et al. 2009, 2011; Ghisellini & Tavecchio 2008), as wellwith other intrinsic parameters such as the spin of the blackhole (Meier 2002; Tchekhovskoy et al. 2010; Garofalo 2019).A main parameter in these models is the accretion rate (cid:219) M onto the black hole. Let us introduce here the Eddingtonrate (cid:219) m ≡ (cid:219) M (cid:219) M Edd , (1)where (cid:219) M Edd is the Eddington mass accretion rate (seeSec. 2.1). Therefore, (cid:219) m gives a measure of the accretion rateof the AGN as a fraction of the Eddington rate. In this workwe will use (cid:219) m to differentiate BL Lacs from FSRQs, so thatBL Lac objects would be those blazars with low (cid:219) m , whileFSRQs those with high (cid:219) m .The so called blazar sequence (Padovani 2007) has beenof strong observational and theoretical focus since the firstmultiwavelenght spectral energy distributions (SEDs) of dif-ferent objects were compared (Fossati et al. 1998; Ghiselliniet al. 1998). Evolutionary scenarios have been proposed inthe past decades which connect both kinds of objects interms of accretion efficiency and the jet formation (B¨ottcher& Dermer 2002; Maraschi & Tavecchio 2003; Celotti & Ghis-ellini 2008; Ghisellini et al. 2011). Thanks to Fermi -LAT ob-servations the view of the blazar sequence has evolved andmore sophisticated trends have been proposed since its intro-duction (e.g. Meyer et al. 2011; Finke 2013; Ajello et al. 2014;Rueda-Becerril et al. 2014). Furthermore, recent works havequestioned if those trends correspond to continuum transi-tion between the two kinds of blazars (Padovani et al. 2019;Keenan et al. 2020).On the theoretical front, AGN jets are believed to belaunched magnetically dominated in the vicinity of a rotat-ing black hole (Blandford & Znajek 1977). Magnetohydrody-namic (MHD) simulations of jet acceleration predict that thebulk acceleration of the jet takes place at the expense of itsmagnetization, i.e., while the bulk Lorentz factor Γ of the jetincreases, its magnetization σ (defined as the Poynting fluxto the total energy flux ratio of the jet) decreases (Komis-sarov et al. 2007, 2009; Tchekhovskoy et al. 2008). Accordingto observations, FSRQs appear with a bulk Lorentz factor Γ of a few tens, in contrast to the slower BL Lacs (see, e.g.,Homan et al. (2009)). This means that FSRQs appear to beassociated with more efficiently accelerated jets, leaving alow energy budget per baryon in the emission region. Thisis in contrast to BL Lacs which do not reach as large of abulk Lorenz factor but, as a result, have an emission regionof high magnetization.It is clear from observations that AGN jets may propa-gate as far as several kpc to a few Mpc from the central en-gine. Relativistic hydrodynamic and MHD simulations haveshown that it is highly probable that instabilities may de-velop in relativistic jets (Perucho et al. 2006; L´opez-C´amaraet al. 2013; Matsumoto & Masada 2013; Komissarov et al.2019; Tchekhovskoy & Bromberg 2016). Instabilities maytranslate into dissipation of energy. In particular, if kink in- stabilities develop in the jet, this could translate into a tan-gled magnetic field in the jet (Tchekhovskoy & Bromberg2016; Barniol Duran et al. 2017). This could in turn in-duce the formation of current sheets, allowing to triggermagnetic reconnection. The theory of magnetic reconnec-tion in the context of blazar flares has been explored in thepast several years (Giannios et al. 2009; Nalewajko et al.2011; Sironi et al. 2015; Petropoulou et al. 2016; Christieet al. 2019), showing that it may be the process respon-sible for the non-thermal particle acceleration and radia-tion (Spruit et al. 2001; Giannios & Spruit 2006; Sironi &Spitkovsky 2014; Barniol Duran et al. 2017). In recent years,first-principle particle in cell (PIC) simulations have demon-strated that magnetic reconnection can account for manyof the extreme spectral and temporal properties of blazars(Sironi & Spitkovsky 2014; Sironi et al. 2015; Petropoulouet al. 2016; Christie et al. 2019). Interestingly, these simu-lations have shown that the crucial parameter that controlsthe distribution of accelerated particles is the jet magnetiza-tion σ . Even for a modest increase in σ of the plasma, mag-netic reconnection results in much harder particle distribu-tions, and, as a result, harder emission spectra (Petropoulouet al. 2016, 2019).In this work we will not focus on the details of the struc-tures that form in the current sheet but only on the globalproperties of the emission region. To determine the fractionof magnetic energy that is dissipated in the reconnection re-gion and the resulting particle distributions, we will exploitthe findings of Sironi et al. (2015) and subsequent work.These studies provide specific predictions for the distribu-tion of the accelerated particles as a function of the jet mag-netization σ . The clear trend is that for σ (cid:46) , the resultingparticle spectra are described by a steep power-law distri-bution function γ (cid:48)− p , where the slope p (cid:38) . A soft particleenergy distribution results in low energy peaks for charac-teristic emission bumps as well as softer resulting spectra.This scenario would correspond to FSRQs which, as we havementioned before, have a modest magnetization at the emis-sion region. On the other hand, a strongly magnetized jetsuch as a BL Lac ( σ (cid:38) ) would be characterized by a hardspectrum of accelerated particles with (cid:46) p (cid:46) .The setup of our model is described in Sec. 2, alongwith its most relevant parameters, and a brief descriptionof the numerical code employed. In Sec. 3 we present anddescribe the results obtained out of our simulations. Finally,in Sec. 4 we discuss the model, the results, its implications,and in Sec. 5 we make the final conclusions from this study. According to MHD theory of relativistic jets, a quantitywhich is conserved along magnetic field lines is the total en-ergy flux per unit rest-mass energy flux µ (see Komissarovet al. 2007; Tchekhovskoy et al. 2009), also known as thebaryon loading. For a cold plasma flow: µ = Γ ( + σ ) , (2)where Γ and σ are the flow bulk Lorentz factor and mag-netization, respectively. The magnetization σ is defined asthe ratio between the Poynting flux and the hydrodynamic MNRAS , 1–11 (2015) aryon loading of blazar jets energy flux. σ = B (cid:48) π ρ (cid:48) c , (3)where B (cid:48) and ρ (cid:48) are the magnetic field strength and the massdensity of the plasma .In this section we will describe a simple model fromwhich we are capable of accounting for the blazar sequence by just considering a simple relation between the jet powerand bulk Lorentz factor Γ , where more powerful jets arethe fastest. We assume that both the jet luminosity L j andthe bulk Lorentz factor Γ depend only on the accretion rateparameter (cid:219) m , keeping the baryon loading µ as a free param-eter. This setup strongly constrains/binds the magnetic andkinetic properties of the emission region. We will quantita-tively test this picture and show that the blazar sequencecan be simply understood in a scenario where µ changeslittle among different objects. Let us define the radiative efficiency of the disk η d ≡ L d / (cid:219) Mc (e.g., Davis & Laor 2011), where c is the speed of light, and L d the disc luminosity. From this parameter let us define theEddington mass accretion rate as follows: (cid:219) M Edd ≡ L Edd η d c , (4)where L Edd ≈ . × ( M / M (cid:12) ) erg s − . The jet luminos-ity L j is related to the accretion power by (e.g., Celotti &Ghisellini 2008) L j = η j (cid:219) Mc (5)where η j is the jet production efficiency. From equations (5)and (4) we get that L j = η j η d L Edd (cid:219) m . (6)According to radio observations there seems to be acorrelation between the bulk Lorentz factor of the emissionregion and the jet power (Lister et al. 2009; Homan et al.2009), or (cid:219) m for this effect, according to Eq. (6). Out of theseempirical relation we make the following ansatz: (cid:219) m = (cid:18) ΓΓ (cid:19) s . (7)It is worth noting here that the parameter Γ has no par-ticular physical meaning. This parameter results from theproportionality relation between (cid:219) m and Γ . In other words,the bulk Lorentz factor of the jet is regulated by the Ed-dington ratio. From empirical results we have set s ∼ and Γ ∼ , so (cid:219) m ≈ . × − Γ . (8)The main spectral features derived in this study remain sim-ilar for . (cid:46) s (cid:46) . and when modestly varying Γ . Forfurther details on this regard see Appendix A. Quantities measured in the comoving frame of the fluid will bedenoted with a prime sign (’), unless noted otherwise. Quantitiesmeasured by a cosmologically distant observer will be denotedwith the subscript ‘obs’. Quantities measured in the laboratoryframe will remain unprimed.
According to the standard model of AGNs (Urry & Padovani1995), the material pumped into the jet will often movethrough an external radiation field produced by the BroadLine Region (BLR). The BLR is believed to be reprocessedradiation from the accretion disk (Sikora et al. 1997; Tavec-chio & Ghisellini 2008). The radius, size and geometry of theBLR are still a topic of debate, although it has been thor-oughly studied over the last decades (e.g., Kaspi et al. 2005,2007; Gaskell 2009, and references therein). As mentionedabove, BL Lacs are considered to have a low-Eddington ac-creting black hole, which translates into a faint BLR radia-tion field; opposed to FSRQs whose black hole is consideredto be accreting at higher rates, and therefore a larger densityof reprocessed photons in the BLR.Different models locate the dissipation either below theBLR (Tavecchio & Ghisellini 2008) or outside the BLR(Marscher & Gear 1985). In the present study we will as-sume that energy dissipation takes place within the BLR(e.g., Sikora et al. 1997; Georganopoulos et al. 2005). In ourmodel we will assume that the emission region is immersed inan isotropic and monochromatic radiation field. The energydensity of the external BLR radiation can be parametrizedas follows (Ghisellini & Tavecchio 2008): u BLR = η BLR L d π cR (9)where R BLR (cid:39) L / , cm is the radius of the BLR, η BLR the covering factor, and L d , = L d /( erg s − ) . Fi-nally, we will consider the radiation field in this regionto be monochromatic with frequency ν BLR . In the comov-ing frame of the plasma flow, ν (cid:48) BLR = ν BLR Γ and u (cid:48) BLR = u BLR Γ ( + β / ) , where β ≡ √ − Γ − is the bulk speed ofthe flow in units of the speed of light. Let us consider an electron-proton jet. According to MHDtheory, instabilities in a Poynting flux dominated flow (i.e.,with σ (cid:38) ) lead to the formation of current sheets, wheremagnetic reconnection is triggered (see Eichler 1993; Begel-man 1998; Giannios & Spruit 2006). In the last decade greatprogress has been made on the understanding of relativisticreconnection trough PIC simulations (Sironi & Spitkovsky2014; Sironi et al. 2015; Petropoulou et al. 2016), show-ing that instabilities develop magnetic islands (plasmoids)in which particles accelerate to ultra-high energies due tomagnetic energy dissipation (see Kagan et al. 2015, for areview).The magnetization of a relativistic jet is defined as theratio of the magnetic energy flux to the matter energy flux(e.g., Janiak et al. 2015) σ = L B L kin = L B L j − L B . (10)By solving the above equation for the Poynting flux lumi-nosity we get that L B = σ + σ L j , (11) MNRAS , 1–11 (2015)
J.M. Rueda-Becerril, A.O. Harrison and D. Giannios which in turn we use to calculate the magnetic energy den-sity of the emitting blob in the comoving frame: u (cid:48) B = L B π R (cid:48) c β Γ , (12)where R (cid:48) b is the size of the emission region or blazar zone ,assumed to be comparable to the cross section of the jet.We also assume that, over a dynamical time t dyn ∼ R (cid:48) b / c , afraction f rec of the magnetic energy in the blob is transferredto the electrons in the system in the form of kinetic energy.In other words, from Eq. (12) we get that the luminosity ofthe electrons in the comoving frame of the blob reads: L (cid:48) e = f rec L B β Γ (13) In blazar jets, magnetic reconnection is believed to takeplace far from the central engine, but at sub-parsec scales(e.g., Petropoulou et al. 2016; Christie et al. 2019). We callsuch place the emission region , which we will assume is ata distance R em from the central engine, and to be a spher-ical blob in the comoving frame of the fluid, covering thecross-sectional area of the jet. We will also assume that theemission region is located close to the outer edge of the BLR,e.g., R em = . R BLR (see Padovani et al. 2019). We can esti-mate the radius of the emitting blob, in the comoving frameof the flow, as follows: R (cid:48) b ≈ R em θ j , (14)where θ j ≈ / Γ is the half-opening angle of the conical jet.Let us take now a distant observer whose line of sightmakes an angle θ obs with respect to the direction of motionof the emitting blob. Assuming that the blob emits isotrop-ically (Gould 1979) ν L ν = f ( τ (cid:48) ν (cid:48) ) τ (cid:48) ν (cid:48) D V (cid:48) ν (cid:48) j (cid:48) ν (cid:48) , (15)where τ (cid:48) ν (cid:48) ≡ R (cid:48) b κ ν (cid:48) , j (cid:48) ν (cid:48) and κ (cid:48) ν (cid:48) are the synchrotron emissiv-ity and self-absorption, respectively (Rybicki & Lightman1979), and f ( τ ) ≡ + exp (− τ ) τ − − exp (− τ ) τ , (16)is the optical depth function for a spherical blob (Gould1979; Dermer & Menon 2009). The transformation fromthe comoving frame of the blob to the central enginereference frame is given by the Doppler factor: D ≡[ Γ ( − β cos θ obs )] − . The magnetization of the plasma undergoing magnetic re-connection in the context of blazars has been studied thor-oughly through PIC simulations in recent years (e.g. Sironiet al. 2015, 2016; Petropoulou et al. 2016). As these sim-ulations have shown, the energy distribution of acceleratedelectrons follows a power-law (non-thermal) profile: Q (cid:48) ( γ (cid:48) ) = Q γ (cid:48)− p H [ γ (cid:48) ; γ (cid:48) min , γ (cid:48) max ] (17) where γ (cid:48) is the electrons Lorentz factor in the comovingframe, H [ x ] the Heaviside function, and γ (cid:48) min and γ (cid:48) max arethe minimum and maximum Lorentz factors of the distri-bution of accelerated electrons. The normalization factor Q can be estimated by calculating the power of these electronsfrom Eq. (17), i.e., L (cid:48) e = V (cid:48) Q m e c ∫ γ (cid:48) max γ min d γγ −( p − ) = V (cid:48) Q m e c γ (cid:48) − p min P( γ (cid:48) max / γ (cid:48) min , p − ) , (18)where V (cid:48) = ( / ) π R (cid:48) b3 is the volume of the emission region,and P( a , s ) : = ∫ a d x x − s (19)is the power-law integral function, numerically computedas in Rueda-Becerril (2017). Finally, from equations (13)and (18) we get that Q = f rec L B β Γ V (cid:48) m e c γ (cid:48) − p min P( γ (cid:48) max / γ (cid:48) min , p − ) . (20)In the reconnection region we have that the magneticenergy available per electron in an electron-proton jet is ∼ σ m p c . As we have mentioned, after reconnection takesplace, a fraction of this energy f rec goes into accelerated elec-trons. This fraction is model dependent as has been shownin Sironi et al. (2015). Additionally, the average energy perinjected electron is f rec σ m p c , which means that the averageLorentz factor of the injected electron is (e.g., Petropoulouet al. 2016) (cid:104) γ (cid:105) ∼ f rec σ m p m e . (21) From the average energy and average Lorentz factor of theinjected electrons one finds that γ (cid:48) min = f rec σ m p m e (cid:18) p − p − (cid:19) . (22)The above result holds for p > and γ (cid:48) max (cid:29) γ (cid:48) min . Onthe other hand, if the distribution has a power-law indexof < p < we can make use of the result found in Sironi &Spitkovsky (2014). In that work it was estimated that themean energy per particle cannot exceed ( σ + ) m p c . Fromthis it is deduced that the maximum Lorentz factor is givenby γ (cid:48) max = (cid:18) f rec ( σ + ) m p m e − pp − (cid:19) /( − p ) γ (cid:48) min − p − p . (23)The minimum and maximum Lorentz factors, γ (cid:48) min and γ (cid:48) max , are set separately for high and low magnetized modelsas described in Sec. 2.3.2. Regarding the value of γ (cid:48) max for thecases with p > is estimated by equating the accelerationrate of the electrons to the synchrotron cooling rate (Dermer& Menon 2009), i.e., γ (cid:48) max = (cid:18) π e (cid:15) acc σ T B (cid:48) (cid:19) / , (24)where the parameter (cid:15) acc could be interpreted as the numberof gyrations the electron experience before it is injected intothe system as part of the non-thermal distribution. MNRAS000
J.M. Rueda-Becerril, A.O. Harrison and D. Giannios which in turn we use to calculate the magnetic energy den-sity of the emitting blob in the comoving frame: u (cid:48) B = L B π R (cid:48) c β Γ , (12)where R (cid:48) b is the size of the emission region or blazar zone ,assumed to be comparable to the cross section of the jet.We also assume that, over a dynamical time t dyn ∼ R (cid:48) b / c , afraction f rec of the magnetic energy in the blob is transferredto the electrons in the system in the form of kinetic energy.In other words, from Eq. (12) we get that the luminosity ofthe electrons in the comoving frame of the blob reads: L (cid:48) e = f rec L B β Γ (13) In blazar jets, magnetic reconnection is believed to takeplace far from the central engine, but at sub-parsec scales(e.g., Petropoulou et al. 2016; Christie et al. 2019). We callsuch place the emission region , which we will assume is ata distance R em from the central engine, and to be a spher-ical blob in the comoving frame of the fluid, covering thecross-sectional area of the jet. We will also assume that theemission region is located close to the outer edge of the BLR,e.g., R em = . R BLR (see Padovani et al. 2019). We can esti-mate the radius of the emitting blob, in the comoving frameof the flow, as follows: R (cid:48) b ≈ R em θ j , (14)where θ j ≈ / Γ is the half-opening angle of the conical jet.Let us take now a distant observer whose line of sightmakes an angle θ obs with respect to the direction of motionof the emitting blob. Assuming that the blob emits isotrop-ically (Gould 1979) ν L ν = f ( τ (cid:48) ν (cid:48) ) τ (cid:48) ν (cid:48) D V (cid:48) ν (cid:48) j (cid:48) ν (cid:48) , (15)where τ (cid:48) ν (cid:48) ≡ R (cid:48) b κ ν (cid:48) , j (cid:48) ν (cid:48) and κ (cid:48) ν (cid:48) are the synchrotron emissiv-ity and self-absorption, respectively (Rybicki & Lightman1979), and f ( τ ) ≡ + exp (− τ ) τ − − exp (− τ ) τ , (16)is the optical depth function for a spherical blob (Gould1979; Dermer & Menon 2009). The transformation fromthe comoving frame of the blob to the central enginereference frame is given by the Doppler factor: D ≡[ Γ ( − β cos θ obs )] − . The magnetization of the plasma undergoing magnetic re-connection in the context of blazars has been studied thor-oughly through PIC simulations in recent years (e.g. Sironiet al. 2015, 2016; Petropoulou et al. 2016). As these sim-ulations have shown, the energy distribution of acceleratedelectrons follows a power-law (non-thermal) profile: Q (cid:48) ( γ (cid:48) ) = Q γ (cid:48)− p H [ γ (cid:48) ; γ (cid:48) min , γ (cid:48) max ] (17) where γ (cid:48) is the electrons Lorentz factor in the comovingframe, H [ x ] the Heaviside function, and γ (cid:48) min and γ (cid:48) max arethe minimum and maximum Lorentz factors of the distri-bution of accelerated electrons. The normalization factor Q can be estimated by calculating the power of these electronsfrom Eq. (17), i.e., L (cid:48) e = V (cid:48) Q m e c ∫ γ (cid:48) max γ min d γγ −( p − ) = V (cid:48) Q m e c γ (cid:48) − p min P( γ (cid:48) max / γ (cid:48) min , p − ) , (18)where V (cid:48) = ( / ) π R (cid:48) b3 is the volume of the emission region,and P( a , s ) : = ∫ a d x x − s (19)is the power-law integral function, numerically computedas in Rueda-Becerril (2017). Finally, from equations (13)and (18) we get that Q = f rec L B β Γ V (cid:48) m e c γ (cid:48) − p min P( γ (cid:48) max / γ (cid:48) min , p − ) . (20)In the reconnection region we have that the magneticenergy available per electron in an electron-proton jet is ∼ σ m p c . As we have mentioned, after reconnection takesplace, a fraction of this energy f rec goes into accelerated elec-trons. This fraction is model dependent as has been shownin Sironi et al. (2015). Additionally, the average energy perinjected electron is f rec σ m p c , which means that the averageLorentz factor of the injected electron is (e.g., Petropoulouet al. 2016) (cid:104) γ (cid:105) ∼ f rec σ m p m e . (21) From the average energy and average Lorentz factor of theinjected electrons one finds that γ (cid:48) min = f rec σ m p m e (cid:18) p − p − (cid:19) . (22)The above result holds for p > and γ (cid:48) max (cid:29) γ (cid:48) min . Onthe other hand, if the distribution has a power-law indexof < p < we can make use of the result found in Sironi &Spitkovsky (2014). In that work it was estimated that themean energy per particle cannot exceed ( σ + ) m p c . Fromthis it is deduced that the maximum Lorentz factor is givenby γ (cid:48) max = (cid:18) f rec ( σ + ) m p m e − pp − (cid:19) /( − p ) γ (cid:48) min − p − p . (23)The minimum and maximum Lorentz factors, γ (cid:48) min and γ (cid:48) max , are set separately for high and low magnetized modelsas described in Sec. 2.3.2. Regarding the value of γ (cid:48) max for thecases with p > is estimated by equating the accelerationrate of the electrons to the synchrotron cooling rate (Dermer& Menon 2009), i.e., γ (cid:48) max = (cid:18) π e (cid:15) acc σ T B (cid:48) (cid:19) / , (24)where the parameter (cid:15) acc could be interpreted as the numberof gyrations the electron experience before it is injected intothe system as part of the non-thermal distribution. MNRAS000 , 1–11 (2015) aryon loading of blazar jets Parameter Value θ obs ◦ M bh M (cid:12) η j η d η BLR ν BLR hf rec s Γ µ
50, 70, 90( σ, p ) (1, 3.0), (3, 2.5), (10, 2.2), (15, 1.5), (20, 1.2) Table 1.
Parameters of the present model. See text for a descrip-tion of each of them.
We will consider a one-zone model in which the emissionregion is a spherical blob of radius R (cid:48) b (see Eq. (14)) whichmoves with constant bulk Lorentz factor Γ for a dynami-cal time. We assume that the accelerated particles radiateisotropically in this region. We perform our simulations us-ing the numerical code Paramo (Rueda-Becerril 2020). Thiscode solves the Fokker-Planck equation using a robust im-plicit method (see Chang & Cooper 1970; Park & Petrosian1996), and for each time-step of the simulation the syn-chrotron, synchrotron self-absorption and inverse Comptonemission (both synchrotron self-Compton, SSC, and externalCompton, EIC) are computed with sophisticated numericaltechniques (Mimica & Aloy 2012; Rueda-Becerril et al. 2017;Rueda-Becerril 2017).For the present work we will focus on solving the Fokker-Planck equation without diffusion terms, i.e., ∂ n (cid:48) ( γ (cid:48) , t (cid:48) ) ∂ t (cid:48) + ∂∂γ (cid:48) (cid:2) (cid:219) γ (cid:48) ( γ (cid:48) , t (cid:48) ) n (cid:48) ( γ (cid:48) , t (cid:48) ) (cid:3) = Q ( γ (cid:48) , t (cid:48) )− n (cid:48) ( γ (cid:48) , t (cid:48) ) t esc , (25)where n (cid:48) is the electrons energy distribution (EED) in theflow comoving frame, Q is a source term (see Eq. (17)), and t esc = t dyn is the electrons the average escape time. The elec-trons radiative energy losses are accounted for with the co-efficient (Rybicki & Lightman 1979): − (cid:219) γ (cid:48) = c σ T m e c β (cid:48) γ (cid:48) ( u (cid:48) B + u (cid:48) BLR ) , (26)where β (cid:48) e is the speed of the electron, in units of c , in thecomoving frame. In this section, we describe the results obtained from oursimulations for different values of the parameters of themodel. In our model, described in the previous section, weaccomplished to reduce parameter space. In Tab. 1 we sum-marize the parameters and values employed in the presentwork. As discussed below, the value of most of these param-eters is constrained by either observations or theory.The accretion disk and jet are parametrized by the blackhole mass M bh , the radiative efficiency of the accretion disk η d , and the jet production efficiency η j . The values for theseparameters were motivated by observations, theory and sim-ulations. For instance, measurements of Bian & Zhao (2003) and Davis & Laor (2011) agree that, for quasars, η d ∼ . .Meanwhile, simulations by Tchekhovskoy et al. (2012) showthat η j may vary between 0.3 and 0.9, depending on thespin of the black hole. Nevertheless, we studied the effectof changing η j in our simulations. We observed that thisparameter controls the luminosity of the synchrotron peakand, to a lesser extent, the luminosity of the EIC peak. With η j = . the bumps increase slightly, while for η j = . theobjects are less luminosity, keeping qualitatively the samespectral features. The BLR is modeled by the covering fac-tor η BLR , and the frequency of the external radiation field ν BLR . Following the formulation by Ghisellini & Tavecchio(2008), we set η BLR = . , while h ν BLR = eV, which is anarbitrary value chosen between the characteristic hydrogenionization frequencies H α and Ly- α .The magnetic reconnection dissipation factor f rec hasbeen set to 0.15, following Petropoulou et al. (2019). Thepower-law index p of the injected particles, Eq. (17), hasbeen estimated by Sironi et al. (2015), and more recently byPetropoulou et al. (2019). Those works report that highlymagnetized flows ( σ (cid:38) ) accelerate electrons with power-law indices in the range (cid:46) p (cid:46) , while mildly magnetizedmodels ( σ (cid:46) ) show electrons distributions with p (cid:38) .A highly magnetized jet will be associated with BL Lac ob-jects, whereas the mildly magnetized to FSRQ jets. Finally,the extrema of the injected particle distribution, γ (cid:48) min and γ (cid:48) max , for FSRQs are given by equations (22) and (24), re-spectively, assuming that the most energetic electrons un-dergo ≈ gyrations before they are injected into the sys-tem (the exact choice for this parameter does not have animportant effect on the results as long as γ (cid:48) max (cid:29) γ (cid:48) min ).Meanwhile, we know that the synchrotron peak of BL Lac-like simulations is given by γ (cid:48) max , which is calculated usingEq. (23). If we take a small value of γ (cid:48) min , the synchrotronpeak will shift to larger frequencies, some of them unrealis-tic, and not shown here. Using the synchrotron peak fromradio observations as a guide, it is therefore possible to con-strain γ (cid:48) min to a reasonable value of ∼ for the injecteddistribution of particles in BL Lacs-like models.Radio observations have shown that the bulk Lorentzfactor of blazar jets ranges from a few to no more than 40(e.g., Lister 2016, found that sources with Γ > are ex-tremely rare). Assuming that blazar jets are ejected withsimilar baryon loading, a jet with µ of a few would implythat the jet will not be able to reach high magnetization andits Lorentz factor will be of order unity. Therefore, we esti-mate that a jet consistent with observations and simulationsshould have a baryon loading µ > .As we have mentioned in the previous section, our modelresides on the hypothesis that all blazars are launched withsimilar baryon loading. In Fig. 1 we show the sequence ofSEDs for three different values of µ . The solid, dashed, dot-dashed, dot-dot-dashed and dotted lines correspond to mag-netization σ = , , , , , respectively. FSRQs are thebrightest of all blazars in all frequencies, their inverse Comp-ton (IC) component tends to be louder than the synchrotronone, and ν syn falls in the infra-red. These features also appearin our simulations with the lowest magnetization, which weassumed as FSRQ-like. On the other hand, the main featuresobserved in SEDs of BL Lac objects are a quieter IC compo-nent, ν syn in the UV–X-rays, and a harder spectral index inthe γ -rays. We find that this is also the case for the highly MNRAS , 1–11 (2015)
J.M. Rueda-Becerril, A.O. Harrison and D. Giannios
12 14 16 18 20 22 24 26 log ν [Hz] l og ν L ν [ e r g s − ] µ = σ = σ = σ = σ = σ =
12 14 16 18 20 22 24 26 log ν [Hz] µ =
70 12 14 16 18 20 22 24 26 log ν [Hz] µ = Figure 1.
Sequence of blazar SEDs for varying model parameters. From left to right, each panel shows the averaged SEDs for differentbaryon loading µ = , , , respectively. Solid, dashed, dot-dashed, dot-dot-dashed and dotted lines correspond to those simulationswith σ = , , , , , respectively. The SEDs were averaged over 1 day after particles start being injected in the emission region.
43 44 45 46 47 48 49 50 log L γ [erg s − ] . . . . . α γ FSRQsBL Lacs − − − − − log L γ / L Edd − − − − − − l og L B L R / L E dd µ = µ = µ = σ = σ = σ = σ = σ = Figure 2. γ -ray spectral index α γ , γ -ray luminosity L γ , and BLR luminosity L BLR . Observational data from Ghisellini et al. (2011) isshown as dark and light gray crosses. Squares, circles and triangles depict the models with baryon loading µ = , , , respectively.Blue, orange, green, red and purple colors show the simulation results with magnetization σ = , , , and 20, respectively. Left panel : γ -ray energy spectral index α γ as a function of the γ -ray luminosity L γ . Right panel : Luminosity of the BLR L BLR as a function of L γ ,both in units of the Eddington luminosity L Edd . MNRAS000
43 44 45 46 47 48 49 50 log L γ [erg s − ] . . . . . α γ FSRQsBL Lacs − − − − − log L γ / L Edd − − − − − − l og L B L R / L E dd µ = µ = µ = σ = σ = σ = σ = σ = Figure 2. γ -ray spectral index α γ , γ -ray luminosity L γ , and BLR luminosity L BLR . Observational data from Ghisellini et al. (2011) isshown as dark and light gray crosses. Squares, circles and triangles depict the models with baryon loading µ = , , , respectively.Blue, orange, green, red and purple colors show the simulation results with magnetization σ = , , , and 20, respectively. Left panel : γ -ray energy spectral index α γ as a function of the γ -ray luminosity L γ . Right panel : Luminosity of the BLR L BLR as a function of L γ ,both in units of the Eddington luminosity L Edd . MNRAS000 , 1–11 (2015) aryon loading of blazar jets magnetized cases. Finally, by contrasting all frames in Fig. 1we can see the blazar sequence trend (cf. Fossati et al. 1998,Fig. 12) is favored for µ > . The jets with larger baryonloading correspond to those sources with larger bulk Lorentzfactor. From Eq. (7), these sources correspond to the mostefficient accretion disks which in turn correspond to thosewith most powerful jets (see Eq. (6)). This effect is moreevident for the highly magnetized cases, whose luminosityincreases for almost two orders of matnitude.In Fig. 2 we present our simulations with σ = , , , , in blue, orange, green, red and purple points,respectively. Those simulations with baryon loading µ = , , are depicted in squares, circles and triangles, re-spectively. Observation data from Ghisellini et al. (2011) isseen in light and dark gray crosses. On the left panel, weshow the spectral index α γ as a function of the bolometricluminosity L γ in the band 0.1-10 GeV (cf. Fig. 1 in Ghiselliniet al. 2011). Observations here are presented in the 1LACcatalogue and range from γ -ray luminosity of 0.1 to 10 GeVand have known redshift. Ghisellini et al. (2011) note thatthe division between BL Lacs and FSRQs is usually around ergs s − , interpreted as a shift from an efficient ac-cretion disk to a relatively inefficient disk. Our simulationsshow a similar trend: efficiently accreting sources with pow-erful jets (FSRQ-like) inhabit the area with L γ (cid:38) erg s − and softer γ -rays spectral index. Mild and highly magnetizedsimulations fall in the area of BL Lac objects with low γ -raysluminosity.On the right panel of Fig. 2 we show the BLR lumi-nosity, L BLR , as a function of L γ , both in units of the Ed-dington luminosity L Edd , together with observational datapoints from Fig. 3 in Ghisellini et al. (2011). According toGhisellini et al. (2011), those sources with a stronger emis-sion lines, i.e., showing a more luminous BLR, appear louderin the γ -ray band. The latter being FSQRs. In our simula-tions, the corresponding ones with a more luminous BLRare those with larger Γ . Our model states that these objectshave larger Eddington ratio (see Eq. (7)), i.e., that wouldcorrespond to highly efficient accretion objects.In the same manner, in Fig. 3 we present our simulationswith σ = , , , , in blue, orange, green, red and purplepoints, respectively. Baryon loadings µ = , , are shownin squares, circles and triangles, respectively. Light and darkgray crosses correspond to BL Lacs and FSRQs sources, re-spectively. On the left panel we show the apparent velocity of our synthetic objects. The observational data correspondto the data in the MOJAVE survey, reported in Lister et al.(2019). A translucent gray arrow draws the trend of incre-ment of the jet luminosity. In this plot we can appreciatehow the synchrotron peak ν syn of our simulations is simi-lar for each magnetization. The apparent velocity is bulkLorentz factor dependent due to relativistic boosting. Thiseffect is clear for those objects with larger Γ (blue and or-ange points), which correspond to those simulations withmore powerful jets. Our simulations with powerful jets con- The apparent velocity is calculated according to the followingexpression: v app = v sin θ obs − vc sin θ obs , where v = β c is the bulk speed of the flow. cur with FSRQs as assumed. This is the case as well withhighly magnetized objects. These objects represent the lesspowerful jets, and fall well in the region of BL Lacs.In the leptonic model of blazars, the Compton domi-nance is defined as the ratio of luminosities between the ICand the synchrotron components of their SED. On the rightpanel of Fig. 3 we contrast the Compton dominance and ν syn of our synthetic sources with the observational data reportedin Finke (2013), depicted as gray crosses. These sources arepresented in the 2LAC clean sample where all had knownredshift and could clearly be classified. In that same work,sources with unknown redshift were also taken into account,finding that the relation between Compton dominance andsynchrotron peak frequency have a physical origin ratherthan it being a redshift selection effect. Regarding our sim-ulations, we can observe that all our simulations fall withinthe observational points. The gray transparent arrow showsthe trend of increment of the jet luminosity. Our simulationsshow that, keeping µ constant, changing the magnetizationwill give the transition from synchrotron-dominant (highlymagnetized) to Compton-dominant and γ -ray loud sources. According to our model, BL Lacs are those blazars withlargest magnetization ( σ (cid:38) ) at the dissipation region.FSRQs, on the other hand, are those with powerful jets butwith low/mild magnetization ( σ (cid:46) ) at the blazar zone.In Fig. 4, it is shown the relation between the main param-eters of our study: the magnetization σ , the bulk Lorentzfactor Γ , and the baryon loading µ , as prescribed by the µσ Γ relation (2). In color gradient we have included thecorresponding jet luminosity L j , in units of the Eddingtonluminosity L Edd (see Eq. (6)). The µσ Γ relation constrainsthese objects to have a mild baryon loading since our modelstands on the assumption that blazars are launched withsimilar baryon loading. Jets launched with µ > wouldgive values of Γ way beyond those inferred from radio obser-vations, for those cases with low magnetization. If blazarswere launched with too low baryon loading, the resulting Γ ∼ would contradict both simulations and observations.These scenarios have been discarded from our analysis. BLLac objects, as blazars with low jet luminosity, fall in theblue–gray region with (cid:46) − L Edd . According to our results(described in Sec. 3), this same region corresponds to oursimulations with high magnetization. FSRQs, the most pow-erful of observed blazars, fall in the the gray–red region. Jetswith super-Eddington power, i.e., those cases with (cid:219) m ∼ , be-long to the orange region in upper-left corner (see App. A).Mildly magnetized blazars, e.g., σ = , develop aparticular behavior. These models have an Eddington rate L j / L Edd ∼ . , synchrotron peak ν syn (cid:38) Hz, like someFSRQs. However, their IC component is less ( µ = ) orsimilar ( µ = ) in luminosity to the synchrotron compo-nent, and the γ -ray spectral index is harder; characteris-tics of BL Lac objects. According to Padovani et al. (2019),the object TXS 0506+056, a “masquerading” BL Lac ob-ject, shows properties like (cid:46) L γ /( erg s − ) (cid:46) and (cid:46) ν syn / Hz (cid:46) . According to our simulations, mildlymagnetized ones (dot-dashed lines in Fig. 1, and green dotsin Figures 2 and 3) have also these features. Moreover, in MNRAS , 1–11 (2015)
J.M. Rueda-Becerril, A.O. Harrison and D. Giannios
Figure 3.
Apparent velocity, Compton dominance and synchrotron peak. Similar to Fig. 2, squares, circles and triangles depict the modelswith baryon loading µ = , , , respectively. Blue, orange, green, red and purple colors show the simulation results with magnetization σ = , , , , respectively. The gray transparent arrow shows the increasing trend of the jet luminosity, L j . Left panel : We show theapparent velocity as a function of the synchrotron peak frequency ν syn . Observational data from Lister et al. (2019). Right panel : Weshow the Compton dominance as a function of the synchrotron peak. In red dashed vertical lines we separate the LBL ( (cid:46) Hz), IBL( (cid:38) Hz and (cid:46) Hz) and HBL ( (cid:38) Hz) regions. Observational data from Finke (2013).
Fig. 4 we can place our mildly magnetized model in the re-gion L j / L Edd ≈ . , which would correspond to an Eddingtonratio (cid:219) m (cid:38) . .In this work we have associated the most extreme ac-cretion systems with blazar jets with large Γ . The brightaccretion disk may dominate the ionizing flux received bythe gas clouds living in the BLR, obscuring the central en-gine and populating that space with a denser photon fieldfrom the reprocessed disk radiation. A denser photon field,in conjuction with the larger blulk Γ , translate into a moreluminous EIC component of the blazar SED. A denser exter-nal radiation field would also mean a strong cooling factor (cid:219) γ , steepening the EED. This agrees with recent findings byKeenan et al. (2020). They agree with the scenario in whichpowerful blazar have a broad-emitting gas surrounding thecore. This also agrees with recent findings of Zhang et al.(2020), regarding the jet properties of other kind of γ -rayemitting AGNs known as Compact Steep-spectrum Sources.Regarding the core surrounding environment, accordingto Ghisellini et al. (2011), there is a clear division betweenFSRQs and BL Lacs in the L BLR – L γ plane at L BLR / L Edd = × − (although for these results fewer sources are pre-sented). According to our model, this divide is not so clear.As we have mentioned before, mildly magnetized simulationshave been setup as FSRQ-like, however, comparing with ob-servables, these show BL Lac features as well. It may be thecase that there is not such a sharp divide between BL Lacsand FSRQs.Looking back into the Compton dominance plot (rightpanel of Fig. 3), if we focus on a particular value of µ , e.g., triangles, one can move through all the observational regionby increasing the jet luminosity, following the gray translu-cent arrow. In other words, blazar jets may indeed launchwith similar baryon loading. In low (cid:219) m systems, a fainter ac-cretion disc means low density external photon field sur-rounding the emission region and a less powerful jet (blueregion in the lower right region of Fig. 4). Jets in low (cid:219) m sys-tems have mostly the synchrotron photons produced in situ as seed photons for upscattering, showing a dim EIC, justlike BL Lac objects whose inner core shows no significantsign of a broad-emitting gas. The SSC component is there-fore dominant in these sources, although not expected to beas γ -ray loud as the EIC of powerful jets, where Dopplerboosting plays a leading role in enhancing the EIC compo-nent.In both panels of Fig. 3, the synchrotron peak ν syn ofpowerful jet simulations corresponds to the synchrotron fre-quency of the cooling break of the EED . It turns out that,because of strong radiative losses, powerful jets have a syn-chrotron peak deeper into the far infrared (bellow these fre-quencies, most of the synchrotron emission is self-absorbed).On the oposite side, ν syn of our simulated BL Lac objects(i.e., simulations with high magnetization and < p < ), The cooling break of a particle energy distribution correspondsto the energy at which the distribution changes slope and is givenby the cooling factor (cid:219) γ in the kinetic equation (25). This pointdepends on how fast particles are being cooled down. The analysisof the cooling stages of the EED in the emission region of ourblazar model here presented is beyond the scope of this work.MNRAS000
Fig. 4 we can place our mildly magnetized model in the re-gion L j / L Edd ≈ . , which would correspond to an Eddingtonratio (cid:219) m (cid:38) . .In this work we have associated the most extreme ac-cretion systems with blazar jets with large Γ . The brightaccretion disk may dominate the ionizing flux received bythe gas clouds living in the BLR, obscuring the central en-gine and populating that space with a denser photon fieldfrom the reprocessed disk radiation. A denser photon field,in conjuction with the larger blulk Γ , translate into a moreluminous EIC component of the blazar SED. A denser exter-nal radiation field would also mean a strong cooling factor (cid:219) γ , steepening the EED. This agrees with recent findings byKeenan et al. (2020). They agree with the scenario in whichpowerful blazar have a broad-emitting gas surrounding thecore. This also agrees with recent findings of Zhang et al.(2020), regarding the jet properties of other kind of γ -rayemitting AGNs known as Compact Steep-spectrum Sources.Regarding the core surrounding environment, accordingto Ghisellini et al. (2011), there is a clear division betweenFSRQs and BL Lacs in the L BLR – L γ plane at L BLR / L Edd = × − (although for these results fewer sources are pre-sented). According to our model, this divide is not so clear.As we have mentioned before, mildly magnetized simulationshave been setup as FSRQ-like, however, comparing with ob-servables, these show BL Lac features as well. It may be thecase that there is not such a sharp divide between BL Lacsand FSRQs.Looking back into the Compton dominance plot (rightpanel of Fig. 3), if we focus on a particular value of µ , e.g., triangles, one can move through all the observational regionby increasing the jet luminosity, following the gray translu-cent arrow. In other words, blazar jets may indeed launchwith similar baryon loading. In low (cid:219) m systems, a fainter ac-cretion disc means low density external photon field sur-rounding the emission region and a less powerful jet (blueregion in the lower right region of Fig. 4). Jets in low (cid:219) m sys-tems have mostly the synchrotron photons produced in situ as seed photons for upscattering, showing a dim EIC, justlike BL Lac objects whose inner core shows no significantsign of a broad-emitting gas. The SSC component is there-fore dominant in these sources, although not expected to beas γ -ray loud as the EIC of powerful jets, where Dopplerboosting plays a leading role in enhancing the EIC compo-nent.In both panels of Fig. 3, the synchrotron peak ν syn ofpowerful jet simulations corresponds to the synchrotron fre-quency of the cooling break of the EED . It turns out that,because of strong radiative losses, powerful jets have a syn-chrotron peak deeper into the far infrared (bellow these fre-quencies, most of the synchrotron emission is self-absorbed).On the oposite side, ν syn of our simulated BL Lac objects(i.e., simulations with high magnetization and < p < ), The cooling break of a particle energy distribution correspondsto the energy at which the distribution changes slope and is givenby the cooling factor (cid:219) γ in the kinetic equation (25). This pointdepends on how fast particles are being cooled down. The analysisof the cooling stages of the EED in the emission region of ourblazar model here presented is beyond the scope of this work.MNRAS000 , 1–11 (2015) aryon loading of blazar jets σ µ − − L j / L E dd Figure 4.
The baryon loading. The color gradient shows the areawhere, following relation (2), the bulk Lorentz factor ≤ Γ ≤ .The gray area depicts the low magnetization, σ < , region. corresponds to the synchrotron frequency of the maximumLorentz factor of the EED, γ (cid:48) max , given by Eq. (23). For thesecases, in contrast with simulations with low magnetization, γ (cid:48) max is highly dependent on γ (cid:48) min . This setup of the EEDsin our model doesn’t give any restriction or upper limit forthe synchrotron peak ν syn (see Keenan et al. 2020). However,from Eq. (23), observations can indeed constrain the valueof γ (cid:48) min . In the present work, we have applied a simple idea thataccounts for the blazar sequence and several observable fea-tures of the blazar population. This model relies on the ideathat all jets are launched with similar energy per baryon,independently of their power. FSRQs, those with the mostpowerful jets, manage to accelerate to high bulk Lorentzfactor and have luminosities (cid:38) . L Edd . FSRQ-like simula-tions were set to have a rather modest magnetization in theemission region and a steep particle energy distribution. Ourpredicted SEDs of these models show similar features as ac-tual FSRQs observations: peak synchrotron ν syn (cid:46) Hz,Compton dominance, soft spectra in the γ -rays, and are γ -ray louder. In the case of BL Lacs, the jet does not achievea very high bulk Lorentz factor, leading to more magneticenergy available for non-thermal particle acceleration. Ac-cording to our model (see Sec. 2), these sources develop highsynchrotron peak, weaker Compton component, and harderemission spectra at frequencies (cid:38) GeV.With our model and simulations reported in this work,we were able to recover observables of blazars. Namely, the blazar sequence was (qualitatively) reproduced, in a simi-lar manner as it was first reported by Fossati et al. (1998),for those models with mild baryon loading. This result con-strains the energy per baryon of blazar jets to (cid:46) µ (cid:46) .The L-like region observed for the apparent velocity andCompton dominance as functions of ν syn was also recovered by changing L j , assuming that it tracks (cid:219) m . With our simplemodel we are also able to show that the brightness of theBLR scales linearly with the γ -rays loudness of the source.Finally, we propose an indirect method to estimate γ (cid:48) min for BL Lacs. From the value of ν syn given by observationswe can directly calculate γ (cid:48) max . Following Eq. (23) we aretherefore able to calculate γ (cid:48) min . PIC simulations of magneticreconnection may be able to test whether our adopted valuesare reasonable.It is worth highlighting the particular case in which anFSRQ-like simulation (green points in figures 2 and 3), isin fact γ -ray quieter. This object would in principle have amild Eddington rate (cid:219) m , and a mildly luminous BLR. How-ever it is not powerful enough to develop an IC componentlouder than its synchrotron component. Additionally, it hasa harder spectral index α γ , and emits close the TeV band,just like BL Lacs. Similar “contradicting” properties havealso been observed in objects like TXS 0506+056.In summary, our model assumes that all jets are injectedwith energy per baryon in a narrow range (cid:46) µ (cid:46) andthat the jet bulk Lorentz factor and power scale positivelywith the accretion rate, and can account for or predict: • That (cid:219) m controls many of the observable features ofblazars such as the high-energy spectral index and lumi-nosity, the brightness of the BLR, the apparent speed, andthe synchrotron spectrum and synchrotron peak frequency. • Sources that are γ -ray brighter have softer γ -ray spec-tral index α γ . Lower values of α γ (i.e., harder spectra) werefound for the γ -ray quieter sources. • The BLR luminosity L BLR scales linearly with the γ -rayluminosity of the object. • Fastest objects have low-frequency synchrotron peak ν syn while objects with intermediate-to-high synchrotronpeak move rather slow. • Low jet luminosity sources are non-Compton dominantbut high synchrotron-peaked, whereas those with higherCompton dominance have a ν syn (cid:46) Hz.
ACKNOWLEDGEMENTS
The research was partly supported by Fermi Cycle 12 GuestInvestigator Program
REFERENCES
Ackermann M., et al., 2011, ApJ, 743, 171Ajello M., et al., 2014, ApJ, 780, 73Barniol Duran R., Tchekhovskoy A., Giannios D., 2017, MNRAS,469, 4957Begelman M. C., 1998, ApJ, 493, 291Bian W.-H., Zhao Y.-H., 2003, PASJ, 55, 599Blandford R. D., Znajek R. L., 1977, MNRAS, 179, 433Blinov D., et al., 2015, Mon. Not. R. Astron. Soc., 453, 1669MNRAS , 1–11 (2015) J.M. Rueda-Becerril, A.O. Harrison and D. Giannios
B¨ottcher M., Dermer C. D., 2002, ApJ, 564, 86Celotti A., Ghisellini G., 2008, MNRAS, 385, 283Chang J. S., Cooper G., 1970, J. Comp. Phys, 6, 1Christie I. M., Petropoulou M., Sironi L., Giannios D., 2019, MN-RAS, 482, 65Davis S. W., Laor A., 2011, ApJ, 728, 98Dermer C. D., Menon G., 2009, High Energy Radiation fromBlack Holes. Gamma Rays, Cosmic Rays, and Neutrinos.Princeton Series in Astrophysics, Princeton University Press,PrincetonEichler D., 1993, ApJ, 419, 111Finke J. D., 2013, ApJ, 763, 134Fossati G., Maraschi L., Celotti A., Comastri A., Ghisellini G.,1998, MNRAS, 299, 433Garofalo D., 2019, ApJ, 876, L20Gaskell C. M., 2009, New Astron. Rev., 53, 140Georganopoulos M., Kazanas D., Perlman E., Stecker F. W.,2005, ApJ, 625, 656Ghisellini G., Celotti A., 2001, A&A, 379, L1Ghisellini G., Tavecchio F., 2008, MNRAS, 387, 1669Ghisellini G., Celotti A., Fossati G., Maraschi L., Comastri A.,1998, MNRAS, 301, 451Ghisellini G., Maraschi L., Tavecchio F., 2009, MNRAS, 396, L105Ghisellini G., Tavecchio F., Foschini L., Ghirlanda G., MaraschiL., Celotti A., 2010, MNRAS, 402, 497Ghisellini G., Tavecchio F., Foschini L., Ghirlanda G., 2011, MN-RAS, 414, 2674Giannios D., Spruit H. C., 2006, A&A, 450, 887Giannios D., Uzdensky D. A., Begelman M. C., 2009, MNRAS,395, L29Giommi P., Padovani P., Polenta G., Turriziani S., D’Elia V.,Piranomonte S., 2012, MNRAS, 420, 2899Giommi P., Padovani P., Polenta G., 2013, MNRAS, 431, 1914Giustini M., Proga D., 2019, A&A, 630, A94Gould R. J., 1979, A&A, 76, 306Homan D. C., Kadler M., Kellermann K. I., Kovalev Y. Y., ListerM. L., Ros E., Savolainen T., Zensus J. A., 2009, Astrophys.J., 706, 1253Janiak M., Sikora M., Moderski R., 2015, MNRAS, 449, 431Jorstad S., Marscher A., 2016, Galaxies, 4, 47Kagan D., Sironi L., Cerutti B., Giannios D., 2015, Space ScienceReviews, 191, 545Kaspi S., Maoz D., Netzer H., Peterson B. M., Vestergaard M.,Jannuzi B. T., 2005, ApJ, 629, 61Kaspi S., Brandt W. N., Maoz D., Netzer H., Schneider D. P.,Shemmer O., 2007, ApJ, 659, 997Keenan M., Meyer E. T., Georganopoulos M., Reddy K., FrenchO. J., 2020, arXiv e-prints, p. arXiv:2007.12661Komissarov S. S., Barkov M. V., Vlahakis N., K¨onigl A., 2007,MNRAS, 380, 51Komissarov S. S., Vlahakis N., K¨onigl A., Barkov M. V., 2009,MNRAS, 394, 1182Komissarov S. S., Gourgouliatos K. N., Matsumoto J., 2019, MN-RAS, 488, 4061Kovalev Y. Y., et al., 2009, ApJ, 696, L17Lister M., 2016, Galaxies, 4, 29Lister M. L., et al., 2009, Astron. J., 138, 1874Lister M. L., et al., 2011, ApJ, 742, 27Lister M. L., et al., 2019, ApJ, 874, 43L´opez-C´amara D., Morsony B. J., Begelman M. C., Lazzati D.,2013, ApJ, 767, 19Maraschi L., Tavecchio F., 2003, ApJ, 593, 667Marscher A. P., Gear W. K., 1985, ApJ, 298, 114Matsumoto J., Masada Y., 2013, ApJ, 772, L1Meier D. L., 2002, New Astron. Rev., 46, 247Meyer E. T., Fossati G., Georganopoulos M., Lister M. L., 2011,Astron. J., 740, 98Mimica P., Aloy M. A., 2012, MNRAS, 421, 2635 Nalewajko K., Giannios D., Begelman M. C., Uzdensky D. A.,Sikora M., 2011, MNRAS, 413, 333Padovani P., 2007, Ap&SS, 309, 63Padovani P., Oikonomou F., Petropoulou M., Giommi P., ResconiE., 2019, MNRAS, 484, L104Park B. T., Petrosian V., 1996, ApJS, 103, 255Perucho M., Lobanov A. P., Mart´ı J. M., Hardee P. E., 2006,A&A, 456, 493Petropoulou M., Giannios D., Sironi L., 2016, MNRAS, 462, 3325Petropoulou M., Sironi L., Spitkovsky A., Giannios D., 2019, ApJ,880, 37Rani P., Stalin C. S., Rakshit S., 2017, Mon. Not. R. Astron. Soc.,466, 3309Rueda-Becerril J. M., 2017, PhD thesis, Universitat de Val`encia,Valencia, Spain, http://roderic.uv.es/handle/10550/60003
Rueda-Becerril J. M., 2020, Numerical approach for radiativecooling in the Klein-Nishina regime, In preparation.Rueda-Becerril J. M., Mimica P., Aloy M. A., 2014, MNRAS, 438,1856Rueda-Becerril J. M., Mimica P., Aloy M. A., 2017, MNRAS, 468,1169Rybicki G. B., Lightman A. P., 1979, Radiative processes in as-trophysics. Wiley-Interscience, New YorkSikora M., Madejski G., Moderski R., Poutanen J., 1997, ApJ,484, 108Sironi L., Spitkovsky A., 2014, ApJ, 783, L21Sironi L., Petropoulou M., Giannios D., 2015, MNRAS, 450, 183Sironi L., Giannios D., Petropoulou M., 2016, MNRAS, 462, 48Spruit H. C., Daigne F., Drenkhahn G., 2001, A&A, 369, 694Tavecchio F., Ghisellini G., 2008, MNRAS, 386, 945Tchekhovskoy A., Bromberg O., 2016, MNRAS, 461, L46Tchekhovskoy A., McKinney J. C., Narayan R., 2008, MNRAS,388, 551Tchekhovskoy A., McKinney J. C., Narayan R., 2009, ApJ, 699,1789Tchekhovskoy A., Narayan R., McKinney J. C., 2010, ApJ, 711,50Tchekhovskoy A., McKinney J. C., Narayan R., 2012, in Journalof Physics Conference Series. p. 012040 ( arXiv:1202.2864 ),doi:10.1088/1742-6596/372/1/012040Urry C. M., Padovani P., 1995, PASP, 107, 803Zhang J., Zhang H.-M., Gan Y.-Y., Yi T.-F., Wang J.-F., LiangE.-W., 2020, Astrophys. J., 899, 2
APPENDIX A: ACCRETION RATE
In the present section we will describe the parametriza-tion of our model. In our formulation, the accretion rateparameter is given by Eq. (7). The main effects of chang-ing values of the accretion index s are shown in Fig. A1.There we show the averaged SEDs from simulations with µ = , , (columns from left to right, respectively), ( σ, p ) = ( , . ) , ( , . ) , ( , . ) , ( , . ) , and ( , . ) (rowsfrom top to bottom, respectively), and s = . , . , . , . (blue, orange, green and red lines, respectively). The syn-chrotron, SSC, EIC and total fluxes are depicted in dashed,dot-dashed, dot-dot-dashed and solid lines, respectively. Forall simulations we set Γ = , and in each panel it isnoted the corresponding bulk Lorentz factor, Γ , accordingto Eq. (2).The first three rows (top to bottom) correspond to mod-els setup FSRQ-like, i.e., with low-to-mild magnetizationand p > . The first two are the brightest and the mostCompton dominant. In fact, the EIC component is the dom-inant radiative process in all this set of simulations. Not so MNRAS000
In the present section we will describe the parametriza-tion of our model. In our formulation, the accretion rateparameter is given by Eq. (7). The main effects of chang-ing values of the accretion index s are shown in Fig. A1.There we show the averaged SEDs from simulations with µ = , , (columns from left to right, respectively), ( σ, p ) = ( , . ) , ( , . ) , ( , . ) , ( , . ) , and ( , . ) (rowsfrom top to bottom, respectively), and s = . , . , . , . (blue, orange, green and red lines, respectively). The syn-chrotron, SSC, EIC and total fluxes are depicted in dashed,dot-dashed, dot-dot-dashed and solid lines, respectively. Forall simulations we set Γ = , and in each panel it isnoted the corresponding bulk Lorentz factor, Γ , accordingto Eq. (2).The first three rows (top to bottom) correspond to mod-els setup FSRQ-like, i.e., with low-to-mild magnetizationand p > . The first two are the brightest and the mostCompton dominant. In fact, the EIC component is the dom-inant radiative process in all this set of simulations. Not so MNRAS000 , 1–11 (2015) aryon loading of blazar jets the middle-row ones, which show an EIC component withsimilar brightness, or dimmer, than the synchrotron com-ponent. BL Lac-like models are those with higher magne-tization and lower Γ (last two rows from top to bottom).These simulations show synchrotron, SSC and EIC compo-nents with similar luminosities.The main effect that the normalization bulk Lorentzfactor Γ has on our simulations is the overall increase in lu-minosity. In the same manner, we noticed in the SEDs thatby increasing the accretion index s , overall brightness de-creases, but the overall spectral structure remains the same.Furthermore, this effect occurs regardless of the magnetiza-tion and baryon loading. From these results we can concludethat (cid:219) m regulates the intensity of the SEDs without changingany local nor broadband spectral feature. This was expectedaccording to Eq. (6), which tells us that (cid:219) m is a measure of L j .The cases that have reached the super-Eddington limit, i.e.,those models with (cid:219) m ≥ , appear in uppermost right panel.In our setup, this frontier is set by the parameter Γ . This paper has been typeset from a TEX/L A TEX file prepared bythe author.MNRAS , 1–11 (2015) J.M. Rueda-Becerril, A.O. Harrison and D. Giannios Γ = Γ = Γ = Γ = Γ = Γ = l og ν L ν [ e r g s − ] Γ = Γ = Γ = Γ = Γ = Γ =
12 14 16 18 20 22 24 26 log ν [Hz] Γ =
12 14 16 18 20 22 24 26 log ν [Hz] Γ =
12 14 16 18 20 22 24 26 log ν [Hz] Γ = Figure A1.
In this figure we show the averaged SEDs of our simulations with µ = , and 90 in the left, middle and right columns,respectively. Simulations with ( σ, p ) = ( , . ) , ( , . ) , ( , . ) , ( , . ) , and ( , . ) are shown from top to bottom, respectively. Thesolid, dashed, dot-dashed and dot-dot-dashed lines correspond to the total, synchrotron, SSC and EIC components, respectively. In blue,orange, green and red are depicted the simulations with accretion index s = . , . , . , . , respectively. The normalization bulk Lorentzfactor in Eq. (7) is set to Γ = . The spectra are averaged over 1 dy since particles start being injected into the emitting blob. Thevalue of the bulk Lorentz factor Γ shown in each panel is given by Eq. (2). MNRAS000