Introducing Aurora: A New Black Hole Seeding Model With Explicit Environmental Coupling

abstract · abstract

We present Aurora, a black hole seeding module for the Mimic semi-analytic framework, in which the seed mass of each newly resolved halo is set by its environment. Motivated by evidence that early quiescent galaxies and overmassive black holes assemble in dense environments, Aurora replaces the fixed seed masses of conventional prescriptions with a continuous spectrum, coupling the seed mass to two environmental parameters, the local overdensity and the radial skewness of the surrounding haloes, and to the halo mass at seeding. Requiring only halo masses and positions, it applies to dark-matter-only simulations. Calibrated against Sage on the Millennium merger trees, this prescription leaves a lasting signal. Massive black holes are assembled from the seeds of many progenitors, seeds supplying 29 per cent of all black hole mass formed by $z=0$. A seed-dominated growth pathway, distinct from the accretion-dominated majority, emerges among galaxies assembled rapidly in dense, asymmetric environments. Environment-dependent seeding raises the quiescent fraction of seeded galaxies in massive haloes at $z \simeq 3.1$ by up to 70 per cent, with the excess concentrated in overdense environments, and a deficit in voids. Scrambling the seeding environments while preserving the mass spectrum erases this gradient. The two environmental couplings imprint separable signatures, broadening the scatter of the galaxy-black hole mass relation while preserving the stellar-halo mass relation and other Sage summary statistics. We discuss applying Aurora at gigaparsec scales and higher resolution, to understand the origins of massive black hole systems in contemporary deep surveys.

S6: 6 Conclusions / block 1 · html

A seed-dominated growth pathway emerges, distinct from the accretion-dominated majority ( Section 4.2.1 ). A minority of massive centrals retain more than half of their black hole mass as their original seed; these were seeded in markedly denser and more asymmetric environments than the accretion-dominated majority ( δ = + 1.6 \delta=+1.6 against + 0.06 +0.06 at z ≃ 2.1 z\simeq 2.1 , and + 2.3 +2.3 against + 0.19 +0.19 at z ≃ 3.1 z\simeq 3.1 ), and respond to the couplings by up to + 0.55 +0.55 dex in black hole mass where the majority do not shift noticeably. This channel corresponds closely to the dense environment assembly identified for the earliest massive quiescent galaxies in the Thesan weak feedback experiment ( Chittenden et al., 2026 ) , weak because of an implementation error in the public release ( Garaldi et al., 2026 ) , and to the earlier seeding and faster black hole growth which distinguish them in Colibre ( Chandro-Gómez et al., 2026a ) , recovered here from dark matter information alone. That a gravity-only model reproduces a channel found in two hydrodynamical simulations under contrasting feedback regimes is consistent with the conclusion of Chittenden et al. (2026) that its conditions for black hole synthesis are set by assembly history rather than by the subgrid physics, even where the black hole growth those conditions enable is not.

S2.SS2.SSS1.Px2: Radial Skewness: / block 1 · html

Figure 1 shows the seed mass spectrum which results from the configuration optimised to the stock Sage16 model (see section 3 ), as a joint function of the two environmental metrics: seed mass increases with overdensity across some four decades in 1 + δ 1+\delta , with the radial skewness modulating the relation at fixed density. Without the overdensity coupling, the median seed mass varies by less than 0.1 0.1 dex across the full range of 1 + δ 1+\delta , showing that the trend is not inherited from the host haloes. Figure 1 does not correspond to a snapshot; it pools every seeding event in the volume across all 64 64 snapshots, with each halo entering exactly once, at the moment of seeding. The quantity plotted is the mass assigned at that instant, from which the black hole subsequently grows under Sage ’s accretion physics. The continuous, unimodal spectrum this produces replaces the discrete channel masses of conventional prescriptions. The remainder of this section will illustrate how this relation is established.

S4.SS2.SSS1: 4.2.1 A seed-dominated growth pathway / block 1 · html

