New X-ray bound on density of primordial black holes

TL;DR

This paper establishes new upper limits on primordial black hole abundance using X-ray data, constraining their density across a broad mass range and assessing their potential role in explaining LIGO black hole detections.

Contribution

It provides the first X-ray based upper bounds on primordial black hole density over a wide mass spectrum, refining previous constraints.

Findings

PBH density constrained across mass range from a few solar masses to 2×10^7 solar masses.

X-ray data limits PBH contribution to observed compact objects in galaxies.

PBH density compatible with LIGO black hole origins is marginally allowed.

Abstract

We set a new upper limit on the abundance of primordial black holes (PBH) based on existing X-ray data. PBH interactions with interstellar medium should result in significant fluxes of X-ray photons, which would contribute to the observed number density of compact X-ray objects in galaxies. The data constrain PBH number density in the mass range from a few $M_\odot$ to $2\times 10^7 M_\odot$. PBH density needed to account for the origin of black holes detected by LIGO is marginally allowed.

Equations

$$\dot{M}$$

$$\times$$

$$L$$

$$\times$$

$$\frac{r_{\mathrm{disk}}}{r_s}$$

$$\simeq$$

$$\frac{dN}{d{L}_X}=\frac{dN}{d\dot{M}}\frac{d\dot{M}}{d{L}_X},$$

$$\frac{dN}{d\dot{M}}$$

$$\frac{d{f}_n}{dn}=f_0(n/n_{\min }{)}^{-\beta },(n_{\min }\le n\le n_{\max }).$$

$$\frac{d{f}_v}{dv}=\sqrt{\frac{2}{\pi }}\frac{v^2}{{\sigma }_v^3}\exp \left(-\frac{v^2}{2{\sigma }_v^2}\right).$$

$$F_{\mathrm{df}}=-\frac{4\pi G^2{M}_{\mathrm{BH}}^{2}n\mu m_p}{v^2}\ln \left(\frac{r_{\max }}{r_{\min }}\right),$$

$$\frac{dN}{d\dot{M}}=N_{\mathrm{PBH},\mathrm{disk}}{\int }_{n_{\min }}^{n_{\max }}dn\frac{d{f}_n}{dn}\frac{d{f}_v}{dv}(v_0)\frac{{\widetilde{v^{\prime }}}_{v=v_0}^2{v_0^{\prime }}^2}{3v_0^3\dot{M}},$$

Full text

1 Introduction

Primordial black holes (PBH) could form in the early universe in a variety of plausible scenarios [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17]. They can account for all dark matter (DM) in a narrow mass window around $10^{20}$g, but even if PBH make up only a small fraction of DM, they can play a role in $r$-process nucleosynthesis [18], provide seeds for supermassive black holes [19, 8], or contribute to gravitational wave signals [20, 21, 22, 23, 24].

The mass of PBHs surviving until present range from $10^{15}$g [2, 3] to as high as $10^{5}M_{\odot}$ [8]. To explore this wide range of masses, various observational techniques and theoretical ideas have been employed [10]. These include extragalactic gamma rays from PBH evaporation [6], femtolensing of gamma-ray bursts [25], r-process nucleosynthesis [18], white-dwarf explosions [26], Subaru/HSC microlensing of stars in M31 [27], Kepler microlensing of stars [28], MACHO/EROS/OGLE mcirolensing of stars [29] and quasar microlensing [30], survival of a star cluster in Eridanus H [31], wide binary disruption [32], millilensing of quasars [33], generation of large-scale structure through Poission fluctuations [34], dynamical friction on halo objects [35], non-detection of gravitational wave background [36], and accretion effects on the cosmic microwave background [37, 38].

Recently, gravitational waves from binary black hole mergers have been detected by LIGO [39]. The origin of those binary black hole systems is under debate. PBHs are among the possible progenitors [20, 21, 40, 22, 9, 12]. Some published constraints in the stellar and intermediate mass range have been recently questioned [22], and there is a strong interest in finding additional reliable constraints at this range of masses.

Freely-floating black holes (BHs), such as kicked stellar remnant BHs, are expected to interact with ambient interstellar medium (ISM) gas through Bondi-Hoyle-Lyttleton accretion [e.g. 41, 42, 43, 44, 45, 46, 47]. PBH in the stellar or intermediate mass ranges should also emit significant fluxes of X-ray photons through this process [48, 49]. By comparing electromagnetic wave data with expected radiation signals from floating PBHs, one can set constraints on the abundance of PBHs [48, 49, 50].

