Primordial black hole constraints for extended mass functions

TL;DR

This paper reevaluates constraints on primordial black holes as dark matter candidates, considering extended mass functions and showing that constraints are generally more stringent than for monochromatic cases, limiting their contribution to dark matter.

Contribution

It introduces a comprehensive analysis of extended mass functions for PBHs and compares constraints with monochromatic models, highlighting the impact on dark matter viability.

Findings

Three mass windows where PBHs could make up all dark matter.

Including all bounds reduces PBHs' contribution to about 10%.

Extended mass functions lead to more stringent constraints than monochromatic ones.

Abstract

We revisit the cosmological and astrophysical constraints on the fraction of the dark matter in primordial black holes (PBHs) with an extended mass function. We consider a variety of mass functions, all of which are described by three parameters: a characteristic mass and width and a dark matter fraction. Various observations then impose constraints on the dark matter fraction as a function of the first two parameters. We show how these constraints relate to those for a monochromatic mass function, demonstrating that they usually become more stringent in the extended case than the monochromatic one. Considering only the well-established bounds, and neglecting the ones that depend on additional astrophysical assumptions, we find that there are three mass windows, around $4\times 10^{-17}M_\odot,$ $2\times 10^{-14}M_\odot$ and $25-100M_\odot$, where PBHs can constitute all dark matter. However, if one includes all the bounds, PBHs can only constitute of order $10\%$ of the dark matter.

Figures

Equations

$$\psi (M)\propto M\frac{\mathrm{d}n}{\mathrm{d}M},$$

$$f_{\mathrm{PBH}}\equiv \frac{{\Omega }_{\mathrm{PBH}}}{{\Omega }_{\mathrm{DM}}}=\int \mathrm{d}M\psi (M),$$

$$\psi (M)=\frac{f_{\mathrm{PBH}}}{\sqrt{2\pi }\sigma M}\exp \left(-\frac{{\log }^2(M/M_c)}{2{\sigma }^2}\right),$$

$$\psi (M)\propto M^{\gamma -1}(M_{\min }<M<M_{\max }).$$

$$\psi (M)\propto M^{2.85}\exp (-(M/M_f{)}^{2.85}),$$

$$\log M_c\equiv ⟨\log M{⟩}_{\psi },{\sigma }^2\equiv ⟨{\log }^{2}M{⟩}_{\psi }-⟨\log M{⟩}_{\psi }^2,$$

$$M_c=M_{\mathrm{cut}}{e}^{-\frac{1}{\gamma }},\sigma =\frac{1}{|\gamma |},$$

$$A[\psi (M)]=A_0+\int \mathrm{d}M\psi (M)K_1(M)$$

$$+\int \mathrm{d}{M}_1\mathrm{d}{M}_2\psi (M_1)\psi (M_2)K_2(M_1,M_2)+\dots ,$$

$$A[\psi (M)]\le A_{\exp },$$

$${\psi }_{\mathrm{mon}}(M)\equiv f_{\mathrm{PBH}}(M_c)\delta (M-M_c),$$

$$f_{\mathrm{PBH}}(M_c)\le \frac{A_{\exp }-A_0}{K_1(M_c)}\equiv f_{\max }(M_c).$$

$$\int \mathrm{d}M\frac{\psi (M)}{f_{\max }(M)}\le 1.$$

$$\sum _{j=1}^N{\left(\int \mathrm{d}M\frac{\psi (M)}{f_{\max ,j}(M)}\right)}^2\le 1,$$

$$A[\psi (M)]\approx A_0+\int \mathrm{d}M{K}_1(M){\psi }_{\mathrm{eff}}(M).$$

$${\psi }_{\mathrm{eff}}(M)=\sum _n{\alpha }_n{\psi }_n(M),$$

$${\psi }_n(M)\equiv \int \prod _{i=1}^n\mathrm{d}{M}_i\psi (M_i)\delta (M-\Sigma M_i),$$

$${\chi }^2(f_{\mathrm{PBH}})=\sum _{j=1}^N\frac{(A_j(f_{\mathrm{PBH}})-{\mu }_j{)}^2}{{\sigma }_j^2}.$$