Figure 7 shows the resulting M BH − M ∗ M_{\text{BH}}-M_{*} plane at z = 4.179 z=4.179 , an epoch at which black holes overmassive relative to their hosts are reported in abundance ( Harikane et al., 2023 ; Maiolino et al., 2024 ; Übler et al., 2024 ; Kocevski et al., 2025 ) , many of them little red dots. Above 10 10 ​ M ⊙ 10^{10}\,\text{M}_{\odot} the median relation is indistinguishable from that of the stock model, agreeing to within 0.03 0.03 dex in log 10 ⁡ ( M BH / M ∗ ) \log_{10}(M_{\text{BH}}/M_{*}) ; below 10 9.5 ​ M ⊙ 10^{9.5}\,\text{M}_{\odot} it is raised by up to 0.40 0.40 dex. Below 10 10 ​ M ⊙ 10^{10}\,\text{M}_{\odot} , where every observed object above the stock model’s 99 99 th percentile lies, the high ratio tail is extended at every mass. Among the galaxies above 10 9 ​ M ⊙ 10^{9}\,\text{M}_{\odot} , the median ratio moves by only + 0.05 +0.05 dex, while the fraction exceeding 10 − 2 10^{-2} rises from 1.9 1.9 to 3.7 3.7 per cent, the fraction exceeding 3 × 10 − 2 3\times 10^{-2} rises from 0.04 0.04 to 0.22 0.22 per cent, and a population above 10 − 1 10^{-1} appears where the stock model has almost none.

S5.SS1: 5.1 A persistent and verifiable seeding signal / block 1 · html

The N N -body simulation used in this work demonstrates Aurora ’s key capabilities and the sensitivity of its parameters, and the scientific results of the Aurora programme will come from its application to high resolution, large volume simulations. What this paper establishes is the property that makes such applications worthwhile: an environmentally graded seeding prescription leaves a real, measurable signal in the galaxy population long after the seeding events themselves, which may be used in future to constrain black hole seeding mechanisms from observations. Self-regulated accretion drives individual black hole masses towards an attractor set by the host, yet the population level budget of Table 3 shows the seeds persisting in bulk, because what one lineage loses to accretion the population recovers through the coalescence of seeds along its merger tree ( Section 4.2 ). Furthermore, the near-perfect degeneracy of the z = 0 z=0 stellar-halo mass relation and median galaxy-black hole mass relation across the sweep ( Section 3.2 ), which the black hole mass function does not share ( Figure 3 ), does not extend to earlier epochs. The quiescent fraction of seeded galaxies in massive haloes carries an enhancement which grows with the overdensity and skew couplings, reaching Θ = + 0.28 \Theta=+0.28 at α ˇ = 2.6 \check{\alpha}=2.6 and z ≃ 3.1 z\simeq 3.1 , and an environmental gradient across overdensity terciles which grows from Δ ​ Θ = + 0.13 \Delta\Theta=+0.13 in the calibrated configuration to + 0.36 +0.36 at α ˇ = 2.6 \check{\alpha}=2.6 ( Section 4.3 ).

S5.SS1: 5.1 A persistent and verifiable seeding signal / block 1 · html

The fact that the discriminating information should lie in configuration space rather than in the global mass function is a conclusion others are reaching independently. Wang et al. (2026) find the candidate origins of the compact red population to separate in their large scale clustering, with the direct collapse route predicting the strongest bias and a primordial one predicting none above that of the matter field, whereas their mass functions do not separate at all. Shen et al. (2026) reach the same point from the other direction. In Lumina ( Zier and others, 2026 ) , a ( 500 ​ cMpc ) 3 (500\,\text{cMpc})^{3} radiation-hydrodynamical simulation whose seeding depends on halo mass alone, the cross-correlation of little red dots with star forming galaxies is almost identical to the galaxy auto-correlation, and matches the observed amplitude for host haloes slightly above 10 11 ​ M ⊙ 10^{11}\,\text{M}_{\odot} . The authors caution that a formation channel tied to environment or assembly history could change that clustering at fixed halo mass. Environmental seeding is one such channel, and Aurora adds it to the halo mass threshold on which the seeding of Lumina rests.

