Evidence for cosmic evolution in the spin of the most massive black holes

We use results from simulations of the production of magnetohydrodynamic jets around black holes to derive the cosmic spin history of the most massive black holes. We assume that the efficiency of jet production is a monotonic function of spin a, as given by the simulations, and that the accretion flow geometry is similarly thick for quasars accreting close to the Eddington ratio and for low-excitation radio galaxies accreting at very small Eddington rates. We use the ratio of the comoving densities of the jet power and the radiated accretion power associated with supermassive black holes with Mbh>~10^8 Msol to estimate the cosmic history of the characteristic spin a. The evolution of this ratio, which increases with decreasing z, is consistent with a picture where the z~0 active galactic nuclei have typically higher spins than those at z~2 (with typical values a~0.35-0.95 and a~0.0-0.25 respectively). We discuss the implications in terms of the relative importance of accretion and mergers in the growth of supermassive black holes with Mbh>~10^8 Msol.


INTRODUCTION
Astrophysical black holes are described by two parameters, mass m• and spinâ, defining the structure of space-time within regions close to the event horizon. Every massive galaxy has a supermassive black hole (SMBH) at its centre, and the combination of accretion and coalescence of SMBHs of similar mass, major merging, can lead to a wide spread in spins (e.g. Berti & Volonteri 2008). A possible constraint on the cosmic spin history is the ratio of kinetic to radiated outputs in active galactic nuclei (AGNs), since these are powered by SMBHs.
During accretion onto black holes, the amount of energy available for radiation and for the production of jets is determined by the mass and spin of the black hole, as well as the rate of accretion of matter and the geometry of the accretion flow. More massive black holes can accrete more matter, but those that are spinning rapidly can extract energy more efficiently from the infalling material (e.g. Novikov & Thorne 1973;Blandford & Znajek 1977). In addition, geometricallythick accretion flows can power jets more effectively than thin ones because they are able to sustain more powerful poloidal magnetic fields (e.g. Meier 2001).
In galactic black holes (GBHs), the presence of steady jets occurs typically during low accretion rates (where by low we mean an Eddington ratio λ≡ṁ/ṁ Edd∼ < 10 −2 ). During transitions to moderate accretion rates, λ ∼ > 10 −2 the jet ⋆ alejo.martinez-sansigre@port.ac.uk is temporarily enhanced in power, and can then become steady, transient or totally absent (Fender et al. 2004). This is interpreted as a transition from advection-dominated accretion flows (ADAFs Narayan & Yi 1995) that are geometrically thick and can sustain significant jets, to gas-pressuredominated accretion flows which are thin, with aspect ratios (the ratio of height to radius) ∼0.01 (Novikov & Thorne 1973).
Quasars and high-excitation radio galaxies are SMBHs accreting at a high fraction of their Eddington limiting luminosity, so that the accretion flow is geometrically thick (Novikov & Thorne 1973, note GBHs with λ∼1 also produce continous jets, e.g. GRS 1915+105, Fender et al. 2004). Low-excitation radio galaxies are SMBHs with very low accretion rates, ADAFs, that have similarly thick accretion flows. The geometry is thus similar at both high and low accretion rates.
SMBHs are known to display a wide range of jet powers with no observable differences in the accretion flow, suggesting that a hidden variable must be controlling the jet power. The 'spin paradigm' assumes the black hole spin to be the physical parameter controlling the kinetic output in the form of jets (e.g. Blandford & Znajek 1977;Sikora et al. 2007).
GBHs show a similar wide spread in jet powers, but these are not found to correlate with the published estimates of the spins (Fender et al. 2010). If the spin estimates are correct, this provides very stringent evidence against the spin paradigm for GBHs. However, due to the uncertainties in the measurements of both the black c RAS hole spin and the total jet power, the results found by Fender et al. (2010) do not provide robust evidence against the spin paradigm. The reader is referred to Section 7.2 of Martínez-Sansigre & Rawlings (2011, hereafter MSR11) for a detailed discussion.
In this letter we make the assumption that spin is indeed an important factor in powering jets, following the results of general-relativistic magnetohydrodynamic simulations (e.g. Hawley & Krolik 2006). Hence, assuming the spin paradigm to be correct, we present an observational constraint on the evolution of the mean spin of the most massive SMBHs as inferred from a measured cosmic change in the ratio of the power output in jets to the power radiated by active galactic nuclei (AGNs). We adopt a ΛCDM cosmology with the following parameters: h = H0/(100 km s −1 Mpc −1 ) = 0.7; Ωm = 0.3; ΩΛ = 0.7.

