# Gamma-ray Diagnostics of r-process Nucleosynthesis in the Remnants of
Galactic Binary Neutron-Star Mergers

Yukikatsu Terada Graduate School of Science and Engineering, Saitama
University, 255 Sakura-ku, Saitama-shi, Saitama 338-8570, Japan  Institute of
Space and Astronautical Science, Japan Aerospace Exploration Agency,
Sagamihara, Kanagawa, Japan  Yuya Miwa  Graduate School of Science and
Engineering, Saitama University, 255 Sakura-ku, Saitama-shi, Saitama 338-8570,
Japan  Hayato Ohsumi  Graduate School of Science and Engineering, Saitama
University, 255 Sakura-ku, Saitama-shi, Saitama 338-8570, Japan  Shin-ichiro
Fujimoto  National Institute of Technology Kumamoto College, Kumamoto
861-1102, Japan  Satoru Katsuda  Graduate School of Science and Engineering,
Saitama University, 255 Sakura-ku, Saitama-shi, Saitama 338-8570, Japan  Aya
Bamba  Department of Physics, The University of Tokyo, 7-3-1 Hongo, Bunkyo,
Tokyo 113-0033, Japan  Research Center for the Early Universe, School of
Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo, Tokyo 113-0033, Japan
Ryo Yamazaki  Department of Physical Sciences, Aoyama Gakuin University,
5-10-1 Fuchinobe, Sagamihara, Kanagawa 252-5258, Japan  Institute of Laser
Engineering, Osaka University, 2-6, Yamadaoka,Suita, Osaka 565-0871, Japan

(Received Mar 03, 2022; Revised Apr 26, 2022; Accepted May 10, 2022)

###### Abstract

We perform a full nuclear-network numerical calculation of the $r$-process
nuclei in binary neutron-star mergers (NSMs), with the aim of estimating
$\gamma$-ray emissions from the remnants of Galactic NSMs up to $10^{6}$ years
old. The nucleosynthesis calculation of 4,070 nuclei is adopted to provide the
elemental composition ratios of nuclei with an electron fraction $Y_{\rm e}$
between 0.10 and 0.45 . The decay processes of 3,237 unstable nuclei are
simulated to extract the $\gamma$-ray spectra. As a result, the NSMs have
different spectral color in $\gamma$-ray band from various other astronomical
objects at less than $10^{5}$ years old. In addition, we propose a new line-
diagnostic method for $Y_{\rm e}$ that uses the line ratios of either
137mBa/85K or 243Am/60mCo, which become larger than unity for young and old
$r$-process sites, respectively, with a low $Y_{\rm e}$ environment. From an
estimation of the distance limit for $\gamma$-ray observations as a function
of the age, the high sensitivity in the sub-MeV band, at approximately
$10^{-9}$ photons s-1 cm-2 or $10^{-15}$ erg s-1 cm-2, is required to cover
all the NSM remnants in our Galaxy if we assume that the population of NSMs by
Wu et al. (2019). A $\gamma$-ray survey with sensitivities of
$10^{-8}$–$10^{-7}$ photons s-1 cm-2 or $10^{-14}$–$10^{-13}$ erg s-1 cm-2 in
the 70–4000 keV band is expected to find emissions from at least one NSM
remnant under the assumption of NSM rate of 30 Myr-1. The feasibility of
$\gamma$-ray missions to observe Galactic NSMs are also studied.

nuclear reactions, nucleosynthesis, abundances – gamma rays: general –
gravitational waves

††journal: ApJ††facilities: INTEGRAL, CGRO, NuSTAR, Hitomi, Fermi, SMILE,
GRAMS, AMEGO, e-ASTROGAM, FORCE

##  1 Introduction

Elements heavier than Bi exist in our universe, but their origin remains a
mystery. Most cosmic isotopes heavier than the iron group are expected to be
created by the rapid-neutron capture process, also known as the r-process
(Cameron, 1957; Burbidge et al., 1957; Cowan et al., 1991; Wanajo & Ishimaru,
2006; Qian & Wasserburg, 2007; Arnould et al., 2007), but the actual
nucleosynthesis sites capable of achieving such neutron-rich environments
remains a matter of debate. Before the discovery of binary neutron-star
mergers (NSM) observed as gravitational-wave objects like GW170817 (Abbott et
al., 2017), NSMs were considered to be more promising as r-process
nucleosynthesis sites than other primary candidates, such as core-collapse
supernovae (SNe), because NSMs could achieve more neutron-rich (lower electron
fraction $Y_{\rm e}$) environments (Wanajo et al., 2011; Lattimer & Schramm,
1974; Metzger et al., 2010). The event rate of NSMs is much lower than that of
SNe, but the yield of $r$-process nuclei in one event is expected to be very
high (Wallner et al., 2015; Hotokezaka et al., 2015). Observational evidence
of the existence of $r$-process nuclei has already been obtained by infrared
observations of kilonovae (also called macronovae or $r$-process novae) in
some short gamma-ray bursts, such as GRB 130603B (Tanvir et al., 2013) and the
gravitational wave event GW170817 (Villar et al., 2017). However, the infrared
radiation from NSMs is, in principle, the result of indirect emissions from
unstable $r$-process nuclei, and any hint of elements heavier than the
lanthanoids is still missing from the infrared information. Given that the
nuclear levels of nuclei are in the MeV energy range, the $\gamma$ rays from
$r$-process nuclei should be the best probe for searching for $r$-process
sites in the universe.

According to theoretical estimates of the $\gamma$-ray flux from binary NSMs
(Hotokezaka et al., 2016), the $\gamma$-ray radiation immediately following a
merging event is very dim at about $10^{-8}$–$10^{-7}$ photons s-1 cm-2 keV-1,
even at an extremely close distance $d$ of 3 Mpc. This flux is comparable with
or below what the sensitivities of current and near-future MeV missions can
detect. The precise measurements of photon energies are, in principle, rather
difficult in the MeV band, where Compton scattering dominates over the photon-
absorption process. Therefore, the ability to detect $\gamma$ rays from NSMs
by an immediate follow-up observation (a Target-of-opportunity observation;
ToO) would be limited by the sensitivity of the $\gamma$-ray instruments.
Instead, a non-ToO observation of $\gamma$ rays from long-lived nuclei in NSMs
would be an alternative way to survey $r$-process sites, and this has been
proposed by Wu et al. (2019) and Wang et al. (2020). The $\gamma$-ray
luminosity from nuclei with long lifetimes, on the order of $10^{3}$–$10^{6}$
years, becomes much lower than that from short-lived nuclei, but if we limit
the survey area within our Galaxy ($d\lesssim$ 10 kpc), then the $\gamma$-ray
flux in non-ToO observations is expected to become comparable with what is
required for ToO observations. Therefore, non-ToO observations should provide
more sensitive $\gamma$-ray surveys of NSMs, because the exposure time (the
accumulation time of signals) is not limited like it is in ToO observations.
Another benefit from performing a non-ToO survey is the better identification
of $\gamma$-ray lines; we expect the effect of Doppler broadening to be
smaller for older NSM remnants than for very young NSMs.

Here we focus on the non-ToO survey of $\gamma$ rays from $r$-process nuclei
in a possible Galactic NSM remnant. In this paper, we estimate $\gamma$-ray
emissions from Galactic NSM remnants in an older age range than in previous
work (Hotokezaka et al., 2016; Wang et al., 2020) by using nuclear-network
numerical calculations with a complete nuclear database. This paper also aims
to provide $\gamma$-ray diagnostic methods for NSMs, showing the required
sensitivities for future $\gamma$-ray observatories. In our study, we assume
that the $\gamma$-ray instruments have a wider field-of-view (FOV) than the
object size of the NSM remnants, which are larger than early NSMs in a ToO
observation. We also assume that the instruments accumulate all of the
$\gamma$-ray emissions from the NSM remnants, even though the nuclei may mix
with the circumstellar medium (CSM) during the evolution of the remnants.

The rest of this paper is organized as follows. In Section 2, we summarize our
environments and procedures for the nuclear-network numerical calculation and
show the results for $\gamma$-ray emissions from NSM remnants. In Section 3,
we present the $\gamma$-ray diagnostics that utilize spectral color to
identify NSM remnants and provide the line properties for estimating the age
$t$ and $Y_{\rm e}$. In Section 4, we discuss the survey distance and coverage
in our Galaxy permitted by the instrument sensitivities, the corresponding
limitation of the NSM rate in our Galaxy, and expectations for future
missions.

##  2 Numerical estimation of Gamma Rays from Neutron-Star Merger Remnants

###  2.1 Overview of Numerical Calculation

To estimate the $\gamma$-ray emissions from binary NSM remnants of various
ages, we performed a numerical simulation comprising the following three
steps: 1) calculation of the mass distribution of $r$-process nuclei for NSMs
at $t=1$ year, 2) calculation of the decay processes of unstable nuclei
emitting $\gamma$ rays, and 3) a simple calculation of the radiation transfer
of $\gamma$ rays from NSMs.

