The Mass Function of Supermassive Black Holes in the Direct-collapse Scenario

TL;DR

This paper models the mass distribution of supermassive black holes formed via direct collapse, revealing a broken power-law shape influenced by growth rates and accretion duration, linking early black hole formation to observed quasar luminosities.

Contribution

It introduces a theoretical model describing the SMBH mass function as a broken power-law based on formation and accretion rates in the direct-collapse scenario.

Findings

Mass function follows a broken power-law with a characteristic break.

The intermediate mass slope depends on the ratio of formation to growth rates.

The high-mass decline is related to the duration of the growth period.

Abstract

One of the ideas to explain the existence of supermassive black holes (SMBH) that are in place by z~7 is that there was an earlier phase of very rapid accretion onto direct collapse black holes (DCBH) that started their lives with masses ~ 10^4-10^5 M_solar. Working in this scenario, we show that the mass function of SMBH after such a limited time period with growing formation rate paired with super-Eddington accretion can be described as a broken power-law with two characteristic features. There is a power-law at intermediate masses whose index is the dimensionless ratio {\alpha} = {\lambda}/{\gamma}, where {\lambda} is the growth rate of the number of DCBH during their formation era, and {\gamma} is the growth rate of DCBH masses by super-Eddington accretion during the DCBH growth era. A second feature is a break in the power law profile at high masses, above which the mass function declines rapidly. The location of the break is related to the dimensionless number \b{eta} = {\gamma} T, where T is the duration of the period of DCBH growth. If the SMBH continue to grow at later times at an Eddington-limited accretion rate, then the observed quasar luminosity function can be directly related to the tapered power-law function derived in this paper.

Equations

$$L_{\mathrm{E}}=\frac{4\pi cG{m}_{\mathrm{p}}M}{{\sigma }_T}.$$

$$\frac{dM}{dt}={\gamma }_{0}M\Rightarrow M(t)=M_0\exp ({\gamma }_{0}t),$$

$$\frac{dn}{d\ln M_0}=\frac{1}{\sqrt{2\pi }{\sigma }_0}\exp \left(-\frac{(\ln M_0-{\mu }_0{)}^2}{2{\sigma }_0^2}\right).$$

$$\ln M(t)=\ln M_0+\gamma t,$$

$$\frac{dn}{d\ln M}=\frac{1}{\sqrt{2\pi }{\sigma }_0}\exp \left(-\frac{(\ln M-{\mu }_0-\gamma t{)}^2}{2{\sigma }_0^2}\right).$$

$${\int }_0^T\frac{1}{\sqrt{2\pi }{\sigma }_0}\exp \left(-\frac{(\ln M-{\mu }_0-\gamma t^{\prime }{)}^2}{2{\sigma }_0^2}\right)f(t^{\prime })d{t}^{\prime }.$$

$$f(t^{\prime })=\frac{\lambda \exp (-\lambda t^{\prime })}{[1-\exp (-\lambda T)]}.$$

$$\int \exp [-(a{x}^2+bx+c)]dx$$

$$f(\ln M)$$

$$\alpha$$

$$\beta$$

Full text

accretion — black hole physics — galaxies: high-redshift — galaxies: luminosity function — quasars: general — quasars: supermassive black holes

1 Introduction

A key challenge to the theory of the formation of supermassive black holes (SMBH) in the early universe is the observation of very massive ($M\approx 10^{9}M_{\odot}$) and luminous ($L\gtrsim 10^{13}\,L_{\odot}$) quasars already in place by $z\sim 7$, when the universe is just $\sim 800$ Myr old (e.g., Fan et al., 2006; Mortlock et al., 2011; Wu et al., 2015; Bañados et al., 2018, ; see also the review by Woods et al. (2018)). Objects that accumulate at least a billion $M_{\odot}$ in less than a billion years after the Big Bang put a strain on the normal ideas of Eddington-limited growth of black hole seeds that originate from Population III stellar remnants. Starting from a seed mass $M_{0}$, Eddington-limited growth leads to a mass $M(t)=M_{0}\exp[\epsilon^{-1}(1-\epsilon)t/t_{\rm E}]$, where $t_{\rm E}\approx 450$ Myr and $\epsilon$ ($\approx 0.1$) is a radiative efficiency factor. Population III stars are thought to have masses $\lesssim 40M_{\odot}$ (Hosokawa et al., 2011), and their remnants would be less massive, so that there is apparently not enough time available to reach $M\sim 10^{9}M_{\odot}$. These constraints show that a combination of both more massive initial seeds and a super-Eddington growth rate may be necessary to account for the observed SMBH at $z\sim 7$.