In this paper, we set new constraints on PBH abundance by utilizing the observed number density of X-ray binaries (XRBs) including ultra-luminous X-ray sources (ULXs). X-ray observations have already identified many X-ray emitting compact extragalactic objects, mostly XRBs. Moreover, we take into account a thick DM disk, so-called dark disk, which is suggested by numerical simulations [51, 52]. In Sec. 2, we introduce accretion processes on to floating PBHs. Constraints on PBH abundances are presented in Sec. 3. Discussion and Conclusions are given in Sec. 4.

2 Accretion onto Floating PBHs

The accretion rate onto a floating PBH from the ISM is given by the Bondi-Hoyle-Lyttleton rate [53, 54, 55]:

[TABLE]

where $M_{\rm BH}$ is the PBH mass, $n$ is the ISM gas density, $m_{p}$ is the proton mass, $\mu$ is the mean molecular weight $\mu=2.72$, $L_{\rm Edd}$ is the Eddington luminosity, the Bondi radius $r_{B}=GM_{\rm BH}/\tilde{v}^{2}$, and $\tilde{v}\equiv(v^{2}+c_{s}^{2})^{1/2}$. $v$ is the PBH velocity relative to the ISM gas cloud and $c_{s}$ is the sound speed which depends on the gas temperature. The effective value of the speed of sound for ISM gas is $c_{s}\sim 10~{}{\rm km~{}s^{-1}}$ because ISM gas phases have a turbulent velocity of $\sim 10~{}{\rm km~{}s^{-1}}$ in approximate pressure balance with each other. The mass accretion processes on to floating BHs is extensively investigated by [44]. We basically follow their arguments hereafter.

The luminosity of a floating PBH interacting with ISM gas is estimated as

[TABLE]

where $\epsilon$ is the disk radiative efficiency. In the radiatively-inefficient accretion flow (RIAF) regime $L\lesssim 0.01L_{\rm Edd}$, the luminosity is known to scale with $L\propto\dot{M}^{2}$ rather than $L\propto\dot{M}$ as in the standard accretion disk regime [56]. Therefore, we set the radiative efficiency as $\epsilon(\dot{M})=\epsilon_{0}/[1.0+({\dot{M}}/{0.01\dot{M}_{\rm Edd}})^{-1}]$, where $\dot{M}_{\rm Edd}\equiv L_{\rm Edd}/\epsilon_{0}c^{2}$ is the Eddington mass accretion rate at the efficiency of $\epsilon_{0}=0.1$. We note that $\epsilon_{0}$ can be from 0.057 for a Schwarzschild BH to 0.42 for an extreme Kerr BH [e.g. 56].

Does gas form an accretion disk? Small perturbations in the density or velocity of the accreting gas lead to an angular momentum large enough to form a disk [57]. The accreted angular momentum is $l=(1/4)({\Delta\rho}/{\rho})\tilde{v}r_{B}$, where $\Delta\rho$ is the difference in density between the top and bottom of the accretion cylinder. The ISM density has turbulence with a Kolomogorov spectrum of fluctuations $\delta\rho/\rho\sim(L/10^{18}~{}{\rm cm})^{1/3}$ extending down to $\sim 10^{8}~{}{\rm cm}$ [58]. We can find the radius of the accretion disk $r_{\rm disk}$ by setting $\Delta\rho/\rho=\delta\rho/\rho|_{L=2r_{B}}$ and equating the angular momentum of the gas with the Keplerian angular momentum $l_{K}=\sqrt{GM_{\rm BH}r_{\rm disk}}$. This gives

[TABLE]

where $r_{s}$ is the Schwarzschild radius $2GM_{\rm BH}/c^{2}$. Thus, an accretion disk can be formed around a PBH by the Bondi-Hoyle-Lyttleton accretion.

We estimate the accretion timescale as $t_{\rm acc}\approx L_{\rm MC}/v\sim 2.4\times 10^{6}(n/100~{}{\rm cm^{-3}})(v/10~{}{\rm km~{}s^{-1}})^{-1}~{}{\rm yr}$ which is longer than the observation periods of various X-ray sources since the dawn of X-ray astronomy. $L_{\rm MC}\sim 24.7(n/100~{}{\rm cm}^{-3})^{-0.9}~{}{\rm pc}$ is the size of a typical molecular cloud [59].

The number of PBHs interacting with the ISM gas can be expressed in terms of the luminosity function:

[TABLE]

where $dN/\dot{M}$ is the number of PBHs accreting at $\dot{M}$. Function $d\dot{M}/dL_{X}$ can be derived from Eq.(2), and $dN/\dot{M}$ is given by [44, 45]

