Constraints on Cosmological Coupling from the Accretion History of Supermassive Black Holes

Coupling of black hole mass to the cosmic expansion has been suggested as a possible path to understanding the dark energy content of the Universe. We test this hypothesis by comparing the supermassive black hole (SMBH) mass density at z = 0 to the total mass accreted in active galactic nuclei (AGN) since z = 6, to constrain how much of the SMBH mass density can arise from cosmologically coupled growth, as opposed to growth by accretion. Using an estimate of the local SMBH mass density of ≈1.0 × 106 M ⊙ Mpc−1, a radiative accretion efficiency, η, in the range 0.05 < η < 0.3, and the observed AGN luminosity density at z ≈ 4, we constrain the value of the coupling constant between the scale size of the Universe and the black hole mass, k, to lie in the range 0 < k ≲ 2, below the value of k = 3 needed for black holes to be the source term for dark energy. Initial estimates of the gravitational-wave background (GWB) using pulsar timing arrays, however, favor a higher SMBH mass density at z = 0. We show that if we adopt such a mass density at z = 0 of ≈7.4 × 106 M ⊙ Mpc−1, this makes k = 3 viable even for low radiative efficiencies, and may exclude nonzero cosmological coupling. We conclude that, although current estimates of the SMBH mass density based on the black hole mass–bulge mass relation probably exclude k = 3, the possibility remains open that, if the GWB is due to SMBH mergers, k > 2 is preferred.


INTRODUCTION 1.The Soltan Argument
There exists an elegant argument for constraining the growth of supermassive black hole (SMBH) mass in galaxies via accretion.Soltan (1982) showed that the integral of the AGN luminosity density over cosmic history gives an estimate of the total accreted mass onto SMBHs, which can be compared to the local mass density of SMBHs to constrain how they grow.This relation, between two independently derived quantities, can give excellent constraints on both AGN physics and the nature of SMBHs themselves.For example, using then available data, Soltan (1982) showed that the history of accreted mass is close to the local SMBH mass density, if a radiative efficiency of η ≈ 0.1 is assumed.This range in η is in good agreement with other estimates.
Theoretical arguments give a range of 0.05 ≲ η ≲ 0.43, depending on the spin of the black hole (BH) (Bardeen 1970;Thorne 1974).This range is itself in good agreement with most other studies, which find values in the range 0.1 ≲ η ≲ 0.3 (e.g.Yu & Tremaine 2002;Elvis et al. 2002;Hopkins et al. 2007;Shankar et al. 2009;Lacy et al. 2015;Zhang & Lu 2017;Shankar et al. 2020;Shen et al. 2020;Farrah et al. 2022).This consistency argues for the increase in SMBH mass density with cosmic time to be driven by accretion.

Cosmologically-coupled black holes
For several decades it has been recognized that BHs with non-singular interiors (including 'dark energy'-like interiors) may have exteriors identical to BHs with interior singularities (e.g.Gliner 1966;Dymnikova 1992;Faraoni & Jacques 2007;Beltracchi & Gondolo 2019;Cadoni et al. 2023a).Recently, however, it was realized that the masses of non-singular BHs could be affected by the cosmic expansion (Croker & Weiner 2019).Indeed, in a recent paper, Cadoni et al. (2023b) show that such Lacy, et al. a coupling is necessary in General Relativity (GR) if BHs do not contain singularities.Thus the detection of a non-zero cosmological coupling would have important implications for our understanding of BHs.
The mass coupling between BHs and the cosmic expansion can be parameterized via: (1) (Croker et al. 2020(Croker et al. , 2021) ) in which M • (a) and M • (a i ) are the BH mass at some later and earlier times respectively, a and a i are the scale factors at those times (a = (1 + z) −1 at redshift z), and k is the cosmological coupling strength, with −3 ≤ k ≤ 3 from a causality constraint.Croker & Weiner (2019) show that the value of k is related to the equation of state of the medium internal to the BH.Expressing the equation of state as p = wρ, where p is the pressure and ρ the density, Croker et al. (2021) show that k = −3w.Thus, for vacuum energy where p = −ρ, k = 3.In this case, the mass density of BHs is conserved as the Universe expands: their number density decreases with a −3 but their mass increases as a 3 and so BHs could account for the cosmological constant (Croker & Weiner 2019;Farrah et al. 2023a) (see also Prescod-Weinstein et al. 2009, for an alternative way BHs could be the origin of dark energy).In contrast, GR solutions for non-singular BHs prefer k = 1 (Cadoni et al. 2023b).In this case, the BH mass growth would not compensate for cosmological dilution of their density, so they would not contribute to dark energy.
Observationally, the evidence is mixed.Farrah et al. (2023b,a) find, based on an analysis of SMBH versus host galaxy masses out to z ≈ 2, that a value of k ≈ 3 is suggested and k = 0 can be ruled out at ≈ 3σ.On the other hand, Lei et al. (2023) argue, based on a much smaller sample of AGN at z > 4.5 with much lower stellar masses, that k = 3 can be ruled out at ≈ 2σ.For stellar mass black holes, Gao & Li (2023) argue that cosmological coupling of black holes can explain the distribution of compact object masses in low mass X-ray binaries, in particular the "mass gap" between the highest mass neutron stars and lowest mass BHs.However, Rodriguez (2023a) disfavor coupling based on observations of a BH binary in the globular cluster NGC 3201 and Andrae & El-Badry (2023) argue that two BH binaries found by Gaia are likely to have formed with BH masses that were too small if k = 3. 1 Also, Ghodla et al. (2023) argue that k = 3 would result in too many BH mergers compared to observations from gravitational wave detectors.