One promising pathway is that of direct collapse black holes (DCBH) (Bromm & Loeb, 2003). The idea is that Lyman-Werner (LW) photons (having energies 11.2 to 13.6 eV) from the first Population III stars can propagate far from their sources and dissociate H2 in other primordial gas clouds. Without H2 cooling these gas clouds equilibrate to temperatures $T\sim 8000$ K set by atomic cooling, which means that the Jeans mass is $\sim 10^{5}M_{\odot}$ at a number density $n=10^{4}$ cm*-3*, as opposed to $\sim 10^{3}M_{\odot}$ in a normal Population III star formation environment or $\sim 1M_{\odot}$ in present day star formation. Due to their large mass these collapsing cloud fragments may be able to collapse directly into black holes, after a brief period as a supermassive star (Bromm & Loeb, 2003) or quasi-star (Begelman et al., 2006, 2008), if the infalling matter can overcome the angular momentum barrier and disruptive effects of radiative feedback. An interesting joint solution to these barriers is proposed by Sakurai et al. (2016) based on the episodic accretion scenario of Vorobyov et al. (2013) that is powered by gravitational instability in a circumstellar disk. In this model the episodic accretion results in a lower surface temperature of a supermassive star, thereby also reducing the effect of radiative feedback that can limit mass accumulation in the case of normal Population III star formation (Hosokawa et al., 2011). The DCBH model has been extensively developed in the context of galaxy formation models, resulting in a scenario where the formation of atomic cooling halos is seeded by the first stars, and the subsequent DCBH produce LW radiation that triggers the formation of other atomic cooling halos and DCBH in a kind of chain reaction process (Yue et al., 2014). A rapid period of growth of atomic cooling halos and therefore DCBH formation ensues, with the growth rate at any time related to the instantaneous number of DCBH. The rapidly growing phase of DCBH creation is also a period of possible rapid mass growth through super-Eddington accretion (Inayoshi & Haiman, 2016; Pacucci et al., 2017). The whole process comes to a rapid halt however, when the gas in the atomic cooling halos is photoevaporated by the ambient radiation field. According to Yue et al. (2014) the DCBH era lasts from $z\approx 20$ to $z\approx 13$, or a time period $T\approx 150$ Myr, after which DCBH formation is completely suppressed. Here we adopt the picture emerging from their semi-analytic model, however we note that numerical simulations (e.g., Agarwal et al., 2012; Chon et al., 2016; Habouzit et al., 2016) have not reached a consensus on the DCBH formation rate or the termination time of their formation.

In this Letter, we seek a simple model of the growth of DCBH in the early universe that captures just the essential features of the scenario described above, in order to reach an analytic understanding of the mass and luminosity functions of observable quasars that form through the DCBH scenario.

2 Background

Yue et al. (2014) estimate that the rapid formation of DCBH occurs between $z\approx 20$ and $z\approx 13$, after which it ends abruptly. The semi-analytic model of Dijkstra et al. (2014) is broadly consistent with this and shows that the number density of DCBH $n_{\rm DCBH}$ grows rapidly during a similar interval.

Figure 1 shows the results of Dijkstra et al. (2014) plotted against time after the Big Bang, for two possible values of the critical flux of LW photons $J_{\rm crit}$ (written in units of $10^{-21}$ erg s*-1* cm*-2* Hz*-1* sr*-1*) that is needed to create the atomic cooling halos that are DCBH progenitors. The calculation of $J_{\rm crit}$ (e.g., Shang et al., 2010) includes the additional effect of near-infrared radiation with energy above 0.76 eV that can also inhibit the formation of H2 by dissociating the H*-* ion that is an intermediary in the H2 formation chain. Although significant uncertainties exist in the actual number of DCBH created, due to uncertainties in the value of $J_{\rm crit}$ that depends on the source density and spectra, escape fraction from the host halos, etc., the slope of growth, $\lambda(t)\equiv d\ln n_{\rm DCBH}/dt$, is similar in all their modeled cases. For the two cases shown here, we find an average value $\lambda$ measured between data points from $z=20$ to $z=12$, which corresponds roughly to the period of rapid DCBH formation. For their canonical model $J_{\rm crit}=300$ the best fit line yields $\lambda=27.7$ Gyr*-1*, and for $J_{\rm crit}=100$ the best fit is $\lambda=17.7$ Gyr*-1*. Here we use the canonical model and adopt $\lambda=28.0$ Gyr*-1*, slightly steeper than the best fit. This in the interest of rounding off and also because there is evidence that the DCBH growth era ended just prior to this data point, at $z\approx 13$ (Yue et al., 2014), which would tend to drive the slope to a slightly greater value.