Figure 1: The r-process nuclei in the NSM case at $t=1.0$ year on the table of
nuclides. The magic numbers for neutrons and protons are indicated by red
dashed lines Figure 2: Mass number distribution of the abundance for the NSM
case at $t=1.0$ year, with $Y_{\rm e}\sim$ 0.1, 0.2, 0.3, and 0.4, shown in
red, magenta, blue, and green, respectively.

For the first step, we adopted the nucleosynthesis calculation for around
4,070 nuclei performed by Fujimoto et al. (2007), which was cooled using the
adiabatic expansion modeled from Freiburghaus et al. (1999) to provide the
elemental composition ratios of nuclei for $Y_{\rm
e}=0.10,0.15,0.20,0.25,0.30,0.35,0.40$, and $0.45$. This estimation assumes
that the initial environment has a temperature of $9\times 10^{9}$ K, radius
of 100 km, entropy per baryon of 10$k_{\rm B}$, where $k_{\rm B}$ is the
Boltzmann constant, and velocity of $2\times 10^{9}$ cm s-1, along with the
initial abundances of the 4,070 nuclei in nuclear statistical equilibrium. As
a result, the calculation provides the mass fractions at t = 1 year evaluated
with the nuclear reaction network (network A in Fujimoto et al. 2007), by
using $Y_{\rm e}=0.10$–$0.45$ in steps of 0.05. To set up the mass
distribution of nuclei for the NSMs at $t=1$ year, we blended the nuclei with
the mass fraction using the $Y_{\rm e}$ provided in Wanajo et al. (2014).
Specifically, the fractions are 4.54%, 4.85%, 14.6%, 29.7%, 10.3%, 25.1%,
10.5%, and 0.33% for $Y_{\rm e}=0.10,0.15,0.20,0.25,0.30,0.35,0.40$, and
$0.45$, respectively. Note that this $Y_{\rm e}$-fraction model by Wanajo et
al. (2014) describes slightly-less neutron-rich environment than those by the
recent dynamical-ejecta models after the kilonova observations of the
gravitational event GW170817, such as four models in Kullmann et al. (2022)
under two kinds of equation-of-states, density dependent 2 (DD2) (Hempel &
Schaffner-Bielich, 2010; Typel et al., 2010) and SFHo (Steiner et al., 2013).
In this paper, we adopted the first one by Wanajo et al. (2014) as a
pessimistic case for the $r$-process site, but changing the $Y_{\rm
e}$-fraction models does not change the conclusions from the $\gamma$-ray
spectra as tested in the later section 3. Figure 1 shows the mass fraction of
multiple nuclei at $t=1$ year generated in an NSM case, information that is
given in the table of nuclides (neutron number $N$ vs atomic number $Z$).
Using the same data set, Figure 2 summarizes the distribution of nuclei with
mass number $A$ at $t=1$ year, showing the contributions of $Y_{\rm e}$. This
plot demonstrates that the environment with lower $Y_{\rm e}$ contributes to
the generation of heavier elements.

For the second step, we simulated the decay processes of unstable nuclei,
starting from the mass distribution at $t=1$ year calculated in the first
step. We used the Decay Data File 2015 (DDF-2015) (Katakura & Minato, 2016) in
the Japanese Evaluated Nuclear Data Library (JENDL) (Katakura, 2012), which
provides the decay profiles of 3,237 nuclei up to Z=104 (Rf). Here we applied
a correction to the $\gamma$-ray information for 241Am; this error was
reported from our study and was fixed in the next version of the database. The
originality of this study lies in the comprehensiveness of nuclei treated in
the calculation. In the nuclear-decay calculation, we adopted the
$\alpha$-decay, $\beta^{-}$-decay, $\beta^{+}$-decay, electron capture,
isomeric transition, and $\gamma$-decay processes. In our calculation, the
internal conversion process is ignored, which emits soft X-rays and makes a
negligible contribution to the $\gamma$-ray band. The neutron- and proton-
emission processes are also ignored because they contribute to the very early
phase, which is out of the scope of this study. The spontaneous fission
process may occur on 257Es, 256Cf, 254Cf, and 250Cm, but its contribution is
negligible. Thus, this process is also excluded from our calculation. In
addition, we do not calculate the $\gamma$-ray emission from secondary
electrons (electrons from $\beta$-decay, $\delta$-rays, and so on) after the
decay of unstable nuclei.

Figure 3: The relative abundance of 195Pt in the NSM case at $t=$ 1 thousand,
1 million, and 1 billion years, shown in green, blue, and thick black,
respectively, compared with the semi-empirical abundance distribution of Beer
et al. (1997).

To verify the calculations in the second step, we refer to Figure 3, which
represents the relative abundances of nuclei at $t=$ 1 thousand, 1 million,
and 1 billion years as a function of $A$. The distribution in $A$ does not
change dramatically after $t\sim$ 1 thousand years, except for the
disappearance of the small dips at the magic numbers. The distribution for
$t>1$ million years becomes roughly consistent with the semi-empirical
abundance distribution of cosmic $r$-process nuclei in Beer et al. (1997),
which is close to the solar abundance distribution.

For the outputs of the second step, we obtain the $\gamma$-ray flux
$F_{\gamma,i}$ from the NSM at distance $d$. Using the nuclear $\gamma$-ray
intensity $I_{\gamma,i}$ of the $i$-th element with the mass number $A_{i}$,
the mass $M_{i}$, and the half-life $T_{1/2}$, the $F_{\gamma,i}$ in a small
time interval $dt$ is described as

$F_{\gamma,i}=\frac{N_{\rm A}}{4\pi
d}\frac{M_{i}}{A_{i}}\frac{I_{\gamma,i}}{T_{1/2}}\left(\frac{1}{2}\right)^{-dt/T_{1/2}},$
(1)

where $N_{\rm A}$ is Avogadro’s number.

Finally, the third step is to calculate the transfer of $\gamma$ rays through
the NSM ejecta. However, we omitted the detailed Monte-Carlo calculation of
the radiation transfer because the optical depth decreases rapidly after the
merger event, by roughly $\propto\sqrt{t}$ (Li, 2019); within the scope of our
study at $t\gg 1$ year, the optical depth is thin and negligible. Therefore,
the degradation of the line profiles by Compton scattering is not included in
our calculation, which would be dominant in only the very early phase. Note
that the detail calculations of the MeV $\gamma$-ray spectra from Galactic
NSMs in the initial phase were performed by Wang et al. (2020);
2021ApJ...923..219W. In this step, we apply only the bulk Doppler-broadening
effect caused by the expansion velocity $v(t)$. The thermal Doppler-broadening
effect is ignored in this calculation because it is two or three orders of
magnitude smaller than that from the expansion motion of the heavy elements in
the $A=$ 50–200 range. In reality, the line profile from the bulk Doppler
effect becomes complicated due to the complex contributions of various
velocity components, as has been observed in the X-ray lines from heavy
elements in SN remnants (Grefenstette et al., 2017; Kasuga et al., 2018). For
simplicity, we applied the Gaussian distribution function for the line profile
in the calculation. Of the various velocity elements in the remnant, we
applied only single Gaussian broadening to the maximum velocity component,
which we assume to be the forward shock motion. This was done to simulate the
most robust case for considering the $\gamma$-ray sensitivity. In our
assumption, $v(t)$ starts from the initial value $v(0)=0.3$ $c$, where $c$ is
the speed of light, and evolves at a constant rate during the free expansion
phase. During the Sedov-Taylor phase, $v(t)$ evolves as $v(t)\propto
t^{-(3/5)}$ (Taylor, 1950), and then as $v(t)\propto t^{-(0.7)}$ during the
pressure-driven snowplow (PDS) phase (McKee & Ostriker, 1977). We assume that
the free expansion, Sedov-Taylor, and PDF phases end at $t=$ 10, $4.7\times
10^{4}$, and $1.65\times 10^{6}$ years, respectively. The ages of these phase
transitions may change by about one order of magnitude due to differences in
density of the CSM, but this modification affects only the Doppler-broadening
effect. It becomes negligible when compared with the typical energy
resolutions of $\gamma$-ray instruments for ages older than $t\sim 1,000$
years, the range that lies within the scope of this study. Note that even at
$t=10^{6}$ years, $v(t)$ approaches $\sim 20$ km s-1, which is about double
the speed of sound for a typical CSM density of 0.01 cm-3. The radius becomes
$\sim 100$ pc. Finally, we get the $\gamma$-ray spectra for NSMs at $t$,
accumulated from all of the $r$-process nuclei in the ejecta.

###  2.2 Gamma-ray Emission and Evolution