S6: 6 Conclusions / block 1 · html

In this work we have presented Aurora , a black hole seeding module for the Mimic semi-analytic framework, in which the seed mass assigned to a newly eligible halo is graded by the environment in which that halo forms. The motivation is the persistent failure of conventional models to reproduce the abundance of massive quiescent galaxies and overmassive black holes now observed at high redshift ( Carnall et al., 2023b ; Valentino et al., 2023 ; De Graaff et al., 2024 ; Maiolino et al., 2024 ; Nanayakkara et al., 2024 ; Labbé et al., 2025 ; Weibel et al., 2025 ; Ji et al., 2026b ; Zhang et al., 2026 ) , together with the finding, in Thesan ( Chittenden et al., 2026 ) and Colibre ( Chandro-Gómez et al., 2026a ) alike, that the earliest such systems assemble in dense environments and grow their black holes early and rapidly. This work tests the environmental side of that picture, asking whether environmental seeding leaves a measurable imprint on the galaxy population; the abundance of massive quiescent galaxies is left to an application at higher resolution ( Section 5.3 ). Rather than assigning every halo a fixed seed, Aurora couples the seed mass in log space to the local overdensity δ \delta , the mass-weighted radial skewness μ 3 \mu_{3} of the surrounding halo distribution introduced by Chittenden and Tojeiro (2023) , and the halo mass at seeding, through the three indices of Equation 5 ( Section 2.2.3 ). The environmental metrics are computed once per snapshot and cached ( Sections 2.2.1 and 2.2.2 ), so the prescription requires only the dark matter component of the simulation, and can be applied to volumes far beyond the reach of hydrodynamics.

S6: 6 Conclusions / block 1 · html

We calibrated the prescription against the Sage16 semi-analytic model ( Croton et al., 2006 ; Croton et al., 2016 ) run on the Millennium merger trees, matching the z = 0 z=0 black hole mass function and the count of black hole hosts ( Section 3.1 ), and explored the behaviour of each coupling through one-dimensional sweeps about the calibrated configuration ( Section 3.2 ). We then followed the consequences forward in time, from the seed mass spectrum to the assembly of massive black holes and the quiescent galaxy population ( Section 4 ), using the metrics defined in Table 3 .

S6: 6 Conclusions / block 1 · html

Our findings can be summarised as follows:

S6: 6 Conclusions / block 1 · html

The three couplings act in structurally different ways, and the difference follows from the distributions of the quantities which they multiply, rather than from the couplings themselves ( Section 4.1 ). Because δ \delta is normalised to the geometric mean density, ln ⁡ ( 1 + δ ) \ln(1+\delta) has a median of + 0.002 +0.002 across seeding events and the overdensity coupling pivots the seed spectrum, boosting seeds in nodes and suppressing them in voids. The skewness μ 3 \mu_{3} is positive for 72 per cent of events, so the skew coupling instead raises seed masses almost everywhere. The two are nearly orthogonal, with a Spearman rank correlation of − 0.038 -0.038 .

S3.SS1: 3.1 Redshift evolution of the black hole mass function / block 1 · html

