Cosmological simulations of massive black hole seeds: predictions for next generation electromagnetic and gravitational wave observations

TL;DR

This study uses cosmological simulations to compare different models of supermassive black hole seed formation, predicting their observable signatures in electromagnetic and gravitational wave data to constrain early Universe black hole origins.

Contribution

A novel method to analyze black hole growth histories within the Illustris simulation framework for various seed formation models, highlighting their observational differences.

Findings

Seed models differ in black hole number densities and halo occupation fractions.

Faint-end quasar luminosity function is sensitive to seed formation prescriptions.

Black hole merger rates vary significantly across models, affecting gravitational wave detection predictions.

Abstract

We study how statistical properties of supermassive black holes depend on the frequency and conditions for massive seed formation in cosmological simulations of structure formation. We develop a novel method to recalculate detailed growth histories and merger trees of black holes within the framework of the Illustris simulation for several seed formation models, including a physically motivated model where black hole seeds only form in progenitor galaxies that conform to the conditions for direct collapse black hole formation. While all seed models considered here are in a broad agreement with present observational constraints on black hole populations from optical, UV and X-ray studies, we find that they lead to widely different black hole number densities and halo occupation fractions which are currently observationally unconstrained. In terms of future electromagnetic spectrum observations, the faint-end quasar luminosity function and the low-mass-end black hole-host galaxy scaling relations are very sensitive to the specific massive seed prescription. Specifically, the direct collapse model exhibits a seeding efficiency which decreases rapidly with cosmic time and produces much fewer black holes in low mass galaxies, in contrast to the original Illustris simulation. We further find that the total black hole merger rate varies by more than one order of magnitude for different seed models, with the redshift evolution of the chirp mass changing as well. Supermassive black hole merger detections with LISA and International Pulsar Timing Array may hence provide the most direct means of constraining massive black hole seed formation in the early Universe.

Equations

$${\dot{M}}_{\mathrm{BH}}=\frac{4\pi \alpha G^2{M}_{\mathrm{BH}}^2\rho }{c_s^3}$$

$$\lambda <\frac{m_d^2{Q}_c{T}_{\mathrm{vir}}}{8j_d{T}_{\mathrm{gas}}},$$

$$M_{\mathrm{chirp}}=\frac{(M_1{M}_2{)}^{3/5}}{(M_1+M_2{)}^{1/5}}.$$

Full text

††pubyear: 2019

1 Introduction

It is widely understood that supermassive black holes (SMBHs) are found at the centre of massive galaxies (Kormendy & Richstone, 1995), and that properties of the host galaxy strongly correlate with black hole mass, suggesting a causal link between them (e.g. Magorrian et al., 1998; Gebhardt et al., 2000; Graham et al., 2001; Ferrarese, 2002; Tremaine et al., 2002; Häring & Rix, 2004; Gültekin et al., 2009; McConnell & Ma, 2013; Kormendy & Ho, 2013).

Observations have also shown that high-redshift quasars are found at $z\sim 6$ (e.g. Fan et al., 2006; Jiang et al., 2009) and even $z>7$, with black hole masses estimated at $\sim 10^{9}\mathrm{M_{\odot}}$ (Mortlock et al., 2011; Bañados et al., 2018). Reaching such high masses at such an early cosmic time provides a challenge to the theoretical models of supermassive black hole growth, and implies the possibility of the formation of sufficiently massive seeds followed by growth capable of reaching $10^{9}\mathrm{M_{\odot}}$ in less than $10^{9}$ yrs of cosmic time. Cosmological simulations have shown that extended periods of Eddington growth can be sustained by high density gas flows in very over-dense regions of the Universe, providing an almost constant fuel source for the accreting black hole, such that $\sim 10^{9}\mathrm{M_{\odot}}$ black holes are produced by $z\sim 7$ when starting with seeds of $\sim 10^{5}\mathrm{M_{\odot}}$ (e.g. Sijacki et al., 2009; Di Matteo et al., 2012; Costa et al., 2014; Curtis & Sijacki, 2016).

The initial formation of supermassive black holes seeds, however, is an open question, with many different mechanisms proposed which may lead to the formation of light (stellar), intermediate and massive seeds (see the seminal paper by Rees, 1984). Light seed formation has been proposed in the context of remnants of Population III (PopIII) stars (e.g. Madau & Rees, 2001; Volonteri et al., 2003). After the first generation of stars forms in pristine haloes at very high redshift, stars with masses from $\sim 25-140\mathrm{M_{\odot}}$ and $\gtrsim 260\mathrm{M_{\odot}}$ are expected to collapse to form a black hole (Heger et al., 2003; Heger & Woosley, 2010; Whalen & Fryer, 2012; Karlsson et al., 2013), producing light seeds with $M_{\rm{BH,seed}}$ on the order of half the stellar mass (though uncertainty in the PopIII initial mass function could produce $M_{\rm{BH,seed}}\sim 10-1000\>\mathrm{M_{\odot}}$, e.g. Fryer et al., 2001; Hirano et al., 2014). While seed formation from PopIII stellar remnants is potentially common (though with some uncertainty in the mass scales), the low-mass nature of these seeds would require substantially more growth and possibly super-Eddington accretion to reach $\sim 10^{9}\>\mathrm{M_{\odot}}$ by $z\sim 7$ making them potentially problematic as the progenitors of bright z $\gtrsim 6$ quasars.

At an intermediate seed mass scale, it is proposed that black holes may form in dense nuclear star clusters (NSC), where runaway stellar collisions can produce a very massive star which can then collapse to an intermediate mass black hole, on the order of $M_{\rm{BH,seed}}\sim 10^{3}\>\mathrm{M_{\odot}}$ (e.g. Begelman & Rees, 1978; Portegies Zwart & McMillan, 2002; Freitag et al., 2006b, a; Omukai et al., 2008; Devecchi & Volonteri, 2009; Katz et al., 2015). The higher initial mass compared to PopIII remnants makes the prospects of reaching bright quasar masses at $z\sim 6$ more viable but not without a challenge of sustained and very high growth rate.

Massive seed formation is generally proposed to occur in haloes with $T_{\rm{vir}}>10^{4}$ K containing (almost) pristine gas which can potentially collapse directly to a black hole without significant fragmentation. Such a collapse requires a massive gas cloud without significant rotation (which would collapse to a rotationally supported disk rather than a massive black hole), and without significant cooling from metals or molecular hydrogen, either of which would cause the cloud to fragment rather than collapse to a massive black hole. The gas content in primordial galaxies at sufficiently high redshift tends to be pristine, thus lacking the metals needed for metal-line cooling, and a nearby source of photo-dissociating Lyman-Werner radiation can disrupt and prevent formation of molecular hydrogen (see e.g. recent work by Regan et al., 2017). The net result is a massive cloud capable of collapsing directly to a black hole seed (referred to as a ‘Direct Collapse Black Hole’, or DCBH) with mass $M_{\rm{BH,seed}}\sim 10^{4}-10^{6}\>\mathrm{M_{\odot}}$ (e.g. Haehnelt & Rees, 1993; Loeb & Rasio, 1994a; Begelman et al., 2006; Regan & Haehnelt, 2009). Moreover, high baryonic streaming velocities with respect to the dark matter may promote the growth of DCBHs (Hirano et al., 2017; Schauer et al., 2017) with alternatively, or additionally sufficiently high gas turbulent velocities preventing excessive fragmentation as well. While DCBH seeds can be very massive (and thus easier to grow to the high masses observed at $z\sim 6-7$), the strict criteria for formation suggests they may be relatively rare to form (see e.g. Habouzit et al. 2016, Wise et al. 2019 and recent reviews by Latif et al. 2016 and Woods et al. 2019).

Given the uncertainty in formation mechanisms and resolution limitations, cosmological simulations often rely upon a simple prescription to insert a black hole seed of $\sim 10^{5}\>\mathrm{M_{\odot}}$ into sufficiently well-resolved haloes (e.g. Sijacki et al., 2009; Di Matteo et al., 2012; DeGraf et al., 2012; Dubois et al., 2014; Hirschmann et al., 2014; Nelson et al., 2015; Schaye et al., 2015; Feng et al., 2016). This method may be considered to broadly encompass each scenario described above: it is comparable to a recently formed DCBH, or could be a PopIII or NSC remnant which formed earlier and grew by several orders of magnitude before being added to the simulation. However the number densities of black hole and host halo occupation fractions could differ largely. Furthermore, the final masses of the largest black holes have been shown not to be strongly sensitive to the seed mass provided it is rather large ($\sim 10^{5}\mathrm{M_{\odot}}$ or above), and this model has been shown to well reproduce observed populations of black holes (e.g. Springel et al., 2005; Di Matteo et al., 2005; Sijacki et al., 2007; DeGraf et al., 2012; Hirschmann et al., 2014; Sijacki et al., 2015; DeGraf et al., 2015; Beckmann et al., 2017; DeGraf & Sijacki, 2017; Di Matteo et al., 2017; McAlpine et al., 2017; Weinberger et al., 2018). However, by construction these models are unable to directly investigate early black hole formation and growth. Furthermore, they treat each black hole’s seed equivalently, and so cannot address how sensitive individual black hole histories or time-dependent statistical properties of the whole population are on the frequency or method by which supermassive black hole seeds form. In particular, they typically rely upon two hand-selected mass scales: the threshold halo mass above which black hole seeds are inserted, and the mass of those seeds. These simple assumptions may not be realistic, and suggest that more advanced models for seed formation should be explored, and there have been several recent analyses investigating this, using both simulations and semi-analytic models (see, e.g. Habouzit et al., 2017; Ricarte & Natarajan, 2018; Wang et al., 2019). This is especially timely in view of the flood of next generation observational data, from, for example, JWST, WFIRST, eROSITA, Athena, LSST, IPTA, EPTA, NANOGrav, PPTA, and LISA which will push into the presently inaccessible regime of lower mass black holes and higher redshifts, together with entirely new constrains from gravitational wave observations. Hence, in this paper we address these questions by studying how changing the seed frequency and criteria for seed formation can impact individual black hole growth and the overall populations of black holes across cosmic time.