From the numerical calculation described in Section 2.1, the $\gamma$-ray
spectra from $t=3$ to 1 million years, under the assumption that the ejecta
mass is $M_{\rm ej}=0.01M_{\odot}$ at $d=10$ kpc, are summarized in Figures 4
and 5. As described in Section 1, we assume that all of the emissions from the
NSM remnants are observable within the wider FOV; this is assumed to be larger
than the object size, which becomes around 10 pc at $t=1,000$ years and
expands into around 100 pc at $t>10^{6}$ years. The spectra contain many
nuclear lines broadened by the Doppler effect. They appear to form a
continuous spectrum in the early phase, but they become separated at ages
older than 1,000 years. Note that the $\gamma$-ray data without the Doppler-
broadening effect (the outputs from the second step of the calculation in
Section 2.1) is provided as the numerical model for the XSPEC tool (Arnaud,
1996) in the HEAsoft package (Appendix A).

Figure 4: Gamma-ray spectra simulated for the NSM case at $t=$ 3 - 1,000
years, assuming a distance of 10 kpc with an initial velocity of 0.3 $c$ (see
the text). The red and black plots represent the spectra with and without the
Doppler-broadening effect, respectively. The time since the merging event and
the velocities are shown in each top label.

Figure 5: Same as Figure 4 but for $t=$ 3,162–1,000,000 years. Energy (keV) | Nuclei | Major $Y_{\rm e}$ | Half-life
---|---|---|---
40.0 $\ddagger$ |  225Ra | 0.10–0.20 |  $\sim 7.9\times 10^{3}$ yr,
|  |  |  $1.6\times 10^{5}$ yr,
|  |  |  $2.1\times 10^{8}$ yr
58.6 $\ddagger$ |  60mCo | 0.35–0.40 |  $\sim 2.6\times 10^{6}$ yr
59.5 |  241Am | 0.10–0.20 |  $\sim 4.3\times 10^{2}$ yr,
|  |  |  $8.4\times 10^{3}$ yr
74.7 $\ddagger$ |  243Am | 0.15–0.20 |  $\sim 7.4\times 10^{3}$ yr
87.6 $\dagger$ |  126Sn | 0.10–0.30 |  $\sim 2.3\times 10^{5}$ yr
106.1 $\ddagger$ |  239Np | 0.10–0.20 |  $\sim 4.5\times 10^{9}$ yr
236.0 |  227Th | 0.10 |  $\sim 21$ yr,
|  |  |  $3.2\times 10^{4}$ yr,
|  |  |  $7.0\times 10^{8}$ yr
276.0 |  81Kr | 0.45 |  $\sim 2.2\times 10^{5}$ yr
311.9 |  233Pa | 0.10–0.20 |  $\sim 2.1\times 10^{6}$ yr
328.4 |  194Ir | 0.10–0.15 |  $\sim 6$ yr
427.9 $\dagger$ |  125Sb | 0.25–0.30 |  $\sim 2.7$ yr
440.5 $\ddagger$ |  213Bi | 0.10–0.20 |  $\sim 7.3\times 10^{3}$ yr
511.9 |  106Rh | 0.30–0.35 |  $\sim 1.1$ yr
513.9 $\dagger$ |  85Kr | 0.30–0.40 |  $\sim 10$ yr
561.1, 834.5 |  92Nb | 0.45 |  $\sim 3.4\times 10^{7}$ yr
609.3 $\ddagger$ |  214Bi | 0.10–0.20 |  $\sim 1.6\times 10^{3}$ yr,
|  |  |  $\sim 7.5\times 10^{4}$ yr
661.7 $\dagger$ |  137mBa | 0.10–0.25 |  $\sim 30$ yr
765.8 |  95Nb | 0.30–0.35 |  $\sim 34$ d, 64 d
871.1, 702.6 |  94Nb | 0.45 |  $\sim 2.0\times 10^{4}$ yr
1115.5 |  65Zn | 0.45 |  $\sim 243$ d
1157.0 |  44Sc | 0.45 |  $\sim 60$ yr
1332.5 |  60Co | 0.45 |  $\sim 5.2$ yr
Table 1: List of nuclei and their most prominent nuclear $\gamma$-ray lines,
arranged by $Y_{\rm e}$. The $\dagger$ and $\ddagger$ marks indicate the lines
used in the $Y_{\rm e}$ diagnostics in Section 3.3 for young ($t<100$ years)
and old ($t>100$ years) cases, respectively. Note that the ’m’ in 137mBa
indicates that it is meta-stable.

To identify the $\gamma$-lines in the spectra, we checked the most prominent
lines in the $\gamma$-ray spectra generated by a single $Y_{\rm e}$ condition.
Table 1 lists the brightest lines shown for the nuclei in each $Y_{\rm e}$.
Roughly speaking, the bright lines seen for objects of a younger age lie in
the higher-energy $\gamma$-ray band of the spectrum. Further details of the
diagnostics will be discussed in Section 3.

##  3 Gamma-ray Diagnostics of Neutron-Star Merger Remnants

###  3.1 Spectral Color Changes of Neutron-Star Merger Remnants

Using the energy spectra of the NSM remnants (Figures 4 and 5), we first
checked the properties of the spectral shapes from the hard X-ray to the soft
$\gamma$-ray bands. As shown in the normalized spectra plotted in Figure 6,
the energy spectra roughly evolve from hard to soft slopes. $\gamma$-ray
emission decreases rapidly leaving the hard X-ray emission in old age, as is
indicated in Table1. This phenomenon of $\gamma$ rays is equivalent to the
Sargent law for $\beta$ decay.

Figure 6: The energy spectra with Doppler broadening, same as the red plots in
Figures 4 and 5 but normalized to the 11.5 keV flux at 30, 100, 300, 1000,
3000, 10000, 30000, 100000, 300000, and 1000000 years, which are shown in
black, red, green, blue, light blue, magenta, yellow, orange, yellow green,
and olive green, respectively.

To see the evolution of the shape of the $\gamma$-ray spectra more
quantitatively, we plotted the light curves of the $\gamma$-ray flux in three
bands: 70 to 200 keV, 200 to 500 keV, and 500 to 3,000 keV, which cover
multiple lines around 100 keV and 300 keV, and a prominent line around 700
keV, respectively. As indicated in the top panel of Figure 7, the flux in the
higher-energy bands decreases more quickly than that in the low-energy bands.
A decaying trend is also seen in the time dependency of the hardness ratio
among these bands, as is indicated in the lower panel of Figure 7. The ratio
drops dramatically at around 200–300 years, indicating that the $\gamma$-ray
flux above the 500 keV band quickly decreases at this age. This phenomenon is
primarily due to the decay of 125Sb and 137mBa listed in Table 1. Note that
this result does not change even if we adopt other $Y_{\rm e}$-fraction models
of DD2-125145, DD2-135135, SFHo-125145, and SFHo-135135 in Kullmann et al.
(2022), as shown in Figure 7.

Figure 7: The top panel represents the light curve of the $\gamma$-ray spectra
simulated for the case of an NSM at a distance of 10 kpc, with the Doppler-
broadening effect and an initial velocity of 0.3 $c$. The energy bands of
70–200 keV, 200–500 keV, and 500–3000 keV, are shown in red, blue, and yellow,
respectively. The bottom panel represents the ratios among the three bands in
the top panel: the ratio between the 70–200 keV and 200–500 keV bands, that
between the 200–500 keV and 500–3000 keV bands, and that between the 70–200
keV and 500–3000 keV bands are shown in purple, green, and orange,
respectively. The dim colors represent the results using the $Y_{\rm
e}$-fraction by Wanajo et al. (2014) as the nomical condition of this paper,
and the lighter color lines represent the results by the other models,
DD2-125145, DD2-135135, SFHo-125145, and SFHo-135135 in Kullmann et al.
(2022).

