Stirring the cosmic pot: how black hole feedback shapes the matter power spectrum in the fable simulations

ABSTRACT

Understanding the impact of baryonic physics on cosmic structure formation is crucial for accurate cosmological predictions, especially as we usher in the era of large galaxy surveys with the Rubin Observatory as well as the Euclid and Roman Space Telescopes. A key process that can redistribute matter across a large range of scales is feedback from accreting supermassive black holes. How exactly these active galactic nuclei (AGNs) operate from sub-parsec to Mega-parsec scales however remains largely unknown. To understand this, we investigate how different AGN feedback models in the fable simulation suite affect the cosmic evolution of the matter power spectrum (MPS). Our analysis reveals that AGN feedback significantly suppresses clustering at scales |$k \sim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, with the strongest effect at redshift |$z = 0$| causing a reduction of |$\sim 10~{{\ \rm per\ cent}}$| with respect to the dark matter-only simulation. This is due to the efficient feedback in both radio (low Eddington ratio) and quasar (high Eddington ratio) modes in our fiducial fable model. We find that variations of the quasar and radio mode feedback with respect to the fiducial fable model have distinct effects on the MPS redshift evolution, with radio mode being more effective on larger scales and later epochs. Furthermore, MPS suppression is dominated by AGN feedback effects inside haloes at |$z = 0$|⁠, while for |$z \gtrsim 1$| the matter distribution both inside and outside of haloes shapes the MPS suppression. Hence, future observations probing earlier cosmic times beyond |$z \sim 1$| will be instrumental in constraining the nature of AGN feedback.

1 INTRODUCTION

The underlying cosmology of our Universe dictates the properties and evolution of cosmic structure. One of these is the distribution of mass in our Universe, which has been mapped through both structure formation and late-time surveys (e.g. Heymans et al. 2021; Abbott et al. 2022; Qu et al. 2024), and early Universe observations of the cosmic microwave background (Hinshaw et al. 2013; Ade et al. 2016). Numerical and observational studies have shown that baryonic physics, specifically feedback processes from stars and black holes, may significantly impact the distribution of matter across the cosmic web (e.g. Seljak 2000; van Daalen et al. 2011; Chisari et al. 2018; Secco et al. 2022). With programmes like the Vera C. Rubin Observatory (Ivezić et al. 2019), the Euclid Space Telescope (Laureijs et al. 2011), and the Roman Space Telescope (Spergel et al. 2013) preparing to map the Large-Scale Structure (LSS) with unprecedented accuracy, precise theoretical model predictions are urgently required to understand the processes that shape the distribution of galaxies and the underlying matter across cosmic time.

Baryonic feedback plays a crucial role in the formation of individual galaxies, where supernovae and AGN activity have been identified as some of the key processes (e.g. White & Frenk 1991; Sijacki et al. 2007; Hopkins et al. 2014; Habouzit, Volonteri & Dubois 2017; Rosdahl et al. 2018; Trebitsch, Volonteri & Dubois 2020). Such feedback is required to reconcile local observations with small-scale challenges to our lambda cold dark matter (⁠|$\Lambda$|CDM) model (Bullock & Boylan-Kolchin 2017) as well as to produce a realistic global star formation history (Madau & Dickinson 2014), or massive quenched elliptical galaxies and the brightest cluster galaxies (McNamara & Nulsen 2007; Fabian 2012). Through galactic outflows and AGN-driven winds, these feedback processes also provide a channel for galaxies to interact with their larger-scale environment and the local distribution of matter, and influence statistics such as the matter power spectrum (MPS) (Daalen et al. 2011; Chisari et al. 2019). Moreover, a wide range of not-so-well-understood baryonic feedback processes such as cosmic rays or magnetism are gaining popularity in the realistic modelling of outflows from galaxies (e.g. Pakmor et al. 2016; Girichidis et al. 2018; Hopkins et al. 2020; Martin-Alvarez et al. 2020, 2023; Beckmann et al. 2022; Farcy et al. 2022; Rodríguez Montero et al. 2024), and may significantly affect how galaxies shape the local distribution of matter.

AGN feedback is the main process regulating the evolution of the most massive galaxies, galaxy groups, and galaxy clusters (Sijacki et al. 2007; Cavagnolo et al. 2010; Bourne & Sijacki 2017; Chisari et al. 2018; Beckmann, Devriendt & Slyz 2019; Bourne & Yang 2023). The outflows driven by AGN are extremely energetic and can reach scales up to |$\sim$|Mpc, making this form of feedback the most important for the cosmic distribution of matter (e.g. Daalen et al. 2011; Mead et al. 2015; McCarthy et al. 2018; Chisari et al. 2019), primarily through redistribution of matter within and beyond the largest haloes (e.g. van Daalen & Schaye 2015; van Daalen, McCarthy & Schaye 2020; van Loon & van Daalen 2023). While the influence of AGN feedback can be captured through simple halo models (Seljak 2000; Mead et al. 2021), due to the complex relationship between the small-scale regulation of accretion onto supermassive black holes (SMBHs), galaxy formation physics, and the large-scale effects of AGN feedback, cosmological simulations are required to understand its effect on the MPS.

Multiple studies employing some of the largest and most sophisticated cosmological simulations to date (e.g. Vogelsberger et al. 2014c; Hellwing et al. 2016; Chisari et al. 2018; Springel et al. 2018; Daalen et al. 2020; Sorini et al. 2022; Loon & van Daalen 2023; Schaye et al. 2023; Gebhardt et al. 2024) have established that AGN feedback affects the MPS at scales |$k \gtrsim 0.5 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. The resulting power suppression with respect to the dark matter non-linear prediction in these models reaches up to |$\sim 20~{{\ \rm per\ cent}}$| (Chisari et al. 2019). While these different simulations display similar qualitative behaviour, quantitative differences across results are significant, emerging from different feedback implementation strategies and configurations, as well as from different model resolutions and numerical solvers.

To better comprehend the discrepancies between different simulations, a more detailed understanding of how different AGN feedback models affect the MPS is required, and how this impact emerges around different galaxies and environments. AGN feedback is an inherently multiscale phenomenon, spanning from event horizon and accretion disc scales at which the feedback (in the form of radiation, winds, and jets) is produced, out to scales beyond the host galaxy itself. As such, modelling this process in full is virtually impossible within a single simulation. Instead, cosmological simulations have to employ sub-grid models that can capture the effects of AGN feedback and how it couples to baryons at resolvable scales. These models can vary in their sophistication and their made assumptions. The simplest approach is direct thermal energy injection into cells or particles close to the black hole (Springel, Di Matteo & Hernquist 2005; Booth & Schaye 2009; Schaye et al. 2015; Tremmel et al. 2019), often combined with numerically motivated modifications, such as minimum heating temperatures (Booth & Schaye 2009; Schaye et al. 2015; McCarthy et al. 2017), fixed duty cycles (Henden et al. 2018; Koudmani, Sijacki & Smith 2022) or artificial prevention of radiative cooling (Tremmel et al. 2017, 2019) in order to avoid overcooling (see discussions in Bourne, Zubovas & Nayakshin 2015; Crain et al. 2015; Schaye et al. 2015; Zubovas, Bourne & Nayakshin 2016). Other models inject momentum to surrounding gas as bipolar wind or jet-like outflows (Dubois et al. 2014; Weinberger et al. 2018; Davé et al. 2019), with several works including separate quasar and radio mode phases that use different energy injection schemes for each (Sijacki et al. 2007, 2015; Dubois et al. 2014, 2021; Henden et al. 2018). Simulations are additionally performed over a wide range of resolutions, which itself can impact the range of gas phases captured and how feedback couples to these different phases (e.g. Bourne et al. 2015; Beckmann et al. 2019; Koudmani et al. 2019; Talbot, Sijacki & Bourne 2024; Hopkins et al. 2024a). Taking this into account, as well as the use of different codes to perform cosmological simulations, model parameters are typically calibrated to match low-redshift observables such as the galaxy stellar mass function and BH scaling relations (Dubois et al. 2014; Schaye et al. 2015; Sijacki et al. 2015; Pillepich et al. 2018) meaning that different feedback models, in different codes and at different resolutions can result in comparable galaxy populations. As such it is the galaxy properties to which simulation parameters are not tuned that can be used to differentiate between models.