To avoid the prohibitive computational expense of repeatedly re-running a large cosmological simulation with many varying parameters for black hole seeding, we developed a post-processing method of re-calculating the growth of black holes, which we apply to the Illustris simulation (Nelson et al., 2015). Using the full output of black holes and their surrounding gas properties at every timestep, we are able to calculate the efficiency with which black holes should grow, even if the frequency and conditions for seed formation are changed. We apply this method to a variety of prescriptions for seed formation, including a less frequent stochastic seed model, a gas-spin model, and a physically motivated model based on the mass, metallicity, and gas spin of progenitor galaxies based on the DCBH mechanism. For each of these models, we study which galaxies we expect to form black holes, and how the growth history of each individual black hole from the original simulation is changed by each of these models. Using the full set of re-computed black hole histories, we characterize the impact on the black hole mass function and quasar luminosity function, late-time accretion efficiencies, black hole-host scaling relations, and black hole merger rates, quantifying how the frequency and mass ratios of black hole mergers are impacted by each seed model.

The outline of the paper is as follows. In Section 2 we discuss the Illustris simulation (Section 2.1) and associated construction of galaxy merger trees (Section 2.2), our post-processing method of calculating black hole evolution (Section 2.3), and the seeding models used in our analysis (Section 2.4). In Section 3 we discuss the results of our post-processing analysis. Section 3.1 shows how the seed model affects seed formation and the galaxies in which seeds form; Section 3.2 shows the impact on global black hole statistics over cosmic time; and Section 3.3 shows the impact on black hole mergers. Finally, we summarize our conclusions in Section 4.

2 Method

2.1 Illustris Simulation

In this work we analyse the Illustris111http://www.illustris-project.org; Nelson et al. (2015). suite of simulations with a cosmological box $106.5\>\rm{Mpc}$ on a side, performed with the moving mesh code AREPO (Springel, 2010). We only consider the highest resolution simulation incorporating the full physics model (for further details see Vogelsberger et al., 2014a; Genel et al., 2014; Sijacki et al., 2015) which has target gas cell mass $m_{\rm gas}=1.26\times 10^{6}\>\mathrm{M_{\odot}}$ and dark matter particle mass $m_{\rm DM}=6.26\times 10^{6}\>\mathrm{M_{\odot}}$. This run uses a standard $\Lambda$CDM cosmology, with $\Omega_{m,0}=0.2726$, $\Omega_{\Lambda,0}=0.7274$, $\Omega_{b,0}=0.0456$, $\sigma_{8}=0.809$, $n_{s}=0.963$, $H_{0}=70.4\,\rm{km}\,\rm{s}^{-1}\,\rm{Mpc}^{-1}$ (consistent with Hinshaw et al., 2013).

For the purpose of this study the modelling of black hole seeding, growth and feedback is most relevant, so we briefly summarize it here. Black holes are treated as collisionless sink particles. In Illustris, their seed mass is $10^{5}\>\mathrm{h^{-1}\>M_{\odot}}$ and they are seeded in all haloes with mass above $5\times 10^{10}\>\mathrm{h^{-1}\>M_{\odot}}$ if void of a black hole particle. This seeding model is loosely motivated by the direct collapse scenario (see e.g. Haehnelt & Rees, 1993; Loeb & Rasio, 1994b; Bromm & Loeb, 2003; Begelman et al., 2006; Regan & Haehnelt, 2009), but may remain broadly consistent with lighter seed formation models provided that their occupation fraction is similar and that they grow efficiently enough in the early Universe. After seeding, black holes may grow in mass either through gas accretion, which is modelled assuming a Bondi-Hoyle-like rate (Bondi & Hoyle, 1944; Bondi, 1952) capped at the Eddington limit, or by merging with other black holes which happen to be within each others’ smoothing lengths. In the absence of feedback such a black hole growth prescription would result in unrealistically large black hole masses. It is hence necessary to incorporate feedback prescriptions which in Illustris is accomplished by invoking three separate modes: “quasar”, “radio” and “radiative”, which depend on the accretion efficiency (for more details, see Vogelsberger et al., 2014b; Sijacki et al., 2015).

We note that while this black hole model leads to a wide range of black hole properties in good agreement with observational constraints (Sijacki et al., 2015), the seeding prescription is rather simplistic. Thus, the aim here is to explore different assumptions regarding seeding to pin down how much they can affect the properties of a representative cosmological black hole population, and in turn how this can inform us about the likely seeding scenario.

2.2 Constructing galaxy merger trees

The Illustris simulation identifies dark matter haloes using a friends-of-friends (FOF) halo-finder (Davis et al., 1985) with linking-length of 0.2 and baryonic matter assigned to the same group as the nearest dark matter particle. Following this, the SUBFIND algorithm (Springel et al., 2001; Nelson et al., 2015) is run to produce a catalog of gravitationally self-bound substructures within each halo, which we use to determine galaxy properties. Furthermore, to track galaxies over time (in particular to identify the progenitors of a given galaxy), we use the SubLink catalog produced by Rodriguez-Gomez et al. (2015). This catalog identifies descendant subhaloes across snapshots by summing a merit function (based on particle binding energies) for each particle common to a given pair of subhaloes in different snapshots. The descendant is then defined as the subhalo with the highest summed merit score. In this manner, we have access to a complete tree of all resolved subhaloes, such that we can identify galaxy mergers and the progenitors of any given galaxy, a key point when investigating if a given galaxy ever satisfied the conditions for SMBH formation at any earlier time.

2.3 Post-processing black hole growth

One of the challenges in investigating black hole seeding mechanisms with cosmological simulations is the large computational expense to run such simulations. Ideally, these simulations must have high resolution to model the regions in which black hole seeding occurs, large volume so as to include rare objects such as massive high-redshift quasars (see e.g. Mortlock et al., 2011), and be run to low-redshift to provide information on the SMBH population in the local Universe where observational constraints are strongest (both for predictive purposes, and also to confirm the accuracy of the simulation). Furthermore, the parameter space to investigate is large, spanning seeding environment (where and when black hole seeds form) and the seed black hole itself222Note that these considerations assume a fixed black hole accretion and feedback model, and that relaxing this would further enlarge the parameter space for exploration.. Hence rather than re-running the entire Illustris simulation for each variation in these parameters (at a cost of $\sim$16 million CPU hours, see Vogelsberger et al., 2014b), which would be numerically prohibitive, we take a post-processing approach in which we re-calculate the entire growth history of the full black hole population, using the high temporal resolution output from the original run, which contains local gas properties of each black hole at every timestep.

Black hole accretion within the simulation is given by

[TABLE]

(see Bondi & Hoyle, 1944; Bondi, 1952), with an imposed upper limit of the Eddington rate [$\dot{M}_{\rm{Edd}\ }=(4\pi GM_{\rm{BH}}m_{p})/(\epsilon_{r}\sigma_{T}c)$]. Thus the accretion rate depends upon black hole mass ($M_{\rm{BH}}$) and the density ($\rho$) and sound speed ($c_{s}$) of the nearby gas, with dimensionless $\alpha=100$, radiative efficiency $\epsilon_{r}=0.2$, gravitational constant $G$, speed of light $c$, proton mass $m_{p}$ and Thompson cross-section $\sigma_{T}$. Since the simulation saves each of these properties every time the black hole properties are calculated (with typical timesteps on the order of $2\times 10^{4}$ yrs), all the information necessary to re-calculate accretion rates is saved.

We perform this growth recalculation by taking the necessary properties ($M_{\rm{BH}}$, $\rho$, $c_{s}$) to calculate the expected accretion rate (Equation 1, with imposed Eddington limit). We repeat this calculation at each timestep to find the expected mass at the following timestep (i.e. $M_{\rm{BH},i+1}=M_{\rm{BH},i}+\dot{M}_{\rm{BH,i}}\Delta t$). In this way we are able to reproduce the full accretion history of the entire black hole population, but with the freedom to change certain parameters (e.g. seeding a subset of the black hole population or changing the initial seed mass). We also incorporate black hole mergers in our new calculation. Because we are using the gas properties local to black holes in the original simulation, our black holes are necessarily subsets of the original Illustris black holes (though often with different masses). Thus the mergers for our new black hole populations are also subsamples of the mergers occurring in the original simulation (and for each merger, the simulation saves the merger time and the ID and mass of each black hole). We determine mergers in our recalculated history by extracting mergers in the original simulation which involve only black holes seeded in our new model, at which point we combine the masses of the merging black hole pair (using our newly-calculated mass histories).

We note that black holes in the Illustris simulation merge as soon as they are within the black hole smoothing length of one another, without considering the hardening/merging timescale, which may delay when a given black hole pair should merge. In this analysis we use the merging times of the original simulation, but note that an extended hardening times may delay black hole mergers and thus impact our predicted black hole merger rate. However, an initial investigation finds that adding a delay of 1Gyr has minimal impact on the overall merger rates, except at the highest redshifts where the first black hole mergers occur, consistent with results from Salcido et al. (2016), who found minimal impact for delay timescales of 5Gyr in gas-poor galaxies.

We note that the primary limitation of our model is the inability to incorporate self-consistent feedback and self-regulation into post-processing calculations. In the original simulation, AGN feedback is deposited into the surrounding gas as the black hole grows; in our post-processing analysis we are able to modify the accretion history of the black hole, but the feedback energy deposited remains unchanged. At low black hole masses we expect this to have a minimal impact, as the black holes are generally too small to have a strong impact on the gas supply of their host halo. At high black hole mass, we know that quasar feedback plays an important role, and acts as a catalyst for self-regulation of further growth (e.g. Sijacki et al., 2007; DeGraf et al., 2012). Thus our recalculated growth histories may underestimate the final black hole mass, but this will only be relevant for black holes whose recalculated mass is substantially below the original mass at the time when feedback-induced self-regulation begins (we further discuss this issue and potential implications in Section 3.2).

We further note that the Illustris simulation output is incomplete for black holes below $z\sim 0.38$. This missing data necessarily prevents recalculating growth and mergers at lower redshifts, so we limit our analysis to $z\geq 0.5$. This limitation has a minimal impact on our results, however, as it is after the majority of black hole growth. We do not expect any qualitative changes at lower redshifts.

