Hadronic supercriticality & GRB emission

Hadronic supercriticality as a trigger for GRB emission

M. Petropoulou , S. Dimitrakoudis, A. Mastichiadis , D. Giannios
Department of Physics and Astronomy, Purdue University, 525 Northwestern Avenue, West Lafayette, IN, 47907, USA
NASA Einstein Postdoctoral Fellow
Institute for Astronomy, Astrophysics, Space Applications & Remote Sensing, National Observatory of Athens, 15 236 Penteli, Greece
Department of Physics, University of Athens, Panepistimiopolis, GR 15783 Zografos, Greece
E-mail:mpetropo@purdue.edu (MP)E-mail: sdimis@noa.gr (SD)E-mail: amastich@phys.uoa.gr (AM)E-mail: dgiannio@purdue.edu (DG)

We explore a one-zone hadronic model that may be able to reproduce -ray burst (GRB) prompt emission with a minimum of free parameters. Assuming only that GRBs are efficient high-energy proton accelerators and without the presence of an ab initio photon field, we investigate the conditions under which the system becomes supercritical, i.e. there is a fast, non-linear transfer of energy from protons to secondary particles initiated by the spontaneous quenching of proton-produced -rays. We first show analytically that the transition to supercriticality occurs whenever the proton injection compactness exceeds a critical value, which favours high proton injection luminosities and a wide range of bulk Lorentz factors. The properties of supercriticality are then studied with a time-dependent numerical code that solves concurrently the coupled equations of proton, photon, electron, neutron and neutrino distributions. For conditions that drive the system deep into the supercriticality we find that the photon spectra obtain a Band-like shape due to Comptonization by cooled pairs and that the energy transfer efficiency from protons to -rays and neutrinos is high reaching . Although some questions concerning its full adaptability to the GRB prompt emission remain open, supercriticality is found to be a promising process in that regard.

astroparticle physics – instabilities – radiation mechanisms: non-thermal –gamma ray burst: general
pagerange: Hadronic supercriticality as a trigger for GRB emissionLABEL:lastpagepubyear: 2014

1 Introduction

The prompt emission of gamma-ray bursts (GRBs) is observed in the  keV- MeV energy band (Preece et al., 2000; von Kienlin et al., 2014), with fluences that generally range between  erg cm and  erg cm, where the lower limit does not necessarily reflect an intrinsic property of GRBs but depends on the sensitivity of the detectors. The fluence distribution of bursts detected by the Fermi/GBM peaks at  erg cm (von Kienlin et al., 2014), which is considered as the typical GRB fluence. GRB light curves are highly variable and have a complex structure: they consist of several pulses (10-100) each of them having typical width ms-s (Norris et al., 1996; Nakar & Piran, 2002), thus making the total duration of the burst longer, e.g.  s. The GRB spectra can, in most cases, be described by a smoothly connected broken power-law (Band et al., 1993, 2009) with values of the break energy clustered around  MeV in the observer’s frame, while even higher peak energies ( MeV) have been detected (e.g. Goldstein et al. 2012). The typical photon indices of the spectrum below and above the peak are and , respectively (Preece et al., 2000; Goldstein et al., 2012).

The origin of the GRB emission is still an open issue, although various models have been proposed during the past decades trying to address all or most of the above properties. Borrowing the terminology from the field of blazar modelling (see e.g. Böttcher 2007; Boettcher 2010, for reviews) GRB emission models can be divided in two categories, namely leptonic and hadronic, according to the species of the radiating particles.

The former try to attribute the gamma-ray emission by employing radiation processes of relativistic electrons (and/or positrons). The classical scenario of the prompt emission, which belongs to the first category, is the optically thin synchrotron model (e.g. Katz 1994; Sari et al. 1996; Tavani 1996; Chiang & Dermer 1999), where the kinetic energy of the flow is dissipated via shocks and the prompt emission is the result of synchrotron radiation of relativistic electrons. The difficulties that this scenario has in dealing with several issues, such as the low-energy photon index (Crider et al., 1997; Preece et al., 1998, 2000) and the fine tuning of parameters required to explain the MeV peak of GRB spectra (see also Beloborodov 2010 for relevant discussion), motivated works on alternative scenarios. For example, variants of non-thermal emission models were discussed in order to overcome the ‘line-of-death’ problem: (i) effects of adiabatic and/or inverse Compton (IC) cooling in the Klein-Nishina regime on the low energy part of the synchrotron spectrum (Derishev et al., 2001; Wang et al., 2009; Daigne et al., 2011); (ii) synchrotron emission from an electron distribution with a smooth low energy cutoff and an anistropic distribution in pitch angles (Lloyd & Petrosian, 2000); (iii) the emission observed in the BATSE energy band 20 keV-1 MeV was interpeted as the result of IC scattering of slow cooling electrons on the self-absorbed part of the synchrotron spectrum (Panaitescu & Mészáros, 2000); (iv) synchrotron self-Compton emission under the assumption of continuous electron acceleration (Stern & Poutanen, 2004); (v) jitter radiation emitted by relativistic electrons moving in non uniform small scale magnetic fields (Medvedev, 2000); (vi) synchrotron electron cooling in a decaying magnetic field (Pe’er & Zhang, 2006); (vii) gamma-ray emission through the Compton-drag process (Lazzati et al., 2000). The so-called photospheric models, where the radiation is released when the outflow becomes transparent, constitute an interesting alternative to the non-thermal ones. If the energy is dissipated at the very inner parts of the outflow, it thermalizes and the radiation that escapes from the GRB photosphere has a quasi-thermal spectrum that peaks at  MeV (e.g. Goodman 1986; Thompson 1994; Beloborodov 2010). In the presence of continuous energy dissipation, however, the resulting spectra may obtain a non-thermal appearance via Comptonization of the quasi-thermal emission by thermal electrons (Mészáros & Rees, 2000; Pe’er et al., 2006; Giannios, 2006, 2012).

Hadronic models for the GRB prompt emission constitute a viable alternative. They are built upon the common basis that the (sub)MeV -ray emission, which serves as the target field for photopion interactions, is not of hadronic origin but it is either the synchrotron radiation of primary electrons (e.g. Dermer & Atoyan 2003; Asano & Inoue 2007; Murase 2008) or the emission from the photosphere itself (e.g. Gao et al. 2012; Asano & Mészáros 2013). The high-energy part of the gamma-ray spectrum (100 MeV) is typically explained by relativistic proton synchrotron radiation (Vietri, 1997; Totani, 1998) or by proton-induced cascades (Dermer & Atoyan, 2006; Asano & Inoue, 2007; Asano et al., 2009). The latter scenario has been applied to explain the underlying power-law components seen in some bright Fermi bursts (e.g. GRB 090902B (Abdo et al., 2009); GRB 080319B (Racusin et al., 2008)), which extend from the hard X-rays up to GeV energies and do not agree with simple extrapolations of the MeV spectrum (Asano et al., 2010). In any case, the suggestion that GRBs are the sources of ultra-high energy cosmic rays (UHECRs) (Waxman, 1995b; Vietri, 1995; Murase et al., 2008) makes hadronic models attractive. Moreover, the associated high-energy neutrino emission has been calculated in various studies (Waxman & Bahcall, 1997; Murase, 2008; Gao et al., 2012; He et al., 2012; Zhang & Kumar, 2013; Asano & Meszaros, 2014; Baerwald et al., 2014; Reynoso, 2014; Petropoulou et al., 2014) and now starts becoming testable by ongoing observations (IceCube Collaboration, 2013; Aartsen et al., 2014b, a).

