Hadronic origin of the TeV flare of M87 in April 2010

M87 is a giant radio galaxy with FR-I morphology. It underwent three episodes of TeV flaring in recent years with the strongest one in April 2010 which was jointly monitored by MAGIC, VERITAS and H.E.S.S. We explain its spectral energy distribution in the energy range 0.3–5 TeV by assuming that the flaring occurs in the innermost region of the jet. In this region the low energy SSC photons serve as the target for the Fermi-accelerated high energy protons of energy ≤30\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\le }30$$\end{document} TeV to form a delta resonance. The TeV photons are produced from the subsequent decay of the delta resonance to neutral pions. In this scenario the observed TeV flux of April 2010 flare is fitted very well.


Introduction
M87 is a giant radio galaxy in the Virgo cluster at a luminosity distance of 16.7±0.2 Mpc [1] and a redshift of z = 0.00436. The mass of the central supermassive black hole (SMBH) is estimated to be M BH = (3-6) × 10 9 M [2]. Based on its radio morphology it is classified as FR-I galaxy [3]. The radio images and modeling of its interaction with the surrounding environment suggests that the jet is misaligned with respect to the line of sight [3,4]. The substructures in the plasma jet originated from the center of M87 is resolved at different wavelengths (radio [5], optical [3] and X-ray [6]). Due to the harbouring of SMBH in the center and the presence of the jet, M87 was considered as a potential candidate for TeVemission. The evidence for very high energy (VHE) γ -rays (E γ > 100 GeV) emission from M87 was reported by the HEGRA Collaboration in 2003 [7] and was later confirmed by H.E.S.S., VERITAS [8,9] and MAGIC. The AGN M87 is normally a weak VHE source, but it shows strong variability at VHE with time scale of the order of days, which indicates a compact emission region < 5 × 10 15 D cm (where D is the a e-mail: sarira@nucleares.unam.mx Doppler factor of the emitting plasma), corresponding to only a few Schwarzschild radii (R s = 2GM BH /c 2 10 15 cm).
So far, there are three episodes of enhanced VHE γ -ray emission observed from the AGN M87 in the years 2005 [10,11], 2008 [8] and 2010 [8,12]. The latest one of April 2010, is the strongest TeV γ -ray flare ever detected from the AGN M87 with a peak flux of (2.7 ± 0.68) × 10 −11 erg cm −2 s −1 for E γ > 350 GeV [8,9,13]. The detected single isolated flare is well described by two sided exponential functions with significantly different flux rise and decay times [8]. The rising (5-8 of April), the peak (9-10 of April), and the falling (11)(12)(13)(14)(15) of April) parts of the flux during the flare are consistent with power-law behavior. This flare was detected simultaneously by VERITAS, MAGIC, and H.E.S.S. [9,12], and triggered further multi-wavelength observations in the radio, optical, and X-ray ranges. This was also observed by Fermi-LAT at MeV-GeV energies but it could not observe day-scale variability [13].
Different theoretical models have been proposed to explain the flaring in M87. Wagner et al. [14] have complied the multi-wavelength data sets spanning almost all the energy range and presented a spectral energy distribution (SED) of M87 along with leptonic and hadronic models predictions. The hadronic synchrotron-proton blazar model [15] suggests the emission of synchrotron photons from protons, charged pions and muons in the jet magnetic field. However, the SED produced using the archival data before 2004 shows a steep drop-off at TeV energies and to explain above TeV energy a strong Doppler boosting is needed which is not the case in M87. So this model is not compatible with any of the VHE spectral measurements after 2004. The leptonic decelerating inner jet model by Georganopoulos et al. [16] does not describe the hard TeV spectra well as it has a strong cut-off. The multi-blob synchrotron self Compton (SSC) model by Lenain et al. [7] needs a low magnetic field in the VHE emitting region, which is unlikely because this region is of the order of the Schwarzschild radius and is expected to have a strong field. Thus the so called one-zone homogeneous leptonic models of Georganopoulos et al. [16] and Lenain et al. [7] are very unlikely to reproduce the observed VHE spectrum. The lepto-hadronic model [17] fits to the low energy γ -ray spectrum by Fermi/LAT and HESS low state but not the flaring state. However, this model has been extended by Reynoso et al. [18] to study the production of very high energy γ -rays in blazars by introducing a two-components jet. The spine-sheath model by Tavecchio and Ghiselline [19] has difficulties to achieve a harder spectrum in the VHE range due to strong absorption of the TeV photons from interactions with the optical-infra red (IR) photons from the spine. In the jet-in-jet model of Giannios et al. [20] minijets are formed within the jet due to flow instabilities and these minijets move relativistically with respect to the main jet flow. The interaction of the daughter jets with the main jet are responsible for the production of VHE gamma rays. While the minijets are aligned with our line of sight, the VHE gamma rays are beamed with a large Doppler factor. This scenario can explain the 2010 flare but does not provide a quantitative prediction of the light curve of the flare. Similarly, the magnetosphere model [21][22][23] can explain the hard TeV spectrum but in this case also there is no detailed quantitative predication for the VHE light curve. The work by Cui et al. [24] explains the VHE gamma ray flare in an external inverse Compton model with a very wide jet. In another picture, the dense and cold cloud called broad line region (BLR) located in a small region surrounding the central black hole may interact with the hadronic component of the AGN jet and produce γ -rays as well as neutrinos [25]. Barkov et al. [26] have proposed a scenario where a red giant star with a loosely bound envelope of mass ∼10 29 g interacts with the base of the M87 jet. The VHE emission is produced near the SMBH due to the interaction of the cosmic ray protons emerging from the jet with the disrupted dense cloud of the red giant through proton-proton interaction. This model reproduced well the light curve and the energy spectrum of the April 2010 flare. An alternative mechanism is proposed to explain the TeV flaring of the blazar 1ES 1959+650 [27] and the multi-TeV emission from Centaurus A [28]. In this mechanism, it is assumed that the flaring occurs within a compact and confined region which lies inside a blob of radius R b . In this region the photon density is very high compared to the rest of the jet. The Fermi-accelerated high energy protons collide with the high density target photons to produce -resonance. These high density photons have the same energies as the low energy tail of the SSC photons (1ES 1959+650) or the SSC peak photons (Centaurus A). The -resonance subsequently decays to high energy photons and neutrinos. Here we further assume that the ratio of photon densities at different energies in the flaring and non-flaring states are the same. This scenario neither needs any intervening foreign object nor any special jet cloud geometry [29] to produce the high energy photons.
The plan of the paper is as follows: In Sect. 2 we review in detail the flaring model and the kinematical condition for the production of -resonance. The results are discussed in Sect. 3 and we briefly conclude in Sect. 4.