2.4 Seeding models

2.4.1 Illustris seeding

Within our framework, the simplest seeding prescription is to re-use the data from seeding black holes in the original Illustris simulation, i.e. when a halo crosses a specified threshold of $5\times 10^{10}\>\mathrm{h^{-1}\,M_{\odot}}$. By only seeding black holes that were originally seeded in Illustris, we explicitly have the local gas properties necessary for re-calculating growth, and the black hole mergers in this new case will simply match the black hole merger trees from Illustris. If using the same seed mass as Illustris ($10^{5}\>\mathrm{h^{-1}\,M_{\odot}}$), this will mean that each black hole should be ‘re-grown’ in a way that matches the data from the original simulation, which we have explicitly confirmed.

2.4.2 Stochastic seeding

An obvious variation on using the Illustris seed criterion is to simply seed a fixed fraction ($f_{\rm{seed}}$) of Illustris-seeded black holes, accomplished using a random number generator to determine whether a given Illustris seed should be considered ‘seedable’ in our re-calculation. As this model continues to only seed subset of black holes which have been seeded by Illustris, we again have the complete set of local gas properties. Furthermore, since all black holes are chosen from the original Illustris catalog, all black hole mergers can be extracted from the original Illustris merger trees.

We note that although we seed using a fixed probability, the total number of black holes is not the same fraction of original black holes (i.e. $N_{\rm{BH,f_{seed}}}\neq f_{\rm{seed}}\times N_{\rm{BH,Illustris}}$), due to black holes which involve multiple black hole seeds that eventually merge together. Furthermore, this is a mass-dependent effect, since high-$M_{\rm{BH,Illustris}}$ black holes tend to have undergone the largest number of mergers and thus have the largest number of chances to have at least one seed in our new calculation (e.g. for $f_{\rm{seed}}=0.1$, at z = 0.5 $18\%$ of Illustris black holes above $10^{6}\>\mathrm{M_{\odot}}$ have been reseeded in our new model, compared to $58\%$ of Illustris black holes above $10^{8}\>\mathrm{M_{\odot}}$. This is discussed further in Section 3.2 and especially Figure 8). This different impact on high- and low-mass black holes means that the impact on overall black hole populations, such as the black hole mass function and luminosity function as functions of cosmic time, is more complicated than a simple decrease in total black hole number (for further details, see Section 3.2).

2.4.3 Spin- and metallicity-based seeding

In Figure 1 we investigate the galaxies seeded by black holes in the Illustris simulation, showing the distribution in both spin-redshift space and metallicity-redshift space (both important criteria in direct-collapse formation models), where $\lambda$ and $Z$ are computed for the gas within the gas half-mass radius of the galaxy. This distribution shows that gas spin in the center of newly black hole-seeded galaxies is essentially redshift independent, and metallicity only shows a slight increase at low redshift (below $z\sim 1.5$). However, this is explicitly the case for haloes which have just recently crossed the $5\times 10^{10}\>\mathrm{h^{-1}\,M_{\odot}}$ threshold; we consider the impact on other galaxies in more detail in Section 3.1.

As a variant on the stochastic seeding model, here we seed a subset of Illustris black holes based on the dimensionless spin parameter of gas within the gas half-mass radius of the host galaxy at the time when the Illustris black hole is seeded (as low-spin may be needed for gas to collapse directly to a black hole; see Section 2.4.4 and Lodato & Natarajan, 2006; Natarajan & Volonteri, 2012, for more details). Specifically, we impose a cut on spin such that only black holes in galaxies with $\lambda<\lambda_{\rm{crit}}$ are seeded, where we select $\lambda_{\rm{crit}}$ to correspond to a given fraction of galaxies (i.e. to be comparable to a given value of $f_{\rm{seed}}$). Since $\lambda$ is approximately redshift-independent, seeding based on gas spin in this manner is qualitatively similar to the random seed selection model (note that the 50th-, 25th-, and 10th-percentiles are marked in Figure 1, corresponding to $f_{\rm{seed}}=0.5,0.25,0.1$ in Section 2.4.2). As a further test, we imposed an additional criterion based on gas metallicity (see also Section 2.4.4). The metallicity distribution is also roughly redshift-independent for $z>2$, so this additional criterion tends to only affect a small number of low-redshift, low-mass black holes, and otherwise produces qualitatively similar results.

2.4.4 Progenitor-based direct collapse seeding

The models described in Sections 2.4.2 and 2.4.3 are based solely upon instantaneous properties when a black hole particle was first inserted into the Illustris simulation. However, we note that the insertion of a sink particle in Illustris (as in other cosmological simulations) is not necessarily intended to represent the actual formation of a seed black hole, but rather a point at which we would expect a seed to have previously formed and grown to a “resolvable” mass. A more physically meaningful approach is to consider which black holes seeded in Illustris have a progenitor galaxy which satisfies the necessary criteria for the formation of a black hole seed, in this case by direct collapse. We are able to do this by considering the full history of all Illustris galaxies (and the merger trees involving them; see Section 2.2), checking for a set of conditions necessary for DCBHs where we impose three criteria:

• •

Galaxy mass: We impose that $M_{\rm{subhalo}}>3\times 10^{9}\>\mathrm{M_{\odot}}$ (where $M_{\rm{subhalo}}$ is the total mass contained within the SUBFIND-defined subhalo, as discussed in Section 2.2). This is largely a resolution-based criterion (corresponding to $\sim 500$ dark matter particles) such that we have sufficient gas particles to reasonably resolve the remaining criteria.

• •

Gas spin: Sufficiently high gas spin implies a more gravitationally stable structure which may not lead to the further collapse necessary to form a DCBH seed. A rough approximation for this criterion is

[TABLE]

where $m_{d}$ is the fraction of gas involved in infall, $Q_{c}$ is the Toomre stability parameter, $j_{d}$ is the angular momentum in the disc, $T_{\rm{vir}}$ is the virial temperature of the halo, and $T_{\rm{gas}}$ is the gas temperature (Lodato & Natarajan, 2006). We calculate $T_{\rm{vir}}$ and $T_{\rm{gas}}$ directly from the simulation, noting that $T_{\rm{gas}}$ is typically on the order of $1-3\times 10^{4}$K. The remaining parameters ($m_{d},j_{d}$, and $Q_{c}$) are much more sensitive to the structure of the disk which is less well-resolved in our simulation, especially in the low-mass galaxies at high-$z$ relevant for the DCBH formation. As such, we follow Natarajan & Volonteri (2012) and use $m_{d}=j_{d}=0.05$ and $Q_{c}=2$ rather than use the resolution-sensitive values from the simulation.

We show an example of these constraints using 1D histograms in Figure 1. The solid (dashed) black histogram shows the distribution of $\lambda$ ($\lambda_{\rm{max}}$) of galaxies newly-seeded in the original Illustris simulation, where there is a minor overlap (i.e. $\lambda<\lambda_{\rm{max}}$). For comparison, the green histograms show the distribution for $\lambda$ and $\lambda_{\rm{max}}$ for all galaxies above $3\times 10^{9}\mathrm{M_{\odot}}$ (our resolution-based mass threshold for DCBH seeding) at $z=4$. Here we see that both $\lambda$ and $\lambda_{\rm{max}}$ are higher than for the newly-seeded haloes in Illustris, but again there is a region of overlap where in principal DCBH formation can occur.

• •

Gas metallicity: In addition to rotational support, one of the main obstacles for massive DCBH seed formation is fragmentation of the collapsing gas cloud. Specifically, metals within the gas can facilitate cooling, which leads to fragmentation of the cloud prior to collapse into a single supermassive black hole seed. To account for this, we impose a metallicity threshold of $Z<10^{-5}Z_{\odot}$ for the galactic gas (consistent with Devecchi & Volonteri, 2009). We note that our results are not sensitive to the exact threshold, since Illustris galaxies rarely have non-zero metallicity near this value (e.g., at $z>5$, $<0.05\%$ of galaxies with non-zero metallicity have metallicity below $10^{-3}Z_{\odot}$). This is also shown in Figure 1, where the green 1D histogram shows the distribution of metallicities for galaxies above $3\times 10^{9}\>\mathrm{M_{\odot}}$ at $z=4$: the peak of non-zero metallicity galaxies occurs at $Z\sim 10^{-3.5}\>Z_{\odot}$, and nearly one-third of galaxies have zero metallicity.

We note however that this is sensitive to the metallicity enrichment model of the simulation. In the Illustris simulation, star particles are treated as single-age stellar populations which deposit metals into the neighboring gas cells based on the stellar evolution model, which includes supernovae (including core collapse, see Portinari et al. 1998 and Woosley & Weaver 1995, and Type Ia, see Thielemann et al. 2003) and AGB winds (see Karakas, 2010) [note that the simulation directly models nine distinct elements, though we only consider total metallicity in this analysis]. The metals are also distributed by stellar feedback via metal-loaded winds, and metal mixing is modeled purely by advection (for complete details of the models used, see Vogelsberger et al., 2013).

Our results may be sensitive to this model as a less efficient metal mixing model could lead to more haloes satisfying the metallicity criteria; conversely, a full incorporation of PopIII star formation may contribute more high-redshift metallicity, thus decreasing the haloes satisfying our metal criterion as there is no high redshift metallicity floor in Illustris.

These criteria are used to determine which Illustris black holes have progenitor galaxies with the necessary conditions for DCBH seed formation. Specifically, for the host subhalo of every newly seeded black hole in the Illustris simulation we check every progenitor galaxy in the SubLink galaxy merger tree (see Section 2.2 and Rodriguez-Gomez et al., 2015). Any Illustris black hole which has at least one progenitor galaxy simultaneously satisfying all three criteria is considered a ‘seedable’ black hole according to our DCBH model. We note here that our method of post-processing black hole growth relies on the output data provided for each Illustris black hole. As such, we are unable to directly calculate growth for black holes which do not exist in the original simulation. Thus our seeding model here does not seed black holes as soon as a galaxy satisfies the DCBH criteria, but rather when the original simulation seeds a black hole which had a progenitor galaxy which satisfies the criteria.