$$n^2\ge {\chi }^2-{\chi }_{\min }^2=\sum _{j=1}^N{\left(\int \mathrm{d}M\psi (M)\frac{K_{1,j}(M)}{{\sigma }_j}\right)}^2,$$

$$f_{\mathrm{PBH}}(M)\le \frac{n{\sigma }_j}{K_{1,j}(M)}\equiv f_{\max ,j}(M).$$

$$\sum _{j=1}^N{\left(\int \mathrm{d}M\frac{\psi (M)}{f_{\max ,j}(M)}\right)}^2\le 1.$$

$$f_{\max }(M)={\left(\sum _{j=1}^N{f}_{\max ,j}(M{)}^{-2}\right)}^{-1/2}.$$

$$\sum _{j=1}^N{\left(\int \mathrm{d}M\frac{\psi (M)}{f_{\max ,j}(M)}\right)}^2\le {\left(\int \mathrm{d}M\frac{\psi (M)}{f_{\max }(M)}\right)}^2.$$

Full text

I Introduction

Besides its gravitational interaction, little is known about the nature of dark matter (DM) except that it is dynamically “cold”. Although the cold dark matter (CDM) is usually assumed to be some form of elementary particle Jungman et al. (1996); Bertone et al. (2005), there is still no evidence for this and PBHs which are too large to have evaporated by now are a possible alternative Chapline and Frampton (2016). Because they form when the baryons only comprise a small fraction of the total cosmological density Hawking (1971); Carr and Hawking (1974); Carr (1975), they are exempt from the Big Bang nucleosynthesis limits on the baryonic density Chapline (1975). Being much more massive than elementary particles, they could also have a greater variety of observational consequences. Indeed the PBH scenario is already severely constrained by cosmological and astrophysical observations Carr et al. (2010, 2016a).

The recent detection of gravitational waves from merging black holes with mass ${\mathcal{O}}(10)M_{\odot}$ by LIGO Abbott et al. (2016a, b) has revived interest in the possibility of PBH DM Kashlinsky (2016); Bird et al. (2016); Clesse and García-Bellido (2016); Sasaki et al. (2016). Although the PBH coalescence rate depends on very uncertain astrophysical parameters, explaining the observed event rate would require the PBHs to contain at least a substantial fraction of the DM.

This has led to a reassessment of the existing PBH bounds in two directions. First, it has been argued that the existing constraints on PBHs with monochromatic mass functions can be relaxed by invoking extended mass functions Clesse and García-Bellido (2015); Inomata et al. (2017), the latter arising naturally if the PBHs are created from inflationary fluctuations  Carr et al. (1994); Garcia-Bellido et al. (1996); Kawasaki et al. (1998); Yokoyama (1998a); Kohri et al. (2008); Frampton et al. (2010); Drees and Erfani (2011); Clesse and García-Bellido (2015); Garcia-Bellido and Ruiz Morales (2017); Orlofsky et al. (2017); Inomata et al. (2017); Kannike et al. (2017) or some form of cosmological phase transition Crawford and Schramm (1982); Hawking et al. (1982); Khlopov (2010); Belotsky et al. (2014). However, there is still no rigorous treatment of how to apply the PBH bounds for an extended mass function and different analyses have led to different conclusions. For example, Ref. Carr et al. (2016a) concludes that intermediate mass PBHs could provide the DM, whereas Refs. Green (2016) and Kuhnel and Freese (2017) reach the opposite conclusion. Second, there have been revisions to the constraints themselves, with some previous bounds being weakened (e.g. those associated with accretion Ricotti et al. (2008); Chen et al. (2016); Ali-Haïmoud and Kamionkowski (2016); Aloni et al. (2016)) and some new bounds being added. Indeed, the constraints are being constantly revised and the recent review of Ref. Carr et al. (2016a) already needs to be updated.

Currently there is no comprehensive study which combines these two approaches. In this paper, we fill this gap by presenting a general method for analysing the latest PBH constraints over the broad mass range $10^{-18}-10^{4}M_{\odot}$ and applying them to an extended PBH mass function.

It should be stressed that PBHs could play an important cosmological role even if they have much less than the DM density. For example, they could be useful in explaining the rapid structure formation at small cosmological scales, provide seeds for supermassive black holes or galaxies and explain other unsolved astrophysical and cosmological puzzles (García-Bellido, 2017; Carr and Silk, 2017). This underlines the importance of knowing how the PBH density is distributed between different masses.

