Prepared for submission to JCAP

Constraining decaying very heavy

dark matter from galaxy clusters
with 14 year Fermi-LAT data
arXiv:2308.00589v2 [astro-ph.HE] 7 Feb 2024

Deheng Songa Kohta Muraseb,c,d,e,a Ali Kheirandishf,g

a Center for Gravitational Physics and Quantum Information, Yukawa Institute for Theoret-
ical Physics, Kyoto University, Kyoto 606-8502, Japan
b Department of Physics, The Pennsylvania State University, University Park, Pennsylvania

16802, USA
c Department of Astronomy and Astrophysics, The Pennsylvania State University, University

Park, Pennsylvania 16802, USA

d Center for Multimessenger Astrophysics, The Pennsylvania State University, University

Park, Pennsylvania 16802, USA

e School of Natural Sciences, Institute for Advanced Study, Princeton, New Jersey 08540,

USA
f Department of Physics & Astronomy, University of Nevada, Las Vegas, NV, 89154, USA
g Nevada Center for Astrophysics, University of Nevada, Las Vegas, NV 89154, USA

E-mail: songdeheng@yukawa.kyoto-u.ac.jp, murase@psu.edu,

ali.kheirandish@unlv.edu

Abstract. Galaxy clusters are promising targets for indirect detection of dark matter
thanks to the large dark matter content. Using 14 years of Fermi-LAT data from seven
nearby galaxy clusters, we obtain constraints on the lifetime of decaying very heavy dark
matter particles with masses ranging from 103 GeV to 1016 GeV. We consider a variety of
decaying channels and calculate prompt gamma rays and electrons/positrons from the dark
matter. Furthermore, we take into account electromagnetic cascades induced by the primary
gamma rays and electrons/positrons, and search for the resulting gamma-ray signals from
the directions of the galaxy clusters. We adopt a Navarro-Frenk-White profile of the dark
matter halos, and use the profile likelihood method to set lower limits on the dark matter
lifetime at a 95% confidence level. Our results are competitive with those obtained through
other gamma-ray observations of galaxy clusters and provide complementary constraints to
existing indirect searches for decaying very heavy dark matter.
Contents

1 Introduction 1

2 Signals of VHDM decay 2

3 Samples of galaxy clusters 4

4 Data analysis 5

5 Results and discussion 7

5.1 Combined analysis 8
5.2 Systematic uncertainties 10

6 Summary 11

A Electromagnetic cascades inside galaxy clusters 13

1 Introduction

Dark matter constitutes the majority of matter in the Universe [1]. Nevertheless, we are
still unaware of the particle nature of dark matter. Weakly Interacting Massive Particles
(WIMPs) are among the most widely accepted candidates for dark matter [2]. However,
direct detection experiments have placed significant constraints on WIMPs across a range of
masses and they are particularly sensitive to those with masses between approximately 0.01-1
TeV. For WIMPs with masses above 1 TeV, the parameter space remains less constrained,
due in part to the lower expected interaction rates. For thermal relics, the unitarity bounds
on the WIMP annihilation cross section set an upper limit of dark matter mass at around 100
TeV [3]. However, more massive dark matter models are allowed especially if the dark matter
is produced non-thermally, for instance, when a period of early matter domination occurs
before the Big Bang Nucleosynthesis [4–6], or when dark matter is produced in a hidden
sector [7–12].
Indirect detection, which looks for standard model particles originating from dark matter
in the cosmos, is a viable approach for searching dark matter. Thermally produced dark
matter particles can still self-annihilate to Standard Model particles in the Universe today.
They can also have finite lifetimes, and decay into Standard Model products. High-energy
gamma-ray and neutrino telescopes can detect gamma rays and neutrinos produced through
dark matter annihilation or decay, which can travel long distances in the Universe. The Fermi
Large Area Telescope (Fermi-LAT) has been extensively used to search for indirect signals
from dark matter. Towards the Galactic center, Fermi-LAT has identified an excess in the
GeV energies, which may be caused by WIMP dark matter annihilation or a population of
unresolved millisecond pulsars [13–21]. Although the origin of the gamma-ray excess is still
debated, the Galactic center region still places meaningful constraints on the dark matter
annihilation cross section [22]. Fermi-LAT is also used to search for dark matter in the range
of ∼ GeV to TeV from dwarf galaxies [23–32], galaxy clusters [33–40], and the diffuse isotropic
gamma-ray background [41–47]. Galaxy clusters are the largest gravitationally-bound systems

–1–
in the Universe and host enormous dark matter halos. TeV gamma-ray telescopes, such
as the Very Energetic Radiation Imaging Telescope Array System (VERITAS), the Major
Atmospheric Gamma Imaging Cherenkov (MAGIC) and the High Altitude Water Cherenkov
(HAWC) observatory, have searched for dark matter signals from different galaxy clusters
and have placed stringent limits on dark matter annihilation and decay [48–53]. Future TeV
gamma-ray telescopes are expected to extend the search for heavy dark matter. [54, 55].
Dedicated analyses of multi-messenger astrophysical data have set strongest constraints
on the lifetime of very heavy dark matter (VHDM) with masses up to EeV or even the GUT
(Grand Unified Theory) scale (∼1016 GeV) [43, 56–60]. The IceCube Neutrino Observa-
tory [61–63] and the large high altitude air shower observatory (LHAASO) [64] also searched
for decaying dark matter with masses above PeV and set competitive limits. Future neutrino
experiments will extend the search for VHDM to ≳ 1012 GeV [65].
In this work, we search for signals from decaying VHDM in galaxy clusters using 14
year gamma-ray data from Fermi-LAT. We select seven nearby galaxy clusters based on their
dark matter content [66]. We examine a wide range of dark matter masses, from ∼ 103
GeV to the GUT scale (∼ 1016 GeV). As Fermi-LAT is only sensitive to gamma rays with
energies up to ∼ TeV, it cannot detect prompt gamma rays from most of the mass range
we consider. However, the high-energy gamma rays and electrons/positrons (e± ) produced
by dark matter decay interact with the radiation and magnetic fields inside galaxy clusters,
causing electromagnetic cascades. The secondary gamma rays generated through inverse-
Compton and synchrotron processes are cascaded down to the energy range that Fermi-LAT
can detect [36, 43, 67–69]. We search for gamma-ray signals from dark matter decay in galaxy
clusters taking into account the electromagnetic cascades and place limits on the lifetime of
dark matter.
The rest of this paper is organized as follows. In section 2, we describe the expected
signals from VHDM decay, including the electromagnetic cascades. In section 3, we summarize
the selection of galaxy clusters. We provide the details about our data analysis procedures in
section 4. Subsequently, we present and discuss our results in section 5 and summarize the
paper in section 6.

2 Signals of VHDM decay

The expected intensity of generated gamma rays from decaying dark matter for given gamma-
ray energy Eγ = Eγ′ /(1 + z) and direction Ω on the sky is
1 1 dSγ
Z
(gen)
Iγ (Eγ , Ω) = dlρχ (r), (2.1)
4π τχ mχ dEγ′
where τχ and mχ are the lifetime and mass of dark matter particle χ, and dSγ /dEγ′ is the
gamma-ray spectrum from dark matter decay to Standard Model particles, which is normal-
ized via
dSi
Z
Σi dEi′ Ei′ ′ = mχ , (2.2)
dEi
where the summation is for all particle species. Here, ρχ is the smooth density distribution
of dark matter, which is integrated over the solid angle Ω and the comoving line-of-sight
distance l. For dark matter halos of galaxy clusters, we use a Navarro–Frenk–White (NFW)
profile [70] in comoving coordinates as the following
ρs
ρχ (r) = , (2.3)
(r/rs )(1 + r/rs )2