One such quantity is the baryon content of groups and clusters, which has been suggested as a proxy for the expected suppression in the MPS (Semboloni et al. 2011; Semboloni, Hoekstra & Schaye 2013; McCarthy et al. 2018; Schneider et al. 2019; Daalen et al. 2020; Debackere, Schaye & Hoekstra 2020; Salcido et al. 2023). The AGN model in the original Illustris suite of simulations was too effective at expelling gas from groups and low-mass clusters (Genel et al. 2014), and indeed, the MPS suppression found in Illustris is more extreme than that found in other simulations that retain higher baryon fractions (Chisari et al. 2019; Daalen et al. 2020). The fable simulation suite remedied this problem by modifying the feedback models employed in Illustris, making the quasar mode more effective and the radio mode less explosive (Henden et al. 2018) to achieve a better agreement to observations of group and cluster baryon content. In determining their fiducial AGN model, other variations were performed with a total of four presented in appendix A of Henden et al. (2018), which result in different present-day stellar and gas fractions in groups and clusters. These variations provide an ideal testbed to study the effect of different AGN feedback models on the MPS, which provides a key motivation for the work presented here.

We describe the fable simulations in Section 2.1, and our procedure to MPS in Section 2.2. Our main results are explored in Section 3, where we compare various AGN feedback models (Section 3.2), comparing its effect on fable with previous simulations (Section 3.3). We explore in more detail how feedback effects vary around galaxies under different selections (halo mass, stellar mass, and black hole mass) in Section 3.4. Section 3.5 briefly reviews how different halo mass components trace the MPS suppression from AGN at different scales and times. In Section 4, we review the main caveats of our study, mainly stemming from the size of the fable computational box employed. Finally, we conclude this manuscript in Section 5 with a summary of our work.

2 NUMERICAL METHODS

2.1 fable simulations

In this section, we provide a brief summary of the fable simulation suite (Henden et al. 2018; Henden, Puchwein & Sijacki 2019, 2020), which we employ for our investigation into the impact of AGN feedback models on the MPS and the galaxy bias. For a detailed description of the fable set-up and the calibration of the simulations see Henden et al. (2018).

2.1.1 Basic simulation properties

The fable simulations were performed with the arepo code (Springel 2010), where the equations of hydrodynamics are solved on a moving unstructured mesh defined by the Voronoi tessellation of a set of discrete points which (approximately) move with the velocity of the local flow. The gravitational interactions are modelled via the TreePM method with stars and DM represented by collisionless particles.

The fable simulation suite comprises cosmological volumes as well as zoom-in simulations of groups and clusters. Here, we focus on the cosmological volume simulations to investigate the clustering of matter at large scales (rather than examining individual haloes). These 40 |$h^{-1} \, \mathrm{Mpc}$| (⁠|$h = 0.679$|⁠) boxes are evolved using initial conditions for a uniformly sampled cosmological volume based on the Planck cosmology (Planck Collaboration XIII 2016) with |$512^3$| DM particles, yielding a resolution of |$m_\mathrm{DM} = 3.4 \times 10^7 ~ h^{-1}\, \mathrm{{\rm M}_{\odot }}$|⁠, and initially |$512^3$| gas elements with target gas mass resolution |$\overline{m}_\mathrm{gas} = 6.4 \times 10^{6} ~ h^{-1}\, \mathrm{{\rm M}_{\odot }}$|⁠. The gravitational softening is set to |$2.393 ~ h^{-1} \, \mathrm{kpc}$| in physical coordinates below |$z=5$| and held fixed in comoving coordinates at higher redshifts. Notably, this leads to a suite of high-resolution simulations, albeit with a comparatively small cosmological volume when addressing cosmological stastics such as the MPS. We discuss the main caveats resulting from this in Section 4.

The fable galaxy formation model is based on Illustris (Vogelsberger et al. 2013; Genel et al. 2014; Torrey et al. 2014; Vogelsberger et al. 2014a; Sijacki et al. 2015), with the prescriptions for radiative cooling (Katz, Weinberg & Hernquist 1996; Wiersma, Schaye & Smith 2009a), uniform ultraviolet background (Faucher-Giguère et al. 2009), chemical enrichment (Wiersma et al. 2009b), and star formation (Springel & Hernquist 2003) unchanged from the Illustris model. The stellar and AGN feedback prescriptions, on the other hand, are modified to improve agreement with the present-day galaxy stellar mass function and to match the gas mass fractions in observed massive haloes.

2.1.2 Stellar feedback

In the Illustris galactic wind model (Vogelsberger et al. 2013), wind particles are launched from star-forming regions driven by the available energy from core-collapse SNe.

This model is also adopted in fable with a few modifications to the parameters that govern the wind energetics. Specifically, the wind energy factor |$\epsilon _\mathrm{W,SN}$|⁠, which gives the fraction of energy available from each core collapse supernova, is increased to |$\epsilon _\mathrm{W,SN} = 1.5$| in fable compared to the Illustris value of |$\epsilon _\mathrm{W,SN} = 1.09$|⁠. Furthermore, one-third of the wind energy is injected as thermal energy in fable, whilst in Illustris, the stellar-feedback-driven winds are purely kinetic. Overall, this leads to more energetic stellar feedback which more efficiently dissipates the released energy to the gas, and somewhat more effectively regulating star formation in low-mass haloes [see Henden et al. 2018 for details; the same method is used by Marinacci, Pakmor & Springel (2014)].

2.1.3 Black hole seeding and growth

BHs are modelled as collisionless particles and are seeded into DM haloes above a mass threshold of |$5 \times 10^{10} \ h^{-1} \, \mathrm{{\rm M}_{\odot }}$| with a seed mass of |$M_\mathrm{BH,seed}= 10^{5} \ h^{-1}\, \mathrm{{\rm M}_{\odot }}$|⁠.

Subsequently, these BHs may grow via BH–BH mergers and gas accretion following the Eddington-limited Bondi–Hoyle–Lyttleton accretion rate with boost factor |$\alpha = 100$| (Hoyle & Lyttleton 1939; Bondi & Hoyle 1944; Springel et al. 2005). For all AGN models the radiative efficiency is set to a constant |$\epsilon _\mathrm{r} = 0.1$| and |$(1-\epsilon _\mathrm{r})$| of the accreted mass is added to the BH particle mass at each timestep.

Lastly, we note that the BHs are pinned to the potential minimum of their host halo to prevent spurious BH movement due to numerical heating (see Sijacki et al. 2007; Vogelsberger et al. 2013, for details on the BH seeding and growth models).

2.1.4 AGN feedback

Analogously to Illustris, the AGN feedback in fable is based on a two-mode model, with the quasar mode operating at high Eddington ratios (see Di Matteo, Springel & Hernquist 2005; Springel et al. 2005) and the radio mode being activated at low Eddington ratios (see Sijacki et al. 2007). For the fiducial fable simulation set-up, this switch occurs at an Eddington ratio of |$f_\mathrm{Edd,QM} = 0.01$| (compared to |$f_\mathrm{Edd,QM} = 0.05$| in Illustris).

In the quasar mode, a fraction |$\epsilon _\mathrm{f}=0.1$| of the AGN luminosity is isotropically injected as thermal energy. In Illustris, this thermal energy injection happens continuously, which can lead to artificial overcooling as small amounts of energy are distributed preferentially into the densest material over a large gas mass due to the limited gas mass resolution. In fable, this issue is alleviated by introducing a duty cycle with an approach similar to that of Booth & Schaye (2009), whereby thermal energy is accumulated over |$\delta t_\mathrm{QM} = 25$| Myr before being released in a single event, allowing high feedback temperatures, and hence longer cooling times, to be reached. Such a feedback cycle is also, at least qualitatively, consistent with episodic accretion observed in high-resolution simulations (Ciotti, Ostriker & Proga 2010; Torrey et al. 2017; Costa et al. 2018).

In the radio mode, the feedback energy is coupled to the gas as hot buoyantly rising bubbles to mimic those inflated by jets (McNamara & Nulsen 2007; Fabian 2012; Bourne & Yang 2023), with the duty cycle of these bubble injections set by the fractional BH mass growth |$\delta _\mathrm{BH} = \delta M_\mathrm{BH}/M_\mathrm{BH}$|⁠. In fable, this threshold is set to |$\delta _\mathrm{BH} = 0.01$| – much smaller than the Illustris value of |$\delta _\mathrm{BH} = 0.15$|⁠. The bubble energy content is determined as |$\epsilon _\mathrm{m} \epsilon _\mathrm{r} \mathrm{c^{2}} \delta M_\mathrm{BH}$| with the radio mode coupling efficiency set to |$\epsilon _\mathrm{m} = 0.8$| in the fiducial fable model. This yields a similar effective radio mode efficiency |$\epsilon _\mathrm{m} \epsilon _\mathrm{r} = 0.08$| as in the Illustris model (where the effective radio mode efficiency is set to 7 per cent). The lower |$\delta _\mathrm{BH}$| then results in more frequent and less energetic bubbles in fable compared to the Illustris set-up.