However, one-zone hadronic models are inherently more complex than pure leptonic ones, since they require modelling of the coupled emission and energy loss processes (see e.g. Dimitrakoudis et al. 2012a) between various species in order to track the evolution of the different components (protons, neutrons, pairs, mesons, neutrinos, photons). Several of these feedback processes were proven to give rise to radiative instabilities (e.g. Stern & Svensson 1991; Kirk & Mastichiadis 1992; Mastichiadis et al. 2005) that share a common feature: the abrupt, i.e. in a few dynamical times, release of energy that is initially stored in protons and is subsequently transfered to photons. The proton synchrotron pair-production instability for example was proposed to give rise to gamma-ray emission peaking at  MeV (e.g. Kazanas et al. 2002; Mastichiadis & Kazanas 2006) offering, at the same time, a physical connection between the prompt and afterglow phases (Mastichiadis & Kazanas, 2009; Sultana et al., 2013). Moreover, possible implications of the automatic -ray quenching instability (Stawarz & Kirk, 2007; Petropoulou & Mastichiadis, 2011) were studied in Petropoulou & Mastichiadis (2012a, b) in the context of hadronic blazar emission.

In the present work we extend this analysis by exploring its emission signatures in gamma-rays, high-energy neutrinos and cosmic-rays for parameters relevant to GRB sources. We begin with the sole assumption of a source that is an efficient ultra high-energy (UHE) proton accelerator and sufficiently magnetized in order to confine the accelerated protons. Contrary to the majority of GRB hadronic models we do not assume an external source of photons. We show that for low values of the proton injection luminosity, the only photon field present in the source is the one emitted by the proton component mainly through synchrotron radiation. In this case, the emitted spectrum cannot be assigned to that of a typical GRB because of its spectral shape and its low luminosity. We follow the evolution of protons by balancing their losses to the respective gains of their secondaries as done in Mastichiadis & Kirk (1995) and show that if their injection luminosity exceeds a critical value, the system undergoes a transition that is triggered by the instability of automatic -ray quenching. The transition is easily identified by an abrupt increase of the photon luminosity that causes the source to enter in a high photon compactness state making the energy exchange between leptons and photons dominant. This results in photon spectral shapes that in general resemble GRBs in the sense that (i) they match the required luminosity and (ii) they can be fit by a Band function. It is this self-consistently produced radiation field that becomes the target for photopion interactions and efficiently drains energy from UHE protons, part of which is transfered to high-energy electron and muon neutrinos produced through the charged pion decay. This constitutes one of the fundamental differences between the present study and others where an ad hoc Band spectrum is assumed (see also Murase et al. 2008 for a similar approach to ours). Finally, relativistic neutrons produced from the same photopion interactions provide an effective means for UHECR escape from the source, since they are not magnetically confined and their decay time is long enough as to allow them to escape freely before converting into protons through decay. We show that for modest values of the bulk Lorentz factor () all three components, namely neutrinos, UHECRs and photons, are energetically similar.

The present work is structured as follows. In §2 we present the physical conditions of the model and derive an analytical expression for the critical proton injection luminosity. In §3.1 we present the numerical code and in §3.2 we continue with a presentation of the photon, UHECR and neutrino emission spectra; we also discuss the effects of a lower value for the high-energy cutoff of the proton distribution. We discuss our results in §4 and conclude with a summary in §5. Throughout this study we use  km Mpc s, , and as an indicative value for the redshift of the gamma-ray source. We also introduce the notation .

2 Analytical approach

2.1 Physical conditions

The GRB lightcurve in the soft -ray band (10 MeV) is variable and consists of several pulses with durations in the range  s (Norris et al., 1996; Nakar & Piran, 2002). Here we adopt the internal shock scenario (see Piran (2004); Zhang & Mészáros (2004), for reviews) according to which the emitting region that corresponds to each individual pulse is modelled as a homogeneous shell with Lorentz factor that forms at a distance from the central engine. In the comoving frame the width of the shell is


As long as the beaming angle is smaller than the opening angle of the jet, which holds during the internal shock phase, we may treat the emission region as a spherical blob of radius .

Since the typical isotropic energy emitted in -rays is  erg (Bloom et al., 2003; Kocevski & Butler, 2008), the isotropic -ray luminosity defined as ranges between and  erg/s, for a fiducial burst duration  s. A minimal requirement is that , where where is the total power of the jet and equals to with , and being the kinetic, Poynting and proton111 refers to the luminosity of a power-law proton distribution. luminosities, respectively. Although electron acceleration at high energies is expected to take place too, here, in our attempt to keep the number of free parameters as low as possible, we assume that the injection luminosity of primary relativistic electrons is much lower than that of protons, making their contribution to the energetics and the overall spectra negligible.

We introduce next the parameters and . We use throughout the text and as indicative values, unless stated otherwise. Using the definition of the Poynting luminosity


we may write the magnetic field strength measured in the comoving frame as follows:


As we show next, all physical quantities in addition to and may be expressed through five essential variables: , and .

We proceed with the derivation of the respective expressions for the two basic quantities that describe the proton distribution, namely its injection compactness and its high-energy cutoff. The former, is defined as


and using eq. (1) it is also written as


We assume that protons are being injected into the blob after having been accelerated into a power-law distribution with index starting from up to a Lorentz factor which is usually determined by the balance between the acceleration and the energy loss processes. Because protons do not, in principal, suffer as severe radiative losses as electrons do, the acceleration of protons to large in GRBs is possible (e.g. Asano et al. 2009; Murase et al. 2012). For the lower cutoff of the proton distribution and the power-law index we use the indicative values and respectively. We note that the exact value of does not alter the main results of this work, while steeper proton spectra would increase the energy demands; for this reason we will not consider such cases here. Assuming that the synchrotron losses are the dominant energy loss process for high energy protons (see also Petropoulou et al. 2014) and that the acceleration process operates close to the Bohm diffusion limit (e.g. Giannios 2010), the typical energy loss and acceleration timescales are given respectively by and , where and . Using the above expressions and eq. (3) we find


Thus, protons can in principle be accelerated to ultra-high energies (UHE), i.e.  eV in the comoving frame. If the gyroradius of these highly energetic protons is, however, larger than the typical size of the emission region , they cannot be confined and escape from it. Thus, the maximum energy of protons in the emission region is given by , where is derived using the Hillas criterion (Hillas, 1984), namely . Using eqs. (1) and (3) this is also written as


Combining eqs. (6) and (7) we find that , unless


For the fiducial parameter values used here and for the purposes of the analytical treatment presented in §2.2. it is sufficient to define the maximum proton Lorentz factor through eq. (7). However, in §3 where we study the problem numerically for different parameter sets we use the appropriate expression for .

2.2 Transition to supercriticality

The injected protons will emit synchrotron radiation that peaks in the comoving frame at , where ,  G, and ; here, is the typical Lorentz factor of protons that cool due to synchrotron losses within the dynamical timescale and it is given by


Using eqs. (3) and (7) we find that


where for simplicity we assumed (for the validity of the assumption, see Appendix A). Thus, for the fiducial parameter values used here the peak of the proton synchrotron spectrum falls well within the -ray energy band. As long as


is satisfied222For the derivation we used the feedback criterion of automatic photon quenching – see Petropoulou & Mastichiadis (2011) and eq. (2) therein., photons at the high-energy part of the proton synchrotron spectrum can be successible to the instability of spontaneous gamma-ray quenching (Stawarz & Kirk, 2007; Petropoulou & Mastichiadis, 2011). According to this, the -ray compactness () produced in a highly magnetized region cannot become arbitrarily high. Whenever it exceeds a critical value () it is spontaneously absorbed producing relativistic pairs that cool by emitting a large333An electron with energy will emit photons, where . number of synchrotron photons, thus providing more targets for further absorption. As the -ray compactess in our framework is related to the proton injection compactess, the existence of an upper limit to is also translated to a limiting value for the proton injection compactess (). This has been already pointed out by Petropoulou & Mastichiadis (2012a), although in a different context.

We refer the reader to Appendix A for the derivation of and here we present the final result


where the first branch is relevant for and the second otherwise. From this point on we will refer to cases with as ‘subcritical’ and ‘supercritical’ otherwise. We will also drop the logarithmic dependance and set instead the numerical factor in the bracket of eq. (13) to 20, i.e. .

Figure 1: Proton injection compactness (solid lines) and proton critical compactness (dashed lines) as a function of the Lorentz factor for three values of marked on the plot. Thick lines denote , while the loci of points with and 10 are plotted with red lines. Other parameters used are: , and  s.