II Constraints on monochromatic PBH mass function

The main constraints on a PBH population derive from PBH evaporations, various gravitational lensing experiments, neutron star capture, numerous dynamical effects, and PBH accretion. The form of these constraints for a monochromatic PBH mass function is indicated in the upper left panel of Fig. 1, together with relevant references. It must be stressed that these constraints depend on various cosmological and astrophysical assumptions, as well as unknown black hole physics. We therefore list these uncertainties explicitly.

The constraints on PBH evaporation via Hawking radiation Hawking (1974) depend on the observed extragalactic photon flux intensity, $I\propto E^{-1-\epsilon},$ where $E$ is the photon energy and $\epsilon$ parametrizes the spectral tilt Carr et al. (2016a). There is some uncertainty in this parameter, so we present our results for the two extreme cases, $\epsilon=0.4$ (solid purple line) Strong et al. (2004) and $\epsilon=0.1$ (dotted purple line) Sreekumar et al. (1998).

The Cosmic Microwave Background (CMB) anisotropy constraints on PBH accretion are subject to uncertainties in the accretion process and its effect on the thermal history of the universe at early times. To account for this, we show the bounds for both collisional ionisation (solid dark blue line) and photoionisation (dotted dark blue line) Ali-Haïmoud and Kamionkowski (2016). Recently, another sort of accretion limit has been obtained in the mass range from a few to $10^{7}M_{\odot}$ on the grounds that PBH accretion from the interstellar medium should result in a significant population of X-ray sources Inoue and Kusenko (2017). Indeed, several earlier papers have considered such a limit Carr (1979); Gaggero et al. (2016). However, all these limits are very dependent on the accretion scenario and are therefore not shown.

Lensing is the only phenomenon which has been claimed to provide positive evidence for PBHs. For example, the results of the MACHO project – searching for microlensing of stars in the Magellanic clouds – originally suggested halo DM in the form of $0.5M_{\odot}$ objects Alcock et al. (1997) and these could plausibly be PBHs formed at the quark-hadron phase transition at $10^{-5}$s. However, the DM fraction was later reduced to 20% Alcock et al. (2000). The interpretation of the MACHO results – and also the EROS and OGLE results – is very sensitive to the properties of the Milky Way halo. In particular, it has been argued that the recent low-mass Galactic halo models would relax the constraints and allow the halo to consist entirely of solar mass PBHs Hawkins (2015). Where only a constraint is claimed, rather than a positive detection, it is important to specify the associated confidence level (CL). For all lensing constraints shown in Fig. 1, we use the 95% CL constraint given in Refs. Niikura et al. (2017); Griest et al. (2014); Tisserand et al. (2007); Allsman et al. (2001).

There is also positive evidence for microlensing from quasar observations, these indicating that 20% of the total mass is in compact objects in the mass range $0.05-0.45M_{\odot}$ Mediavilla et al. (2017). This is compatible with the expected characteristics of the stellar component and the observations may also exclude a significant population of PBHs outside this mass range. However, this constraint is not used in our analysis because it is difficult to express the result of Ref. Mediavilla et al. (2017) as a quantitative upper bound on the PBH mass fraction. This conclusion may also conflict with recent claim that long-term radio variability in the light-curves of active galactic nuclei (AGN) arises from gravitational millilensing of features in AGN jets Vedantham et al. (2017). This claim could imply that the DM is either individual black holes of mass $10^{3}-10^{6}M_{\odot}$ or clusters of this mass comprising smaller black holes.

In the latter context, additional relaxing of constraints would apply if the PBHs were spatially clustered into sub-haloes. As stressed in Ref. Clesse and García-Bellido (2016), this is expected if PBHs are part of a larger-scale overdense region. However, this effect depends on details of small-scale structure formation which are not fully understood, so we simply adopt the results presented in the current literature.

Observations of neutron stars limits the PBH abundance and indeed it has been claimed that this excludes PBH DM over a wide range of masses. However, these limits are dependent on the DM density in the cores of globular clusters, which is very uncertain. Following Ref. Capela et al. (2013), the neutron star capture constraint is presented for three values of this density (dashed and dot-dashed yellow lines).