To compare the spectral shape of NSM remnants with other astronomical objects,
we plotted the color-color diagrams in the hard X-ray band (10–500 keV) and in
the hard X-ray to $\gamma$-ray band (70–3000 keV), in the top and bottom of
Figure 8, respectively. We divided the energy band-pass for these spectra into
three ranges: 10–25, 25–70, and 70–500 keV for the hard X-ray band (top of
Figure 8), and 70–500, 500–1,000, and 1,000–3,000 keV for the hard X-ray to
$\gamma$-ray band (Figure 8 bottom). Note that the divisions of the energy
bands are defined so that they follow the energy band-pass of current
$\gamma$-ray instruments on board NuSTAR (Harrison et al., 2013), INTEGRAL
(Winkler et al., 2003), and other observatories. For comparison, the spectral
colors of other astronomical objects, calculated using the INTEGRAL catalog
version 0043 111https://www.isdc.unige.ch/integral/science/catalogue, are also
plotted in the same figures. In the hard X-ray band (the 10–500 keV band in
the top of Figure 8), the spectra of NSM remnants older than $t\sim$ 1,000
years have spectral colors similar to those of supernova remnants or active
galactic nuclei, but NSM remnants younger than $t\sim$ 1,000 years can be
distinguished from other known objects by their hard X-ray colors. In other
words, the spectral color in the hard X-ray band below 500 keV is a good
indicator of young NSM remnants. Furthermore, this differentiation from known
objects becomes more prominent when we include the higher-energy band covering
the MeV portion of the spectrum, as is clearly indicated in the bottom of
Figure 8. Note that this result does not change even if we adopt other $Y_{\rm
e}$-fraction models of DD2-125145, DD2-135135, SFHo-125145, and SFHo-135135 in
Kullmann et al. (2022), as shown in Figure 8. Therefore, NSM remnants have
unique spectral colors in the hard X-ray to $\gamma$-ray bands. This
observation is one of the important conclusions from our calculation. Note
that the spectral models in the INTEGRAL catalog are simple enough that the
colors of known objects in the gamma-ray band (bottom of Figure 8) are less
scattered than those in the hard X-ray band (top of Figure 8). The spectral
separation between NSM remnants and other objects in the bottom of Figure 8
does not change dramatically, even if we lower the low-energy threshold (70
keV in the bottom of Figure 8) to cover 20 keV, for example. However, it
becomes worse if we set it higher so that everything up to a certain point,
200 keV, for example, is ignored. This implies that hard X-rays around 100 keV
provide the key information for distinguishing NSM remnants from other
objects. Note that these results are based on the pure-nuclear $\gamma$ rays
from $r$-process nuclei in NSMs, and thus the synchrotron radiation from
electrons that are accelerated by the shocks may contaminate the hard X-ray
band for young remnants. Additionally, when taking actual observations, we
must be careful to isolate the contamination of the hard X-ray spectrum that
arises from other objects located behind the NSM, such as active galactic
nuclei within the FOV.

Figure 8: (Top) A color-color diagram of the flux ratio between the 10–25 keV
and 25–70 keV bands versus that between the 25–70 keV and 70–500 keV bands is
shown using simulated $\gamma$-ray spectra for the case of an NSM at a
distance of 10 kpc, with the Doppler-broadening effect and an initial velocity
of 0.3 $c$. The spectral evaluations of the NSM both with and without K X-ray
emission are shown by magenta and gray lines, respectively, and those by the
other $Y_{\rm e}$-fraction models, DD2-125145, DD2-135135, SFHo-125145, and
SFHo-135135 in Kullmann et al. (2022), are shown in the light-pink lines, for
reference. The color-color diagram of the X-ray objects listed in the INTEGRAL
catalog (version 0043) are also plotted: X-ray binaries with the purple
squares, galactic compact objects, such as cataclysmic variables, with the
green squares, supernova remnants with the cyan circles, active galactic
nuclei with the orange circles, and unidentified objects with the red
triangles. (Bottom) Same plot as the top panel but in the 70–500, 500–1,000,
and 1,000–3,000 keV bands.

###  3.2 Nuclear Line Emissions from Older Neutron-Star Merger Remnants

For ages older than $t>$ 3,000 years, nuclear lines are clearly seen in the
$\gamma$-ray spectra of NSM remnants due to the minimal Doppler-broadening
effect, as is shown in Figures 4 and 5. Using the $\gamma$-ray spectra of NSM
remnants that were shown in Section 2 (i.e., the $Y_{\rm e}$ distribution for
the NSM case with $M_{\rm ej}=0.01M_{\odot}$ at $d=10$ kpc), we selected the
brightest nuclear lines in each energy band, 3–75 keV, 75–500 keV, and
500–4000 keV, for at least one epoch in the age range spanning $t=$10 to
$4\times 10^{6}$ years. Note that these energy bands are defined in such that
they simulate the energy bands that are observable by current and near-future
instruments aboard satellites, such as the hard X-ray focusing missions NuSTAR
(Harrison et al., 2013) and FORCE (Nakazawa et al., 2018), and $\gamma$-ray
missions like INTEGRAL (Winkler et al., 2003), e-ASTROGAM (De Angelis et al.,
2017; de Angelis et al., 2018), AMEGO (Kierans, 2020), and GRAMS (Aramaki et
al., 2020).

Figure 9 presents the time evolution of the brightest nuclear $\gamma$-ray
lines in these energy bands. To account for the reduction in the line
sensitivities as a result of the Doppler-broadening effect, we accumulated the
photons that were within the energy resolution of $\Delta E=3\pm 1$ keV from
the center energy of their associated lines. This chosen value for the energy
resolution is typical for semiconductor $\gamma$-ray detectors. For reference,
the evolution of lines without Doppler broadening is also shown in the figure
as dashed lines. As indicated in Figure 9, the Doppler-broadening effect
becomes less dominant in the hard X-ray band after a few hundred years, but it
is still present until about $t=10^{3}$ and $10^{4}$ years in the soft
$\gamma$-ray and the hard $\gamma$-ray bands, respectively. Note that the
reason why several lines, such as those of 126mSb and 239Np, increase as $t$
approaches $10^{3}$–$10^{5}$ years is that the number of parent nuclei
increases in these phases.

Figure 9: The time evolution of the brightest $\gamma$-ray lines in the energy
bands of 3–75 keV, 75–500 keV, and 500–4000 keV are shown in the top, middle,
and bottom panels, respectively. The line fluxes with and without the Doppler
effect are shown as thick and dashed lines, respectively. The thick lines plot
the $\gamma$-ray fluxes within the $\Delta E=3$ keV width, and the upper and
lower boundaries of the corresponding hatched areas show the flux within the
$\Delta E=2$ keV and 4 keV widths, respectively. The nuclei and energies are
indicated by the same color as the lines and hatches.

From Figure 9, we can identify the nuclear lines that are useful as indicators
for the ages of NSMs. The ages can be categorized into three epochs: $t<100$
years, $t\sim 10^{3}$– $10^{4}$ years, and $t>10^{4}$ years. In summary, if we
detect the lines from 125Sb, 194Os, 227Th, or 194Ir, then we can determine the
age of the NSM to be very young at $t<100$ years. Similarly, lines from 137mBa
in the $\gamma$-ray band indicate that the age is around $t\sim 10^{2}$ years.
In the age range spanning $t\sim 10^{3}$– $10^{4}$ years, nuclear lines will
be detected from 241Am, 243Am, 214Pb, 239Np, and/or 214Bi. A nuclear line from
126mSb indicates that the NSM is very old at $t>10^{4}$ years. In the wide age
range from $t=400$ to $10^{5}$ years, the line from 126Sn stays almost
constant at $10^{-9}$ photons s-1 cm-2 for a distance of $d=10$ kpc, and thus
it can be used as a standard candle for measuring $d$.

###  3.3 Line Diagnostics for the electron fraction

In addition to the spectral colors (Section 3.1), nuclear lines can be used to
identify NSM remnants among astronomical objects, especially when the remnants
are of an older age. Since the NSMs are thought to have both a more neutron-
rich environment and a lower $Y_{\rm e}$ condition than SNs (Wanajo et al.,
2011; Lattimer & Schramm, 1974; Metzger et al., 2010), a new line-diagnostic
method utilizing $Y_{\rm e}$ values will be useful for distinguishing NSMs
from SNe. In this subsection, we search for $\gamma$-ray line diagnostics for
$Y_{\rm e}$. We use the $\gamma$-ray spectra calculated under the pure $Y_{\rm
e}$ conditions in the $Y_{\rm e}=0.10$–0.45 range, whereas in the previous
sections we used the mixed $Y_{\rm e}$ condition for NSMs.

To identify the best candidates among the nuclear $\gamma$-ray lines for the
identification of $Y_{\rm e}$, we first selected the five brightest lines for
each age, $t=100,1000,10^{4},10^{5}$, and $10^{6}$ years, and for each $Y_{\rm
e}$ ($=0.10,0.15,0.20,0.25,0.30,0.35,0.40$, and $0.45$). Then among these $5$
(ranks) $\times 5$ ($t$) $\times 7$ ($Y_{\rm e}$) lines, we selected the
nuclear lines which appeared in two or more of the conditions for $t$ and
$Y_{\rm e}$. In total, ten $\gamma$-ray lines are selected and are marked as
$\dagger$ and $\ddagger$ in Table 1 for $t<$ 100 years and $t>$ 100 years,
respectively. Therefore, the lines from 137mBa (661.66 keV), 85Kr (513.9 keV),
and 125Sb (427.87 keV) in ages below 100 years are indicators of low, middle,
and high $Y_{\rm e}$ environments, respectively. Here, low, middle, and high
are numerically defined as $Y_{\rm e}\sim 0.10$–$0.20$, $0.20$–$0.35$, and
$0.35$–$0.45$, respectively. In ages older than $t=$100 years, the lines from
225Ra (40.0 keV), 243Am (74.66 keV), 239Np (106.1 keV), 213Bi (440.5 keV), and
214Bi (609.3 keV) are emitted from a low $Y_{\rm e}$ environment, whereas the
lines from 60mCo (58.6 keV) and 126Sn (87.6 keV) become bright in the middle
and high $Y_{\rm e}$ environments, respectively.