The flaring model
In a recent paper Sahu et al. [27] have explained the orphan TeV flare of 4th June, 2002, from the blazar 1ES1959+650 through hadronic model. In this work the standard interpretation of the leptonic model is used to explain both, low and high energy peaks, by synchrotron and SSC photons, respectively, as in the case of any other AGNs and blazars. Thereafter, it is proposed that the low energy tail of the SSC photons in the blazar jet serves as the target for the Fermiaccelerated high energy protons of energy ≤100 TeV, within the jet to produce TeV photons through the decay of neutral pions from the -resonance. This model explains very nicely the observed TeV flux from the orphan flare. Also it is interesting to note that this scenario is self sufficient and does not need any external medium for the production of gamma rays. As discussed, the flaring occurs within a compact and confined volume of radius R f inside the blob of radius R b (R f < R b ) which is shown in Figure 1 of Ref. [27]. Both the internal and the external jets are moving with the same bulk Lorentz factor and the Doppler factor D. In normal situation within the jet, we consider the injected spectrum of the Fermi-accelerated charged particles having a power-law spectrum dN /dE ∝ E −α with the power index α ≥ 2. But in the flaring region we assume that the Fermiaccelerated charged particles follow a power-law spectrum with an exponential cut-off [27,30] and is given as where the high energy proton has the cut-off energy E p,c and again the spectral index has the restriction α > 2. Also probably due to the copious annihilation of electron positron pairs, splitting of photons in the magnetic field and enhance IC photons around the base of the jet [20,31], the comoving photon density n γ,f (flaring) in the flaring region is much higher than the rest of the blob n γ (non-flaring) i.e. n γ,f ( γ ) n γ ( γ ). Here we assume that the ratio of photon densities at two different background energies γ 1 and γ 2 in flaring and nonflaring states remains the same, that is, In general, in the leptonic one-zone synchrotron and SSC jet model the emitting region is a blob with comoving radius R b moving with a velocity β c corresponding to a bulk Lorentz factor and seen at an angle θ ob by an observer which results with a Doppler factor D = [ (1 − β c cos θ ob )] −1 . The emitting region is filled with an isotropic electron population and a randomly oriented magnetic field B . The electrons have a power-law spectrum. The energy spectrum of the Fermiaccelerated protons in the blazar jet is also assumed to be of power-law. Due to high radiative losses, electron acceleration is limited, however, protons and heavy nuclei can reach UHE through the same acceleration mechanism.
Here, we assume that, due to photohadronic interaction in the jet, the pions are produced through the intermediate -resonance given by which has a cross section σ ∼ 5×10 −28 cm 2 . Subsequently, the neutral and charged pions will decay through π 0 → γ γ and π + → e + ν e ν μνμ , respectively. The produced neutrinos and photons are in the GeV-TeV range energy. For the production of the -resonance, the kinematical condition is where E p and ε γ are, respectively, the proton and the background photon energies in the comoving frame of the jet. We define the quantities with a prime in the comoving frame and without prime in the observer frame. For high energy protons we assume β p 1. Since in the comoving frame the protons collide with the SSC photons from all directions, in our calculation we consider an average value (1 − cos θ) ∼ 1 (θ in the range of 0 and π ). Going from comoving frame to observer frame, the proton and photons energies can be written as, respectively, and the kinematical condition given in Eq. (4) can be written in the observer frame as In the jet comoving frame, each pion carries ∼0.2 of the proton energy while 50 % of the π 0 energy will be given to each γ . So the relationship between high energy γ -ray and the E p is E γ D E p /10. From these relations we can express the -resonance kinematical condition in terms of photon energies (target photon energy ε γ and the observed photon energy E γ ) as The optical depth to produce the -resonance is given as The comoving photon number density within the confined volume can be given in terms of the luminosity L γ as with κ ∼ (0-1) (depending on whether the jet is continuous or discrete) and η ∼ 1. Here in this work we consider κ = 0. For κ = 1, the photon density will be reduced by a of factor D −1 in the discrete jet as compared to continuous one. The relationship between observed γ -ray flux F γ , high energy proton flux, and the background SSC photon density in the flaring region is given as [27] Then the observed high energy γ -ray flux at two different energies will scale as where E γ 1,2 are two different γ -ray energies and correspondingly the proton energies are E p 1,2 . In Eq. (12) we have used the relations in Eq. (2) and For a self consistent treatment, in principle we should use the photon density n γ, f in the hidden internal jet, rather than using the photon density n γ in the external jet. Unfortunately the only qualitative knowledge we have as regards the hidden internal jet is that the photon density, the energy density in photon, and the magnetic field are higher than the outer jet so that pγ interaction process will take place to produce multi-TeV γ -rays and neutrinos. For this reason we assume the scaling behavior of the photon densities in different background energies as shown in Eq. (2). This approximation has been used to interpret the orphan flaring of 1ES 1959+650 which explains very well the observation [27]. By taking into account an observed flux at a given energy, we can calculate the fluxes at other energies by using Eq. (12). A better approach is to solve self consistently the coupled transport equations for leptons and photons along the jet by taking into account their respective cooling mechanisms as well as the injection spectrum of the primary particles [18]. To avoid this complication we assume the scaling behavior of photon number densities in different regions as discussed above and then use the one-zone leptonic model fit to the SED of M87 for the calculation of the flaring events. By doing so we could simplify the problem to a great extent.
The optical depth τ pγ implies that out of τ −1 pγ many protons, one will interact with the SSC background photons to produce -resonance. In this case the fluxes of the TeV photons and the Fermi-accelerated high energy protons F p , are related through Like photons, the proton fluxes at different energies E p 1 and E p 2 , scale as By using the above relation we can calculate the proton fluxes at different energies.
The Fermi-accelerated protons could also interact inelastically with the non-relativistic protons present in the jet and produce γ -rays and neutrinos. But the efficiency of the pp process depends on the density of the background protons. In this process, to extract efficiently energy from the accelerated protons, the mean energy loss time scale in the comoving frame t pp = (K p σ pp n p ) −1 has to be shorter than the dynamical time scale t dyn = R f /c. This imposes a lower bound on the cold proton density n p (σ pp K p t dyn ) −1 . For R f ∼ 5 × 10 15 cm, we obtain n p 10 10 cm −3 . On top of this, the mass loaded jet also needs greater kinetic energy in the blob. Thus, from the energetic arguments, it is shown that high proton densities are unlikely in the γ -ray emission region of the blazars [32,33], so the pp process can be neglected.

