On the Prevalence of Super-Massive Black Holes over Cosmic Time

TL;DR

This study models the abundance and distribution of super-massive black hole seeds over cosmic time, constraining their formation and evolution through simulations and observations, and introduces a new diagnostic for gravitational wave detection.

Contribution

It provides a comprehensive, seed-mechanism-independent analysis of SMBH prevalence, linking halo mass, seeding probability, and observable signatures across cosmic history.

Findings

Seeding efficiency depends on critical halo mass M_c.

Constraints on SMBH density from local and high-redshift observations.

Predictions for SMBH merger rates and gravitational wave signals.

Abstract

We investigate the abundance of Super-Massive Black Hole (SMBH) seeds in primordial galaxy halos. We explore the assumption that dark matter halos outgrowing a critical halo mass M_c have some probability p of having spawned a SMBH seed. Current observations of local, intermediate-mass galaxies constrain these parameters: For $M_c=10^{11}M_\odot$, all halos must be seeded, but when adopting smaller M_c masses the seeding can be much less efficient. The constraints also put lower limits on the number density of black holes in the local and high-redshift Universe. Reproducing z~6 quasar space densities depends on their typical halo mass, which can be constrained by counting nearby Lyman Break Galaxies and Lyman Alpha Emitters. For both observables, our simulations demonstrate that single-field predictions are too diverse to make definitive statements, in agreement with mixed claims in the literature. If quasars are not limited to the most massive host halos, they may represent a tiny fraction (~10^-5) of the SMBH population. Finally, we produce a wide range of predictions for gravitational events from SMBH mergers. We define a new diagnostic diagram for LISA to measure both SMBH space density and the typical delay between halo merger and black hole merger. While previous works have explored specific scenarios, our results hold independent of the seed mechanism, seed mass, obscuration, fueling methods and duty cycle.

Figures

Equations

$$N_{z=0}=\frac{1}{2}(M_{z=0}/M_c{)}^{3/4}$$

$$P(\mathrm{BH}|N)=1-{\left(1-p\right)}^N$$

$$p\text{apprge}(M_c/10^{11}{M}_{\odot }{)}^{3/4}.$$

$$\bar{A}:=\mu -4\sigma .$$

Full text

quasars: supermassive black holes, galaxies: halos, galaxies: high-redshift, galaxies: evolution

1 Introduction

Supermassive Black Holes (SMBHs) are ubiquitous in local, massive galaxies (Kormendy & Ho, 2013; Graham, 2016). This raises the question of how these black holes came into existance, and how they evolved over cosmic time. Mass growth of black holes over cosmic time is thought to be traced by active galactic nuclei (AGN; e.g., Soltan, 1982; Yu & Tremaine, 2002). The physical process creating black holes in the first place remains unknown, although multiple scenarios have been proposed, including collapses in the cores of massive stars, dense nuclear star clusters or pristine gas clouds (Rees, 1984; Latif & Ferrara, 2016). These lead to different initial mass regimes. Much research has been dedicated to develop detailed physically meaningful seed mechanisms and timely mass growth (see e.g., Volonteri, 2010; Reines & Comastri, 2016). This work however addresses a simpler problem: How frequently do massive black holes need to arise to explain current observations?

Modern cosmological simulations now routinely incorporate SMBHs to quench star formation in massive halos (see e.g. Croton et al., 2006). For computational convenience, the most common method of creating SMBHs is to seed dark matter halos once they reach a certain critical mass, $M_{c}$, e.g., $M_{h}>M_{c}=10^{10-11}M_{\odot}$ (Somerville & Davé, 2015). However, massive galaxies in the local Universe have been assembled from many such progenitor halos. Menou et al. (2001) computed that one can achieve high SMBH occupation fractions at the massive end by only seeding 3% of $M_{c}\approx 10^{9}M_{\odot}$ halos at redshift $z=5$. They demonstrated this by letting the seeds trickle down simulated merger trees, the representation of the hierarchical galaxy evolution. They then made predictions about the rate of SMBH births over cosmic time, merger rates (relevant for space-based long-baseline gravitational wave experiments). Since then, our understanding of the cosmology has improved, numerical simulations have become more sophisticated and local observations give a better idea down to which mass limit SMBHs occupy galaxies.

While dedicated simulations tend to be highly specific to a particular seeding and feeding scenario (see, e.g., Volonteri, 2010; Naab & Ostriker, 2017), occupation calculations under this framework are independent of the mass of the seeding process and the mass growth of these seeds. All that matters is whether a massive black hole is present or not.