Note that the above expression is valid as long as the gamma-rays that are spontaneously absorbed and initiate the instability are the result of proton synchrotron radiation. If the parameters are such as to push the peak of the proton synchrotron radiation to GeV energies, e.g. , then the emission from pairs produced by Bethe-Heitler and/or photopion interactions of protons with their own synchrotron radiation dominates in TeV energies. It can be shown that even in such cases the instability can still operate (Petropoulou & Mastichiadis, 2012a). Although the derivation of an expression similar to eq. (13) in this case is out of the scope of the present study, we will present a detailed numerical example in §3.

The absolute value of is not so important by itself. More important for the evolution of the system is the ratio which measures how deep in the supercritical regime the system is driven for given physical conditions. Using eqs. (5) and (13) this is written


where we used the first branch of eq. (13) for simplifying reasons. In any case, both expressions of are similar. Interestingly, eq. (14) shows that the condition is satisfied for a wide range of parameter values, with lower values of being prefered. Because of the strong dependance of the ratio on , slightly different values of the bulk Lorentz factor lead to very different photon and neutrino spectra as we show in §3. It is useful, therefore, to define a ‘critical’ value of the Lorentz factor too. Setting we find


which has very weak dependance on the parameters.

The above are summarized in Fig. 1 where we plot (solid lines) and (dashed lines) for three values of marked on the plot. Other parameters used are: , and  s. The dashed lines are plotted up to (see eq. (12)), as the derived expression for the critical compactness is not relevant for larger values of . For each value of we plotted with thick lines the part of the curve that lies above . Finally, the red lines denote the loci of points with and 10. The parameter space that leads to supercriticality becomes wider as increases. For high enough values, e.g.  erg/s, supercriticality is ensured for almost all values of relevant to GRBs. Note also that low values of the bulk Lorentz factor not only favour the transition to supercriticality but also correspond to .

The condition implies that the system lies deep in the supercritical regime and, as we will show in the next section with detailed numerical examples, this has the following implications: (i) the neutrino production efficiency is high, (ii) UHE protons cool down effectively through photopair and photopion processes, and (iii) the photon spectrum may be adequately described by a Band function (Band et al., 1993, 2009).

2.3 Comparison with the WB model

Waxman and Bahcall derived an elegant expression for the energy lost by protons through pion production within a dynamical timescale (see eq. (4) in Waxman & Bahcall 1997), which simply depends on the observed peak energy and -ray luminosity as well as on and . Since proton cooling due to photopion interactions becomes more inefficient for larger values of , their expression can also be seen as an upper limit for


In other words, in the WB model, pion and neutrino production is efficient for , which requires either low values of the Lorentz factor if  s or modest values of , e.g. , if the -ray variability is extremely fast, i.e.  ms (see e.g. Waxman & Bahcall 1997; Guetta et al. 2004; Abbasi et al. 2010). In our framework, however, efficient pion production is ensured, even without the requirement of an ab initio target photon field, for parameters leading to the supercritical regime (see e.g. Fig. 9 in §3). In the previous paragraph we showed by analytical menas that the transition to supercriticality occurs for , where is defined in eq. (15). A comparison between the two Lorentz factors is shown in Fig. 2, where and are plotted with dashed and solid lines, respectively for  erg/s,  MeV,  erg/s, and .

The grey colored area denotes the parameter space where the transition to supercriticality occurs due to the non-linear feedback loops leading to efficient pion production. As we show in §3, the self-consistently produced gamma-ray spectrum starts resembling a typical GRB one, for or . Thus, our model is equivalent to the WB model in the sense that efficient pion production on a Band-like gamma-ray spectrum is ensured for roughly similar parameter values, but with a fundamental difference: here the gamma-ray spectrum is not assumed a priori but is produced self-consistently by a series of processes, which we describe in detail in §3.2.

Figure 2: - plane and the two characteristic Lorentz factors derived by our analysis (solid line) and the analysis of Waxman & Bahcall 1997 (dashed line). The grey colored region denotes additional parameter space with respect to the WB model, where efficient pion production can occur due to non-linear feedback loops. For the rest of the parameters used, see text.

3 Numerical investigation

3.1 Numerical code

In the previous section we showed by analytical means that the transition to supercriticality is ensured for a wide range of parameter values. To verify this we employ the time-dependent numerical code as presented in Dimitrakoudis et al. (2012a) – hereafter DMPR12, that follows the evolution of protons, neutrons, secondary pairs, photons and neutrinos by solving the coupled differential equations that describe the various distributions. The coupling of energy losses and injection introduces a self-consistency in this approach that allows the study of the system at various conditions, e.g. in the presence of non-linear electromagnetic (EM) cascades and other feedback loops (see also Petropoulou & Mastichiadis 2012b for a relevant discussion). As a word of caution we stress that the aforementioned loops can be fully understood only if the coupled kinetic equation approach is used. The often used Monte Carlo techniques are intrinsically linear and fail to capture complex, non-linear effects such as this. While Monte Carlo codes are an excellent tool for the description of the system in the subcritical regime, should the choice of parameters drive the system into the supercritical regime the results of such codes may be in error.

We assume that protons are being injected in the source at a constant rate given by


where , , , and is the time measured in the comoving frame in units. Protons, as well as secondary particles, are allowed to leave the emission region in an average time . This may account in an approximate way for the expansion of the source, since the steady-state particle distributions derived by solving a kinetic equation containing a physical escape term or an adiabatic loss term are similar. All particles, primary and secondary, lose energy through various processes. Although details can be found in DMPR12, for the sake of completeness, we summarize here the physical processes that are included in the code:

  • proton-photon pair production (photopair)

  • proton-photon pion production (photopion)

  • neutron-photon pion production

  • proton synchrotron radiation

  • pion, kaon, muon and electron synchrotron radiation

  • synchrotron self-absorption

  • electron inverse Compton scattering

  • photon-photon pair production

  • electron-positron pair annihilation

  • Compton scattering of photons by cooled pairs

Photohadronic interactions are modelled using the results of Monte Carlo simulations. In particular, for Bethe-Heitler pair production the Monte Carlo results by Protheroe & Johnson (1996) were used (see also Mastichiadis et al. 2005). Photo-pion interactions were incorporated in the time-dependent code by using the results of the Monte Carlo event generator SOPHIA (Mücke et al., 2000). Synchrotron radiation of charged pions and muons was not included in the version of the code presented in DMPR12 and for the exact treatment we refer the reader to Dimitrakoudis et al. (2014). Pairs that cool down to Lorentz factors contribute to the Thomson depth and they are treated as a separate population. Following Lightman & Zdziarski (1987), we assume that this population thermalizes at a temperature , where . For the pair annihilation and photon downscattering processes we followed Coppi & Blandford (1990) and Lightman & Zdziarski (1987), respectively (for more details see Mastichiadis & Kirk 1995).

The photon escape timescale , in particular, is modelled as


where is the photon energy in units and is a function that equals unity for , it decreases as for and it becomes zero for . Moreover, and is the Klein-Nishina cross section. On the one hand, the above expression takes into account in an approximate way the fact that photons may be ‘trapped’ in the source for longer than one crossing time, i.e. for . On the other hand, eq. (18) does not take into account the effect of the expansion of the source during the photon escape. If the source expands on a dynamical time, the photon escape time is found to be (Giannios, 2006). Note, however, that in the examples shown here , making thus the effect of expansion rather modest.

3.2 Results

We investigated in total 10 parameter sets that were divided in two groups, namely A and B. All simulations in groups A are obtained for  erg/s, , and  s, whereas for group B  s. In each group we performed 5 simulations with Lorentz factors varying between and with a logarithmic step of 0.1. All other parameters that are used as an input for the numerical code, i.e. , , , and , are then derived using eqs. (1), (3), (6)-(7), and (5) respectively. These are summarized in Table 1. In all cases we let the system reach a steady-state, where the photon and neutrino emission as well as the proton and neutron energy distributions were then calculated. Although, in most cases, a steady-state is achieved in , a time-dependent treatment of the GRB emission that is intrinsically variable, is more adequate and it will be the subject of a future work.