[TABLE]

where the number of PBHs passing through the Galactic disk region is given by $N_{\rm PBH,disk}$ and $v_{0}^{2}=(4\pi G^{2}M^{2}n\mu m_{p}/\dot{M})^{2/3}-c_{s}^{2}$. $df_{n}/dn$ and $df_{v}/dv$ is the probability distribution of the number density of the ISM gas and the velocity distribution of PBHs, respectively. We can not take all values of $n$ between $n_{\rm min}$ and $n_{\rm max}$ since $v_{0}^{2}>0$ and $n$ must be greater than $c_{s}^{3}/(4\pi G^{2}M^{2}\mu m_{p}/\dot{M})$.

Using the Navarro-Frenk-White (NFW) profile [60] with parameters given in [61] for the Milky way, the DM mass included within 15 kpc radius is about $\sim$7% of $M_{\rm halo}=10^{12}M_{\odot}$. $N_{\rm PBH,disk}$ can be estimated as $0.07(\pi r^{2}2H_{g})/(4\pi r^{3}/3)N_{\rm BH,tot}\sim 0.11H_{g}/(15~{}{\rm kpc})N_{\rm BH,tot}$, where $N_{\rm PBH,tot}$ is the total number of PBHs in our galaxy $(\Omega_{\rm PBH}/{\Omega_{M}})(M_{\rm halo}/M_{\rm BH})$. $H_{g}$ is the disk scale height of the ISM gas. Recently, star formation activity in the extreme outer Galaxy at 22 kpc away from the Galactic center has been reported [62]. Since the radial distribution of such clouds is not well understood yet, we do not take into account such distant clouds. If we extrapolate the ISM gas distribution out to 22 kpc, we would expect a 10 % increase in the X-ray emitting PBH number density.

The probability distribution of the number density of the ISM gas is defined as [44, 45]

[TABLE]

The normalization factor is defined as $f_{0}\equiv{(\beta-2)\langle\Sigma\rangle}/{2\mu m_{p}n_{\rm min}^{2}H_{g}}$, where $\langle\Sigma\rangle$ is the mean surface mass density of the ISM gas in the Galactic disk, which is averaged for PBHs passing through the disk by $\langle\Sigma\rangle=\int 2\pi r\Sigma(r)\rho_{\rm DM}(r)dr/\int 2\pi r\rho_{\rm DM}(r)dr$, where $r$ is the Galactocentric radius and $\rho_{\rm DM}(r)$ is the dark halo density given by the NFW profile.

For the ISM gas, we consider three gas phases, galactic disk molecular clouds, cold HI clouds, and the central molecular zone (CMZ). We do not take into account warm or hot ISM gas, which would host only low luminosity objects below X-ray detection thresholds. We take $\Sigma(r)$ of [63] and [64], for the molecular clouds and the HI clouds, respectively. We find $\langle\Sigma\rangle=27M_{\odot}~{}{\rm pc}^{-2}$ for the molecular clouds and $2.8M_{\odot}~{}{\rm pc}^{-2}$ for the cold HI clouds. The mass of CMZ is $\sim 3\times 10^{7}M_{\odot}$ and its radius is $\sim 150~{}{\rm pc}$ with a scale height of $H_{g}=30~{}{\rm pc}$ [see 65, and references therein]. In this paper, we assume a constant radial distribution of the surface mass density leading $\langle\Sigma\rangle=420M_{\odot}~{}{\rm pc}^{-2}$ for the CMZ region. Although the mass function of molecular cloud cores has been recently studied [65], the density probability distribution function is not well constrained. Since the CMZ is composed of molecular clouds, we assume the same number density distribution function shape as that of molecular clouds. For the scale height of the molecular clouds and the HI clouds, we take $H_{g}=75$ pc and 150 pc, respectively, following [44]. For the probability distribution function of the ISM gas, we take $\beta=2.8$, $n_{\rm min}=10^{2}~{}{\rm cm^{-3}}$, and $n_{\rm max}=10^{5}~{}{\rm cm^{-3}}$ for the molecular gas and the CMZ, while $\beta=3.8$, $n_{\rm min}=10~{}{\rm cm^{-3}}$, and $n_{\rm max}=10^{2}~{}{\rm cm^{-3}}$ for the HI atomic gas [44].

The velocity distribution of PBHs can be described by a Maxwell-Boltzmann distribution,

[TABLE]