–2–
where the scale density ρs and scale radius rs vary from cluster to cluster.
The generated gamma-ray energy flux within a cone with a sold angle of Ω is given by

1 1 ′ 2 dSγ 4πd2c dec

Z
(gen)
Eγ FEγ (Ω) = dΩEγ2 Iγ(gen) (Ω) = E × J , (2.4)
4πd2L τχ mχ γ dEγ′ 4π

where dL and dc are luminosity distance and comoving distance, respectively, and the angle-
integrated J-factor J dec (Ω) for decaying dark matter is given by [66]
Z Z Z
J (Ω) = dΩ dlρχ (r) = dΩJ dec (Ω),
dec
(2.5)

where J dec is the angle-dependent J-factor of the cluster. In this work, we use the public
software CLUMPY [71–73] to calculate J dec following the NFW profile. Note that the integration
over the virial radius gives J dec ≈ Mvir /d2c , where Mvir is the virial radius, and the normalized
J-factor becomes [36]
2π dθθEγ2 Iγ
R
d2
Jˆ (θ) =
dec
≈ c J dec (θ), (2.6)
Eγ FEγ Mvir
where θ = tan−1 (r/dc ).
In this study, we consider VHDM in the mass range between 103 GeV to 1016 GeV
assuming that each dark matter particle decays into a pair of Standard Model particles. We
consider seven decaying channels: bb̄, tt̄, e+ e− , µ+ µ− , τ + τ − , Z 0 Z 0 , and W + W − . To calculate
the prompt spectra, we use HDMSpectra [74], which incorporates all relevant electroweak
interactions. This contribution is crucial for VHDM.
Equation 2.4 is for the flux of generated gamma rays, but for gamma rays from VHDM
we must further take into account contributions from electromagnetic cascades that are de-
veloped inside galaxy clusters following the treatments described in Refs. [36, 43]. We solve
the Boltzmann equations of photons and e± considering e± pair creation, inverse-Compton
scattering and synchrotron radiation. Prompt e± from VHDM decay are trapped in the
galaxy clusters and lose energies by interacting with the cosmic microwave background (CMB)
and extragalactic background light (EBL) as well as intra-cluster magnetic fields, leading to
inverse-Compton and synchrotron emission. When dark matter is sufficiently heavy (mχ ≳
10 TeV), prompt gamma rays are attenuated by radiation fields in the galaxy clusters and
secondary e± are generated. These secondary e± also produce inverse-Compton and syn-
chrotron emissions. The final gamma-ray spectra from dark matter decay in galaxy clusters
combine the remaining prompt and the total secondary gamma rays. We include the details
about intra-cluster electromagnetic cascades in Appendix A. After the gamma rays leave the
galaxy clusters and enter the intergalactic space, they are further attenuated by the EBL. We
calculate such attenuation, but ignore intergalactic cascade emission as in previous works [36].
This is because they are expected to form a giant pair halo with a size of λγγ /d [75–77], where
λγγ is the mean free path to the two-photon annihilation process and d is the source distance.
The mean free path of 100 TeV gamma rays is around λγγ ∼ 3 Mpc [78], and higher-energy
gamma rays are cascaded inside clusters. The mean free path of 10 TeV gamma rays is
λγγ ∼ 100 − 200 Mpc, which is larger than the distance to nearby clusters. Therefore most of
the intergalactic cascade emission contributes to diffuse or even quasi-isotropic emission, and
it is reasonable and more conservative to ignore the contribution for the present analyses. See
Ref. [79] for the total contribution from a source and extended/diffuse intergalactic cascade
emission.

–3–
 10−8 10−8
bb̄ τ +τ − bb̄ τ +τ −

Ee± FE ± (Ω) [GeV cm−2 s−1 ]

Z 0Z 0 Z 0Z 0
Eγ FEγ (Ω) [GeV cm−2 s−1 ]
tt̄ tt̄
e+ e− W +W − e+ e− W +W −
10−9 10−9
µ+ µ− µ+ µ−

mχ = 104 GeV mχ = 104 GeV

−10 −10
10 10

(gen)
(gen)

e
10−11 10−11

10−12 −1 10−12 −1
10 100 101 102 103 104 10 100 101 102 103 104
Eγ [GeV] Ee± [GeV]

Figure 1. Prompt spectra of gamma rays from Figure 2. Same as figure 1, but of e± .
dark matter decay. Here, we show the case of
mχ = 104 GeV and show seven decaying channels:
bb̄, tt̄, e+ e− , µ+ µ− , τ + τ − , Z 0 Z 0 , and W + W − .
The grey area shows the Fermi-LAT energy range.
The fluxes are normalized by J dec (Ω) of Virgo.

Figure 1 and figure 2 show the prompt gamma-ray and e± spectra resulting from dark
matter decay, respectively. The gamma-ray and e± fluxes are normalized by assuming τχ =
1028 s and J dec (Ω) of the Virgo cluster. The grey area denotes the energy range to which
Fermi-LAT is sensitive. These figures assume mχ = 104 GeV and include seven channels. For
the bb̄, tt̄, Z 0 Z 0 , and W + W − channels, the decaying final states are largely hadronic, resulting
in similar gamma-ray spectra at low energies. The gamma-ray spectra of the e+ e− and µ+ µ−
channel is very hard and peaks at ∼ mχ /2 as they are dominated by the final-state radiation.
Therefore, including electromagnetic cascades is more critical for the e+ e− and µ+ µ− channel.
The τ + τ − channel is intermediate since the leptonic and hadronic decays are comparable [66].
Figure 3 shows the expected gamma-ray spectra from dark matter decay when electromagnetic
cascades are considered, for mχ = 104 GeV. We assume that the magnetic field in the cluster
is Bcl = 0.3 µG [80] and the size of the emission region is set to Rcl = 3 Mpc. We discuss the
impact of these parameters later in section 5. Comparing to figure 1, the spectra including
electromagnetic cascades extend to lower energies due to contributions from inverse-Compton
and synchrotron radiation processes. The spectra are also undulating when transitioning from
prompt to secondary gamma rays. When dark matter is sufficiently heavy, Fermi-LAT will
only be able to detect secondary gamma rays, as demonstrated in figure 4 for mχ = 1012 GeV.
In this case, gamma-ray spectra in the Fermi-LAT energy range are dominated by synchrotron
radiation and are broad and smooth.

3 Samples of galaxy clusters

We choose the galaxy clusters from the catalog presented in Ref. [66]. This catalog was
initially developed based on [81, 82]. Here, we select seven clusters with the largest J dec. .
Table 1 summarizes the properties of the clusters, including locations, redshifts (z), and dark

–4–
 10−8 10−8
bb̄ τ +τ − bb̄ τ +τ −
tt̄ Z 0Z 0 tt̄ Z 0Z 0
Eγ FEγ (Ω) [GeV cm−2 s−1 ]

Eγ FEγ (Ω) [GeV cm−2 s−1 ]

e+ e− W +W − e+ e− W +W −
10−9 10−9
µ+ µ− µ+ µ−

mχ = 104 GeV mχ = 1012 GeV

10−10 10−10

10−11 10−11

10−12 −1 10−12 −1
10 100 101 102 103 104 10 100 101 102 103 104
Eγ [GeV] Eγ [GeV]

Figure 3. Gamma-ray spectra from dark matter Figure 4. Same as figure 3, but for mχ = 1012
decay when electromagnetic cascades are consid- GeV.
ered.

Cluster l b z × 103 cvir rvir θvir log10 Mvir log10 J dec.