It must be stressed that the constraints in Fig. 1 have varying degrees of certainty and they all come with caveats. For some, the observations are well understood (e.g. the CMB and gravitational lensing data) but there are uncertainties in the black hole physics. For others, the observations themselves are not fully understood or depend upon additional astrophysical assumptions. To address the associated uncertainties in a systematic way, we split the constraints into two classes. The first class, presented in Fig. 1 by solid lines, are relatively robust, while the second class, presented by dashed lines, are somewhat less firm and depend upon astrophysical parameters. In particular, this applies to most of the dynamical and accretion constraints (e.g. those associated with dwarf galaxies, wide binaries and neutron stars). However, we stress that this division is not completely clear-cut. In the following, we present our results for the two classes of constraints both separately and together.

III Constraints on extended PBH mass function

If the PBHs span an extended range of masses, the mass function is usually written as ${\rm d}n/{\rm d}M$ where ${\rm d}n$ is the number density of PBHs in the mass range $(M,M+{\rm d}M)$. For our purposes it is more convenient to introduce the function

[TABLE]

normalised so that the fraction of the DM in PBHs is

[TABLE]

where $\Omega_{\rm PBH}$ and $\Omega_{\rm DM}\approx 0.25$ are the PBH and DM densities in units of the critical density. The lower cut-off in the mass integral necessarily exceeds $M_{*}\approx 4\times 10^{14}$g, the mass of the PBHs evaporating at the present epoch Carr et al. (2010). Note that $\psi(M)$ is the distribution function of $\log{M}$ and has units [mass]-1.

In this paper we consider three types of mass function.

A lognormal mass function of the form:

[TABLE]

where $M_{c}$ is the mass at which the function $M\psi(M)$ peaks and $\sigma$ is the width of the spectrum. This was first suggested in Ref. Dolgov and Silk (1993) and is often a good approximation if the PBHs result from a smooth symmetric peak in the inflationary power spectrum. This was demonstrated numerically in Ref. Green (2016) and analytically in Ref. Kannike et al. (2017) for the case in which the slow-roll approximation holds. It is therefore representative of a large class of extended mass functions. Note that the lognormal mass function used in Refs. Green (2016); Horowitz (2016); Kuhnel and Freese (2017) omits the $M^{-1}$ term in Eq. (3). In this case, the position of the peak of $M\psi(M)$ is no longer $M_{c}$ but $e^{\sigma^{2}}M_{c}$. The form (3) is more useful for our purposes because $M\psi(M)$ relates to the DM fraction in PBHs of mass $M$. 2. 2.

A power-law mass function of the the form

[TABLE]

For $\gamma\neq 0$, either the lower or upper cut-off can be neglected if $M_{\rm min}\ll M_{\rm max}$, so this scenario is effectively described by two parameters. Only in the $\gamma=0$ case are both cut-offs necessary. For example, a mass function of this form arises naturally if the PBHs form from scale-invariant density fluctuations or from the collapse of cosmic strings. In both cases, $\gamma=-2w/(1+w)$, where $w$ specifies the equation of state, $p=w\rho$, when the PBHs form Carr (1975). In a non-inflationary universe, $w\in(-1/3,1)$ and so the natural range of the mass function exponent is $\gamma\in(-1,1)$. Equation (4) is not applicable for $w\in(-1,-1/3)$, corresponding to $\gamma\in(1,\infty)$, because PBHs do not form during inflation but only after it as a result of inflation-generated density fluctuations. Special consideration is also required in the $w=0$ (matter-dominated) case Khlopov and Polnarev (1980); Polnarev and Khlopov (1985), because then both cut-offs in (4) can be relevant and this is discussed elsewhere Carr et al. (2017). In the following analysis we will consider both positive and negative values for $\gamma$ but not zero. 3. 3.

A critical collapse mass function Yokoyama (1998b); Niemeyer and Jedamzik (1999); Musco and Miller (2013); Carr et al. (2016b):

[TABLE]