Figure 10: (Top) The $Y_{\rm e}$ dependencies of the line intensities of the
nuclei listed with $\dagger$ marks in Table 1 are shown in the upper panel.
The reduction due to the Doppler-broadening effect is not considered. The line
intensities of 137mBa, 85Kr, and 125Sb are shown in red, black, and blue,
respectively. Thick, dotted, and dashed lines represent the intensities at
$t=$ 10, 100, and 1000 years. The ratios among the lines are plotted in the
lower panel; the ratios between 137mBa and 85Kr, 125Sb and 137mBa, and 85Kr
and 125Sb are shown in purple, orange, and cyan, respectively. (Middle) Same
plot as the top panel but for older ages ($t>100$ years, $\ddagger$ marks in
Table 1). The dependencies of 243Am (and 239Np), 214Bi (and 225Ra, 213Bi),
60mCo, and 126Sn are shown in red, green, black, and blue, respectively, and
the thick, dotted, dashed, and dotted dash lines represent the intensity at
$t=$ 1000, $10^{4}$, $10^{5}$, and $10^{6}$ years, respectively, in the upper
panel. The ratios in the lower panel between 243Am and 60mCo, between 214Bi
and 60mCo, between 126Sn and 243Am, between 126Sn and 214Bi, and between 60mCo
and 126Sn are shown in purple, brown, orange, dark-yellow, and cyan,
respectively. (Bottom) The time dependencies of the line ratios of
137mBa/85Kr, 125Sb/137mBa, 85Kr/125Sb, 243Am/60mCo, 126Sn/243Am, and
60mCo/126Sn are shown in the thin purple, thin orange, thin cyan, thick
purple, thick orange, and thick cyan lines, respectively. The NSM case by the
mass fraction of Wanajo et al. (2014) and the case of $Y_{\rm e}=0.45$ are
shown in the straight and dotted lines, respectively. The results by the other
NSM models, DD2-125145, DD2-135135, SFHo-125145, and SFHo-135135 in Kullmann
et al. (2022), are shown in the lighter colors.

Since the absolute flux of a single line changes with respect to $t$ and $d$,
the ratio between two or more lines should be a good indicator for $Y_{\rm
e}$. Figure 10 summarizes the line intensities and their ratios using the ten
nuclei selected above. For simplicity, the plots for young and old ages, $t=$
10–$10^{3}$ years and $10^{3}$–$10^{6}$ years, respectively, are shown
separately. In the young age range (top of Figure 10), the ratios of
137mBa/85Kr, 125Sb/137mBa, and 85Kr/125Sb become larger than unity in the low,
middle, and high $Y_{\rm e}$ environments, respectively. These low and high
$Y_{\rm e}$ indicators (i.e., 137mBa/85Kr and 85Kr/125Sb, respectively)
exhibit more prominent ratios over time since 85Kr decays slower than 125Sb
and faster than 137mBa, whereas the middle indicator (125Sb/137mBa) becomes
dim after 100 years. Note that the plots use the incident line fluxes
calculated in step 2 of Section 2.1, and the reduction due to the Doppler-
broadening effect is not considered. The Doppler effect is particularly
significant in plots for the young age range (top of Figure 10).
Quantitatively, the ratios of 137mBa/85Kr, 125Sb/137mBa, and 85Kr/125Sb change
by factors of $1.57,0.82$, and $0.77$, respectively, for $t<1,000$ years. In
the old age range (middle of Figure 10), the ratios of 243Am/60mCo,
126Sn/243Am, and 60mCo/126Sn indicate the low, middle, and high $Y_{\rm e}$
environments, respectively. The line from 239Np has the same flux and time
evolution as that from 243Am (red plots), because they are in the same decay
chain. Similarly, the lines from 214Bi, 225Ra, and 213Bi (green plots) follow
almost the same trend as those from 243Am and 239Np (red plots). Among them,
the low $Y_{\rm e}$ indicator (243Am/60mCo) in the old age range is valid up
to $t=1$ million years, and the middle $Y_{\rm e}$ indicator (126Sn/243Am)
shows more significant ratios with $t>10^{3}$ years. On the other hand, lines
for the high $Y_{\rm e}$ indicator 60mCo/126Sn decay quickly and become
unavailable after $10^{3}$ years; that is, if the ratio 60mCo/126Sn is larger
than unity, then the object is in a high $Y_{\rm e}$ environment with an age
of $t\sim 10^{3}$ years. In summary, using these indicators, which become
larger than unity in specific $Y_{\rm e}$ conditions, we can estimate the
$Y_{\rm e}$ environment independently from the spectral-color diagnostics
shown in Section 3.1.

Finally, we checked the line ratios blended by $Y_{\rm e}$ distributions of
the NSM cases. The time evolution is plotted in Figure 10 bottom. The
difference of the $Y_{\rm e}$-fraction models between Wanajo et al. (2014) and
Kullmann et al. (2022) does not affect the trend of the NSMs so much. These
$\gamma$-ray lines are also expected to be observed from the remnants of core-
collapse SNe, which are considered to be less neutron-rich environment at
$Y_{\rm e}\sim 0.5$ (Andrews et al., 2020) than NSMs. However, our numerical-
calculation model in this paper has limitations in estimating $\gamma$-ray
radiation from the core-collapse SNe, because the mass fractions of the
$r$-process nuclei in the ejecta are different between the maximum $Y_{\rm e}$
condition of our calculation (i.e., $Y_{\rm e}=0.45$) and the SNe case
($Y_{\rm e}\sim 0.50$), and the neutron-rich nuclei in the nominal core-
collapse SNe are predominantly generated via the $s$-process rather than the
$r$-process. For reference, we plotted the time evolution of the line ratios
in the $\gamma$-ray spectra of $Y_{\rm e}=0.45$, which should still reproduce
well an environment with almost-equal numbers of neutrons and protons.
According to Figure 10 bottom, we expect the low-$Y_{\rm e}$ indicators (i.e.,
137mBa/85Kr and 243Am/60mCo) in the NSM case become many orders-of-magnitude
larger than those in the SNe case. In the core-collapse SNe where neutron-rich
nuclei are generated via the $s$-process, relatively large amount of 85Kr and
almost no 243Am are synthesized. Therefore, the difference of these
low-$Y_{\rm e}$ indicators between the NSMs and SNe cases are expected to
become larger than Figure 10 bottom. As for the middle-$Y_{\rm e}$ indicators
(125Sb/137mBa and 126Sn/243Am), they may not be useful to distinguish
$\gamma$-rays from NSMs and SNe according to Figure 10 bottom. In the
$s$-process environment of core-collapse SNe, almost no 137mBa and 243Am are
synthesized and thus theose middle-$Y_{\rm e}$ indicators can be larger than
the values in the figure. Finally, the high-$Y_{\rm e}$ indicators (85Kr/125Sb
and 60mCo/126Sn) can discard the NSMs from the SNe cases as indicated by
Figure 10 bottom. In summary, the new line-diagnostic method for $Y_{\rm e}$
provides a tool for distinguishing between NSMs and SNe.

##  4 Discussion

In Section 2, we presented a nuclear-decay simulation using a large nuclear
database, the goal of which was to estimate the $\gamma$-ray spectra of NSMs
up to the age of $t=10^{6}$ years. We have identified many nuclear lines,
listed in Table 1, that can be used for identifying the nucleosynthesis
environments of NSMs, even with the Doppler-broadening effect altering the
profiles of these lines in the young age range. In Section 3, we numerically
analyzed the simulated $\gamma$-ray spectra from NSMs and found that the
spectral slope in the soft $\gamma$-ray band above 500 keV changes at around
$t=$ 200–300 years. We also found that the spectral colors of NSMs in the hard
X-ray to soft $\gamma$-ray bands differ from those of other astronomical
objects up to $t=10^{5}$ years old. Consequently, we can identify a
$\gamma$-ray object as an NSM remnant using the $\gamma$-ray spectral colors
(Section 3.1). Among the many nuclear lines in the spectra, we identified that
the nuclear lines from 241Am, 243Am, 214Pb, 239Np, and 214Bi are prominent for
$t=10^{3}$–$10^{4}$ years, and that the lines from 126Sn and 126mSb are
prominent for $t>10^{4}$ years (Section 3.2). In addition, we proposed a new
line-diagnostic method for distinguishing $Y_{\rm e}$ environments that uses
the line ratios of 137mBa/85K and 243Am/60mCo, which become larger than unity
for low $Y_{\rm e}$ objects with young and old ages, respectively (Section
3.3). This diagnostic method distinguishes NSMs from SNe. In the next section,
we focus on the sensitivities in the $\gamma$-ray band that are required for
current and future MeV $\gamma$-ray missions that aim to search for Galactic
NSM remnants.