By z ∼ 7 z\sim 7 , growing black holes in the heavy-seed astrophysical models ( hs-edd1 and hs-edd5 ) generally accrete a fraction of the available gas mass (since this is lower than the fractional Eddington mass). This results in a convergence of the BHMFs - for example, for a black hole mass of 10 7 ​ M ⊙ 10^{7}\mathrm{M}_{\odot} , the number densities predicted by these two models are only different by a factor 2. The black holes in these models grow to masses of about 10 10 ​ M ⊙ 10^{10}\mathrm{M}_{\odot} for number densities of about 10 − 8 ​ c ​ M ​ p ​ c − 3 10^{-8}{\mathrm{c}Mpc^{-3}} . While black holes with M B ​ H ∼ < 10 7 ​ M ⊙ \mathrm{M}_{BH}\lower 2.15277pt\hbox{$\;\mathrel{\mathop{\kern 0.0pt\sim}\limits^{<}}\;$}10^{7}\mathrm{M}_{\odot} still grow slower in the ls-edd5 model, for higher masses, that BHMF lies between the heavy Eddington-limited and super-Eddington models and reaches out to 10 10 ​ M ⊙ 10^{10}\mathrm{M}_{\odot} . Again, the ls-edd1 model shows the lower limit to the growth of light black holes, only reaching maximal masses of about 10 7.5 ​ M ⊙ 10^{7.5}\mathrm{M}_{\odot} for a number density of 10 − 8 ​ c ​ M ​ p ​ c − 3 10^{-8}{\mathrm{c}Mpc^{-3}} . Interestingly, the PBH BHMF shows a clear excess - by about two orders of magnitude - of intermediate-mass black holes at M B ​ H ∼ 10 3 − 4 ​ M ⊙ \mathrm{M}_{BH}\sim 10^{3-4}\mathrm{M}_{\odot} compared to any astrophysical seeding model; at M B ​ H ∼ 10 5 − 8 ​ M ⊙ \mathrm{M}_{BH}\sim 10^{5-8}\mathrm{M}_{\odot} , it predicts results that are very similar to the heavy seeding models. We also find both the heavy seeding ( hs-edd1 and hs-edd5 ) and PBH models to be in excellent agreement with the observationally-inferred BHMF at z ∼ 7 z\sim 7 ( Fei et al., 2026 ) . Light seeding models predict a much flatter slope than the observations although the ls-edd5 model is in accord with data for black holes with M B ​ H ∼ > 10 6.5 ​ M ⊙ \mathrm{M}_{BH}\lower 2.15277pt\hbox{$\;\mathrel{\mathop{\kern 0.0pt\sim}\limits^{>}}\;$}10^{6.5}\mathrm{M}_{\odot} ; given the slower growth of light-seeded Eddington limited black holes, these observations can already be used to rule out the ls-edd1 model. We

S3.SS1: 3.1 Redshift evolution of the black hole mass function / block 2 · html

briefly note that given our “heavy seeds” have a range of masses between 10 3 − 5 ​ M ⊙ 10^{3-5}\mathrm{M}_{\odot} , we do not find any specific peak in the BHMF; this is in contrast to the results of Jeon et al. (2025) who find a well-defined peak at 10 5 ​ M ⊙ 10^{5}\mathrm{M}_{\odot} , possibly since all of their heavy seeds are initiated with this fixed mass.

S3: III Results and Discussions / block 1 · html

Based on the initial conditions and numerical setup described above, we now present the simulation results and discuss their physical implications. To highlight the environmental dependence of PBH-seeded evolution, we organize this section around four connected aspects. We first examine how the large-scale environment drives distinct BH–host evolutionary pathways in Section III.1 . Building on these trends, we identify the main evolutionary phases and possible transitions toward normal galaxy–AGN configurations in Section III.2 . We then connect these regimes to stellar growth and metal enrichment in Section III.3 before translating the resulting evolutionary stages into observable signatures for current and future surveys in Section III.4 .

S3.SS2: III.2 Road to Normalcy / block 1 · html