Based on reliable halo trees from dark-matter simulations presented in §2, this work makes three main contributions: 1) Section 3 explores how efficient seeding mechanisms need to be to explain current black hole occupation observations. 2) Section 4 explores the evolution of the required black hole population. In §5 we make predictions for the future gravitational wave experiment LISA, and demonstrate how the observed events can be used to infer both the space density of SMBHs and merger inspiral delays. 3) Section 8 focuses on the environment of high-redshift quasars, and clarifies how nearby galaxies trace the halo mass of quasar hosts.

2 Simulations

The premise of our analysis is that galaxy halos that grow to a critical mass $M_{c}$ (e.g., $10^{10}M_{\odot}$) have a chance probability $p$ to have produced a seed black hole by this point. For an individual galaxy the likelihood of forming a seed black hole will in reality depend in some unknown manner on other properties (e.g. environment, radiation fields, gas metallicity, evolutionary history). Considering however the halo population crossing a given mass threshold, we can define an effective occurance rate $p$ up to that evolutionary stage. If the future evolution of the halo and the seeding probability correlate only weakly, this formalism allows an estimate of the efficiency of the seeding process. We explore the implications of such seeding recipes and viable parameter ranges for $M_{c}$ and $p$.

Connecting (potentially early) black hole seeding with observations in the local universe requires simulations that follow the mass evolution of the Universe. High mass resolution is necessary to follow sites of emerging proto-galaxies. Additionally, because seeds may be rare, large cosmological volumes need to be probed. To this end, we use the MultiDark cosmological dark matter N-body simulations (Prada et al., 2012; Klypin et al., 2016). These assume a Planck cosmology ($h=0.6777$, $\Omega_{\Lambda}=0.693$, $\Omega_{m}=0.307$, Planck Collaboration et al., 2014), which we use throughout111We explicitly multiply out $h$ e.g. in distances, but for comparison with literature we keep halo masses in units of $h^{-1}$.. The simulation evolves an initial dark matter density distribution over cosmic time under gravity. This encompasses the gravity-dominated collapse into sheets, filaments and finally halos, wherein galaxies should reside, as well as the mergers of structures. To study the evolution of halos, merger trees222Constructed with ROCKSTAR and ConsistentTrees (Behroozi et al., 2013a, b). summarise the merging of (sub-)halos as well as the dark matter halo mass at each simulation snapshot. In particular we focus on the highest-resolution, Small MultiDark simulation (SMDPL) with a box size of $400\mathrm{Mpc/h}$, populated with $3840^{3}$ dark matter particles of mass $10^{8}M_{\odot}/h$, which resolves well halos of masses down to $10^{10}M_{\odot}$. When smaller halo resolution is needed, we use the $40\mathrm{Mpc}/h$ “Following ORbits of Satellites” (FORS; González & Padilla, 2016) simulation, whose small particle masses ($\approx 4\times 10^{6}M_{\odot}/h$) resolve halos down to $10^{8}M_{\odot}$ well. The FORS cosmology is very similar to the above-mentioned values. This work adopts virial masses in $M_{\odot}h^{-1}$ throughout333To convert into physical units, our reported masses need to be multiplied by $h$ and divided by $\sim 0.75$ to account for the impact of baryons (see e.g., Sawala et al., 2013). and physical units for all other quantities.

Throughout, we refer to the created black holes as seeds or (S)MBHs. However, because we only count the halo occupation, our model does not need to assign any specific mass to them. Our calculations are thus independent of seed mass and mass evolution. Consequently, we refrain from exploring correlations of galaxy properties with black hole mass. By-products of the seeding process that do not become MBHs are not considered here, but may additionally be present in the Universe.

3 SMBHs in local galaxies

Owing to the hierarchical growth of structures, a high number of halos go into building massive galaxies. When halos merge, occupant black holes are generally passed on to the merger product (this assumption is discussed below). This implies that the fraction of halos that need to be seeded could be very low (Menou et al., 2001). As a starting point, we simply count for $z=0$ halos, from how many $10^{10}M_{\odot}$ halos they were built. Specifically, we count in each merger tree the leaves with $M>M_{c}$ . This is shown in Figure 1. If galaxies only were built by mergers of existing halos, we would expect a linear relation, $N_{z=0}=M_{z=0}/M_{c}$ (gray line). However, growth through accretion of non-collapsed structures (e.g., filaments) is also important. A better empirical description is a sub-linear powerlaw relation with sub-unity normalisation:

[TABLE]

This is plotted as a red dashed curve in Figure 1, and is an appropriate approximation as we will show below. Around this mean number of constituents, Figure 1 shows modest scatter of approximately $0.2\mathrm{dex}$. We verified that this relation also holds for other values of $M_{c}$ ($10^{8.5-11}M_{\odot}$).

We then compute the black hole occupation fraction. If a halo is made from $N$ building blocks, each with equal chance $p$ to contain a seed black hole, it contains a black hole with probability:

[TABLE]

This is derived by considering the probability that none of the $N$ halos have a MBH, and taking the complement. A subtle point here is that we do not need to keep track of the number of black holes inside a halo and whether/when they merge or get ejected from the system. Equation 2 only assumes that if multiple halos with black holes merge, at least one of these black holes remains in the merged halo. The average of $P(BH|N)$ over a halo population yields the expected occupation fraction. Figure 2 plots the occupation fraction as a function of halo mass and seeding fraction $p$, for the case of $M_{c}=10^{10}M_{\odot}$. We now constrain $p$ using observations.

We focus on a recent measurement of the fraction of local galaxies containing a massive black hole as a function of mass by Miller et al. (2015, M15). They surveyed local early-type galaxies for central X-ray point sources, a telltale sign of accretion onto SMBHs, and carefully corrected for flux limits. Correcting for the X-ray inactive fraction of black holes is systematically uncertain. The correction anchors on the X-ray active fraction at the higher masses, where the occupation fraction is thought to be $\sim 100\%$ (Greene, 2012). M15 used an advanced methodology to incorporate that the X-ray luminosity distribution changes with host mass. They detect sources down to and below host galaxy stellar masses of $10^{9}M_{\odot}$ and their inferred lower limits on the black hole occupation fraction are shown in the bottom panel of Figure 3 as a thick black dashed line. Greene (2012) previously reviewed observational constraints and also considered the late-type spiral sample of Desroches & Ho (2009), obtaining slightly elevated lower limits (gray connected points in Figure 3).

To compare our simulation to these observational constraints, we convert from halo masses to stellar masses. We test several methods: Firstly, we assign stellar masses from the Moster et al. (2010) conditional stellar mass function $P(M_{\star}|M_{h})$444We only consider the central galaxies because we work with a subhalo catalogue and satellites are unimportant in the mass range of interest. derived from matching the local stellar mass function and galaxy clustering to simulated dark matter halos. Importantly, this empirical method does not assume any galaxy evolution physics and is observational. Secondly, we try assigning stellar masses according to the distribution produced by the SAGE semi-analytic model for the MDPL2 simulation (Croton et al., 2016) and that produced by the EAGLE hydro-radiative simulation. The latter (Schaye et al., 2015, version Recal-L025N0752) reproduces the galaxy mass function and sizes very well, in particular in this mass regime. However, populating halos based only on a $M_{\star}-M_{h}$ distribution neglects that the number of progenitors $N$ could influence $M_{\star}$ at a given $M_{h}$, inducing a correlation between occupation probability and $M_{\star}$. To test this, we also derive the number of progenitors $N$ and $M_{\star}$ from the merger trees of EAGLE. The groups of curves in Figure 3 demonstrate that for a given $p$ the results from these four methods are consistent. Thus, the choice of $M_{h}-M_{\star}$ conversion method does not materially affect our conclusions.

Figure 3 shows the occupation fraction predictions for $M_{c}=10^{10}M_{\odot}$. The comparison with the observational lower limit implies $p\gtrsim 30\%$, with lower seed probabilities ruled out. For higher masses, even higher seed probabilities are necessary, e.g., $M_{c}=10^{11}M_{\odot}$ matches observations only with full occupation ($p=100\%$). We explore the allowed seeding fractions $p$ by extrapolating equations 1 and 2 into lower mass regimes. We compute curves like in Figure 3 with the mean $M_{\star}/M_{h}$ ratio of Moster et al. (2010) and accept those that pass the M15 limits. Figure 4 shows the allowed parameter space. Because the number of halo building blocks has a $3/4$ powerlaw exponent (see Figure 1 and eq. 1), we find the constraint:

[TABLE]

In the following, we consider only models satisfying this constraint.

4 The SMBH population over cosmic time