To illustrate our DCBH seeding method, in Figure 2 we show the full galaxy progenitor trees of three $z=5$ galaxies with newly-seeded black holes in Illustris, colour-coded by gas spin and symbols distinguishing between galaxies with $M_{*}=0$ and $M_{*}>0$ (vertical lines represent individual galaxy evolution with time; horizontal lines represent galaxy mergers). The first is a representative moderate-sized galaxy tree with no seedable progenitors; the second is a very large merger tree with a single seedable progenitor galaxy (circled in red), and the third shows a representative small merger tree containing a single seedable progenitor galaxy (circled in red).

As seen in Figure 2, the conditions for DCBH seed formation typically occur significantly earlier than the Illustris seeding. As such, although we apply our post-processing analysis to a given black hole with an initial seed mass of $10^{5}\rm{M}_{\odot}$ starting when Illustris seeds, there is a substantial period of time since actual seeding (when DCBH conditions are satisfied) during which the seed could grow. For high-redshift seeds, black holes tend to have several hundred Myr in which to grow, introducing two additional parameters to characterize this early growth: the seed mass when the initial DCBH forms, and the typical growth rate from that time until Illustris seeds the BH. Among $z>5$ seeds, we find that if the initial seed mass is $10^{4}\rm{M}_{\odot}$, 96% of DCBHs would be able to reach $10^{5}\rm{M}_{\odot}$ by the time they are included in Illustris if limited to sub-Eddington rates, and 83% would be able to do so even if initially formed at only $10^{3}\rm{M}_{\odot}$. As such, we find that the seed mass used in Illustris is fully consistent with a wide range of lower formation masses, followed by a reasonable period of sub-Eddington accretion. Recent work by Wang et al. (2019) investigated this early growth directly, using GADGET-3 simulations which seed at $M_{\rm{BH,seed}}=10^{3}\>\mathrm{M_{\odot}}$ based on local gas properties, and found that the mean mass for black holes in haloes at the threshold for seeding in Illustris is $\sim 10^{5.26}\>h^{-1}\>\mathrm{M_{\odot}}$, confirming that a seed mass of $\sim 10^{5}\>\mathrm{M_{\odot}}$ is reasonable, and that we can expect our DCBH seeds to reach this mass by the time the Illustris black holes are included.

We further note that an additional constraint for DCBH formation is a nearby source of photo-dissociating Lyman Werner radiation, capable of dissociating molecular hydrogen and thus prevent fragmentation resulting from molecular cooling. However, the Lyman-Werner intensity threshold needed remains very poorly constrained (see, e.g. Agarwal et al., 2016, who find $J_{\rm{crit}}$ ranges from $5\times 10^{-22}$ to $10^{-18}$ erg s*-1* cm*-2* sr*-1* Hz*-1*). Furthermore, recently Wise et al. (2019) proposed that lower Lyman-Werner fluxes can still be sufficient, with the primary driver for DCBH formation being high gas inflow rates onto the halo. Hence estimating for the incoming Lyman-Werner flux in the Illustris simulation would be highly uncertain, given the resolution limitations and the lack of PopIII stars. Given these caveats, we do not include any Lyman-Werner constraint in our analysis here. Note that we expect our seed model to over-estimate the DCBH formation rate compared to a simulation which fully incorporates this additional criterion, and thus may be considered an upper limit on the formation efficiency.

This DCBH seeding model is the most physically meaningful of all the models we consider, since we only take into account black holes in galaxies which have, at some point in the past, satisfied the necessary criteria for DCBH formation. Because our post-processing growth calculations rely upon the local gas properties which are only saved around Illustris black holes, this model by construction only includes black holes which are also seeded in the original Illustris simulation. As such, we again have the complete history of any mergers involving any combination of black holes seeded according to this new model. Note that this approach ignores galaxies which satisfy the criteria for DCBH seed formation but which never have a black hole form in Illustris. However, the only cases in which this can occur are galaxies which satisfy the above criteria but whose host haloes do not cross the $5\times 10^{10}\>\mathrm{h^{-1}\,M_{\odot}}$ mass threshold for seeding. These black holes would likely not grow into the regime in which we consider our results well-resolved ($M_{\rm{BH}}>10^{6}\>\mathrm{M_{\odot}}$), and so will not largely impact our results.

We note that we are only able to check progenitor galaxy conditions in the saved snapshots, which sets the time-resolution for this model. This timescale ($\sim 40$ Myr for $z\sim 5-10$) is longer than the typical timescale for DCBH formation, such that it is possible to miss a galaxy which satisfies the criteria during a short time period between snapshots. The typical variation of $\lambda$ is relatively slow compared to the snapshot timescales (see the generally gradual colour variation in individual vertical lines in Figure 2), so we would expect this to be a sub-dominant effect. The metallicity evolution can be fairly rapid, however, suggesting more frequent snapshots might show additional galaxies which satisfy the conditions for DCBH formation. By sub-sampling the available snapshots and extrapolating to shorter inter-snapshot timescales, if we assume that the collapse time is on the order of $\sim 5$ Myr, it is possible we might underestimate the seed frequency by up to 50%. We expect this is a conservative upper limit based on significant extrapolation (and the unresolved metal enrichment); nonetheless, we would expect that even an increase at this upper limit should not qualitatively change the results presented here. A more precise treatment would require a full re-simulation which checks for seed conditions on-the-fly, which is beyond the scope of this analysis.

2.4.5 Individual black hole growth

In Figure 3 we illustrate our post-processing technique, showing the full history of a SMBH and all its progenitors using the four different seeding models described above. Each panel shows the mass history for the most massive black hole from the original Illustris simulation as well as all black holes it merged with, having re-calculated the history according to each seed model. Over-plotted in colour on each panel is the most massive progenitor history for each of the four models, determined by tracking the final black hole back in time, and at each black hole merger selecting the more massive progenitor black hole at the time of the merger (and as such the coloured lines are identical in all four panels). There are significant differences in the total number of seeded black holes and the masses such black holes reach at both early and late times.

The top panel shows the history for the $f_{\rm{seed}}=1$ model, i.e. the seeding conditions as in the original Illustris simulation. This seed model involves a large number of black hole mergers, with black holes being seeded continually throughout cosmic time. The second panel shows the $f_{\rm{seed}}=0.1$ model. By construction this model exhibits an equivalent seeding pattern as the original simulation, merely more rare. This shows that forming only a subset of the original seed population can have a significant impact on the overall black hole growth history, affecting not only the seed number but also how efficiently the most massive progenitor grows (comparing the dark green and blue curves), and even the final mass (in this case resulting in a black hole $\sim$0.5 dex less massive). In the third panel, we show seeding using a cut on the instantaneous gas spin and metallicity when Illustris first seeds black holes, as described in Section 2.4.3. Again we see a moderate impact on the final black hole mass, as well as a significant impact on the growth histories reaching that point.

Finally, in the bottom panel of Figure 3 we show the analogous history according to the DCBH seeding model. Again we see a moderate impact on the final $M_{\rm{BH}}$ (though weaker than seeding with $\sim 10\%$ efficiency). We furthermore note a dramatic impact on when DCBH seeds form: rather than forming continuously throughout the simulation, the seeding occurs preferentially at high redshifts ($z>4$) when gas clouds are more likely to consist of pristine gas capable of collapsing directly to a massive black hole without fragmenting. This example growth history demonstrates the impact the seed formation can have on black hole populations, which we investigate quantitatively throughout the remainder of the paper.

3 Results

3.1 Black hole seeds and their host galaxies

Having performed the full recalculated growth histories of black holes according to the models described in Section 2.4, we consider the galaxies in which seeds form in each model. In Figure 4 we show the distribution of dark matter mass ($M_{\rm{DM}}$), gas mass ($M_{\rm{gas}}$), stellar mass ($M_{*}$), and baryon fraction ($f_{\rm{baryon}}$) of galaxies with newly-seeded black holes as a function of redshift, for the original Illustris seeding model (right panels) and our DCBH-based seeding (left panels). The most obvious differences between the seeding models is that the original Illustris seeding tends to take place in galaxies with higher masses and at lower redshifts than our DCBH-based seeding. This is to be expected, since the idea behind Illustris seeding (and that of most cosmological simulations) is to insert a black hole roughly consistent with the black hole - host galaxy scaling relation possibly well after the actual seed should have formed. In contrast, our new DCBH model finds galaxies in which the conditions for DCBH formation are satisfied, which tend to occur at earlier times and in lower-mass galaxies (similar to previous analyses, e.g. Habouzit et al., 2017).

Furthermore, we note that the masses and baryon fractions of galaxies hosting Illustris seeds tend to be relatively redshift independent (though $M_{*}$ does increase slightly with time), since Illustris always seeds haloes as they cross a fixed halo mass threshold ($M_{\rm{halo}}>5\times 10^{10}\>\mathrm{h^{-1}M_{\odot}}$). In contrast, in our DCBH-based seeding we find that at low redshift (below $z\sim 3$), only galaxies with relatively low baryon fractions are seeded, as they are more likely to have the low-metallicity necessary to prevent excessive fragmentation of gas.

The top left-hand panel of Figure 5 shows the total number of galaxies and haloes above several mass thresholds as a function of redshift (galaxies - solid lines; haloes - dot-dashed lines), showing that the number of high-mass galaxies grows monotonically with time, and the number of low-mass galaxies peaks at moderate redshift (below which we find a slight decrease due to increasing number of galaxy mergers). In particular, we note that the dot-dashed blue line ($>5\times 10^{10}\>\mathrm{h^{-1}M_{\odot}}$) is the number of haloes above the mass threshold for black hole seeding in the original Illustris simulation, and thus corresponds roughly to the number of black holes seeded in Illustris (neglecting mergers). We compare these total halo numbers to the number of galaxies satisfying our criteria for DCBH seed formation at redshift $z$ (dotted line), and the cumulative number of galaxies satisfying the same seedability criterion at any earlier time (dashed line). This shows that the total number of seeded black holes by $z=0.5$ is comparable between the original Illustris simulation (blue line) and our DCBH-based model (dashed black line); however, where the original simulation has a monotonically increasing seeded number, our DCBH seed model peaks at $z\sim 8$ where the conditions for direct collapse are most common. Furthermore, we note that our total number of DCBH seeds is actually higher than the original simulation. This increase is due to the lower host galaxy mass (see also Figure 4), and thus represents a population of black holes in lower-mass galaxies than included in the Illustris black hole population. Given how little the scaling relation is affected, this would represent a correspondingly low-$M_{\rm{BH}}$ population and would thus have minimal impact on the population of resolved black holes we consider in this analysis (see also discussion in Section 3.2).