Each DCBH can be modeled as growing in mass at an exponential rate, but the starting times of the accretion process will be spread throughout the DCBH formation era. However, the super-Eddington growth will cease for all DCBH at about the same time, so that there will be a distribution of accretion times among the population of DCBH. This is a key part of our model as developed in Section 3.

The growth of an individual DCBH is thought to proceed by default at an Eddington-limited rate, but periods of super-Eddington growth are also possible (Pacucci et al., 2017). The Eddington luminosity is

[TABLE]

where $M$ is the black hole mass, $m_{\rm p}$ is the proton mass, and $\sigma_{T}=(8\pi/3)(e^{2}/m_{\rm e}c^{2})^{2}$ is the Thomson cross section in which $m_{\rm e}$ is the electron mass. At this luminosity the radiation pressure can balance the gravitational pressure. The accretion of mass to very small radii comparable to the Schwarzschild radius will release a significant portion of the rest mass energy, hence the luminosity is normally estimated as $L_{\rm acc}=\epsilon\dot{M}_{\rm acc}c^{2}$, where $\dot{M}_{\rm acc}$ is the mass accretion rate and $\epsilon$ is the radiative efficiency, typically set to 0.1. Since the accretor will gain rest mass at the rate $\dot{M}=(1-\epsilon)\dot{M}_{\rm acc}$, we equate $L_{\rm acc}$ with $L_{\rm E}$ to find that

[TABLE]

where $\gamma_{0}=(1-\epsilon)/(\epsilon t_{\rm E})$ and $t_{\rm E}=2e^{4}/(3Gm_{\rm p}m_{\rm e}^{2}c^{3})=450$ Myr is the Eddington time. Here we follow Pacucci et al. (2017) in accounting for the idea that accretion (especially of the super-Eddington variety) may be episodic, by identifying the duty cycle $\mathcal{D}$ as the fraction of time spent accreting, and the Eddington ratio $f_{\rm Edd}$ that is $=1$ for Eddington-limited accretion but can be $<1$ for sub-Eddington accretion and $>1$ for super-Eddington accretion. We use a generalized accretion rate $\gamma=\chi\gamma_{0}$ where $\chi=\mathcal{D}f_{\rm Edd}$ is a correction factor to account for the fact that the accretion rate could be super-Eddington for some periods of time. The quantitatively relevant parameter is $\chi=\mathcal{D}f_{\rm Edd}$. Pacucci et al. (2017) find that objects with $M\gtrsim 10^{5}M_{\odot}$ can have high efficiency accretion, $0.5\leq\mathcal{D}\leq 1$ and $1\leq f_{\rm Edd}\leq 100$, but objects with $M\lesssim 10^{5}M_{\odot}$ have low efficiency accretion, $0\leq\mathcal{D}\leq 0.5$ and $0\leq f_{\rm Edd}\leq 1$. The simplest assumption is that $\chi=1$ for Eddington-limited growth, but our model allows for the putative super-Eddington growth in the DCBH formation era.

3 Mass Function

We assume that the distribution of initial black hole masses is lognormal, i.e., the differential number density per logarithmic mass bin is distributed normally:

[TABLE]

Here $\mu_{0}$ and $\sigma_{0}$ are the mean and standard deviation of the distribution of $\ln M_{0}$, respectively. A lognormal distribution for the birth mass function of DCBH seeds is consistent with the results of Ferrara et al. (2014) for intermediate masses ($4.75<\log(M/M_{\odot})<6.25$), and we fit those results with $\mu_{0}=11.7$ (corresponding to a peak at $\log M/M_{\odot}=5.1$) and $\sigma_{0}=1.0$.

Since the growth law implies that

[TABLE]

we can write the mass function at a later time as

[TABLE]

The accretion time $t$ may not be a fixed constant that applies to all objects, therefore we can integrate over a function $f(t)$ (which has units of inverse time) that describes the distribution of accretion times. In this case the final observed mass function $f(\ln M)\equiv dn/d\ln M$ is

