Neutrinos and Gamma Rays from Galaxy Clusters Constrained by the Upper Limits of IceCube

Clusters of galaxies possess the capability to accelerate cosmic rays (CRs) to very high energy up to ∼1018 eV due to their large size and magnetic field strength, which favor CR confinement for cosmological times. During their confinement, they can produce neutrinos and γ-rays out of interactions with the background gas and photon fields. In recent work, Hussain et al. have conducted three-dimensional cosmological magnetohydrodynamical simulations of the turbulent intracluster medium combined with multidimensional Monte Carlo simulations of CR propagation for redshifts ranging from z ∼ 5 to z = 0 to study the multimessenger emission from these sources. They found that when CRs with a spectral index in the range 1.5–2.5 and cutoff energy Emax=1016–1017 eV are injected into the system, they make significant contributions to the diffuse background emission of both neutrinos and γ-rays. In this work, we have revisited this model and undertaken further constraints on the parametric space. This was achieved by incorporating the recently established upper limits on neutrino emission from galaxy clusters, as obtained by the IceCube experiment. We find that for CRs injected with spectral indices in the range 2.0–2.5, cutoff energy Emax=1016–1017 eV, and power corresponding to (0.1–1)% of the cluster luminosity, our neutrino flux aligns with the upper limits estimated by IceCube. Additionally, the resulting contribution from clusters to the diffuse γ-ray background remains significant with values of the order of ∼10−5 MeV cm−2 s−1 sr−1 at energies above 500 GeV.