#  (G)  (cm) a
Group A:  s
Group B:  s
  • It is calculated using the first branch of eq. (13) without the logarithmic dependance.

Table 1: Parameter values used for the calculation of the photon, neutrino and proton energy spectra shown in Figs. 3 and 7. Other parameters used are:  erg/s, and .

3.2.1 Emission spectra

Figure 3: Observed spectra for a single GRB pulse with duration  s at redshift . Photon, total neutrino and proton escaping luminosities are shown with solid, dotted and dashed lines, respectively, for (black lines), (blue lines) and (grey lines). The first two cases are supercritical with and 7, whereas the last one is subcritical with . Other parameters used are:  erg/s, , .

Figure 3 shows the observed multiwavelength photon spectra (solid lines) obtained from a single GRB pulse for  erg/s, , and  s (Group A) and for three indicative values of , i.e. (black lines), (blue lines) and (black lines). For comparison reasons, the total neutrino444We refer to the sum of electron/muon neutrino and antineutrino fluxes as the total neutrino flux. (dotted lines) and proton escaping (dashed lines) energy spectra are overplotted. Since is kept fixed, the proton injection compactness increases for decreasing (see also eq. (5)). In particular, it increases gradually from (lower curve) to (upper curve) with increaments of 0.5 in logarithm. The ratio calculated using eq. (14) for each case is given in label of Fig. 3.

Spectra shown with grey lines are obtained for , which for the particular choice of parameters, leads to low and . This is a typical example of emission signatures obtained when the system operates in the subcritical regime. The photon emission is characterized by the following: the MW spectra are dominated by the proton synchrotron component, i.e. synchrotron radiation is the dominant energy loss mechanism for UHE protons, and the radiative efficiency, which is defined as , is low in agreement with typical proton synchrotron emission models (see e.g. Mücke & Protheroe (2001)). In the subcritical regime, the proton spectra at steady state are the same as at injection because of negligible energy losses. Note that the sharp cutoff of the proton spectra at reflects the injection spectrum . The high energy cutoff would have been smoother if we were to use a more physically motivated injection spectrum, e.g. . The produced neutrinos in this regime are the result of photopion interaction of protons with their own emitted synchrotron radiation (see also DMPR12). However, the photopion energy loss rate of protons is very low and the neutrino emission is supressed. Note that the peak neutrino luminosity is .

As the injection compactness increases, at some point it will exceed the critical value (see eq. (13)) and the system will undergo a phase transition, which can be easily identified by a radical change of the spectral shape and an abrupt increase of the emitted luminosity (see black and blue lines in Fig. 3). Note that all spectra obtained in the supercritical regime are also characterized by in agreement with the analysis of §2.2. The underlying reason for this abrupt transition is the instability of spontaneous -ray quenching that redistributes the energy from the -ray energy band to lower energy parts of the spectrum through an EM cascade (Stawarz & Kirk, 2007; Petropoulou & Mastichiadis, 2011). Because of the increased production of low energy photons in the source, it is the photopion/photopair energy loss channels that are favoured over the the proton synchrotron one while in the supercritical regime. Thus, the EM cascade initiated by the instability of spontaneous -ray quenching is further supported by the injection of secondary pairs, which are highly relativistic and produced by the photopair and photopion interactions of UHE protons with the automatically generated soft photons. In other words, once proton-produced rays reach a certain compactness they are spontanteously absorbed, giving rise to electron-positron pairs and radiation, which causes more proton cooling via photopair and photopion processes; and eventually more -rays, thus sustaining the loop, which is illustrated in Fig. 4.

Figure 4: Schematic diagram of the feedback loop that is formed whenever spontaneous quenching of proton-produced -rays takes place (dashed lines).
Figure 5: Plot of the Thomson optical depth as a function of in logarithmic scale for the same parameters as in Fig. 3. The values corresponding to the injection and to the steady state, as derived from the numerical code, are shown with circles and squares, respectively. The analytical value of for three values of is also shown.

All the processes discussed above lead naturally to a large number density of cooled pairs () that is related to the source’s Thomson depth () through . Since the number of secondaries produced in the EM cascade is related to the proton injection compactness, we expect that larger values of lead to larger optical depths. This trend is exemplified in Fig. 5 where we plot the Thomson optical depth versus at the injection (circles) and at the steady state (squares) for the same parameters as in Fig. 3. The hatched area corresponds to the supercritical regime, while indicative values of (eq. 14) are marked on the plot for comparison reasons. At injection we assumed that any electrons present in the source are these related to the injected protons for conservation of neutrality and, hence, . At the steady state we find for , and this is another manifestation of the increased secondary pair injection while the system lies in the supercritical regime.

If the value of is such as to result in (blue lines in Fig. 3), then the process that determines the spectral shape is photon Comptonization by cooled electrons. In particular, our numerical simulations indicate that the peak energy of the photon spectrum scales approximately as , which is in good agreement with the dependance found by the solving the Kompaneets equation (Kompaneets, 1956). Although similar work can be found in the literature (e.g. Arons (1971); Illarionov & Syunyaev (1972); Lightman et al. (1981)), we present in Appendix B an analytic derivation for the case of continuous power-law photon injection that better describes the physical problem under investigation, and further supports our numerical results. Note that if we were to neglect the effects of photon downscattering, then all spectra obtained in the supercritical regime would have had the universal shape of a spectrum that peaks at  GeV in the observer’s frame (see e.g. blue line in Fig. 3), which reminds of the power-law underlying component seen in several bright GRBs (e.g. 080319B, 090902B, 090926A) detected with Fermi (Racusin et al., 2008; Abdo et al., 2009; Ackermann et al., 2011). Cascade emission produced through proton interactions with the MeV photons of the GRB were suggested as an alternative explanation for this emission (e.g. Asano et al. 2010). Interestingly, we derive similar photon spectral shapes which are also the result of a cascade but with the key difference that a photon field is not required ab initio; the targets are provided by the quenching loop. In the succession of spectra shown in Fig. 3 only those obtained for correspond to large enough and to start resembling to a Band-like photon spectrum (Band et al., 1993). The exact spectral shapes, however, should be considered with caution, since a better description for the cooling of pairs below is required and which we plan to address in the future.

As the system is driven deeper to the supercritical regime the photopion energy losses become gradually more significant and the energy drained from UHE protons is transfered to photons and secondary particles, such as pairs and neutrinos. The abrupt increase of the neutrino fluence is an additional sign of the transition to the supercritical regime – see grey and black lines in Fig. 3. As the proton injection compactness progressively increases (bottom to top) the neutrino spectrum becomes flatter in units, and it extends to lower energies tracing the evolution of the peak photon luminosity. The neutrino spectra for share many common features with those obtained in studies where the photon target field is modelled by an ab initio Band spectrum (e.g. Murase 2008; Baerwald et al. 2011, 2013; Petropoulou 2014).

Figure 6: Observed spectra for protons (solid lines) and neutrons (dashed lines) for the same parameters as in Fig. 3.
Figure 7: Comparison of observed photon, neutrino and proton energy spectra for  s and  s shown with black and blue lines, respectively. Spectra in panels (a) and (b) are obtained for and , respectively. All other parameters are same as in Fig. 3.
Figure 8: Gamma-ray spectra obtained for  s (black lines) and  s (grey lines). Other parameters used are:  erg/s, , and . The bowties show the range of observed values for the GRB photon indices, namely and (Preece et al. 2000), and are plotted for guiding the eye.