2.1.5 AGN model variations

In addition to the fiducial fable model, Henden et al. (2018) explore three additional AGN feedback parametrizations:

1. The RadioStrong set-up, which has the same radio mode parameters as the fiducial run but no quasar duty cycle.

2. The Quasar set-up, which employs a quasar duty cycle but has significantly weaker radio mode feedback with a lower threshold for bubble injections (⁠|$\delta _\mathrm{BH} = 0.001$|⁠) and a lower coupling efficiency (⁠|$\epsilon _\mathrm{m} = 0.4$|⁠).

3. The RadioWeak set-up, which does not have a quasar duty cycle and employs the weaker radio mode feedback.

Together with the fiducial run, these three alternative AGN set-ups then allow us to isolate the impact of the quasar duty cycle and increasing the strength of the radio mode feedback. Note that all of the additional runs also have a higher critical Eddington fraction for switching to the quasar mode (⁠|$f_\mathrm{Edd,QM} = 0.05$|⁠, as in Illustris).

Furthermore, we also analyse the results from an additional fable model variation, NoAGN, which was performed without seeding any black holes therefore providing a useful reference run without any AGN feedback.

The four AGN runs and the no-AGN run form the core of our analysis and the corresponding AGN parameters for these five set-ups are listed in Table 1. For reference, the corresponding parameters for the original Illustris simulation set-up are also given.

2.1.6 Halo and galaxy identification

For our analysis, we identify DM haloes (‘groups’) and galaxies (‘subhaloes’) via the friends-of-friends (FoF) and subfind algorithms (Davis et al. 1985; Springel et al. 2001; Dolag et al. 2009), respectively. The FoF search linking length is set to 0.2 times the mean particle separation. Within the FoF groups gravitationally bound systems are identified as subhaloes, as found by subfind. The central subhalo corresponds to the subhalo at the minimum potential of the FoF group whilst all other subhaloes in the same group are categorized as satellites.

We characterize the galaxy properties employing the total stellar mass of each subhalo as the stellar mass of each galaxy |$M_*$|⁠, the central black hole mass in each subhalo as |$M_\text{BH}$|⁠, and the subhalo total mass as |$M_\text{halo}$|⁠. For the NoAGN and DMO models, we estimate the |$M_\text{BH}$| and |$M_{*}$| of each galaxy interpolating its |$M_\text{halo}$| in the Fiducial model relationships. These relationships are presented and described in Section 3.4.

2.2 Power spectra

We focus our investigation on the matter power spectrum, which provides information of the matter clustering at different scales, studied here in Fourier space and characterized by a wavelength k in units of |$\,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. In order to extract from the fable simulations the studied matter power spectra and cross-correlations between quantities, we make use of the fftw library.¹ We project the entire computational domain of the simulation onto a uniform grid of |$1024^3$| cells. Hence, the computational domain with a physical size of |$L_\text{FFT} \sim 40\, \text{h}^{-1}\, \text{cMpc}$| (⁠|$59\, \text{cMpc}$|⁠) is resolved down to |$dx_\text{FFT} \sim 39 \, \text{h}^{-1}\, \text{ckpc}$| (⁠|$58\, \text{ckpc}$|⁠). Consequently, each of our spectra spans from |$k_\text{min} \sim 0.2\,\,\text{h}\,\,\text{cMpc}^{-1}$| to |$k_\text{max} \sim 74\,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, although we note that our results are affected by the limited simulation volume, outlined in Section 4, particularly at |$k_\text{max} \lesssim 1\,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. To obtain the MPS, we project onto a 3D grid all the particles included in fable (i.e. dark matter, stars, gas, and black holes), employing a simple nearest grid point (NGP) interpolation. Finally, whenever performing a cross-correlation between two scalar fields, we compute this through their multiplication in Fourier space, finally averaging onto a 1D k-space binning.

3 RESULTS

3.1 A qualitative comparison of AGN impact around massive galaxies

As the hot gas ejected by AGN feedback escapes from galaxies and expands against the circumgalactic and intergalactic medium, it leads to the ejection of gas from the densest regions of the cosmic web, reducing the amount of clustering on the smallest scales. To provide a qualitative visualization of this effect, we display in Fig. 1 overdensity and density contrast maps, with each of the studied simulations corresponding to a different column. The first row of panels shows RGB projections where the overdensity at different scales is represented in colours: large (red; |$k < 11 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠), intermediate (green; |$11 \,\,\text{h}\,\,\text{cMpc}^{-1}< k < 22 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠), and small (blue; |$k > 22 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠). Such separation is generated in 3D Fourier space, through a band-pass filter that isolates a specific range of scales applied to the entire simulated domain. The resulting remaining power is then converted back to configuration space, and the corresponding overdensity field employed to generate the projections. As the baryonic impact amounts to proportions of no more than |$\sim 20~{{\ \rm per\ cent}}$| in most simulations (e.g. the set compiled by Chisari et al. 2019), the overdensity projections in this first row show only subtle differences across simulations.

To explore the impact of different AGN feedback models, we focus on relative changes between simulations, which are normalized with respect to the dark matter only (DMO) fable model (Fig. 1). Accordingly, the bottom row of Fig. 1 shows the overdensity contrast for the gas mass with respect to the DMO model, |$\mathcal {D}_\text{gas, DMO}$|⁠, calculated as:

where |$\delta _\text{component}$| is the overdensity of a given component, computed for a given simulation sim with respect to a reference model ref. In Fig. 1, positive values (shown in blue) indicate that the baryonic gas has a higher overdensity than the (total) overdensity of the DMO case, whereas negative values (shown in red) indicate a lower overdensity instead. The colour scales are fixed equally for all panels in each row, and we separate large scales in the bottom row, following the Fourier low-k-pass filter method outlined above. We include the underlying gas overdensity distribution in grey to guide the eye. All models show some increase of power with respect to the DMO simulation within densest nodes of the cosmic web. Such denser structures are primarily driven by baryonic cooling. The NoAGN panel (rightmost column; note that NoAGN does not include AGN feedback but still has supernova stellar feedback), illustrates how baryonic cooling and SN feedback suffice to drive some local mild power suppression with respect to the DMO scenario. However, models with strong AGN feedback show clear circular red structures around massive galaxies. These are associated with a considerable reduction of power at large scales, driven by AGN activity evacuating matter towards larger radii. The Quasar model appears to mostly enhance events of bi-directional power suppression. By combining the duty cycle (Quasar) and increased radio mode strength (RadioStrong) modifications to the AGN model, the Fiducial simulation has an enhanced suppression of power, where both the characteristic large-scale ring-like structures from a stronger AGN and the bi-channel ejection of the duty-cycle are intensified. We note that the isotropic or anisotropic impact of the feedback is also driven by the environment impacted, the scales reached by the AGN feedback, and even the redshift when the effect takes place (see Section 3.5).

Overall, the presence of ring-like structures for the efficient AGN models in such large-scale projections (e.g. central object in top panels), illustrates how AGN feedback is responsible for re-shaping the distribution of matter around the clusters and galaxies of the cosmic web. These structures suggest an approximately isotropic displacement of gas for the strongest AGN models (Fiducial and RadioStrong), and more anisotropic effects for the duty-cycle model (Quasar). We will show below how simulations with these spherical matter ejections display a larger power suppression at intermediate cosmological scales, whereas our duty cycle model is more efficient in driving small-scale effects, and has a more complex redshift evolution. The lack of significant power suppression in the NoAGN and RadioWeak models suggests that in the absence of efficient AGN models, other baryonic physics may only have marginal effects.

3.2 The impact of baryons on the matter power spectrum of the fable simulations

The main statistic of interest to understand how baryonic feedback affects the distribution of matter is the MPS, which we compute as described in Section 2.2. The resulting MPS for some representative fable models at |$z = 0$| are presented in the top panel of Fig. 2, which includes the camb predictions for this redshift, employing the Boltzmann solver (with |$k_\text{h min} = 10^{-4} \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, |$k_\text{h max} = 7 \cdot 10^{1} \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, using 1000 points) and for the employed fable cosmology (see Section 2.1). For the non-linear prediction, we employ the Halofit model version by Mead et al. (2021). The full-physics fable MPS significantly deviates from the non-linear prediction at scales |$k \gtrsim 30 \,\,\text{h}\,\,\text{cMpc}^{-1}$| due to clustering from the collapse of baryons into haloes.