###  4.1 Detectable distance to Galactic Neutron-Star Merger Remnants

A $\gamma$-ray flux from the brightest line in a particular age range can be
used to estimate the distance that a virtual $\gamma$-ray instrument with a
certain line sensitivity will be able to detect. If we assume the age of the
NSM remnants in Figure 9, then we can estimate the limit of the distance that
can be detected with the line sensitivity of a specific instrument. The top of
Figure 11 summarizes the achievable limit of $d$ for NSM remnants as a
function of $t$ and is calculated for the three energy bands 3–75 keV, 75–500
keV, and 500–4,000 keV. For example, instruments with a line sensitivity of
$10^{-7}$ photons s-1 cm-2 in the 3–75 keV band (red line in the top of Figure
11), such as Hitomi HXI (Takahashi et al., 2014) and NuSTAR (Harrison et al.,
2013), can observe the brightest lines from NSM remnants with $t<10^{3}$ years
at $d=$ 3 kpc. We also checked the degradation of the distance limit due to
the Doppler-broadening effect, as shown in the middle of Figure 11, but the
results do not dramatically change. If the line sensitivities are the same
among the three energy bands, then hard X-rays (thick lines) will be a
powerful tool in the search for NSMs that are younger than $t<10^{3}$ years,
but the $\gamma$-ray observations (dotted or dashed lines) are better for
surveying NSMs that are older than $t>10^{4}$ years.

Figure 11: (Top, Middle) The distance limit of NSMs as a function of $t$ are
shown with the line sensitivities given in units of photons s-1 cm-2; see keys
for details. The thick, dotted, and dashed lines represent the results for the
3–75 keV, 75–500 keV, and 500–4,000 keV bands, respectively. The top and
middle panels show the results without and with the Doppler-broadening effect,
respectively. The line photons are accumulated within $\Delta E=3$ keV in the
middle panel. Note that the small jump in the 75–500 keV data at $t\sim
10^{3}$ years in the middle panel is due to the interaction between the two
brightest lines when the Doppler-broadening effect is applied. (Bottom)
Probabilities for the existence of NSMs, which we took as an assumption in the
calculation for Figure 12, in the $t$–$d$ plane.

The G4.8+6.2 associated with AD 1163 is one example, from the middle of Figure
11, that provides the requirement for the $\gamma$-ray sensitivity needed to
observe an NSM remnant with a known distance and age. The object is reported
to be a young kilonova remnant with $t\sim 860$ years (Liu et al., 2019). If
it is an NSM remnant at $d\sim 10$ kpc, then a sensitivity of $10^{-8}$ and
$10^{-9}$ photons s-1 cm-2 is required to observe G4.8+6.2 in the hard X-ray
and $\gamma$-ray bands, respectively. This sensitivity is roughly one or two
(or more) orders of magnitude deeper than that of INTEGRAL IBIS (Winkler et
al., 2003). If the distance is closer at $d\sim 3$ kpc, then hard X-ray
instruments with a sensitivity of $10^{-7}$ photons s-1 cm-2 in the 3–75 keV
band, such as Hitomi HXI (Takahashi et al., 2014) and NuSTAR (Harrison et al.,
2013), are expected to be able to observe the emissions from the object.

###  4.2 Direct estimation of local NSM rates using gamma rays

To estimate the coverage of Galactic NSM remnants observable for specific
$\gamma$-ray sensitivities, we first prepare a probability map for the
existence of Galactic NSMs. This is given in the same plane as the top and
middle of Figure 11 (the $t$–$d$ plane). Since NSMs are not uniformly
distributed in our Galaxy, we apply the probabilities for NSMs in the $d$ and
$t$ spaces given by Wu et al. (2019) and multiply them to get the plot shown
in the bottom of Figure 11. We assume that the NSMs are primarily concentrated
around the Galactic plane. Most of the NSMs are expected to exist at around
$d\sim 8$ kpc and $t\sim 10^{4}$–$10^{6}$ years, as has already been described
in Wu et al. (2019).

We then accumulate the probabilities for the existence of NSMs (bottom of
Figure 11) within the distance-limit curves (top and middle of Figure 11). As
a result, we obtained the coverage of Galactic NSMs as a function of the line
sensitivity; this is shown in the top of Figure 12. For example, if we survey
Galactic NSM remnants with an instrument having a line sensitivity of
$10^{-8}$ photons s-1 cm-2 in the 3–75 keV band, then we expect to observe
about 3% of the NSMs in our Galaxy with $M_{\rm ej}=0.01M_{\odot}$. This value
corresponds to about one object if we assume an NSM rate in our Galaxy of 30
per million years. In addition, we performed the same procedure to estimate
the NSM coverage in units of erg s-1 cm-2, which requires a sensitivity that
is $E_{\rm\gamma}^{1}$ times higher than that required for units of photons
s-1 cm-2. $E_{\rm\gamma}$ here is the photon energy (the energy of the
$\gamma$-ray line). The results are shown in the bottom of Figure 12.
Therefore, instruments that can achieve a sensitivity of $10^{-14}$ erg s-1
cm-2 in the 75–500 keV or the 500–4000 keV bands are expected to be able to
observe one NSM remnant with $M_{\rm ej}=0.01M_{\odot}$ in our Galaxy under
the same assumption of the NSM rate mentioned above. Similarly, a sensitivity
of $10^{-15}$ erg s-1 cm-2 is required in the hard X-ray band to observe one
object with the same $M_{\rm ej}$.

Figure 12: Coverage of NSM remnants in our Galaxy as a function of the narrow-
band sensitivity is shown. In the top panel, the results are given in units of
photons s-1 cm-2. The magenta and red plots assume $M_{\rm ej}=0.01M_{\odot}$
with and without the Doppler effect, respectively. Similarly, the cyan and
blue plots assume $M_{\rm ej}=0.05M_{\odot}$ with and without the Doppler
effect. The thick, dotted, and dashed lines represent the coverage values in
the 3–75 keV, 75–500 keV, and 500–4,000 keV bands, respectively. The top axis
(shown in green) represents the detectable number of NSMs under the assumption
of an NSM rate in our Galaxy of 30 per million years. The bottom panel shows
the sensitivity results given in units of erg s-1 cm-2, but omitting the
results without the Doppler effect.

The NSM rates from previous studies are summarized in Figure 13. The NSM rates
are estimated using several methods, and even though the values approach each
other recently, they still have non-negligible uncertainties or systematic
errors that are dependent on the methods used. According to Figure 12,
instruments with higher sensitivities can cover more than 10% of NSMs and
should be able to observe multiple Galactic NSM remnants (meaning a
sensitivity of $10^{-9.5}$–$10^{-8.5}$ photons s-1 cm-2 or
$10^{-16.5}$–$10^{-14.5}$ erg s-1 cm-2 in the hard X-ray to $\gamma$-ray
bands). The actual numbers observed by future NSM surveys with highly
sensitive instruments will provide direct information for the NSM rate in the
local universe.

Figure 13: The NSM rate in our Galaxy estimated by previous studies as a
function of the published year. The red, magenta, green, and blue plots
represent the estimation from the binary or double pulsar population
(Belczynski et al., 2002; Voss & Tauris, 2003; Kim et al., 2006, 2015;
Chruslinska et al., 2017; Pol et al., 2019; Grunthal et al., 2021), star
evolution (Belczynski et al., 2010; Dominik et al., 2012; Mennekens &
Vanbeveren, 2014; Artale et al., 2019; Olejak et al., 2020; Chu et al., 2022),
short gamma-ray bursts (Petrillo et al., 2013; Jin et al., 2015), and
gravitational wave events (Abbott et al., 2019, 2020, 2021), respectively. The
values from the gravitational wave events in Gpc-3 yr-1 are converted into the
Myr-1 unit under the assumption of a galactic density of 0.01 galaxy Mpc-3.

###  4.3 Sensitivity requirements for future missions

To assess the feasibility of detecting Galactic NSM remnants using past,
current, and future $\gamma$-ray missions, the $\gamma$-ray spectra expected
from NSMs (Figures 4 and 5) are compared with the sensitivities of the
instruments for these missions in Figure 14. For the MeV bands, we expect that
in the 2030s sensitivities will be achieved that are one or two orders of
magnitude higher than those of current missions. Consequently, we conclude
that future missions, such as e-ASTROGAM, AMEGO, and GRAMS, have the potential
to detect MeV emissions from young NSM remnants in the age range of $t=$
100–1000 years old at 10 kpc, with a sensitivity of approximately
$10^{-13}$–$10^{-11}$ erg cm-2 s-1. Furthermore, the hard X-ray band below 100
keV is also useful in searching for NSM remnants. NuSTAR data may be able to
indicate very young NSM remnants at about $t=$ 100 years, and FORCE may be
able to detect emissions from an NSM older than $t<10^{4}$ years at 10 kpc.

