Abstract
The Universe is filled with a diffuse background of MeV gamma-rays and PeV neutrinos, whose origins are unknown. Here, we propose a scenario that can account for both backgrounds simultaneously. Low-luminosity active galactic nuclei have hot accretion flows where thermal electrons naturally emit soft gamma rays via Comptonization of their synchrotron photons. Protons there can be accelerated via turbulence or reconnection, producing high-energy neutrinos via hadronic interactions. We demonstrate that our model can reproduce the gamma-ray and neutrino data. Combined with a contribution by hot coronae in luminous active galactic nuclei, these accretion flows can explain the keV – MeV photon and TeV – PeV neutrino backgrounds. This scenario can account for the MeV background without non-thermal electrons, suggesting a higher transition energy from the thermal to nonthermal Universe than expected. Our model is consistent with X-ray data of nearby objects, and testable by future MeV gamma-ray and high-energy neutrino detectors.
Similar content being viewed by others
Introduction
The Universe is filled with high-energy radiation including X-rays1, MeV–TeV gamma-rays2,3, and TeV–PeV neutrinos4. It is widely believed that the cosmic X-ray background (CXB) is predominantly produced by radio-quiet active galactic nuclei (AGN) and star-forming galaxies5,6, and that radio-loud AGN including blazars and star-forming galaxies provide the bulk of the GeV–TeV gamma-ray backgrounds7,8. However, the origins of the soft gamma-ray and high-energy neutrino backgrounds remain unknown9,10. Type-Ia supernovae emit MeV gamma-rays through nuclear decay line emission11,12. However, the event rate estimated by recent observations revealed that the expected signal is below the measured background13,14. Kilonovae, transients powered by r-process nuclei ejected by neutron star mergers, also emit MeV gamma-rays through nuclear decay line emission, but the predicted flux is not enough to explain the data15. Radio-loud AGN16 including blazars17 are possible candidate sources, where nonthermal electrons accelerated in their jets emit MeV gamma-rays. However, recent population studies using Swift Burst Alert Telescope (BAT) data suggest that they can contribute to only a part of the measured MeV background18,19,20. Alternatively, Refs. 21,22 suggested that nonthermal electrons accelerated in the coronae of radio-quiet AGN can account for the MeV gamma-ray background. However, the existence of nonthermal electrons is in question, because electrons are rapidly thermalized by Coulomb collisions in typical coronae. On the other hand, proton acceleration by magnetic reconnection and plasma turbulence may occur due to their slower Coulomb relaxation, from which hadronic gamma-rays are unavoidable through proton-induced electromagnetic cascades inside the magnetized coronae23. In this case, radio-quiet AGN can explain at least ~ 20% of the observed MeV gamma-ray background, if AGN coronae account for the neutrino data in the 10–100 TeV range.
Here, we propose a scenario that naturally explains the soft gamma-ray background in the MeV range without relying on nonthermal mechanisms. Hot accretion flows, or radiatively inefficient accretion flows (RIAFs)24,25, are widely expected in low-luminosity AGN (LLAGN), where the feature of an optically thick accretion disk, namely the big-blue bump, is not seen in their optical/UV spectra26. Synchrotron radiation from hot thermal electrons explains the observed infrared and radio emission, and the associated soft gamma-rays from Comptonization naturally make a significant contribution to the gamma-ray background up to ~3 MeV. In addition, protons in RIAFs may be accelerated and efficiently emit neutrinos via hadronic interactions27, which implies that our RIAF model can account for the MeV gamma-ray and neutrino backgrounds simultaneously. A combination of nonthermal protons and thermal electrons in RIAFs is naturally expected. Protons are collisionless in the sense that the relaxation timescale is longer than the dissipation timescale, whereas electrons are rapidly thermalized within the dissipation timescale via synchrotron-self absorption and Coulomb collisions28.
Results
Guaranteed soft gamma rays from RIAFs
First, we describe the properties of the thermal plasma in RIAFs (see subsection Emission from thermal electrons in RIAFs in Methods and Ref. 29 for technical details). Thermal electrons in RIAFs emit broadband photons through synchrotron radiation, bremsstrahlung, and inverse Compton scattering. For the parameter set shown in Table 1, values of the magnetic field, B, the Thomson optical depth, τT = npσTR (np is number density, σT is Thomson cross section, and R is the size of accreting plasma), and the normalized electron temperature, Θe = kBTe/(mec2) (me is electron mass, kB is Boltzmann constant, and c is speed of light), are tabulated in Table 2. RIAFs are optically thick to synchrotron self-absorption at the synchrotron characteristic energy, \({\varepsilon }_{{{{{{{{{{\rm{syn,ch}}}}}}}}}}}\approx 3{h}_{p}eB{{{\Theta }}}_{e}^{2}/(4\pi {m}_{e}c)\), where hp is Planck constant and e is the elementary charge. The absorption energy above which RIAFs become optically thin is determined by equating the synchrotron emissivity to the blackbody radiation: \(\varepsilon_{\rm syn,ab}\approx 3x_{M}h_{p}{eB{\Theta}}_{e}^{2}/(4\pi m_{e}c )\simeq 7.0 \times 10^{-3} ( B/0.18\,{{{\rm{kG}}}} ) ({\Theta}_{e}/1.5)^2 (x_{M}/10^3)\,{{{\rm{eV}}}}\), where xM = εsyn,ab/εsyn,ch ~ 103 represents the correction from εsyn,ch to εsyn,ab. See Ref. 30 for the values and the method to estimate xM. The synchrotron luminosity from an optically thick source is given by
where M is the mass of the supermassive black hole (SMBH), \({{{{{\mathcal{R}}}}}}=R/R_{S}\), and \({R}_{S}=2GM/c^{2}\) is the Schwarzschild radius. The thermal electrons up-scatter the synchrotron photons, and the resulting spectrum can be approximated by a power-law form, \({\varepsilon }_{\gamma }{L}_{{\varepsilon }_{\gamma }}\propto {\varepsilon }_{\gamma }^{1-{\alpha }_{{{{{{{{{{\rm{IC}}}}}}}}}}}}\), where \({\alpha }_{{{{{{{{{{\rm{IC}}}}}}}}}}}=-{{{{{{{{\mathrm{ln}}}}}}}}}\,{\tau }_{T}/{{{{{{{{\mathrm{ln}}}}}}}}}\,{A}_{{{{{{{{{{\rm{IC}}}}}}}}}}}\) and \({A}_{{{{{{{{{{\rm{IC}}}}}}}}}}}=4{{{\Theta }}}_{e}+16{{{\Theta }}}_{e}^{2}\)31. If αIC < 1, the Comptonization dominates over other cooling processes, which is satisfied for \(\dot{m}\,\gtrsim\, 2\times 1{0}^{-3}\) as seen in Table 2. In this case, the luminosity of the Comptonized photons is
Owing to the strong Θe dependence, our model predicts RIAF electron temperatures which lie in a narrow range Θe ~ 1 − 3, unavoidably leading to a peak in the MeV range for \(\dot{m}\,\gtrsim\, 2\times 1{0}^{-3}\). The normalization is determined by the balance with the Coulomb heating:
where \(\dot{m}\) is the normalized accretion rate, ηrad,sd is radiation efficiency for a standard disk, and α is the viscous parameter (See subsection Emission from thermal electrons in RIAFs in Methods for the Coulomb heating rate and definition of some parameters).
The broadband photon spectra from thermal electrons are shown in Fig. 1 for various values of \(\dot{m}\). The synchrotron emission produces a peak at 0.001–0.01 eV depending on \(\dot{m}\), and the Comptonization of the synchrotron photons creates higher-energy photons up to 1–10 MeV. Cases with higher \(\dot{m}\) have harder spectra because of their higher Thomson optical depths (see also Refs. 32,33), making their spectral peaks in the MeV range. These features are quantitatively consistent with the analytic estimates in Equations (1) and (3).
We use the parameter set for model A (reference model) with M = 108M⊙ and various \(\dot{m}\). The solid, long-dashed, short-dashed, dotted, and dotted-dashed lines are for \(\dot{m}=0.03,\,1.1\times 1{0}^{-2},\,3.7\times 1{0}^{-3},\,1.3\times 1{0}^{-3},\,4.6\times 1{0}^{-4}\), respectively. The photons of energies below the vertical dotted line are mainly emitted by the synchrotron process, while the photons above the energy are produced by the Comptonization process.
Our model is consistent with observations of nearby LLAGN. Ref. 34 reported a softening feature in the hard X-ray band in NGC 3998, from which they claimed that the electron temperature is ≃ 30–40 keV. Our RIAF model can reproduce the softening feature in the NuSTAR band, as well as the Swift BAT data shown in Fig. 2, despite a higher electron temperature. Ref. 34 also provided the X-ray spectrum for NGC 4579, which has a higher \(\dot{m}\) and does not show any softening feature. Our model also produces a hard power-law spectrum consistent with the NuSTAR data (see Fig. 2). In our RIAF model, the resulting spectra for NGC 3998 and NGC 4579 are relatively hard, and well below the longer wavelength data (radio, infrared/optical/ultraviolet, and soft X-rays). These should be attributed to other emission components, such as compact jets or outer accretion disks35. Indeed, radio jets are observed in both objects36,37.
The data by XMM-Newton & NuSTAR (orange regions; with a systematic error of 10%), Swift BAT (pink regions with 90% confidence levels), and Fermi LAT (downward arrows; upper limits with 95% confidence levels) are obtained from Ref. 34, Ref. 112, and Ref. 40, respectively. a Spectra for photons from thermal electrons (dashed lines), nonthermal protons (dotted-dashed), total neutrinos (thick-solid), pγ neutrinos (thin-solid), and photons by electromagnetic cascades (thick dotted) for NGC 3998. We use M = 8.1 × 108M⊙113, \(\dot{m}=2.1\times 1{0}^{-3}\), and DL = 14.1 Mpc. The thin-dotted line is the sensitivity curve of e-ASTROGAM with 1-yr integration41. b Same as (a), but for NGC 4579. We use M = 7.2 × 107M⊙101, \(\dot{m}=8.0\times 1{0}^{-3}\), and DL = 16.4 Mpc. The NuSTAR data is not smoothly connected to the BAT data, and given the huge statistical error bars in the BAT data, we ignore the BAT data.
The Event Horizon Telescope Collaboration (EHT) reported a horizon-scale image of the SMBH in M8738. Its brightness temperature is ~5 × 109 K (equivalent to Θe ~ 0.8), while the real temperature should be Θe ≃ 3–10, because the image is beam-smeared and the RIAF is likely optically thin at the observed frequency. Our model predicts Θe ≃ 3.5 with the parameters appropriate for M87 \((M=6.3\times 1{0}^{9}{{{{{{{{{{\rm{M}}}}}}}}}}}_{\odot },\,\dot{m}=6.1\times 1{0}^{-4})\), which matches the expected temperature. The electron temperature in RIAFs also affects the interpretation of the photon ring observed by EHT38,39. The emission region of the photon ring is determined by the electron temperature and magnetization, which should be clarified through the future multi-wavelength modeling of nearby LLAGN.
Nonthermal particles in RIAFs
Protons in RIAFs are accelerated by magnetohydrodynamic (MHD) turbulence and/or magnetic reconnection generated by the magnetorotational instability (MRI). Here, we focus on the stochastic proton acceleration mechanism, in which nonthermal particles randomly gain or lose their energy via interactions with turbulent MHD waves. The accelerated protons, or cosmic-rays (CRs), produce neutrinos and gamma-rays via hadronuclear and photohadronic interactions with thermal protons and photons inside the RIAFs. Neutrinos freely escape from the system, whereas the gamma-rays create electron/positron pairs via γγ → e+e−, which initiates proton-induced electromagnetic cascades. Figure 2 shows the resulting proton, neutrino, and proton-induced cascade gamma-ray spectra for NGC 3998 and NGC4579. The proton spectrum is hard because of the stochastic acceleration and has a cutoff around 10–100 PeV due to photohadronic interactions. The neutrinos are mainly produced by pp interactions for εν ≲ 105 − 106 GeV, where \({\varepsilon}_{\nu}\) is the neutrino energy, while pγ interactions are more efficient around the cutoff energy (See subsection Nonthermal particles in RIAFs in Methods for details). The resulting cascade gamma-ray spectrum is flat for εγ < εγγ, where \({\varepsilon}_{\gamma}\) is the gamma-ray energy and \({\varepsilon}_{\gamma\gamma}\) is the energy above which gamma-rays are efficiently attenuated. The gamma-ray flux decreases rapidly above the energy. This feature is commonly seen in photon spectra from well-developed electromagnetic cascades23. The cascade gamma-ray spectra are well below the upper limit from the Fermi Large Area Telescope (LAT)40 and the design sensitivity of future MeV satellites41,42,43,44.
Cosmic gamma-ray and neutrino backgrounds
We calculate the gamma-ray and neutrino background intensities from the RIAFs in LLAGN (see subsection Cumulative background intensities in Methods for details). Fig. 3 shows the resulting gamma-ray and neutrino intensities from LLAGN. In our RIAF model, the Comptonized emission from the thermal electrons naturally accounts for the soft gamma-ray background in the 0.3–3 MeV range, below which canonical AGN coronae explain the X-ray background. Furthermore, nonthermal protons in RIAFs can simultaneously reproduce the IceCube data above ~0.3 PeV, below which hot coronae in luminous AGN can account for the neutrino data in the 10–100 TeV range. Hence, our synthesized AGN core scenario provides an attractive unified explanation for a wide energy range of the cosmic keV–MeV photon and TeV–PeV neutrino backgrounds, as demonstrated in Fig. 3. Notably, AGN accretion flows can do this using a reasonable set of plasma parameters of the nonthermal proton component: β ~ 1–10, ηtur ~ 10–20, and PCR/Pth ~ 0.0123 (see Table 2). This value of PCR/Pth is reasonable in the sense that the CR energy density is lower than the magnetic field energy density. Also, the parameters related to nonthermal proton production, such as ηtur and PCR/Pth, are expected to be similar in both AGN coronae and RIAFs, because they share the same CR acceleration and turbulence generation mechanisms. In this sense, the parameters of the nonthermal particles in our RIAF model are effectively calibrated by the 10–100 TeV neutrino data and the AGN corona model.
Data points are provided by Swift BAT115 (brown-circle; 1-σ errors for systematic and statistical errors), SMM116 (brown-triangle; definition of errors are not available), Nagoya-baloon117 (purple-star; definition of errors are not available), COMPTEL2 (purple-square; 2-σ error bars for the linear sum of the systematic and statistical errors), Fermi LAT3 (black-plus; 1-σ errors with systematic and statistical uncertainties), and IceCube4,118 (lightgrey and magenta regions; 68% confidence levels). The red, blue, and green lines indicate background intensities for photons by thermal electrons, neutrinos by CR protons, and gamma-rays by proton-induced electromagnetic cascades, respectively. The thick, thin-dotted, and thin-solid lines are contributions from RIAFs in LLAGN (this work), coronae in luminous AGN23, and sum of these, respectively. Thick-dashed and thick-dotted lines are drawn with the luminosity functions of Ref. 45 and Ref. 53, respectively, and the regions between the two lines are filled with the corresponding colors. The yellow vertical bands represent the energy bands where LLAGN provide the dominant contribution. The Fermi data should be reproduced by another source class, such as radio-loud AGN7,52,114.
As shown in panel (a) of Fig. 4, it is the relatively more luminous LLAGN that mainly contribute to the MeV intensity. In particular, LLAGN with LHα ~ Lcrit provide the dominant contribution, where Lcrit is the Hα luminosity for \(\dot{m}={\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}\), where \({\dot{m}}_{{{{{\rm{crit}}}}}}\) is the critical mass accretion rate above which RIAFs no longer exist (See subsection Emission from thermal electrons in RIAFs in Methods for the value of \({\dot{m}}_{{{{{\rm{crit}}}}}}\)). Thus, we can analytically estimate the gamma-ray background intensity to be23
where H0 ~ 70 km s−1 Mpc−3 is the Hubble constant, L*, ρ*, and s1 are the break luminosity, break density, and power-law index for the luminosity function45 (see subsection Cumulative background intensities in Methods), respectively, and ξz is the redshift evolution factor of the luminosity density. This is consistent with our numerical results and the observed MeV background from COMPTEL. We stress that the resulting MeV gamma-ray intensity does not depend strongly on the parameters related to the RIAF, such as α, B, M, \({\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}\), nor on the electron heating prescription, as long as we can use the luminosity function and bolometric correlation. This is because LLAGN of \(\dot{m} \sim {\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}\) provide the dominant contribution, and the energy budget of such LLAGN does not strongly depend on the value of \({\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}\).
On the other hand, the relatively faint LLAGN mainly contribute to the neutrino background (see panel (b) of Fig. 4). The high-energy gamma-rays accompanying the neutrinos are considerably attenuated inside RIAFs (see Table 2 for values of the break energy due to γγ → e+e−, εγγ), making the GeV gamma-ray intensity well below the Fermi data. Thus, LLAGN can be regarded as gamma-ray hidden neutrino sources46.
Our RIAF model provides a ~5–10% contribution to the cosmic soft X-ray (0.5–8 keV) background, which is dominated by canonical AGN and star-forming galaxies6,47. In observations of distant LLAGN, fainter ones are likely to be classified as star-forming galaxies, while relatively luminous ones will be indistinguishable from the faint end of canonical AGN. Thus, it is difficult to determine the contribution by LLAGN accurately. Here, we make a rough estimate of the LLAGN contribution. The luminous objects of LX > 3.2 × 1042erg s−1 provides 63% of the CXB flux. Our LLAGN should have X-ray luminosities below this range, and thus, 37% contribution can be regarded as an upper limit for the LLAGN contribution. Ref. 26 reported that 43% of nearby galaxies have signatures of AGN activity. Then, the LLAGN contribution can be as high as 0.43 × 0.37 ≃ 0.16. Our model predicts a lower LLAGN contribution than our rough estimate. Future X-ray and optical spectroscopic surveys with better sensitivities are necessary to unravel the X-ray contribution from LLAGN.
a Diffuse soft gamma-ray intensities from relatively luminous (dotted) and faint (dashed) LLAGN. The solid line shows the sum of these. Data points are provided by Swift BAT115 (brown-circle; 1-σ errors for systematic and statistical errors), SMM116 (brown-triangle; definition of errors are not available), Nagoya-baloon117 (purple-star; definition of errors are not available), COMPTEL2 (purple-square; 2-σ error bars for the linear sum of the systematic and statistical errors). b Same as (a) but for diffuse neutrino intensities. Data are provided by IceCube4,118 (lightgrey and magenta regions; 68% confidence levels).
Tests by future observations
High-energy multimessenger observations are essential for identifying the origin of the gamma-ray and neutrino backgrounds. We have estimated the detectability of neutrinos from RIAFs (see subsection Neutrino detectability from nearby LLAGN in Methods for details). The expected number of through-going muon track events from nearby LLAGN, \({{{{{{{{{{\mathcal{N}}}}}}}}}}}_{\mu }\), is shown in Fig. 5. Current facilities are not sufficient to detect the signal, even using stacking techniques. However, stacking ~10 LLAGN will enable the planned IceCube-Gen2 facility to detect a few neutrino events by means of its larger effective area. Its better angular resolution will reduce the atmospheric background, making individual source detections possible.
We plot the expected number of through-going muon track events above a given muon energy for a 10-yr operation with IceCube (thin-blue) and IceCube-Gen2-like detector (thick-red). The solid lines are for the case that stacks the ten brightest LLAGN, and the dashed lines are for the brightest LLAGN, NGC 4565. The dotted-dashed lines show the expected number of the atmospheric background for the case that stacks the ten LLAGN. The dotted horizontal line indicates \({{{{{{{{{{\mathcal{N}}}}}}}}}}}_{\mu }=1\) for comparison. Stacking more LLAGN does not help improving the detectability because the background also increases with the number of the stacked LLAGN29.
Future MeV satellites, such as e-ASTROGAM41, the All-sky Medium Energy Gamma-ray Observatory (AMEGO)48, and Gamma-Ray and AntiMatter Survey (GRAMS)43, will be able to measure the electron temperature in RIAFs by detecting a clear cutoff feature (see Fig. 2), which will shed light on the electron heating mechanisms in the collisionless plasma. In addition, anisotropy tests for the MeV gamma-ray background are promising44. Our RIAF model predicts a smaller anisotropy than those from Seyfert and blazar models owing to the higher source number density of LLAGN.
Possible effects of another acceleration mechanism
If protons are accelerated by another mechanism, such as magnetic reconnection, the resulting proton spectrum is usually described by a power-law form. In the upper panel of Fig. 6, we plot the diffuse neutrino and gamma-ray intensities for power-law injection models with the parameter sets tabulated in Table 3 (see subsection Power-law injection models for CRs in Methods). The power-law injection models can also account for the PeV neutrino data with PCR/Pth ~ 0.01 for model B and PCR/Pth ~ 0.1 for model C. A higher PCR/Pth is demanded for a higher sinj, which may lead some feedback to the MHD turbulence (see discussion in the next subsection). The cascade gamma-ray emissions are well below the Fermi data owing to the low γγ break energy in the RIAFs. Note that the resulting diffuse MeV gamma-ray intensities for models B and C are the same as that for model A, because the thermal electrons in the RIAFs are identical.
Data are provided by Fermi LAT3 (black-plus; 1-σ errors with systematic and statistical uncertainties), and IceCube4,118 (lightgrey and magenta regions; 68% confidence levels). a Power-law injection models that can account for the PeV neutrino data. The blue and green lines are for models B (sinj = 1) and C (sinj = 2) in Table 3, respectively. The solid and dotted-dashed lines indicate the neutrino spectra, and the dashed and dotted lines represent the gamma-ray spectra by proton-induced electromagnetic cascades. b Models for 10–100 TeV neutrinos. The thick-blue and thin-green lines are for neutrinos and proton-induced cascade gamma-rays, respectively. The dashed, dotted, dotted-dashed lines are for models D, E, F in Table 3, respectively.
Medium-energy neutrinos in the 10–100 TeV range
For the sake of completeness, we also investigate the possibility of whether our model can reproduce the 10–100 TeV neutrino background. The observed intensity of the 10–100 TeV neutrinos is higher than that of the 100 GeV gamma-rays. Hence, the emission region of origin of the 10–100 TeV neutrinos is expected to be opaque to GeV–TeV gamma-rays46.
The lower panel of Fig. 6 depicts the diffuse gamma-ray and neutrino intensities for models D (stochastic acceleration), E (power-law injection with sinj = 1), and F (power-law injection with sinj = 2), whose parameters are tabulated in Table 3. Since the thermal component is identical to that for model A, we only show the nonthermal components.
The neutrino intensity is normalized so that it matches the observed 10–100 TeV neutrinos, which requires higher values of ηCR and ηtur or ηacc. The gamma-ray intensity from proton-induced cascades is also higher than that of PeV neutrinos, but still considerably lower than the Fermi data. Therefore, our LLAGN model could in principle be the source of the mysterious 10–100 TeV neutrinos, although the contribution from Seyferts and quasars is typically more important in view of the energetics.
The resulting pressure ratio of the CRs to the thermal component is about 10% for models D and E and 50% for model F. These are much higher than those for models for PeV neutrinos. Nevertheless, even with such a high value of PCR/Pth, the global dynamical structure of RIAFs are very similar as long as CRs are confined inside the flow as shown in Ref. 49. However, such a high value of PCR/Pth is at odds with the picture of turbulent acceleration, and the feedback from CRs to MHD turbulence may be significant. In Ref. 29, we estimated the point-source detectability of the neutrinos from nearby LLAGN using the method in Ref. 50 with the parameter sets similar to those for models D, E, and F. The diffuse neutrino spectra for the models in Ref. 29 are almost identical to those for models D, E, and F. Ref. 29 showed that the TeV–PeV neutrinos and MeV gamma-rays are detectable with the planned neutrino and gamma-ray experiments, IceCube-Gen2 and AMEGO/e-ASTROGAM/GRAMS, respectively.
Discussion
The present work differs from earlier works in a number of respects, as described below. Ref. 23 calculated neutrino and gamma-ray emission from hot coronae surrounding the accretion disks in luminous AGN, i.e., Seyfert galaxies and quasars. The stochastic proton acceleration was considered, and proton-induced electromagnetic cascades are fully taken into account. It was shown that the high-energy emission from coronae can account for the X-ray and 10–100 TeV neutrino backgrounds that require gamma-ray hidden sources. Note that the structure of the systems considered in the present work and Ref. 23 are different. The present work considers the RIAFs in LLAGN, where CRs can be accelerated in the bulk of the accretion flows. A broadband of target photons are provided by the thermal electrons inside the RIAFs. Ref. 23 used the disk-corona paradigm in luminous AGN, where protons are accelerated at the hot corona above the geometrically thin, optically thick standard disk. The disks provide thermal UV photons, while the coronae supplies X-rays by upscattering the UV photons. Strictly speaking, the prescription of the nonthermal particle injection for our RIAF model is also different from that in Ref. 23. In the RIAF model, we use a constant ηCR parameter, i.e., the fraction of CR luminosity to the accretion luminosity is independent of \(\dot{m}\). On the other hand, in the AGN corona model, we considered that \({\dot{{{{{{{{{{\mathcal{F}}}}}}}}}}}}_{p,{{{{{{{{{\rm{inj}}}}}}}}}}}\) is proportional to LX, where ηCR depends on LX and \(\dot{m}\). The injection process of nonthermal CR acceleration remains an open problem, so we examined effects of both prescriptions. We find that our RIAF and the AGN corona models can reproduce 0.1–1 PeV and 10–100 TeV neutrino datas, respectively, with either of the prescriptions, with a similar value of PCR/Pth ~ 0.01. Our conclusions are nearly independent of these prescriptions about the injection. Ref. 27 calculated the neutrino background from RIAFs, considering stochastic particle acceleration by MHD turbulence. It considered neither gamma-ray emission nor neutrino point-source detectability, but studied the neutrino background for either 10–100 TeV or 0.1–1 PeV components. Ref. 29 focused on neutrino and gamma-ray emissions from nearby LLAGN with RIAFs, and discussed the neutrino and gamma-ray detectabilities in future experiments. In addition to the stochastic acceleration model, the power-law injection models were discussed. Ref. 51 considered gamma-ray and neutrino emission from magnetically arrested disks (MADs) in radio galaxies, motivated by Fermi detection of nearby radio galaxies. They assume an efficient acceleration close to the theoretical limit, which could be possible by magnetic reconnections. The maximum energy of the accelerated particles is much higher than in the other papers on luminous AGN and LLAGN, leading to efficient GeV gamma-ray production by the proton synchrotron process.
The Hα luminosity function has uncertainties. Ref. 52 proposed a lower value of Hα luminosity function. Then, an analytic estimate using Eq. (4) gives \({E}_{\gamma }^{2}{{{\Phi }}}_{\gamma } \sim 7\times 1{0}^{-7}\,{{{{{{{{{\rm{GeV}}}}}}}}}}\,{{{{{{{{{{\rm{s}}}}}}}}}}}^{-1}\,{{{{{{{{{{\rm{cm}}}}}}}}}}}^{-2}\,{{{{{{{{{{\rm{sr}}}}}}}}}}}^{-1}\), which is an order of magnitude lower than the MeV gamma-ray data as shown in Fig. 3. Even with the lower value of the luminosity function, the PeV neutrino data may be reproduced if we use ηCR ~ 0.02, which leads to a high CR energy density, PCR/Pth ~ 0.5. However, such a situation is at odds with the assumption that the turbulence energy is the source of the CRs. The uncertainty in the luminosity function should be palliated by future surveys using a line-sensitive optical instrument or sensitive X-ray satellites, such as Subaru Prime Focus Spectrograph53 or extended Roentgen Survey with an Imaging Telescope Array (eROSITA)54 and Focusing On Relativistic universe and Cosmic Evolution (FORCE)55.
The e+e− pair production processes are inefficient in RIAFs. The most efficient process is the two-photon (γγ) interaction. The optical depth for MeV photons to the e+e− pair production is estimated to be \({\tau }_{\gamma \gamma }\approx {n}_{t}{\sigma }_{\gamma \gamma }R\simeq 1.4\times 1{0}^{-3}\,{L}_{{{{{{{{{{\rm{IC,42}}}}}}}}}}}{R}_{14.5}^{-1}\), where σγγ ~ 0.2σT is the cross section for two-photon interactions, \({n}_{t} \sim {L}_{{{{{{{{{{\rm{IC}}}}}}}}}}}/(4\pi {R}^{2}{m}_{e}{c}^{3}) \sim 3.2\times 1{0}^{7}\,{L}_{{{{{{{{{{\rm{IC,42}}}}}}}}}}}{R}_{14.5}^{-2}\,{{{{{{{{{{\rm{cm}}}}}}}}}}}^{-3}\) is the target photon density at MeV energies, and we use convention of QX = Q/10X in cgs unit. Equating the e+e− pair production rate and advective escape rate, the number density of the electron-positron pairs is estimated to be51 \({n}_{\pm }\approx {n}_{t}{\tau }_{\gamma \gamma }(c/{V}_{R}) \sim 4\times 1{0}^{6}\,{L}_{{{{{{{{{{\rm{IC}}}}}}}}}},{{{{{{{{{\rm{42}}}}}}}}}}}^{2}{R}_{14.5}^{-5/2}{\alpha }_{-1}^{-1}\,{{{{{{{{{{\rm{cm}}}}}}}}}}}^{-3}\), which is a few orders of magnitude smaller than the electron density given by Eq. (6). Electron-electron (ee) and electron-proton (ep) interactions can also produce e+e− pairs, whose timescales are roughly approximated to be \({t}_{ep,\pm } \sim {t}_{ee,\pm } \sim 1/({n}_{p}{\alpha }_{{{{{{{{{{\rm{em}}}}}}}}}}}^{2}{\sigma }_{T}c)\), where αem is the fine structure constant. Then, the ratio of the e+e− pair production timescale to the infall timescale is estimated to be \({t}_{{{{{{{{{{\rm{fall}}}}}}}}}}}/{t}_{ee,\pm }\sim 4\times 1{0}^{-4}{\alpha }_{-1}^{-2}{\dot{m}}_{-2}\), and thus, ee and ep interactions cannot provide sufficient amount of pairs. Photon-proton (γp) and photon-electron (γe) interactions also produce pairs, but the cross sections for these processes at εγ ~ 3–4 MeV are similar to those for ee and ep interactions, which leads to pair production rates similar to those by ee and ep interactions. Photons with a higher energy have larger γp and γe cross sections, but the number density of such photons is too small to produce the pairs efficiently. Therefore, the RIAF plasma is unlikely to reach the pair equilibrium, in which the density of the electron-positron pairs is much lower than the proton density. Ref. 56 demonstrated this conclusion by detailed calculations.
While this work demonstrates that LLAGN can significantly contribute to the higher-energy part of IceCube neutrinos, contrary to the guaranteed MeV gamma-ray emission, their non-thermal contribution might in principle be much smaller if the CR acceleration in the RIAF disks is more inefficient than in the coronae of radio-quiet AGN. This may be the case if the CR acceleration occurs predominantly in a low-β plasma, because β ~ 3–30 in a RIAF disk is higher than β ≲ 1 in AGN coronae. In this case, we would need other models that can explain the neutrino background in the PeV range57.
In conclusion, we proposed RIAFs in LLAGN as a promising origin for the soft gamma-ray background. We constructed a one-zone model that can reproduce the observed X-ray features of LLAGN, and demonstrated that LLAGN can also simultaneously account for the high-energy neutrino background. In the RIAFs, electrons are thermalized and emit soft gamma rays through Comptonization. The protons there are naturally accelerated by reconnection or turbulence due to their longer thermalization timescale, and produce high-energy neutrinos efficiently. The accompanying gamma-rays are significantly attenuated by two-photon interactions, resulting in a gamma-ray intensity well below the Fermi data. Since hot coronae in luminous AGN can produce 10–100 TeV neutrinos through the same mechanism23, accretion flows in AGN can account for a wide range of high-energy photon (keV–MeV) and neutrino (TeV–PeV) backgrounds. This scenario does not require any nonthermal electron population to account for the MeV background, which implies that the transition energy from the thermal to nonthermal Universe is higher than previously expected.
Methods
Emission from thermal electrons in RIAFs
Here, we describe the properties of thermal plasma in RIAFs24,58 in detail. Hereafter, we use the notation of QX = Q/10X in cgs unit unless otherwise noted. We consider an accreting plasma of size R around a supermassive black hole (SMBH) of mass M with an accretion rate \(\dot{M}\). We use the normalized radius and mass accretion rate, \({{{{{{{{{\mathcal{R}}}}}}}}}}=R/{R}_{S}\) and \(\dot{m}=\dot{M}{c}^{2}/{L}_{{{{{{{{{{\rm{Edd}}}}}}}}}}}\), where RS is the Schwarzschild radius and LEdd is the Eddington luminosity. Plasma quantities in RIAFs are described by two parameters, the viscosity parameter, α, and the pressure ratio of gas to magnetic field, β. Based on recent numerical simulations (e.g., Refs. 59,60,61,62,63,64), the radial velocity, number density, proton thermal temperature, magnetic field, Thomson optical depth, and Alfvén velocity in the RIAF are analytically approximated to be
where \({V}_{K}=\sqrt{GM/R}\) is the Keplerian velocity and H ≈ R/2 is the scale height. Observations of X-ray binaries and AGN demand α ~ 0.1–165, while the global MHD simulations result in α ~ 0.01–0.1 and β ~ 3–3061,66. Hence, we set α = 0.1 and β = 7.0 as their reference values.
Thermal electrons in RIAFs emit broadband photons through synchrotron radiation, bremsstrahlung, and inverse Compton scattering. The electron temperature is determined so that the resulting photon luminosity is equal to the bolometric luminosity estimated by \(\dot{m}\). Assuming that thermal electrons are heated by Coulomb collisions with protons, the bolometric luminosity is estimated to be \({L}_{{{{{{{{{{\rm{bol}}}}}}}}}}}\approx {\eta }_{{{{{{{{{{\rm{rad}}}}}}}}}},{{{{{{{{{\rm{sd}}}}}}}}}}}{\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}{L}_{{{{{{{{{{\rm{Edd}}}}}}}}}}}{(\dot{m}/{\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}})}^{2}\), where ηrad,sd ~ 0.1 is the radiation efficiency for the standard disk, and \({\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}\approx 0.03{(\alpha /0.1)}^{2}\) is the critical mass accretion rate above which RIAFs no longer exist30. We calculate the photon spectra by synchrotron radiation, bremsstrahlung, and inverse Compton scattering by thermal electrons with the steady-state and one-zone approximations. The synchrotron and bremsstrahlung spectra are calculated by the method given in Appendix of Ref. 27, where we use the fitting formulae for the emissivity of these processes and Eddington approximation to take into account the effects of the radiative transfer. For the inverse Compton scattering spectrum, we utilize the corrected delta-function method given in Ref. 67. In this method, the distribution of the scattered photon energy is approximated to be a delta function, but the mean energy of the scattered photon is calculated using the exact kernel. This method approximately takes into account the electron recoil effect for εγ ≳ kBTe. The error of the method is about 50%. This is sufficiently accurate for our purpose, considering significant uncertainty in the MeV gamma-ray data. Given the uncertainty, our calculation results are consistent with those by the Monte Carlo simulations68. The spectral decline due to the cutoff in this method is somewhat stronger than that in the exact method, implying that our results on the MeV fluxes are regarded as conservative. In this work, we assume that the Coulomb heating is the dominant electron heating mechanism, and \({L}_{{{{{{{{{{\rm{bol}}}}}}}}}}}\propto {\dot{m}}^{2}\) is used. This treatment is qualitatively different from the previous work27, where \({L}_{{{{{{{{{{\rm{bol}}}}}}}}}}}\propto \dot{m}\) is assumed. Such a treatment may be more appropriate if the electrons are directly heated by the plasma dissipation process69,70,71,72. Nevertheless, these details will not change our conclusions that LLAGN are bright in soft MeV gamma-rays. Also, the electron heating prescription do not strongly affect the critical mass accretion rate above which RIAFs no longer exist28,30,73,74.
Nonthermal particles in RIAFs
Here, we describe the details of the stochastic acceleration model, where protons are accelerated by MRI turbulence64,75. To obtain the CR spectrum, we solve the transport equation for CR protons, which is a diffusion equation in the momentum space76:
where \({{{{{{{{{{\mathcal{F}}}}}}}}}}}_{p}\) is the momentum distribution function for protons (\(dN/d{\varepsilon }_{p}=4\pi {p}^{2}{{{{{{{{{{\mathcal{F}}}}}}}}}}}_{p}/c\)), \({D}_{{\varepsilon }_{p}}\) is the diffusion coefficient that mimics the stochastic particle acceleration, tcool is the cooling time for protons, tesc is the escape time, \({\dot{{{{{{{{{{\mathcal{F}}}}}}}}}}}}_{p,{{{{{{{{{\rm{inj}}}}}}}}}}}={\dot{{{{{{{{{{\mathcal{F}}}}}}}}}}}}_{0}\delta ({\varepsilon }_{p}-{\varepsilon }_{p,{{{{{{{{{\rm{inj}}}}}}}}}}})\) is the injection function, and εp,inj is the initial energy of the particles injected to the stochastic acceleration process. We assume that the particles are injected to the stochastic acceleration process by fast acceleration processes such as magnetic reconnections, which are induced by MRI77,78,79. We consider a delta-function injection term with εp,inj much higher than the thermal proton energy, which may mimic the injection by the reconnection. The value of εp,inj has no influence on the resulting spectrum as long as εp,inj is much lower than the cutoff energy. We consider resonant scatterings between MHD waves and CR particles, where CR particles interact with the turbulent eddy of their gyration radius, rL = εp/(eB). Then, diffusion coefficient in energy space can be written as
where ηtur = B2/(8π∫Pkdk) is the turbulence parameter (Pk is the turbulence power spectrum), q is the power-law index of Pk, and we set the scale height, H, to the injection scale of the MHD turbulence. We assume a Kolmogorov turbulence, q = 5/3. The acceleration time is given by \({t}_{{{{{{{{{{\rm{acc}}}}}}}}}}}\approx {\varepsilon }_{p}^{2}/{D}_{{\varepsilon }_{p}}\approx {\eta }_{{{{{{{{{{\rm{tur}}}}}}}}}}}{\beta }_{A}^{-2}(H/c){[{\varepsilon }_{p}/(eBH)]}^{1/3}\). ηtur is lower for a stronger turbulence, and a low value of ηtur results in a higher maximum energy of the CR protons. Since RIAFs are expected to be turbulent due to MRI80,81, the value of ηtur should be small. The turbulent strength is also related to the value of α, as stronger turbulence leads to a higher α. Our fiducial value of ηtur ~ 15 is reasonable in the sense that ηtur is close to α−1. The amount of CRs is determined so that \(\int {L}_{{\varepsilon }_{p}}d{\varepsilon }_{p}={\eta }_{{{{{{{{{{\rm{CR}}}}}}}}}}}\dot{m}{L}_{{{{{{{{{{\rm{Edd}}}}}}}}}}}\) is satisfied, where \({L}_{{\varepsilon }_{p}}={t}_{{{{{{{{{{\rm{loss}}}}}}}}}}}^{-1}{\varepsilon }_{p}d{N}_{p}/d{\varepsilon }_{p}\) is the differential proton luminosity23, ηCR is the CR production efficiency, and \({t}_{{{{{{{{{{\rm{loss}}}}}}}}}}}^{-1}={t}_{{{{{{{{{{\rm{cool}}}}}}}}}}}^{-1}+{t}_{{{{{{{{{{\rm{esc}}}}}}}}}}}^{-1}\) is the total energy loss rate including cooling and escape processes.
We solve the transport equation until the steady-state is reached using the Chang-Cooper method82,83. As the proton cooling mechanism, we consider the proton synchrotron, Bethe-Heitler (p + γ → p + e+ + e−), photomeson (p + γ → p + π), and pp inelastic collision (p + p → p + p + π) processes, The timescale of the pp inelastic collisions is given as \({t}_{pp}^{-1}\approx {n}_{p}{\sigma }_{pp}{\kappa }_{pp}c\), where σpp and κpp are the cross section and inelasticity of pp interactions84. The photomeson production and Bethe-Heitler cooling timescales are estimated to be
where γp = εp/(mpc2), \({\overline{\varepsilon }}_{\gamma }\) is the photon energy in the proton rest frame, and σi and κi are the cross section and inelasticity for the process (i = pγ for photomeson production85 and i = BH for Bethe-Heitler process86,87). The cooling time by the proton synchrotron is \({t}_{p,{{{{{{{{{\rm{syn}}}}}}}}}}}\approx 6\pi {m}_{p}^{4}{c}^{3}/({\sigma }_{T}{B}^{2}{\varepsilon }_{p})\). The total cooling rate is then given by \({t}_{{{{{{{{{{\rm{cool}}}}}}}}}}}^{-1}={t}_{pp}^{-1}+{t}_{p\gamma }^{-1}+{t}_{{{{{{{{{{\rm{BH}}}}}}}}}}}^{-1}+{t}_{p,{{{{{{{{{\rm{syn}}}}}}}}}}}^{-1}\). Regarding the escape process, we only consider the advective escape, i.e., infall to the SMBH. We write the escape timescale as tesc = tadv ≈ R/VR. We consider that the CR component has the same bulk velocity as that for the thermal component owing to efficient interactions through the turbulent magnetic field. The diffusive escape may be inefficient because the high-energy protons tend to move in the azimuthal direction, which is the direction of the background magnetic field in differentially rotating accretion flows64,75. See also Refs. 23,29 for technical details of the calculation methods for the cooling and escape timescales.
In Fig. 7, we plot the acceleration and loss rates of CRs for model A (reference model) for NGC 4579, whose parameters are shown in Table 1 and caption. The dominant loss process is the advective escape for εp ≲ 1 × 108 GeV. Although the acceleration time is longer than the advective escape time for εp ≳ 2 × 105 GeV, the CR proton spectrum continues to a higher-energy because the weak εp dependence of the advective escape leads to a very gradual cutoff in the proton spectrum27,88. At εp ~ 1 × 107 GeV, the photomeson production becomes more efficient than the acceleration, which makes a sharp cutoff owing to a strong εp dependence (see Fig. 2).
The thick-solid, thick-long-dashed, thick-dotted, thick-short-dashed, and thick-dotted-dashed lines are energy loss rates by photomeson production, Bethe-Heitler, synchrotron, pp inelastic collision, and infall processes, respectively. The thin-dotted line indicates the acceleration rate. We use M = 7.2 × 107M⊙ and \(\dot{m}=8.0\times 1{0}^{-3}\) with a parameter set for model A (reference model: see Table 1).
The CRs produce pions via both pp and pγ interactions, and pions decay to gamma-rays, electron/positron pairs, and neutrinos (π0 → 2γ; π± → e± + 3ν). These neutrinos are believed to explain IceCube neutrinos89,90,91. We calculate neutrino spectra by pp collisions using the formalism given by Ref. 92. For the neutrinos by pγ interactions, we use a semi-analytic prescription given in Refs. 93,94 (see, e.g., Ref. 85 for numerical results). Since the effect of meson cooling is negligible in the RIAFs, neutrino flavor ratio is νe:νμ:ντ = 1:2:0 at the source, which becomes ~ 1:1:1 on Earth after the flavor mixing. As shown in Fig. 2, neutrinos are mainly produced by the pp collisions for εν ≲ 4 × 105 GeV, and the photomeson production is effective above the energy. The Bethe-Heitler and proton synchrotron processes are subdominant in the range we investigated. For higher \(\dot{m}\) cases, the photomeson production is more efficient, and hence, the peak energy of the proton spectrum is lower. In our model, the neutrino background is dominated by relatively faint LLAGN, and their target photon spectra are not hard. Then, the multi-pion production channel is subdominant, and our approximation provides reasonably accurate results.
The hadronic interactions also produce gamma-rays and electron/positron pairs. The gamma-rays are absorbed by two-photon annihilation, and create electron/positron pairs. These pairs also emit gamma-rays, and electromagnetic cascades are initiated. We calculate the cascade emission by solving the kinetic equations of electron/positron pairs and photons95,96:
where \({n}_{{\varepsilon }_{i}}\) is the differential number density (i = e or γ), \({\dot{n}}_{{\varepsilon }_{i}}^{(xx)}\) is the particle source term from the process xx (xx = IC (inverse Compton scattering), γγ (e+e− pair production by γγ interactions), syn (synchrotron), or ff (bremsstrahlung)), \({\dot{n}}_{{\varepsilon }_{i}}^{{{{{{{{{{\rm{inj}}}}}}}}}}}\) is the injection term from the hadronic interaction, and Pyy is the energy loss rate for the electrons from the process yy (yy = IC (inverse Compton scattering), syn (synchrotron), ff (bremsstrahlung), or Cou (Coulomb collision)). See also Refs. 95,96 for technical details. We approximately treat the pair injection processes by the Bethe-Heitler process and photomeson production as in Refs. 23,29.
In our RIAF models, secondary pairs do not contribute to the the high-energy gamma-ray background, even for the case that the secondary pairs are re-energized by MHD turbulence. The secondary pairs suffer from a strong cooling by the synchrotron emission, whose timescale is estimated to be \({t}_{e,{{{{{{{{{\rm{syn}}}}}}}}}}}=6\pi {m}_{e}^{2}{c}^{3}/({\sigma }_{T}{B}^{2}{\varepsilon }_{e})\). When their energy becomes sufficiently low, they may be re-energized by MHD turbulence, as suggested by Ref. 23. The energization timescale is the same with the proton acceleration timescale: \({t}_{e,{{{{{{{{{\rm{acc}}}}}}}}}}}\approx {\eta }_{{{{{{{{{{\rm{tur}}}}}}}}}}}{\beta }_{A}^{-2}(H/c){[{\varepsilon }_{e}/(eBH)]}^{2-q}\). Equating these two timescales, we obtain the critical energy at which the secondary pairs piles up:
where we use q = 5/3 for the last equation. The turbulence power below the mean thermal proton energy, \(3{k}_{B}{T}_{p} \sim 35{{{{{{{{{{\mathcal{R}}}}}}}}}}}^{-1}\) MeV, is significantly reduced due to dissipation by plasma kinetic effects, and the turbulent re-energization is not expected when εe,crit < 3kBTp. In our model A (fiducial parameters), we cannot expect re-energization of secondary pairs even for the case with \({L}_{{{{{{{{{{\rm{H}}}}}}}}}}\alpha }={L}_{\min }\). For a further lower \(\dot{m}\) case, re-energization may occur. In such a case, we can ignore the inverse Compton component by the secondary pairs owing to its lower photon energy density (\({B}^{2}\propto \dot{m}\), \({L}_{{{{{{{{{{\rm{bol}}}}}}}}}}}\propto {\dot{m}}^{2}\)). Then, the synchrotron peak energy for the re-energized pairs are estimated to be
Hence, we conclude that the re-energized pairs cannot contribute to the MeV gamma-ray background for all the \(\dot{m}\) range in our model. Also, primary electrons are not expected to be produced efficiently in RIAFs, because of their rapid thermalization in the range of our interest28. Thus, they do not contribute to the MeV gamma-ray background.
Cumulative background intensities
Here we describe the method to obtain the background intensities. Since the Hα luminosity functions include much fainter sources than the X-ray luminosity functions, we use the luminosity function for type-1 Seyfert galaxies provided by Ref. 45: \(d\rho /d{L}_{{{{{{{{{{\rm{H}}}}}}}}}}\alpha }\approx ({\rho }_{* }/{L}_{* })/[{({L}_{{{{{{{{{{\rm{H}}}}}}}}}}\alpha }/{L}_{* })}^{{s}_{1}}+{({L}_{{{{{{{{{{\rm{H}}}}}}}}}}\alpha }/{L}_{* })}^{{s}_{2}}]\), where ρ* ≈ 1.2 × 10−6Mpc−3, L* = 3.7 × 1042erg s−1, s1 = 2.05, and s2 = 5.12 (the values of L* and ρ* are corrected using Hubble constant of H0 = 70 km s−1Mpc−1). Type-1 Seyfert galaxies exhibit broad emission lines that are unique to AGN. We extrapolate this luminosity function to \({L}_{\min }=1{0}^{38}\,{{{{{{{{{\rm{erg}}}}}}}}}}\,{{{{{{{{{{\rm{s}}}}}}}}}}}^{-1}\), below which the Palomar survey finds a hint of a break26. The survey also indicates a correlation between X-ray luminosity, LX, and LHα for LLAGN. The ratio, κX/Hα = LX/LHα ranges 5 ≲ κX/Hα ≲ 7 in the luminosity range of our interest for type-1 AGN. We use κX/Hα = 6.0 for simplicity26, but the difference from the cases with κX/Hα = 5 or κX/Hα = 7 is less than a factor of 1.2.
Observationally, the X-ray luminosity at the 2–10 keV band can be converted to the bolometric luminosity using the bolometric correction factor, κbol/X = Lbol/LX ≃ 15 for LLAGN97,98,99,100. Using the two correction factor, κbol/X and κX/Hα, we can convert \(\dot{m}\) to LHα if we fix a SMBH mass, M. Ref. 101 provided a sample of LLAGN, and the mean and median values of \({{{{{{{{\mathrm{log}}}}}}}}}\,({M}/{{{{{{{{{{\rm{M}}}}}}}}}}}_{\odot })\) are 8.0 and 8.1, respectively. Also, the X-ray luminosity density is dominated by AGN with M ~ 108 − 3 × 108M⊙ if the Eddington ratio function is independent of the SMBH mass97,102,103. Thus, we set M = 108M⊙ as a reference value. With this prescription, the critical Hα luminosity above which RIAFs no longer exist is estimated to be \({L}_{{{{{{{{{{\rm{crit}}}}}}}}}}}={\eta }_{{{{{{{{{{\rm{rad}}}}}}}}}},{{{{{{{{{\rm{sd}}}}}}}}}}}{\dot{m}}_{{{{{{{{{{\rm{crit}}}}}}}}}}}{L}_{{{{{{{{{{\rm{Edd}}}}}}}}}}}/({\kappa }_{X/{{{{{{{{{\rm{H}}}}}}}}}}\alpha }{\kappa }_{{{{{{{{{{\rm{bol}}}}}}}}}}/X})\simeq 4.2\times 1{0}^{41}\,{{{{{{{{{\rm{erg}}}}}}}}}}\,{{{{{{{{{{\rm{s}}}}}}}}}}}^{-1}\). Finally, we integrate the gamma-ray and neutrino fluxes over the range of \({L}_{\min }\le {L}_{{{{{{{{{{\rm{H}}}}}}}}}}\alpha }\le {L}_{{{{{{{{{{\rm{crit}}}}}}}}}}}\), which corresponds to \(4.6\times 1{0}^{-4} \, < \,\dot{m} \, < \, 0.03 \,{{{{{\rm{for}}}}}}\,\alpha=0.1\). The background intensity is written as104
where we use the cosmological parameters of Ωm ≈ 0.3 and ΩΛ ≈ 0.7. Since dimmer AGN tend to have weaker redshift evolution based on the X-ray, gamma-ray, and radio luminosity functions103,105,106, no redshift evolution of the Hα luminosity function is used. With this treatment, LLAGN at z < 0.5 contribute about a half of the diffuse intensity, objects at 0.5 < z < 2 provide the other half, and the contribution from z > 2 is negligible. We consider the gamma-ray attenuation by the extragalactic background light, and include the exponential suppression with the optical depth given in Ref. 107. The two-photon pair annihilation initiates intergalactic electromagnetic cascades108, but it has little influence on the GeV gamma-ray intensity because TeV gamma-rays cannot escape from the RIAFs due to two-photon annihilation inside the RIAFs (see Table 2 for the value of εγγ, the γγ break energy above which two-photon annihilation attenuates the gamma-rays inside the RIAFs).
Ref. 45 also provided the luminosity function for type-2 AGN that exhibit only a narrow emission line feature, but we find that their contributions to the MeV gamma-ray and PeV neutrino backgrounds are negligible because of two reasons. One is the value of κX/Hα. Type-2 AGN have κX/Hα ~ 126, and neutrino and gamma-ray luminosity scales with \({L}_{\nu }\propto {\dot{m}}^{2}\propto {L}_{{{{{{{{{{\rm{bol}}}}}}}}}}}\propto {\kappa }_{X/{{{{{{{{{\rm{H}}}}}}}}}}\alpha }\). The other is the shape of the luminosity function. The break luminosity for type-2 AGN is lower than the critical luminosity, L* ≈ 2.8 × 1040erg s−1 < Lcrit, resulting in a lower MeV gamma-ray intensity. The slope in the low-luminosity range is lower than 2, s1 = 1.77, which makes the contribution by relatively faint LLAGN smaller, leading to a lower PeV neutrino intensity. In addition, type-2 objects may be contaminated by non-AGN objects, the fraction of which is still uncertain. Therefore, we focus on contribution by type-1 AGN in this study.
Note that the Hα luminosity function does not match the X-ray luminosity function with a constant κX/Hα. The conversion factor depends on the luminosity, and there is some degree of uncertainty in the conversion factor (see Ref. 52). The Hα luminosity function by Ref. 52 matches the observed X-ray luminosity function at a high luminosity range, while it underestimates the number density for LX ≲ 1042erg s−1. This range is the most relevant to the MeV gamma-ray and PeV neutrino backgrounds, and the luminosity function by Ref. 45 gives a higher number density at the range. From this reason, we use one by Ref. 45 as a reference case. Future optical spectroscopic and hard X-ray surveys are necessary to reduce this uncertainty.
Neutrino detectability from nearby LLAGN
The number of through-going muon track events is calculated by50,109
where Eν is the incoming neutrino energy, Eμ is the muon energy, \({\phi }_{\nu }={{\Delta }}t{L}_{{E}_{\nu }}/(4\pi {d}_{L}^{2}{E}_{\nu })\) is the neutrino fluence from a LLAGN for a time interval of Δt, \({{{{{{{{{{\mathcal{N}}}}}}}}}}}_{A}\) is the Avogadro Number, σCC is the charged-current cross section, τνN is the optical depth for the neutrino scattering in the Earth, and αμ + βμEμ represents the energy loss rate of muons. We use Δt = 10 yr and the list of LLAGN given by Ref. 101 (see Ref. 29 for the brightest 10 LLAGN in the list). This method can reproduce the effective area by Ref. 110. For IceCube-Gen2, we use a 102/3 times bigger \({A}_{\det }\) than that for IceCube111. The main background of the astrophysical neutrino is atmospheric background, whose intensity depends on the declination. We appropriately take into account both conventional and prompt atmospheric muon neutrinos as well as the attenuation in the Earth50.
Power-law injection models for CRs
The CR spectrum inside the RIAFs can be obtained by solving the transport equation with a power-law injection term:
where εp,cut is the cutoff energy for the injected protons and \({\dot{N}}_{0}\) is the normalization factor. The injection is normalized such that \(\int {\varepsilon }_{p}{\dot{N}}_{{\varepsilon }_{p},{{{{{{{{{\rm{inj}}}}}}}}}}}d{\varepsilon }_{p}={\eta }_{{{{{{{{{{\rm{CR}}}}}}}}}}}\dot{m}{L}_{{{{{{{{{{\rm{Edd}}}}}}}}}}}\) is satisfied. We estimate \({\varepsilon }_{p,\max }\) by equating the infall timescale to the acceleration timescale, tacc. We set tacc = ηaccrL/c, where ηacc is the acceleration efficiency parameter. The value of ηacc is quite uncertain, but ηacc ≫ 1 is expected in the RIAFs because ηacc ~ 1 is achieved only for relativistic shocks or relativistic reconnections (i.e., magnetic energy density is higher than the rest mass energy). Here, we tune it so that our RIAF models can explain IceCube neutrino data. Equation (20) has an analytic solution:
The cooling and loss processes are the same with the stochastic acceleration model. See Ref. 29 for details.
Data availability
The data generated in this study have been deposited in https://doi.org/10.6084/m9.figshare.15170790. The observational data are obtained from the tables and fitting results provided in the references shown in Figure captions. The extracted data files are available upon reasonable request.
Code availability
The codes used for this study are available from S.S.K. and K.M. upon reasonable request.
References
Fabian, A. C. & Barcons, X. The origin of the X-ray background. Annu. Rev. Astron. Astrophys. 30, 429–456 (1992).
Weidenspointner, G. et al. The cosmic diffuse gamma-ray background measured with COMPTEL. Vol. 510 (eds McConnell, M. L. & Ryan, J. M.) 467–470 (American Institute of Physics Conference Series, 2000).
Ackermann, M. et al. The spectrum of isotropic diffuse gamma-ray emission between 100 MeV and 820 GeV. Astrophys. J. 799, 86 (2015).
Aartsen, M. et al. Characteristics of the diffuse astrophysical electron and tau neutrino flux with six years of IceCube high energy cascade data. Phys. Rev. Lett. 125, 121104 (2020).
Ueda, Y., Akiyama, M., Ohta, K. & Miyaji, T. Cosmological Evolution of the Hard X-Ray Active Galactic Nucleus Luminosity Function and the Origin of the Hard X-Ray Background. Apj 598, 886–908 (2003).
Luo, B. et al. The Chandra Deep Field-South Survey: 7 Ms Source Catalogs. ApJS 228, 2 (2017).
Fornasa, M. & Sánchez-Conde, M. A. The nature of the Diffuse Gamma-Ray Background. Phys. Rep. 598, 1–58 (2015).
Ajello, M. et al. The Origin of the Extragalactic Gamma-Ray Background and Implications for Dark Matter Annihilation. Apj 800, L27 (2015).
Inoue, Y. Cosmic Gamma-ray Background Radiation. arXiv e-prints arXiv:1412.3886 (2014).
Ahlers, M. & Halzen, F. IceCube: neutrinos and multimessenger astronomy. Prog. Theor. Exp. Phys. 2017, 12A105 (2017).
Clayton, D. D. & Ward, R. A. On emission lines in the cosmic gamma-ray background. Apj 198, 241–244 (1975).
Watanabe, K., Hartmann, D. H., Leising, M. D. & The, L. S. The Diffuse Gamma-Ray Background from Supernovae. Apj 516, 285–296 (1999).
Strigari, L. E., Beacom, J. F., Walker, T. P. & Zhang, P. The concordance cosmic star formation rate: implications from and for the supernova neutrino and gamma ray backgrounds. J. Cosmol. Astropart. Phys. 2005, 017 (2005).
Ruiz-Lapuente, P. et al. The Origin of the Cosmic Gamma-ray Background in the MeV Range. Apj 820, 142 (2016).
Ruiz-Lapuente, P. & Korobkin, O. Gamma-Rays from Kilonovae and the Cosmic Gamma-Ray Background. Apj 892, 45 (2020).
Strong, A. W., Wolfendale, A. W. & Worrall, D. M. Diffuse γ rays from discrete extragalactic sources. J. Phys. A Math. Gen. 9, 1553–1566 (1976).
Ajello, M. et al. The Evolution of Swift/BAT Blazars and the Origin of the MeV Background. Apj 699, 603–625 (2009).
Inoue, Y. Contribution of Gamma-Ray-loud Radio Galaxies’ Core Emissions to the Cosmic MeV and GeV Gamma-Ray Background Radiation. Apj 733, 66 (2011).
Ajello, M. et al. The Luminosity Function of Fermi-detected Flat-spectrum Radio Quasars. Apj 751, 108 (2012).
Toda, K., Fukazawa, Y. & Inoue, Y. Cosmological Evolution of Flat-spectrum Radio Quasars Based on the Swift/BAT 105 Month Catalog and Their Contribution to the Cosmic MeV Gamma-Ray Background Radiation. Apj 896, 172 (2020).
Stecker, F. W., Salamon, M. H. & Done, C. On the Origin of the MeV Gamma-Ray Background. arXiv e-prints astro–ph/9912106 (1999).
Inoue, Y., Totani, T. & Ueda, Y. The Cosmic MeV Gamma-Ray Background and Hard X-Ray Spectra of Active Galactic Nuclei: implications for the Origin of Hot AGN Coronae. Apj 672, L5 (2008).
Murase, K., Kimura, S. S. & Meszaros, P. Hidden cores of active galactic nuclei as the origin of medium-energy neutrinos: critical tests with the MeV gamma-ray connection. Phys. Rev. Lett. 125, 011101 (2020).
Narayan, R. & Yi, I. Advection-dominated accretion: a self-similar solution. Apj 428, L13–L16 (1994).
Yuan, F. & Narayan, R. Hot Accretion Flows Around Black Holes. ARAA 52, 529–588 (2014).
Ho, L. C. Nuclear Activity in Nearby Galaxies. ARAA 46, 475–539 (2008).
Kimura, S. S., Murase, K. & Toma, K. Neutrino and Cosmic-Ray Emission and Cumulative Background from Radiatively Inefficient Accretion Flows in Low-luminosity Active Galactic Nuclei. Apj 806, 159 (2015).
Mahadevan, R. & Quataert, E. Are Particles in Advection-dominated Accretion Flows Thermal? Apj 490, 605–618 (1997).
Kimura, S. S., Murase, K. & Mészáros, P. Multimessenger tests of cosmic-ray acceleration in radiatively inefficient accretion flows. Phys. Rev. D. 100, 083014 (2019).
Mahadevan, R. Scaling laws for advection-dominated flows: applications to low-luminosity galactic nuclei. Apj 477, 585–601 (1997).
Rybicki, G. B. & Lightman, A. P. Radiative Processes In Astrophysics (Wiley-VCH, Weinheim, Germany, 1979).
Esin, A. A., McClintock, J. E. & Narayan, R. Advection-Dominated Accretion and the Spectral States of Black Hole X-Ray Binaries: application to Nova Muscae 1991. Apj 489, 865–889 (1997).
Narayan, R. Advection-dominated Models of Luminous Accreting Black Holes. Apj 462, 136 (1996).
Younes, G. et al. NuStar Hard X-Ray View of Low-luminosity Active Galactic Nuclei: High-energy Cutoff and Truncated Thin Disk. Apj 870, 73 (2019).
Nemmen, R. S., Storchi-Bergmann, T. & Eracleous, M. Spectral models for low-luminosity active galactic nuclei in LINERs: the role of advection-dominated accretion and jets. MNRAS 438, 2804–2827 (2014).
Contini, M. The complex structure of low-luminosity active galactic nuclei: NGC 4579. MNRAS 354, 675–683 (2004).
Sridhar, S. S. et al. LOFAR view of NGC 3998, a sputtering AGN. AA 634, A108 (2020).
Akiyama, K. et al. First M87 Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole. Astrophys. J. 875, L1 (2019).
Akiyama, K. et al. First M87 Event Horizon Telescope Results. V. Physical Origin of the Asymmetric Ring. Astrophys. J. 875, L5 (2019).
de Menezes, R., Nemmen, R., Finke, J. D., Almeida, I. & Rani, B. Gamma-ray observations of low-luminosity active galactic nuclei. MNRAS 492, 4120–4130 (2020).
De Angelis, A. et al. The e-ASTROGAM mission. Exploring the extreme Universe with gamma rays in the MeV - GeV range. Exp. Astron. 44, 25–82 (2017).
Moiseev, A. & Amego Team. All-Sky Medium Energy Gamma-ray Observatory (AMEGO). International Cosmic Ray Conference (POS ICRC 2017). 301, 798 (2017).
Aramaki, T., Adrian, P. O. H., Karagiorgi, G. & Odaka, H. Dual MeV gamma-ray and dark matter observatory—GRAMS Project. Astropart. Phys. 114, 107–114 (2020).
Inoue, Y., Murase, K., Madejski, G. M. & Uchiyama, Y. Probing the Cosmic X-Ray and MeV Gamma-Ray Background Radiation through the Anisotropy. Apj 776, 33 (2013).
Hao, L. et al. Active Galactic Nuclei in the Sloan Digital Sky Survey. II. Emission-Line Luminosity Function. AJ 129, 1795–1808 (2005).
Murase, K., Guetta, D. & Ahlers, M. Hidden Cosmic-Ray Accelerators as an Origin of TeV-PeV Cosmic Neutrinos. Phys. Rev. Lett. 116, 071101 (2016).
Bauer, F. E. et al. The Fall of Active Galactic Nuclei and the Rise of Star-forming Galaxies: a Close Look at the Chandra Deep Field X-Ray Number Counts. AJ 128, 2048–2065 (2004).
McEnery, J. et al. All-sky Medium Energy Gamma-ray Observatory: Exploring the Extreme Multimessenger Universe. Bullet. Am. Astron. Soc. 51, 245 (2019).
Kimura, S. S., Toma, K. & Takahara, F. Effects of High-energy Particles on Accretion Flows onto a Supermassive Black Hole. Apj 791, 100 (2014).
Murase, K. & Waxman, E. Constraining High-Energy Cosmic Neutrino Sources: implications and prospects. Phys. Rev. D94, 103006 (2016).
Kimura, S. S. & Toma, K. Hadronic High-energy Emission from Magnetically Arrested Disks in Radio Galaxies. Apj 905, 178 (2020).
Greene, J. E. & Ho, L. C. The Mass Function of Active Black Holes in the Local Universe. Apj 667, 131–148 (2007).
Takada, M. et al. Extragalactic science, cosmology, and Galactic archaeology with the Subaru Prime Focus Spectrograph. PASJ 66, R1 (2014).
Merloni, A. et al. eROSITA Science Book: Mapping the Structure of the Energetic Universe. ArXiv e-prints (2012).1209.3114.
Mori, K. et al. A broadband x-ray imaging spectroscopy with high-angular resolution: the FORCE mission. Society of Photo-Optical Instrumentation Engineers (SPIE). Conference Series 9905, 99051O (2016).
Kusunose, M. & Mineshige, S. Effects of Electron-Positron Pairs in Advection-dominated Disks. Apj 468, 330 (1996).
Fang, K. & Murase, K. Linking high-energy cosmic particles by black-hole jets embedded in large-scale structures. Nat. Phys. 14, 396–398 (2018).
Blandford, R. D. & Begelman, M. C. On the fate of gas accreting at a low rate on to a black hole. MNRAS 303, L1–L5 (1999).
Ohsuga, K. & Mineshige, S. Global Structure of Three Distinct Accretion Flows and Outflows around Black Holes from Two-dimensional Radiation-magnetohydrodynamic Simulations. Apj 736, 2 (2011).
Narayan, R., SÄdowski, A., Penna, R. F. & Kulkarni, A. K. GRMHD simulations of magnetized advection-dominated accretion on a non-spinning black hole: role of outflows. MNRAS 426, 3241–3259 (2012).
Penna, R. F., Sadowski, A., Kulkarni, A. K. & Narayan, R. The Shakura-Sunyaev viscosity prescription with variable α (r). MNRAS 428, 2255–2274 (2013).
Ressler, S. M., Tchekhovskoy, A., Quataert, E. & Gammie, C. F. The disc-jet symbiosis emerges: modelling the emission of Sagittarius A* with electron thermodynamics. MNRAS 467, 3604–3619 (2017).
Chael, A., Rowan, M., Narayan, R., Johnson, M. & Sironi, L. The role of electron heating physics in images and variability of the Galactic Centre black hole Sagittarius A*. MNRAS 478, 5209–5229 (2018).
Kimura, S. S., Tomida, K. & Murase, K. Acceleration and escape processes of high-energy particles in turbulence inside hot accretion flows. MNRAS 485, 163–178 (2019).
Martin, R. G., Nixon, C. J., Pringle, J. E. & Livio, M. On the physical nature of accretion disc viscosity. N. A 70, 7–11 (2019).
Chael, A., Narayan, R. & Johnson, M. D. Two-temperature, Magnetically Arrested Disc simulations of the jet from the supermassive black hole in M87. MNRAS 486, 2873–2895 (2019).
Coppi, P. S. & Blandford, R. D. Reaction rates and energy distributions for elementary processes in relativistic pair plasmas. MNRAS 245, 453–507 (1990).
Pozdnyakov, L. A., Sobol, I. M. & Syunyaev, R. A. Comptonization and the shaping of X-ray source spectra - Monte Carlo calculations. Astrophys. Space Phys. Res. 2, 189–331 (1983).
Sironi, L. & Narayan, R. Electron Heating by the Ion Cyclotron Instability in Collisionless Accretion Flows. I. Compression-driven Instabilities and the Electron Heating Mechanism. Apj 800, 88 (2015).
Sironi, L. Electron Heating by the Ion Cyclotron Instability in Collisionless Accretion Flows. II. Electron Heating Efficiency as a Function of Flow Conditions. Apj 800, 89 (2015).
Zhdankin, V., Uzdensky, D. A., Werner, G. R. & Begelman, M. C. Electron and Ion Energization in Relativistic Plasma Turbulence. Phys. Rev. Lett. 122, 055101 (2019).
Kawazura, Y., Barnes, M. & Schekochihin, A. A. Thermal disequilibration of ions and electrons by collisionless plasma turbulence. Proc. Natl Acad. Sci. 116, 771–776 (2019).
Narayan, R. & Yi, I. Advection-dominated accretion: underfed black holes and neutron stars. Apj 452, 710 (1995).
Xie, F.-G. & Yuan, F. Radiative efficiency of hot accretion flows. MNRAS 427, 1580–1586 (2012).
Kimura, S. S., Toma, K., Suzuki, T. K. & Inutsuka, S.-i. Stochastic Particle Acceleration in Turbulence Generated by Magnetorotational Instability. Apj 822, 88 (2016).
Petrosian, V. Stochastic Acceleration by Turbulence. Space Sci. Rev. 173, 535–556 (2012).
Hoshino, M. Angular Momentum Transport and Particle Acceleration During Magnetorotational Instability in a Kinetic Accretion Disk. Phys. Rev. Lett. 114, 061101 (2015).
Kunz, M. W., Stone, J. M. & Quataert, E. Magnetorotational Turbulence and Dynamo in a Collisionless Plasma. Phys. Rev. Lett. 117, 235101 (2016).
Comisso, L. & Sironi, L. Particle Acceleration in Relativistic Plasma Turbulence. Phys. Rev. Lett. 121, 255101 (2018).
Balbus, S. A. & Hawley, J. F. A powerful local shear instability in weakly magnetized disks. I - Linear analysis. II - Nonlinear evolution. Apj 376, 214–233 (1991).
Balbus, S. A. & Hawley, J. F. Instability, turbulence, and enhanced transport in accretion disks. Rev. Mod. Phys. 70, 1–53 (1998).
Chang, J. S. & Cooper, G. A Practical Difference Scheme for Fokker-Planck Equations. J. Comput. Phys. 6, 1–16 (1970).
Park, B. T. & Petrosian, V. Fokker-Planck equations of stochastic acceleration: a study of numerical methods. Astrophys. J. Suppl. Ser. 103, 255 (1996).
Kafexhiu, E., Aharonian, F., Taylor, A. M. & Vila, G. S. Parametrization of gamma-ray production cross sections for p p interactions in a broad proton energy range from the kinematic threshold to PeV energies. Phys. Rev. D. 90, 123014 (2014).
Murase, K. & Nagataki, S. High energy neutrino emission and neutrino background from gamma-ray bursts in the internal shock model. Phys. Rev. D. 73, 063002 (2006).
Stepney, S. & Guilbert, P. W. Numerical FITS to important rates in high temperature astrophysical plasmas. MNRAS 204, 1269–1277 (1983).
Chodorowski, M. J., Zdziarski, A. A. & Sikora, M. Reaction rate and energy-loss rate for photopair production by relativistic nuclei. Apj 400, 181–185 (1992).
Becker, P. A., Le, T. & Dermer, C. D. Time-dependent Stochastic Particle Acceleration in Astrophysical Plasmas: exact solutions including Momentum-dependent Escape. Apj 647, 539–551 (2006).
Aartsen, M. G. et al. First Observation of PeV-Energy Neutrinos with IceCube. Phys. Rev. Lett. 111, 021103 (2013).
IceCube Collaboration. Evidence for High-Energy Extraterrestrial Neutrinos at the IceCube Detector. Science 342, 1242856 (2013).
Aartsen, M. G. et al. Evidence for Astrophysical Muon Neutrinos from the Northern Sky with IceCube. Phys. Rev. Lett. 115, 081102 (2015).
Kelner, S. R., Aharonian, F. A. & Bugayov, V. V. Energy spectra of gamma rays, electrons, and neutrinos produced at proton-proton interactions in the very high energy regime. Phys. Rev. D. 74, 034018 (2006).
Kimura, S. S., Murase, K., Mészáros, P. & Kiuchi, K. High-energy Neutrino Emission from Short Gamma-Ray Bursts: prospects for coincident detection with gravitational waves. Apj 848, L4 (2017).
Kimura, S. S. et al. Transejecta high-energy neutrino emission from binary neutron star mergers. Phys. Rev. D. 98, 043020 (2018).
Murase, K. New prospects for detecting high-energy neutrinos from nearby supernovae. Phys. Rev. D. 97, 081301 (2018).
Murase, K., Franckowiak, A., Maeda, K., Margutti, R. & Beacom, J. F. High-energy Emission from Interacting Supernovae: new Constraints on Cosmic-Ray Acceleration in Dense Circumstellar Environments. Apj 874, 80 (2019).
Marconi, A. et al. Local supermassive black holes, relics of active galactic nuclei and the X-ray background. MNRAS 351, 169–185 (2004).
Hopkins, P. F., Richards, G. T. & Hernquist, L. An Observational Determination of the Bolometric Quasar Luminosity Function. Apj 654, 731–753 (2007).
Lusso, E. et al. Bolometric luminosities and Eddington ratios of X-ray selected active galactic nuclei in the XMM-COSMOS survey. MNRAS 425, 623–640 (2012).
Liu, Z. et al. X-ray spectral properties of the AGN sample in the northern XMM-XXL field. MNRAS 459, 1602–1625 (2016).
Saikia, P. et al. 15-GHz radio emission from nearby low-luminosity active galactic nuclei. AA 616, A152 (2018).
Li, Y.-R., Ho, L. C. & Wang, J.-M. Cosmological Evolution of Supermassive Black Holes. I. Mass Function at 0 < z < 2. Apj 742, 33 (2011).
Ueda, Y., Akiyama, M., Hasinger, G., Miyaji, T. & Watson, M. G. Toward the Standard Population Synthesis Model of the X-Ray Background: Evolution of X-Ray Luminosity and Absorption Functions of Active Galactic Nuclei Including Compton-thick Populations. Apj 786, 104 (2014).
Murase, K., Inoue, Y. & Dermer, C. D. Diffuse neutrino intensity from the inner jets of active galactic nuclei: Impacts of external photon fields and the blazar sequence. Phys. Rev. D. 90, 023007 (2014).
Padovani, P. et al. The VLA Survey of Chandra Deep Field South. V. Evolution and Luminosity Functions of Sub-millijansky Radio Sources and the Issue of Radio Emission in Radio-quiet Active Galactic Nuclei. Apj 740, 20 (2011).
Ajello, M. et al. The Cosmic Evolution of Fermi BL Lacertae Objects. Apj 780, 73 (2014).
Franceschini, A. & Rodighiero, G. The extragalactic background light revisited and the cosmic photon-photon opacity. AA 603, A34 (2017).
Murase, K., Beacom, J. F. & Takami, H. Gamma-ray and neutrino backgrounds as probes of the high-energy universe: hints of cascades, general constraints, and implications for TeV searches. J. Cosmol. Astro-Part. Phys. 2012, 030 (2012).
Laha, R., Beacom, J. F., Dasgupta, B., Horiuchi, S. & Murase, K. Demystifying the PeV Cascades in IceCube: less (Energy) is More (Events). Phys. Rev. D88, 043009 (2013).
Aartsen, M. G. et al. Observation and Characterization of a Cosmic Muon Neutrino Flux from the Northern Hemisphere Using Six Years of IceCube Data. Apj 835, 3 (2016).
Aartsen, M. G. et al. IceCube-Gen2: a vision for the future of neutrino astronomy in Antarctica (2014).1412.5106.
Oh, K. et al. The 105-Month Swift-BAT All-sky Hard X-Ray Survey. ApJS 235, 4 (2018).
Walsh, J. L., van den Bosch, R. C. E., Barth, A. J. & Sarzi, M. A Stellar Dynamical Mass Measurement of the Black Hole in NGC 3998 from Keck Adaptive Optics Observations. Apj 753, 79 (2012).
Stecker, F. W., Shrader, C. R. & Malkan, M. A. The Extragalactic Gamma-Ray Background from Core-dominated Radio Galaxies. Apj 879, 68 (2019).
Ajello, M. et al. Cosmic X-Ray Background and Earth Albedo Spectra with Swift BAT. Apj 689, 666–677 (2008).
Watanabe, K. et al. The Cosmic γ-ray Background from supernovae. Vol. 410 (eds Dermer, C. D., Strickman, M. S. & Kurfess, J. D.) Proceedings of the Fourth Compton Symposium. 1223–1227 (American Institute of Physics Conference Series, 1997).
Fukada, Y. et al. Energy spectrum of diffuse component of cosmic soft gamma rays. Nature 254, 398 (1975).
Stettner, J. Measurement of the diffuse astrophysical muon-neutrino spectrum with ten years of IceCube data. In 36th International Cosmic Ray Conference (ICRC2019), Vol. 36 of International Cosmic Ray Conference, 1017 (2019).1908.09551.
Acknowledgements
This work is supported in part by the IGC postdoctoral fellowship program, JSPS Oversea Research Fellowship, JSPS Research Fellowship, KAKENHI No. 19J00198 (S.S.K.), the Alfred P. Sloan Foundation, NSF Grant No. AST-1908689, and KAKENHI No. 20H01901 and No. 20H05852 (K.M.), and the Eberly Foundation (P.M.).
Author information
Authors and Affiliations
Contributions
K.M. provided the basic idea of this study, and S.S.K. constructed the model and performed all the calculations. K.M. developed the calculation code for the proton-induced cascade emission. All the authors (S.S.K., K.M., and P.M.) discussed the implications of the results, and contributed to the paper.
Corresponding author
Ethics declarations
Competing interests
The authors declare no competing interests.
Additional information
Peer review information Nature Communications thanks Ramesh Narayan and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Supplementary information
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Kimura, S.S., Murase, K. & Mészáros, P. Soft gamma rays from low accreting supermassive black holes and connection to energetic neutrinos. Nat Commun 12, 5615 (2021). https://doi.org/10.1038/s41467-021-25111-7
Received:
Accepted:
Published:
DOI: https://doi.org/10.1038/s41467-021-25111-7
This article is cited by
-
Sagittarius \(\textrm{A}^{\star }\) as a Plausible Source Candidate for PeV Neutrinos
Brazilian Journal of Physics (2024)
Comments
By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.