Scope of this paper
As the observational picture seems far from settled, we decided to investigate an independent constraint on BH growth.For cosmologically-coupled mass growth of BHs to be viable, the mass increase must fit within the constraint provided by the Soltan argument -that is, the local SMBH mass density must accommodate both the mass increase due to accretion, and the mass increase due to cosmological coupling.In this paper, we examine whether cosmologically-coupled BH growth can be accommodated within the constraints from the Soltan argument, and, if so, what constraints this can set on the value of k.We adopt the latest measures of both the AGN luminosity function (LF) over cosmic history, and of the local SMBH mass density, and explore values for the cosmological coupling strength in the range 0 < k < 3. We discuss uncertainties on our results arising from uncertainties in both η and the Eddington ratio λ. §2 presents our methods and §3 our sources of data.§4 presents our results and §5 discusses these results in the context of uncertainties on observed parameters.§6 summarizes our conclusions.We assume a cosmology with H 0 = 70 kms −1 Mpc, Ω M = 0.3 and Ω Λ = 0.7.

Generalization of the Soltan Argument for cosmologically-coupled black holes
The SMBH mass and luminosity are connected for each AGN through where L is the bolometric luminosity of the AGN, η is the radiative efficiency, and Ṁacc is the mass accretion rate of the BH.Mass/energy that is not radiated is accreted onto the BH, so the rate of mass increase of the BH, Ṁ In the absence of cosmological coupling (and ignoring the mass density of any high redshift seed population, which we assume to be negligible throughout this paper, the mass density in SMBHs today is equal to the integral of the AGN luminosity density, U (z) over Cosmic Time, modulated by η: where U : where and H 0 is the Hubble Constant at the current epoch.If, in addition, BH mass is increasing with the scale factor of the Universe as per Equation 1, then the contribution to the BH mass density at redshift z ′ , δρ and the version of Equation 3 with cosmologicallycoupled BHs is thus:

Constraints from the AGN luminosity function at high redshifts
The high redshift AGN luminosity function also provides a constraint on SMBH growth.If the growth of SMBHs by accretion at high redshifts is too slow, there is insufficient mass density in SMBHs to power the observed AGN luminosity density at early times (z > ∼ 2), when a large fraction of the SMBH population was actively accreting.We relate the luminosity density in AGN to the mass density in SMBHs via a mean Eddington ratio, λ = L/L Edd , where L Edd is the Eddington luminosity.(We note that both this analysis and the Soltan analysis above could be made more constraining by considering the observed distribution of SMBH masses and AGN luminosities and using a continuityequation approach (e.g.Raimundo et al. 2012;Tucci & Volonteri 2017).However, given the uncertainties in our understanding of the possible mass/luminosity dependencies of variables such as λ and η, and uncertainties in the exact forms of the mass and luminosity functions and SMBH merger histories, we retain an integral approach for this study.) 3. DATA Two measurements are needed to use the formalism in §2.First is a measurement of the local SMBH mass density.Second is a measure of the accretion history of SMBHs since the first quasars, via a bolometric AGN luminosity function.We describe both these sources of data in this section.

The local mass density of SMBHs
Graham & Driver (2007, hereafter GD07) compared several estimates of the local SMBH mass density available in the literature at the time(see Figure 1, left-hand panel).These were derived from the SMBH mass -bulge velocity dispersion (M • −σ) and the SMBH mass -bulge mass (M • −M bulge ) relations, and adjusted to a common H 0 = 70kms −1 Mpc −1 .GD07 determined a local SMBH mass density in the range 4.4-5.9×10 5 M ⊙ Mpc −3 .Subsequently, Kormendy & Ho (2013) presented a recalibration of the M • − σ and M • − M bulge relations by excluding pseudobulges and galaxy-galaxy mergers from the fitting procedure.This resulted in an increase of the M • /M bulge ratio from ≈ 0.2% to ≈ 0.5%.Thus, the SMBH mass density calculated from a stellar mass function and the Kormendy & Ho (2013) M • − M bulge relation results in a local SMBH mass density substantially higher than the Graham & Driver ( 2007) estimate (we note also that the work of McConnell & Ma (2013) also suggests a relatively high normalization, ≈ 0.3%, however, to be conservative, we adopt the higher Kormendy & Ho (2013) normalization).Furthermore, the difference in the shapes of the stellar mass functions for disks and bulges means that the assumption of a single bulge fraction independent of stellar mass is not strictly correct.
We therefore construct an updated range for the local SMBH mass density as follows.We use the study of Thanjavur et al. (2016), who split the stellar mass function into spheroid and disk components.Using their spheroid mass function and the M • − M bulge relation of Kormendy & Ho (2013) results in an estimate of ρ • (0) = 1.0 × 10 6 M ⊙ Mpc −3 , higher than that of Graham & Driver (2007), but lower than that obtained by assuming a fixed bulge fraction and using an all-galaxy stellar mass function.
The uncertainty in this estimate is substantial.It may be too low -there is evidence that a significant fraction of SMBH mass may be missed by using spheroid mass of velocity dispersion as a proxy.For example, Voggel et al. (2019) argue for the existence of a population of 'stripped' galaxy nuclei in the haloes of giant ellipticals.These objects have modest stellar masses, but can harbor SMBHs of up to ∼ 10 7 M ⊙ .Voggel et al. (2019) fur-Lacy, et al. ther argue that this population may increase the total SMBH mass density by up to ∼ 30%.At high velocity dispersions (or equivalently, high stellar masses) the M • − σ relation may also be biased low, as there are several examples of 'overmassive' SMBHs for their hosts (e.g.Dullo et al. 2021;Nightingale et al. 2023).Another plausible source of error is ejected SMBHs (Postman et al. 2012;Chu et al. 2023).Conversely, it is possible our estimate is too high.Selection biases (Schulze & Wisotzki 2011;Shankar et al. 2016) could lead to an overestimate of the SMBH mass density by a factor ≈ 3, though the existence of this bias has not been confirmed observationally.Also, the spheroid mass function still includes pseudobulges, which, as Kormendy & Ho (2013) emphasize, typically fall below the M • − M bulge relation (Hopkins et al. 2008, though pseudobulges are typically less massive than classical bulges, making up < ∼ 10% of the stellar mass density in spheroids).It is beyond the scope of this paper to synthesize these disparate results into a firm range for ρ • .However, a maximum value of ≈ 1.6 × 10 6 M ⊙ Mpc −1 does not appear to be ruled out (assuming ≈ 30% of the mass density is stripped and ≈ 30% ejected), and a minimum value of 3 × 10 5 M ⊙ Mpc −1 also seems plausible, based on the most extreme estimates of bias.We therefore adopt an estimate and range of log 10 (ρ Finally, we note that the recent detection of a gravitational wave background (GWB) by pulsar timing arrays, if ascribed purely to the mergers of SMBHs, results in a much higher SMBH mass density estimate than those discussed above.Using the results of Agazie et al. (2023) we find that the peaks of their posterior distributions of the most relevant parameters from the fits based on their Phenom+Astro priors (the galaxy stellar mass function characteristic density and cutoff mass, the M • /M bulge ratio, and the logarithmic width of the M • /M bulge distribution) lead to an estimate of 7.4 × 10 6 M ⊙ Mpc −3 .Using the posterior distribution of ρ • constructed from their Monte-Carlo simulations results in a 95% confidence range of 2 − 1400 × 10 6 M ⊙ Mpc −3 and an average of ≈ 1 × 10 8 M ⊙ Mpc −3 (we note also that Casey-Clyde et al. 2022, found a similar result using the earlier 12.5 yr NanoGrav data release).Given the large uncertainties in the SMBH merger estimates it is probably too early to be concerned with the tension between the GWB versus the M • − M bulge relation (see Afzal et al. 2023, for a discussion of the strength of the tension between the SMBH binary model and the measured GWB in the context of new physics), nevertheless we consider the implications of a GWB-based estimate of ρ • (0) = 7.4×10 6 M ⊙ Mpc −3 in our analysis as well.

The bolometric AGN luminosity function
To compute the comoving AGN luminosity density, we need an estimate of the bolometric AGN luminosity function as a function of redshift.Integrating under this luminosity function as a function of redshift then gives the comoving AGN luminosity density.To perform this calculation, we adopt the Shen et al. (2020) (their "global fit B") luminosity function (Equation 10), which agrees well with, but is somewhat more refined than other estimates (e.g.Hopkins et al. 2007;Lacy et al. 2015;Runburg et al. 2022): where the redshift dependencies of ϕ * , L * , γ 1 and γ 2 are detailed in Shen et al. (2020).

RESULTS
Figure 1 (middle and right panels) plots the SMBH mass densities from §3.1 with the result of integrating Equation 9 from z = 6 to z = 0 for fixed η values of 0.1 and 0.3 (Thorne 1974).We show curves for k = 0 (no cosmological coupling), k = 1 (e.g. the exact solutions of Faraoni & Jacques 2007) and k = 3 (Farrah et al. 2023a).This indicates that, if we assume the GD07 value for the local SMBH mass density and ⟨η⟩ ≃ 0.3, then an upper limit can be set of k < ∼ 2. Alternatively, using the GWB value for the local SMBH mass density, then a non-zero value of k is required for any reasonable value of η.
To generalize this argument, we note that, for a given value of k, η can uniquely determined by matching the observed local black hole mass density, so if we assume a value of the Eddington ratio, λ between the minimum value needed to explain the AGN luminosity function at high redshift (assuming all SMBHs are activelyaccreting) and λ = 1 then we can define an allowed region in AGN accretion parameter space as a function of k.In Figure 2, we use the Salpeter time: as a proxy for accretion rate as a function of η, and plot k vs t S to define an allowed region of parameter space in which quasars can (on average) lie (assuming no evolution in the mean η or λ and that the increase in the SMBH mass density is dominated by accretion).Here, a similar result is seen.With the GD07 local density we constrain k to < 2, with our derived density range then 0 < k < 3 can (just) be accommodated (though only with very high values of η ≈ 0.55 and λ ≈ 1), but with the GWB local density range then an approximate lower limit can be set of k > 2, assuming standard accretion disk dominated growth.

DISCUSSION
Conventional assumptions about SMBH accretion, combined with estimates of the local SMBH mass density based on the M • − M bulge relation, imply an allowed range of 0 < k < ∼ 2 (Figure 2).This range is in slight tension (∼ 90% confidence interval) with the results in Farrah et al. (2023a), though it is consistent with some of the recent constraints discussed in §1.2 (e.g.Lei et al. 2023;Rodriguez 2023b;Andrae & El-Badry 2023;Ghodla et al. 2023).
Our constraints on k from the Soltan argument are, however, dependent on four variables: the local SMBH mass density, the integrated luminosity density from the AGN luminosity function and the mean values of λ and η, which regulate the total accretion derived from integrating the AGN luminosity function.We first briefly review the uncertainties on the luminosity density, η and λ (the uncertainties on the local SMBH density are discussed in §3.1), before discussing the effect of these uncertainties on our limit on k.

Uncertainties in the luminosity density
The bolometric luminosity function of AGN is still uncertain, but most of the uncertainty lies at the faint end at high redshifts, whereas the dominant contribution to the luminosity density arises around the break in the luminosity function, where AGN surveys are relatively complete.
As briefly discussed in §3.2, we adopt the luminosity function of Shen et al. (2020), however, if, for example, we use the mid-infrared luminosity function of Glikman et al. (2018) and convert to the bolometric luminosity  1) and the Salpeter time (Equation 11) tS for quasars assuming accretion-dominated growth of the supermassive black hole density.Left: assuming the GD07 estimate of ρ•(0), middle assuming the best estimate from this paper and right assuming the value from the peaks of the GWB posterior fits.The hatched blue regions represent the parameter space where the Salpeter time is too long to grow enough black holes at z = 4 to explain the observed luminosity function and the red regions are those where the Eddington ratio needs to exceed unity to do this.The yellow hatched regions are excluded based on the plausible range of 0.05 < η < 0.3.function using the prescription of Lacy et al. (2015) we recover a black hole mass density 60% higher than that from the Shen et al. ( 2020) luminosity function.(Using this value would be even more constraining on the value of k.) Based on this small difference between two independently-derived luminosity density estimates we believe that the uncertainties from the AGN luminosity function, although not negligible, are smaller than those from the local black hole mass density.

Uncertainties in accretion efficiency
The typical accretion efficiency of an SMBH is poorly constrained (e.g.Raimundo et al. 2012).The canonical value of η ≈ 0.1 (e.g.Soltan 1982) is not measured directly, but obtained by matching early estimates of accretion luminosity to local SMBH mass density.The theoretical limit is thought to be η ≈ 0.3 (Thorne 1974) Observational estimates of η vary from < ∼ 0.1 to potentially > 0.5 (e.g.Trakhtenbrot 2014;Jana et al. 2021;Farrah et al. 2022;Lai et al. 2023).Several factors are thought to affect η; for example, magnetized discs may have η > 0.3 (Avara et al. 2016;Kinch et al. 2021), but there is no consensus on how these factors interrelate.For Figure 2 we have assumed the conventional range of 0.05 < η < 0.3, but a larger range is certainly not excluded.We note though that high values of η reduce the available time for black hole growth by accretion in the early Universe, which is already problematic (e.g.Bañados et al. 2018).

Uncertainties in Eddington ratio
Directly observed Eddington ratios in quasars are consistent with λ ∼ 0.05, with a broad tail towards higher Eddington ratios (e.g.Farrah et al. 2023b).This ratio may however be biased low.Most SMBH mass is likely accreted in obscured phases (Martínez-Sansigre et al. 2005) during which Eddington ratios can be higher (Teng et al. 2014;Farrah et al. 2022).Furthermore, quasar surveys may find AGN whose SMBH masses and Eddington ratios are biased low and high, respectively (e.g.Shen et al. 2008).

CONCLUSIONS
By analysing the growth of the mass density in SMBHs via accretion we show that high values of cosmological coupling (k ≈ 3) needed by Farrah et al. (2023a) to account for dark energy via BHs is disfavored if conventional assumptions regarding the local SMBH mass density and accretion mechanisms are made.However, we cannot rule out scenarios in which k = 3 is possible, and indeed, recent estimates of the GWB from pulsar timing arrays favor a much higher mass density of SMBHs than estimates based on the M • − M bulge relation and would comfortably allow k = 3 (Figure 2, right panel).Better constraints on the local mass function of SMBHs are thus essential both to our understanding of the GWB and for more reliably constraining any cosmological coupling of black holes.

Figure 1 .
Figure 1.Middle and right panels: constraints on the value of the cosmological coupling strength k (Equation 1) based on comparing the growth of SMBH density via accretion to the local density of SMBHs.Middle: assuming η = 0.1 (which agrees well with the GD07 estimate of ρ•(0) if k = 0) and right: η=0.3, the maximum theoretical value for a maximally-spinning BH (Thorne 1974).In both, the red triangle and error bar represents our estimate of the local SMBH mass density from the M• − M bulge relation and the blue square and error bar (which extends well above the plotted range) the preliminary estimate from the GWB, as discussed in §3.2.Also in both, the black lines represent the minimum mass density in SMBHs needed to account for the observed luminosity function for different mean values of λ and the cyan-solid, magenta-dotted and red-dashed lines the growth in BH mass density for different values of k assuming the AGN luminosity function of Shen et al. (2020).The shaded region around the k = 0 (cyan-solid) line corresponds to the uncertainty in Shen et al. (2020) (uncertainty ranges for other values of k are similar).Far left a compilation of estimates of the local mass density in black holes from the literature: S04: Shankar et al. (2004), M04: Marconi et al. (2004), H07: Hopkins et al. (2007), G07: Graham et al. (2007), GD07: Graham & Driver (2007), YL08: Yu & Lu (2008), S09: Shankar et al. (2009), with the cyan shaded area indicating the range of uncertainty at z = 0 for the Shen et al. (2020) luminosity function, assuming η = 0.1.

Figure 2 .
Figure 2. The allowed values of the cosmological coupling strength k (Equation1) and the Salpeter time (Equation11) tS for quasars assuming accretion-dominated growth of the supermassive black hole density.Left: assuming the GD07 estimate of ρ•(0), middle assuming the best estimate from this paper and right assuming the value from the peaks of the GWB posterior fits.The hatched blue regions represent the parameter space where the Salpeter time is too long to grow enough black holes at z = 4 to explain the observed luminosity function and the red regions are those where the Eddington ratio needs to exceed unity to do this.The yellow hatched regions are excluded based on the plausible range of 0.05 < η < 0.3.