To explore how the SMBH population evolves over cosmic time, we populate the halo tree according to our recipe. When a halo for the first time reaches a mass $M_{c}$, with chance $p$ it is marked as occupied by a black hole. The occupation is inherited by the merger tree descendant. As a consequence, in our seeding framework black holes are not only “born” in the high-redshift Universe. To show this, the left panel of Figure 5 presents birth rates over cosmic time. Seeding is most frequent at early times ($z>4$). Depending on the model parameters, the median seeding time lies in the second, third or forth billion years. At $z<6$, the model curve shapes are almost independent of the input parameters. Regarding the normalisation, we note that a $100\%$ seeding above mass $M_{c}$ creates a black hole population at a similar rate as $10\%$ seeding at mass $10\%\times M_{c}$ (left panel of Figure 5).

As halos containing black holes merge, interactions of multiple SMBHs are possible. The middle panel of Figure 5 shows mergers of halos hosting black holes. The normalisation again scales as just described for the birth rate. Changes in the critical mass shift the peak slightly. We explore the possible merging of the black holes themselves in the next section.

The right panel of Figure 5 presents the total number of massive black holes over cosmic time. For $p=1$, this is just the number of halos above the mass threshold. Overall, we find a present-day SMBH space density of $>0.01\,\mathrm{Mpc^{-3}}$ (with models not ruled out in Figure 4). In the local Universe, the fraction of occupied halos is mass-dependent in the way presented in Figure 3. The comparison between the number of SMBHs to the observed number of AGN provides the fraction of actively accreting black holes in excess of that luminosity threshold. Modern surveys detect hard X-ray emission even in heavily obscured AGN, and recover the intrinsic accretion luminosity using X-ray spectra (e.g., Buchner et al., 2015). In Figure 5 we include the space density of AGN at $L(2-10\mathrm{keV})>3\times 10^{42}\mathrm{erg/s}$ and $>10^{44}\mathrm{erg/s}$, which correspond approximately to accretion rates of $5\times 10^{6}$ and $5\times 10^{8}M_{\odot}/\mathrm{Gyr}$ (Marconi et al., 2004), respectively, assuming 10% radiative efficiency. The AGN density inferred from observations is $1-3$ orders of magnitude below the total SMBH population in Figure 5. The comparison indicates that the vast majority of the SMBH population is dormant or accretes at low levels.

5 Gravitational Waves from SMBH mergers

The merging of halos in principle brings their SMBHs together as well. The merging of SMBHs produces gravitational waves, detectable with the proposed Laser Interferometer Space Antenna (LISA, Amaro-Seoane et al., 2017). Salcido et al. (2016) predicted in detail waveforms of SMBH mergers from the EAGLE cosmological simulation. Because of the sensitivity down to $10^{4}M_{\odot}$ (Amaro-Seoane et al., 2017), Salcido et al. (2016) find that essentially all mergers are within the sensitivity of LISA, with their first mergers of the low-mass SMBH seeds dominating. Therefore, one can neglect mass and sensitivity considerations and focus on the occurance of SMBH merger events. We predict the number and redshifts of observable GW events from our scenarios. In this computation, the rest-frame black hole mergers at a given redshift $r(z)$ within a redshift interval $dz$ in a comoving volume $dV_{c}$ are converted into observable rates across the entire sky using $\frac{dz}{dt}\times\frac{dV_{c}}{dz}\times\frac{1}{1+z}$ to account for the corresponding rest-frame time, comoving volume on the sky and time dilation (e.g., Sesana et al., 2004; Salcido et al., 2016).

However, the time delay between halo merger and black hole merger is uncertain. Kelley et al. (2017) explored in detail the physical processes as SMBH binaries overcome nine orders of magnitude in separation from kilo- to micro-parsecs. They also present several physically reasonable delay functions, quantifying the fraction of SMBH binaries have merged after a given time. These are shown in the inset of Figure 6. These delay functions have the effect of changing the observed merger redshift distribution in a characteristic fashion. For example, the dashed curve in the main panel of Figure 6 is assuming no delay (see also Salcido et al., 2016) between halo merger and black hole merger, while each solid curve applies a different delay function. These generally move the peak to lower redshifts, narrow the distribution and can suppress the rate of black hole mergers substantially. Overall, all of the distributions appear approximately Gaussian in the scale factor, $\log a=-\log(1+z)$.

We define two LISA observables: The number of GW events observable per year, $N$, and, $\bar{A}$, which describes the location and shape of the distribution and is defined as:

[TABLE]

The $\bar{A}$ statistic combines the mean $\mu$ and standard deviation $\sigma$ of the $\log(1+z)$ distribution (see for example Figure 6). Through experimentation, we found that an $\bar{A}$ centered on $\mu$ but skewed 4$\sigma$ to the left of the distribution captures well both the shift and narrowing of the redshift distribution caused by merger delays. It can be readily computed from detected gravitational wave events, assuming a cosmology.

These two observables diagnose the underlying SMBH population and differentiate delay functions. Figure 7 presents our LISA diagnostic diagram, populated by all combinations of the seeding prescriptions and delay functions. The number of events predicted (y-axis) generally reflects the overall space density of SMBHs (purple contour curves and arrow), and is thus primarily driven by the chosen seeding prescription. It is notable that a wide range of predictions are possible, from few to thousands of events per year. However, even in the most unfavorable scenario with slow inspirals and a high mass threshold ($10^{11}M_{\odot}$), a few detections per year are predicted. The color-coding of the LISA model points in Figure 7 indicates the median delay time corresponding to a delay function in the inset of Figure 6. Because the colors approximately change along the x-axis, the $\bar{A}$ statistic (x-axis) is a proxy for the delay function, and can distinguish slow from quick inspirals (green to light blue points, from left to right). To quantify $\bar{A}$ with an uncertainty of $\pm 0.05$ requires approximately 200 GW detections.

We verify that we reproduce the same results as Salcido et al. (2016) under their assumptions ($M_{c}=10^{10}M_{\odot}$, $p=1$, no or short delays). In other words, their GW event rates are special cases of our framework. As noted there, however, Sesana et al. (2007) obtained substantially different results, as in their models seeding peaks at $z\sim 10$ and ends by $z\sim 6$. In that case, $\bar{A}>0$555Specifically, $(N,\bar{A})$ for their models are: VHM: (10.67, 0.143), KBD: (79.75, 0.26) BVRhf: (0.97, 0.076) BVRlf: (4.93,0.015)., because the average $\log(1+z)$ is higher ($\sim 1$ instead of $\sim 0.5$). Nevertheless also in those models, merger delays will shift the GW events to lower redshifts, and thus the LISA diagnostic diagram can diagnose the shift within a particular model. As discussed in Salcido et al. (2016) and Barausse (2012), the number of GW events and their redshift distribution are powerful diagnostic to compare seeding models. The LISA diagnostic diagram, based on $(N,\bar{A})$ or just $(N,\mu)$, provides a useful visual summary of model predictions.

6 SMBHs active as quasars at $z\sim 6$

The earliest census of the SMBH population is available from very high-redshift quasar surveys. The Sloan Digital Sky Survey revealed several dozens of high-redshift optical quasars at $z\sim 6-8$ (e.g., Fan et al., 2001; Jiang et al., 2016). Because quasars are interpreted as accreting SMBHs, the quasar number density measurement places a lower limit666It is a lower limit because this selection misses obscured, faint and dormant black holes, and possibly, by virtue of the Eddington limit, also those with the lowest masses. As shown before, the fraction of dormant SMBH can be very high at least in the local Universe, and so is the fraction of obscured AGN (e.g., 75% in Buchner et al., 2015). on the black hole volume density. An assumption made by some previous works is that these quasars live in the most massive halos at that time (e.g., Li et al., 2007; Volonteri & Rees, 2006; Wyithe & Padmanabhan, 2006; Romano-Díaz et al., 2011). This makes it easier for modelers to explain high black hole masses. In this and the next section we explore the consequences of this halo density matching assumption and relax it.

Figure 8 shows the observational lower limit as a dotted horizontal line for quasars (Willott et al., 2010b; Jiang et al., 2016) and AGN (Onoue et al., 2017). Well above this, the top-most curve shows the cumulative halo mass function, i.e., the number of halos above a given mass. If quasars populated only the most massive halos, e.g., $M>10^{12.5}M_{\odot}$, the space density would match observations. If instead quasars are permitted to inhabit lower masses, e.g., $M_{h}\sim 10^{11}M_{\odot}$, only a very small fraction ($f\approx 10^{-6}$, bottom panel) of halos need be quasars. According to our seeding prescription, we expect all of these halos to have SMBHs, but their triggering as quasar requires closer investigation.

6.1 Environments of $z\sim 6$ quasars