which may apply generically if the PBHs form from density fluctuations with a $\delta$-function power spectrum. In this case, the mass spectrum extends down to arbitrarily low masses but there is an exponential upper cut-off at a mass-scale $M_{f}$ which corresponds roughly to the horizon mass at the collapse epoch. If the density fluctuations are themselves extended, as expected in the inflationary scenario, then Eq. (5) must be modified Carr et al. (2016a). Indeed, the lognormal distribution may then be appropriate. So although the mass function (5) is described by a single parameter, two may be required in the more realistic critical collapse situation.

To compare with the lognormal case, we describe the mass function in the last two cases by the mean and variance of the $\log M$ distribution:

[TABLE]

where $\langle X\rangle_{\psi}\equiv f_{\rm PBH}^{-1}\int{\rm d}M\,\psi(M)X(M)$. For a power-law distribution these are

[TABLE]

where $M_{\rm cut}$ stands for ${\rm max}(M_{\rm min},M_{*})$ if $\gamma<0$ or $M_{\rm max}$ if $\gamma>0$. For the critical-collapse distribution (5), the exponential cut-off is very sharp, so the mass function is well approximated by a power law distribution with $\gamma=3.85$ and $M_{\rm max}\approx M_{f}$. As it is relatively narrow, Eq. (7) implying $\sigma=0.26$, even the monochromatic mass function provides a good fit. Since critical collapse should be a fairly generic feature of PBH formation, $\sigma=0.26$ will usually provide a lower limit to the width of the mass function. However, critical collapse may not be relevant in all cases, for example in the cosmic string or matter-dominated ($w=0$) scenarios.

It should be stressed that two parameters should always suffice to describe the PBH mass function locally (i.e. close to a peak) since this just corresponds to the first two terms in a Taylor expansion. However, in principle the mass function could be more complicated than this. For example, depending on the form of the inflaton potential, it could have several distinct peaks. Indeed, with a sufficiently contrived form, these peaks could be tuned to exactly match all the constraint windows.

The existing constraints on the allowed fraction of PBH DM are commonly presented assuming a monochromatic mass function (presented in the upper panel of Fig. 1). In the following we introduce a simple method for generalising these results to arbitrary mass functions. For this purpose, consider an astrophysical observable $A[\psi(M)]$ depending on the PBH abundance (e.g. the number of microlensing events of given duration in a given time interval). It can generally be expanded as

[TABLE]

where $A_{0}$ is the background contribution and the functions $K_{j}$ depend on the details of the underlying physics and the nature of the observation. If PBHs of different mass contribute independently to the observable, only the first two terms in Eq. (8) need to be considered. Explicit expressions are given for lensing and survival of stars in Ref. Green (2016), evaporation in Ref. Carr et al. (2016a), and neutron star capture and accretion in Ref. Kuhnel and Freese (2017). In this case, if a measurement puts an upper bound on the observable,

[TABLE]

then for a monochromatic mass function with $M=M_{c}$,

[TABLE]

this translates to

[TABLE]

The function $f_{\rm max}(M)$ corresponds to the maximum observationally allowed fraction of DM in PBHs for a monochromatic mass distribution. Combining Eqs. (8)–(11) then yields

[TABLE]

Once $f_{\rm max}$ is known, it is possible to apply Eq. (12) for an arbitrary mass function $\psi(M)$ to obtain the constraints equivalent to those for a monochromatic mass function.

In detail the procedure is as follows. We first integrate Eq. (12) over the mass range ($M_{1},M_{2}$) for which the constraint applies, assuming a particular function $\psi(M;f_{\rm PBH},M_{c},\sigma)$. Once we have specified $M_{1}$ and $M_{2}$, this constrains $f_{\rm PBH}$ as a function of $M_{c}$ and $\sigma$. (In all cases except lensing, we take the integral limits to be the values of $M$ for which $f_{\rm max}=100$.) The last three panels in Fig. 1 are then derived by assuming $\sigma=2$ for the lognormal mass function (upper right panel) and $\gamma=\pm 1$ for the power law mass function (lower panels).

The procedure must be implemented separately for each observable. As shown in the Appendix, different constraints can be combined by using the relation

[TABLE]