From this figure, we find a very large number of galaxies which could potentially form DCBH seeds. Based on our model criteria (see Section 2.4.4), we predict $\sim 10^{-2}\>\>\mathrm{cMpc^{-3}}$ possible DCBH seed candidates by $z\sim 6$, far above the number density of high-redshift quasars (e.g. Fan et al., 2004; Jiang et al., 2016). However, we note that the majority of these objects will not grow significantly. In fact, at $z\sim 5$, we find $\sim 9\times 10^{-4}\mathrm{cMpc^{-3}}$ possible DCBH seeds whose host halo will have reached the mass threshold for seeding in Illustris, and our regrowth calculations suggest only $\sim 10^{-5}\mathrm{cMpc^{-3}}$ will have reached $M_{\rm{BH}}=10^{7}\mathrm{M_{\odot}}$. This is still well above the observed quasar number density, so incorporating additional seed criteria such as a Lyman-Werner radiation threshold or a gas inflow cut (e.g. Wise et al., 2019) would be a worthwhile exercise, bearing in mind a possibility that only a fraction of DCBH seeds grows sufficiently to power high redshift quasars.

In the lower left-hand panel of Figure 5 we plot the fraction of galaxies at redshift $z$ (with $M_{\rm{gal}}>3\times 10^{9}\>\mathrm{M_{\odot}}$) which would have formed a DCBH seed at some earlier time, according to our DCBH model. We find that the seeded fraction reaches a peak of $\sim 25\%$ at $z\sim 8$, while at lower redshifts the fraction drops, suggesting only $\sim 10\%$ of local galaxies above $3\times 10^{9}\mathrm{M_{\odot}}$ should have a black hole formed through direct collapse. We also note a slight upturn in the seeded fraction below $z\sim 3$. However, this increase is not due to more direct collapse seeds; rather it is a result of the decrease in $N_{\rm{gal}}(M>3\times 10^{9}\>\mathrm{M_{\odot}})$ caused by low redshift galaxy mergers (seen in the red curve of the upper panel).

In the right-hand panel of Figure 5 we plot the fraction of original Illustris black holes which are seeded according to several alternate seeding models (i.e. $N_{\rm{BH,new}}/N_{\rm{BH,Illustris}}$). In particular, we consider seeding based on gas spin only (see Section 2.4.3), having selected a maximum gas spin value corresponding to approximately 50, 25, and 10 per cent of Illustris black hole seeds ($\lambda_{\rm{g},50}$, $\lambda_{\rm{g},25}$, and $\lambda_{\rm{g},10}$, respectively), and for our DCBH seed model. We also show the total number of black holes seeded in Illustris for reference. Consistent with Figure 1 which showed $\lambda$ had minimal z-dependence, seeding according to instantaneous galaxy spin alone is roughly redshift-independent, similar to the simpler $f_{\rm{seed}}$-based seeding (not plotted, but would correspond to horizontal lines by construction). In contrast to this, DCBH seeding grows with time until $z\sim 4$. Below this the fraction of seeded black holes decreases as few new DCBH seeds form and existing black holes proceed to merge together.

This peak at $z\sim 4$ is qualitatively different from the left-hand panel of Figure 5, which demonstrated that the peak in DCBH seeding occurs at $z\sim 8$. This difference is in the definition of what is being plotted: left-hand panel shows the fraction of galaxies which satisfy the criteria for seeding a black hole at redshift $z$, while the right-hand panel shows the fraction of Illustris black holes which have at least one progenitor which satisfied our seeding criteria at some earlier time. At lower-redshift, a given Illustris black hole has a larger number of progenitor galaxies than at higher redshift, and is thus more likely to have at least one progenitor satisfying the necessary criteria, hence the increase in seedable fraction. At lower redshifts ($z\lesssim 4$), however, most of the progenitor galaxies capable of forming a DCBH have already merged into a large enough halo for Illustris to insert a black hole; thus most black holes inserted into Illustris at at $z\lesssim 4$ tend not to have any DCBH-seedable progenitor galaxies making low-$z$ DCBH formation very rare.

In Figure 6 we show the distribution of black hole mass at $z=0.5$ as a function of redshift at which the black hole is seeded ($z_{\rm{seed}}$) for the original Illustris simulation (top) and our DCBH seed model (middle). We define $z_{\rm{seed}}$ as the seed redshift of the most massive progenitor within the given black hole merger tree, rather than the earliest seed redshift within the tree, noting that the most massive progenitor is often not the first to be seeded (e.g. see Figure 3), but is the most relevant quantity to consider here.

We find a direct correlation between $M_{\rm{BH}}(z=0.5)$ and $z_{\rm{seed}}$, consistent between the original Illustris simulation, our DCBH seeding, and our $f_{\rm{seed}}$-based models (not plotted). As previously discussed, there are few low-$z_{\rm{seed}}$ black holes in the DCBH model. This is further quantified in the bottom panel of Figure 6 where we plot the fraction of original Illustris black holes which should be seeded according to our DCBH model. This shows a rapid decline in the fraction of Illustris black holes with $z_{\rm{seed}}\lesssim 4$ which should form via our DCBH model.

3.2 Global black hole properties

Having considered the galaxies in which black holes should form, we use the full recalculated growth history of all black holes to investigate the impact the seeding model has on measurements of the full black hole population. In particular, in Figure 7 we show the black hole mass function (BHMF) and quasar luminosity function (QLF) from the original Illustris simulation, for our fractional seed models and for our DCBH-based seeding model. We also show observational data from Hopkins et al. (2007, black points) and Aird et al. (2015, orange triangles), and the sensitivity limits of the Athena (Nandra et al., 2013) and eROSITA (Merloni et al., 2012) surveys. For the Aird et al. (2015) data, the Athena and eROSITA limits, we convert from X-ray to bolometric luminosities using the SED of (Hopkins et al., 2007).

The BHMF for the fractional seed models shows generally expected behavior: as $f_{\rm{seed}}$ decreases, the BHMF decreases. At high redshift, this shift in BHMF is proportional to the decrease in $f_{\rm{seed}}$, as the lower $f_{\rm{seed}}$ decreases the number of seeded black holes without affecting the growth of those black holes which are seeded. At low redshift, however, black hole mergers become significant for high-mass black holes, and thus the post-seed growth behaviour depends on seeding probability as well. Since massive, low-redshift black holes tend to have a large number of progenitor black holes, most of the high-$M_{\rm{BH}}$ black holes are seeded regardless of $f_{\rm{seed}}$.

We show this explicitly in Figure 8, the upper panel of which shows the correlation between the mass of each black hole at $z=0.5$, $M_{\rm{BH}}(z=0.5)$, and the total number of progenitor seeds it had, $N_{\rm{prog}}$, while the lower panel shows the fraction of black holes which have at least 2, 4, or 10 progenitors (corresponding to $f_{\rm{seed}}$ = 0.5, 0.25, 0.1, respectively). This shows that $\left<N_{\rm{prog}}\right>$ grows with $M_{\rm{BH}}$, and the mass scale above which we would expect most black holes to have at least one progenitor seed would be $M_{\rm{BH}}(z=0.5)\sim 1.5\times 10^{7}\>\mathrm{M_{\odot}}$, $6\times 10^{7}\>\mathrm{M_{\odot}}$, and $4\times 10^{8}\>\mathrm{M_{\odot}}$ for $f_{\rm{seed}}=0.5,0.25$, and $0.1$, respectively. We further note that these scales are redshift-dependent, as higher redshift black holes tend to have fewer progenitors. In addition to the original Illustris data, we show the number of progenitor black holes in our DCBH seeding model using dashed lines. As expected, $N_{\rm{prog}}$ is lower in this model, as fewer black holes are seeded.

The decrease in $M_{\rm{BH}}$ among high-mass black holes is due to the decreased mass gained via mergers (since fewer progenitors were seeded). This is further compounded by a lower accretion rate, both because $M_{\rm{BH}}\propto M_{\rm{BH}}^{2}$ (Equation 1), and the higher accretion rate in the original simulation produced stronger feedback and thus lower gas density and higher temperature than we should have in the new seed model, which each contribute to a lower $\dot{\mathrm{M}}$. Furthermore, the most massive black holes require a relatively early seed in a host with sufficient gas density to fuel an extended phase of near-Eddington growth (e.g. the most massive progenitor curves in Figure 3). When the seed fraction is decreased, there are fewer chances and thus a generally delayed onset of Eddington growth, producing a lower black hole mass. The net effect is that the number of high-mass black holes tends to be roughly unchanged but with slightly lower mass (as seen in Figure 3), producing a slightly steeper slope in the BHMF. We further note that this also means that although the high-end mass function is lower in our DCBH model, the black hole occupation fraction is still high in massive galaxies (rapidly approaching 1 for total galaxy masses $\sim\rm{few}\times 10^{12}M_{\odot}$). Conversely, low mass black holes with single progenitors are less likely to be seeded, but those that do have unchanged masses. Thus the decrease in the low-end BHMF does represent a comparable decrease in the black hole occupation fraction among low-mass galaxies as has also been found to be the case when seeding at lower masses based on local gas properties (see Habouzit et al., 2017).