If quasars at $z\sim 6$ indeed live in very massive halos, galaxy over-densities should be observable around them. Identifying over-densities has however produced contradictory results in the literature. In this section, we use mock observations to understand the difficulties.

From the simulations, we can predict the number of galaxies near quasars. For this, we consider three cases for the minimal quasar halo masses in the $3\times 10^{11}-3\times 10^{12}M_{\odot}$ range (squares in the lower panel of Figure 8), and draw 1000 halos randomly above that threshold. For a mock observation, we then search for surrounding galaxy halos, taking into account the angular separation and redshift selection window. We assign galaxies to halos to mimic observations.

Lyman break galaxies (LBGs) can be found at these redshifts through filter drop-outs (e.g., Stanway et al., 2003; Zheng et al., 2006; Bouwens et al., 2011; Ota et al., 2018). This selects a relatively broad redshift range (we adopt here $z=5.5-7$; e.g., Ota et al., 2018) to a typical detection limit of 26th magnitude AB (e.g., Utsumi et al., 2010; Morselli et al., 2014; Ota et al., 2018). Following the LBG clustering analysis in the UltraVISTA and UDS fields of Hatfield et al. (2017), we use the simple prescription that LBGs occupy halos of $M>10^{11.2}M_{\odot}$, and verify that this matches the number detected in those (non-QSO) fields. Selecting a cylinder as high as the redshift selection ($\sim 500\mathrm{Mpc}$) and projecting it, we take into account chance alignments. For each quasar we compute the number of LBGs enclosed within an angular separation of $\theta$. This is shown in the top panel of Figure 9. The mean number is generally higher in the high-mass case (red), but there is substantial overlap in the 95% confidence intervals. This indicates that a single quasar field is not powerful enough to discriminate halo mass occupations, and wide variations of numbers are expected.

This is indeed what is observed: Some studies claim overdensities (e.g., Zheng et al., 2006), while others find densities comparable to control fields (e.g., Bañados et al., 2013; Morselli et al., 2014; Mazzucchelli et al., 2017; Ota et al., 2018). In the middle panel of Figure 9 we compiled some studies777The luminosities of the targets are consistent with the quasar definition used in Figure 8. The LBG magnitude detection limits in these studies are comparable. and present as thick curves the enclosed number of LBGs as a function of angular distance. Overall, the data are consistent with all three model ranges at all radii. The Zheng et al. (2006) data set888We have considered their C-complex as a single source and chose those objects only $1\sigma$ above their color-cut as less secure candidates. could be an exception just where their field-of-view ends, as well as one of the Morselli et al. (2014) fields at the low end. Clearly, more observations are needed to reduce the variations. The bottom panel of Figure 9 shows model predictions for 20 quasar fields. At this point the most massive scenario can be unambiguously distinguished.

Some of the prediction variance is likely due to the broad redshift range of the LBG selection ($\Delta z\approx 1.5$ corresponding to $\sim 500\mathrm{Mpc}$). To focus on the environment near the quasar, studies have advocated the use of Lyman Alpha emitters (LAEs) to complement LBGs studies. To detect the Ly$\alpha$ line near the target quasar redshift requires narrow (and sometimes custom-made) filters, but has the benefit of probing a narrow redshift range ($\Delta z\sim 0.1$). We again adopt a simple prescription, following Kovač et al. (2007) and Sobacchi & Mesinger (2015), and assign LAEs to halos above a halo mass limit of $M>10^{10.6}M_{\odot}$. We verify that the numbers in non-QSO fields, SXDS (Ouchi et al., 2010) and SDF (Ota et al., 2018), are reproduced. These can also be reproduced by choosing a lower mass cut and a duty cycle, but we find this does not change our conclusions significantly. We again compute the enclosed number of neighbors to the quasar at various angular separations. The top panel of Figure 10 shows the predictions, which are indistinguishable for a single quasar field. Observations (middle panel; Bañados et al., 2013; Mazzucchelli et al., 2017; Goto et al., 2017; Ota et al., 2018) fall within the predicted ranges, which are extremely wide, in part due to Poisson statistics. The only measurement falling outside one of the predicted ranges is that of Goto et al. (2017), which did not detect any sources999This is the same field as the Utsumi et al. (2010) LBG measurement.. This may indicate that our LAE prescription is too abundant. The bottom panel shows predictions when 10 quasar fields have been observed. In that case, the predictions marginally separate.