Summarizing, conditions leading to high photon and neutrino luminosities cause unavoidable cooling of UHE protons (see also Asano 2005) that obtain a steeper power-law spectrum than that of injection, i.e. , where . Figure 6 compares the observed differential luminosity of escaping protons and neutrons for the same parameters as in Fig. 3. In all three cases, we find that protons with energies  eV are unaffected by cooling due to photohadronic processes and actually serve as a large energy reservoir. However, their contribution to the spectral and temporal properties of the gamma-ray spectra is negligible both in the subcritical and supercritical regimes, and this can be understood as follows: protons with  eV emit mainly through synchrotron at low photon energies, e.g. and with a much lower luminosity than the gamma-ray one555The gamma-ray luminosity in the subcritical regime is given by the peak luminosity of the proton synchrotron component, while in the supercritical regime is given by the peak of the cascade emission component (see Fig. 3).. In this context it is the high energy part of the proton distribution that is ‘active’. The luminosity carried by neutrons, which are produced via the photopion channel , increases as the cooling of UHE protons becomes progressively more significant, i.e. as the conditions lead the system deeper into the supercritical regime. At the same time, the peak of the neutron energy spectrum moves towards lower energies, similarly to the cooling break energy of the proton energy spectrum, and the neutron distribution can be described by the same power-law as cooled protons, i.e., with index . In the optically thin limit for photopion interactions the produced neutron distribution would follow the proton injection spectrum. Here, the steepening of the neutron spectrum is the result of efficient neutron cooling through before escape from the emission region. Contrary to protons, neutrons are not confined by magnetic fields and their decay time is large enough to allow them to escape freely before converting through -decay into protons. These will propagate as UHECRs into the intergalactic medium having an injection spectrum similar to that of their parent population (Kirk & Mastichiadis, 1989; Begelman et al., 1990; Giovanoni & Kazanas, 1990; Mannheim et al., 2001; Atoyan & Dermer, 2003). Other possible escape mechanisms of UHECRs are discussed in Baerwald et al. (2013); Asano & Meszaros (2014).

Having explained the basic features of the photon, neutrino and proton spectra as the system is driven progressively from the subcritical to the supercritical regime for a particular parameter set (Group A), we proceed to investigate the role of other parameters, such as – see Fig. 7. If all other parameters are kept fixed, faster variability is translated to higher proton injection compactness (), higher magnetic field strength (), smaller emission region () but approximately constant , as the only dependance on comes through the logarithmic term (see eq. (13)). Thus, for the same but smaller the system is driven more deep to the supercritical regime (panel (a) in Fig. 7). The effect of the stronger magnetic field for  s is also imprinted on the cutoff energy of the neutrino spectrum, which moves from  PeV to 10 PeV (panel (b)). The gamma-ray emission produced for in both cases is shown in Fig. 8. The bowties show the range of observed values for the GRB photon indices, namely and (Preece et al., 2000), and are plotted for guiding the eye. Note the change of the gamma-ray spectrum below the peak, which becomes harder as .

Figure 9: Log-log plot of the ratio as a function of the Lorentz factor for parameter sets from Group A and B shown in panels (a) and (b), respectively. The subscript accounts for photons, electron and muon neutrinos, neutrons, protons, and UHE ( eV) protons. All parameters are the same as in Figs. 3 and 7.

3.2.2 Efficiency

A robust sign of the transition to supercriticality is the abrupt increase of the photon, neutron and neutrino luminosities. The efficient conversion of energy originally stored in relativistic protons into radiation can be the result of other underlying feedback loops, such as the ‘PPS-loop’ (Kirk & Mastichiadis, 1992), which was also studied in the framework of GRB prompt emission (Mastichiadis & Kazanas, 2006, 2009). Here we define the efficiency of the -th component as , where , (photons), (neutrinos), (neutrons), (protons) and (protons with energies  eV). The efficiency as a function of the Lorentz factor is shown in Fig. 9 with panels (a) and (b) corresponding to Groups A and B, respectively. Few things are worth commenting:

  • The luminosity of the proton component is the dominant one. The photon luminosity becomes comparable to the luminosity carried by the proton component only for parameters that drive the system deep into the supercritical regime (see e.g. Fig. 1), where we typically find that or .

  • The abrupt increase of , and is seen for (see also eq. (15)) and it marks the transition to supercriticality.

  • In panel (b) the abrupt increase of the neutrino, neutron and photon luminosities is not evident, as our simulations do not extend above , which in this case is .

  • The steep decrease of ( orders of magnitude for a 0.3 order of magnitude change in ) indicates the significant cooling of UHE protons.

  • For we find . For on the other hand, decreases faster than , since photopion interactions are the sole source of neutrinos contrary to photons, which are produced mainly via synchrotron radiation in the subcritical regime.

  • The efficiencies in neutrons and neutrinos are, generally, of the same order of magnitude. They show, however, different trends: remains approximately constant (panel (a)) or decreases (panel (b)) for smaller values of , whereas increases. This implies that neutrons interact with photons before escaping from the source and contribute to the neutrino production through the process .

3.2.3 Maximum proton energy

We continue our study on the emission features while the system is in the supercritical regime by considering a fiducial case where protons are not accelerated up to UHE, i.e. the maximum proton energy in the comoving frame is  eV. One could imagine a scenario where the magnetization of the burst is very low (e.g. ) or the size of the emission region is small enough to confine higher energy protons.

Figure 10: Examples of multiwavelgth spectra obtained for a lower value of the proton injection energy () in the subcritical (panel a) and supercritical (panel b) regime. Panel (a): The total photon spectrum when all processes are taken into account is plotted with thick black line, whereas when absorption is omitted the result is shown with thin black line. The contribution to the total flux of various components is also shown: proton synchrotron (dotted line), gamma-rays from decay (dashed line) and synchrotron from pairs from BH and photopion processes (dashed-dotted line). The neutrino spectra (grey line) are also overplotted. The grey colored region marks the 0.1-10 TeV energy band. Panel (b): The total photon (solid line), neutrino (dotted line) and proton (dashed line) energy spectra are shown . The typical low- and high-energy photon indices are also shown. Other parameters used are:  erg/s, ,  s, ,  G and  cm; the proton injection luminosity used in panels (a) and (b) is  erg/s and  erg/s, respectively.

To exemplify the above we adopt the following parameters: ,  erg/s, ,  s, and two values of the proton injection luminosity,  erg/s and  erg/s, which correspond to the subcritical and supercritical regime, respectively. Other parameters used are:  G and  cm.

The results for the subcritical and supercritical cases are summarized in panels (a) and (b) of Fig. 10, respectively. In the subcritical regime the various components of the overall photon spectrum can be identified, in contrast to the supercritical regime where the formation of the EM cascade blurs the emission signatures from individual processes. In panel (a) we plot the various contributions in order to demonstrate that in the TeV energy band (grey colored region) the gamma-ray emission is no longer dominated by the proton synchrotron component (for comparison see grey line in Fig. 3). Instead, it is the synchrotron radiation from secondary pairs produced through the Bethe-Heitler and photopion processes that is being emitted as Tev gamma-rays and that it is going to initiate the instability of automatic photon quenching (for relevant discussion see also Petropoulou & Mastichiadis 2012b).

Note that the emission signatures of hadronic plasmas in the subcritical regime may differ significantly for various parameter sets, as we exemplified in Figs. 3 and 10. However, the photon and neutrino emission produced when the system is driven to the supercritical regime, and in particular for parameters which ensure , is rather ‘universal’. Besides the overall energetics, the photon SEDs shown in Figs. 3 (black line), 8 and 10 have similar features.

We also found that for lower values of , a higher value is generally required to enter the supercritical regime. Given that in this regime approximetaly goes to gamma-rays and neutrinos (see also Fig. 9), the above conditions lead inevitably to bright photon and neutrino bursts. For the case shown in Fig. 10 for example, the gamma-ray fluence emitted in the supercritical regime would be erg/cm, if we were to assume a duration of  s. This would place such an event at the high-fluence tail of the Fermi distribution (von Kienlin et al., 2014).

4 Discussion

Hadronic models have served for a long time as alternatives to the more popular leptonic ones for AGN and GRB high-energy emission. One of their unique features is the prediction of copious neutrino emission which can be produced alongside the photon spectrum and is particularly attractive because it allows for hadronic models to be further tested. These suffer, however, from low efficiencies, as hadrons have in general long cooling timescales.