[deg] [deg] [kpc] [deg] [M⊙ ] [GeV cm−2 sr]
Virgo 283.94 74.52 3.58 6.36 2023.32 7.28 14.66 20.41
Centaurus 302.22 21.65 8.44 6.40 1971.53 3.02 14.62 19.63
Norma 325.29 −7.21 17.07 5.60 2842.43 2.16 15.10 19.51
Perseus 150.58 −13.26 17.62 5.63 2796.52 2.06 15.08 19.46
Coma 57.20 87.89 24.45 5.40 3132.32 1.67 15.23 19.33
Hydra 269.55 26.41 10.87 6.74 1687.23 2.01 14.42 19.21
Fornax 239.98 −56.68 4.17 8.47 873.63 2.71 13.56 19.17

Table 1. Sample of galaxy clusters. The parameters are taken from Ref. [66].

matter halo parameters. The selected clusters have virial masses of Mvir ∼ 1014 − 1015 M⊙ .
The virial radius rvir and the concentration parameter cvir are related to the scale radius of
the dark matter halo, rs = rvir /cvir . We also report the angular size of the dark matter halo
at the virial radius, which is given by

θvir = tan−1 (rvir /dc ). (3.1)

The Virgo cluster has the largest J dec. , which is ∼ 1020 GeV cm−2 sr, and is the closest
cluster with a redshift z = 3.58 × 10−3 . Therefore, the dark matter halo of Virgo extends
to θvir = 7.28◦ . The remaining clusters are: Centaurus, Norma, Perseus, Coma, Hydra, and
Fornax. They typically have J dec. ∼ 1019 GeV cm−2 sr, and θvir ∼ 2◦ − 3◦ .

4 Data analysis
We analyze 14 year Fermi-LAT Pass8 data collected between Aug 4 2008 and Aug 4 2022.
We select photon events from the ULTRACLEANVETO class and include both FRONT and BACK
types. We apply a quality filter “DATA_QUAL>0 && LAT_CONFIG==1” and limit the maximum

–5–
zenith angle to 90◦ . For each galaxy cluster, the region of interest (ROI) is a 20◦ by 20◦
square centered at the cluster’s location with an angular resolution of 0.1◦ . The energy range
of the gamma rays is from 100 MeV to 1 TeV, with 5 logarithmic energy bins per decade.
We fit the Fermi data using the open-source Python package fermipy [83]. For each
cluster, we include the Galactic interstellar emission model gll_iem_v07.fits and isotropic
diffuse template iso_P8R3_ULTRACLEANVETO_V3_v1.txt as the background models. Also in-
cluded are the point sources resolved in the LAT 12-year source catalog (4FGL-DR3) [84, 85]
that are up to 5◦ away from the boundary of the ROI. Since the 4FGL-DR3 only contains
point sources resolved with 12 years of Fermi data, new point sources may emerge in the 14
year data we use. These new point sources, if they exist, must be included in the background
model. Otherwise, their photons may be attributed to dark matter models. We use the
find_sources algorithm provided by fermipy to search for new point sources in each ROI
and include them in the background model. We find seven new point sources in the ROIs for
four of the galaxy clusters in our list. Table 2 summarizes these new point source locations,
offsets from the galaxy clusters, and test statistic (TS)1 . There sources are included in the
background models and their normalizations are allowed to vary.
We consider the gamma-ray signal from dark matter as an extended emission in the ROI
up to the virial radius of the dark matter halo. We use the public software CLUMPY [71–73] to
calculate the angular-dependent J-factor dJ dec. /dΩ(ϕ) of the cluster following an NFW profile
and use it as a diffuse template in the analysis. CLUMPY has been widely used to generate
templates of dark matter halos that can be directly used in the analysis of Fermi-LAT data.
In each fit, we fix mχ and the decaying channel of dark matter and the only free parameter
of the dark matter model is the normalization of the gamma-ray flux, which is inversely
proportional to τχ . The best-fitting model is found by maximizing the Poissonian likelihood
function summed over each pixel i and energy bin j for the k-th cluster, which is given by
k k
Y µkij (θ)nij e−µij (θ)
k
L (θ) = . (4.1)
ij
nkij !

Here, nkij represents the observed photon counts, and µkij (θ) represents the predicted
P photon
counts at the pixel i and energy bin j of the k-th cluster. The parameter θ = {τχ , l λl } is a
set of parameters includes the dark matter lifetime τχ and the parameters of the astrophysical
backgrounds λl . For the backgrounds, we vary the normalizations and spectral shapes of the
diffuse emissions and the point sources within 10◦ from the centers of the ROIs.
We obtain the limit on τχ at the 95% confidence interval (CL) using the profile likelihood
method [86, 87]. For each individual cluster, we scan τχ and refit the background parameters
λl to find the change in the likelihood function from the best-fitting model. The 95% CL
limit on the dark matter lifetime is set by finding the set of nuisance parameters for which
 
−2∆ log Lk = −2 log Lk (θ̂) − log Lk (θ) = −2.71, (4.2)

where θ̂ is the best fit parameters we found for cluster.

1
We note that the find_sources algorithm provided by fermipy only offers a straightforward routine to
search for potential new sources after fitting the known background sources. We include this procedure to
optimize our background models for the search of VHDM signals. The true nature of these sources requires
further investigation.

–6–
 Cluster Name l b offset TS
[deg] [deg] [deg]
Centaurus PS J1320.6−4524 308.29 17.15 7.28 68.12
Centaurus PS J1158.2−4419 293.00 17.51 9.62 30.01
Norma PS J1524.6−6320 319.16 −5.40 6.36 33.34
Perseus PS J0254.8+3934 147.40 −17.36 5.12 38.28
Perseus PS J0400.1+3710 159.61 −11.91 8.92 35.55
Coma PS J1323.3+3212 71.07 81.42 6.56 45.37
Coma PS J1316.6+2056 347.17 81.55 7.98 32.06

Table 2. Additional point sources included in the ROIs.

Virgo χ → bb̄

1027 1027

1026 1026
τχ [s]

τχ [s]

1025 1025
bb̄ τ +τ − Virgo Coma
10 24 tt̄ Z 0Z 0 1024 Centaurus Hydra
e+ e− W +W − Norma Fornax
µ+ µ− Perseus
1023 3 1023 3
10 105 107 109 1011 1013 1015 1016 10 105 107 109 1011 1013 1015 1016
mχ [GeV] mχ [GeV]

Figure 5. Limits on τχ for five decaying channels Figure 6. Limits on τχ for the bb̄ channel from
from the Virgo cluster. seven galaxy clusters.

5 Results and discussion

In figure 5, we present the 95% CL lower limits on τχ for five decaying channels from the Virgo
cluster for 103 GeV < mχ < 1016 GeV. Overall, the most stringent limits have τχ ≳ 1027 s
at mχ ∼ TeV. The limits become less stringent when mχ increases since the expected signals
change from prompt dominated region to secondary (inverse-Compton) dominated region,
which gradually becomes less intense. The limits weaken to τχ ≳ 1026 s when mχ is around
EeV. For mχ ≳ EeV, synchrotron emission from prompt and secondary e± becomes important
in the Fermi energy range, and the limits tend to recover from mχ ∼ EeV. For bb̄ and tt̄ and
Z 0 Z 0 and W + W − channels, the limits gradually increase and eventually reach to τχ ≳ 1027
at mχ ∼ 100 EeV. Above mχ ∼ 100 EeV, the limits gradually become weaker as the peak of
the synchrotron emission pass the Fermi energy range. For e+ e− , µ+ µ− and τ + τ − channels,
the behavior of the limits is different since these channels are largely leptonic. The limits
recover early and quickly to τχ ≳ 1027 s at mχ ∼ EeV. The limits also show a second peak at
1014 GeV since the inverse-Compton photons are attenuated and pair-create e± , generating
a secondary synchrotron component.

–7–
 Figure 6 shows the limits for the bb̄ channel from seven galaxy clusters. Generally, larger