The impact of baryonic feedback on the MPS is typically concentrated on relatively small cosmological scales |$k \sim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| (Chisari et al. 2019). In order to study such impact in fable in more detail, we show in Fig. 2 the |$z = 0$| ratio of the different studied models with respect to the DMO case (central panel) and with respect to the NoAGN model (bottom panel). We find the baryonic feedback in fable to only have a significant effect for scales |$k > 1 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, (see also Section 4). The clustering effect of baryons dominates for scales |$k \gtrsim 30 \,\,\text{h}\,\,\text{cMpc}^{-1}$| in all our models. Amongst all the simulations with AGN feedback, the RadioWeak model shows the lowest suppression of power, with an MPS almost equivalent to that of the NoAGN model. The Quasar and RadioStrong models both significantly suppress the amount of baryonic clustering. This leads to a large deviation from the NoAGN simulation, reaching an MPS comparable to the DMO case at scales |$k \lesssim 20\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. These two models also show an interesting power cross-over at scales of approximately |$k \sim 12 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, where the Quasar model yields a lower power suppression at scales larger than this cross-over. The AGN duty cycle in this model concentrates feedback into periodic bursts, and AGN activity appears particularly efficient at scales |$k \gtrsim 5 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. On the other hand, increasing the radio mode feedback strength (RadioStrong model) smoothly suppresses the MPS on scales of |$k \sim 2\!-\! 10 \,\text{h}\,\,\text{cMpc}^{-1}$|⁠. When both a stronger radio mode and quasar duty cycle are combined in the Fiducial model, a maximal power suppression of |$\sim 10~{{\ \rm per\ cent}}$| below the DMO scenario is reached, with the suppression peaking at scales |$k \sim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. As expected, the largest relative difference in power between two models with AGN feedback at this suppression peak occurs between the RadioWeak and Fiducial .

Interestingly, the Fiducial power suppression has features characteristic from both models, such as the plateau in power for |$k \in \left[3, 7\right]\,\,\text{h}\,\,\text{cMpc}^{-1}$| from RadioStrong and a dip in power at |$k \sim 15\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| from Quasar . To understand whether the modifications with respect to RadioWeak in each of the two models can be directly combined to recover the suppression observed in Fiducial, we include in the bottom panel of Fig. 2(a) line for RadioStrong |$*$| Quasar (pink dashed line). This is computed by removing from the RadioWeak |$P_\text{mm}$| the differences in power between both RadioStrong and RadioWeak, and Quasar and RadioWeak . Overall, ‘RadioStrong |$*$| Quasar’ traces the Fiducial model well, with only a slight underestimation of the suppression towards large scales (⁠|$k < 5\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠), and a suppression overestimate at |$k \gtrsim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. This hints towards an independent impact of applying a quasar duty cycle and increasing the radio mode efficiency, at least within the variation studied in fable. We will show below that this is due to the two modes being active and effective at different redshifts.

To better understand how AGN feedback, through the two separate modes, progressively carves its impact on the MPS, Fig. 3 displays the redshift evolution of this quantity between |$z = 8$| and |$z = 0$| in the four studied AGN models. At high redshifts, |$P_\text{mm}$| remains similar to |$P_\text{mm, DMO}$|⁠, but as the simulations evolve, clustering is increased at small scales and AGN feedback progressively leads to power suppression at the intermediate scales. In the two simulations without the quasar duty cycle (i.e. RadioWeak and RadioStrong; right column), baryonic cooling leads to higher clustering than in the DMO case down to |$z \sim 0.5$|⁠, with the models being comparable at |$z \gtrsim 1$|⁠. After this redshift, the RadioStrong model undergoes a significant power suppression, particularly prominent during the |$z \in \left(0.5, 0.0\right]$| interval. The RadioWeak simulation evolves to closely resemble the NoAGN case at |$z = 0$|⁠, with its maximal relative deviation from it at |$z \sim 0.5$|⁠, and with only a mild suppression afterwards.

On the other hand, both models with the quasar duty cycle (i.e. Quasar and Fiducial; left column) display an early suppression of power from |$z \sim 3$| onwards. The early impact takes place at larger comoving scales and progressively shifts towards smaller scales as the simulation evolves, to eventually reach the peak of suppression observed at |$k \sim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$| at |$z = 0$|⁠. Despite this, the AGN feedback in the Quasar model is unable to maintain suppression over the DMO case after |$z \sim 1$|⁠, and develops the noticeable peaks in power (⁠|$k \sim 7 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) and suppression (⁠|$k \sim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) after |$z \sim 0.5$|⁠. The Fiducial model continuously builds up relative suppression with respect to the DMO model during the |$z \sim 3$| to |$z = 0$| interval, with most of the deviation taking place after |$z \sim 0.5$| as observed in the RadioStrong case.

The differences in the redshift evolution of the matter power spectrum for these different models are primarily driven by the temporal evolution of the quasar mode and radio mode fractions (AGN mode fractions are shown in Fig. A1, discussed in Appendix  A), with the quasar mode dominating at high redshifts and the radio mode becoming increasingly more important towards low redshift. The quasar duty cycle in combination with strong radio mode feedback, as in the Fiducial model, then ensures efficient AGN feedback injection throughout cosmic history. Consequently, the effective quasar mode will lead to an earlier power suppression (⁠|$z > 1$|⁠) whereas the radio mode will be important at late times (⁠|$z < 2$|⁠). These introduce an evolutionary degeneracy that should be addressed as upcoming observatories such as the Simons Observatory probe |$z \gtrsim 1$| (Ade et al. 2019).

Along these considerations, while matter clustering in the RadioWeak model is always above the DMO case, and comparable to the NoAGN simulation, all the other models have |$P_\text{mm}/P_\text{mm, DMO}< 1$| at some point during their evolution in the |$3 \ge z \ge 0$| interval. The largest power suppression for both Fiducial and RadioStrong models takes place at |$z = 0$|⁠, whereas it takes place during |$z \sim 2\!-\!1$| for the Quasar simulation. When comparing the strongest and weakest AGN feedback models (top row; Fiducial versus RadioWeak) their relative suppression of power at |$k \sim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$| is of the order of 40 per cent, and primarily develops at |$z \lesssim 1$|⁠. The bottom row illustrates how, at small scales (⁠|$k > 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) and after |$z \sim 3$|⁠, the Quasar AGN model has a considerably higher suppression of power than the RadioStrong model (bottom row). Finally, all AGN simulations, excepts perhaps RadioWeak, have a complex evolution of their MPS with respect to the DMO case.

3.3 Comparing fable with other simulations

In Fig. 4, we compare the relative impact of baryons on the MPS (through |$P_\text{mm}/P_\text{mm, DMO}$|⁠) of fable with their relative impact in other well-known cosmological galaxy formation simulations.² For each of the models, we include in parenthesis the cosmological box size, with different lines coloured according to the maximum spatial resolution of the respective simulation. Overall, the suppression of the relative MPS in simulations due to baryonic physics, primarily due to AGN feedback, occurs at scales of |$k \sim 5\!-\! 20 \,\text{h}\,\text{cMpc}^{-1}$|⁠. Focusing on the fable Fiducial model, the scale of maximal power suppression |$k_\text{peak} \sim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$| is comparable to most other simulations (typically |$k_\text{peak} \sim 10\!-\!20 \,\text{h}\,\,\text{cMpc}^{-1}$|⁠) and the overall shape of the relative MPS lies within the bulk of the outcomes from other simulation projects.

We first compare fable with the HorizonAGN and IllustrisTNG-100 simulations, as they have comparable order of magnitude resolutions, finding that they all have similar maximum power suppression of |$\sim 10\!-\!20~{{\ \rm per\ cent}}$|⁠. Interestingly, fable has a larger impact at scales |$1 < k/\,\,\text{h}\,\,\text{cMpc}^{-1}< 5$| than IllustrisTNG-100 and especially HorizonAGN, with the latter finding their relative power suppression to be concentrated at |$k \gtrsim 5 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. Instead, AGN feedback in fable has a shape of the |$P_\text{mm}/P_\text{mm, DMO}$| curve that resembles that of the OWLS or BAHAMAS simulations, but with a lower suppression magnitude. This resemblance is possibly the result of all three of these simulations employing episodic quasar mode feedback.

When comparing alternative fable AGN feedback physics with other simulations, the NoAGN and the weakest AGN (RadioWeak) models resemble the behaviour of the HorizonNoAGN simulation, although the fable cases have a less pronounced clustering at |$k \gtrsim 4 \,\,\text{h}\,\,\text{cMpc}^{-1}$| probably due to stronger SN feedback. This illustrates the weak impact of RadioWeak, despite featuring relatively unchanged galaxy populations (see Henden et al. 2018; or Section 3.4). Instead, RadioWeak and NoAGN models have higher gas mass fractions within the largest haloes, in some tension with observations. The Fiducial and RadioStrong cases are in good agreement with observations, (see fig. A2 in Henden et al. 2018), and suggest radio mode AGN feedback in fable also leads to the discussed correlation between cluster gas mass fractions and MPS power suppression (Daalen et al. 2020).