ASSUMPTIONS
High accretion rates will lead to efficient radiation of the thermal energy originating from viscous forces caused by magnetohydrodynamic turbulence. A fraction ǫ of the accreted energy will be radiated away, and the remaining energy will be advected into the black hole. The bolometric power available for radiation can therefore be described by: where L bol is the bolometric luminosity due to radiation, m• is the rate of accretion of mass onto the SMBH, and the term ǫ is the radiative efficiency. The shearing of magnetic fields frozen into the accreting plasma will allow the extraction of electromagnetic energy to power jets. Close to a rotating black hole, the dragging of inertial frames will lead to a higher angular velocity (and hence larger shear) and additional amplification of the magnetic fields occurs by extraction of energy from the rotating black hole. The power available for the production of jets can be described as: where Qjet is the jet power and η is the jet efficiency (e.g. Hawley & Krolik 2006, this is also the same as ǫ kin in Merloni & Heinz 2008). Here it is assumed that η is a function of spin so that rapidly spinning black holes will have a higher jet efficiency. The AGNs concerned in this letter are powered by SMBHs accreting either at very low or at very high fractions of their Eddington limiting luminosity (so either λ ∼ < 10 −2 or λ∼1 Willott et al. 2001;Smolčić et al. 2009). In these two limits of low-and high-accretion rates the geometry of the accretion flows will be similar, with aspect ratios ∼ > 0.1. Hence, we assume that η is a monotonic function ofâ and that variations due to the accretion flow thickness are a secondary effect (see also Section 2.2 of MSR11 for more details).

THE COMOVING JET POWER AND RADIATED ENERGY
To infer the cosmic spin history of SMBHs, we will use the mean comoving kinetic jet power density from AGNs pro-ducing jets, ρQ jet , and the mean comoving jet power from radiating AGNs, ρL bol . These quantities can both be measured using: where φQ jet (Qjet, z) = d 2 N/dQjetdz and φL bol (L bol , z) = d 2 N/dL bol dz are the luminosity functions for Qjet and L bol respectively, representing the comoving number density of AGNs with a given kinetic jet power or bolometric luminosity, and their evolution with cosmological redshift z The jet power can be estimated from the low-frequency radio luminosity density due to synchrotron radiation (e.g. Willott et al. 1999), assuming that the jet output results in energy stored in radio source lobes together with associated and ineluctable work done on the source environment. Adopting standard assumptions for the way this energy is shared between magnetic fields and particles, the work done inflating these lobes as well as typical jet advance speeds and a typical density profile for the intergalactic medium, ρQ jet can be estimated from the radio luminosity function φL ν (Lν , z): where Lν is the luminosity density at 151 MHz in W Hz −1 sr −1 . The term f represents the combination of several uncertainty terms when estimating Qjet from Lν (for example plasma volume filling factor or the fraction of energy in non-radiating particles). An alternative method to estimating jet powers is to estimate the energy stored in cavities in the X-ray emitting intra-cluster gas, and both methods agree best for f ∼20 (Cavagnolo et al. 2010). We therefore initially take f =20 but in Section 4 consider deviations from this assumption.
The radio luminosity function includes two components: the luminosity function for the high radio luminosity classical double sources (Willott et al. 2001) and the luminosity function for the lower radio luminosity AGNs (Smolčić et al. 2009) 1 . For the lower radio luminosity AGNs, the luminosity function extends to z=1.4 only, and we have extrapolated it to z=2.0 using the same functional form for the evolution. The luminosity functions of other AGNs suggest a turnover around z=2, so no further extrapolation is attempted. We neglect any errors in the quantity inside the integral in Equation 5 which is reasonable because it is essentially the luminosity density, a quantity dominated by objects near the break in the radio luminosity function: the normalization at this break is directly measured to z ∼ 2 (Willott et al. 2001;Vardoulaki et al. 2010).
In AGNs the hard (high-energy) X-ray luminosity originates from inverse-Compton scattering of photons from the accretion flow by a corona of hot electrons and it is a good tracer of the bolometric luminosity. The X-ray luminosity function φL X (LX, z ) (Silverman et al. 2008) can be used to estimate ρL bol : The factor FCT accounts for the fraction of luminous Compton-thick AGNs missed by the hard X-ray surveys: (FCT≈0.5, e.g. Martínez-Sansigre et al. 2007;Gilli et al. 2007), and we adopt an uncertainty (from Table 8 of Fiore et al. 2009). so that FCT = 0.5 ± 0.1. The bolometric correction CX converts from a monochromatic luminosity to a bolometric luminosity, L bol = CXLX: we use the values given by (Hopkins et al. 2007) with, following Martínez-Sansigre & Taylor (2009), an assumed 8 per cent uncertainty.
The value of ρQ jet is always dominated by the output from the most massive SMBHs, with m• ∼ > 10 8 M⊙ at all redshifts Smolčić et al. 2009). To ensure that ρL bol is also dominated by the most massive black holes, we only integrate above X-ray luminosities that can be achieved by SMBHs with m• 10 8 M⊙ accreting at 25% of the Eddington accretion rate. This is a compromise between completeness to SMBHs with m• 10 8 M⊙ and contamination from lower-mass SMBHs. This value securely avoids serious contamination from SMBHs with masses <10 8 M⊙. From the distribution in λ found by Heckman et al. (2004, their Figure 3) we estimate that it leads to ρL bol at low-z being underestimated by a factor ∼3. At high redshift we expect it to underestimate ρL bol by a negligible amount.
This evolution of the ratio ρQ jet /ρL bol with redshift z for the most massive SMBHs (with m• ∼ > 10 8 M⊙) is shown in Figure 1. There is likely some exaggeration of the rise in ρQ jet /ρL bol towards low z, because of the cut-off assumed in X-ray luminosity, but only by a factor ∼0.5 (base-10 logarithm), smaller than the inferred cosmic rise towards low z. Note also that the absolute value of the height of the curve will scale with f as f 3/2 .