In earlier works (Mastichiadis et al., 2005; Dimitrakoudis et al., 2012b; Petropoulou & Mastichiadis, 2012a, b) that were targeting AGN high-energy emission, we showed that we can divide the parameter space of hadronic plasmas into two regimes. In the first one, which we shall call ‘subcritical’, protons carry the majority of the energy while a small amount is radiated away mainly by proton synchrotron radiation; photon induced processes like photopair and photopion carry even less luminosity – this is especially true if one assumes that there are no ambient photons illuminating the source. In this regime the system is inefficient, i.e. protons lose a very small part of their energy to secondaries. The other regime, the ‘supercritical’ one, is separated by a sharp boundary in phase space from the subcritical one and is characterised by exactly the opposite trend. Here, processes like photopair and photopion production dominate the losses and essentially drain the protons of their stored energy giving it to secondaries, thus increasing the efficiency to high values. Because of the radical change not only in the efficiency but also in the emission signatures of such plasmas, the transition from the sub-to the super-critical regime can be characterised as a ‘phase’ transition, with the underlying reason being the existence of non-linear feeback loops.

In the present study we re-examined the radiative signatures of hadronic plasmas in the context of GRBs. The key question we wanted to answer was whether the energetics required for the supercriticality are compatible with typical GRB parameters. For this, we adopted the usual GRB hadronic picture used in the literature, i.e. we assumed the injection of high energy protons having a power law distribution in a source of a given size that contains a certain magnetic field and adopted parameter values relevant to GRB sources. However, departing from the usual GRB assumptions, neither external photons nor an extra population of accelerated relativistic electrons were considered. This choice minimized the number of free parameters to five: , , , and ; all other quantities, such as the proton injection luminosity () and the maximum energy () of their distribution, can be expressed in terms of these parameters (see eqs. (1), (3), (5), (6), (7)). For the derivation of we assumed that proton acceleration is fast, i.e. , where is the proton gyration timescale. For example, stochastic Fermi acceleration due to magnetic turbulence (e.g. Waxman 1995a; Dermer & Humi 2001) as well as relativistic magnetic reconnection (e.g. Giannios 2010) have been suggested as viable processes for fast, UHE proton acceleration. Whether or not the accelerated protons can form a power-law distribution remains to be shown, although there are some indications for the formation of a high-energy proton power-law tail in simulations of relativistic reconnection in electron-ion plasmas (private communication with Dr. L. Sironi). Note that the same kind of numerical studies performed for pair plasmas clearly show the formation of a power-law electron distribution (Sironi & Spitkovsky, 2014).

As a tool for our study we employed a recently developed numerical code (Dimitrakoudis et al., 2012a) for solving the system of spatially averaged kinetic equations, that describes the coupling between protons and their stable by-products, namely photons, electrons (and positrons), neutrons and neutrinos. This, in contrast to most of the work performed on hadronic plasmas thus far, has allowed us to make a self-consistent study of the evolution of the system by keeping track of the energy lost and gained by the various species. We have confirmed that this coupling between the species is the key to understanding its behaviour as only with this scheme can one follow the development and growth of non-linear loops that eventually lead the system to supercriticality.

Altough there are various loops that may cause the phase transition from sub- to super-criticality (e.g. Kirk & Mastichiadis 1992; Kazanas et al. 2002; Mastichiadis et al. 2005), our analysis presented in §2 and in Appendix A shows that ray quenching (Stawarz & Kirk, 2007; Petropoulou & Mastichiadis, 2011) is the leading one for the parameters used. This means that once proton-produced rays reach a certain compactness they are spontanteously absorbed, giving rise to electron-positron pairs and radiation, which causes more proton cooling via photopair and photopion processes; and eventually more -rays, thus sustaining the loop. The cycle will continue until protons are drained of their energy and this, depending on the initial conditions, can lead the system either to steady state or to a limit cycle behaviour – in the latter case energy is gradually built into protons and is abruptly released in a few crossing times. Note that all the examples presented in §3 were obtained for parameters that led quickly (in 1-2 dynamical times) to a steady state.

As we have shown in §2, it is the existence of a critical -ray compactness that makes the proton compactness the most important parameter in our study. On the one hand, if we use a similar definition to the radiation compactness and relate the observed proton luminosity, which is assumed to be a fraction of the jet’s kinetic luminosity , to the one measured in the comoving frame of the flow through a Lorentz transformation, we find . On the other hand, we showed that the critical -ray compactness is translated to a critical value for the proton one, which roughly scales as . The combination of the and dependances favours supercriticality for most parameter values, except for flows with high bulk Lorentz factors. These results are summarized in Fig. 1, which answers in the most satisfactory manner our earlier posed key question.

Apart from the above, there are some far reaching consequences of this model when applied to GRBs. First, the efficiency of photons and neutrinos becomes quite high, reaching 0.1-0.4 of the total available proton luminosity. Moreover, the high number of electron-positron pairs created as secondaries cool to low energies producing a high Thomson optical depth (ranging from a few up to ) which can downscatter high energy photons producing a bump at   in the rest frame of the flow (see Appendix B). It is noteworthy that a similar idea for explaining the peak of the GRB emission was proposed by Brainerd (1994), although in a different context. Moreover, the boundary in the phase space between the two regimes is very sharp, in the sense that a small perturbation in one of the proton injection parameters, while the system is still in the subcritical regime, can push it over to the supercritical one. The transition is in most cases very abrupt and it manifests itself with a photon flare which lasts a few crossing times (Dimitrakoudis et al., 2012a; Petropoulou & Mastichiadis, 2012b).

Clearly there are open questions that one needs to address before the present model can successfully explain the GRB phenomenology, such as gamma-ray spectra and time-variability, which we plan to investigate in a forthcoming publication. For the former, we plan to include an additional equation for low-energy electrons () and follow the cooling/heating due to Compton process in more detail. Regarding the latter, one of the central issues that has to be addressed is whether or not our hadronic model can produce the fast variability observed in -rays. Preliminary results of photon lightcurves derived in the supercritical regime assuming a variable proton injection rate show that rapid variability (of the order of the source crossing time) can be obtained. Note that even for a constant proton injection rate there are regimes in the parameter space that are relevant to GRBs and lead to a limit cycle behaviour (Petropoulou & Mastichiadis, 2012b). Thus, in cases where the proton injection is variable we find that the structure of the gamma-ray light curves is complex due to the superposition of the intrinsic periodicity of the hadronic plasma and of the variability pattern of the ‘external’ source of proton injection.

5 Summary

We showed that the injection of high-energy protons, e.g.  eV, with luminosities  erg/s in a region which is part of a GRB-like flow with bulk Lorentz factor can lead to an abrupt energy transfer from protons to photons and neutrinos, which may carry approximately of the injected proton luminosity. We also found that the proton injection luminosity that is required for triggering the efficient cooling of protons is higher, e.g. reaching  erg/s, for proton distributions extending to lower energies, e.g. . In this case, our model has the testable prediction of a contemporaneous bright burst in -rays and high-energy neutrinos. If the neutrino non-detection from Fermi bright GRBs with IceCube will be established (e.g. Abbasi et al. 2012; He et al. 2012; Liu & Wang 2013), in the context of our model, it will mean that the conditions in sources where proton acceleration to UHE is not possible, should be such as not to drive the system to supercriticality, e.g.  erg/s.

We showed that in our framework the gamma-ray photon spectrum is self-consistently determined. Its shape is sensitive on how deep in the supercritical regime the system is driven. When the proton luminosity is marginally above the critical value, the photon spectrum peaks at  GeV, whereas for higher proton luminosities it is modified due to Comptonization and appears more as a Band spectrum. Not only different values of the proton luminosity but also small changes in one of the other parameters, and in particular of the bulk Lorentz factor, may lead to photon spectra ranging between the two shapes described above.