INTRODUCTION
The diffuse neutrino (Aartsen et al. 2015a) and gamma ray (Ackermann et al. 2015) backgrounds provide a unique prospective of the high-energy Universe but their origin is still debated.It is closely related with ultra-high-energy cosmic rays (UHECRs) (Abu-Zayyad et al. 2013;Aab et al. 2017).The most plausible scenario is that the neutrino and gamma-rays are produced by interactions of UHECRs with the background gas and photon fields in astrophysical environments.Several sources embedded in galaxy clusters arise as candidates for the production of very high energy CRs including active galactic nuclei and starburst galaxies.The turbulent intracluster medium (ICM) is also believed to be particularly suitable to accelerate and confine CRs up to energies ∼ 10 18 eV (Inoue et al. 2007;Wiener & Zweibel 2019;Alves Batista et al. 2019).
Recently, the IceCube collaboration (Abbasi et al. 2022) performed a stacking analysis for 1094 galaxy clusters of masses ≳ 10 14 M ⊙ and up to redshift z ≲ 1.0 taken from the 2015 PLANCK survey (Ade et al. 2016).To complete the catalogue, they calculated the distribution of galaxy clusters by drawing ∼ 10 5 samples with mass range 10 14 ≲ M/M ⊙ ≲ 10 15 extending to redshift z ≲ 2.0, using the Tinker−2010 halo mass function (Tinker et al. 2010).Based on these samples, they estimated upper limits according to which the contribution from the clusters to the observed diffuse neutrino background is at most ∼ 4.6% at 100 TeV for the most realistic scenario (see for details Raghunathan et al. 2022).
Several recent studies (see e.g., Murase et al. 2013;Fang & Olinto 2016;Fang & Murase 2018;Hussain et al. 2021Hussain et al. , 2023) ) have predicted that clusters of galaxies can contribute to a sizeable percentage to the diffuse neutrino and γ−ray background.In particular, Hussain et al. (2021Hussain et al. ( , 2023) used the most detailed numerical approach to date combining three-dimensional (3D) cosmological magneto-hydro-dynamic (MHD) simulations with Monte-Carlo simulations of CR propagation and cascading and predicted that clusters can potentially contribute to a fairly large fraction to the diffuse neutrino and gamma-ray background observed by IceCube (Aartsen et al. 2015a) and Fermi-LAT (Ackermann et al. 2015), respectively.However, the aforementioned upper limit reported by the IceCube (Abbasi et al. 2022) excludes part of the parametric space they considered in their analysis of the neutrino flux.
In this work, we constrain the parametric space considered in Hussain et al. (2021Hussain et al. ( , 2023) ) employing the latest upper limits reported by the IceCube for the neutrino flux of galaxy clusters and derive new limits to the diffuse γ−ray flux.We find that these new constraints eliminate the harder CR spectral indices α ≲ 2.0, but still predict a substantial contribution from the clusters to the very high energy range of the γ−ray flux.
In section 2 we describe the methodology for the calculation of the fluxes of neutrinos and γ−rays in galaxy clusters, in section 3 we present the results, and section 4 is dedicated to the discussion of our results and the conclusions.

METHODOLOGY
We have employed here the same set of simulations as in Hussain et al. (2021Hussain et al. ( , 2023)).The ICM was modeled with 3D-MHD cosmological simulations (Dolag et al. 2005) up to redshift z ≲ 5.0 taking into account the non-uniform distribution evolution of the magnetic field, gas density, and temperature.CRs were proppagated in the ICM and intergalactic medium (IGM) to produce neutrinos and γ−rays, employing the Monte Carlo code CRpropa (Batista et al. 2016).To calculate the neutrino flux CRs were injected in the energy range 10 14 ≤ E/eV ≤ 10 19 because we are interested in neutrino energies above TeV.On the other hand, for γ−rays, since we are interested in energies above 10 GeV, CRs were injected with energies 10 11 ≤ E/eV ≤ 10 19 .However, in order to normalize the total energy to the luminosity of the clusters, it was considered the whole energy range starting from 1 GeV.Also, in order to account for different sources of acceleration, CRs were injected at different locations within the clusters: in the center, at a radius ∼ 300 kpc, and in the outskirt ∼ 1 Mpc.Obviously, the dominant contribution to γ−ray and neutrino fluxes come from CRs sources located in the center (see Hussain et al. 2021Hussain et al. , 2023, for details), for details).In these previous studies, it has been considered that 1% of the cluster luminosity (L C ) goes to the CRs.
The simulations have two steps.In the first step, CRs are propagated inside the clusters and γ−rays and neutrinos are collected at the edge of them.All relevant CR interactions were considered during their propagation inside clusters, namely, proton-proton (pp) interactions, photopion production, Bethe-Heitler pair production, pair production, and inverse Compton scattering (ICS).Energy losses due to the expansion of the universe and the synchrotron emission were also considered, but the photons produced in these processes are below the energy range of interest of this work.In the second step, the γ−rays collected at the boundary of the clusters were propagated through the intergalactic medium (IGM) across the redshift interval.During this propagation, the electromagnetic cascade processes including (single, double, and triplet) pair production and ICS were also accounted for.The effect of the IGM magnetic field was neglected in this propagation step since it does not significantly affect the γ−ray flux above 10 GeV (Hussain et al. 2023).
As in Hussain et al. (2021) and Hussain et al. (2023), the integrated neutrino and γ−ray fluxes (Φ) from all clusters in the mass range 10 12 ≲ M/M ⊙ ≲ 2 × 10 15 and z ≤ 5.0, are obtained from: where dN/dM is the number of clusters per mass interval calculated from the MHD simulation, g(E obs , E, z) accounts for the interactions of gamma rays in the ICM and the IGM, ψ ev (z) is a function that describes the cosmological evolution of the emissivity of the CR sources (see e.g., Alves Batista et al. 2019), the quantity E 2 d Ṅ /dE denotes the neutrino or γ−ray power spectrum obtained from the simulations, d L is the luminosity distance, and f (M ) is a factor that accounts for stellar and AGN feedback (Planelles et al. 2014).For detailed calculations, we refer to Hussain et al. (2023) and Hussain et al. (2021).

