Primordial Black Holes and $r$-Process Nucleosynthesis

TL;DR

This paper proposes that interactions between primordial black holes and neutron stars can produce r-process elements and observable electromagnetic transients, offering a novel explanation for certain astrophysical phenomena and the origin of heavy elements.

Contribution

It introduces a new mechanism where primordial black holes induce neutron star destruction, leading to r-process nucleosynthesis and unique electromagnetic signals, distinct from traditional merger models.

Findings

Primordial black hole interactions can produce r-process elements in neutron stars.

Ejected neutron-rich matter can generate kilonova-like transients and radio bursts.

The scenario aligns with pulsar statistics, dark matter distribution, and observed 511-keV emission.

Abstract

We show that some or all of the inventory of $r$-process nucleosynthesis can be produced in interactions of primordial black holes (PBHs) with neutron stars (NSs) if PBHs with masses ${10}^{-14}\,{\rm M}_\odot < {\rm M}_{\rm PBH} < {10}^{-8}\,{\rm M}_\odot$ make up a few percent or more of the dark matter. A PBH captured by a neutron star (NS) sinks to the center of the NS and consumes it from the inside. When this occurs in a rotating millisecond-period NS, the resulting spin-up ejects $\sim 0.1-0.5\,{\rm M}_{\odot}$ of relatively cold neutron-rich material. This ejection process and the accompanying decompression and decay of nuclear matter can produce electromagnetic transients, such as a kilonova-type afterglow and fast radio bursts. These transients are not accompanied by significant gravitational radiation or neutrinos, allowing such events to be differentiated from compact object mergers occurring within the distance sensitivity limits of gravitational wave observatories. The PBH-NS destruction scenario is consistent with pulsar and NS statistics, the dark matter content and spatial distributions in the Galaxy and Ultra Faint Dwarfs (UFD), as well as with the $r$-process content and evolution histories in these sites. Ejected matter is heated by beta decay, which leads to emission of positrons in an amount consistent with the observed 511-keV line from the Galactic Center.

Equations

$$F_{0}=\sqrt{6\pi}\dfrac{\rho_{\text{DM}}}{m_{\rm PBH}}\Big{[}\dfrac{R_{\rm NS}R_{\rm S}}{\overline{v}(1-R_{\rm S}/R_{\rm NS})}\Big{]}\Big{(}1-e^{-3E_{\rm loss}/(m_{\rm PBH}\overline{v}^{2})}\Big{)}~{},$$