The LAE and LBG populations are both promising approaches to constrain the halo mass of quasars, and thereby the fraction of black holes active as quasars (Figure 8). However, because dark matter halo neighborhoods are diverse, observations of many quasar fields are needed to make definitive statements.

7 Discussion & Conclusion

We made an analysis of SMBH seeding that is independent of the seed mass, growth mechanism and feedback processes, to estimate the number of massive black holes across cosmic time. Our assumption, following Menou et al. (2001), is simply that by the time a dark matter halo reaches a critical mass $M_{c}$, some process has seeded it with a black hole with efficiency $p$. Details of the seeding process involved are not needed in our analysis.

The high black hole occupation fraction observed in local low-mass galaxies constrains these parameters. The seeding fraction $p$ has to exceed $p\apprge(M_{c}/10^{11}M_{\odot})^{3/4}$ (see Figure 4). This simultaneously constrains the underlying population, for example the local SMBH space density is above $>0.01/\mathrm{Mpc}^{3}$ in all our models. In our framework, black hole seeding happens continuously throughout cosmic time, but is generally most frequent at $z>4$. If we require seeding to only occur before the epoch of reionization ($z>7$), generally lower mass thresholds are required to increase the population (dotted line in Figure 4). However, the merger distributions shift only slightly to higher redshifts. Indeed, when seeding only a small fraction of halos (e.g., $p=1\%$), a redshift cut-off has virtually no effect on the shape of the host merger history.

Our framework can approximate the behavior of physically-motivated seeding mechanisms. Some examples are shown in Figure 4 for comparison. The light seed scenario of Volonteri et al. (2003) considers the end-products of Population III stars. This mechanism effectively populates all $M_{c}\gtrsim 10^{7}M_{\odot}$ halos at $z\sim 20$ with SMBH seeds, and as such taken as a special case of our formalism. Heavy seed scenarios operate at higher masses $M_{c}\gtrsim 3\times 10^{7}M_{\odot}$ halos at $z\sim 10$, but are less efficient ($p=3-30$, e.g., Volonteri et al. 2008; Agarwal et al. 2012). These require UV radiation from neighboring galaxies, which usually merge after the host halo creates the seed (by z 6, e.g., Agarwal et al., 2014). Their effective behavior may thus perhaps be better approximated with a higher mass threshold with lower efficiencies down to lower redshift. In practice however, the differences in the model predictions at lower redshift ($>1\mathrm{Gyr}$ after seeding) are negligible. If we apply a $z>15$ constraint to our seeding prescription this yields very stringent constraints on the parameter space, ruling out $M_{c}\gtrsim 10^{9}M_{\odot}$, and giving $p\gtrsim 20\%$ for $M_{c}=10^{8}M_{\odot}$, the smallest halo masses our simulations can reliably probe. These efficiency constraints are also consistent with Greene (2012), which showed that the aforementioned models are near current observational constraints. In any case however, all surviving models have populated all $M>10^{11}M_{\odot}$ halos with SMBH seeds at $z\sim 6$.

The quasar population at $z\sim 6$ can also be connected to our prescription. Willott et al. (2010a) argued that the observed $z\sim 6-7$ quasars are the $M_{\mathrm{BH}}>10^{7}M_{\odot}$ population accreting at the Eddington-limit, while the remaining black holes are pristine seeds without substantial accretion. Under this interpretation we consider $f$, the success fraction of turning a seed into a quasar at that redshift. Indeed, the above constraints indicate that at $z\sim 6$, there is an abundance of black holes. Fewer than $\sim 10^{-6}$ of the seeds are quasars at that time. If seeds become SMBH and potential quasars only above a certain halo masses, the active fraction can still be as low as $f\sim 10^{-5}$, indicating that we may only see very “lucky” seeds.

Thus, our constraints show that physical seeding mechanisms and the activation of quasars can be highly inefficient. For example, when studying a comoving cosmological volume of $20\,\mathrm{Mpc}$ side length, the theorist needs to create $\sim 80$ seeds by $z\approx 6$, but only one in a million of those need become a quasar with a SMBH. It is thus encouraging to consider seeding and feeding mechanisms requiring rare conditions, such as special configurations of halos or multiple, complex mergers (see e.g., Agarwal et al., 2012; Mayer & Bonoli, 2018; Inayoshi et al., 2018), rather than focusing on massive halos.