Results
With the homogeneous leptonic one-zone synchrotron and SSC jet model [34] the SED is fitted assuming the viewing angle 10 • and bulk Lorentz factor = 2.3. This corresponds to a Doppler factor D = 3.9, which is shown in Figure 4 of Ref. [10]. Based on the multi-band correlations detected in 2005 and 2008 flaring events of M87, the core and the HST-1 are favored as the emitting regions. But during the VHE flare of 2008 and 2010, Chandra detected an enhanced X-ray flux from the core region which are the two highest measurements since the start of its observation in 2002. During this time HST-1 remained in a low state [6]. During the 2005 VHE flaring episode no enhanced X-ray emission from the core was detected. On the other hand, at that time, HST-1 was more than 30 times brighter than the core region in X-rays Table 1 These parameters (up to B ) are taken from the one-zone synchrotron model of Ref. [10] which are used to fit the SED of M87. The last three parameters are obtained from the best fit to the observed flare data in our model leading to uncertainty in the flux estimation of the core [35]. The coincidence in X-ray and VHE emission as well as the observed timescales of short variability (∼ day) at VHE/Xray suggests that the size of the emitting region is compact and leads one to believe that the 2010 VHE flare probably originates in the innermost region of the jet. Hence, here we assume that the flaring occurs within the confined volume of radius R s < R f < R b , which is in the core region. To fit the multi-wavelength SED of M87, the blob radius is taken to be R b = 1.4 × 10 16 cm and the magnetic field is B = 55 mG [10]. The above blob radius is consistent with the few days timescale variability in TeV energy. Due to our ignorance about the flaring region we assume the scaling behavior shown in Eq. (2) and our results depend on blob radius R b but not on the flaring radius R f . Although we do not know exactly how the magnetic field behaves as a function of distance r in the jet, a radial dependence of the form B (r ) = B 0 (r 0 /r ) n can be assumed, where n = 1 or 2 [17] and r 0 is few times the Schwarzschild radius. It is reasonable to take B 0 ∼ 6 G so that the magnetic field in the inner jet is ∼0.5 G, which is stronger than the outer region. In this work we use the SED and parameters of the one-zone synchrotron model given in Ref. [10]. In Table 1 we have summarized these parameters. During the flaring in April 2010, the high energy γ -ray flux was observed in the energy range 0.3 TeV ≤ E γ ≤ 5 TeV. Also the flaring spectra had distinct rise time and fall time. The rising, the peak, and the falling of the flux are fitted with power-law with different flux normalization and the spectral index α [8]. In the hadronic model, the above E γ range corresponds to the proton energy in the range 1.9 TeV ≤ E p ≤ 30 TeV. Protons in this energy range will collide with the background photons in the energy range 1.5 MeV(3.7 × 10 20 Hz) ≥ ε γ ≥ 0.1 MeV(2.3 × 10 19 Hz) to produce the -resonance and the subsequent decay of it   Fig. 2 The rising, the peak, and the falling parts of the TeV flare are fitted with the power-law exponential SED. For all these we use the same spectral index α = 2.83 and E c 12 TeV. The rising part is fitted with two different normalized flux will produce both γ -rays and neutrinos through neutral and charged pion decay, respectively. We can observe that the above energies ε γ s lie in the rising part of the SSC photons, shown as a shaded region in Fig. 1. The number den- As discussed in Ref. [27], for the calculation of the TeV flux, first we take into account one of the observed flaring fluxes with its corresponding energy for normalization e.g. F γ (E γ 2 = 3.18 TeV) 3.8 × 10 −12 TeV cm −2 s −1 and n γ (ε γ 2 = 0.15MeV) 387 cm −3 . Using it, we calculate the flux for other energies with the help of Eq. (12). We have done this for different observed fluxes for a better fit. The spectral index α and the cut-off energy E c are the free parameters in the model and the best fit is obtained for α = 2.83 and E c 12 TeV. The γ -ray cut-off energy of ∼12 TeV corresponds to E p,c 71 TeV and above the cut-off energy the flux decreases rapidly which is clearly shown in Fig. 1. With the same α = 2.83 and E c 12 TeV but different normalized fluxes we fitted the rising, the peak, and the falling fluxes which are shown in Fig. 2. The rising flux is fitted with two different normalization to have a better picture of it. In Fig. 3, the fitting of the peak flux in our model is compared with the jet cloud interaction hadronic model [26]. It shows that in the high energy limit (above ∼5 TeV), the flux in our model falls faster than the jet cloud interaction scenario while in the low energy limit (below ∼1 TeV) it is the opposite. It is interesting to note that the spectral indices α fitted to the TeV flaring SEDs of M87 and the blazar 1ES 1959+650 have the same 1e+18 1e+19 1e+20 1e+21 Fig. 4 The ratio of photon densities n γ (ε γ1 )/n γ (ε γ2 ) for a given value n γ (ε γ2 = 3.7 × 10 19 Hz) is plotted as a function of SSC photon energy. The points are the density ratios for observed data points. The curve is fitted with a straight line value, 2.83, which probably hints to a common mechanism of particle acceleration during the flaring [27].
We have plotted the ratio n γ (ε γ 1 )/n γ (ε γ 2 ) of Eq. (12 ) for a given value n γ (ε γ 2 = 3.7 × 10 19 Hz) 386 cm −3 in Fig. 4 as a function of SSC photon energy. It shows that the density ratio is almost a linear function of energy. We have specifically chosen the energy range in the vicinity of the shaded region of Fig. 1, which is responsible for the TeV spectra. The standard power-law fitting with an exponential cut-off to the SED [30] is expressed as where F 0 is a constant. But here F 0 is replaced by energy dependent photon density of the background and due to this energy dependent coefficient, fitting to SED in this model is different from the standard one. During the flaring period, not only protons but also electrons are Fermi accelerated in the inner jet with the same energy as the protons. The e + produced during the π + decay has energy in the range 0.095 TeV ≤ E γ ≤ 1.5 TeV. These electrons and positrons will produce synchrotron radiation in the jet magnetic field. While the Fermi-accelerated electrons will emit synchrotron photons in the frequency band 2.5 × 10 18 Hz ≤ ε γ ≤ 6.3 × 10 20 Hz, the positrons will radiate in the frequency band 6.3 × 10 15 Hz ≤ ε γ ≤ 1.6 × 10 18 Hz. So the flaring in the TeV energy should be accompanied by a simultaneous enhanced synchrotron emission in the frequency band 6.3 × 10 15 Hz ≤ ε γ ≤ 6.3 × 10 20 Hz. In principle, the interaction of TeV γ -rays with the extragalactic background photons can produce electron-positron pairs and reduce the multi-TeV flux from M87. But for the energy range 0.3 T eV ≤ E γ ≤ 5 T eV , the optical depth for the pair production process is almost constant, and the spectral shape remains nearly unchanged.

Conclusions
The strongest TeV flaring of the radio galaxy M87 in April 2010 is explained by assuming it to be due to the photohadronic interaction of the Fermi-accelerated protons of energy ≤30 TeV with the SSC photons in the energy range ∼ (0.1-1.5) MeV. In this scenario the proton spectrum is a power-law with an exponential cut-off. For the fitting of the rising, the peak, and the falling parts of the TeV flare we use the same spectral index α = 2.83 and the γ -ray cut-off energy E c 12 TeV but different normalization. Our results fit very well to these distinct phases of the flare.