$$\left\{\begin{array}[]{ll}F_{0}^{\rm MW}=1.5\times 10^{-11}{\rm/yr}\\ F_{0}^{\rm UFD}=6.0\times 10^{-10}{\rm/yr}~{}.\\ \end{array}\right.$$

$$\langle M_{\rm ej}\rangle=\dfrac{\int_{P_{\min}}^{\infty}\Big{(}\dfrac{dN_{\rm MSP}}{dP}\Big{)}M_{\rm ej}(P)dP}{\int_{P_{\min}}^{\infty}\Big{(}\dfrac{dN_{\rm MSP}}{dP}\Big{)}dP}~{},$$

$$R_{e^{+}}=N_{e^{+}}/{\tau }_d\sim 10^{50}{\mathrm{yr}}^{-1}.$$

$$F_0^{\mathrm{MW}}=\int d^3{v}_n{f}_n^{3D}({\vec{v}}_n,{\overline{v}}_n)\int d^3{v}_d{f}_d^{3D}({\vec{v}}_d,{\overline{v}}_d)F_0(|{\vec{v}}_n-{\vec{v}}_d|),$$

$$f^{1D}(v,\overline{v})=\Big{(}\dfrac{3}{2\pi\overline{v}^{2}}\Big{)}^{3/2}e^{-3v^{2}/(2\overline{v}^{2})}$$

$$M_{\rm B}=\dfrac{c_{s}^{2}}{G}\Big{(}\dfrac{2\sqrt{3}\nu^{2}}{\omega_{0}}\Big{)}^{1/3}=1.8\times 10^{-3}~{}(P_{1}^{1/3}T_{5}^{-4/3})~{}M_{\odot}~{},$$

$$M_{\mathrm{crit}}=$$

$$=$$

$$t_{\rm viol}=\Big{(}\dfrac{M_{\rm crit}}{M_{B}}\Big{)}^{2}\tau=400~{}\Big{(}\dfrac{P}{1\,\text{ms}}\Big{)}^{-7}\Big{(}\dfrac{T}{10^{7}~{}\text{K}}\Big{)}^{4}~{}\text{yr}~{},$$

$$\frac{m_{\mathrm{crit}}^{8/3}}{\log m_{\mathrm{crit}}}<$$

$$\displaystyle\times\Big{(}\dfrac{1\,\text{ms}}{P}\Big{)}\Big{(}\dfrac{B_{0}}{\text{G}}\Big{)}^{2}~{},$$

$${\rho }_c=M_{\mathrm{NS}}/(4\pi R_{\mathrm{NS}}^3)\times |{\xi }_1/{\theta }^{\prime }({\xi }_1)|$$

$$m_{\mathrm{BH}}(r_b)=4\pi {\alpha }^3{\rho }_c{\xi }_b^2|{\theta }^{\prime }({\xi }_b)|.$$

$$R_p(r_b)=R_{\mathrm{NS}}+R_{\mathrm{S}}(r_b)-r_b,$$

$$J_0=\frac{2}{5}{M}_{\mathrm{NS}}{R}_{\mathrm{NS}}^2{\Omega }_0=$$

$$+J_{\mathrm{BH}}+J_{\mathrm{tran}}.$$

$$J_{tran}$$

$$\displaystyle=J_{0}-\dfrac{2}{5}(M_{\rm NS}-m_{\rm BH})\Big{(}\dfrac{R_{p}^{5}-R_{\rm S}^{5}}{R_{p}^{3}-R_{\rm S}^{3}}\Big{)}\Omega~{},$$

$$v_{\rm eq}^{mod}=\dfrac{R_{\rm eq}M_{\rm NS}R_{\rm NS}^{2}}{M_{\rm NS}-m_{\rm BH}}\Big{(}\dfrac{R_{p}^{3}-R_{\rm S}^{3}}{R_{p}^{5}-R_{\rm S}^{5}}\Big{)}\Omega_{0}~{}.$$

$$R_{\mathrm{p},\mathrm{ej}}=$$

$$R_{\mathrm{eq},\mathrm{ej}}=$$

$$M_{\mathrm{ej}}=$$

$$v_{\mathrm{esc},\mathrm{ej}}=$$

$$v_{\mathrm{eq},\mathrm{ej}}=$$

Full text

††preprint: IPMU17-0048

Primordial black holes (PBHs) can account for all or part of the dark matter (DM) Zel’dovich and Novikov (1967); Hawking (1971); Carr and Hawking (1974); Garcia-Bellido et al. (1996); Khlopov (2010); Frampton et al. (2010); Kawasaki et al. (2012, 2016); Cotner and Kusenko (2016); Carr et al. (2016); Inomata et al. (2016, 2017); Georg and Watson (2017). If a PBH is captured by a neutron star (NS), it settles into the center and grows until the supply of nuclear matter is exhausted by accretion and ejection.

In this Letter we show that NS disruptions by PBHs in DM-rich environments, such as the Galactic center (GC) and dwarf spheroidal galaxies, provide a viable site for r-process nucleosynthesis, thus offering a solution to a long-standing puzzle Burbidge et al. (1957); Arnould et al. (2007); Cameron (1957); Truran et al. (2002); Qian and Wasserburg (2007). The transients accompanying NS disruption events and the positrons produced in these events are consistent with present observations, and they offer a way of testing the NS–PBH scenario in the future.

We will demonstrate that, when a PBH accretes matter inside a rapidly rotating millisecond pulsar (MSP), the resulting pulsar spin-up causes $\sim 0.1-0.5\,M_{\odot}$ of neutron-rich material to be ejected without significant heating and only modest neutrino emission. This provides a favorable setting for $r$-process nucleosynthesis, occuring on the Galactic time-scales, which can evade several problems that have challenged the leading proposed $r$-process production sites, such as neutrino-heated winds from core collapse supernovae or binary compact object mergers (COM) Argast et al. (2004); Freiburghaus et al. (1999). The unusual distribution of $r$-process abundances within the ultra-faint dwarf spheroidal galaxies Ji et al. (2016); Beniamini et al. (2016) (UFDs) is naturally explained by the rates of PBH capture in these systems. The rates are also consistent with the paucity of pulsars in the GC Dexter and O’Leary (2014). A similar distribution of $r$-process material in UFDs can be expected from NS disruptions due to black holes produced in the NS interiors by accretion of particle dark matter onto the NS Bramante and Linden (2016), although the rates and the implications for dark matter properties are, of course, different. The probability of PBH capture depends on both the PBH and the NS densities. DM-rich environments, such as the GC and dwarf spheroidal galaxies, are not known to host NSs. An exception is the young magnetar Mori et al. (2013); Kennea et al. (2013), whose age is small compared to the time scales of PBH capture. NSs are found in the disk and the halo, as well as in the globular clusters, where the dark matter density Bradford et al. (2011); Ibata et al. (2013) is too low to cause a substantial decrease in the pulsar population. Positrons emitted from the heated neutron-rich ejecta can account for the observed 511-keV line from the GC Prantzos et al. (2011); Lingenfelter et al. (2009). The final stages of neutron star demise can be the origin Fuller and Ott (2015); Shand et al. (2016) of some of the recently observed Lorimer et al. (2007) fast radio bursts (FRBs), as well as X-ray and $\gamma$-ray transients. A kilonova-type Kasen et al. (2015); Li and Paczynski (1998); Metzger et al. (2010); Roberts et al. (2011); Hotokezaka and Piran (2015); Piran et al. (2013); Barnes and Kasen (2013); Martin et al. (2015) afterglow can accompany the decompressing nuclear matter ejecta, but unlike COM, these events are not associated with a significant release of neutrinos or gravitational radiation. Therefore, future observations of gravitational waves and kilonovae will be able to distinguish between $r$-process scenarios.

Millisecond pulsars are responsible for the predominant contribution to the nucleosynthesis initiated by PBH-induced centrifugal ejection of neutron-rich material, since MSPs have the highest angular velocities at the time of PBH capture. The most prominent sites of $r$-process production must have a high density of MSPs as well as PBHs. The latter trace the DM spatial distribution. The DM density is high in the GC and in the UFDs. On the other hand, the MSP density is high in molecular clouds, including the Central Molecular Zone (CMZ), and the globular clusters. While the CMZ is located within the GC with an extremely high DM density, observations imply that the DM content of globular clusters is fairly low Bradford et al. (2011); Ibata et al. (2013). The product of the DM density and the MSP density is still sufficient to allow for some $r$-process nucleosynthesis in both the UFDs and the globular clusters, but we estimate that CMZ accounts for 10% to 50% of the total Galactic production. We include contributions from the GC (CMZ) and the rest of the halo (which may be comparable, within uncertainties) Read et al. (2009, 2008).

The CMZ has an approximate size of $\sim 200$ pc and is located near the GC, where supernova rates are the highest. Since the DM density peaks at the GC, the CMZ is a site of frequent PBH-MSP interactions. The pulsar formation rate O’Leary et al. (2016) in the CMZ is $1.5\times 10^{-3}$yr*-1* ($7\%$ of the Galactic formation rate), consistent with the GeV $\gamma$-ray flux observed Ajello et al. (2016) from the GC by the Fermi Large Area Telescope. Hence, we expect $N_{p}^{\text{GC}}\approx 1.5\times 10^{7}$ NSs to be produced at the GC during the lifetime of the Galaxy, $t_{G}\sim 10^{10}$ yr. Roughly, 30%-50% of these NSs become MSPs Lorimer (2013), and the number of MSPs with a particular rotation period can then be estimated from a population model Cordes and Chernoff (1997); Faucher-Giguere and Kaspi (2006); Lorimer (2013).

Simulations and observations of UFDs imply Ji et al. (2016) that $\sim 2000$ core-collapse supernovae have occurred in 10 UFDs during $t_{\text{UFD}}\sim 5\times 10^{8}$ yr. Hence, we expect that $N_{p}^{\text{UFD}}\sim 10^{2}$ pulsars have been produced in each of these systems, on average, and we estimate the fraction of fast-rotating MSPs using a population model Cordes and Chernoff (1997).

The black hole capture rate can be calculated using the initial PBH mass $(m_{\rm PBH})$, DM density and velocity dispersion (assuming a Maxwellian distribution). With the base Milky Way (MW) and UFD capture rates denoted as $F_{0}^{\rm MW}$ and $F_{0}^{\rm UFD}$, one obtains Capela et al. (2013) the full PBH-NS capture rates $F=(\Omega_{\rm PBH}/\Omega_{\rm DM})F_{0}^{\rm MW}$ and $F=(\Omega_{\rm PBH}/\Omega_{\rm DM})F_{0}^{\rm UFD}$ in the MW and UFD, respectively. The number of PBHs captured by a NS is $F\cdot t$, where the time $t$ is either $t_{G}$ or $t_{\text{UFD}}$ for the MW and the UFD, respectively. For our analysis we consider a typical NS with mass $M_{NS}=1.5M_{\odot}$ and radius $R_{NS}=10$ km. The base NS-PBH capture rate $F_{0}$ is given by  Capela et al. (2013)

[TABLE]

where $\overline{v}$ is the DM velocity dispersion, $R_{\rm S}=2GM_{\rm NS}$ and $R_{\rm NS}$ are the BH Schwarzschild radius and the NS radius, respectively. $E_{\rm loss}$ is the energy loss associated with PBH-NS interaction. BH capture can occur when $E_{\rm loss}>m_{\rm PBH}v_{0}^{2}/2$, with $v_{0}$ being the asymptotic velocity of the PBH. Taking a uniform flux of PBHs across the star, the average energy loss for a typical NS is found to be $E_{\rm loss}\simeq 58.8G^{2}m_{\rm PBH}^{2}M_{\rm NS}/R_{\rm NS}^{2}$. Since MSPs originate from binaries, a higher binary gravitational potential causes an increase in the capture rate. We, therefore, assume that the capture rate for MSPs is a factor 2 higher than for isolated NSs. For example, for typical values of parameters and $m_{\rm PBH}=10^{19}$ g, one obtains

[TABLE]

For the MW we have used a velocity dispersion of 48 km/s and 105 km/s for NS and DM, respectively, as well as DM density $8.8\times 10^{2}$ GeV/cm3. The pulsar and DM velocity dispersions are simultaneously taken into account for the MW as described below. For UFD we have used DM velocity dispersion 2.5 km/s and DM density 10 GeV/cm3.

A PBH could also be captured by a NS progenitor prior to supernova core collapse Capela et al. (2014), but this does not increase the capture rates significantly.

Natal pulsar kicks can enable pulsars to escape from the region of interest. We include this effect in our calculations (see Supplemental Material sup ).

Pulsar lifetimes in the presence of PBHs with a given number density can be estimated as $\langle t_{\rm NS}\rangle=1/F+t_{\rm loss}+t_{\rm con}$, where the first term describes the mean BH capture time, $t_{\rm loss}$ is the time for the PBH to be brought within the NS once it is gravitationally captured, and $t_{\rm con}$ is the time for the black hole to consume the NS. For a typical NS one finds Capela et al. (2013) that $t_{\rm loss}\simeq 4.1\times 10^{4}\,(m_{\rm PBH}/10^{22}\,\text{g})^{-3/2}\,{\rm yr}$. The spherical accretion rate of NS matter onto the PBH is described by the Bondi equation $dm_{\rm BH}/dt=4\pi\lambda_{s}G^{2}m_{\rm BH}^{2}\rho_{c}/v_{s}^{3}=C_{0}m_{\rm BH}^{2}$, where $m_{\rm BH}(t)$ is the time-dependent mass of the central black hole, $v_{s}$ is the sound speed, $\rho_{c}$ is the central density and $\lambda_{s}$ is a density profile parameter. For typical NS values of Shapiro and Teukolsky (1983) $v_{s}=0.17,\rho_{c}=10^{15}$ g/cm3 and $\lambda_{s}=0.707$ (for a star described by an $n=3$ polytrope) we obtain that $t_{\rm con}=10~{}(10^{19}~{}\text{g}/m_{\rm PBH})$ yr. If PBHs make up all of the DM, we calculate that $\langle t_{\rm NS}\rangle<10^{12}$ yr for $10^{17}~{}\text{g}<m_{\rm PBH}<10^{25}\,\text{g}$, implying that a $\mathcal{O}(1-10)\%$ fraction of pulsars should have been consumed in the age of the Galaxy. This is consistent with observations Dexter and O’Leary (2014) suggesting an under abundance of MSPs near the central Galactic black hole, Sgr A*∗*. A recently discovered young, $4\times 10^{4}$ yr old, magnetar J1745-2900 located just 0.1 pc from the GC Mori et al. (2013); Kennea et al. (2013) is also consistent with our results, since this magnetar’s age is shorter than $\langle t_{\rm NS}\rangle$. The unusual surface temperature Zelati et al. (2015) and X-ray luminosity of J1745-2900 warrants scrutiny, as this activity might be consistent with PBH destruction in progress.

Angular momentum transfer determines the dynamics of a NS spin-up. As the captured PBH starts to grow and consume the spinning pulsar from the inside, the radius of the neutron star decreases and angular momentum conservation forces a spin-up. As the star contracts, the fractional change in radius is greater for accreted matter in the inner regions than it is for material further out. This could lead to differential rotation. However, if angular momentum can be efficiently transferred outward, the star can maintain rigid-body rotation. Viscosity Kouvaris and Tinyakov (2014); Markovic (1995) and magnetic stresses Markovic (1995) can prevent differential rotation from developing. It can be shown (see Supplemental Material  sup ) that angular momentum is transferred efficiently on the relevant time scales and that Bondi accretion proceeds nearly uninterrupted throughout the BH evolution.

Ejected mass originates from the star spin-up when matter at the equator exceeds the escape velocity. Using polytropic NS density runs with different indices Shapiro and Teukolsky (1983), we have calculated analytically (see end of “Ejected Mass” section in Supplemental Material sup ) the amount of ejected material. The results are shown in Figure 1 for NS victims with a range of initial rotation periods. Based on our estimates discussed in Supplemental Materials, NS with periods of a few ms can eject more than $10^{-1}M_{\odot}$ of material. A detailed calculation taking into account general relativistic effects Andersson and Kokkotas (2001); Owen et al. (1998) is needed to improve understanding of the ejected mass.

The number of MSPs in the disk with periods greater than $P$ is described by a power-law distributed population model Cordes and Chernoff (1997). We assume that the distribution in the CMZ is the same, and we normalize the total to the number of neutron stars produced in supernova explosions. According to the population model Cordes and Chernoff (1997), N_{\rm MSP}\simeq 1.6\times 10^{4}\Big{(}1.56\,{\rm ms}/P\Big{)}. Using the differential distribution $d(N_{\rm MSP})/dP$, we obtain the population-averaged ejected mass:

[TABLE]

where $P_{\min}$ is the minimal MSP period in the population, and $M_{\rm ej}(P)$ is the ejected mass function interpolated from the distribution shown in Figure 1. We find that the population-averaged ejected mass is $\langle M_{\rm ej}\rangle=0.18M_{\odot}$ and $0.1M_{\odot}$ if we take the shortest period to be $P_{\min}=0.7$ ms (theoretically predicted) and $P_{\min}=1.56$ ms (observed), respectively. Since realistic nuclear matter equations of state suggest flatter NS density profiles than our polytropic approximations, our estimate is conservative and $\langle M_{\rm ej}\rangle$ can be up to a factor of few larger. Alternative population models Lorimer (2013), such as those based on Faucher-Giguere and Kaspi (2006), do not significantly alter the results.

Nucleosynthesis takes place in the ejecta. Heating accompanying the growth of a BH inside a NS results in a temperature increase near the event horizon that is only a factor of few higher than the NS surface temperature Kouvaris and Tinyakov (2014): $T_{\rm h}/T_{\rm surf}\sim 3$. Consequently, neutrino emission is negligible and ejected material does not suffer significant heating or exposure to neutrinos.

Decompression of the centrifugally-ejected, relatively low entropy and very low electron fraction nuclear matter in this scenario could be expected to result in a significant mass fraction of this material participating in $r$-process nucleosynthesis Lattimer et al. (1977); Meyer (1989); Lippuner and Roberts (2015); Jaikumar et al. (2007); Goriely et al. (2011a, 2016); Eichler et al. (2015). The large neutron excess in this scenario, relatively unmolested by neutrino charged current capture-induced reprocessing of the neutron-to-proton ratio, could lead to fission cycling Eichler et al. (2015); Mumpower et al. (2016), thereby tying together the nuclear mass number $A=130$ and $A=195$ $r$-process abundance peaks. Unlike COM $r$-process ejecta, which will have a wide range of neutrino exposures, entropy, and electron fraction and thereby can reproduce the solar system $r$-process abundance pattern Wanajo et al. (2014), the PBH scenario may be challenged in producing the low mass, $A<100$, $r$-process material.

The material ejected in the PBH-NS destruction process is heated by beta decay and fission, resulting in thermodynamic conditions and abundances closely akin to those in the COM-induced “tidal tail” nuclear matter decompression that gives rise to kilonova-like electromagnetic signatures Kasen et al. (2015); Li and Paczynski (1998); Metzger et al. (2010); Roberts et al. (2011); Hotokezaka and Piran (2015); Piran et al. (2013); Barnes and Kasen (2013); Martin et al. (2015). This could be a more luminous and longer duration transient compared to the classic COM-generated kilonovae, as the ejecta in the PBH scenario can have more mass than the tidal tails of COM.

The total amount of ejected $r$-process material in the PBH-NS destruction process can be estimated via $M_{\rm tot}^{r}=FtN_{\rm MSP}\langle M_{\rm ej}\rangle$, assuming that the bulk of the ejecta undergoes $r$-process nucleosynthesis. The overall mass of $r$-process material in the Galaxy is $M_{\rm tot}^{r,{\rm MW}}\sim{10}^{4}\,{\rm M}_{\odot}$. The required fraction of dark matter in the form of PBHs is $(\Omega_{\rm PBH}/\Omega_{DM})=M_{\rm tot}^{r,{\rm MW}}/(F_{0}^{{\rm MW}}t_{G}N_{\rm MSP}^{{\rm GC}}\langle M_{\rm ej}\rangle)$. If the mass of ejected r-process material in a single event is $0.1-0.5\,{\rm M}_{\odot}$, the PBH capture rate ${10}^{-5}-{10}^{-6}\,{\rm Mpc}^{-3}\,{\rm yr}^{-1}$ can account for all of the $r$-process in the Galaxy. At this rate, $10^{5}$ NS disruption events have occurred in the lifetime of the Galaxy.

This rate of NS disruptions in UFDs is also consistent with the observationally inferred UFD $r$-process content and with the uneven distribution of this material among the observed UFD. Observations imply that one in ten of UFDs have been a host to r-process nucleosynthesis events, which must, therefore, be rare Ji et al. (2016); Hirai et al. (2015); Beniamini et al. (2016). The rate $F_{0}^{\rm UFD}$ implies that the probability of a NS disruption in a single UFD is about 0.1, which explains the uneven distribution. The amount of $r$-process material supplied by a single event, $\sim 0.1\,M_{\odot}$, is more than sufficient to explain the observations Ji et al. (2016); Hirai et al. (2015); Beniamini et al. (2016). Only a small fraction $\sim(\bar{v}^{\rm UFD}/v_{\infty})\sim 10^{-3}$ of the produced r-process material is likely to remain in the shallow gravitational potential well of a UFD because it is produced with a velocity $v_{\infty}\sim 0.1\,v_{\rm esc}$. The observations are consistent with this: the required $10^{-4}M_{\odot}$ of r-process material is consistent with a 0.1% fraction of the $0.1M_{\odot}$ produced in a single event.

We have separately fit to the $r$-process abundances for the MW and UFD, accounting for uncertainties in various quantities as described below. The combined requirements result in the allowed region of parameter space shown on Figure 2, along with the current constraints for PBH contribution to the DM abundance. The region denoted “all $r$-process” shows parameter space for which $r$-process observations are fully explained simultaneously in the MW and in UFDs. For our fit we have varied the input parameters over a broad range, covering significant parameter space (see Supplemental Material  sup ). The enclosed region above the line can be interpreted as a constraint of $r$-process material over-production from PBH-NS interactions, subject to large uncertainties in astrophysical input parameters. Energy losses and capture rates for black holes with masses below $\sim 10^{18}$ g are not well understood, and there is an uncertainty in the range of parameters for small masses.

We note that COM-produced $r$-process, with an event rate of ${10}^{-4}-{10}^{-5}~{}{\rm Mpc}^{-3}\,{\rm yr}^{-1}$, could also be consistent with this analysis Goriely et al. (2011b). However, COM simulations suggest an ejecta mass of $\sim 0.01\,{\rm M}_{\odot}$. This would imply a COM rate near the upper end of the allowed range, if COM are to explain all of the $r$-process. Such a rate is still marginally consistent with the current Advanced LIGO (aLIGO) limits, but readily verifiable or refutable when aLIGO reaches its design sensitivity Cote et al. (2017) in a few years. Sensitivity similar to aLIGO is expected in the upcoming Advanced Virgo Acernese et al. (2015) (aVirgo) and KAGRA Aso et al. (2013) experiments. A recent analysis of kilonova Jin et al. (2016) also exhibits tension with observations and highlights the need for an extremely efficient ejection of $r$-process material in COM scenario.

Positron emission from ejecta can explain the observed 511 keV emission line from the Galactic central region Jean et al. (2003), which is consistent with the $e^{+}e^{-}$ annihilation line via positronium formation. The origin of the positrons remains unknown Prantzos et al. (2011). The 511 keV line flux in the bulge component is Knodlseder et al. (2003) $\sim 10^{-3}$ photons cm*-2* s*-1*. The line can be explained through electron–positron annihilations that occur at a rate of $\Gamma(e^{+}e^{-}\rightarrow\gamma\gamma)\sim 10^{50}\,{\rm yr}^{-1}$. Ejected cold nuclear matter expands on a dynamical time scale of $\tau_{e}\sim\alpha/\sqrt{G_{N}\rho}=446\alpha(\rho/{\rm g\,cm}^{-3})^{-1/2}{\rm s}~{},$ where $\alpha=0.01-10$ is a model-dependent parameter Lattimer et al. (1977); Meyer (1989); Goriely et al. (2011b). At the same time, beta decays and fission raise the temperature to $T\sim 0.1$ MeV Lattimer et al. (1977); Meyer (1989); Goriely et al. (2011b). This temperature is high enough to generate a sizable equilibrium density of positrons, which leak through the surface of each clump. Taking the radius of each clump as $R\sim 0.1$ km and the density as $\rho\sim 10^{8}$g/cm3, the total surface area of $0.1\,M_{\odot}$ of ejected material is $A\sim 4\pi R^{2}(0.1M_{\odot}/\rho)/(4\pi R^{3}/3)\sim 10^{20}{\rm cm}^{2}$. The number of positrons emitted in a single event, during the time $\tau_{e}$, while the temperature $T\sim 0.1$ MeV is maintained, can be estimated as $N_{e^{+}}\sim Av\tau_{e}\times 2(m_{e}T/2\pi)^{3/2}\exp(-m_{e}/T)$, where $v$ is the average speed of positrons emitted with a relativistic $\gamma$ factor $\gamma\sim(3T/m_{e})$. If the neutron star disruption events occur in the GC on the time scale of $\tau_{d}\sim 10^{5}$ yr, the average rate of positron production is

[TABLE]

Since the average positron energy $E_{e^{+}}\approx 3T$ is below 3 MeV, the positrons do not annihilate in flight in the interstellar medium Beacom and Yuksel (2006).

Fast radio bursts, kilonovae and other signatures are expected from the PBH capture-induced NS demise. During the final stages of the event, described by dynamical time scales of the order of a few to tens of milliseconds, $10^{41}-10^{43}$ ergs of energy stored in the magnetic field are released. Inside the cold NS, at temperatures below 0.4 MeV, the nuclear matter is a Type II superconductor and magnetic field is concentrated in flux tubes. A consequence of the rapid rearrangement of nuclear matter accompanying ejection is a prodigious release of electromagnetic radiation from magnetic field reconnection and decay. The resulting bursts of radio waves Fuller and Ott (2015); Shand et al. (2016) with duration of a few milliseconds can account for some of the observed Lorimer et al. (2007) FRBs. One FRB is known to be a repeater, while the others appear to be one-time events. The FRB energy of $10^{41}$ erg is consistent with observations Dolag et al. (2015). If 1–10% of the magnetic field energy is converted to radio waves, an FRB could accompany an NS destruction event. The rapidly evolving magnetic field can also accelerate charged particles leading to X-ray and $\gamma$-ray emission.

Detection of an “orphaned” kilonova (macronova) within the aLIGO, aVirgo and KAGRA sensitivity distance ($\sim 200\,{\rm Mpc}$) that is not accompanied by a binary compact object in-spiral gravitational wave signal or a short $\gamma$-ray burst, but possibly associated with an FRB, would constitute an indirect argument that NS disruptions via PBH capture occur and could account for a significant fraction of the $r$-process. Sophisticated numerical simulations of PBH-induced NS collapse and of the accompanying nucleosynthesis and electromagnetic emission (including FRB) could help enable feasible observational search strategies. The search can be further assisted by detailed mapping of chemical abundances that will be made possible by the future Hitomi-2 detector. The stakes are high, as finding evidence for PBH-NS destruction could have profound implications for our understanding of the origin of the heavy elements and for the source and composition of dark matter.

Acknowledgments. We thank B. Carr, Y. Inoue, S. Nagataki, R. Rothschild, and E. Wright for helpful discussions. The work of G.M.F. was supported in part by National Science Foundation Grants No. PHY-1307372 and PHY-1614864. The work of A.K. and V.T. was supported by the U.S. Department of Energy Grant No. DE-SC0009937. A.K. was also supported by the World Premier International Research Center Initiative (WPI), MEXT, Japan.

I Natal Pulsar Kicks

Only $\sim 2\%$ of NS from UFDs remain in the central region on the relevant time scales Bramante and Linden (2016). On the other hand, since some MSPs reside in binaries, their kick velocities are retarded by gravitational binding to the companion and, consequently, the number of resident neutron stars can be higher, $\sim 10\%$. We modify the number $N_{p}^{\rm UFD}$ accordingly.

To account for the pulsar velocity dispersion in the MW, we replace the star capture rate $F_{0}^{\rm MW}$ for stationary NS with an effective rate, combining the Maxwellian velocity distributions of DM and pulsars with a pulsar velocity dispersion Cordes and Chernoff (1997) of $\sim 48$ km s*-1*. Including the velocity dispersions of PBHs and pulsars, the modified rate is

[TABLE]

where the subscripts $n,d$ refer to NS and DM/PBH, respectively, while $f^{3D}(\vec{v},\overline{v})$ denotes the appropriate 3-D Maxwellian velocity ($\vec{v}$) distribution with a mean of $\overline{v}$. The 3-D integrals are decomposed by switching to spherical coordinates, where $d^{3}v_{n}=v_{n}^{2}dv_{n}d(\cos\chi)d\nu$, $d^{3}v_{d}=v_{d}^{2}dv_{d}d(\cos\eta)d\beta$, with the $d(\cos\chi),d(\cos\eta)$ integrals evaluated on the interval $[-1,1]$ and $d\beta,d\nu$ on $[0,2\pi]$. The vector difference of velocities decomposes as $|\vec{v}_{n}-\vec{v}_{d}|=(v_{n}^{2}+v_{d}^{2}-2v_{n}v_{d}\cos\chi)^{1/2}$, where $\chi$ is the lab frame angle between the two velocity vectors. The velocity integrand is then weighted by the 1-D Maxwellian distribution

[TABLE]

and is evaluated over the interval $[0,\infty)$.

II Viscosity and Magnetic Differential Rotation Breaking

Viscosity efficiently breaks differential rotation Kouvaris and Tinyakov (2014); Markovic (1995) and spherical Bondi accretion can be maintained up to BH mass

[TABLE]

where $\nu=2\times 10^{11}T_{5}^{-2}~{}\text{cm}^{2}/\text{s}$ is the kinematic sheer viscosity of the neutron superfluid inside the NS, $\omega_{0}=2\pi/P$ is the angular velocity, $P_{1}=(P/1\,\text{s})$ is the period and $T_{5}=(T/10^{5}~{}\text{K})$ is the temperature. On the other hand, the angular momentum of infalling matter can’t stall Bondi accretion above a critical BH mass

[TABLE]

where $\psi$ quantifies the BH spin and is 1 for a Schwarzchild BH and 1/3 for an extreme Kerr BH Markovic (1995).

Bondi accretion can be violated for central BH masses in the regime $M_{B}<m_{\rm BH}<M_{\rm crit}$. This could slow down the growth of the central BH substantially. The maximal time $t_{\rm viol}$ that a BH can spend violating Bondi accretion is Kouvaris and Tinyakov (2014)

[TABLE]

where $\tau=1/(C_{0}M_{B})$ and where we have used the values relevant for MSPs: $P=1$ ms and $T=10^{7}$ K, where the temperature reflects Bondi accretion heating. Since $t_{\rm viol}$ is short on the relevant time scales ($t_{\rm G}$, $t_{\rm UFD}$), Bondi accretion continues effectively uninterrupted.

The Bondi-violating mass window $M_{B}<m_{\rm BH}<M_{\rm crit}$ may be closed in the presence of magnetic fields, since magnetic torques can facilitate efficient angular momentum transfer. The initial seed field $B_{0}$ can be amplified via a dynamo effect. For the resulting magnetic torque per unit area to exceed the angular momentum current density and ensure rigid rotation, the initial seed field should satisfy the condition Markovic (1995)

[TABLE]

where we take $m_{\rm crit}=M_{\rm crit}/M_{B}$ and $x\sim 1$ is a parameter describing the fraction of the NS mass that has been consumed. Here, $C_{0}$ is a numerical pre-factor of the Bondi accretion as defined in the main text. The resulting condition gives $B_{0}>10^{11}$ G, which is consistent with magnetic fields associated with pulsars.

III Black hole inside a rotating star, ejected mass

As the star spins-up, the matter at the equator can exceed the escape velocity and be ejected. For the purpose of making analytic estimates, we take the star’s density profile to follow a polytropic relation Shapiro and Teukolsky (1983) with index $n=1.5-3$. This range for $n$ subsumes various NS models, from less to more centrally condensed. Some stiffer NS equations-of-state suggest a nearly flat density profiles and, if that is the case, our estimates will be under-estimates of the amount of ejected material. A polytropic relation connects pressure $P$ and density $\rho$ as $P=K\rho^{(1+1/n)}$, where $K$ is a constant. After a change of variables to $r=\alpha\xi$, where $\alpha$ is a constant, the density is given by $\rho(r)=\rho_{c}\theta^{n}(\xi)$, where $\theta(\xi)$ is a solution to the Lane-Emden equation and $\rho_{c}$ is the central density. The first zero at $\xi_{1}$ corresponds to the star’s radius. For a given mass $M_{\rm NS}$ and radius $R_{\rm NS}$, one obtains $\alpha=R_{\rm NS}/\xi_{1}$,

[TABLE]

with the derivative being with respect to $\xi$. As the BH inside the NS consumes the star up to radius $r_{b}=\alpha\xi_{b}$, the BH mass becomes

[TABLE]

The new NS polar radius $R_{p}$ is

[TABLE]

where $R_{S}(r_{b})=2Gm_{\rm BH}(r_{b})$ is the BH Schwarzchild radius.

For a star rotating near the mass-shedding limit, assuming a Roche lobe model description for the extended matter envelope (now a rotationally-squashed spheroid), the equatorial radius is related to the polar radius as Shapiro and Teukolsky (1983); Shapiro and Shibata (2002) $R_{\rm eq}(r_{b})=(3/2)R_{p}(r_{b})$. The equatorial escape velocity is $v_{\rm esc}(R_{\rm eq})=\sqrt{2GM_{\rm NS}/R_{\rm eq}}$ and is $v_{\rm esc}\sim 0.7$ at the onset of BH growth. Initially, before the BH size is appreciable, the equatorial velocity is $v_{\rm eq}^{0}=\Omega_{0}R_{\rm eq}^{0}=(3/2)R_{\rm NS}\Omega^{0}$ for a pulsar rotating at the initial angular velocity $\Omega_{0}$. As the BH grows, the conservation of angular momentum $J=mvr$ for an increment of mass $m$ at the equator determines the equatorial velocity at a later time as $v_{\rm eq}(r_{b})=v_{\rm eq}^{0}R_{\rm eq}^{0}/R_{\rm eq}(r_{b})$.

Conservation of angular momentum demands that the total initial star angular momentum $J_{0}$ is distributed between the black hole and the remainder of the neutron star above the Schwarzschild radius:

[TABLE]

Here, the first term on the right describes the angular momentum of the star’s solid spherical shell $(R_{p}-R_{\rm S})$ in terms of $J_{\rm sh}=I_{\rm sh}\Omega$, where $I_{\rm sh}$ is the shell’s moment of inertia. The last two terms correspond to the angular momentum acquired by the BH itself and the angular momentum transferred out to the remainder of the star. $\Omega(t)$ is the instantaneous angular velocity. Assuming that the infalling angular momentum is transferred outward efficiently, and that the BH spin contribution is negligible Kouvaris and Tinyakov (2014), the outer spherical shell acquires an additional angular momentum contribution

[TABLE]

where $\Omega_{\rm sh}$ is the additional angular velocity acquired by the outer shell. For $\Omega(r_{b})=v_{\rm eq}(r_{b})/R_{\rm eq}(r_{b})$, the spin-up yields a new equatorial speed $v_{\rm eq}^{mod}=R_{\rm eq}(\Omega+\Omega_{\rm sh})$. As long as the star is a rigid rotator,

[TABLE]

When the equatorial velocity exceeds the escape velocity, mass is ejected. For total ejected mass $M_{\rm ej}$, the radius of the rigidly rotating NS is $R_{\rm ej}<R_{\rm NS}$. Hence, the previous equations are modified:

[TABLE]

The modified velocity at the new equator after mass ejection, including the central angular momentum transfer, is $v_{\rm eq,ej}^{\rm mod}=v_{\rm eq}^{\rm mod}(R_{\rm eq,ej}/R_{\rm eq})$.

The condition $v_{\rm esc,ej}(r_{b},R_{\rm ej})\leq v_{\rm eq,ej}^{\rm mod}(r_{b},R_{\rm ej})$ signals mass ejection. The maximal possible amount of ejected material corresponds to $R_{\rm ej}=r_{b}$. For a pulsar of period $P=2\pi/\Omega_{0}$, we solve for $r_{b}$ that satisfies the above condition with $R_{\rm ej}=r_{b}$. Having determined $R_{\rm ej}$, we then numerically integrate the polytropic equation from $R_{\rm ej}$ to $R_{\rm NS}$ to determine the maximal ejected mass $M_{\rm ej}$ from Eq. (19).

While we have neglected BH spin, it is possible for the spin of a slowly-spinning BH to grow to a critical value and thus alter accretion and momentum transfer. However, as found in Ref. Kouvaris and Tinyakov (2014), this occurs in the regime $M_{B}<m_{BH}<M_{\rm crit}$ where spherical Bondi accretion is violated. Since such regime is short and occurs during the late stages of BH growth, this consideration will not significantly affect our conclusions.

IV Allowed Parameter Ranges

Assuming conservatively that the overall amount of DM should not be larger than the baryonic content, which has a mass of $10^{8}M_{\odot}$ in the inner 0.1 kpc volume around the GC, we take the volume-averaged DM density in the CMZ to be $8.8\times 10^{2}$ GeV/cm3. For velocity dispersion we take Kaplinghat et al. (2014) 105 km/s. In a UFD, at 1 pc from the center, the density is Navarro et al. (1997) $\rho=1.5\times 10^{2}$ GeV/cm3 and velocity dispersion is Koposov et al. (2015) $\overline{v}=2.5$ km/s.

To map out the parameter space consistent with the $r$-process abundances of the MW and UFD (Figure 2), we vary the fit input parameters as follows. For the MW, we vary the DM density in the range $50~{}\text{GeV}/\text{cm}^{3}<\rho_{\text{DM}}^{\text{MW}}<8.8\times 10^{2}$ GeV/cm3. The lower bound corresponds to the “flat-core” Burkert profile Burkert (1996) with uniform density in the central kpc. The upper bound is the volume-averaged maximum allowed mass of DM within 0.1 kpc of the GC. For UFD we vary the DM density in the range $0.3~{}\text{GeV}/\text{cm}^{3}<\rho_{\text{DM}}^{\text{UFD}}<15$ GeV/cm3, which corresponds to the Navarro-Frenk-White Navarro et al. (1997) (NFW) profile evaluated in the $1-50$ pc range from the galactic center, respectively. We consider MW DM velocity dispersion values in the range $50~{}\text{km/s}<\overline{v}^{\text{MW}}<200~{}\text{km/s}$, where the lower limit corresponds to possible DM disk within the halo  Read et al. (2008, 2009) and the upper limit corresponds to NFW DM density profile without adiabatic contraction Kaplinghat et al. (2014) at 0.1 kpc from the GC. Additionally, we take values of the pulsar velocity dispersion in the MW to be in the range of Cordes and Chernoff (1997) 48 km/s to Lyne et al. (1998) 80 km/s. The UFD DM velocity dispersion is varied in the range Koposov et al. (2015) $2.7~{}\text{km/s}<\overline{v}^{\text{UFD}}<10~{}\text{km/s}$. Similarly, we also vary the fraction of MSPs not kicked out of the inner $\sim{\rm pc}$ region of the UFD Bramante and Linden (2016) from 2% to 25%. The population-averaged ejected MSP mass is varied in the range $0.09~{}M_{\odot}<\langle M_{\rm ej}\rangle<0.32~{}M_{\odot}$, corresponding to NS density profile that is more centrally condensed for $n=3$, or less so for $n=1.5$ polytropic index, respectively. This provides a conservative estimate of ejected mass, as realistic nuclear equations of state suggest even flatter NS density profiles. In the estimates of ejected mass, we have also included the 1-$\sigma$ uncertainty inherent in the pulsar period-population distribution  Cordes and Chernoff (1997). We find that (see main text), within these ranges of parameters, the UFD and MW $r$-process abundances are explained simultaneously if PBHs contribute at least a few per cent to the overall DM density. In Figure 3 we display the parameter space span for MW and UFD where PBHs can account $r$-process element production.