The BHMF for the DCBH model is roughly consistent with the $f_{\rm{seed}}\sim 0.25$ model, except for the low-redshift, low-$M_{\rm{BH}}$ end, where the DCBH seed model falls off drastically (e.g. below $M_{\rm{BH}}\sim 10^{7}\mathrm{M_{\odot}}$ at $z=2$) due to the difference in redshifts at which black hole seeds form. In particular, most seeds in the DCBH model form at an early time (see Figure 5), in contrast to $f_{\rm{seed}}$ based models which continue seeding black holes at all redshifts as is the case in the original Illustris model. Since the DCBH model seeds very few black holes at low redshift and $M_{\rm{BH}}(z=0.5)$ is directly correlated with $z_{\rm{seed}}$, nearly all black holes have had sufficient time to reach at least $\sim 10^{6.5}\>\mathrm{M_{\odot}}$ (given our accretion prescription), and thus the mass function decreases rapidly below this scale. Thus our model predicts very few low redshift, low mass black holes which formed via direct collapse route. However, we note that this does not necessarily mean that there will be a turnover in the overall observed BHMF. Our model only considers DCBHs, but we also expect black holes may form via NSCs or from PopIII remnants. Black holes formed through these mechanisms may contribute significantly to the mass function, providing the potential to partially or even completely overcome the turnover found in the DCBH population. Alternatively, or additionally, it is possible that some of the direct collapse seeds will not grow by accretion as efficiently as in our model which may lead to a less significant decrease in the low mass-end BHMF shape. Furthermore, our analysis does not include black holes which may have been seeded by direct-collapse but whose host halo never reaches the threshold for seeding in the original simulation ($5\times 10^{10}\>h^{-1}\>\mathrm{M_{\odot}}$). As discussed in Section 3.1, this will represent a low-$M_{\rm{BH}}$ population, and thus should have minimal impact on the results shown here. In particular, we note that the mean mass for black holes in haloes below $10^{11}\>\mathrm{M_{\odot}}$ in the original simulation is only $2.6\times 10^{5}\>\mathrm{M_{\odot}}$ at z = 0.5 and less than 3% have $M_{\rm{BH}}>10^{6}\>\mathrm{M_{\odot}}$, suggesting these black holes will have minimal impact on the results considered here, where we generally limit our analysis to $M_{\rm{BH}}>10^{6}\>\mathrm{M_{\odot}}$. This is supported by results of Wang et al. (2019), who found that seeding black holes at earlier times with $M_{\rm{BH,seed}}=10^{3}\>h^{-1}\>\mathrm{M_{\odot}}$ resulted in comparable mass black holes when reaching the Illustris-seed threshold (for haloes at $M_{\rm{DM}}=10^{10}\>h^{-1}\>\mathrm{M_{\odot}}$, $\langle log(M_{\rm{BH}})\rangle\sim 5.18$ with standard deviation $\sim 0.55$). Habouzit et al. (2017) also looked at low-$M_{\rm{BH,seed}}$ simulations, and found that the black holes would not grow to $10^{6}\>\mathrm{M_{\odot}}$ until the host galaxy reached a stellar mass of a few $\times 10^{8}\>\mathrm{M_{\odot}}$ (and a few $\times 10^{9}$ in simulations which used a delayed cooling model for SN feedback), broadly consistent with Illustris’ halo mass threshold.

Similar to the BHMF, the QLF tends to exhibit a simple normalization shift equivalent to $f_{\rm{seed}}$ at either high redshift or low-$L_{\rm{BH}}$, the regimes in which most black holes have undergone few or no mergers and thus modifying the seed prescription changes which black holes form but does not affect their subsequent growth. As with the high-mass BHMF, the low redshift, high-luminosity black holes tend to have undergone numerous mergers, and thus the number is generally unaffected by the seed model. However, we find that the low-redshift, high-end of the QLF is somewhat less affected by the seed prescription than the BHMF is, such that current observational data cannot differentiate between our seed formation models. This is broadly consistent with the semi-analytic model of Volonteri et al. (2008), who also found the mass function to be more sensitive to seed model than the luminosity function.

We note that we find good agreement with observational data across a wide range of luminosities and redshifts. At low redshift and high luminosities, all of our seed models reproduce the observed data well (though noting that the simulations are volume-limited, and thus do not produce the brightest quasars). At the lowest luminosities the rarer seed models under-predict the luminosity function, but additional seed models which provide lower-mass and lower-redshift seeds would be expected to fill this gap in. At moderate redshifts ($z\sim 2$) we find that the rarest seed models (including the DCBH model) match the faint end most effectively. Although the brightest objects are underestimated, we expect this to be at least partially due to the limited boxsize and somewhat overestimating the impact of feedback in our DCBH seed model (see the discussion of the scaling relation in Figure 10, below). At the highest redshifts ($z\sim 4$), all our models remain broadly consistent with observations, though this is largely due to the uncertainty in the observed high-redshift faint-end luminosity function. We note that the upcoming eROSITA data will not help differentiate between the different black hole seeding models (see redshift dependent eROSITA limits), which will only become possible (both at high and low redshifts) with the Wide Field Imager of the Athena mission.

In Figure 9 we show the total comoving number density of AGN binned by bolometric luminosity to investigate cosmic downsizing, comparing the original Illustris simulation (dashed lines) to the DCBH model (solid lines). We also show observational constraints in the form of shaded regions, showing the areas spanned by the data from (Aird et al., 2015) and a range of models fitted to observational data (the ‘Full’, ‘PLE’, ‘Scatter’, ‘Schechter’, and ‘LDDE’ fitting models from Hopkins et al., 2007). We find that Illustris overpredictss low-luminosity AGN ($10^{43}-10^{44}{\rm erg\,s^{-1}}$, blue) at high-redshift, whereas the DCBH model more closely matches the observed data. Similarly, at intermediate luminosities ($10^{44}-10^{45}{\rm erg\,s^{-1}}$, green) the DCBH model matches observations very closely, where the original simulation over-predicted the number density, especially at high-redshift. In the brightest sample ($10^{45}-10^{46}{\rm erg\,s^{-1}}$, red) the DCBH under-predicts the number density at all but the highest redshifts. However, we note that the brightest AGN are least well-constrained in our model, as they are the objects where the AGN feedback (which is unaffected by our DCBH model) plays the strongest role, which we can better understand by looking at the relation between black holes and their host galaxies. Furthermore, we emphasize that the comoving number densities from Hopkins et al. (2007) are based on analytic fits to the observed data. The observed data includes significant scatter (see Figure 7) and can vary significantly depending on the fitting model used, especially when extrapolating toward faint luminosities.

In Figure 10 we show the impact the DCBH seeding method has on the scaling relation between the stellar mass of a galaxy ($M_{\rm{b,*}}$; defined as the stellar mass within the stellar half-mass radius of the subhalo) and the corresponding $M_{\rm{BH}}$ (the mass of the central black hole within the galaxy) for the original Illustris simulation (contours) and our new DCBH model (colourmaps). At high redshift we see excellent agreement between the two, with the exception of some low-$M_{\rm{BH}}$ outliers in the DCBH model. These outliers are caused when the black hole which undergoes the fastest growth in the original simulation is not seeded. Instead there is a lower-mass black hole that grows up toward the main relation, but will not reach it until a lower redshift.

At lower redshift, we note that the DCBH model only fills the high-mass end of the scaling relation, without any low-mass objects (similar to Volonteri et al., 2008, who also found that stricter seed models result in a significant fraction of low-mass galaxies which do not host black holes). This lack of a low-mass end is due to the low-redshift, low-mass galaxies not satisfying the conditions for direct-collapse seed formation while those that are seeded did so at sufficiently high-redshift, providing enough time for both the black hole and the galaxy to grow to at least a moderate mass scale333We stress again this result may also be sensitive to the accretion model. If a less-efficient accretion model were adopted, the slower growth of black holes seeded at high-redshift would also lead to a larger population of low-mass black holes at low redshift.. We emphasize that this suggests that low-redshift, low-mass galaxies rarely host massive black holes formed via direct collapse; however they may still host SMBHs formed through other channels, e.g. from NSCs or PopIII stellar remnants. We also note that the threshold for seeding in these low-mass galaxies may be sensitive to the metal enrichment model in Illustris (see Section 2.4.4 and Vogelsberger et al., 2013). Thus the precise mass scale below which we do not find DCBHs may be model-dependent, but the qualitative result should hold.

In addition, we note that below $z\sim 2$ the high-mass end in the DCBH model lies below the relation from the original simulation. As discussed earlier in this section, black holes in the DCBH model often have lower masses than in the original simulation (see e.g. the mass function in Figure 7). However, we note that one caveat of our model is that the feedback from the original simulation remains unchanged, and thus the growth of massive black holes in the DCBH model may be underestimated, as the gas supply around black holes which reach the self-regulated regime in the original simulation has been impacted by a stronger feedback than it would be with the lower-mass DCBH. As such, while we do expect the DCBH model to have a lower high-$M_{*}$ end, our results provide a lower bound rather than a direct prediction.

We also consider the global black hole mass density as a function of redshift for a range of $M_{\rm{BH}}$-bins in the upper panel of Figure 11. As noted above, high-$M_{\rm{BH}}$ are likely somewhat under-massive due to overestimation of the effect of feedback. In this figure we hence compare the original Illustris mass densities (dashed lines) to a range of possible values in our DCBH model (shaded region). The shaded areas show the region bounded by two constraints: the explicit calculation of $M_{\rm{BH}}$ from our DCBH model, and an adjusted value to account for this overestimated feedback. For the adjusted mass, we compare the best-fitting relation for the scaling relation (see Figure 10) of the original simulation and the DCBH model, and at the high-mass end where the DCBH model is lower, we increase $M_{\rm{BH}}$ of the DCBH by the difference between the two relations. The lower-$M_{\rm{BH}}$ bound is thus expected to be a lower limit as the black hole growth has been underestimated. Similarly, the higher-$M_{\rm{BH}}$ bound is likely the upper limit, as we expect the DCBHs to have a somewhat lower mass than the original simulation due to the less frequent seeding. As such, the shaded regions represent a well-defined constraint for the mass densities. The width of the shaded regions is clearly mass dependent, with only the most massive black holes above $10^{9}\,{\rm M_{\odot}}$ being affected significantly by the overestimated feedback in our DCBH model.