RESULTS
In Fig. 1, we show the flux of neutrinos from clusters of galaxies, considering that their contribution is constrained by the upper limits recently reported by the IceCube (Abbasi et al. 2022).As stressed, in the previous works Hussain et al. (2021Hussain et al. ( , 2023) ) had assumed that 1% of each cluster luminosity (L C ) goes to CRs, and considered a range for the CR spectral index 1.5 ≤ α ≤ 2.5 and maximum energy 10 16 ≤ E max /eV ≤ 5 × 10 17 .Here, we find that the CR parameters which are more suitable to fulfill the new IceCube limits are the following: 2.0 ≤ α ≤ 2.5, E max = 10 16 − 10 17 eV, and luminosities in the range (0.1 − 1.0)% L C .The figure shows three bands all constrained by this interval of α values, but with three different luminosities 1%L C , 0.5%L C , and 0.1%L C , from light to dark-blue, respectively.Also, while the light-blue and blue bands have E max = 10 17 eV, the dark-blue band shows the flux for E max in the range 10 16 −10 17 eV.In each of these bands, the larger the value of the index α, the smaller the flux is.We see that the band with 0.1% L C is the most constrained one by the IceCube limits.In particular, for α ≥ 2.3, this band falls entirely below the IceCube limits.Also, decreasing the value of E max results in larger flux at smaller neutrino energies.This explains why the light-blue and blue bands produce fluxes around or below the IceCube limits only for E max ≃ 10 17 eV.For E max ≃ 10 16 eV, the peak of the flux in the figure increases by almost an order of magnitude at neutrino energies ∼ 10 13 eV.
Fig. 2 summarizes our results showing both the γ−ray and the neutrino fluxes for the same parameters of Fig. 1.Besides the IceCube upper limits for clusters and diffuse neutrino background, it also depicts the diffuse γ−ray background (DGRB) data from Fermi-LAT (Ackermann et al. 2015) and the upper limits from HAWC (Albert et al. 2022).The γ−ray flux is given by the light-green, dark-green and olive-green bands, which are the counterparts of the light-blue, blue and dark-blue neutrino bands, respectively.We see that the olive-green, which is the most constrained band by the upper limits of the IceCube, falls below the Fermi data for γ−rays.
In Figure 3 we compare the results obtained in Fig. 2 with the earlier results of (Hussain et al. 2021) and (Hussain et al. 2023).We note that the parametric space considered in these works cannot be excluded entirely by the upper limits of the IceCube (Abbasi et al. 2022).10 12 10 13 10 14 10 15 10 16 10 17  The green bands give the diffuse flux of γ−ray obtained for the same CR parametric space as in Fig. 1 that suits the IceCube upper limits.The γ−ray fluxes in the light-green, dark-green and olive-green bands, have the same parameters as the neutrino fluxes in the light-blue, blue and dark-blue, respectively.The diffuse neutrino background upper limits reported by the IceCube (Aartsen et al. 2015a,b), the DGRB observed by Fermi-LAT (Ackermann et al. 2015), and the upper limits for the DGRB from HAWC (Albert et al. 2022) are also depicted.

DISCUSSION AND CONCLUSIONS
We have computed the flux of neutrinos and γ−rays from the entire population of galaxy clusters using the most detailed numerical simulations to date, considering 3D-MHD cosmological simulations of galaxy clusters up to redshift z ≲ 5.0, as in Hussain et al. (2021Hussain et al. ( , 2023)).According to these authors, clusters can contribute to a sizeable percentage of up to 100% to diffuse neutrino background, depending on the parameters adopted for the CR spectrum.However, the IceCube collaboration (Abbasi et al. 2022) has reported that this contribution cannot exceed (9−13)%.Evaluating upper limits for the neutrino flux of the clusters, this collaboration concluded that these new constraints would exclude the Hussain et al. ( 2021) models for hard CR spectral indices α < 2.0.Our present results indicate that, in fact, the new Icecube limits point to harder spectral indices α ≥ 2.0.Nevertheless, these new results are entirely compatible with the parametric space explored in Hussain et al. (2021), which also included values of α ≥ 2.0 (see Figure 3).
In particular, for α = 2.5 and E max = 10 17 eV the neutrino flux obtained in (Hussain et al. 2021) is below the upperlimits of the IceCube and decreases approximately by an order of magnitude if we assume CR luminosity 0.1% L C instead of 1% L C .
We have also computed the γ−ray flux constrained by these IceCube upper limits and found that it decreases by roughly an order of magnitude in comparison with Hussain et al. (2023), for CR luminosity ∼ 0.1% L C .Despite that, the contribution of clusters to the DGRB is still substantial above 500 GeV (Figure 3).Moreover, the flux falls within the sensitivity ranges of the HAWC, LHAASO, and the upcoming CTA observatories, suggesting the possibility for direct detections of γ−rays from clusters of galaxies (Figure 4).As in the case of the neutrinos, the γ−ray flux goes for values lower than the Fermi-LAT observations (Ackermann et al. 2015) for α ∼ 2.5, and reduces even more, if we consider a CR luminosity 0.1% L C , rather than 1% L C .While discrete source categories such as AGNs (Di Mauro 10 10 10 11 10 12 10 13 10 14 10 15 10 16 10 17   (Hussain et al. 2023(Hussain et al. , 2021)).These are compared with the results of Figure 2, i.e., the flux of neutrinos (blue bands) obtained with the new parametric space constrained by the upper limits estimated by the IceCube (Abbasi et al. 2022).The corresponding γ−ray flux (green bands) for the same parameters is also shown.The diffuse neutrino background reported earlier by the IceCube (Aartsen et al. 2015a,b) is also depicted.The DGRB observed by Fermi-LAT (Ackermann et al. 2015) and the upper limits for the DGRB from HAWC (Albert et al. 2022) and the CASA-MIA (Chantell et al. 1997a) are also depicted.et al. 2013;Ajello et al. 2015) and star-forming galaxies (Roth et al. 2021) can make an important contribution to the DGRB for energy levels below the TeV range, our results highlight that the combined gamma-ray flux originating from clusters can surpass the combined impact of individual classes of unresolved sources for energies exceeding 500 GeV.This aligns with the findings of Hussain et al. (2023).
We should emphasize that our aim here was not to fit either IceCube upper limits (Abbasi et al. 2022) or Fermi-LAT data (Ackermann et al. 2015).Instead, we have only compared our evaluation of the integrated flux of neutrinos and γ− rays for the entire population of clusters with those observations.Our estimations are dependent on the parametric space and so does the IceCube results.To obtain their upper-limits, Abbasi et al. (2022) considered sources with masses in the range 10 14 ≤ M/M ⊙ ≤ 10 15 up to redshift z ≤ 2.0.Our analysis, on the other hand, has considered the entire mass range of clusters 10 12 ≤ M/M ⊙ ≤ 2 × 10 15 , and redshifts z ≤ 5.0.Though major contribution to neutrino and γ−ray background comes from the nearby sources (z ≤ 1), and more massive clusters are more frequent, the contribution from clusters in the mass range 10 13 ≤ M/M ⊙ ≤ 10 14 is not negligible (Hussain et al. 2021(Hussain et al. , 2023)).Therefore, including this mass interval might change the upper limits estimated by IceCube (Abbasi et al. 2022).
Finally, we would like to stress that some specific parameters may have the potential to influence our simulations.For instance, a mixed composition of cosmic rays (CRs) and the distribution of CR sources within the clusters could lead to variations in our results.If we were to consider a CR composition involving heavy elements like iron (Fe), it could potentially alter our conclusions.Similarly, slight adjustments in the assumed distributions of CR sources within the structures (see Section 2) might also have an impact on our outcomes.Exploring further these effects is a direction we plan to pursue in the future.Furthermore, it is worth noting that we have not accounted for the influence of the uncertain diffuse magnetic field outside the clusters, which could introduce further minor variations in the gamma-ray observations.Flux of γ−rays from the entire population of galaxy clusters as in Figures 2 and 3.It is compared with sensitivity curves from different experiments, namely, from Hawc for point-like sources (Abeysekara et al. 2013), LHAASO (Lhaaso Collaboration et al. 2016), andCTA (CTA Consortium et al. 2018).Also depicted are the upper limits for DGRB from HAWC (Albert et al. 2022) and CASA-MIA (Chantell et al. 1997b).