We take the velocity dispersion $\sigma_{v}=150~{}{\rm km~{}s^{-1}}$ [66], because PBHs are expected to make some fraction of DM particles.

Dynamical friction from gas reduces the speed of PBHs in molecular clouds [67, 45]. The dynamical friction force is given by [67]

[TABLE]

where $r_{\rm max}$ and $r_{\rm min}$ are the sizes of the surrounding medium and the object receiving the force, respectively. We take $r_{\max}=L_{\rm MC}$ and $r_{\rm min}=r_{s}$. Then, by taking into account dynamical friction during the crossing time, one can rewrite Eq.(2.8) as

[TABLE]

where $\tilde{v^{\prime}}^{2}=({v^{\prime}_{0}}^{2}+c_{2}^{2})^{1/2}$, ${v^{\prime}_{0}}^{2}=(4\pi G^{2}M^{2}n\mu m_{p}/\dot{M})^{2/3}-c_{s}^{2}$, and $v_{0}^{4}={v^{\prime}_{0}}^{4}+\\ 16\pi G^{2}M_{\rm BH}n\mu m_{p}L_{\rm MC}\ln(L_{\rm MC}/r_{s})$. Hereafter, we take into account the effect of the dynamical friction.

Numerical simulations suggest that a thick DM disk, so-called dark disk, forms in a galactic halo after a merger at $z<2$ due to existence of stellar/gas disk [51, 52]. The dark disk corotates with the stellar disk with a scale height of $\sim 1.5~{}{\rm kpc}$. The dark disk can contribute to $f_{\rm DM,disk}\sim 0.25$-1.5 times of the non-rotating halo DM density at the solar position. The lower bound is the conservative limit since it takes into account only one merger with a Large Magellanic Cloud mass. The velocity dispersion of the dark disk is $50~{}{\rm km~{}s^{-1}}$ [51, 52, 68, 69]. This low velocity would boost up interaction rates. We assume a radially constant $f_{\rm DM,disk}$ in the disk, which is reasonable assumption for our order-of-magnitude calculation [52].

The luminosity functions of XRBs have been established in literature and they are also known to be correlated with star formation rate (SFR) of their host galaxies [e.g. 70, 71]. We adopt the SFR normalized luminosity function in [71], where the SFR is determined for the stellar mass range of $0.1-100~{}M_{\odot}$ with a Salpeter initial mass function. We take the SFR of the Milky way as $\sim 1~{}M_{\odot}~{}{\rm yr}^{-1}$ [e.g., 72].

3 Results

The results are shown in Fig. 1 for $M_{\rm BH}=10^{3}M_{\odot}$ assuming $\Omega_{\rm PBH}=\Omega_{\rm DM}$. It is clear that $\Omega_{\rm PBH}$ need to be less than $\Omega_{\rm DM}$ even for $f_{\rm DM,disk}=0.0$, otherwise it will violate observations of XRBs. By taking into account the dark disk component, the luminosity function is enhanced because of its low velocity dispersion. Interaction with molecular clouds dominate the luminosity function. Other ISM components will not contribute significantly.

The PBH density is constrained by data points having X-ray luminosities above $10^{39}~{}{\rm erg~{}s^{-1}}$ corresponding to ULXs. Those luminous objects would prevent efficient mass accretion as given in Equation 2.2 due to radiation feedback effects. Mass accretion rate is known to decrease by a factor of $(1-L/L_{\rm Edd})^{2}$ in the Bondi-Hoyle-Lyttleton accretion [73]. The radiation feedback effect is also shown in Figure. 1. Although it will decrease the XLF at high luminosity end, the constraint will not be significantly different.

Constraints on $f_{\rm PBH}\equiv\Omega_{\rm PBH}/\Omega_{\rm DM}$ can be set based on the requirement that the predicted luminosity function does not violate the observed luminosity function of XRBs at any luminosities. In Fig. 2, we show the upper bound on the PBH mass fraction to DM based on the XRB luminosity function together with other constraints. For simplicity, we assume a monochromatic mass distribution for PBHs. We show three limits. One is our fiducial model with $f_{\rm DM,disk}=0.25$. The others are the cases with $f_{\rm DM,disk}=0.0$ and $1.0$.

If we do not take into account the effect of dynamical friction, the constraints can become $\sim 3$% tighter at $10~{}M_{\odot}\lesssim M_{\rm PBH}\lesssim 1000~{}M_{\odot}$. This is because dynamical friction changes the velocity by about 1 ${\rm km~{}s}^{-1}$ in our parameter space.