In the lower panel of Figure 11 we show the ratio of the DCBH mass density to that of the original Illustris simulation. Here we see that the density of low-mass black holes ($<10^{7}\>\rm{M}_{\odot}$) tends to drop off below z$\sim$3-4, as new seeds rarely form and thus black holes tend to grow out of the lowest mass bin. For higher mass black holes, we see that at high-redshift we have much lower mass densities in the DCBH model, but as time passes the more efficient growth starts to increase the mass density up toward that of the original Illustris simulation.

From Figure 11 several interesting features can be deduced. First, consistent with our findings above, black hole mass density of the DCBH model is systematically lower with respect to the Illustris prediction. Second, by splitting the contribution of total black hole mass density into different mass bins, we see that at high redshift the density is dominated by low-mass black holes. Hence, the discrepancy with regard to current high-$z$ observational estimates for black hole mass density is likely due to underestimating this low-mass population, requiring either a slightly higher DCBH seed efficiency (but see our discussion above that our DCBH seeding model is likely rather conservatively high), or other channels of black hole formation contributing to the total mass density (e.g. SMBHs from PopIII or NSC seeds). Furthermore, at low redshift the overall density is dominated by high mass black holes, suggesting that, if the growth of SMBHs from low mass seeds is negligible, matching current observational densities would require more efficient growth than our post-processing analysis finds, possibly up to including some super-Eddington phases. It will be interesting to see if future JWST, Athena, IPTA, and LISA observations will sufficiently inform us about these growth channels to break these degeneracies in theoretical models (e.g., see Sesana et al., 2007; Plowman et al., 2010; Sesana et al., 2011; Pacucci et al., 2015; Klein et al., 2016; Natarajan et al., 2017; Valiante et al., 2018).

3.3 Black hole merger rates

In addition to the impact on the overall black hole population and its cosmic growth, changing the seeding prescription strongly modifies the rate of black hole mergers, as well as the mass ratio of the merging black holes.

Naively one might expect that the number of mergers $N_{\rm{merger}}$ scales roughly with $f_{\rm{seed}}^{2}$, if $f_{\rm{seed}}$ were the probability of any given black hole having formed. However, this ignores the full history of the merger trees: more massive black holes tend to have undergone a large number of mergers (see Figure 8 and the discussion in Section 3.2), and thus have a much higher probability of at least one of those progenitors having been seeded. For this reason, mergers involving more massive black holes will be less sensitive to the seeding prescription with the number of mergers in the DCBH model becoming significantly closer to the Illustris prediction for higher mass black holes.

In Figure 12 we plot the 2-dimensional distribution of the number density of mergers, binned by both $M_{1}$ and $M_{2}$ (the masses of the more massive and less massive black holes involved in the merger, respectively). For fixed $f_{\rm{seed}}$ models, we find that the merger frequency peaks at $M_{1}\sim 10^{6.3}\>\rm{M}_{\rm\odot}$ and $M_{2}\sim 10^{6}\>\rm{M}_{\rm\odot}$. The overall $M_{1}-M_{2}$ distribution appears roughly comparable (barring an overall normalization shift), but we do note that the lower-mass mergers are affected more strongly by a lower $f_{\rm{seed}}$, so lower $f_{\rm{seed}}$ results in a distribution weighted slightly more toward high-mass mergers.

The DCBH based seeding model, however, produces a qualitatively different merger-mass distribution which peaks at $M_{1}\sim 10^{7.75}\>\rm{M}_{\odot}$, $M_{2}\sim 10^{7}\>\rm{M}_{\rm\odot}$ and generally has $M_{2}$ within $\sim 1$ dex of $M_{1}$. We also see this quite clearly in the 1D distributions of $M_{1}$ and $M_{2}$ in Figure 12, which show that mergers in the DCBH model peak at both higher $M_{1}$ and $M_{2}$ (see also Figure 13 and 14).

We also note that an earlier analysis of the EAGLE simulation found merger rates dominated by low-mass black holes, but much more strongly peaked toward seed mass black holes. They also found that decreasing the black hole seed mass to $M_{\rm{BH,seed}}\sim 10^{4}\>\mathrm{M_{\odot}}$ (but maintaining the halo threshold for seed formation) dramatically affected the merger distribution, eliminating nearly all mergers except those with $M_{2}\sim M_{\rm{BH,seed}}$ and a small population of $M_{1},M_{2}\sim 10^{8}\>\mathrm{M_{\odot}}$, reinforcing the importance of improved seed models when looking at merger rates (Salcido et al., 2016).

In Figure 13 we plot 1-dimensional histograms of $M_{2}$ for various bins of $M_{1}$. There is some slight $M_{2}$-dependence, with high-$M_{2}$ values tending to be less sensitive to $f_{\rm{seed}}$ due to the increase in progenitor black holes. Other than this small effect, the distribution is qualitatively similar for all $f_{\rm{seed}}$ models beyond an overall normalization shift. However, as in Figure 12, the DCBH-based seeding model shows a dramatically different distribution, dropping off rapidly at low $M_{2}$. This dropoff is consistent with the dropoff in mass function seen in Figure 7, which shows dramatically fewer low-mass black holes in the DCBH model at low redshifts. Because the DCBH seeding model has most black holes form at relatively high redshift, at lower redshift there are far fewer low-$M_{\rm{BH}}$ black holes available to undergo mergers, and thus we find fewer low-$M_{2}$ mergers. This also means that the low-$M_{2}$ dropoff is partially a redshift-dependent factor, and so high redshift mergers exhibit a less extreme dropoff (see the left- vs. right- panels of Figure 13).

In Figure 14 we show the number density of mergers for original Illustris and our DCBH model as a function of primary mass ($M_{1}$, the more massive black hole involved in the merger; left panel) and secondary mass ($M_{2}$, the less massive black hole; center panel). In addition to a significant decrease in merger number, we again see that the DCBH model tends to have much higher secondary black hole masses, while the primary masses tend to be relatively similar (though without reaching as high- or low- masses). We note that the Illustris simulation predicts that $M_{1}$ should peak near $\sim 10^{7}\>\mathrm{M_{\odot}}$, broadly consistent with the findings from the EAGLE simulation, though EAGLE produced more mergers with $M_{2}\sim M_{\rm{BH,seed}}$ and relatively fewer intermediate-$M_{2}$ mergers (Salcido et al., 2016). In addition, we show (right panel) the merger number as a function of chirp mass, a measure of the effective mass of a binary system affecting the gravitational wave signal of the merger, defined as

[TABLE]

We see that the DCBH model peaks at a higher $M_{\rm{chirp}}$ ($\sim 10^{7}-10^{8}\>\rm{M}_{\rm\odot}$, compared to $\sim 10^{6}-10^{7}\>\rm{M}_{\rm\odot}$ from the original Illustris simulation), and does not extend to as low values (ending just below $10^{6}\>\rm{M}_{\rm\odot}$), suggesting DCBH-DCBH mergers are most common at slightly higher masses then the original Illustris simulation predicts. It is worth emphasizing that the merger density of our DCBH model remains consistent with predicted constraints from Pulsar Timing Array (e.g. Middleton et al., 2016), though the current PTA constraints remain very broad, hence future tighter constraints would help differentiating the models.

In Figure 15 we look at the redshift distribution of mergers, showing the cumulative number of mergers before redshift $z$, binned by $M_{1}$ (left panels) and $M_{2}$ (right panels), for $f_{\rm{seed}}=1$ (top), $f_{\rm{seed}}=0.1$ (middle), and the DCBH model (bottom). Again we find the $f_{\rm{seed}}$ models are qualitatively similar to one another, but the DCBH seed model is significantly different. Mergers with low-$M_{1}$ in the DCBH seed model tend to occur exclusively at high redshifts. Most DCBH mergers with $M_{1}<10^{7}\>\rm{M}_{\rm\odot}$ occur by $z\sim 3.5$, whereas in $f_{\rm{seed}}$-based models $N_{\rm{mergers}}$ grows with time throughout the simulation. Similarly, there are few $M_{2}<10^{6}\>\rm{M}_{\rm\odot}$ DCBH mergers below $z\sim 3.5$, though $10^{6}\>\rm{M}_{\rm\odot}<M_{2}<10^{7}\>\rm{M}_{\rm\odot}$ do continue to lower redshifts.

In Figure 16 we show the relation between chirp mass and merger redshift. In the original Illustris simulation (orange contours), we see that the typical $M_{\rm{chirp}}$ remains roughly constant with time. In these models, black holes are seeded all the way to $z=0$, so at all redshifts we find mergers spanning the full range of masses (see also Figure 15, which shows comparable increases in $N_{\rm{mergers}}$ for all $M_{\rm{BH}}$). In contrast, however, the DCBH mergers (blue contours) show a strong increase in $M_{\rm{chirp}}$ with decreasing redshift, as expected from our earlier analysis: at low-redshift, low mass black holes tend to be much rarer (Figure 7), resulting in far fewer low-mass mergers (the flattening of the $N_{\rm{mergers}}$ for low-$M_{\rm{BH}}$ in Figure 15). Thus the redshift distribution of $M_{\rm{chirp}}$ detections can provide a strong means of discriminating between seed models.

Finally, in Figure 17 we consider the rate at which mergers occur throughout the observable universe, both at a given chirp mass ($dN/dlogMdt$, black curves) and the cumulative number above a given chirp mass ($dN/dt$, red curves). Here we see that the merger rate peaks for a chirp mass of $\sim 3-4\times 10^{6}\>\rm{M}_{\rm\odot}$ (original Illustris) and $\sim 4\times 10^{6}\>\rm{M}_{\rm\odot}$ (DCBH model). For the Illustris simulation, this peak is comparable to the peak seen in Figure 14, as we would expect; although the observable volume is redshift-dependent, $M_{\rm{chirp}}$ is approximately independent of redshift (see Figure 16), so the merger rate is comparable to the merger density. However, for the DCBH model, $M_{\rm{chirp}}$ evolves strongly with redshift, shifting the merger rate toward a lower chirp mass. From the cumulative rate of mergers above a given chirp mass (red curves), we see that the DCBH model predicts $\sim 1\>dex$ fewer mergers, such that we would only expect to detect $\sim 0.1$ supermassive black hole merger per year. However, we note that this rate is only for mergers involving black holes which are well resolved by the simulation, so detectors sensitive to lower-mass ranges can be expected to make more detections than our simulation predicts.