Several works have studied overdensities to measure the halo mass of quasars, with mixed results. We demonstrate that the diversity of single-QSO field observations is expected. This is because even at a given mass, halos have diverse neighbourhoods, with the number of surrounding halos varying by an order of magnitude. Previously, Overzier et al. (2009) came to a similar conclusion considering a mock field observation of i-band dropout galaxies. Observations of several quasar fields are necessary to make definitive statements about the host halo mass of quasars. LBGs and LAEs are both suitable probes. However, given that detection of LAEs requires filters specific to the targeted quasar redshift, observations of LBGs may be more economical. In the clustering predictions we have made highly simplified assumptions, especially in how LAEs and LBGs populate and cluster around halos at high redshift. This is still an open research question also in non-quasar fields.

A promising future probe of the seeding scenarios is the space-based LISA gravitational wave experiment (eLISA Consortium et al., 2013), because of its broad mass and redshift sensitivity. We have presented a new diagnostic diagram for LISA GW events. The number of observed events and their redshift distribution can probe simultaneously both the space density of the SMBH population and the typical delay between halo and black hole mergers. We populate the diagram with various seeding scenarios, taking into account possible delays between halo mergers and black hole coalescence, and find that the local SMBH occupation constraints already imply a lower limit of at least one GW event per year detected by LISA. Current SMBH constraints are lower limits on their occurance, as inactive SMBH can often go undetected. To strongly constrain the parameter space, e.g., to near the dashed line in Figure 4, it would be necessary to establish that 50% of the $M_{\star}\sim 10^{9}M_{\odot}$ galaxies lack a SMBH (under some mass definition of SMBH seeds). This is challenging because with an imaging resolution $\theta$ the black hole sphere of influence is detectable only to distances $D\approx 2.3\mathrm{Mpc}\times(\theta/0.1^{\prime\prime})\times\sqrt{M_{\mathrm{BH}}/10^{6}M_{\odot}}$ (e.g., Do et al., 2014)101010This contains $M_{\mathrm{BH}}\propto\sigma^{4}$, however the relation is uncertain at low masses (see, e.g., Kormendy & Ho, 2013; Graham, 2016)). Thus many local galaxies are beyond the resolution limit of current instruments (Ferrarese & Ford, 2005), but can hopefully be probed with upcoming Extremely Large Telescopes (Do et al., 2014). Thanks to its sensitivity to black hole mergers of various sizes, LISA will provide the most powerful constraints on the abundance of SMBHs. Figure 7 shows that the number of LISA detections per year correlates with the local black hole number density, which relates to a $M_{c}\times p$ combination (right panel of Figure 5).

8 Acknowledgements

JB thanks Roberto E. Gonzalez and Nelson D. Padilla for feedback on the manuscript. We acknowledge support from the CONICYT-Chile grants Basal-CATA PFB-06/2007 and Basal AFB-170002 (JB, FEB, ET), FONDECYT Regular 1141218 (FEB) and 1160999 (ET), FONDECYT Postdoctorados 3160439 (JB), CONICYT PIA Anillo ACT172033 (ET), and the Ministry of Economy, Development, and Tourism’s Millennium Science Initiative through grant IC120009, awarded to The Millennium Institute of Astrophysics, MAS (JB, FEB). LFS and KS acknowledge support from SNSF Grants PP00P2_138979 and PP00P2_166159. Discussion of some aspects of this work was carried out at the Aspen Center for Physics, which is supported by National Science Foundation grant PHY-1607611.

The authors gratefully acknowledge the Gauss Centre for Supercomputing e.V. (www.gauss-centre.eu) and the Partnership for Advanced Supercomputing in Europe (PRACE, www.prace-ri.eu) for funding the MultiDark simulation project by providing computing time on the GCS Supercomputer SuperMUC at Leibniz Supercomputing Centre (LRZ, www.lrz.de). The MultiDark simulations were performed on Pleiades supercomputer at the NASA Ames supercomputer centre. The SMDPL simulation have been performed on SuperMUC at LRZ in Munich within the pr87yi project. The MultiDark Database used in this paper and the web application providing online access to it were constructed as part of the activities of the German Astrophysical Virtual Observatory as result of a collaboration between the Leibniz-Institute for Astrophysics Potsdam (AIP) and the Spanish MultiDark Consolider Project CSD2009-00064. The Geryon cluster at the Centro de Astro-Ingenieria UC was extensively used for the FORS simulation. The Anillo ACT-86, FONDEQUIP AIC-57, and QUIMAL 130008 11 provided funding for several improvements to the Geryon cluster.