Another interesting comparison is that between our Quasar model and MassiveBlack, where both simulations show a somewhat flatter MPS, with a modest peak in clustering at |$k \sim 12 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. MassiveBlack employs a thermal feedback prescription with a constant energy parameter of |$f = 0.05$|⁠, and does not include an alternative feedback injection mechanism in the radio mode regime (Khandai et al. 2015). Consequently, it is possible that the comparable behaviour is the result of feedback being predominantly more effective at higher redshifts, whilst at low redshifts, as the BH accretion rate density decreases, the impact of AGN feedback on the MPS declines significantly. In agreement with Chisari et al. (2018), we attribute this behaviour in Quasar (and possibly in MassiveBlack) to a late-time decrease in AGN regulation (due to the lack of effective’radio-mode’ feedback), where power builds up more rapidly at intermediate scales (⁠|$1 < k/\,\,\text{h}\,\,\text{cMpc}^{-1}< 10$|⁠) when gas is re-accreted into massive haloes (Beckmann et al. 2017; Habouzit et al. 2021).

Finally, we note that the fable AGN model has an impact well below the |$\sim 40~{{\ \rm per\ cent}}$| suppression with respect to the non-linear dark matter-only scenario observed for simulations such as Illustris or OWLS, which likely have ejective feedback that may be too effective (Genel et al. 2014).

3.4 The impact of different galaxies on the matter power spectrum

To understand how different galaxies and AGN feedback from their central SMBHs influences the distribution of matter in our simulations, we explore the variation of the MPS around galaxies separated according to multiple property cuts. We analyse separations according to the halo mass (⁠|$M_\text{halo}$|⁠), stellar mass (⁠|$M_{*}$|⁠) and black hole mass (⁠|$M_\text{BH}$|⁠) of galaxies, computed as described in Section 2.1.6. The baryon fraction |$f_\text{baryon}$| (defined as the fraction of baryonic mass to the total mass inside each halo) is another important quantity to understand how different feedback contributes to MPS power suppression (e.g. Semboloni et al. 2011; Salcido et al. 2023). We show the distribution of fable galaxies across the parameter space of these properties in Fig. 5. From top left to bottom right, each panel displays the |$f_\text{baryon}\!-\!M_\text{halo}$|⁠, |$f_\text{baryon}\!-\!M_\text{BH}$|⁠, |$M_\text{halo}\!-\!M_{*}$|⁠, and |$M_\text{BH}\!-\!M_{*}$| relations, respectively. We rescale our |$f_\text{baryon}$| measurements by |$\Omega _m/\Omega _b$| to display the proportion with respect to the universal baryon fraction in our simulations. Solid lines in each of the panels display the median relation for each of the shown models, whereas the shaded bands encompass the 20–80 per cent quantiles.

Overall, fable runs have relatively similar mean relations for the simulated galaxies, with the most important variations taking place with |$M_\text{BH}$|⁠. Simulations with weaker feedback reach higher |$M_\text{BH}$|⁠, but the effects of the duty cycle and increased radio mode strength affect differently systems across various |$M_\text{BH}$| ranges. This effect appears more prominent on the low |$M_{*}$| boundary of the |$M_\text{BH}\!-\!M_{*}$| relation (also see Koudmani, Henden & Sijacki 2021). We refer the reader to Henden et al. (2018) for further analysis of the galaxy and cluster populations in fable. When reviewing the |$f_\text{baryon}$| changes with respect to |$M_\text{halo}$|⁠, we find the relative baryon content to decrease above masses of |$\sim$||$10^{12}\, \mathrm{{\rm M}_{\odot }}$|⁠, and remain relatively flat for the stronger AGN feedback models. As expected, the relation of this same quantity with respect to the |$M_\text{BH}$| reflects the AGN-driven reduction of baryonic content in these haloes, as we find a clear decrease of |$f_\text{baryon}$| with increasing |$M_\text{BH}$| for increasingly stronger AGN feedback models.

Fig. 5 also includes dashed lines corresponding to the cuts employed to separate our galaxies into different mass ranges. We present only a subset of all the investigated cuts, varying these cuts results only in monotonic and minor variations. Our mass range divisions are first set to separate halo masses (i.e. |$M_\text{halo},_\text{cut} \in \left[7.5 \cdot 10^{10}, 2.5 \cdot 10^{11}, 7.5 \cdot 10^{11}, 2.5 \cdot 10^{12}\right] \mathrm{{\rm M}_{\odot }}$|⁠). These |$M_\text{halo},_\text{cut}$| are then converted to |$M_{*},_\text{cut}$| and |$M_\text{BH},_\text{cut}$| following the population medians of the Fiducial fable scalings. The intersects across such scalings are shown as orange points in Fig. 5, and the resulting range division values correspond to |$M_{*},_\text{cut} = \left[1.7 \cdot 10^8, 2.5\cdot 10^9, 1.5\cdot 10^{10}, 6.0\cdot 10^{10}\right] \mathrm{{\rm M}_{\odot }}$| and |$M_\text{BH},_\text{cut} = \left[0.0,\, 2.5\cdot 10^5, 1.4\cdot 10^6, 1.8\cdot 10^7\right] \mathrm{{\rm M}_{\odot }}$|⁠, where the first cut in |$M_\text{BH}$| separates galaxies with and without an SMBH. These divisions serve to review comparable populations selected under different criteria.

To illustrate the effects of filtering our mass distribution, Fig. 6 presents a comparative analysis of the MPS inside and outside of haloes (⁠|$M_\text{halo}> 7.5 \cdot 10^{10}\, \mathrm{{\rm M}_{\odot }}$|⁠), along with twice their cross-correlation (i.e. |$2\, \mathcal {C}_\text{filt,exc}= P_\text{mm}-(P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}+ P_\text{mm}^\text{excluded})$|⁠). The volume outside of haloes constitutes most of the simulation domain, and altogether with the cross-correlation, dominates the contribution to |$P_\text{mm}$| at |$k < 2\,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. At these large scales, haloes contribution to |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}/P_\text{mm, DMO}$| is of order of |$\sim 0.2$|⁠.

Across different AGN models, variations in |$P_\text{mm}^\text{excluded}$| and |$\mathcal {C}_\text{filt,exc}$| are relatively minor, with only a slightly higher amount of power outside of haloes in Fiducial and RadioStrong . Models with stronger radio mode AGN feedback also display somewhat lower cross-correlation in the |$10 \lesssim k/(\,\,\text{h}\,\,\text{cMpc}^{-1}) \lesssim 4$| range, with this quantity rapidly approaching zero towards smaller-scales. It will be important for the analysis below to emphasize that at scales |$k \gtrsim 7\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, the mass outside of haloes and its cross-correlation with |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| only constitutes a small fraction of |$P_\text{mm}$| (⁠|$\lesssim 0.1 P_\text{mm, DMO}$|⁠). Furthermore, differences across models are smaller than the variations observed in |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| or in Fig. 4. Consequently, any suppression observed in |$P_\text{mm}$| emerges from variations of the power within haloes.

We now focus on the variations of the MPS within haloes, selected according the thresholds describe above. Fig. 7 shows the resulting spectra |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$|⁠, for our lowest and highest threshold selections: |$M_\text{halo}> 7.5 \cdot 10^{10}~\mathrm{{\rm M}_{\odot }}$| and |$M_\text{halo}> 2.5 \cdot 10^{12}~\mathrm{{\rm M}_{\odot }}$| (leftmost column), as well as their corresponding thresholds in |$M_{*}$| and |$M_\text{BH}$|⁠. From left to right, columns correspond to cuts on |$M_\text{halo}$|⁠, |$M_{*}$| and |$M_\text{BH}$|⁠. The top set of panels show the filtered spectra divided by the total MPS of the DMO model (⁠|$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}/P_\text{mm, DMO}$|⁠). For a better view of variations across fable models, the central set of panels shows the same mass thresholds, now divided with respect to the halo mass filtered DMO simulation (⁠|$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}/P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}$|⁠). Consequently, the ratio of all panels in a row is computed with respect to the |$P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| in its leftmost column. This allows for a direct comparison across selection masses and fable models. Finally, and to facilitate further comparison with Fig. 2, we show in the bottom row the ratio of |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| with respect to the |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| of the NoAGN simulation.