J dec. in clusters provide more stringent limits, as expected. Across the entire mass range, the
most stringent limits are set by two clusters: Virgo and Centaurus. The Virgo cluster has a
J dec. nearly one order of magnitude larger than that of Centaurus. However, the limits set
by Centaurus are overall comparable to and often stronger than Virgo. This is largely caused
by the fact that the Messier 87 (M87) galaxy at the center of Virgo cluster is detected as an
active galactic nucleus by Fermi-LAT [88]. Since we refit the background sources when we
derive the limits, a fraction of the gamma-ray flux from M87 is attributed to dark matter
in the Virgo cluster due to their spatial overlapping, leading to weakened limits on τχ . On
the top panel of figure 7, we show TSVirgo , the TS value of the dark matter component in
the Virgo cluster for each mχ and decaying channel. In most cases, TSVirgo are less than 1.
On the bottom panel, we show −∆TSM87 , the reduction of the TS value of M87 when the
corresponding dark matter model is included in the fit. It is obvious that −∆TSM87 correlates
with TSVirgo . This clearly shows how M87 affects the inferred dark matter limits from the
Virgo cluster. It is known that background sources are impactful for dark matter constraints
in general. By refitting background sources in all clusters, we obtain conservative and robust
limits. Figure 8 shows the limits for the remaining decaying channels. Same as the bb̄ channel,
the Virgo and Centaurus clusters together set the most stringent constraints on τχ .
In figure 9, we compare our results with limits set by galaxy clusters from recent gamma-
ray and neutrino observations. For decaying dark matter, MAGIC has reported limits on τχ in
the mass range between 200 GeV and 200 TeV from the Perseus cluster using 400 hours of data
between 2009 and 2017 [52]. HAWC has also performed a search from the Virgo cluster for
dark matter masses between 1 TeV and 100 TeV using 1523 days of observation [53]. Using 6-
year neutrino data, IceCube has constrained decaying dark matter in the mass range between
10 TeV and 10 PeV from observations of three clusters (Virgo, Coma, and Perseus) [63]. In
figure 9, we show the HAWC limits in dashed lines, the MAGIC limits in dashed-dotted lines,
and the IceCube limits, which we have converted from 90% CL to 95% CL according to the χ2
distribution, in dotted lines. We compare these limits with our results from Virgo (blue solid
lines) and Centaurus (orange solid lines) since they set the most stringent limits. We compare
the limits for 103 GeV < mχ < 107 GeV. Our limits are more stringent than MAGIC for
the four available channels (bb̄, µ+ µ− , τ + τ − , and W + W − ). At mχ = 103 GeV, our limits
are stronger than MAGIC’s by about 2 to 3 orders of magnitude since Fermi-LAT can fully
cover the prompt emission from dark matter. HAWC limits are not available for the µ+ µ−
channel. For other channels, HAWC limits are generally stronger than MAGIC by one order
of magnitude. Our limits are usually stronger than HAWC for mχ ≲ 104 GeV. We note that
limits from MAGIC and HAWC do not consider the electromagnetic cascade from VHDM
decay. IceCube limits are available for the bb̄ and τ + τ − channels. Our limits are competitive
with IceCube for mχ ≲ 105 GeV. However, IceCube can access much higher energies, and
therefore has set the strongest limits for mχ > 105 GeV.

5.1 Combined analysis

So far, we have reported limits on τχ from individual clusters. It is also possible to combine
the analyses and define the total likelihood function of all seven clusters:
Y
L(θ) = Lk (θ). (5.1)
k

–8–
 103
bb̄ τ +τ −
tt̄ Z 0Z 0
102
e+ e− W +W −
µ+ µ−
101
TSVirgo

100

10−1

10−2

10−3 3
10 104 105 106 107 108 109 1010 1011 1012 1013 1014 1015 1016
mχ [GeV]

103

102

101
−∆TSM87

100

10−1 bb̄ τ +τ −
tt̄ Z 0Z 0
10−2 e+ e− W +W −
µ+ µ−
10−3 3
10 104 105 106 107 108 109 1010 1011 1012 1013 1014 1015 1016
mχ [GeV]

Figure 7. Top. TS value of the dark matter component in the Virgo cluster for each mχ and decaying
channel. Bottom. Reduction of the TS value of M87 for the corresponding mχ and decaying channel.

We derive the combined limits at 95% CL at

 
−2∆ log L = −2 log L(θ̂) − log L(θ) = −2.71. (5.2)

Naturally, τχ is assumed to be the same in all clusters and the likelihood functions are
combined. The background parameters λl now include all the point sources and diffuse
emissions in the ROIs of the clusters. Although we use the same models for the diffuse
emissions, their parameters are allowed to have different values for different clusters.
Figure 10 shows the combined limits from all seven clusters. Overall, they are slightly
weaker than the strongest limits set by Virgo and Centaurus individually. This is driven by
clusters with small J dec. . Figure 11 shows the contributions to −2∆ log(L) from the seven
clusters when the 95% CL limits on τχ are set for the bb̄ channel and for different mχ . Fornax,
Hydra, Coma, and partially Norma contribute negatively to the −2∆ log(L), meaning that
including them in the combined analyses weakens the limits. This is due to their small
expected signals and potentially more complicated background environments (e.g., Norma is
located close to the Galactic plane). Previous study also finds similar results in combined

–9–
 χ → tt̄ χ → e+ e− χ → µ+ µ−

1027 1027 1027

1026 1026 1026

τχ [s]

τχ [s]

τχ [s]
1025 1025 1025
Virgo Coma Virgo Coma Virgo Coma
1024 Centaurus Hydra 1024 Centaurus Hydra 1024 Centaurus Hydra
Norma Fornax Norma Fornax Norma Fornax
Perseus Perseus Perseus
1023 3 1023 3 1023 3
10 105 107 109 1011 1013 1015 1016 10 105 107 109 1011 1013 1015 1016 10 105 107 109 1011 1013 1015 1016
mχ [GeV] mχ [GeV] mχ [GeV]

χ → τ +τ − χ → Z 0Z 0 χ → W +W −

1027 1027 1027

1026 1026 1026

τχ [s]

τχ [s]

τχ [s]
1025 1025 1025
Virgo Coma Virgo Coma Virgo Coma
1024 Centaurus Hydra 1024 Centaurus Hydra 1024 Centaurus Hydra
Norma Fornax Norma Fornax Norma Fornax
Perseus Perseus Perseus
1023 3 1023 3 1023 3
10 105 107 109 1011 1013 1015 1016 10 105 107 109 1011 1013 1015 1016 10 105 107 109 1011 1013 1015 1016
mχ [GeV] mχ [GeV] mχ [GeV]

Figure 8. Limits on τχ for the other decaying channels: tt̄, e+ e− , µ+ µ− , τ + τ − , Z 0 Z 0 , and W + W − .

likelihood analysis [34]. Therefore, we find that promising individual clusters provide the
strongest probe of VHDM.

5.2 Systematic uncertainties