THE COSMIC SPIN HISTORY
The ratio of ρQ jet and ρL bol , for the most massive black holes ( ∼ > 10 8 M⊙) is then given by that includes the comoving density of accretion of the AGNs contributing to ρQ jet , namely ρṁ •Q jet and of AGNs contributing to ρL bol , namely ρṁ •L bol . Separating the contribution to ρQ jet from ADAFs and quasars, ρQ jet =ηADAFρṁ ADAF + ηQSOρṁ QSO . In the case of ρL bol , although ρṁ ADAF ρṁ QSO , ǫADAF≪ǫQSO, so that ρL bol is totally dominated by quasars and ρL bol =ǫQSOρṁ QSO . Hence, We consider two approaches to using ρQ jet /ρL bol to constrain the cosmic spin history. The first one is to attempt to use the ratio to find a fiducial spin, under some simplifying assumptions. The second is to compare a parametric model to the observations, to see how well it reproduces them. In both approaches, the jet efficiency η is described using the results from a set of three-dimensional generalrelativistic magnetohydrodynamic simulations, which can be approximated as (Hawley & Krolik 2006):

The fiducial cosmic spin
At z ∼ > 1, the density of accretion onto ADAFs is insignificant compared to the density of accretion onto quasars (e.g. see Figure 9 of MSR11 and Figure 4 of Merloni & Heinz 2008), ρṁ ADAF ≪ρṁ QSO so that ρQ jet /ρL bol =ηQSO/ǫQSO. Therefore, at high redshift ρQ jet /ρL bol gives a good estimate of the 2 In Section 5, we will discuss the relevance of chaotic accretion. In this paradigm, an average efficiency of co-and counter-rotating accretion should be used (King et al. 2008). The results of MSR11 used jet efficiencies from the literature for co-rotating accretion, so for consistency we limit ourselves to these. The radiative efficiency is determined for quasars and their spins are close toâ=0.0 (MSR11) so the co-and counter-rotating radiative efficiencies are virtually identical. Therefore, although there is a slight conceptual inconsistency between our model radiative efficiency and our interpretation, in practice it makes no difference to the results. ratio of efficiencies for quasars and, under the spin paradigm, of their typical spin.
At very low redshifts the comoving density of accretion onto ADAFs is comparable to that of quasars, so ρṁ ADAF ≈ ρṁ QSO (Merloni & Heinz 2008, MSR11), in which case ρQ jet /ρL bol ≈ (ηADAF/ηQSO + 1) ηQSO/ǫQSO. It is therefore not straightforward to interpret ρQ jet /ρL bol in terms of typical efficiencies, but some progress can be made under the assumption that ηADAF ∼ > ηQSO, so that ρQ jet /ρL bol ∼ > 2ηQSO/ǫQSO. This shows that at z∼0 the ratio ρQ jet /ρL bol will overestimate the typical efficiencies of quasars and ADAFs. Figure 2 represents a fiducial cosmic spin history of the most massive black holes inferred from ρQ jet /ρL bol . It represents the characteristic spin if at a given redshift all SMBHs have the same spin. It is derived by assuming the first term in Equation 8 is negligible, so we expect this figure to provide an accurate estimate at moderate and high redshifts, but to overpredictâ fid (z) at the low redshifts, marked in gray. Overestimating η/ǫ by a factor of >2 at z∼0 suggests that the typical spin at z∼0 is <0.9 rather than ∼0.95. Consideration of f suggests that lower values of f also lead to lower estimates of the typical low-z spin: re-running our analysis with f = 10 (rather than f = 20) leads to no changes at high z but a typical low-z spin of 0.8. We note that high values of f are supported by X-ray cavity observations undertaken at low radio powers and low redshift so lower values of f for the low-z population are hard to envisage. If we assume that different effects (e.g. f and underestimating ρL bol at low-z) add to an overall overestimate of ρQ jet /ρL bol by a factor of ∼2 (4), this will result in a fiducial spin at z∼0 of a fid ∼0.9 (0.7). In these cases the inferred evolution would be less extreme, but still present.
To test wether this fiducial spin history is consistent with constraints from the radiative efficiency of quasars (the So ltan 1982, argument), we compute the luminosity-weighted spin (see also MSR11): â L ≡ â fid (z)ρL bol dz/ ρL bol dz. The spin history of Fig. 2 yields a luminosity-weighted spin is â L=0.055, which corresponds to a radiative efficiency of ǫ L=0.059. This is in good agreement with observational constraints, which typically yield values ǫ L∼0.05-0.10 (see e.g. Martínez-Sansigre & Taylor 2009).