where $f_{{\rm max},j}(M)$ correspond to the different bounds for a monochromatic mass function, as defined by Eq. (11). Most of the constraints shown in Fig. 1 rely on a single observable. For lensing this is the number of lensing events Barnacka et al. (2012); Niikura et al. (2017); Tisserand et al. (2007); Allsman et al. (2001), for neutron star capture it is the age of neutron stars Capela et al. (2013), and for white dwarfs and wide binaries it is their abundance Graham et al. (2015); Monroy-Rodríguez and Allen (2014).

However, some monochromatic constraints reported in the literature contain contributions from multiple observables. For example, consider the Planck constraint of Ref. Ali-Haïmoud and Kamionkowski (2016). Earlier analyses calculated the optical depth from the CMB data and used that to constrain the PBH abundance Ricotti et al. (2008); Chen et al. (2016). In this case, there is only one observable, the optical depth, and Eq. (12) is applicable. However, the constraint from CMB anisotropies calculated in Ref. Ali-Haïmoud and Kamionkowski (2016) combines experimental data at various multipole moments by performing a $\chi^{2}$ analysis for $f_{\rm PBH}$ using the CMB data. It is shown in the Appendix that using the combined bound $f_{{\rm max},{\rm CMB}}(M)$ in Eq. (12), instead of considering the contributions from different multipoles separately, will result in a more stringent constraint. However, since the mass dependence of the PBH contribution is expected to be roughly proportional for different multipoles Ali-Haïmoud and Kamionkowski (2016), the error should be small compared to the theoretical uncertainties, shown in Fig. 1, associated with PBH accretion physics. To accurately estimate the size of this error, one should repeat the analysis of Ref. Ali-Haïmoud and Kamionkowski (2016) for each mass function separately, which is beyond the scope of this work.

The important qualitative point is that the form of Fig. 1 in the non-monochromatic case is itself dependent on the PBH mass function. One cannot just compare a predicted extended mass function with the monochromatic form of the constraints, as some authors have done. In displaying the constraints, one also needs to select values of the parameters which describe the mass function. In both the lognormal and power-law cases, we have taken these to be $\sigma$ and $M_{c}$. For the critical collapse model, there is only one parameter ($M_{f}$) but this model is practically indistinguishable from the monochromatic one because only a small fraction of the PBH density is associated with the low-mass tail. So this case is not shown explicitly.

We now discuss some caveats that have to be kept in mind when applying Eq. (12). The mass function evolves in time if the PBH merge or if new black holes are created. This can have an important impact on the constraints. For example, if mergers between recombination and the present are significant, the accretion constraints will be relaxed, since the mass function at recombination would have peaked at a lower mass than today. A period of merging after recombination is not implausible, as this would be induced by the small-scale density fluctuations which are likely to accompany PBH production Hawking (1971); Carr and Hawking (1974).

We next discuss the effect of the higher order terms in (8) since these can induce errors in Eq. (12). These terms become relevant if the contribution of a black hole population with a given mass is influenced by the presence of another black hole population with a different mass. For example, the non-detection of a stochastic gravitational wave background from PBH binaries may constrain the PBH abundance in the near future Wang et al. (2016). This constraint depends on the $K_{2}$ term in Eq. (8) and it may also depend on $K_{3}$ if the formation of PBH binaries depends on $3$-body effects Ioka et al. (1998).

In some cases, it is possible to remove the higher order terms by introducing an ‘effective’ mass function. For example, compact gravitationally bound systems (such as binaries) may behave as a single objects in the context of lensing. Consider an idealised scenario in which the observable depends purely on the total mass of the object, so that $K_{n}(M_{1},M_{2},\ldots,M_{n})=K_{n}(M_{1}+M_{2}+\ldots+M_{n})$, but is otherwise independent of the composition, so that $K_{n}\propto K_{1}$. If we additionally assume that the mass function within these compact bound systems follows the overall mass mass function, we obtain

[TABLE]

Here the effective mass function is given by

[TABLE]

where $\alpha_{n}$ relates to the fraction of $n$-body bound objects,

[TABLE]

and the effective mass function $\psi_{\rm eff}$ has to satisfy the normalisation condition $f_{\rm PBH}\leq 1$. The constraints for the general and monochromatic mass functions are still related by (12) but likely overestimate the allowed PBH mass since $\psi_{\rm eff}(M)$ is always shifted towards higher masses. In principle, all the constraints discussed below and shown in our figures relate to the effective mass functions, which can be different for different constraints.