To derive limits on τχ , we take into account the electromagnetic cascades of prompt gamma
rays and e± from dark matter decay. Therefore, we have to consider systematic uncertain-
ties raised from the calculation of the cascade flux. Two factors are mostly influential: the
magnetic field Bcl in the cluster in which e± lose energies and the escaping distance Rcl
at which gamma rays leave the cluster and enter the intergalactic space. We have adopted
Bcl = 0.3 µG and Rcl = 3 Mpc. Typical intracluster magnetic fields in galaxy clusters av-
erage to ∼ 0.1 − 1 µG [80]. For the seven clusters in this work, their virial radii are ∼ 1
– 3 Mpc (see table 1). In figure 12, we show the limits from Virgo for the bb̄ channel by
using Bcl = 0.1, 0.3, and 1.0 µG. The limits are largely the same at lower masses when the
prompt emission is more important. The limits slightly change shapes for mχ ≳ 105 GeV
but maintain roughly the same level. Realistically, the magnetic fields within the volume of
the galaxy clusters vary, and the spectral features of the VHDM signals would alter corre-
spondingly. Nonetheless, the constraints are not expected to dramatically change, as they
are primarily determined by the predicted fluxes in the Fermi-LAT energy range and electro-
magnetic cascades smear out the spectra. Therefore, the impacts of varying magnetic fields
within the volume of the clusters are expected to lie around the range of the constraints shown
in figure 12. Figure 13 shows the limits from Virgo for the bb̄ channel, using Rcl = 1 and 3
Mpc. For mχ ≳ 1010 GeV, the limits from Rcl = 3 Mpc are slightly stronger since stronger
synchrotron emission is expected from larger Rcl . However, the effect is not significant. We

– 10 –
 χ → bb̄ χ → µ+ µ−

1028 1028

1027 1027
τχ [s]

τχ [s]
1026 1026
Virgo (This Work)
1025 Centaurus (This Work)
1025
HAWC (2021) Virgo (This Work)
1024 MAGIC (2018) 1024 Centaurus (This Work)
IceCube (2021) MAGIC (2018)
1023 3 1023 3
10 104 105 106 107 10 104 105 106 107
mχ [GeV] mχ [GeV]

χ → τ +τ − χ → W +W −

1028 1028

1027 1027
τχ [s]

τχ [s]

1026 1026
Virgo (This Work)
1025 Centaurus (This Work)
1025 Virgo (This Work)
HAWC (2021) Centaurus (This Work)
1024 MAGIC (2018) 1024 HAWC (2021)
IceCube (2021) MAGIC (2018)
23 23
10 10
103 104 105 106 107 103 104 105 106 107
mχ [GeV] mχ [GeV]

Figure 9. Comparisons between our limits and recent gamma-ray (MAGIC [52] and HAWC [53])
and neutrino (IceCube [63]) observations on galaxy clusters.

expect similar uncertainties for different clusters and decay channels. Overall, our limits are
robust when Bcl and Rcl are in reasonable ranges.

6 Summary

We searched for gamma-ray signals from very heavy dark matter decay from nearby galaxy
clusters using 14 year Fermi-LAT data and set lower limits on the lifetime of decaying dark
matter as a function of dark matter mass. Notably, we considered gamma rays from both the
prompt emission and the electromagnetic cascades in the clusters. The cascade component is
essential for constraints at large masses. Seven decaying channels are included: bb̄, tt̄, e+ e− ,
µ+ µ− , τ + τ − , Z 0 Z 0 , and W + W − . We use the profile likelihood method to obtain the limits

– 11 –
 Combined χ → bb̄
Virgo
1027 10.0 Centaurus
Norma
Perseus
7.5
Coma
1026

−2∆ log(L)
Hydra
5.0
τχ [s]

bb̄ Fornax
tt̄
1025 2.5
e+ e−
µ+ µ−
1024 τ +τ − 0.0
Z 0Z 0
W +W − −2.5
23
10
103 105 107 109 1011 1013 1015 1016 103 105 107 109 1011 1013 1015 1016
mχ [GeV] mχ [GeV]

Figure 10. Limits on τχ for five decaying chan- Figure 11. Contributions to −2∆ log(L) at the
nels from combined analyese of seven clusters. 95% CL limits for the bb̄ channel from seven clus-
ters in the combined analyses.

Virgo, χ → bb̄ Virgo, χ → bb̄

Bcl = 0.1 µG Rcl = 1 Mpc
Bcl = 0.3 µG Rcl = 3 Mpc
Bcl = 1.0 µG

1027 1027
τχ [s]

τχ [s]

1026 3 1026 3
10 105 107 109 1011 1013 1015 1016 10 105 107 109 1011 1013 1015 1016
mχ [GeV] mχ [GeV]

Figure 12. Limits on τχ from Virgo for the bb̄ Figure 13. Limits on τχ from Virgo for the bb̄
channel assuming different Bcl . channel assuming different Rcl .

for seven galaxy clusters: Virgo, Centaurus, Norma, Perseus, Coma, Hydra, and Fornax. We
found that gamma-ray observations in the direction of Virgo and Centaurus clusters provide
the most stringent limits on dark matter lifetime. We also reported the limits from combined
analyses of seven clusters. Our limits remain strong and robust when we consider different
magnetic field strengths and sizes of the clusters. Our results are competitive with limits
from other gamma-ray observations on galaxy clusters and are complementary to previous
gamma-ray and neutrino constraints on decaying very heavy dark matter.

– 12 –
Acknowledgments

The authors thank the anonymous referee for their valuable comments and suggestions. We
thank Saikat Das, Nagisa Hiroshima, and B. Thedore Zhang for fruitful discussions. D.S. is
supported by KAKENHI No. 20H05852. The work of K.M. was supported by the NSF Grants
No. AST-2108466 and No. AST-2108467, and KAKENHI No. 20H01901 and No. 20H05852.
A.K. is supported by the NASA grant 80NSSC23M0104. This work uses the computational
resources provided by the super computer Yukawa-21 at Yukawa Institute for Theoretical
Physics.

A Electromagnetic cascades inside galaxy clusters

Following Ref. [36] (see also Refs. [79, 89, 90] in the context of cosmic rays) , our numerical
code solves the following Boltzmann equations in a time-dependent manner,
∂Nγ (Eγ ) dn dµ Nγ
Z Z
= − Nγ dε (1 − µ)cσγγ (ε, µ) − (A.1)
∂t dε 2 tesc
dn dµ dσIC
Z Z Z
+ dE ′ Ne E ′ ε, µ, E ′



dε (1 − µ)c
dε 2 dEγ
syn
∂Nγ
+ + Qinj
γ ,
∂t
∂Ne (Ee ) dn dµ
Z Z
= − Ne dε (1 − µ)cσIC (ε, µ) (A.2)
∂t dε 2
dn dµ dσγγ
Z Z Z
+ dE ′ Nγ E ′ ε, µ, E ′



dε (1 − µ)c
dε 2 dEe
dn dµ dσIC
Z Z Z
+ dE ′ Ne E ′ ε, µ, E ′



dε (1 − µ)c
dε 2 dEe
∂
− [Psyn Ne ] + Qinj
e .
∂E
The background photon number density is dn/dε at the energy ε. For infrared and optical
photons in galaxy clusters, we adopt the low-IR model of the EBL light from Ref. [91] with
a 10 times enhancement to account for contributions from the galaxies in the clusters [36,
92, 93]. This is reasonable within uncertainties (see Refs. [90, 93, 94]). The CMB is also
included. We have implemented the accurate cross sections for pair production (σγγ ) and
inverse-Compton scattering (σIC ) including the Klein-Nishina effect, as in calculations of
intergalactic propagation of gamma rays and cosmic rays [77, 78, 95, 96]. Here, Psyn is the
energy loss rate due to the synchrotron radiation. The differential synchrotron radiation
spectrum ∂Nγsyn /∂t is
√ 3
∂Nγsyn (Eγ ) 3e B
Z
′ ′


= dE Ne E 2
G(x), (A.3)
∂t me c 2πℏEγ
1.81e−x
G(x) ≈
1/2 , (A.4)
x−2/3 + (3.62/π)2

where x = Eγ /Ec and Ec is the critical energy. The function form G(x) is an accurate fit
to the synchrotron radiation spectrum [97, 98]. The diffusion of electrons and positrons is
negligible for high-energy radiation [78]. For example, at TeV energies, the e± cooling time

– 13 –
   