Can we explain ρQ jet /ρL bol ?
In MSR11, we inferred a moderate evolution amongst the SMBHs with m• 10 8 M⊙. This evolution was really driven by the change in relative space densities of two populations: the high Eddington-rate SMBHs, or quasars, which dominate at high redshifts, were found to have spin distributions centred around low spin,â∼0. The low Eddington-rate SMBHs, or ADAFs, become significant at low redshifts, and these were found to have a bimodal spin distribution, with a peak centred atâ∼0 and another atâ∼1. These spin distributions allowed us to explain the local radio luminosity function, as well as explaining its evolution up to z∼1.
The results of MSR11 suggest that the evolution of the mean spin of SMBHs is driven by a transition from a high-accretion rate, low spin population, to a low-accretion rate population showing a bimodality in spins (and hence a slightly higher mean spin). This evolution was qualita- Figure 2. Fiducial spin of the most massive black holes as a function of redshift. It shows the posterior probability distribution forâ fid (z) given the ratio ρ Q jet /ρ L bol plotted in Figure 1, given η(|â|) and ǫ(|â|). The posterior probability has been normalised at each redshift. The underlying assumption is that the aspect ratio of the accretion flow is similar for all the AGNs contributing to ρ Q jet /ρ L bol , so that geometry has a negligible effect in the evolution of ρ Q jet /ρ L bol , which is determined by the evolution in cosmic spin. The colours and contours are the same as those in Figure 1. At z ∼ > 0.9, the maximum posterior probability line as well as the -1 and -2σ contours go along theâ=0 axis. The grey shades mark the regions where the assumption that ρṁ ADAF ≪ ρṁ QSO breaks down. In the light gray region, we estimate that ρṁ ADAF ≈ 0.1ρṁ QSO while in the dark gray region, ρṁ ADAF ≈ 0.25ρṁ QSO (MSR11). Hence in these regions we overestimateâ fid (z).
tively similar to that seen in Fig. 2, but it was significantly weaker, changing from â ∼0.25 at high redshift to ∼0.35 at low redshift.
We can test whether the spin distributions we inferred in MSR11 can reproduce the observed evolution in ρQ jet /ρL bol . We use these spin distributions to predict how the ratio ρQ jet /ρL bol evolves with redshift using ρQ jet ρL bol pred = ηADAF ǫQSO Note that the angled brackets indicate expectation values: ηADAF = η(â)PADAF(â)dâ, ηQSO = η(â)PQSO(â)dâ and ǫQSO = ǫ(â)PQSO(â)dâ. The variance of the ADAF term is large, due to the bimodal nature, while the variance of the QSO term is small. As mentioned earlier, at z ∼ > 1 ρṁ QSO ≫ρṁ ADAF so that the ηQSO term dominates, but at low redshift ρṁ QSO ∼ρṁ ADAF . Fig. 3 shows the resulting evolution of ρQ jet /ρL bol , from Equation 10 and using the spin distributions from MSR11. At z∼0 the most likely predicted value of ρQ jet /ρL bol is about 0.1 dex higher than the observed value, and at z ∼ > 1 it is about 0.3 dex higher. This can be understood since in MSR11 (Figs. 4 and 7 of that work) we overpredicted the radio luminosity function of high-Eddington rate SMBHs, which dominate ρQ jet at high redshifts.
The expectation value is systematically slightly higher than the most likely value of ρQ jet /ρL bol shown in Fig. 1, but the curvature is similar and the two figures agree very well within their uncertainties. Similar results are found using Figure 3. The evolution of ρ Q jet /ρ L bol as predicted by the spin distributions from MSR11. In MSR11 we predicted spin distributions for high-Eddington rate objects (QSOs) and low-Eddington rate ones (ADAFs), as well as the accretion density onto these two sub-populations (ρṁ QSO and ρṁ ADAF respectively). Following Equations 8 and 10, and the spin distributions from MSR11, we predict an expectation value that behaves very similarly to that of Fig. 1. The contours are wide due to the variance of the spin probability distributions for high-and low-Eddington rate objects. At low redshifts the variance is larger due to the bimodality in the spin distribution of the low-Eddington rate objects.
other jet efficiencies from the literature (see supplementary figures).
Although the amount of evolution in the mean spin is different, Fig. 2 and the mean spin from MSR11, which can reproduce Fig. 1, agree qualitatively: the typical spin was very low at high redshifts, and has increased at low redshifts.