The results of Section III.1 show that the cosmological environment regulates the gas supply and co-evolution of PBH-hosting systems. We now examine how the initial PBH mass and its placement within the assembling overdense region modify this evolution. Here, “normalcy” is defined in an evolutionary, trend-driven sense. It does not require exact agreement with the local BH–galaxy scaling relations. Instead, it refers to a progressive reduction in the relative dominance of the initial seed, such that stellar and halo growth outpace BH growth and drive both M ∙ / M ⋆ M_{\bullet}/M_{\star} and M ∙ / M h , PBH M_{\bullet}/M_{\rm h,PBH} toward the ranges occupied by more conventional high-redshift galaxy–AGN systems.

S3.SS2: III.2 The variation due to feedback models / block 1 · html

To capture this mass dependence, we introduce a critical stellar mass at which feedback strength varies due to the deepening gravitational potential of galaxies. Motivated by cosmological galaxy-formation simulation results, we parameterize the stellar-mass-dependent feedback model using a hyperbolic tangent function,

S4: IV Summary & Discussions / block 1 · html

In this Letter, we present a BH growth model coupled with DM halo assembly to explore the emergence of overmassive BH populations with M BH / M ⋆ > 0.01 M_{\rm BH}/M_{\star}>0.01 , exceeding the values observed in the local universe. We parameterize AGN feedback-driven mass loss as M ˙ ​ ( r ) ∝ r p \dot{M}(r)\propto r^{p} ( 0 < p < 1 0<p<1 ), providing a semi-analytical framework to map evolutionary trajectories in the M BH − M ⋆ M_{\rm BH}-M_{\star} diagram. The growth rates of BHs and galaxies form a converging flow pattern toward a linear attractor set by feedback strength p p and gas supply efficiency to the nuclei. Overmassive seed BHs grow more slowly than their hosts, while undermassive seeds experience rapid growth before aligning with the attractor trajectory.

S5: V Conclusions / block 1 · html

We have proposed a gas-triggered SIDM-accretion model for Little Red Dots. The central idea is that a compact nuclear thick disk can do two things at once: compress the central SIDM enough to ignite rapid dark accretion onto a pre-existing black hole, and provide the high-column gas that later reprocesses AGN radiation into an LRD spectrum. For a fiducial RUBIES-EGS-55604-like halo, the model gives M BH ≃ 9 × 10 7 ​ M ⊙ M_{\rm BH}\simeq 9\times 10^{7}M_{\odot} in 0.32 0.32 Myr for σ / m = 1 ​ cm 2 ​ g − 1 \sigma/m=1\,{\rm cm^{2}\,g^{-1}} , while the subsequent gas-powered luminosity remains sub-Eddington.

S5: V Conclusions / block 1 · html

The model naturally links several observed LRD properties. The black hole can become overmassive without requiring hyper-Eddington radiative growth. The compact disk provides Compton-thick columns, explaining weak escaping X-rays and strong reprocessing. A smooth Balmer-active disk atmosphere produces the Balmer-limit curvature. Orientation determines whether the same event appears as a red LRD or a bluer compact AGN. The LRD phase ends when the high-covering column is consumed, expelled, diluted, or no longer replenished, giving a natural lifetime of tens to hundreds of Myr.

abstract / block 1 · abstract

We present a theoretical framework for seeding and growing supermassive black holes (SMBHs) in dark matter halos along their assembly histories. Seeds are bred out of Pop-III stars formed during the first collapse of pristine gas in mini-halos that have reached the $\rm H_2$-cooling limit, modulated by UV radiation from star formation and dynamical heating from fast halo assembly. Such breeding persists until the enrichment of the intergalactic medium (IGM) enables Pop-II stars to form. Post-seeding growth of black holes (BHs) is driven by distinct channels, starting with episodic super-Eddington accretion associated with nuclear bursts induced by global disturbances of galaxies, followed by sustained sub-Eddington accretion via capturing sub-clouds formed in self-gravitating gas clouds (SGCs) in halos of fast assembly, and ending with merger-dominated, quiescent growth. We implement the model in subhalo merger trees to build a coherent framework to follow SMBH-galaxy-halo co-evolution across the whole history and structural hierarchy. BH seeds are bred with a broad mass spectrum of $M_{\rm BH} = 10 - 10^5\,{\rm M}_\odot$ at $z \approx 20 - 30$ in mini-halos with masses of $10^5 - 10^8\,{\rm M}_\odot$. Nuclear bursts provide the key condition for seeds to grow into SMBHs. The $M_{\rm BH}$-$M_*$ relation is a multi-piece, redshift-dependent function shaped by the interplay among different growth channels. Our model predictions are broadly consistent with existing observations; especially, a population of BHs reminiscent of 'little red dots' (LRDs) discovered by JWST naturally results from the seeding and growing processes. Potential future tests of the model are discussed.