[TABLE]

Here $f(t^{\prime})$ is a normalized distribution of accretion times $t^{\prime}$ and $T$ is the maximum possible accretion time. The function $f(t^{\prime})$ is determined by considering the creation rate of black holes in the DCBH scenario. The number density $n$ of black holes grows in a type of chain reaction (Yue et al., 2014; Dijkstra et al., 2014) with the instantaneous creation rate $dn/dt=\lambda(t)\,n$. The simplest case where $\lambda(t)=\lambda$ has a constant value leads to pure exponential growth. If this growth continues from a time $t=0$ when the first DCBH is created until a time $T$ when the creation of all DCBH is terminated, then each black hole that was created at time $t$ has an accretion lifetime $t^{\prime}=T-t$ in the range $[0,T]$. The normalized distribution of accretion lifetimes $t^{\prime}$ is then

[TABLE]

Using the indefinite integral identity

[TABLE]

valid for $a>0$, we evaluate the integral in Equation (6) using Equation (7) and obtain a full expression

[TABLE]

Here $\alpha\equiv\lambda/\gamma$, the dimensionless ratio of the growth rate of DCBH formation to the growth rate of the mass of individual DCBH, and $\beta\equiv\gamma\,T$, the dimensionless number of DCBH growth times within the DCBH formation era.

In the limit $T\rightarrow\infty$, the function becomes

[TABLE]

which is the modified lognormal power law (MLP) function (Basu et al., 2015). Equation (9) represents a tapered version of the MLP, with the break in the power-law occurring at $\log M\approx(\mu_{0}+\beta)/2.3$, meaning that the peak of the original lognormal is shifted in $\ln M$ by an amount $\beta=\gamma T$. Henceforth, we refer to Equation (9) as the tapered power law (TPL) function.

Figure 2 shows the TPL function using parameter values obtained from models of the DCBH growth era. We pick $\mu_{0}=11.7$ and $\sigma_{0}=1.0$ based on the model of Ferrara et al. (2014). From Dijkstra et al. (2014) (see Figure 1) we adopt $\lambda=28.0$ Gyr*-1* for the era of rapid DCBH formation using their canonical model. The length of the DCBH growth era is $T=0.15$ Gyr (Yue et al., 2014). For accretion growth during this period we expect that super-Eddington growth can occur (Dijkstra et al., 2014) but with a wide range of possible values. We pick a series of values $\chi=[1,2,3]$ covering Eddington-limited growth and two values of super-Eddington growth. Since $\gamma=\chi\gamma_{0}=20\chi$ Gyr*-1*, this leads to $\alpha=[1.4,0.7,0.47]$ and $\beta=[3,6,9]$ for our adopted values of $\lambda$ and $T$. Figure 2 shows that the super-Eddington growth models allow for the development of a mass function that has both a visually evident power-law profile as well as a notable break in the power law at high mass. This break is a marker of the end of the DCBH growth era, since both the creation of new DCBH as well as their super-Eddington growth ceases after the time interval $T$.

4 The Quasar Luminosity Function

Once the DCBH growth era has ended at $z\approx 13$, the population of DCBH may continue to undergo Eddington-limited accretion, and the luminosity function can be estimated using Equation (1). Over time, the mass function $f(\ln M)$ will retain its shape but move to the right since $\ln M$ at the end of the DCBH era will shift by an amount $\gamma\Delta t_{z}$ where $\Delta t_{z}$ is the time interval between the end of the DCBH era ($z\approx 13$) and an observable redshift $z$. However, the duty cycle $\mathcal{D}$ and therefore $\chi=\mathcal{D}f_{\rm Edd}$ may be $\ll 1$ in this era, rendering mass growth to small fractional levels. A random sampling of $\mathcal{D}$ and $f_{\rm Edd}$ for individual object growth after $z\approx 13$ shows that the overall distribution maintains its shape and moves to the right in $\log M$. We also note that the mass growth of SMBH may be quenched above $\sim 10^{10}\,M_{\odot}$ (Inayoshi & Haiman, 2016; Ichikawa & Inayoshi, 2017), in agreement with results of current quasar surveys (Ghisellini et al., 2010; Trakhtenbrot, 2014).