Figure 14: The $3\sigma$ sensitivities of the missions in the hard X-ray to
$\gamma$-ray bands are compared with the $\gamma$ rays expected from NSMs that
are 100, $10^{4}$, and $10^{6}$ years old, with the assumption of an ejecta
mass of 0.05 $M_{\cdot}$ at a 10 kpc distance with $\Delta E=3$ keV. Past and
future missions are shown with the dotted and thick lines, respectively. The
sensitivities are taken from the following references: CGRO/COMPTEL (9 years),
INTEGRAL/SPI (1 year), e-ASTROGAM (1 year), and Fermi LAT (10 years in survey
mode) from de Angelis et al. (2018), INTEGRAL/IBIS ($10^{6}$ sec) and AMEGO (5
years) from Kierans (2020), CGRO/EGRET (9 years) from Si-liang et al. (2021),
NuSTAR ($10^{6}$ sec) from Harrison et al. (2013), FORCE (3.5$\sigma$,
$10^{6}$ sec) from Nakazawa et al. (2018), Hitomi HXI and SGD (100 ksec) from
Takahashi et al. (2014), SMILE3 ($10^{6}$ sec, best condition) from Takada et
al. (2021), and GRAMS (1 year) from Aramaki et al. (2020). For reference, 1.0,
$10^{-3}$, $10^{-6}$, and $10^{-9}$ times the hard X-ray flux from the Crab
Nebula (Kouzu et al., 2013), with a simple extension into the $\gamma$-ray
band with the single power-law spectrum, are shown in dashed lines and are
noted as m-crab, $\mu$-crab, and n-crab, respectively.

This work was supported in part by JSPS KAKENHI [Grant Nos. JP18H04571 and
JP20K04009 (YT), 18H01232 and 22H01251(RY), 20K03957 (SF), 20H00174 (SK),
21H01121 (SK, YT, SF), 19K03908 (AB)]. YT and SK are deeply appreciative of
the Observational Astrophysics Institute at Saitama University for supporting
the research fund, and R.Y. deeply appreciates the Aoyama Gakuin University
Research Institute for helping to fund our research. Finally, we would like to
thank the anonymous referee for his/her careful reading of our manuscript and
helpful comments.

## References

  * Abbott et al. (2017) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017, Phys. Rev. Lett., 119, 161101, doi: 10.1103/PhysRevLett.119.161101
  * Abbott et al. (2019) —. 2019, Physical Review X, 9, 031040, doi: 10.1103/PhysRevX.9.031040
  * Abbott et al. (2020) —. 2020, ApJ, 892, L3, doi: 10.3847/2041-8213/ab75f5
  * Abbott et al. (2021) Abbott, R., Abbott, T. D., Abraham, S., et al. 2021, ApJ, 913, L7, doi: 10.3847/2041-8213/abe949
  * Andrews et al. (2020) Andrews, S., Fryer, C., Even, W., Jones, S., & Pignatari, M. 2020, ApJ, 890, 35, doi: 10.3847/1538-4357/ab64f8
  * Aramaki et al. (2020) Aramaki, T., Adrian, P. O. H., Karagiorgi, G., & Odaka, H. 2020, Astroparticle Physics, 114, 107, doi: 10.1016/j.astropartphys.2019.07.002
  * Arnaud (1996) Arnaud, K. A. 1996, in Astronomical Society of the Pacific Conference Series, Vol. 101, Astronomical Data Analysis Software and Systems V, ed. G. H. Jacoby & J. Barnes, 17
  * Arnould et al. (2007) Arnould, M., Goriely, S., & Takahashi, K. 2007, Phys. Rep., 450, 97, doi: 10.1016/j.physrep.2007.06.002
  * Artale et al. (2019) Artale, M. C., Mapelli, M., Giacobbo, N., et al. 2019, MNRAS, 487, 1675, doi: 10.1093/mnras/stz1382
  * Beer et al. (1997) Beer, H., Corvi, F., & Mutti, P. 1997, ApJ, 474, 843, doi: 10.1086/303480
  * Belczynski et al. (2010) Belczynski, K., Dominik, M., Bulik, T., et al. 2010, ApJ, 715, L138, doi: 10.1088/2041-8205/715/2/L138
  * Belczynski et al. (2002) Belczynski, K., Kalogera, V., & Bulik, T. 2002, ApJ, 572, 407, doi: 10.1086/340304
  * Burbidge et al. (1957) Burbidge, E. M., Burbidge, G. R., Fowler, W. A., & Hoyle, F. 1957, Reviews of Modern Physics, 29, 547, doi: 10.1103/RevModPhys.29.547
  * Cameron (1957) Cameron, A. G. W. 1957, PASP, 69, 201, doi: 10.1086/127051
  * Chruslinska et al. (2017) Chruslinska, M., Belczynski, K., Bulik, T., & Gladysz, W. 2017, Acta Astron., 67, 37, doi: 10.32023/0001-5237/67.1.2
  * Chu et al. (2022) Chu, Q., Yu, S., & Lu, Y. 2022, MNRAS, 509, 1557, doi: 10.1093/mnras/stab2882
  * Cowan et al. (1991) Cowan, J. J., Thielemann, F.-K., & Truran, J. W. 1991, Phys. Rep., 208, 267, doi: 10.1016/0370-1573(91)90070-3
  * De Angelis et al. (2017) De Angelis, A., Tatischeff, V., Tavani, M., et al. 2017, Experimental Astronomy, 44, 25, doi: 10.1007/s10686-017-9533-6
  * de Angelis et al. (2018) de Angelis, A., Tatischeff, V., Grenier, I. A., et al. 2018, Journal of High Energy Astrophysics, 19, 1, doi: 10.1016/j.jheap.2018.07.001
  * Dominik et al. (2012) Dominik, M., Belczynski, K., Fryer, C., et al. 2012, ApJ, 759, 52, doi: 10.1088/0004-637X/759/1/52
  * Freiburghaus et al. (1999) Freiburghaus, C., Rosswog, S., & Thielemann, F. K. 1999, ApJ, 525, L121, doi: 10.1086/312343
  * Fujimoto et al. (2007) Fujimoto, S.-i., Hashimoto, M.-a., Kotake, K., & Yamada, S. 2007, ApJ, 656, 382, doi: 10.1086/509908
  * Grefenstette et al. (2017) Grefenstette, B. W., Fryer, C. L., Harrison, F. A., et al. 2017, ApJ, 834, 19, doi: 10.3847/1538-4357/834/1/19
  * Grunthal et al. (2021) Grunthal, K., Kramer, M., & Desvignes, G. 2021, MNRAS, 507, 5658, doi: 10.1093/mnras/stab2198
  * Hanisch et al. (2001) Hanisch, R. J., Farris, A., Greisen, E. W., et al. 2001, A&A, 376, 359, doi: 10.1051/0004-6361:20010923
  * Harrison et al. (2013) Harrison, F. A., Craig, W. W., Christensen, F. E., et al. 2013, ApJ, 770, 103, doi: 10.1088/0004-637X/770/2/103
  * Hempel & Schaffner-Bielich (2010) Hempel, M., & Schaffner-Bielich, J. 2010, Nucl. Phys. A, 837, 210, doi: 10.1016/j.nuclphysa.2010.02.010
  * Hotokezaka et al. (2015) Hotokezaka, K., Piran, T., & Paul, M. 2015, Nature Physics, 11, 1042, doi: 10.1038/nphys3574
  * Hotokezaka et al. (2016) Hotokezaka, K., Wanajo, S., Tanaka, M., et al. 2016, MNRAS, 459, 35, doi: 10.1093/mnras/stw404
  * Jin et al. (2015) Jin, Z.-P., Li, X., Cano, Z., et al. 2015, ApJ, 811, L22, doi: 10.1088/2041-8205/811/2/L22
  * Kasuga et al. (2018) Kasuga, T., Sato, T., Mori, K., Yamaguchi, H., & Bamba, A. 2018, PASJ, 70, 88, doi: 10.1093/pasj/psy085
  * Katakura (2012) Katakura, J. 2012, JAEA-Data/Code, 2011-025
  * Katakura & Minato (2016) Katakura, J., & Minato, F. 2016, JAEA-Data/Code, 2015-030
  * Kierans (2020) Kierans, C. A. 2020, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 11444, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, 1144431, doi: 10.1117/12.2562352
  * Kim et al. (2006) Kim, C., Kalogera, V., & Lorimer, D. R. 2006, arXiv e-prints, astro.  https://arxiv.org/abs/astro-ph/0608280
  * Kim et al. (2015) Kim, C., Perera, B. B. P., & McLaughlin, M. A. 2015, MNRAS, 448, 928, doi: 10.1093/mnras/stu2729
  * Kouzu et al. (2013) Kouzu, T., Tashiro, M. S., Terada, Y., et al. 2013, PASJ, 65, 74, doi: 10.1093/pasj/65.4.74
  * Kullmann et al. (2022) Kullmann, I., Goriely, S., Just, O., et al. 2022, MNRAS, 510, 2804, doi: 10.1093/mnras/stab3393
  * Lattimer & Schramm (1974) Lattimer, J. M., & Schramm, D. N. 1974, ApJ, 192, L145, doi: 10.1086/181612
  * Li (2019) Li, L.-X. 2019, ApJ, 872, 19, doi: 10.3847/1538-4357/aaf961
  * Liu et al. (2019) Liu, Y., Zou, Y.-C., Jiang, B., et al. 2019, MNRAS, 490, L21, doi: 10.1093/mnrasl/slz141
  * McKee & Ostriker (1977) McKee, C. F., & Ostriker, J. P. 1977, ApJ, 218, 148, doi: 10.1086/155667
  * Mennekens & Vanbeveren (2014) Mennekens, N., & Vanbeveren, D. 2014, A&A, 564, A134, doi: 10.1051/0004-6361/201322198
  * Metzger et al. (2010) Metzger, B. D., Martínez-Pinedo, G., Darbha, S., et al. 2010, MNRAS, 406, 2650, doi: 10.1111/j.1365-2966.2010.16864.x
  * Nakazawa et al. (2018) Nakazawa, K., Mori, K., Tsuru, T. G., et al. 2018, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 10699, Space Telescopes and Instrumentation 2018: Ultraviolet to Gamma Ray, ed. J.-W. A. den Herder, S. Nikzad, & K. Nakazawa, 106992D, doi: 10.1117/12.2309344
  * Olejak et al. (2020) Olejak, A., Belczynski, K., Bulik, T., & Sobolewska, M. 2020, A&A, 638, A94, doi: 10.1051/0004-6361/201936557
  * Petrillo et al. (2013) Petrillo, C. E., Dietz, A., & Cavaglià, M. 2013, ApJ, 767, 140, doi: 10.1088/0004-637X/767/2/140
  * Pol et al. (2019) Pol, N., McLaughlin, M., & Lorimer, D. R. 2019, ApJ, 870, 71, doi: 10.3847/1538-4357/aaf006
  * Qian & Wasserburg (2007) Qian, Y. Z., & Wasserburg, G. J. 2007, Phys. Rep., 442, 237, doi: 10.1016/j.physrep.2007.02.006
  * Si-liang et al. (2021) Si-liang, F., Peng, F., Yi-fan, H., Tian-yu, M., & Yan, X. 2021, Chinese Astron. Astrophys., 45, 281, doi: 10.1016/j.chinastron.2021.08.002
  * Steiner et al. (2013) Steiner, A. W., Hempel, M., & Fischer, T. 2013, ApJ, 774, 17, doi: 10.1088/0004-637X/774/1/17
  * Takada et al. (2021) Takada, A., Takemura, T., Yoshikawa, K., et al. 2021, arXiv e-prints, arXiv:2107.00180.  https://arxiv.org/abs/2107.00180
  * Takahashi et al. (2014) Takahashi, T., Mitsuda, K., Kelley, R., et al. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9144, Space Telescopes and Instrumentation 2014: Ultraviolet to Gamma Ray, ed. T. Takahashi, J.-W. A. den Herder, & M. Bautz, 914425, doi: 10.1117/12.2055681
  * Tanvir et al. (2013) Tanvir, N. R., Levan, A. J., Fruchter, A. S., et al. 2013, Nature, 500, 547, doi: 10.1038/nature12505
  * Taylor (1950) Taylor, G. 1950, Proceedings of the Royal Society of London Series A, 201, 159, doi: 10.1098/rspa.1950.0049
  * Typel et al. (2010) Typel, S., Röpke, G., Klähn, T., Blaschke, D., & Wolter, H. H. 2010, Phys. Rev. C, 81, 015803, doi: 10.1103/PhysRevC.81.015803
  * Villar et al. (2017) Villar, V. A., Guillochon, J., Berger, E., et al. 2017, ApJ, 851, L21, doi: 10.3847/2041-8213/aa9c84
  * Voss & Tauris (2003) Voss, R., & Tauris, T. M. 2003, MNRAS, 342, 1169, doi: 10.1046/j.1365-8711.2003.06616.x
  * Wallner et al. (2015) Wallner, A., Faestermann, T., Feige, J., et al. 2015, Nature Communications, 6, 5956, doi: 10.1038/ncomms6956
  * Wanajo & Ishimaru (2006) Wanajo, S., & Ishimaru, Y. 2006, Nucl. Phys. A, 777, 676, doi: 10.1016/j.nuclphysa.2005.10.012
  * Wanajo et al. (2011) Wanajo, S., Janka, H.-T., & Müller, B. 2011, ApJ, 726, L15, doi: 10.1088/2041-8205/726/2/L15
  * Wanajo et al. (2014) Wanajo, S., Sekiguchi, Y., Nishimura, N., et al. 2014, ApJ, 789, L39, doi: 10.1088/2041-8205/789/2/L39
  * Wang et al. (2020) Wang, X., N3AS Collaboration, Vassh, N., et al. 2020, ApJ, 903, L3, doi: 10.3847/2041-8213/abbe18
  * Winkler et al. (2003) Winkler, C., Courvoisier, T. J. L., Di Cocco, G., et al. 2003, A&A, 411, L1, doi: 10.1051/0004-6361:20031288
  * Wu et al. (2019) Wu, M.-R., Banerjee, P., Metzger, B. D., et al. 2019, ApJ, 880, 23, doi: 10.3847/1538-4357/ab2593