Figure 1 .EFigure 2 .
Figure1.Neutrino flux from the entire population of clusters up to redshift z ≤ 5.0.From top to bottom blue bands correspond to 1%, 0.5%, and 0.1% of the cluster luminosity (LC ), respectively.Light-blue and blue bands represent the flux for 2.0 ≤ α ≤ 2.5 and Emax = 10 17 eV.The dark-blue band shows the flux for the same spectral indices, but Emax ranges from 10 16 − 10 17 eV.The IceCube upper limits for clusters(Abbasi et al. 2022) as well as the diffuse neutrino background(Aartsen  et al. 2015b,a)  are also depicted.

Figure 3 .
Figure 3. Neutrino (brown band) and γ−ray (gray band) fluxes form the entire population of galaxy clusters as obtained in(Hussain et al. 2023(Hussain et al. , 2021)).These are compared with the results of Figure2, i.e., the flux of neutrinos (blue bands) obtained with the new parametric space constrained by the upper limits estimated by the IceCube(Abbasi et al. 2022).The corresponding γ−ray flux (green bands) for the same parameters is also shown.The diffuse neutrino background reported earlier by the IceCube(Aartsen et al. 2015a,b) is also depicted.The DGRB observed by Fermi-LAT(Ackermann et al. 2015) and the upper limits for the DGRB from HAWC(Albert et al. 2022) and the CASA-MIA(Chantell et al. 1997a) are also depicted.
Figure 4.Flux of γ−rays from the entire population of galaxy clusters as in Figures2 and 3.It is compared with sensitivity curves from different experiments, namely, from Hawc for point-like sources(Abeysekara et al. 2013),LHAASO  (Lhaaso Collaboration et al. 2016), andCTA (CTA Consortium et al. 2018).Also depicted are the upper limits for DGRB from HAWC(Albert et al. 2022) and CASA-MIA(Chantell et al. 1997b).