Assuming that observed quasars are undergoing Eddington-limited accretion, we use Equation (1) to transform the mass function into a luminosity function. We expect that the mass of the quasars are not growing substantially during this time, for reasons discussed above, and we are really most interested in fitting the shape of the function, which should remain much the same for a variety of redshifts in the post-DCBH-growth era. Figure 3 shows the inferred quasar luminosity function (QLF) $\phi(L)\propto dn/d\log L$ for a suitable pair of values for $\alpha$ and $\beta$, overlaid on $z=3$ quasar bolometric luminosities compiled by Hopkins et al. (2007). Here we are only interested in fitting the shape of the luminosity function and not the absolute number of sources. The normalization can be scaled to fit the observed number at any redshift.

We hold $(\mu_{0},\sigma_{0})$ fixed at their model-inspired values $(11.7,1.0)$ since they are most important in determining the unobserved low end of the luminosity function. We effectively fit observations with two parameters $(\alpha,\beta)$. This is in contrast to the usual practice of fitting the observed QLF with a double power law (e.g., Hopkins et al., 2007; Masters et al., 2012; Schindler et al., 2019) that requires three parameters: the two power-law indices and a joining point.

In our model, the values of $\alpha,\beta$ that fit the QLF are not merely mathematical parameters. Instead, they reveal the history of the putative DCBH growth era. For the adopted DCBH number growth rate $\lambda=28.0$ Gyr*-1* and a duration $T=0.15$ Gyr, and individual masses growing at a rate $\gamma=\chi\gamma_{0}$, the two fitted parameters are related to the super-Eddington factor $\chi$ by

[TABLE]

We find an excellent fit to the QLF with $[\alpha,\beta]=[0.5,8.4]$. Both the values of $\alpha$ and $\beta$ imply a super-Eddington factor $\chi=2.8$, revealing the self-consistency of our model. In principle the QLF could have been fit with any $\alpha$ and $\beta$ that could individually imply very different values of $\chi$. In that case our underlying model would be inconsistent, or at least need to explore values of $\lambda$ and $T$ that were quite different than those implied by current models of the DCBH growth era.

To elaborate on the above point, our model could in principle also be applied to SMBH formation from alternate scenarios such as Population III remnants (Madau & Rees, 2001; Whalen & Fryer, 2012) or mergers of primordial Population III stars (Boekholt et al., 2018; Reinoso et al., 2018). It could apply as long as the black hole production could be described as growing exponentially at some rate $\lambda$ and for a finite time $T$, during which the individual masses grew at an Eddington-limited or super-Eddington rate.

5 Summary

We have presented an analytic model that captures some essential features of the DCBH growth scenario and uses them to derive an analytic mass function and by implication a luminosity function. A double power-law function has been commonly used in the literature to mathematically fit the quasar luminosity function (QLF). Here, we instead use a physically-motivated formula based on the scenario of the DCBH growth era that has been developed by many researchers. It is not a double-power law at high mass and luminosity, but rather a tapered power-law. We believe that the rapid fall off in the QLF at high luminosity is better modeled as a tapered part of a power law rather than as a second power law. The break point of the power law identifies the end of the era of DCBH creation.

We have fit an observed QLF with a power-law index $\alpha=0.5$ and the break-point related parameter $\beta=8.4$. These are consistent with a period of rapid mass growth of DCBH with super-Eddington factor $\chi=2.8$, for a time period $T=150$ Myr during which the growth rate of the number density $n_{\rm DCBH}$ was $\lambda=28.0$ Gyr*-1*. In principle, the best fit values to QLF data can be used to constrain such theoretical models of DCBH growth.

Our model has two key components. Initially high mass $\sim 10^{5}\,M_{\odot}$ seeds grow rapidly in number during a limited time period in the early universe, since DCBH formation leads to the emission of LW photons that seed the formation of other DCBH. These objects also live within gas-rich halos and undergo super-Eddington mass accretion. Then at some time both the formation of DCBH as well as the super-Eddington accretion of the existing DCBH comes to a rapid halt due to the photoevaporation of the host halos. What remains is a tapered power-law (TPL) distribution of masses and therefore also of luminosity if the observed quasars are undergoing subsequent Eddington-limited accretion. Future modeling can relax some of these assumptions, for example the formation of DCBH may continue long enough to outlive the period of rapid (super-Eddington) mass growth, especially if driven by mechanisms other than the LW flux (Wise et al., 2019), and the super-Eddington accretion may not apply to all objects (Pacucci et al., 2017; Latif et al., 2018).

Acknowledgements

We thank the referee for constructive comments. SB was supported by a Discovery Grant from NSERC.