TeV uCMB
is tcool ≈ tIC ∼ 1 Myr [78]. This is much shorter than the light crossing
Ee  urad
Rcl
time, tesc = Rcl /c ∼ 10 Myr . The diffusion time scale tdiff is always longer than
3 Mpc
the light crossing time tesc , implying tcool ≪ tdiff , which is also the case in the calculations of
cascade gamma rays induced by cosmic rays confined in galaxy clusters [36, 79].
The Boltzmann equations are numerically binned in finite energy and time steps and are
solved iteratively to obtain essentially steady-state spectra for high-energy gamma rays [36].

References

[1] Planck collaboration, Planck 2018 results. VI. Cosmological parameters, Astron. Astrophys.
641 (2020) A6 [1807.06209].
[2] G. Bertone, D. Hooper and J. Silk, Particle dark matter: Evidence, candidates and constraints,
Phys. Rept. 405 (2005) 279 [hep-ph/0404175].
[3] K. Griest and M. Kamionkowski, Unitarity Limits on the Mass and Radius of Dark Matter
Particles, Phys. Rev. Lett. 64 (1990) 615.
[4] J. Bramante and J. Unwin, Superheavy Thermal Dark Matter and Primordial Asymmetries,
JHEP 02 (2017) 119 [1701.05859].
[5] M. Cirelli, Y. Gouttenoire, K. Petraki and F. Sala, Homeopathic Dark Matter, or how diluted
heavy substances produce high energy cosmic rays, JCAP 02 (2019) 014 [1811.03608].
[6] D. Bhatia and S. Mukhopadhyay, Unitarity limits on thermal dark matter in (non-)standard
cosmologies, JHEP 03 (2021) 133 [2010.09762].
[7] V.K. Dubrovich, D. Fargion and M.Y. Khlopov, Primordial bound systems of superheavy
particles as the source of ultrahigh-energy cosmic rays, Astropart. Phys. 22 (2004) 183
[hep-ph/0312105].
[8] J. Hisano, S. Matsumoto, M. Nagai, O. Saito and M. Senami, Non-perturbative effect on
thermal relic abundance of dark matter, Phys. Lett. B 646 (2007) 34 [hep-ph/0610249].
[9] B. von Harling and K. Petraki, Bound-state formation for thermal relic dark matter and
unitarity, JCAP 12 (2014) 033 [1407.7874].
[10] I. Baldes and K. Petraki, Asymmetric thermal-relic dark matter: Sommerfeld-enhanced
freeze-out, annihilation signals and unitarity bounds, JCAP 09 (2017) 028 [1703.00478].
[11] J. Smirnov and J.F. Beacom, TeV-Scale Thermal WIMPs: Unitarity and its Consequences,
Phys. Rev. D 100 (2019) 043029 [1904.11503].
[12] I. Baldes, Y. Gouttenoire, F. Sala and G. Servant, Supercool composite Dark Matter beyond 100
TeV, JHEP 07 (2022) 084 [2110.13926].
[13] L. Goodenough and D. Hooper, Possible Evidence For Dark Matter Annihilation In The Inner
Milky Way From The Fermi Gamma Ray Space Telescope, 0910.2998.
[14] D. Hooper and T. Linden, On The Origin Of The Gamma Rays From The Galactic Center,
Phys. Rev. D 84 (2011) 123005 [1110.0006].
[15] D. Hooper and L. Goodenough, Dark Matter Annihilation in The Galactic Center As Seen by
the Fermi Gamma Ray Space Telescope, Phys. Lett. B 697 (2011) 412 [1010.2752].
[16] K.N. Abazajian and M. Kaplinghat, Detection of a Gamma-Ray Source in the Galactic Center
Consistent with Extended Emission from Dark Matter Annihilation and Concentrated
Astrophysical Emission, Phys. Rev. D 86 (2012) 083511 [1207.6047].

– 14 –
[17] C. Gordon and O. Macias, Dark Matter and Pulsar Model Constraints from Galactic Center
Fermi-LAT Gamma Ray Observations, Phys. Rev. D 88 (2013) 083521 [1306.5725].
[18] K.N. Abazajian, N. Canac, S. Horiuchi and M. Kaplinghat, Astrophysical and Dark Matter
Interpretations of Extended Gamma-Ray Emission from the Galactic Center, Phys. Rev. D 90
(2014) 023526 [1402.4090].
[19] F. Calore, I. Cholis and C. Weniger, Background Model Systematics for the Fermi GeV Excess,
JCAP 03 (2015) 038 [1409.0042].
[20] T. Daylan, D.P. Finkbeiner, D. Hooper, T. Linden, S.K.N. Portillo, N.L. Rodd et al., The
characterization of the gamma-ray signal from the central Milky Way: A case for annihilating
dark matter, Phys. Dark Univ. 12 (2016) 1 [1402.6703].
[21] Fermi-LAT collaboration, Fermi-LAT Observations of High-Energy γ-Ray Emission Toward
the Galactic Center, Astrophys. J. 819 (2016) 44 [1511.02938].
[22] D. Hooper, C. Kelso and F.S. Queiroz, Stringent and Robust Constraints on the Dark Matter
Annihilation Cross Section From the Region of the Galactic Center, Astropart. Phys. 46 (2013)
55 [1209.3015].
[23] I. Cholis and P. Salucci, Extracting limits on Dark Matter annihilation from gamma-ray
observations towards dwarf spheroidal galaxies, Phys. Rev. D 86 (2012) 023528 [1203.2954].
[24] E. Carlson, D. Hooper and T. Linden, Improving the Sensitivity of Gamma-Ray Telescopes to
Dark Matter Annihilation in Dwarf Spheroidal Galaxies, Phys. Rev. D 91 (2015) 061302
[1409.1572].
[25] MAGIC, Fermi-LAT collaboration, Limits to Dark Matter Annihilation Cross-Section from a
Combined Analysis of MAGIC and Fermi-LAT Observations of Dwarf Satellite Galaxies, JCAP
02 (2016) 039 [1601.06590].
[26] A. Lopez, C. Savage, D. Spolyar and D.Q. Adams, Fermi/LAT observations of Dwarf Galaxies
highly constrain a Dark Matter Interpretation of Excess Positrons seen in AMS-02, HEAT, and
PAMELA, JCAP 03 (2016) 033 [1501.01618].
[27] Fermi-LAT collaboration, Searching for Dark Matter Annihilation from Milky Way Dwarf
Spheroidal Galaxies with Six Years of Fermi Large Area Telescope Data, Phys. Rev. Lett. 115
(2015) 231301 [1503.02641].
[28] M.G. Baring, T. Ghosh, F.S. Queiroz and K. Sinha, New Limits on the Dark Matter Lifetime
from Dwarf Spheroidal Galaxies using Fermi-LAT, Phys. Rev. D 93 (2016) 103009
[1510.00389].
[29] F. Calore, P.D. Serpico and B. Zaldivar, Dark matter constraints from dwarf galaxies: a
data-driven analysis, JCAP 10 (2018) 029 [1803.05508].
[30] S. Hoof, A. Geringer-Sameth and R. Trotta, A Global Analysis of Dark Matter Signals from 27
Dwarf Spheroidal Galaxies using 11 Years of Fermi-LAT Observations, JCAP 02 (2020) 012
[1812.06986].
[31] A. Alvarez, F. Calore, A. Genina, J. Read, P.D. Serpico and B. Zaldivar, Dark matter
constraints from dwarf galaxies with data-driven J-factors, JCAP 09 (2020) 004 [2002.01229].
[32] Fermi-LAT collaboration, Dark Matter search in dwarf irregular galaxies with the Fermi Large
Area Telescope, PoS ICRC2021 (2021) 509 [2109.11291].
[33] M. Ackermann, M. Ajello, A. Allafort, L. Baldini, J. Ballet, G. Barbiellini et al., Constraints on
dark matter annihilation in clusters of galaxies with the Fermi large area telescope, J.
Cosmology Astropart. Phys. 2010 (2010) 025 [1002.2239].
[34] X. Huang, G. Vertongen and C. Weniger, Probing Dark Matter Decay and Annihilation with
Fermi LAT Observations of Nearby Galaxy Clusters, JCAP 01 (2012) 042 [1110.1529].