For the ISM gas phases, we assumed $n_{\rm min}=10^{2}~{}{\rm cm^{-3}}$ for the molecular gas and the CMZ and $n_{\rm min}=10~{}{\rm cm^{-3}}$ for the HI atomic gas [44]. Given in Eq. 2, luminosities of accreting massive PBHs become too high for observed XRB luminosities (See. Fig. 1). Thus, we see an upper cutoff at high PBH mass range in our constraints.

4 Discussion and Conclusions

Recently, the constraint by [37] utilizing the cosmic microwave background (CMB) data is revisited by [38] with the latest Planck data by taking account into detailed physical processes. For example, the radiative efficiency is estimated from the first principle in [38], while it is fixed in [37]. As a consequence, the CMB constraints on the PBH density is significantly weakened (Fig. 2). However, we note that [74] obtained a tighter constraint than that by [37] with Planck data following the method in [37].

Our methods constrain the PBH abundance in the mass range from a few $M_{\odot}$ to $2\times 10^{7}M_{\odot}$. The obtained constraint is tighter than the other constraints at $10~{}M_{\odot}\lesssim M_{\rm PBH}\lesssim 1000~{}M_{\odot}$. Since we consider the ISM gas phases of molecular clouds, HI clouds, and CMZ whose density is at $10~{}{\rm cm^{-3}}\leq n\leq 10^{5}~{}{\rm cm^{-3}}$, we can not put tight constraints at low and high mass end of PBHs with X-ray data ranging $35\leq\log L_{X}\leq 41$.

Very recently, [50] has also studied the PBH abundances based on the study of Bondi accretion processes in Galactic PBHs. Although our estimates are based on different reasoning, the results are consistent for the case of $f_{\rm DM,disk}=0.0$. The approach of Ref. [50] is based on comparison with the X-ray and radio source catalogs in the Galactic ridge region, while we use the luminosity function of x-ray binaries. While in Ref. [50] it is assumed that all the sources listed in the X-ray catalogs are BH candidates, it is known that white dwarfs and coronally active stars are the dominant sources in the Galactic ridge region [75]. To make use of the radio catalog data (which we do not use) the authors of Ref. [50] converted x-ray luminosities to radio luminosities using the fundamental plane of black holes [76]. However, the existence of the plane is questioned by recent studies [77].

At $30M_{\odot}$, we find $f_{\rm PBH}\lesssim 4\times 10^{-2}$ with our fiducial model $f_{\rm DM,disk}=0.25$. This bound rules out PBH scenarios requiring $f_{\rm PBH}$ of an order of unity to explain the GW events [22], while our limit still marginally allows the models in which binary PBH systems are efficiently formed in the radiation dominant era [78] requiring $f_{\rm PBH}\sim 2\times 10^{-3}-2\times 10^{-2}$ at $30M_{\odot}$ [40, 79, 80]. Further detailed studies of $f_{\rm DM,disk}$ and the binary black hole merger rate will tighten the current limit on to the PBH scenario as the dominant GW sources.

X-ray emission from free-floating BHs interacting with ISM gas is similar to XRBs [41, 42, 43, 44, 45, 46, 47], but some features may help distinguish such isolated systems from XRB systems. First, x-ray emitting isolated BHs must be associated with molecular clouds or cold HI clouds if they are to emit significant fluxes of x-ray photons (Eq. 2). Second, one will not see any orbital motion, because the BHs are isolated. Lastly, these isolated BHs are expected to be associated with ionized bubbles, because the ISM gas needed for the emission is ionized by the emitted x-rays. Thus, a spatially extended iron K$\alpha$ fluorescent line structure can be found, although one needs to resolve the source of the size $L_{\rm MC}$.

To summarize, we have set a new observational bound on the abundance of PBHs in the mass range from a few $M_{\odot}$ to $2\times 10^{7}M_{\odot}$ utilizing X-ray data by considering PBH interactions with interstellar medium. The density of ULXs gives the most stringent constraint on the PBH abundance. Our constraints are tighter than recent constraints by Planck data at $5M_{\odot}\lesssim M_{\rm PBH}\lesssim 300M_{\odot}$.

Acknowledgments

The authors thank Yutaka Fujita, Kunihito Ioka, Norita Kawanaka, Tomonori Totani, Masato Tsuboi, Masahiro Tsujimoto for comments and discussions. YI is supported by the JAXA international top young fellowship and JSPS KAKENHI Grant Number JP16K13813. The work of A.K. was supported by the U.S. Department of Energy Grant No. DE-SC0009937, and by the World Premier International Research Center Initiative (WPI), MEXT, Japan.