Summarizing, we showed that supercriticalities are a generic feature of hadronic models, as they manifest themselves for a wide range of parameters; from those relevant to AGN high-energy emission to those relevant to GRBs. In the context of the latter, such supercriticalities can lead naturally to bursts of gamma-rays that share several properties with the typical GRB prompt emission, such as the observed photon luminosity, the Band-like spectra and the time duration. These offer also a unique way of transferring energy from protons to photons and neutrinos in a very efficient way, and of overcoming the usual low-efficiency problem of hadronic models.


Support for this work was provided by NASA through Einstein Postdoctoral Fellowship grant number PF 140113 awarded by the Chandra X-ray Center, which is operated by the Smithsonian Astrophysical Observatory for NASA under contract NAS8-03060. DG acknowledges support from the Fermi 6 cycle grant number 61122.


  • Aartsen et al. (2014a) Aartsen M. G. et al., 2014a, PhysRevD, 89, 102001
  • Aartsen et al. (2014b) Aartsen M. G. et al., 2014b, ArXiv e-prints
  • Abbasi et al. (2012) Abbasi R. et al., 2012, Nature, 484, 351
  • Abbasi et al. (2010) Abbasi R. et al., 2010, ApJ, 710, 346
  • Abdo et al. (2009) Abdo A. A. et al., 2009, ApJL, 706, L138
  • Ackermann et al. (2011) Ackermann M. et al., 2011, ApJ, 729, 114
  • Arons (1971) Arons J., 1971, ApJ, 164, 437
  • Asano (2005) Asano K., 2005, ApJ, 623, 967
  • Asano et al. (2009) Asano K., Guiriec S., Mészáros P., 2009, ApJL, 705, L191
  • Asano & Inoue (2007) Asano K., Inoue S., 2007, ApJ, 671, 645
  • Asano et al. (2010) Asano K., Inoue S., Mészáros P., 2010, ApJL, 725, L121
  • Asano & Mészáros (2013) Asano K., Mészáros P., 2013, JCAP, 9, 8
  • Asano & Meszaros (2014) Asano K., Meszaros P., 2014, ArXiv e-prints
  • Atoyan & Dermer (2003) Atoyan A. M., Dermer C. D., 2003, ApJ, 586, 79
  • Baerwald et al. (2013) Baerwald P., Bustamante M., Winter W., 2013, ApJ, 768, 186
  • Baerwald et al. (2014) Baerwald P., Bustamante M., Winter W., 2014, ArXiv e-prints
  • Baerwald et al. (2011) Baerwald P., Hümmer S., Winter W., 2011, PhysRevD, 83, 067303
  • Band et al. (1993) Band D. et al., 1993, ApJ, 413, 281
  • Band et al. (2009) Band D. L. et al., 2009, ApJ, 701, 1673
  • Begelman et al. (1990) Begelman M. C., Rudak B., Sikora M., 1990, ApJ, 362, 38
  • Beloborodov (2010) Beloborodov A. M., 2010, MNRAS, 407, 1033
  • Bloom et al. (2003) Bloom J. S., Frail D. A., Kulkarni S. R., 2003, ApJ, 594, 674
  • Boettcher (2010) Boettcher M., 2010, ArXiv e-prints
  • Böttcher (2007) Böttcher M., 2007, Ap&SS, 309, 95
  • Brainerd (1994) Brainerd J. J., 1994, ApJ, 428, 21
  • Chiang & Dermer (1999) Chiang J., Dermer C. D., 1999, ApJ, 512, 699
  • Coppi & Blandford (1990) Coppi P. S., Blandford R. D., 1990, MNRAS, 245, 453
  • Crider et al. (1997) Crider A. et al., 1997, ApJL, 479, L39
  • Daigne et al. (2011) Daigne F., Bošnjak Ž., Dubus G., 2011, A&A, 526, A110
  • Derishev et al. (2001) Derishev E. V., Kocharovsky V. V., Kocharovsky V. V., 2001, A&A, 372, 1071
  • Dermer & Atoyan (2003) Dermer C. D., Atoyan A., 2003, Physical Review Letters, 91, 071102
  • Dermer & Atoyan (2006) Dermer C. D., Atoyan A., 2006, New Journal of Physics, 8, 122
  • Dermer & Humi (2001) Dermer C. D., Humi M., 2001, ApJ, 556, 479
  • Dimitrakoudis et al. (2012a) Dimitrakoudis S., Mastichiadis A., Protheroe R. J., Reimer A., 2012a, A&A, 546, A120
  • Dimitrakoudis et al. (2012b) Dimitrakoudis S., Petropoulou M., Mastichiadis A., 2012b, International Journal of Modern Physics Conference Series, 8, 19
  • Dimitrakoudis et al. (2014) Dimitrakoudis S., Petropoulou M., Mastichiadis A., 2014, Astroparticle Physics, 54, 61
  • Gao et al. (2012) Gao S., Asano K., Mészáros P., 2012, JCAP, 11, 58
  • Giannios (2006) Giannios D., 2006, A&A, 457, 763
  • Giannios (2010) Giannios D., 2010, MNRAS, 408, L46
  • Giannios (2012) Giannios D., 2012, MNRAS, 422, 3092
  • Giovanoni & Kazanas (1990) Giovanoni P. M., Kazanas D., 1990, Nature, 345, 319
  • Goldstein et al. (2012) Goldstein A. et al., 2012, ApJS, 199, 19
  • Goodman (1986) Goodman J., 1986, ApJL, 308, L47
  • Guetta et al. (2004) Guetta D., Hooper D., Alvarez-Mun~Iz J., Halzen F., Reuveni E., 2004, Astroparticle Physics, 20, 429
  • He et al. (2012) He H.-N., Liu R.-Y., Wang X.-Y., Nagataki S., Murase K., Dai Z.-G., 2012, ApJ, 752, 29
  • Hillas (1984) Hillas A. M., 1984, ARA&A, 22, 425
  • IceCube Collaboration (2013) IceCube Collaboration, 2013, Science, 342
  • Illarionov & Syunyaev (1972) Illarionov A. F., Syunyaev R. A., 1972, SvA, 16, 45
  • Kardashev (1962) Kardashev N. S., 1962, SvA, 6, 317
  • Katz (1994) Katz J. I., 1994, ApJL, 432, L107
  • Kazanas et al. (2002) Kazanas D., Georganopoulos M., Mastichiadis A., 2002, ApJL, 578, L15
  • Kirk & Mastichiadis (1989) Kirk J. G., Mastichiadis A., 1989, A&A, 213, 75
  • Kirk & Mastichiadis (1992) Kirk J. G., Mastichiadis A., 1992, Nature, 360, 135
  • Kocevski & Butler (2008) Kocevski D., Butler N., 2008, ApJ, 680, 531
  • Kompaneets (1956) Kompaneets A. S., 1956, Zh.E.F.T, 31, 876
  • Lazzati et al. (2000) Lazzati D., Ghisellini G., Celotti A., Rees M. J., 2000, ApJL, 529, L17
  • Lightman et al. (1981) Lightman A. P., Lamb D. Q., Rybicki G. B., 1981, ApJ, 248, 738
  • Lightman & Zdziarski (1987) Lightman A. P., Zdziarski A. A., 1987, ApJ, 319, 643
  • Liu & Wang (2013) Liu R.-Y., Wang X.-Y., 2013, ApJ, 766, 73
  • Lloyd & Petrosian (2000) Lloyd N. M., Petrosian V., 2000, ApJ, 543, 722
  • Mannheim et al. (2001) Mannheim K., Protheroe R. J., Rachen J. P., 2001, PhysRevD, 63, 023003
  • Mastichiadis & Kazanas (2006) Mastichiadis A., Kazanas D., 2006, ApJ, 645, 416
  • Mastichiadis & Kazanas (2009) Mastichiadis A., Kazanas D., 2009, ApJL, 694, L54
  • Mastichiadis & Kirk (1995) Mastichiadis A., Kirk J. G., 1995, A&A, 295, 613
  • Mastichiadis et al. (2005) Mastichiadis A., Protheroe R. J., Kirk J. G., 2005, A&A, 433, 765
  • Medvedev (2000) Medvedev M. V., 2000, ApJ, 540, 704
  • Mészáros & Rees (2000) Mészáros P., Rees M. J., 2000, ApJ, 530, 292
  • Mücke et al. (2000) Mücke A., Engel R., Rachen J. P., Protheroe R. J., Stanev T., 2000, Computer Physics Communications, 124, 290
  • Mücke & Protheroe (2001) Mücke A., Protheroe R. J., 2001, Astroparticle Physics, 15, 121
  • Murase (2008) Murase K., 2008, PhysRevD, 78, 101302
  • Murase et al. (2012) Murase K., Asano K., Terasawa T., Mészáros P., 2012, ApJ, 746, 164
  • Murase et al. (2008) Murase K., Ioka K., Nagataki S., Nakamura T., 2008, PhysRevD, 78, 023005
  • Nakar & Piran (2002) Nakar E., Piran T., 2002, MNRAS, 331, 40
  • Norris et al. (1996) Norris J. P., Nemiroff R. J., Bonnell J. T., Scargle J. D., Kouveliotou C., Paciesas W. S., Meegan C. A., Fishman G. J., 1996, ApJ, 459, 393
  • Panaitescu & Mészáros (2000) Panaitescu A., Mészáros P., 2000, ApJL, 544, L17
  • Pe’er et al. (2006) Pe’er A., Mészáros P., Rees M. J., 2006, ApJ, 642, 995
  • Pe’er & Zhang (2006) Pe’er A., Zhang B., 2006, ApJ, 653, 454
  • Petropoulou (2014) Petropoulou M., 2014, submitted to MNRAS
  • Petropoulou et al. (2014) Petropoulou M., Giannios D., Dimitrakoudis S., 2014, ArXiv e-prints
  • Petropoulou & Mastichiadis (2011) Petropoulou M., Mastichiadis A., 2011, A&A, 532, A11
  • Petropoulou & Mastichiadis (2012a) Petropoulou M., Mastichiadis A., 2012a, MNRAS, 426, 462
  • Petropoulou & Mastichiadis (2012b) Petropoulou M., Mastichiadis A., 2012b, MNRAS, 421, 2325
  • Piran (2004) Piran T., 2004, Reviews of Modern Physics, 76, 1143
  • Preece et al. (1998) Preece R. D., Briggs M. S., Mallozzi R. S., Pendleton G. N., Paciesas W. S., Band D. L., 1998, ApJL, 506, L23
  • Preece et al. (2000) Preece R. D., Briggs M. S., Mallozzi R. S., Pendleton G. N., Paciesas W. S., Band D. L., 2000, ApJS, 126, 19
  • Protheroe & Johnson (1996) Protheroe R. J., Johnson P. A., 1996, Astroparticle Physics, 4, 253
  • Racusin et al. (2008) Racusin J. L. et al., 2008, Nature, 455, 183
  • Reynoso (2014) Reynoso M. M., 2014, ArXiv e-prints
  • Sari et al. (1996) Sari R., Narayan R., Piran T., 1996, ApJ, 473, 204
  • Sironi & Spitkovsky (2014) Sironi L., Spitkovsky A., 2014, ApJL, 783, L21
  • Stawarz & Kirk (2007) Stawarz Ł., Kirk J. G., 2007, ApJL, 661, L17
  • Stern & Svensson (1991) Stern B., Svensson R., 1991, in Lecture Notes in Physics, Berlin Springer Verlag, Vol. 391, Relativistic Hadrons in Cosmic Compact Objects, Zdziarski A. A., Sikora M., eds., p. 41
  • Stern & Poutanen (2004) Stern B. E., Poutanen J., 2004, MNRAS, 352, L35
  • Sultana et al. (2013) Sultana J., Kazanas D., Mastichiadis A., 2013, ApJ, 779, 16
  • Tavani (1996) Tavani M., 1996, ApJ, 466, 768
  • Thompson (1994) Thompson C., 1994, MNRAS, 270, 480
  • Totani (1998) Totani T., 1998, ApJL, 509, L81
  • Vietri (1995) Vietri M., 1995, ApJ, 453, 883
  • Vietri (1997) Vietri M., 1997, Physical Review Letters, 78, 4328
  • von Kienlin et al. (2014) von Kienlin A. et al., 2014, ApJS, 211, 13
  • Wang et al. (2009) Wang X.-Y., Li Z., Dai Z.-G., Mészáros P., 2009, ApJL, 698, L98
  • Waxman (1995a) Waxman E., 1995a, Physical Review Letters, 75, 386
  • Waxman (1995b) Waxman E., 1995b, ApJL, 452, L1
  • Waxman & Bahcall (1997) Waxman E., Bahcall J., 1997, Physical Review Letters, 78, 2292
  • Zhang & Kumar (2013) Zhang B., Kumar P., 2013, Physical Review Letters, 110, 121101
  • Zhang & Mészáros (2004) Zhang B., Mészáros P., 2004, International Journal of Modern Physics A, 19, 2385