It is also possible that the mass function is position-dependent. This is expected in dwarf galaxies because mass segregation causes lighter PBHs to migrate outwards, with the heavier ones occupying the central region. This will introduce corrections for constraints arising from the evolution of stars in the Galaxy Koushiappas and Loeb (2017); Brandt (2016). Again, it might be possible to invoke an effective mass function $\psi_{\rm eff}$ that only accounts for the heavier PBHs. However, an estimate of this effect requires detailed numerical simulations which are beyond the scope of this work.

IV Results and discussion

Our main results are presented in Fig. 2, where we show constraints on the maximum allowed fraction of PBH DM, $f_{\rm max}$, in the ($M_{c},\sigma$) plane for lognormal and power-law PBH mass functions. In the upper right panel all the constraints shown in Fig. 1 are considered, using the most restrictive forms for the evaporation, accretion and neutron star constraints, as depicted by the dotted lines. In the other panels only the constraints corresponding to the solid lines are taken into account. We have combined the constraints using Eq. (13). The black lines in Fig. 1 correspond to constant $\sigma$ slices in Fig. 2. The regions where 10%, 20%, 50% and 100% of DM can consist of PBHs are indicated in Fig. 2 by the dotted, dot-dashed, dashed and solid lines, respectively, while less than 0.1% of the DM can be in PBHs in the white region.

The shape of the constraints in Fig. 2 makes it clear that the allowed mass range for fixed $f_{\rm PBH}$ decreases with increasing the width $\sigma$, thus ruling out the possibility of evading the constraints by simply extending the mass function. Moreover, Fig. 2 gives an upper bound $\sigma\lesssim 1$ if all dark matter is in the form of PBHs. This implies $|\gamma|\gtrsim 1$, which effectively rules out PBH DM from the collapse of cosmic strings or scale-invariant density fluctuations.

Our results agree with the conclusions of Refs. Green (2016); Horowitz (2016); Kuhnel and Freese (2017). However, Refs. Green (2016); Horowitz (2016) focus on PBHs in the solar to intermediate mass range, considering microlensing and dynamical constraints from Eridanus II. Ref. Kuhnel and Freese (2017) performs a more comprehensive analysis, covering the mass range $10^{-18}-10^{4}M_{\odot}$, but their study does not include the recent constraint from Subaru Hyper Suprime-Cam Niikura et al. (2017) and they calculate the Planck constraint as in Ref. Chen et al. (2016), resulting in a more stringent constraint than the one from Ref. Ali-Haïmoud and Kamionkowski (2016) used in this work. Also they use the potential SKA pulsar timing constraints Schutz and Liu (2017), even though these are not yet realised. Some of the difference between our figures and those in Refs. Green (2016); Horowitz (2016); Kuhnel and Freese (2017) results from the difference in the definition of $M_{c}$.

The same conclusion can be drawn if one compares the constraints presented in the upper left and right panels of Fig. 1. In the latter case, we show the corresponding $(f_{\rm PBH},M_{c})$ constraints for extended mass functions with fixed width. The effect of the extension is to ‘smooth’ the constraints. Although the most restrictive constraints for the PBH fraction are weakened, it can be seen that the regions allowing a relatively large PBH fraction are reduced. So the constraints become wider, as indicated in Fig. 1. We conclude that previous claims in the literature that wide mass functions allow one to avoid PBH bounds are premature and not supported by our more rigorous computations.

The shape of the colored region of Fig. 2 can be understood as follows: The lognormal mass function is symmetric in the $\log M$ scale, while the power law with $\gamma<0$ has a high-mass tail and $\gamma>0$ is skewed towards low masses. Since the evaporation constraint Carr et al. (2016a) is much stronger than the accretion one Carr et al. (2010), the low-mass tail excludes wider mass functions, whereas $\gamma<0$ allows it.

There are three regions in the upper left panel of Fig. 2 where all DM can consist of PBHs. Two of them are at very low mass, just above the evaporation limit, and the third is in the mass window relevant for the LIGO black hole coalescence events. However, this neglects the dynamical constraints, shown by the dashed lines in Fig. 1. As explained above, this might be justified for reasons associated with the dynamics of the observed astrophysical systems.