DISCUSSION
Using the observed ratio between the comoving densities of jet power and radiation from SMBHs, ρQ jet /ρL bol , we have inferred that the spin of the SMBHs must be evolving. We have done this in two ways. The first method included some simplifying assumptions to allow us to "invert" the problem. This allowed us to infer a fiducial spin for all SMBHs, and we find this fiducial spin evolves strongly. The second method consists of comparing some highly parametrised description of the spin distributions, inferred from fitting the local radio luminosity function in our previous work (MSR11). There are no simplifying assumptions, but rather than infer the evolution of a fiducial spin value, we simply test whether spin distributions previously inferred in MSR11 can explain the observed evolution of ρQ jet /ρL bol : we find that they can explain this evolution very well.
The results from both methods are quantitatively different, in that they differ in the spin values inferred at low redshift. However the two methods agree qualitatively: a lowspin epoch occurs at z ∼ > 1 and coincides with the period during which the black holes with masses m• ∼ > 10 8 M⊙ accreted most of their mass (Yu & Tremaine 2002). This epoch is one characterised by high accretion rates onto low-spin objects, while the present-day epoch (z∼0) is characterised by low accretion rates onto a population of SMBHs where a significant fraction have high spins.
Somewhat counter-intuitively, the epoch of lowest-spin corresponds to the epoch of highest accretion. This, however, is consistent instead with the paradigm of "chaotic accretion" (e.g. King et al. 2008): approximately half of the accretion occurs with the disc counter-rotating with respect to the SMBH, hence decreasing the spin of the hole as opposed to the increase brought by a co-rotating disc.
Chaotic accretion is expected to lead to SMBHs with low spins,â∼0.1 (King et al. 2008), and our result is thus in agreement with black holes being "spun down" at high redshift, when the typical accretion rate was higher, due to the greater available supply of concentrated cold gas in the central regions of the host galaxy (e.g. Obreschkow & Rawlings 2009).
When two black holes of similar mass coalesce, the orbital angular momentum contributes significantly to the final angular momentum of the coalesced black hole, leading to high values of the final spin, typicallyâ∼0.7 (Rezzolla et al. 2008). Hence, mergers of SMBHs of similar mass (major mergers, with mass ratios ∼ < 4:1) are expected to provide the "spin-up" mechanism.
SMBHs of mass ∼ > 10 8 M⊙ are hosted by galaxies with bulge masses ∼ > 10 11 M⊙, typically elliptical galaxies. A significant fraction of these galaxies have undergone major mergers since z ∼ 1 (Robaina et al. 2010), meaning their central SMBHs will have most likely undergone a major merger too.
Due to the steepness of the mass function above the break in the mass function, our result is most likely dominated by the behaviour of SMBHs hosted by galaxies with masses ∼1 × 10 11 M⊙ and reflects their merger history, rather than reflecting that of the rarer, more massive galaxies with >3 × 10 11 M⊙, which seem to have undergone an earlier evolution (Banerji et al. 2010).
The mechanism for spinning-up SMBHs, major merging, is present during the entire epoch of 0 z 2 since the merging of galaxies with similar masses is also a common occurrence at redshifts z∼2. However, the centrally concentrated, typically molecular, cold gas content of the merging galaxies, the fuel for accretion, decreases significantly with cosmic time.
It is therefore the gradual disappearance of the "braking" mechanism, chaotic accretion, which determines the increase in the typical spin of the most massive black holes.