The differences across the various models, cuts, and selection mass types are more prominent at small scales (⁠|$k > 5\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠). When all galaxies are considered (a row panels) a trend towards higher clustering of the hydrodynamical runs with respect to the haloes in the DMO model is observed. Such transition from lower to higher clustering occurs at different scales for different models, depending on the efficiency of AGN feedback (a1 and c1 panels). For NoAGN, RadioWeak, and Quasar, this occurs at |$k \sim 4 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. For RadioStrong and Fiducial, it takes place at smaller scales, with |$k \sim 12 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$| and |$\sim 25 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, respectively. As haloes are progressively discarded by increasing the threshold mass, the amount of small-scale power is progressively reduced until only the largest haloes are considered (as shown in b row panels). Comparing the least and most restrictive threshold, we find a considerable decrease (⁠|$\sim 0.3 P_\text{mm, DMO}$|⁠) of power at scales |$k \gtrsim 20 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, but a negligible reduction in the |$k \in [5\!-\!10]\,\text{h}\,\text{cMpc}^{-1}$| range. The proportional separation between the |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| of the hydrodynamical models and |$P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| varies differently when increasing mass cuts (c versus d rows), with AGN feedback affecting differently the power clustering contribution from different halo masses (Loon & van Daalen 2023). For example, |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| of Fiducial is reduced from about 20 per cent to 10 per cent over |$P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| at scales of |$k \sim 30\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, whereas RadioWeak increases from 55 per cent over to 65 per cent over |$P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}$|⁠. We attribute this to AGN feedback in this model becoming less capable of impacting more massive galaxies (panels c1 versus d1). Despite their differences, the scale at which these weaker AGN models (i.e. NoAGN, RadioWeak, and Quasar) transition from having less power than the DMO case (⁠|$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}/P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}< 1$|⁠; large scales) to more power than the DMO case (⁠|$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}/P_\text{mm, DMO}^{\text{M}_\text{min} < \text{M}_\text{gal}}>1$|⁠; small scales) remains approximately unchanged at a scale of |$k \sim 4 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. This scale is approximately independent of the mass cut employed. For the models with strongest AGN feedback (RadioStrong and Fiducial), this scale remains unchanged for |$M_\text{halo}$| and |$M_{*}$| sample selections, but shifts to almost a factor |$2\times$| higher scales when the galaxies hosting the most massive black holes are selected (panels c3 versus d3). Applying the most stringent mass cut according to |$M_\text{BH}$| rather than |$M_\text{halo}$| (or |$M_{*}$|⁠) enhances differences between models. This |$M_\text{BH}$|-based selection further increases the power suppression observed due to AGN across all scales (panel d3). This is also in agreement with the expectation of higher power suppression with decreasing |$f_\text{baryon}$|⁠, as selecting galaxies in our Fiducial model with increasingly larger |$M_\text{BH}$| leads to a faster decrease of |$f_\text{baryon}$| than when the analogous selection is applied on |$M_\text{halo}$| (Fig. 5, panels a1 and a2).

Such lack of variations depending on whether |$M_\text{halo}$| and |$M_{*}$| sample selection criteria are employed reflects a tighter interrelation between these two quantities at the high mass end of the galaxy population, which dominate the |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$|⁠. The enhanced relative power suppression when massive SMBH are selected confirms the mass of SMBHs as an integrated measure of AGN feedback and power suppression. However, the specific implementation of AGN feedback may significantly affect the amount of the MPS suppression, motivating further exploration of more sophisticated AGN model implementations (e.g. Zubovas et al. 2016; Bourne & Sijacki 2017; Costa, Pakmor & Springel 2020; Talbot, Bourne & Sijacki 2021; Beckmann et al. 2022; Huško et al. 2022; Koudmani et al. 2023; Rennehan et al. 2023).

When reviewing such differences across models, the trends observed in Fig. 4 are clearly reflected in the |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}/P_\text{mm, DMO}$| ratio at scales |$k \gtrsim 7\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. This is true both for the ratio with respect to |$P_\text{mm, DMO}$| (a and b rows) and with respect to the NoAGN model (e row). As discussed above (Fig. 6), haloes dominate MPS power at such scales. Consequently, the suppression observed in |$P_\text{mm}$| emerges necessarily from the effects of AGN inside haloes. The lower power across all scales in Fiducial (c row panels) and RadioStrong to a lower extent, is the result of clustering being reduced at small scales and matter being ejected outside of haloes (Loon & van Daalen 2023). This is in agreement with the slight increase of power outside of haloes (Fig. 6; top panel) at scales of |$k \sim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. As above, whether systems are selected according to their |$M_{*}$| or their |$M_\text{halo}$| does not have a significant effect on these results. Selecting systems according to the SMBH mass accentuates differences between models, with the Fiducial, RadioStrong, and Quasar models exhibiting a lower amount of power at scales |$k \gtrsim 2\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. The resulting selection of systems preferentially focuses on haloes from which a larger proportion of matter has been removed or ejected.³ Finally, the peak of power at |$k \sim 6\!-\!7\, \,\text{h}\,\text{cMpc}^{-1}$| in the Quasar model is somewhat suppressed by the SMBH selection. We attribute the peak to more inefficient feedback around less massive SMBHs: as the Quasar simulation transitions from a more prevalent quasar mode at high redshift to radio mode dominance (Fig. A1), clustering increases around massive haloes with AGN less efficient in mass removal (Beckmann et al. 2017).

To provide further insight into the effect of feedback around smaller galaxies, we repeat a similar analysis now selecting only galaxies within non-intersecting mass intervals, instead of a minimum mass threshold. This isolates their power contribution, as they are secondary when compared with the most massive systems (Daalen & Schaye 2015). The resulting spectra for different haloes are shown in Fig. 8, employing a |$M_\text{halo}$| selection. The spectra of haloes is shown as the ratio of each model with respect to the NoAGN case. Most models do not have any significant effect on |$P_\text{mm}^{\text{M}_\text{min} < \text{M}_\text{gal}}$| around the smallest galaxies (top panel) except Fiducial . In this simulation, a large suppression takes place across all studied scales, and especially at |$k \lesssim 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. Due to these systems being hosts of small |$M_\text{BH}$|⁠, with low integrated power budgets, the observed power suppression is potentially driven by AGN residing in large neighbouring galaxies affecting their satellites and nearby smaller galaxies (Dashyan et al. 2019; Martín-Navarro et al. 2021; Shuntov et al. 2022). The lack of any significant suppression around intermediate mass haloes (⁠|$2.5 \cdot 10^{11}~\mathrm{{\rm M}_{\odot }}< M_\text{halo}< 2.5 \cdot 10^{12}~\mathrm{{\rm M}_{\odot }}$|⁠; central panel) supports a scenario where only the largest SMBH are capable of such clustering suppression. Once again, only the Fiducial simulation experiences some notable clustering reduction, notable at the smallest scales. This lack of suppression around haloes |$M_\text{halo}< 5 \times 10^{12}\,\mathrm{{\rm M}_{\odot }}$| was also found by Loon & van Daalen (2023). This behaviour remains approximately unchanged regardless of whether an |$M_\text{halo}$|⁠, |$M_{*}$|⁠, or |$M_\text{BH}$| selection is employed.

3.5 Tracing AGN power suppression with haloes at different scales and times

To further understand how baryonic physics modifies clustering inside and outside of haloes across cosmic time, we show the evolution of |$P_\text{mm}(k_\text{scale})/P_\text{mm, DMO}$| in Fig. 9. The top panel displays our results at |$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. At early times (⁠|$z \gtrsim 1$|⁠), the power inside of haloes is the primary contribution to the total |$P_\text{mm}$|⁠, but is not fully dominant. The Quasar and RadioStrong models reveal how clustering is sensitive to the efficiency of the radio versus quasar modes, both inside and outside haloes. The second panel shows the ratio of the difference |$\Delta P = P_\text{mm}-P_\text{mm, halos}$|⁠, which highlights how well is the total MPS traced by |$P_\text{mm, halos}$| for each model. By |$z \sim 1$|⁠, haloes constitute approximately 95 per cent of |$P_\text{mm}$|⁠, and despite significant variations in total power within haloes across models, the proportion of power outside of haloes remains similar across simulations.

Hence, in our hydrodynamical simulations, any large deviations of |$P_\text{mm}$| (or |$P_\text{mm,haloes}$|⁠) from |$P_\text{mm}$| in the NoAGN simulation is primarily driven by AGN feedback modifying the clustering of matter. With the Quasar and Fiducial models evolving comparably, the lower |$P_\text{mm}$| at |$z \sim 0$| due to the higher efficiency of the implemented quasar duty cycle is driven by a suppression of power inside haloes. In the Quasar model, the enhanced impact due to the quasar duty cycle at early times preserves its imprint down to |$z = 0$|⁠. After |$z \sim 0.5$|⁠, the enhanced radio feedback power in the Fiducial and RadioStrong simulations leads to a considerable suppression of halo clustering. This drives a late-time decrease of |$P_\text{mm}$| in both models. At |$z \lesssim 1$|⁠, the difference between the total |$P_\text{mm}$| and that of all haloes remains considerably smaller than the deviations across models, and further confirms the impact of baryons at |$k_\text{scale} \gtrsim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$| to be driven by AGN feedback, and constrained to the interior of haloes.