##  Appendix A XSPEC Model for gamma rays from r-process objects

The $\gamma$-ray spectra for NSM remnants in this paper are implemented as the
table spectral model for XSPEC (Arnaud, 1996) version 12 in the HEAsoft
package. The model is provided with the file name of rprocgamma.mod in the
flexible image transport system (FITS) format (Hanisch et al., 2001). The
model parameters and descriptions are summarized below. The model does not
contain the Doppler effect, but instead it gives the output of the second step
described in Section 2.1.

  * •

time : time from the merging event, in units of years.

  * •

Ye$NN$ ($NN=10,15,20,25,30,35,40,45$): the ejecta mass of the $r$-process
nuclei in solar mass units $M_{\odot}$ under the $Y_{\rm e}=0.NN$ environment.

  * •

$z$ : redshift

  * •

norm : normalization in units of photons s-1 cm-2 $d^{-2}$ $(1+z)^{-2}$, where
$d$ is the distance to the object in kpc.

For example, the spectral model of a Galactic NSM remnant at $t=100.0$ years
and $d=10$ kpc, calculated under the assumption that the ejecta masses of
$Y_{\rm e}=$0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40, and 0.45 are 0.0454,
0.0485, 0.146, 0.297, 0.103, 0.251, 0.105, and 0.00326 $M_{\odot}$,
respectively, without the Doppler-broadening effect, is described with the
following parameters:

    Model atable{rprocgamma.mod}<1> Source No.: 1   Active/Off
    Model Model Component  Parameter  Unit     Value
     par  comp
       1    1   rprocgamma time       year     100.000      +/-  0.0
       2    1   rprocgamma Ye10       Msun     4.54000E-02  +/-  0.0
       3    1   rprocgamma Ye15       Msun     4.85000E-02  +/-  0.0
       4    1   rprocgamma Ye20       Msun     0.146000     +/-  0.0
       5    1   rprocgamma Ye25       Msun     0.297000     +/-  0.0
       6    1   rprocgamma Ye30       Msun     0.103000     +/-  0.0
       7    1   rprocgamma Ye35       Msun     0.251000     +/-  0.0
       8    1   rprocgamma Ye40       Msun     0.105000     +/-  0.0
       9    1   rprocgamma Ye45       Msun     3.26000E-03  +/-  0.0
      10    1   rprocgamma z                   0.0          frozen
      11    1   rprocgamma norm                0.01000      +/-  0.0

An example of the XSPEC commands to use this numerical model with the Doppler
broadening effect are followings. The comments shown to the right of %% and
are not executed. For detail, please check the manual for XSPEC.

        XSPEC12> model gsmooth(atable{rprocgamma.mod})   %% set the model and parameters
        XSPEC12> cpd /xw              %% change the plot device using X-Windows system
        XSPEC12> dummyrsp 1 4000 1000 %% set a dummy response from 1 to 4,000 keV in 1,000 bins
        XSPEC12> plot model           %% plot the model on the screen