To clarify what role different constraints play in the regions of interest, we present these regions in detail in Fig. 3 for $\Omega_{\rm PBH}=\Omega_{\rm DM}$. The masses $25-100M_{\odot}$ satisfy the microlensing and accretion constraints but conflict with dynamical constrains from ultra-faint dwarfs and wide binaries. At the lower mass end, there is a narrow window around $4\times 10^{-17}M_{\odot}$ if we assume a conservative bound from evaporations and another window around $2\times 10^{-14}M_{\odot}$ if the dynamical constrains associated with neutron stars and white dwarfs are neglected. Both masses are in the asteroid range. If all the constraints are taken into account, the maximally allowed fraction of PBH DM is 13% in the high mass window and 30% and 15% in the two low mass windows, respectively. Note that whether the DM can be in PBHs in the asteroid window is sensitive to the form of the PBH evaporation limit and this depends on the precise form of the extragalactic $\gamma$-ray background.

V Conclusions

We have studied the constraints on PBH DM with an extended mass function, presenting a general method for extracting these constraints from those for monochromatic PBH mass functions and discussing possible caveats associated with their interpretation. Our computations cover the broad mass range $10^{-18}-10^{4}M_{\odot}$ and show that extended mass functions do not generally alleviate the already existing constraints on the PBH DM fraction, because the allowed fraction decreases with increasing the width of the mass function. We have identified three mass windows where an appreciable fraction of DM can still consist of PBHs: $4\times 10^{-17}M_{\odot}$, $2\times 10^{-14}M_{\odot}$ and $25-100M_{\odot}$. If all the constraints discussed in the literature are taken at face value and treated on an equal footing, then at most ${\cal O}(10\%)$ of DM can be in PBHs. However, if some of the dynamical constraints can be circumvented, then 100% PBH DM might be allowed in these windows. Even ${\cal O}(10\%)$ DM in the ${\cal O}(10)M_{\odot}$ window might suffice to explain the LIGO events.

Acknowledgements.

We thank Juan Garcia-Bellido, Sebastien Clesse, Anne Green, Gert Hütsi, Chloe James-Turner, Luca Marzola, Paul Schechter and Danton Weil for useful discussions. We also thank an anonymous referee for important input. This work was supported by the Estonian Research Council grants IUT23-6 and EU through the ERDF CoE program grant TK133. T.T. is supported by the U.K. Science and Technology Facilities Council grant ST/J001546/1. *

Appendix A Combined constraints

In this appendix we describe how different constraints can be combined. Consider $N$ independent observables $A_{j}(f_{\rm PBH})$ with observed expectation values $\mu_{j}$ and variances $\sigma_{j}^{2}$. The $\chi^{2}$ for these observables is

[TABLE]

We assume that the observables $A_{j}(f_{\rm PBH})$ are linear in the mass function, so that only the first two terms in (8) are relevant, and that the mean values coincide with $f_{\rm PBH}=0$, which implies $\mu_{j}=A_{0,j}$. Since $\psi\propto f_{\rm PBH}$, the $n\sigma$ constraint on $f_{\rm PBH}$ is then

[TABLE]

where $\chi_{\rm min}$ is the minimum of $\chi$. As in Eq. (11), the kernel $K_{1,j}$ can be extracted from the constraint for a monochromatic mass function if $N=1$. This corresponds to

[TABLE]

It follows that Eq. (18) can be recast as

[TABLE]

Note that the upper bound in Eq. (9) is $A_{{\rm exp},j}=A_{0,j}+n\sigma_{j}$, where $n\sigma_{j}$ is the confidence level of the bound.

Consider then the combined constraint for a monochromatic mass function, which from Eq. (20) can be expressed as

[TABLE]

Using this in Eq. (12) will always lead to an overestimation of the actual constraint given by (20) because the triangle inequality implies

[TABLE]

For a more precise estimate, the contribution of each observable has to be included separately for any constraint derived from multiple observables. However, if the constraints $f_{{\rm max},j}(M)$ are proportional to each other, the expressions (20) and (12) are equivalent, so no error is made.