At large scales (⁠|$k_\text{scale} = 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠; third and fourth panel), and subject to the caveats described in Section 4, redshift evolution of all models resembles their |$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| counterparts. Quasar and Fiducial feature some power suppression until |$z \sim 1$|⁠. However, the effect of the enhanced radio mode is already in place by |$z \sim 1$| (instead of |$z \sim 0.5$| for |$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠). It leads to a suppression of power in Fiducial only slightly above 1 per cent. As shown in Fig. 6, and reflected by the bottom panel, power at this scales is dominated by matter outside of haloes and its cross-correlation with galaxies. With variations in the fraction of power outside of haloes (bottom panel) being larger than the actual suppression, understanding the power decrease with respect to the DMO requires constraining power both inside and outside of haloes. Overall, Fig. 9 illustrates how different AGN implementations will not only modify the details of |$P_\text{mm}$| and |$P_\text{mm, halos}$| at |$z = 0$|⁠, but also its redshift evolution. Hence, any precision measurements (⁠|$< 1~{{\ \rm per\ cent}}$|⁠) will require a more detailed understanding of how AGN feedback operates across cosmic time (Semboloni et al. 2011; Chisari et al. 2019; Huang et al. 2019). For the interested reader, the time evolution of |$P_\text{mm, halos}$| across all studied scales is shown in Appendix  B.

We conclude by revisiting in Fig. 10 how the baryon content in the most massive galaxies, as well as that of the gas and stellar mass, correlates with power suppression at each of these two scales in fable simulations. These are selected as |$M_{*}> 6 \cdot 10^{10}\, \mathrm{{\rm M}_{\odot }}$|⁠, to investigate how well the most massive systems trace the power suppression of our different AGN models. Their power suppression is shown as a function of the total mass in the haloes of all selected systems (⁠|$M_{i, \text{halos}}$|⁠), separately for baryons (circles), gas (squares), and stars (star symbol). We normalize |$M_{i, \text{halos}}$| to their values in the NoAGN simulation (⁠|$M_{i, \text{halos}}/M_{\text{NoFb}, \text{halos}}$|⁠). In addition, we calculate power-law best fits to the stellar measurements (yellow band), gas measurements (green band), and total baryonic mass measurement (blue band). Note that due to too high star formation rates in massive haloes in NoAGN model, its gas content is depleted and is somewhat lower than that for RadioWeak and Quasar simulations. The total baryonic component provides an accurate tracer of power suppression at both |$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| and |$k_\text{scale} = 1\, \text{h}\,\text{cMpc}^{-1}$| (coefficient of determination |$R^2 \sim 0.95$|⁠), in agreement with previous work (e.g. Semboloni et al. 2011; McCarthy et al. 2018; Daalen et al. 2020; Debackere et al. 2020; Salcido et al. 2023). Separate baryonic mass components provide a less tight constraint on the amount of suppression explored at these two different scales. Focusing on the specific baryonic components, the gas mass provides a better tracer at |$k_\text{scale} = 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| (⁠|$R^2_\text{gas,1} \sim 0.89$|⁠; |$R^2_\text{stars,1} \sim 0.62$|⁠), whereas the stellar mass performs better at scales of |$k \sim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| (⁠|$R^2_\text{gas,10} \sim 0.73$|⁠; |$R^2_\text{stars,10} \sim 0.77$|⁠).

Comparing now the different AGN models, note that the quasar mode duty cycle provides a mechanism that efficiently suppresses the stellar mass of massive haloes with a lower gas ejection than the enhanced radio mode. Hence, the Quasar simulation has a higher suppression at small scales |$k_\text{scale} \sim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| (closer to RadioStrong), whereas it is closer to the RadioWeak case at large scales (⁠|$k_\text{scale} \sim 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠). As a result, its large suppression of stellar mass provides a better correlation at |$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, where it displays more clustering suppression. RadioStrong and Quasar feedback models displace the correlation in different directions, once again illustrating how different AGN implementations may allow to modify galaxy properties, mass content and the impact of baryons on the |$P_\text{mm}$| separately.

4 CAVEATS AND CONSIDERATIONS

While the studied fable simulations provide a high-resolution model for AGN feedback in galaxies and clusters, and their backreaction on the MPS, the comparatively small box size (i.e. 40 cMpc |$h^{-1}$|⁠) leads to a number of important shortcomings that should be considered.

As a result of simulation boxes having a limited sampling of their largest wavelength modes, domains with modest cosmological sampling feature a low representation of relevant domain modes, which leads to an underestimation of the MPS power at large scales (Daalen & Schaye 2015; Schaye et al. 2023). This shortcoming may be important down to scales of |$k_\text{peak} \lesssim 0.02\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, with MPS power only decreasing at |$k < k_\text{peak}$|⁠. For our simulated domain, this limited sampling of large modes with wavenumber below a few |$\,\,\text{h}\,\,\text{cMpc}^{-1}$| leads to an underrepresentation of massive haloes, as well as the proportion of the total mass that is contained within this systems. For the specific box size of the fable simulations studied here, this corresponds to galaxy clusters with halo masses of |$\gtrsim 10^{13} \mathrm{{\rm M}_{\odot }}$| (Daalen & Schaye 2015). Importantly, at scales |$k \sim 1\!-\!5\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, clusters with such masses contribute most of the MPS power (see e.g. Loon & van Daalen 2023). Previous studies identify such massive clusters, specifically those with |$M_\text{halo}\sim 10^{14}\, \mathrm{{\rm M}_{\odot }}$|⁠, to dominate the suppression of power in the MPS (e.g. Daalen et al. 2020; Salcido et al. 2023). These caveats emerging from the relatively small box size of our models imply that our results in the |$k \lesssim 1 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$| regime should be interpreted with particular care.

Despite these caveats, we note that the primary goal of this investigation is highlighting how variations of AGN feedback modelling in the high resolution regime leads to important changes in the baryonic feedback of galaxies and clusters on the MPS. A broad comparison of our fiducial model with larger simulations further places this variability into context, showing that it follows consistent qualitative trends to models such as BAHAMAS or FLAMINGO; while featuring a higher power suppression at scales of |$k \sim 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| than other models such as IllustrisTNG300 or MilleniumTNG. Finally, Bigwood et al. (2025) find the fable Fiducial model has a comparable impact on the MPS across different cosmological seeds and simulated domain sizes spanning from the studied 40 up to 100 cMpc |$h^{-1}$|⁠, particularly at scales |$k \lesssim 1 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, further supporting the qualitative robustness of our conclusions down to scales of |$\sim 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠. Bigwood et al. (2025) also reveal that at scales of |$k \sim 10 \, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, cosmic variance can lead to variations of relative power suppression between |$0.86-0.91$|⁠, with the Fiducial fable model studied here being representative of the spread and average of their new simulations.

5 CONCLUSIONS

In this work, we study how variations in the radio and quasar mode around the fiducial fable AGN model (Henden et al. 2018) impact the distribution of matter at different cosmic times. The fable simulations are performed with the arepo code (Springel 2010), evolving a uniform cosmological box with |$40\, \text{h}^{-1}\, \text{Mpc}$| on a side and featuring a galaxy formation model following illustris (Genel et al. 2014; Vogelsberger et al. 2014b; Sijacki et al. 2015). In addition to a dark matter-only simulation, the studied suite of 5 models spans: no AGN feedback (NoAGN), weak AGN feedback (RadioWeak), stronger AGN radio mode (RadioStrong), a quasar mode duty cycle (Quasar), and a fiducial model combining the stronger radio mode with the quasar duty cycle (Fiducial).

For each of these models, we investigate the matter power spectrum (MPS) and how different haloes selected accordingly to varying |$M_\text{halo}$|⁠, |$M_{*}$| and |$M_\text{BH}$| thresholds contribute to it. Our main findings are summarized as follows:

• The fable AGN models feature the largest MPS power suppression at scales of |$k \sim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| and at |$z = 0$|⁠, with a reduction of |$\sim 10~{{\ \rm per\ cent}}$| with respect to the DMO scenario. At |$k \sim 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠, the Fiducial model has a clustering suppression of |$\sim 0.012 P_\text{mm, DMO}$|⁠. The impact of baryonic feedback on the MPS in the fable Fiducial simulation is in general comparable to Horizon-AGN and IllustrisTNG-100, but is more similar to the BAHAMAS and Flamingo simulations at scales |$k \lesssim 5\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠.