S2.SS1: II.1 The model of galaxy formation and black hole growth / block 1 · html

A feedback-regulated initial mass function (IMF). Here, the feedback from the first Pop-III star cluster formed in the collapsed gas of a pristine halo is evaluated to determine the resulting mass of individual stars in the cluster.

S3.SS1: III.1 The evolution of number density / block 1 · html

Fig. 2 shows the predicted number density ( n n ) of LRDs as a function of z z based on our default selection (solid red curve). For comparison, we show the observational results from Y. Ma et al. [45] and V. Kokorev et al. [37] , and the empirical prediction by K. Inayoshi [29] . The number density of LRDs predicted by our model recovers the piece-wise evolution found in observations: a plateau at z ≳ 4 z\gtrsim 4 where the number density reaches a few times 10 − 5 ​ cMpc − 3 10^{-5}\,{\rm cMpc}^{-3} and a sharp drop at z ≲ 4 z\lesssim 4 .

S7: VII Conclusion / block 1 · html

We investigate the origin of the high-redshift SMBHs in LRDs observed by JWST from small-scale clustering of PBHs followed by runaway mergers. The derived mass function of LRDs under the small-scale PBH clustering matches the JWST observations for 1 ​ – ​ 10 % 1\textup{--}10\% accretion duty cycle. Other properties of LRDs—including their overmassive nature and compactness—are also explained within our framework when the seed effect and the halo spin distribution are taken into account.

article / block 1 · html

Supermassive black holes (SMBHs) in “little red dots” (LRDs) discovered the James Webb Space Telescope (JWST) may result from runaway mergers of primordial black holes (PBHs) in clusters—through long–short mode coupling on small scales in the early Universe. In this framework, we derive the SMBH mass function, together with the compactness and overmassive features of LRDs. We also estimate that the dense gas residing in PBH clusters is consistent with LRD observations. In addition, SMBHs formed from PBH clusters can help accelerate galaxy formation at high redshifts, thus alleviating tension with Λ \Lambda CDM cosmology.

S2.SS2: II.2 Star Formation and Stellar Feedback / block 1 · html

Table 2: Summary of supernova (SN) feedback parameters from stellar particles. We show values for a single stellar particle of m ⋆ ≃ 600 ​ M ⊙ m_{\star}\simeq 600~\mathrm{M}_{\odot} , based on the average over the stellar IMF. More specifically, we assume that all Pop II stellar particles are identical, whereas Pop III particles contain individual stars sampled from the input IMF [ 47 , 46 ] . Here, E SN E_{\rm SN} is the total SN energy, M Z M_{Z} the mass of the metals released, t ⋆ t_{\star} the typical life time of massive stars, and r final r_{\rm final} the final radius of the SN bubble. One can refer to table 1 in Liu and Bromm [62] for a more detailed description.

S3: III Results / block 1 · html

With the stellar feedback recipes described in Section II , we now analyze the resulting PBH simulation runs, comparing to previous work 4 4 4 Also see the appendix in [ 73 ] . that did not include the full feedback physics. To study the synergy between stellar and BH feedback, metal enrichment and gas outflow, we first track the star formation history (SFH) and the assembly of galaxies in Section III.1 , followed by a discussion of the metal enrichment history in Section III.2 .