– 15 –
[35] S. Ando and D. Nagai, Fermi-LAT constraints on dark matter annihilation cross section from
observations of the Fornax cluster, J. Cosmology Astropart. Phys. 2012 (2012) 017
[1201.0753].
[36] K. Murase and J.F. Beacom, Galaxy Clusters as Reservoirs of Heavy Dark Matter and
High-Energy Cosmic Rays: Constraints from Neutrino Observations, JCAP 02 (2013) 028
[1209.0225].
[37] Fermi-LAT collaboration, Search for extended gamma-ray emission from the Virgo galaxy
cluster with Fermi-LAT, Astrophys. J. 812 (2015) 159 [1510.00004].
[38] X. Tan, M. Colavincenzo and S. Ammazzalorso, Bounds on WIMP dark matter from galaxy
clusters at low redshift, Mon. Not. Roy. Astron. Soc. 495 (2020) 114 [1907.06905].
[39] C. Thorpe-Morgan, D. Malyshev, C.-A. Stegen, A. Santangelo and J. Jochum, Annihilating
dark matter search with 12 yr of Fermi LAT data in nearby galaxy clusters, Mon. Not. Roy.
Astron. Soc. 502 (2021) 4039 [2010.11006].
[40] M. Di Mauro, J. Pérez-Romero, M.A. Sánchez-Conde and N. Fornengo, Constraining the dark
matter contribution of γ rays in clusters of galaxies using Fermi-LAT data, Phys. Rev. D 107
(2023) 083030 [2303.16930].
[41] M. Cirelli, P. Panci and P.D. Serpico, Diffuse gamma ray constraints on annihilating or
decaying Dark Matter after Fermi, Nucl. Phys. B 840 (2010) 284 [0912.0663].
[42] K.N. Abazajian, S. Blanchet and J.P. Harding, Current and Future Constraints on Dark
Matter from Prompt and Inverse-Compton Photon Emission in the Isotropic Diffuse
Gamma-ray Background, Phys. Rev. D 85 (2012) 043509 [1011.5090].
[43] K. Murase and J.F. Beacom, Constraining Very Heavy Dark Matter Using Diffuse Backgrounds
of Neutrinos and Cascaded Gamma Rays, JCAP 10 (2012) 043 [1206.2595].
[44] T. Bringmann, F. Calore, M. Di Mauro and F. Donato, Constraining dark matter annihilation
with the isotropic γ-ray background: updated limits and future potential, Phys. Rev. D 89
(2014) 023012 [1303.3284].
[45] M. Ajello et al., The Origin of the Extragalactic Gamma-Ray Background and Implications for
Dark-Matter Annihilation, Astrophys. J. Lett. 800 (2015) L27 [1501.05301].
[46] M. Di Mauro and F. Donato, Composition of the Fermi-LAT isotropic gamma-ray background
intensity: Emission from extragalactic point sources and dark matter annihilations, Phys. Rev.
D 91 (2015) 123001 [1501.05316].
[47] Fermi-LAT collaboration, Limits on Dark Matter Annihilation Signals from the Fermi LAT
4-year Measurement of the Isotropic Gamma-Ray Background, JCAP 09 (2015) 008
[1501.05464].
[48] J. Aleksić, L.A. Antonelli, P. Antoranz, M. Backes, C. Baixeras, S. Balestra et al., MAGIC
Gamma-ray Telescope Observation of the Perseus Cluster of Galaxies: Implications for Cosmic
Rays, Dark Matter, and NGC 1275, ApJ 710 (2010) 634 [0909.3267].
[49] T. Arlen, T. Aune, M. Beilicke, W. Benbow, A. Bouvier, J.H. Buckley et al., Constraints on
Cosmic Rays, Magnetic Fields, and Dark Matter from Gamma-Ray Observations of the Coma
Cluster of Galaxies with VERITAS and Fermi, ApJ 757 (2012) 123 [1208.0676].
[50] MAGIC collaboration, Constraining the Dark Matter decay lifetime with very deep observations
of the Perseus cluster with the MAGIC telescopes, PoS ICRC2015 (2016) 1204 [1509.03974].
[51] HAWC collaboration, Dark Matter Annihilation and Decay Searches with the High Altitude
Water Cherenkov (HAWC) Observatory, PoS ICRC2015 (2016) 1227 [1508.04352].
[52] MAGIC collaboration, Constraining Dark Matter lifetime with a deep gamma-ray survey of the
Perseus Galaxy Cluster with MAGIC, Phys. Dark Univ. 22 (2018) 38 [1806.11063].

– 16 –
[53] HAWC collaboration, Search for TeV decaying dark matter from the Virgo cluster of galaxies,
PoS ICRC2021 (2021) 530.
[54] T.N. Maity, A.K. Saha, A. Dubey and R. Laha, Search for dark matter using sub-PeV γ-rays
observed by Tibet ASγ, 2105.05680.
[55] D. Tak, M. Baumgart, N.L. Rodd and E. Pueschel, Current and Future γ-Ray Searches for
Dark Matter Annihilation Beyond the Unitarity Limit, Astrophys. J. Lett. 938 (2022) L4
[2208.11740].
[56] K. Ishiwata, O. Macias, S. Ando and M. Arimoto, Probing heavy dark matter decays with
multi-messenger astrophysical data, JCAP 01 (2020) 003 [1907.11671].
[57] M. Chianese, D.F.G. Fiorillo, R. Hajjar, G. Miele and N. Saviano, Constraints on heavy
decaying dark matter with current gamma-ray measurements, JCAP 11 (2021) 035
[2108.01678].
[58] M. Chianese, D.F.G. Fiorillo, R. Hajjar, G. Miele, S. Morisi and N. Saviano, Heavy decaying
dark matter at future neutrino radio telescopes, JCAP 05 (2021) 074 [2103.03254].
[59] S. Das, K. Murase and T. Fujii, Revisiting ultrahigh-energy constraints on decaying superheavy
dark matter, Phys. Rev. D 107 (2023) 103013 [2302.02993].
[60] D.F.G. Fiorillo, V.B. Valera, M. Bustamante and W. Winter, Searches for dark matter decay
with ultra-high-energy neutrinos endure backgrounds, 2307.02538.
[61] IceCube collaboration, Search for neutrinos from decaying dark matter with IceCube, Eur.
Phys. J. C 78 (2018) 831 [1804.03848].
[62] IceCube collaboration, Searches for Connections Between Dark Matter and Neutrinos with the
IceCube High-Energy Starting Event Sample, PoS ICRC2019 (2020) 839 [1907.11193].
[63] IceCube collaboration, A Search for Neutrinos from Decaying Dark Matter in Galaxy Clusters
and Galaxies with IceCube, PoS ICRC2021 (2021) 506 [2107.11527].
[64] LHAASO collaboration, Constraints on Heavy Decaying Dark Matter from 570 Days of
LHAASO Observations, Phys. Rev. Lett. 129 (2022) 261103 [2210.15989].
[65] C.A. Argüelles, D. Delgado, A. Friedlander, A. Kheirandish, I. Safa, A.C. Vincent et al., Dark
Matter decay to neutrinos, 2210.01303.
[66] M. Lisanti, S. Mishra-Sharma, N.L. Rodd and B.R. Safdi, Search for Dark Matter Annihilation
in Galaxy Groups, Phys. Rev. Lett. 120 (2018) 101101 [1708.09385].
[67] K. Murase, R. Laha, S. Ando and M. Ahlers, Testing the Dark Matter Scenario for PeV
Neutrinos Observed in IceCube, Phys. Rev. Lett. 115 (2015) 071301 [1503.04663].
[68] T. Cohen, K. Murase, N.L. Rodd, B.R. Safdi and Y. Soreq, γ -ray Constraints on Decaying
Dark Matter and Implications for IceCube, Phys. Rev. Lett. 119 (2017) 021102 [1612.05638].
[69] C. Blanco and D. Hooper, Constraints on Decaying Dark Matter from the Isotropic
Gamma-Ray Background, JCAP 03 (2019) 019 [1811.05988].
[70] J.F. Navarro, C.S. Frenk and S.D.M. White, A Universal Density Profile from Hierarchical
Clustering, ApJ 490 (1997) 493 [astro-ph/9611107].
[71] A. Charbonnier, C. Combet and D. Maurin, CLUMPY: A code for γ-ray signals from dark
matter structures, Computer Physics Communications 183 (2012) 656 [1201.4728].
[72] V. Bonnivard, M. Hütten, E. Nezri, A. Charbonnier, C. Combet and D. Maurin, CLUMPY:
Jeans analysis, γ-ray and ν fluxes from dark matter (sub-)structures, Computer Physics
Communications 200 (2016) 336 [1506.07628].
[73] M. Hütten, C. Combet and D. Maurin, CLUMPY v3: γ-ray and ν signals from dark matter at
all scales, Computer Physics Communications 235 (2019) 336 [1806.08639].