Appendix A Derivation of critical proton compactness

Here we derive the expression for the critical proton compactness () in the case where the gamma-rays that are spontaneously absorbed are the result of proton synchrotron radiation. We derive the expression for a power-law proton distribution with slope . A similar analysis can be followed for different power-law indices.

The peak of the proton synchrotron spectrum in units appears at


where ,  G, and . Here is the Lorentz factor of protons that cool within the dynamical timescale and it is given by eq. (9). The maximum Lorentz factor of protons is , where and are given by eqs. (6) and (7), respectively.

  • The condition is translated to


    In this regime, we find that if


    where we used eqs. (9) and (7). Both conditions are identical apart from a numerical factor . Thus, as long as condition (20) is satisfied we can assume that synchrotron proton cooling is not significant and set .

  • The condition is satisfied as long as


    Using eqs. (9) and (6) we find that only if


    Inspection of conditions (22) and (23) reveals that in the regime where the high-energy part of the proton distribution is affected by synchrotron losses, at least for most parameter values. For this, we set as long as condition (22) is satisfied.

The peak luminosity of the proton synchrotron spectrum in the comoving frame is then given by




is the total number of protons which is related to the proton injection rate as


where , and . The above expression is exact only if . However, we find that is still a good approximation for and .

Using eq. (26) the proton injection compactness defined by eq. (5) can be written in terms of as


Assuming that , where is the integrated -ray luminosity, and using eqs. (24) and (27), the -ray compactness is written as


As shown by Petropoulou & Mastichiadis (2011), the critical -ray compactness () is a function of the ray photon’s energy (see eq. (34) therein) having a minimum at the energy (in units)


and increasing as for . The minimum value of is found to be . For , the critical -ray compactness may be written as


In general, the peak of the proton synchrotron spectrum () appears at . The transition to supercriticality occurs if . By combining eqs. (3), (19), (28), and (30) we find equivalently that


Appendix B Peak energy - Thomson depth relation

The complete form of the Kompaneets equation (Kompaneets, 1956) is


where , , and is the differential photon number density which is related to the photon occupation number as . In most astrophysically related cases and the induced emission term () in the Kompaneets equation may be safely neglected. Here we consider cases where and , i.e. the Kompaneets equation describes the ‘recoil effect’, where a photon loses energy through multiple scatterings with cold electrons.

By adding a source and an escape term, the Kompaneets equation now reads


where is the time in units and is an approximate photon escape timescale given by


and is given by eq. (21b) in Lightman & Zdziarski (1987). For our purposes, however, it is sufficient to use . Finally, the source term appearing in eq. (33) is given by


We note that the Thomson depth is, in principal, a time-dependent quantity. Using knowledge gained from the numerical study of the problem, according to which scales as with in the supercritical regime, we can consider to be constant.

Equation (33) is solved using the method of characteristics that transforms a partial differential equation (DE) into an ordinary DE along the characteristic curves or surfaces in two- or three- dimensional problems, respectively. In our case, the characteristic curve is given by


same as in the case of synchrotron or/and inverse Compton (in the Thomoson regime) cooling of relativistic electrons (see e.g. Kardashev 1962). First, we find the solution to the homogeneous equation by setting . Along the characteristic curve of eq. (36), the PDE now reads