• Stronger radio mode feedback (RadioStrong) is more effective at suppressing power at large scales (particularly |$k \lesssim 5 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) and at late cosmic times (⁠|$z \lesssim 1$|⁠). The effects of the quasar duty cycle (Quasar) are complementary to this, being more effective at smaller scales (⁠|$k \gtrsim 10 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠), and with their most important impact at early cosmic times (⁠|$3 < z < 1$|⁠). Variations in these two modes allow for comparable impacts at |$z \lesssim 0.5$|⁠, but importantly lead to significantly different redshift evolution up to |$z \sim 3$|⁠, which future observations probing into the high-redshift regime will be able to constrain (Huang et al. 2019).

• Clustering suppression in fable takes place around the most massive galaxies (⁠|$M_\text{halo}> 2.5 \times 10^{12}\mathrm{{\rm M}_{\odot }}$|⁠). This is approximately unchanged whether galaxies are selected employing halo or stellar masses. Smaller galaxies display no significant MPS suppression, except for small haloes (⁠|$M_\text{halo}\sim 10^{11}\mathrm{{\rm M}_{\odot }}$|⁠) in the Fiducial simulation, which are likely satellites or neighbours of massive galaxies hosting large SMBH. Interestingly, selecting haloes above a given central SMBH mass threshold leads to the highest amount of relative MPS power suppression with respect to the DMO scenario, particularly at scales |$k \gtrsim 7 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠.

• The baryonic impact on the MPS at scales of |$k \gtrsim 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$| is primarily due to clustering suppression within haloes at |$z \lesssim 1$|⁠. At higher redshift, and larger scales, power suppression comes from a combination of modifications of the matter distribution both inside and outside of haloes.

• The total baryonic mass content in the most massive haloes of the fable simulation provides an accurate tracer of MPS power suppression both at large (⁠|$k_\text{scale} = 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) and small (⁠|$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) scales (Daalen et al. 2020). However, modifications in the quasar and radio AGN modes drive the correlation of the stellar and gaseous components in opposite directions.

Overall, our results illustrate how different AGN implementations, and especially variations of the quasar and radio mode feedback, have distinct effects on the distribution of matter and the MPS redshift evolution. Our work motivates further delving into resolving the AGN feedback physics occurring on small scales (⁠|$\lesssim 1$| kpc) and capturing it within larger computational domains to better understand how these results vary in more representative cosmological domains (Daalen & Schaye 2015; Salcido et al. 2023; see Section 4). Additionally, it also supports further exploration of more sophisticated and physically motivated feedback models, either through enhanced resolution in the innermost regions of galaxies (e.g. Curtis & Sijacki 2016; Beckmann et al. 2019; Bourne, Sijacki & Puchwein 2019; Bourne & Sijacki 2021; Martin-Alvarez et al. 2022; Hopkins et al. 2024b), more realistic modelling of SMBH accretion and AGN activity (e.g. Bourne & Sijacki 2017; Talbot et al. 2021; Huško et al. 2022; Koudmani et al. 2023; Rennehan et al. 2023), or even through the inclusion of non-thermal components (e.g. Pfrommer et al. 2017; Costa et al. 2018; Martin-Alvarez et al. 2021; Su et al. 2021; Beckmann et al. 2022; Ruszkowski & Pfrommer 2023; Wellons et al. 2023).

ACKNOWLEDGEMENTS

We kindly thank the referee for their insightful comments and suggestions that helped improve the quality of this work. S.M.A. is supported by a Kavli Institute for Particle Astrophysics and Cosmology (KIPAC) Fellowship, and by the NASA/DLR Stratospheric Observatory for Infrared Astronomy (SOFIA) under the 08_0012 Program. SOFIA is jointly operated by the Universities Space Research Association, Inc. (USRA), under NASA contract NNA17BF53C, and the Deutsches SOFIA Institut (DSI) under DLR contract 50OK0901 to the University of Stuttgart. S.M.A acknowledges visitor support from the Kavli Institute for Cosmology, Cambridge, where part of this work was completed. V.I. acknowledges support by the Kavli Foundation and the KICC Fellowship. S.K. is supported by a Flatiron Research Fellowship at the Flatiron Institute, a division of the Simons Foundation and a Junior Research Fellowship from St Catharine’s College, Cambridge. S.K, M.A.B., and D.S. acknowledge support by European Research Council Starting Grant 638707 ‘Black holes and their host galaxies: coevolution across cosmic time’. M.A.B and D.S. additionally acknowledge support from the Science and Technology Facilities Council (STFC; grant number ST/W000997/1). M.A.B. acknowledges support from a UKRI Stephen Hawking Fellowship (EP/X04257X/1). This work made use of the following DiRAC facilities (www.dirac.ac.uk): the Data Analytic system at the University of Cambridge [funded by a BIS National E-infrastructure capital grant (ST/K001590/1), STFC capital grants ST/H008861/1 and ST/H00887X/1, and STFC DiRAC Operations grant ST/K00333X/1] and the COSMA Data Centric system at Durham University (funded by a BIS National E-infrastructure capital grant ST/K00042X/1, STFC capital grant ST/K00087X/1, DiRAC Operations grant ST/K003267/1 and Durham University). DiRAC is part of the National E-Infrastructure.

DATA AVAILABILITY

Data employed in this manuscript will be shared upon reasonable request by contacting the corresponding author.

The raw Fable MPS data for each of the studied models at |$z = \left[0, 0.4, 1.0, 2.0\right]$| is publicly available at:

https://github.com/MartinAlvarezSergio/Fable_MPS

Footnotes

REFERENCES

APPENDIX A: TEMPORAL EVOLUTION OF AGN FEEDBACK MODES

In the main body of this work, we found different AGN feedback models to exhibit distinct impacts on the MPS. Here, we briefly discuss the evolution of AGN feedback modes (quasar versus radio) across the redshift evolution of fable.

The fraction of AGN in each mode are shown in Fig. A1, where solid lines represent the quasar mode and dashed lines denote the radio mode. The quasar mode is predominantly active at higher redshifts. This aligns with the expected behaviour where quasar mode, being associated with high accretion rates, is more prevalent during the early, more chaotic epochs of galaxy formation. Conversely, the radio mode, which is often linked to maintenance feedback in more evolved systems, becomes dominant after |$z \lesssim 4$|⁠. Despite a different feedback switch fraction in the Fiducial model, the measured impact on the MPS appears dominated by the combination of its higher radiative efficiency (akin to RadioStrong) and the quasar duty cycle (also included in the Quasar simulation). At high redshift, the largest SMBHs are well within the quasar mode Eddington fraction regimen. Due to these AGN in these haloes being the main drivers of power suppression, we expect the feedback switch fraction to have relatively little effect on the MPS.⁴ Except for this aspect of the Fiducial model, all simulations exhibit a broadly similar evolution pattern. Such evolutionary trend for the AGN feedback modes provide further context for why our modifications to the quasar mode (Quasar and Fiducial) are particularly important at high redshift (⁠|$z \gtrsim 1$|⁠), whereas the so-called ‘maintenance’ radio mode effects reveal themselves after |$z \lesssim 1$|⁠.

APPENDIX B: REDSHIFT EVOLUTION OF THE HALOES MPS

In Section 3.2, we studied how the total MPS evolved over redshift, whereas Sections 3.4 and 3.5 addressed the effect of halo selection on the MPS, and the evolution of power at a fixed large (⁠|$k_\text{scale} = 1\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) and small (⁠|$k_\text{scale} = 10\, \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠) scale, respectively. To complement these two aspects, here we provide further detail on the evolution of the haloes MPS. Understanding the evolution of the MPS for matter inside haloes will be particularly important as future observatories probe systems at |$3 > z \gtrsim 0.5$| (Ade et al. 2019).

Fig. B1 is complementary to Fig. 3, but now shows the MPS from matter within haloes. It displays the redshift evolution of |$P_\text{mm, halos}$| from |$z = 3$| to |$z = 0$| in the fable AGN simulations: Fiducial (top left), RadioWeak (top right), Quasar (bottom left), and RadioStrong (bottom right). Overall, the trends observed in the main text regarding the temporal and scale-dependent impacts of different AGN feedback models are reproduced in the autopower of matter inside haloes here. At larger scales (⁠|$k \lesssim 5 \,\,\text{h}\,\,\text{cMpc}^{-1}$|⁠), the amount of power residing within haloes is significantly smaller than the total |$P_\text{mm}$|⁠, but features a comparable evolution of all models.