For example, we note that LISA (Amaro-Seoane et al., 2017) has primary sensitivity which extends down to $\sim 10^{4}\>\rm{M}_{\rm\odot}$, well below our limit of $\sim 10^{6}\>\rm{M}_{\rm\odot}$. Given the steepness of the halo mass function (and the associated increase of low mass galaxies, see Figure 5), we note that the actual number of detectable mergers will likely be higher as a result of the lower-mass mergers unresolved by our simulation This is consistent with previous analyses (e.g. Sesana et al., 2004, 2007; Tanaka & Haiman, 2009; Sesana et al., 2011; Ricarte & Natarajan, 2018), which have found that the distribution of merger rates peaks at low black hole masses, tending toward the black hole seed mass used in the model. However, at higher masses where we resolve well the black holes, we expect our estimates to be much closer to complete.

We also note that this predicted rate assumes that black holes form only via direct collapse, and the incorporation of black holes formed via alternative seed pathways (e.g. PopIII remnants or in NSCs) would be expected to significantly increase the overall detection rate. In particular, a recent paper by Dayal et al. (2019) used a semi-analytic model to predict high redshift merger rates, and found that for $z>5$ mergers, LISA detections of mergers between a DCBH and a stellar mass seed would outnumber the DCBH-DCBH merger detections by two orders of magnitude or more. Furthermore, similar to Ricarte & Natarajan (2018), Dayal et al. (2019) found that the merger rate would peak at very low masses, with LISA-detectable mergers (based on SNR $>7$) peaking at $M_{\rm{BH}}(1+z)\sim 10^{4}-10^{5}\>\rm{M}_{\rm\odot}$, and Bellovary et al. (2019) found that black hole mergers in dwarf galaxies (below the masses considered in this work) will also be detectable by LISA suggesting that the predictions we make for the rate of high-mass DCBH-DCBH mergers represents only a small subset of all LISA detections.

As such, our predicted merger rate can be treated as a lower-limit for detections. However, this is dependent on the treatment of merging timescales. In particular, the black hole model in Illustris repositions black holes toward the potential minimum of its host halo to prevent N-body scattering with dark matter or stellar particles, which we note could result in underestimating the migration time for a secondary black hole brought into a halo by a galaxy merger. Furthermore, the simulation merges black holes as soon as they are within a smoothing length of one another (see Section 2.3 and Sijacki et al., 2015), although the actual hardening timescale may be quite significant (reaching up to $\sim 10^{9}$ yr for major mergers with $M\sim 10^{9}\>\rm{M}_{\rm\odot}$, depending on orbital parameters; see, e.g. Kelley et al., 2017).

Adding long hardening timescales would substantially delay the merger of any given black hole pair, decreasing the rate of high-redshift mergers which may be detected. However, it would be expected to produce a corresponding increase in late-time merger rates, and the lower-redshift at which the merger occurs would make them more easily detectable. Thus we expect that a more realistic hardening timescale would not substantially impact the detection rate of black hole mergers, only that it might lead to a distribution weighted more toward low-redshift (and we further note that an initial investigation suggests that incorporating a delay of $\sim$1 Gyr would not significantly impact the overall merger rates, consistent with findings of Salcido et al. 2016).

4 Conclusions

In this work, we have developed a novel method to investigate the impact black hole seeding models have on the cosmic evolution of the entire black hole population. Once a given seed model is assumed, our method recalculates entire black hole growth and merging histories taking advantage of the Illustris cosmological simulation of galaxy formation. Hence our method presents a more self-consistent and comprehensive approach than the traditional semi-analytic models while at the same time allowing us to test a large number of seeding models which would be computationally prohibitively expensive with full hydrodynamical simulations of structure formation. Recall the caveat to this approach is that the impact of AGN feedback on the gas supply around the black hole is unchanged during the post-processing approach; thus for black holes with lower mass compared to the original simulation our model may over-estimate the effect of AGN feedback and thus decrease the predicted accretion rate. However, this effect will generally be significant only for a small fraction of black holes, at the highest mass end where feedback is strongest (for details, see Section 3.2 and Figure 11).

In our analysis we test stochastic and physically-motivated seed models based on the direct collapse black hole (DCBH) formation scenario, find the impact on low- and high-redshift black hole populations and merger rates, and determine which observables provide a means for distinguishing between simulated models for seed formation. From this, we draw the following conclusions:

• •

Our novel tool produces post-processed black hole growth histories based on the simulation data available at each timestep. These histories match the original simulation when using the same seeding parameters, but alternate conditions for seeding can also be easily incorporated at a negligible fraction of the original simulation’s computational expense. In particular, we implement simple stochastic seeding models (our $f_{\rm{seed}}$-based models) and a physically motivated seeding model which uses galaxy merger trees to find progenitor galaxies in atomically cooling haloes with low metallicity and low central gas spin in which may be the forming sites of black hole seeds via direct collapse (our DCBH model).

• •

Galaxy properties of newly-seeded black holes in the original Illustris simulation are largely redshift independent (though with some increase in metallicity below $z\sim 2$), and seeding continues throughout cosmic time. In contrast, our DCBH seed model is biased toward higher redshifts. Furthermore, the low redshift seeds (below $z\sim 3.5$) tend to be found in galaxies with relatively low baryon fractions, which have a tendency for low-metallicity while still having enough gas to form DCBH seeds.

• •

About $10-20\%$ of all simulated galaxies with total mass above $3\times 10^{9}\>\rm{M}_{\odot}$ have a progenitor which satisfies all our criteria for direct collapse formation, suggesting only $\sim 20\%$ of Illustris black holes may have formed in this manner.

• •

All of our seed models (down to $f_{\rm{seed}}=0.1$) produce quasar luminosity functions (QLFs) consistent with current observational constraints. Except at very high redshift (where the observed QLF is poorly constrained), the bright end of the QLF ($L_{\rm{BH}}>10^{43.5}$ erg s*-1*) is similar for all seed models. The faint-end QLF shows significant dependence on seeding, however, particularly at low redshift, suggesting that upcoming surveys such as Athena should be able to put constraints on the seeding model. We note, however, that this is based on a single-formation mechanism via direct collapse. Additional formation pathways, such as PopIII remnants and NSCs, can provide additional black hole populations, particularly at lower masses.

• •

The scaling relation between black hole and host galaxy masses is strongly affected by the seeding model assumed, with the DCBH seed model resulting in very few low black hole mass objects at low redshift, in contrast to the original Illustris simulation. Hence, future observational campaigns that aim to probe the low mass end of the scaling relations may provide much needed insight into the likely nature of black hole seeding.

• •

Total black hole mass density is sensitive to the seeding model adopted, and our results imply that either DCBH seed need to undergo more efficient growth to match current observational constraints (possibly including a phase of super-Eddington accretion) or that the other seed formation channels lead to non-neglible cosmic black hole mass growth. Future JWST, Athena, LISA, and PTA observations may be able to shed light on this degeneracy (see also Sesana et al., 2007; Plowman et al., 2010; Sesana et al., 2011; Klein et al., 2016, for how LISA could differentiate between light and heavy seed models).

• •

Black hole - black hole mergers are dramatically affected by the seeding model. The total number of mergers scales more strongly than linearly with $f_{\rm{seed}}$, with the DCBH model expecting approximately one order of magnitude fewer mergers than the original Illustris simulation. Furthermore, the expected merger frequency is mass-dependent, with mergers involving the lowest masses depending most strongly on the seeding frequency. Thus high mass-ratio mergers are most reduced within our DCBH seed model: black holes tend to reach the highest masses only at low-redshift, and the lack of low redshift seeding means there are very few low-mass DCBHs at late times with which the large black hole can merge.

• •

The redshift evolution of merger chirp mass is strongly dependent on the seed formation model as well. The continued seeding across all redshifts in the original simulation results in a roughly redshift-independent distribution of chirp masses. In contrast, the lack of low-redshift seed formation in our DCBH model results in a typical chirp mass that strongly increases with time, reaching $\sim 10^{8}\,\rm{M}_{\rm\odot}$ at low-$z$, compared to $\sim 10^{6.5}\,M_{\rm\odot}$ in the original Illustris simulation.

In summary, this work demonstrates that seed formation mechanisms can have significant impact on black hole populations across all mass and redshift ranges. We have shown that differing seed models, including a more physically-motivated mechanism for seeding black holes in cosmological simulations based on the direct collapse scenario, are capable of broadly matching current observations across all redshifts. However, different seeding models have distinct predictions for higher redshift and lower mass black holes, which upcoming surveys such as JWST and Athena may be able to constrain. Perhaps most significantly, we have shown that seeding models have a strong effect on black hole merger rates, particularly at low redshift and high mass ratios, such that upcoming gravitational wave detections [e.g. LISA (Amaro-Seoane et al., 2017), International Pulsar Timing Array (Verbiest et al., 2016), European Pulsar Timing Array (Desvignes et al., 2016), Parkes Pulsar Timing Array (Reardon et al., 2016), and NANOGrav (Arzoumanian et al., 2018)] may be able to differentiate between seed formation pathways for supermassive black holes.

Acknowledgments

We would like to thank Martin Haehnelt, Alberto Sesana, Marta Volonteri and Priyamvada Natarajan for useful discussions and comments on this work and the referee for a helpful report. CD and DS acknowledge ERC starting grant 638707 and support from the STFC. This research used: The Cambridge Service for Data Driven Discovery (CSD3), part of which is operated by the University of Cambridge Research Computing on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The DiRAC component of CSD3 was funded by BEIS capital funding via STFC capital grants ST/P002307/1 and ST/R002452/1 and STFC operations grant ST/R00689X/1. DiRAC is part of the National e- Infrastructure. Simulations were run on the Harvard Odyssey and CfA/ITC clusters, the Ranger and Stampede supercomputers at the Texas Advanced Computing Center as part of XSEDE, the Kraken supercomputer at Oak Ridge National Laboratory as part of XSEDE, the CURIE supercomputer at CEA/France as part of PRACE project RA0844, and the SuperMUC computer at the Leibniz Computing Centre, Germany, as part of project pr85je.