– 17 –
[74] C.W. Bauer, N.L. Rodd and B.R. Webber, Dark matter spectra from the electroweak to the
Planck scale, Journal of High Energy Physics 2021 (2021) 121 [2007.15001].
[75] F.A. Aharonian, P.S. Coppi and H.J. Voelk, Very high-energy gamma-rays from AGN:
Cascading on the cosmic background radiation fields and the formation of pair halos, Astrophys.
J. Lett. 423 (1994) L5 [astro-ph/9312045].
[76] K. Murase, K. Takahashi, S. Inoue, K. Ichiki and S. Nagataki, Probing Intergalactic Magnetic
Fields in the GLAST Era through Pair Echo Emission from TeV Blazars, Astrophys. J. Lett.
686 (2008) L67 [0806.2829].
[77] K. Murase, J.F. Beacom and H. Takami, Gamma-Ray and Neutrino Backgrounds as Probes of
the High-Energy Universe: Hints of Cascades, General Constraints, and Implications for TeV
Searches, JCAP 08 (2012) 030 [1205.5755].
[78] K. Murase, High-Energy Emission Induced by Ultra-high-Energy Photons as a Probe of
Ultra-high-Energy Cosmic-Ray Accelerators Embedded in the Cosmic Web, Astrophys. J. Lett.
745 (2012) L16 [1111.0936].
[79] K. Murase, S. Inoue and K. Asano, Cosmic rays above the 2nd knee from clusters of galaxies,
Int. J. Mod. Phys. D 18 (2009) 1609.
[80] F. Govoni and L. Feretti, Magnetic Fields in Clusters of Galaxies, International Journal of
Modern Physics D 13 (2004) 1549 [astro-ph/0410182].
[81] R.B. Tully, Galaxy Groups: A 2MASS Catalog, AJ 149 (2015) 171 [1503.03134].
[82] E. Kourkchi and R.B. Tully, Galaxy Groups Within 3500 km s−1 , ApJ 843 (2017) 16
[1705.08068].
[83] M. Wood, R. Caputo, E. Charles, M. Di Mauro, J. Magill, J.S. Perkins et al., Fermipy: An
open-source Python package for analysis of Fermi-LAT Data, in 35th International Cosmic Ray
Conference (ICRC2017), vol. 301 of International Cosmic Ray Conference, p. 824, July, 2017,
DOI [1707.09551].
[84] S. Abdollahi, F. Acero, M. Ackermann, M. Ajello, W.B. Atwood, M. Axelsson et al., Fermi
Large Area Telescope Fourth Source Catalog, ApJS 247 (2020) 33 [1902.10045].
[85] S. Abdollahi, F. Acero, L. Baldini, J. Ballet, D. Bastieri, R. Bellazzini et al., Incremental Fermi
Large Area Telescope Fourth Source Catalog, ApJS 260 (2022) 53 [2201.11184].
[86] W.A. Rolke, A.M. Lopez and J. Conrad, Limits and confidence intervals in the presence of
nuisance parameters, Nucl. Instrum. Meth. A 551 (2005) 493 [physics/0403059].
[87] M. Ackermann, M. Ajello, W.B. Atwood, L. Baldini, G. Barbiellini, D. Bastieri et al.,
Constraints on the Galactic Halo Dark Matter from Fermi-LAT Diffuse Measurements, ApJ
761 (2012) 91 [1205.6474].
[88] A.A. Abdo, M. Ackermann, M. Ajello, W.B. Atwood, M. Axelsson, L. Baldini et al., Fermi
Large Area Telescope Gamma-Ray Detection of the Radio Galaxy M87, ApJ 707 (2009) 55
[0910.3565].
[89] K. Murase, S. Inoue and S. Nagataki, Cosmic Rays Above the Second Knee from Clusters of
Galaxies and Associated High-Energy Neutrino Emission, Astrophys. J. Lett. 689 (2008) L105
[0805.0104].
[90] K. Kotera, D. Allard, K. Murase, J. Aoi, Y. Dubois, T. Pierog et al., Propagation of ultrahigh
energy nuclei in clusters of galaxies: resulting composition and secondary emissions, Astrophys.
J. 707 (2009) 370 [0907.2433].
[91] T.M. Kneiske, T. Bretz, K. Mannheim and D.H. Hartmann, Implications of cosmological
gamma-ray absorption. 2. Modification of gamma-ray spectra, Astron. Astrophys. 413 (2004)
807 [astro-ph/0309141].

– 18 –
[92] L. Silva, G.L. Granato, A. Bressan and L. Danese, Modeling the Effects of Dust on Galactic
Spectral Energy Distributions from the Ultraviolet to the Millimeter Band, Astrophys. J. 509
(1998) 103.
[93] H. Takami and K. Murase, The Role of Structured Magnetic Fields on Constraining Properties
of Transient Sources of Ultra-high-energy Cosmic Rays, ApJ 748 (2012) 9 [1110.3245].
[94] K. Fang and K. Murase, Linking High-Energy Cosmic Particles by Black Hole Jets Embedded
in Large-Scale Structures, Nature Phys. 14 (2018) 396 [1704.00015].
[95] K. Murase, C.D. Dermer, H. Takami and G. Migliori, Blazars as Ultra-High-Energy
Cosmic-Ray Sources: Implications for TeV Gamma-Ray Observations, Astrophys. J. 749
(2012) 63 [1107.5576].
[96] K. Murase, M. Ahlers and B.C. Lacki, Testing the Hadronuclear Origin of PeV Neutrinos
Observed with IceCube, Phys. Rev. D 88 (2013) 121301 [1306.3417].
[97] B.T. Zhang, K. Murase, P.e. Veres and P.e. Mészáros, External Inverse-Compton Emission
from Low-luminosity Gamma-Ray Bursts: Application to GRB 190829A, Astrophys. J. 920
(2021) 55 [2012.07796].
[98] F.A. Aharonian, S.R. Kelner and A.Y. Prosekin, Angular, spectral, and time distributions of
highest energy protons and associated secondary gamma rays and neutrinos propagating through
extragalactic magnetic and radiation fields, Phys. Rev. D 82 (2010) 043002 [1006.1045].

– 19 –