EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 687 (July 2024) A&A, 687 (2024) A68 Full HTML
Open Access
Issue
A&A
Volume 687, July 2024
Article Number A68
Number of page(s) 17
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202349045
Published online 27 June 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

It is a well-established result that a fraction of the galaxy population at all cosmic epochs are not actively forming stars, and typically have red colours and old stellar populations. This fraction is also known to vary with galaxy environment and cosmic epoch. In the local Universe, where rather accurate estimates of the star formation activity can be obtained using spectroscopic information for large statistical samples of galaxies, they are shown to follow a “bimodal” distribution, with galaxies being either blue and actively forming stars or red and “quenched” (e.g. Blanton et al. 2003; Kauffmann et al. 2003; Baldry et al. 2004, just to cite the milestone papers on the subject based on the Sloan Digital Sky Survey). The “transition” region, sometimes dubbed “green valley”, is relatively empty and its position does not vary significantly in different environments, suggesting that if galaxies move from one region to another, this must occur on a relatively short timescale (Wetzel et al. 2012) and must be driven by a process that is effective in different environments.

At earlier cosmic epochs, the separation between quiescent and active galaxies cannot be achieved using direct measurements of the star formation activity (at least not for large galaxy samples), and it is instead based on colour-colour selections. This, together with errors associated with photometric redshifts, make the separation between passive and active galaxies subject to larger uncertainties than when using direct estimates of star formation rates (SFRs). Deep photometric surveys have enabled studies of the luminosity, mass, and colour distribution of galaxies at higher redshifts, confirming that colour bimodality persists up to z ∼ 4 and beyond (e.g. Brammer et al. 2009; Whitaker et al. 2011; Muzzin et al. 2013; Weaver et al. 2023). The advent of the James Webb Space Telescope (JWST) has allowed us to push these studies to even earlier cosmic epochs, revealing a large population of massive and apparently quiescent galaxies at z > 3, in excess of previous observational inferences (e.g. Carnall et al. 2023; Valentino et al. 2023).

This rich amount of information has triggered numerous theoretical studies and highlighted, over the years, important disagreements between observational measurements and theoretical predictions. These have, in turn, led to significant improvements in published models. Since early renditions of semi-analytic models, it has been realized that an additional source of energy at the centre of massive systems is required to avoid the formation of overly massive and star-forming galaxies. In early models (e.g. Kauffmann et al. 1999; Benson et al. 2003; De Lucia et al. 2004), this was implemented in a very simplistic way by suppressing gas cooling above some critical halo mass or circular velocity. Later theoretical developments, including a physical treatment of the impact of feedback from active galactic nuclei (AGN), have confirmed that this process plays a crucial role in suppressing the cooling flows at the centre of massive clusters, thereby reconciling the massive end of the galaxy stellar mass function with observations and leading to a sequence of red and massive galaxies that is qualitatively consistent with data (e.g. Croton et al. 2006; Bower et al. 2006; Sijacki et al. 2007; Lagos et al. 2008; Dubois et al. 2013; Puchwein & Springel 2013, and many others).

Early detailed comparisons between theoretical predictions and colour distributions of galaxies in the local Universe, split by centrals and satellites, showed a significant excess of passive satellites with respect to observational estimates (Weinmann et al. 2006). It was soon realized that these discrepancies were shared by all state-of-the-art theoretical models of galaxy formation published in those years (Fontanot et al. 2009). A revised treatment of satellite galaxies led to some improvements, but reproducing the colour or SFR distributions observed for satellites and central galaxies in the local Universe has long remained, and to some extent still represents, an important challenge for both semi-analytic and hydro-dynamical simulations of galaxy formation (e.g. Hirschmann et al. 2014; Wang et al. 2018; De Lucia et al. 2019; Xie et al. 2020; Donnari et al. 2021; Bravo et al. 2023, and references therein). The challenges become even greater when trying to reproduce the estimated fractions of quiescent galaxies at earlier cosmic epochs (De Lucia et al. 2019; Donnari et al. 2021; Lustig et al. 2023; Weaver et al. 2023), or when focussing on the cluster environment where the galaxy population is dominated by satellites (Bahé et al. 2017; Kukstas et al. 2023).

In this paper, we present an updated version of our GAlaxy Evolution and Assembly (GAEA) model, that combines our recent improved treatments for AGN feedback and ram-pressure stripping of hot and cold gas from satellite galaxies, and an explicit treatment for partitioning the cold gas in its atomic and molecular components. We make a critical assessment of the predicted trends for passive galaxies, as a function of physical properties (stellar mass primarily), redshift, and environment. In a companion paper (Xie et al. 2024), we discuss the emergence and physical origin of the first quiescent galaxies.

The layout of the paper is as follows: in Sect. 2, we present the model and the simulations used in this paper. In Sect. 3, we compare predictions from our new model with observational estimates of the distribution of the specific star formation rate (sSFR) and passive fractions in the local Universe, while in Sect. 4 we consider observational estimates of passive fractions at earlier cosmic epochs. In Sect. 5, we present a comparison between our model predictions and estimated passive fractions in galaxy clusters at z ∼ 1. Finally, in Sect. 6, we discuss our results and give our conclusions.

2. The simulation and the galaxy formation model

GAEA1 traces the formation and evolution of galaxies in a cosmological framework, using state-of-the-art prescriptions for the evolution of different baryonic components, coupled to substructure based merger trees extracted from high-resolution N-body simulations. The galaxy formation model builds on that presented in the original paper by De Lucia & Blaizot (2007), but it has been updated significantly over the years. In particular, our model includes:

– a detailed treatment of the non-instantaneous recycling of gas, energy, and metals, with an accurate accounting of their timings depending on stellar lifetimes, which enables the tracing of individual metal abundances (De Lucia et al. 2014);

– a parametrization of stellar feedback that is partially based on results from hydro-dynamical simulations, and that we have shown allows us to reproduce the observed galaxy stellar mass function up to z ∼ 3 (Hirschmann et al. 2016);

– an explicit partition of the cold gas in its atomic and molecular components, tuned to reproduce the measured HI and H2 galaxy mass functions in the local Universe (Xie et al. 2017);

– a careful tracing of the angular momentum exchanges between different components, and a treatment for the non-instantaneous stripping of the cold and hot gaseous components associated with galaxies infalling onto larger systems (Xie et al. 2020);

– an improved model for cold gas accretion on super massive black holes (BHs) and AGN driven outflows that, in the framework of our model, reproduces the measured evolution of the AGN luminosity function and the so-called AGN downsizing trend (Fontanot et al. 2020).

The last two implementations enumerated above represent so far two independent branches of our model, that we have now merged into the model presented in this work. Our model also features the possibility to adopt a variable stellar Initial Mass Function (IMF; Fontanot et al. 2017a, 2024), that we do not consider in this work where we only present results based on a Chabrier IMF.

In previous work, we have shown that our reference model (Hirschmann et al. 2016, and later versions) is able to reproduce a large number of important observational constraints. These include the galaxy stellar mass function up to z ∼ 7 and the cosmic star formation rate density up to z ∼ 10 (Fontanot et al. 2017b), the relation between galaxy stellar mass and gas metallicity and its secondary dependence as a function of star formation rate and gas mass (De Lucia et al. 2020), and the observed evolution of the galaxy mass – gas/star metallicity relations (Hirschmann et al. 2016; Fontanot et al. 2021). Our model is of course not without problems. In previous work, we have emphasized in particular the difficulty to reproduce the bimodality in colour or specific star formation rate observed in the local Universe (Hirschmann et al. 2016; Fontanot et al. 2020), and the tendency of the model to underpredict the observed fractions of passive galaxies at higher redshift (De Lucia et al. 2019). As mentioned in Sect. 1, these difficulties are shared to different extent by virtually all recently published theoretical models of galaxy formation. We show, in the following, that these problems are largely overcome by the latest rendition of GAEA presented in this work.

We refer to the original papers cited above (and references therein) for a detailed description of the specific implementations adopted. As for AGN feedback, we use what in Fontanot et al. (2020) is referred to as the F06-GAEA set of parametrizations, that include: (i) a modelling for the inflow of cold gas towards the central massive black holes driven by star formation in the central regions of galaxies; (ii) an accretion rate that is determined by the viscous accretion timescale; (iii) an empirical scaling relation between the mass loading of quasar winds and the bolometric luminosity. In Table 1, we list the values of the relevant parameters adopted in Xie et al. (2020, GAEA2020X hereafter) and Fontanot et al. (2020, GAEA2020F), and the corresponding values adopted in the new model whose results are discussed in this work (GAEA2023 in the following). All parameters that are not included in the table have been left unchanged with respect to the model versions published in our previous work. The primary observables that have been considered to recalibrate our model are: the observed galaxy stellar mass function and its evolution up to z ∼ 3, the HI and H2 mass functions in the local Universe, and the evolution of the AGN luminosity function up to z ∼ 4. When calibrating our model parameters, we also make sure that the mass–metallicity relations predicted at z = 0 have approximately the correct normalization. We come back to the impact of these parameter modifications on our model predictions in the final discussion section. In Appendix A, we show the main figures used to calibrate our model parameters.

Table 1.

Values of the relevant parameters used in Xie et al. (2020) and Fontanot et al. (2020), and the corresponding values adopted in the model considered in this study.

The results presented in this paper are based on the Millennium Simulation (Springel et al. 2005). This simulation follows 21603 dark matter particles in a box of 500 Mpc h−1 on a side, and assumes cosmological parameters consistent with WMAP1 (ΩΛ = 0.75, Ωm = 0.25, Ωb = 0.045, n = 1, σ8 = 0.9, and H0 = 73 km s−1 Mpc−1). While more recent measurements provide slightly different cosmological parameters (in particular, a larger value for Ωm and a lower value for σ8), we do not expect these differences to affect significantly our model predictions, once the model parameters have been returned to reproduce a specific set of observational results (see also Wang et al. 2008; Guo et al. 2013). In a forthcoming paper, we will re-address this issue presenting results from the same model considered here but applied to a simulation with cosmological parameters consistent with a Planck cosmology. To verify convergence of our physical model, we also use in this work the MillenniumII Simulation (Boylan-Kolchin et al. 2009). This assumes the same cosmological model and uses the same number of particles of the Millennium Simulation, but has a smaller volume (one-fifth the size of the Millennium, that is, 100 Mpc h−1), and 125 times better mass resolution (the particle mass is 6.9 × 106 M⊙h−1).

3. The specific star formation rate distribution in the local Universe

We start, in this section, by considering model predictions at z = 0. In particular, we focus here on the distribution of sSFR and on the fraction of passive galaxies, as a function of the galaxy stellar mass. The observational data shown in this section are based on data from the Sloan Digital Sky Survey (SDSS) DR8, cross-correlated with the JHU-MPA DR7 catalogue2 and with the DR7 version of the group catalogue by Yang et al. (2007). Stellar masses included in the JHU-MPA DR7 catalogue are estimated based on fits to the photometry following the philosophy of Kauffmann et al. (2003) and Salim et al. (2007). SFRs are based on Brinchmann et al. (2004). Aperture corrections, that are computed by fitting the photometry of the outer regions of galaxies as detailed in Salim et al. (2007), have been used to convert to total SFR estimates.

Figure 1 shows the sSFR distributions for bins of galaxy stellar mass (different panels). Shaded histograms show the observational estimates, while open histograms correspond to predictions from different versions of GAEA, as indicated in the legend. For the observational sample, we have considered only galaxies in the redshift range 0.025 < z < 0.05, where we expect the sample to be approximately volume complete down to a galaxy stellar mass ∼109 M⊙. Following Fontanot et al. (2020), we have assigned an “upper observable value” for the SFR to each model galaxy that has a theoretical SFR < 10−4 M⊙ yr−1. For these galaxies, we have adopted the following relation, that reproduces the locus of the upper limits of passive galaxies in the SDSS sample:

log [ SFR ] = 0.5 × log [ M star ] − 6.59 . $$ \begin{aligned} \log [\mathrm{SFR}] = 0.5 \times \log [M_{\rm star}] - 6.59. \end{aligned} $$

thumbnail Fig. 1.

Specific star formation rate distribution as predicted by the models indicated in the legend, and compared to observational estimates based on SDSS DR8 (grey shaded histograms – see text for details). Each panel corresponds to a different bin in galaxy stellar mass, as indicated in the top-right legend. The histograms corresponding to the GAEA2020X and GAEA2020F models have been slightly shifted (by 0.02 dex) to the left and right with respect to the GAEA2023 model for clarity.

The above relation has also been perturbed by assuming a lognormal scatter of 0.25 dex.

The figure shows that, as discussed in previous work, both the GAEA2020X and GAEA2020F versions of our model have problems in reproducing the observed sSFR distributions. In particular, GAEA2020X follows relatively well the observational distribution for the lowest galaxy mass bin considered, although the predicted distribution is somewhat skewed towards low values of sSFR with respect to data. However, this is the range where incompleteness in the observed sample and/or our treatment of upper observable values of SFR become relevant. At intermediate masses, the GAEA2020X model does not reproduce well the observed bimodality in sSFR, and for the most massive galaxies the predicted distributions are shifted towards sSFR values larger than observed. We had highlighted similar problems in Hirschmann et al. (2016), where we also noted that we were unable to solve the excess of star formation in massive galaxies with simple modifications of the (radio mode) AGN feedback. The reason for this, as discussed in Fontanot et al. (2020), is that this star formation excess is not related to late gas cooling nor to gas brought in by satellites accretion. We found that the main reason for this excess residual activity was the large gaseous content of the main progenitors of these galaxies at z ∼ 2. Figure 1 shows that the inclusion of AGN driven winds reduces significantly, albeit not removing completely, the excess of activity in massive galaxies. This feedback process also introduces a clear bimodality in the distributions of galaxies at intermediate stellar masses. However, for the lowest mass bins considered, the GAEA2020F model over-predicts significantly the number of galaxies with low sSFR values.

The GAEA2023 model presented in this paper captures well the observed trends, for all galaxy mass bins considered. At the lowest stellar masses, model predictions are very close to those from the GAEA2020X model, with no excess of passive galaxies. At intermediate masses, there is a clear bimodality, with a minimum at log[sSFR] ∼ −11 yr−1, consistent with observations. No significant excess of star formation activity is found for the most massive galaxies. Considering that we are not making any attempt to mimick in detail the observational selection, and taking into account uncertainties in observational estimates of SFRs and aperture corrections, the level of agreement shown in Fig. 1 between observational data and our GAEA2023 model is remarkable.

The improved agreement with data of the GAEA2023 model with respect to previous model renditions can be explained as follows: for low-mass galaxies, the better regulation of star formation is driven by the improved treatment of environmental processes introduced in the GAEA2020X model. In fact, as noted above, predictions from the latter model are very similar to those from GAEA2023 in this mass range. In contrast, satellite galaxies tend to be over-quenched in the GAEA2020F model due to the instantaneous stripping of their hot gas reservoir at the time of accretion, combined with a rather efficient supernovae feedback. The slight modification of parameters between the GAEA2020X and GAEA2023 models also leads to slightly larger gas fractions (and therefore larger values of star formation) for low-mass galaxies, resulting in a slight improvement with respect to the GAEA2020X model.

At the most massive end, the inclusion of an efficient feedback channels associated with quasar winds triggered by mergers and disk instability, leads to a sSFR distribution that is peaked at values comparable to those observed but leaves a non negligible excess of galaxies with sSFR ∼ 10−11 yr−1 in the GAEA2020F model. We find that this excess is removed in the GAEA2023 model thanks to an the overall better improvement of the density threshold for star formation. In fact, the prescriptions adopted in GAEA2020F lead to a very bursty behaviour of the star formation, with relatively large amounts of gas accumulating for some time-steps and then leading to intense short bursts of star formation. The partition of the cold gas in its atomic and molecular gas components, and the assumption that star formation can only take place in regions corresponding to a certain surface density of molecular hydrogen (for details, see Xie et al. 2017), limits significantly this behaviour.

At intermediate galaxy stellar masses, the more pronounced bimodality with respect to the GAEA2020X model, as well as the good agreement between the position of the predicted and observed “passive” peak of the sSFR distribution can also be largely ascribed to efficient ejections of gas due to quasar winds.

Figure 2 shows the fraction of passive galaxies for all galaxies (left panel), central galaxies only (middle panel), and satellites (right panel). Solid lines show predictions from the different models, as indicated in the legend, and all correspond to a sSFR threshold of 10−11 yr−1 to identify passive galaxies (as noted above, this describes well the separation between the two sSFR peaks for the observational data). The dashed line corresponds to the GAEA2020X model and to a threshold of 10−10.6 yr−1. This is the threshold that was adopted in Xie et al. (2020) and better describes the separation between the two sSFR peaks in that model. Symbols correspond to observational estimates: squares with error bars are based on GAMA (Davies et al. 2019) and assume a sSFR threshold equivalent to 10−10.5 yr−1 to select passive galaxies; filled and open circles show estimates based on SDSS data assuming a threshold of 10−11 yr−1 and 10−10.6 yr−1, respectively. The figure shows that the latest version of our model reproduces very nicely the observational estimates in the local Universe. In particular: for central galaxies, the predicted passive fractions decrease monotonically with decreasing galaxy stellar mass, down to masses ∼109 M⊙. The same trend is found for satellite galaxies, also in quite good agreement with observational estimates. Previous versions of our model over-predict the passive fraction of satellites below galaxy stellar masses ∼109.5 − 1010.5 M⊙, and under-predict the fraction of passive massive satellites. This is more significant for the GAEA2020X model, especially when considering the same sSFR threshold in all models. To give an idea of the cosmic variance for our model, we divide the simulated volume of the Millennium in 125 independent box of 100 Mpc h−1 on a side, and compute the expected quiescent fractions in each independent sub-volume. The grey regions shown in each panel represent the regions between the minimum and maximum quiescent fractions obtained. The variance is relatively small at these scales, but not negligible, and it becomes significant at the largest galaxy stellar masses (the rarest systems).

thumbnail Fig. 2.

Passive fractions for all galaxies (left panel), central galaxies (middle panel), and satellites (right panels). Symbols correspond to observational measurements. In particular: squares correspond to estimates from Davies et al. (2019) based on the GAMA survey, and adopt a log(sSFR) threshold of −10.5 yr−1 to identify passive fractions. Filled and empty circles correspond to measurements based on SDSS DR8 and are obtained assuming a threshold of −11 yr−1 and −10.6 yr−1, respectively. Solid lines correspond to predictions from different models, as indicated in the legend, assuming a log(sSFR) threshold of −11.0 yr−1 to select passive galaxies. The grey regions show the minimum-maximum fractions obtained when considering predictions from the GAEA2023 model from 125 independent boxes, each of 100 Mpc h−1 on a side. The dashed line shows predictions from the GAEA2020X model corresponding to a log(sSFR) threshold of −10.6 yr−1, as adopted in Xie et al. (2020).

4. Passive galaxies beyond z = 0

In this section, we present model predictions for passive fractions beyond the local Universe, and compare them with different observational estimates. We emphasize that these are true predictions, because these data have not been considered when tuning the model parameters. We use a sSFR selection for model galaxies. Specifically, we select as passive galaxies all those with sSFR smaller than 0.3 × t Hubble − 1 $ 0.3\times t^{-1}_{Hubble} $ (Franx et al. 2008), which corresponds to ∼10−11 yr−1 at z = 0. In Sect. 6, we discuss the impact of a different selection.

Figure 3 shows the predicted passive fractions in three different redshift bins, for central (solid lines) and satellite galaxies (dashed lines). Symbols with error bars show the corresponding observational estimates by Fossati et al. (2017). These have been computed using spectroscopic and grism redshifts in the five CANDELS/3D-HST fields, and multiwavelength photometry to classify galaxies as passive. The central/satellite information has been obtained by assigning, to each galaxy, probability functions computed using a mock galaxy catalogue that reproduces the 3D-HST sample selection and redshift accuracy. The fraction of passive galaxies tends to increase, both in the models and in the data, with decreasing redshift. For the GAEA2020X model, the fractions of passive centrals and satellites are very similar at all redshifts and for the entire galaxy mass range considered, but for the most massive galaxies where the fraction of passive satellites turns over (as it does at z = 0 – see right panel of Fig. 2). This trend is in contrast with observational estimates that always find a fraction of quenched satellites larger than that of quenched centrals, with the difference decreasing with increasing redshift. In addition, the overall predicted fractions of passive satellites from the GAEA2020X model are lower than estimated, at all redshifts. This is true also for the fraction of passive centrals at the highest redshifts considered. The fractions of passive centrals predicted by the GAEA2020F model are very similar to those predicted by the GAEA2020X model, while the fractions of passive satellites are larger. This brings model predictions in better agreement with observational estimates for the most massive satellites, but also leads to a significant over-prediction of passive satellites with masses below 1010 M⊙, over the entire redshift range considered. Our latest model version exhibits the best agreement with the observational data considered in Fig. 3, with a similar decrease of the difference between satellites and centrals with increasing redshift. The predicted fractions of passive centrals tend to be larger than observed for the most massive galaxies, but the uncertainties here are larger due to small number statistics and larger uncertainties in the classification of centrals and satellites.

thumbnail Fig. 3.

Fraction of passive centrals (solid lines) and satellites (dashed lines) as predicted by the models considered in this study. Symbols with error bars (squares are for centrals and triangles for satellites) show observational estimates by Fossati et al. (2017), based on five CANDELS/3D-HST fields.

Pushing the comparison to earlier cosmic epochs, we consider the observational estimates by Muzzin et al. (2013) and Weaver et al. (2023). The former3 are based on a sample of ∼100 000 Ks selected galaxies in the COSMOS/UltraVISTA field, with photometric redshifts estimated using photometry covering the wavelength range 0.15 − 8.0 μm. The mean estimated uncertainty of photometric redshifts is δz/(1 + z) = 0.013, and the estimated fraction of catastrophic outliers is only ∼1.6%. Stellar masses have been estimated using SED fitting and a Kroupa IMF, consistent with the one adopted in our model. Passive galaxies have been identified using a (UVJ) colour-colour selection. Estimates4 by Weaver et al. (2023) are based on the most recent release of the COSMOS catalogue (the COSMOS2020, Weaver et al. 2022). The sample, including approximately 1 million galaxies, is complete down to 109 M⊙ out to z ∼ 3, and is supported by extensive photometry ranging from the UV to 8 μm and including deep near-infrared UltraVISTA imaging and Subaru Suprime-Cam intermediate bands. The quoted photo-z uncertainties are δz/(1 + z) < 0.01 at i ∼ 20, and < 0.04 at i ∼ 26 AB. Passive galaxies have been selected adopting a (NUVrJ) colour-colour selection, and galaxy stellar masses have been estimated assuming a Chabrier stellar IMF.

Figure 4 shows a comparison between our model predictions and these observational estimates. The red open symbols with error bars and shaded regions show the maximum likelihood measurements and the 1/Vmax measurements by Muzzin et al. (2013). The blue symbols with error bars show the observational estimates by Weaver et al. (2023). Our model predictions are shown by lines of different colours, as indicated in the legend. For this comparison, we have selected all model galaxies from snapshots falling in the redshift bins used by Muzzin et al. (2013). These match well those considered in Weaver et al. (2023) but for the second and last redshift bins shown (see figure legend – the redshift bins adopted by Weaver et al. 2023 are written in blue). However, we have verified that the comparison is not affected significantly by this slight mismatch.

thumbnail Fig. 4.

Passive fraction at different cosmic epochs. Lines correspond to predictions from the models considered in this study. The shaded area and red symbols with error bars show observational estimates from Muzzin et al. (2013). Specifically, the shaded areas show the 1/Vmax measurements, while the open symbols with error bars show the maximum-likelihood measurements (see bottom panels of Fig. 6 in Muzzin et al. 2013). The blue symbols with error bars correspond to the observational estimates by Weaver et al. (2023). A sSFR cut has been employed to select model passive galaxies (see text for details), while the observational samples are based on different colour-colour selections, a SED fitting estimate of the galaxy stellar mass, and photometric redshifts estimates.

The GAEA2020X model exhibits a progressive worsening of the underprediction of passive fractions with increasing redshift. Only about 20% or less of the most massive galaxies are predicted to be passive at 1.5 < z < 2.0, against a fraction that is estimated to be as large as ∼60% by Muzzin et al. (2013). Less that 1% of the galaxies in the entire volume considered from the GAEA2020X run are predicted to be passive at z > 2.5. The inclusion of an explicit feedback mode associated with AGN (in the form of quasar driven winds), mainly relevant for radiatively efficient BH accretion, increases significantly the fraction of passive galaxies in the GAEA2020F model. For galaxies more massive than ∼1010 M⊙, the passive fractions predicted by this model are lower than observational estimates by Muzzin et al. (2013) but in quite nice agreement with estimates by Weaver et al. (2023). For less massive galaxies, this model over-predicts significantly the observational estimates. As shown in Fig. 3, the overprediction of passive fractions for lower mass galaxies is due primarily to an excess of quenched satellites in this model. Predictions from our new model, GAEA2023, are in very good agreement with the observational estimates by Muzzin et al. (2013), over the entire redshift range considered, but over-predict the passive fractions measured by Weaver et al. (2023) at galaxy stellar masses larger than ∼1011 M⊙ at 1 < z < 2. It is interesting that the observational measurements considered differ significantly over this redshift and mass range. It is difficult to identify one single possible reason for this difference: the two studies are based on a different colour-colour selection (UVJ in the case of Muzzin et al. 2013 and NUVrJ in the case of Weaver et al. 2023), but there are other important differences both in the data used and in the data analysis. For example, the deeper NIR selection employed by Weaver et al. (2023) might be more sensitive to dusty star-forming galaxies. In fact, many massive systems that are identified in COSMOS2020 are red systems best-fit by dusty star-forming spectral templates (Weaver, priv. comm.). While the differences between these published estimates can be considered as a rough indication of the systematic uncertainties for quenched fractions at large stellar masses, we note that measurements by Weaver et al. (2023) show a rather strong evolution between z ∼ 1 and z ∼ 1.5. This appears to be inconsistent with stellar population studies of massive quiescent galaxies in the local Universe finding very old stellar ages (and therefore early quenching times). At low galaxy stellar masses, the observational estimates considered are consistent, and our GAEA2023 model is in very nice agreement with the observed trends down to the lowest masses shown in the figure, with some hint for a turnover of passive fractions at the lowest masses accessible through the Millennium Simulation. We come back to this point in Sect. 6.

In addition to the quenched fractions shown above, it is instructive to consider the corresponding mass distributions from which they have been computed, that is, the galaxy stellar mass functions for all and for quiescent galaxies. A comparison between our model predictions and the galaxy mass functions published in Muzzin et al. (2013) and in Weaver et al. (2023) is shown in Figs. 5 and 6. Observational estimates are shown by symbols with error bars, while lines correspond to our model predictions as indicated in the legend. For the GAEA2023 model only, we show both the intrinsic model predictions, as well as the corresponding predictions obtained assuming that stellar masses have uncertainties that are Gaussian distributed with a width of 0.25 dex (dashed lines). We note that, while this is a typical statistical uncertainty quoted, it is well known that stellar mass estimates based on SED fitting depend on a number of assumptions including metallicity, stellar population models, stellar initial mass function, and star formation history (see e.g. Conroy et al. 2009), and that these can result into larger uncertainties and even systematic trends as a function of galaxy stellar mass (Mitchell et al. 2013). As emphasized above, we have used observational estimates of the galaxy stellar mass function, up to z ∼ 3, to constrain our model parameters for the GAEA2023 and GAEA2020X models. Therefore, the relatively good agreement that is shown in Fig. 5, when considering the galaxy mass function for all galaxies predicted by these models, does not come as a surprise. The slight excess of galaxies below the knee of the mass function, that is visible in particular in the top panels, is already visible in our previous papers (see e.g. Fig. A.1 in Xie et al. 2020), but becomes somewhat less evident when considering a compilation of independent observational estimates (see Appendix A). The deficit of massive galaxies that is evident for these two models, particularly in the bottom panels, is more susceptible to the statistical errors associated with the estimates of stellar mass, as shown by the dashed line in each panel. It is interesting to note that the GAEA2020F model, that has been tuned primarily to reproduce the AGN luminosity function rather than the galaxy stellar mass function, predicts both lower number densities of galaxies below the knee and larger number densities for the most massive galaxies. At the highest redshift considered, our GAEA2023 model appears to be in quite good agreement with observational estimates by Muzzin et al. (2013) while underpredicting the number densities of the most massive galaxies when compared to observational measurements by Weaver et al. (2023). In their work, Weaver et al. (2023) comment about the larger number densities of massive galaxies at these redshifts, when compared to earlier studies, and speculate that these might be due to the presence of proto-clusters in the field.

thumbnail Fig. 5.

Galaxy stellar mass function at different cosmic epochs. Symbols with error bars correspond to the observational estimates published in Muzzin et al. (2013, in this case we are showing only the 1/Vmax estimates published in this study) and Weaver et al. (2023). Solid lines correspond to predictions from the models considered in this study, as indicated in the legend. The dashed line in each panel shows, for the GAEA2023 model only, the corresponding predictions obtained assuming that the stellar masses have uncertainties that are Gaussian distributed with a width of 0.25 dex.

As for the galaxy mass function of passive galaxies (Fig. 6), both the GAEA2020X and the GAEA2023 models seem to work reasonably well at the lowest redshift bin considered, with the GAEA2020X model slightly underpredicting the number densities of quiescent galaxies around the knee and showing a turn-over below stellar masses ∼109.5 M⊙. At the same redshift, the GAEA2020F model exhibits a clear excess of passive galaxies below ∼1010 M⊙. These trends become stronger at higher redshifts, with a dramatic steepening of the galaxy mass function of quiescent galaxies below ∼1010 M⊙ for the GAEA2020F model. As we have discussed earlier, this is mainly driven by an excess of passive satellite galaxies. A turn over of the number densities of quiescent galaxies is predicted also by the GAEA2023 and GAEA2020X models, but it is moved towards lower masses and does not appear to be excluded by available data. Deeper observations at intermediate to high redshift can provide important constraints on these models. In Sect. 6, we discuss if and how the turn over at low galaxy stellar masses is affected by the resolution of the adopted simulation. As reflected in the low passive fractions presented above, the GAEA2020X model significantly underpredicts the number densities of quiescent galaxies at high redshift, with the under-prediction becoming more severe with increasing redshift. The latest version of our model predicts number densities of passive galaxies that are in rather good agreement with observational estimates up to z ∼ 3, over the entire range of masses sampled by the data considered in this figure, and no significant under-prediction of the most massive quiescent galaxies when one accounts for the uncertainties in the stellar mass estimates. In fact, as noted above, the typical uncertainties can also be much larger than the fiducial 0.25 dex assumed for the dashed lines shown. At the highest redshift shown, the number densities of massive quiescent galaxies predicted by our model, appear to under-predict the observational measurements by Weaver et al. (2023). As mentioned above, these authors argue that the larger number densities estimated at this redshift (with respect to other literature estimates) might be due to proto-clusters in the field (in one case, they also have spectroscopic confirmation). We get back to a more detailed comparison with observational measurements of massive quiescent fractions at z > 3 later on.

thumbnail Fig. 6.

As in Fig. 5, but for quiescent galaxies only.

5. Passive galaxies in clusters at z ∼ 1

As mentioned in Sect. 1, state-of-the-art models struggle to reproduce in particular the estimated fraction of quiescent galaxies in overdense regions, where the galaxy population is expected to be dominated by satellite galaxies (Bahé et al. 2017; Kukstas et al. 2023, see also Fig. 2 in Xie et al. 2020). In this section, we discuss how our model predictions compare to observational estimates based on a subset of clusters from the GOGREEN survey (Balogh et al. 2017, 2021). This is based on homogeneous deep imaging and spectroscopy of 21 galaxy groups and clusters at 1 < z < 1.5. Below, we restrict the observational sample to the ten massive clusters at 1.5 < z < 1.0 that were analysed in van der Burg et al. (2020). Cluster masses, determined from Jeans modelling as described in detail in Biviano et al. (2021), range between 1014.1 and 1014.8 M⊙. In the Millennium Simulation, there are ∼300 to ∼700 haloes with M200c5 in the same mass and redshift range of the GOGREEN cluster considered. Kukstas et al. (2023) recently compared the same observational measurements with three state-of-the-art suites of hydrodynamical simulations (BAHAMAS+MACSIS, EAGLE+Hydrangea, IllustrisTNG) showing that simulations generally reproduce well the stellar mass function of quenched galaxies in the field, but all struggle to reproduce the corresponding measurements in the cluster environment. In particular, all simulations considered in the work by Kukstas et al. (2023) over-predict the fraction of quenched low-mass satellite galaxies and two out of the three suites considered under-predict the number densities of quenched massive galaxies, probably due to inefficient AGN feedback.

To compare our model predictions with GOGREEN observational estimates, we have selected all simulated haloes with mass larger than 1014 M⊙ from the simulated volume, at all snapshots covering the redshift range of the observed sample. Then, we have selected 25 sub-samples that match both the halo mass and redshift distributions of the ten clusters used in van der Burg et al. (2020). As in this observational work, we also consider as cluster members only galaxies that are within a physical distance of 1 Mpc from the brightest cluster galaxies.

Figure 7 shows a comparison between the observational estimates of the galaxy stellar mass function presented in van der Burg et al. (2020) and our model predictions. Symbols with error bars show the galaxy mass function for star-forming (upper panel) and quiescent galaxies (lower panel), classified on the basis of their location in a UVJ colour-colour diagram. Solid lines in each panel show the corresponding predictions from our different models, as indicated in the legend. As in previous sections, we have used a time-dependent sSFR cut to select quiescent model galaxies. The shaded regions show the area enclosed between the 10th and 90th percentiles of the distribution of number densities obtained in the GAEA2023 model, for each galaxy stellar mass bin considered, when considering the 25 simulated samples discussed above. All model predictions shown in this section assume stellar mass uncertainties of 0.25 dex.

thumbnail Fig. 7.

Galaxy stellar mass function for the star-forming (top panel) and quiescent galaxies (bottom panel) in clusters at 1 < z < 1.4. Solid lines show predictions from the models considered in this study, as indicated in the legend. Symbols with error bars show corresponding observational estimates published in van der Burg et al. (2020). The shaded regions show the area enclosed by the 10th and 90th percentiles of the distributions of number densities obtained, for each galaxy stellar mass bin, when considering 25 samples of simulated haloes that match the mass and redshift distributions of the GOGREEN sample (see text for details). All model lines shown in this figure assume an uncertainty in stellar mass of 0.25 dex.

The figure clearly shows that previous versions of our model fail to reproduce simultaneously the observational measurements for quiescent and for star-forming galaxies: the GAEA2020X model undepredicts the number densities of quiescent galaxies at all stellar masses, particularly below the knee. The predicted mass function for star-forming galaxies is closer to observational data, with a slight excess of galaxies where the data exhibit a “dip” at low stellar masses, and an excess of more massive star-forming galaxies. The GAEA2020F model undepredicts the number densities of quiescent galaxies around the knee, and underpredicts significantly the number densities of star-forming galaxies below the knee. The GAEA2023 model exhibits the best agreement with data, especially when considering the effect of sample variance. However, it still undepredicts slightly the number densities of quiescent galaxies around the knee, and overpredicts slightly the number densities of star-forming galaxies at low-stellar masses.

Figure 8 shows the passive fractions corresponding to the mass functions considered above. The dashed lines (these are shown only for the GAEA2023 model for clarity – the dispersion predicted by the other model versions considered is similarly large) represent the 20th and 80th percentiles of the distributions obtained considering our 25 simulated samples. Consistently with what found in the field, the GAEA2020X model under-predicts the fractions of passive galaxies for the entire galaxy stellar mass range considered, with only few samples predicting a fraction of massive passive galaxies consistent with the data. The GAEA2020F model predicts a quiescent fraction that is rather constant as a function of galaxy stellar mass, consistent with data for the most massive galaxies, and significantly larger than observed for low-mass galaxies. The passive fractions predicted by the GAEA2023 model are slightly below the observational estimates when considering the average results, but consistent with the data within the expected sample variance. Considering that the selection of quiescent galaxies in the models and in the data is based on different criteria, the agreement of our GAEA2023 with data is overall remarkable.

thumbnail Fig. 8.

Passive fraction in the cluster environment at z ∼ 1. Lines show predictions from the models considered in this study, while symbols with error bars show the observational estimates from van der Burg et al. (2020). The dashed lines show the area enclosed by the 20th and 80th percentiles of the quiescent fractions obtained, for each galaxy stellar mass bin, when considering 25 samples of simulated haloes that match the mass and redshift distributions of the GOGREEN sample (see text for details). We show the dashed lines only for the GAEA2023 model for clarity. All model lines shown in this figure assume an uncertainty in stellar mass of 0.25 dex.

6. Discussion

In this paper, we have presented an updated version of our GAEA theoretical model for galaxy formation that combines two major implementations we have published in the last years: (i) an improved treatment for satellite galaxies that features a non-instantaneous stripping of the cold and hot gas components associated with galaxies being accreted onto larger systems, and a careful tracing of the angular momentum exchanges between different baryonic components (Xie et al. 2020; ii) an improved modelling of cold gas accretion on super-massive black holes and an explicit implementation of quasar winds (Fontanot et al. 2020). We have shown that the combination of these implementations leads to a remarkable agreement between our model predictions and the observed distributions of sSFR in the local Universe (see our Fig. 1). This represents a significant improvement over previous versions of our model that were either predicting an excess of low-mass quiescent galaxies (the GAEA2020F model), or a non negligible offset towards larger sSFR values for massive galaxies and the lack of a clear bimodality at intermediate galaxy masses (the GAEA2020X model).

As explained in Sect. 3, the main reasons for the improvement with respect to our previous model renditions are the following: at low galaxy stellar masses, the updated treatment adopted for satellite galaxies (as implemented in Xie et al. 2020) reduces drastically the significant excess of quenched galaxies found in the GAEA2020F model. At intermediate to large galaxy stellar masses, the improvements are mainly driven by the inclusion of AGN feedback in the form of quasar winds, with a non negligible role played by an overall better treatment of the density threshold for star formation.

6.1. Impact of different selections for passive galaxies

Throughout this work, we have used a classification of passive galaxies based on the sSFR. While this is appropriate for the local Universe, the observational selection of quiescent galaxies at z > 0 is typically based on a colour-colour selection. In particular, most previous work relies on the rest-frame U − V versus V − J diagram (e.g. Williams et al. 2009; Whitaker et al. 2011; Muzzin et al. 2013). This preference for a photometric selection is driven by the fact that rest-frame colours can be easily calculated for any galaxy sample, while an accurate estimate of the sSFR has stronger requirements on the depth and wavelength coverage of the data. Given that the division into two populations is motivated by the existence of a clear bimodality in physical and observable properties of galaxies, the expectation is that a different selection for model galaxies does not have a strong impact on the results discussed above. However, it is useful and interesting to quantify explicitly the impact of a different galaxy selection on the results discussed in the previous sections.

Figure 9 shows, in each panel, the distribution of all model galaxies (grey shaded regions) at the redshift indicated in the legend. The dashed lines indicate the box regions of the UVJ diagrams that are typically used to distinguish between star-forming and quiescent galaxies. The black contours show the regions enclosing 90, 60, and 30% of the model galaxies (increasing thickness of the lines). Red contours show the regions occupied by galaxies that are classified as passive on the basis of their sSFR, as used in the previous sections. For this analysis, we have only used a (representative) sub-volume of the Millennium Simulation (about 10%). In each panel, we also give the total number of model galaxies used and the following information:

thumbnail Fig. 9.

U − V versus V − J diagram for model galaxies at different cosmic epochs. Grey shaded regions show the distributions of all model galaxies more massive than ∼109 M⊙ at each redshift. Contours mark the regions enclosing 90, 60, and 30% of the galaxies. The red contours show the corresponding contour levels for galaxies that are classified as quiescent according to their sSFR. Dashed lines indicate the regions typically considered for a photometric selection of passive galaxies, in different studies. The fractions given in each panel indicate the completeness and contamination of the passive and active samples selected using the UVJ diagram (see text for details).

– on the left of each panel, in red, the completeness: the fractions of quiescent galaxies (on the basis of their sSFR) that fall in the passive region of the diagram. We also give the contamination in black, that is, the fraction of star-forming galaxies that fall in the same passive region of the diagram;

– on the right in each panel, we give the same quantities but for the active region of the diagram.

We give three different estimates of completeness and contamination associated with different colour cuts widely adopted in the literature. Although slightly different selections have been proposed in the literature, these do not lead to significant differences in the completeness and contamination that we can estimate using our model. Some differences are found only at the two highest redshifts considered, where the selection suggested by Whitaker et al. (2011) corresponds to systematically larger completeness fraction with respect to the selection used in Muzzin et al. (2013), with no increase in contamination. The completeness of the passive sample selected using the UVJ diagram is very high, up to z ∼ 1.5, decreasing to ∼60% at z ∼ 2 and to only ∼20% at z ∼ 3. We note that these numbers are based on a sample of model galaxies including all galaxies down to ∼109 M⊙, which is typically (much) below the completeness limit of the observational data at intermediate-to-high redshift. We show in Fig. 10 the same U − V versus V − J diagram, but this time considering only model galaxies with galaxy stellar mass larger than ∼1010 M⊙, and for the three highest redshifts considered. When considering only galaxies more massive than ∼1010 M⊙, the completeness of the UVJ classification goes up to ∼70% at z ∼ 2 and 20 − 40% at z ∼ 3 depending on the specific colour cuts adopted. We stress that, at these redshifts, current observational samples are typically complete at even larger galaxy stellar masses.

thumbnail Fig. 10.

As in Fig. 9 but only considering model galaxies more massive than ∼1010 M⊙ and for the three highest redshifts considered.

Considered the numbers given above, and the fact that the colour cuts applied in different observational studies are often slightly modified to account for example for slight zero-point differences of the available photometry, we can conclude that a different selection of quiescent galaxies would not have a significant impact on the agreement between model predictions and observational estimates discussed in previous Sections up to z ∼ 2. At higher redshifts, the UVJ selection fails to identify a significant fraction of galaxies that have very low sSFR values (see also Fig. 1 in Xie et al. 2024). While previous studies have demonstrated that this classical colour-colour selection is indeed incomplete at z > 3, the quoted “failure rate” is smaller than what found for our model. We note that several model ingredients or assumptions might have an impact on synthetic colours, for example assumptions about the stellar initial mass function, stellar evolutionary tracks, and dust. In future work, we plan to carry out a more detailed comparison with observational data by mimicking the same colour-colour selections adopted in observational studies, extending the analysis to alternative selections that have been recently employed at high-z.

6.2. A turn-over of the number densities of low-mass quiescent galaxies

In Sect. 4, we have shown that all models considered in this study predict a turn-over of the quiescent galaxy stellar mass function at low masses (see our Fig. 6). While this appears to be in agreement with recent observational findings (Santini et al. 2022; Weaver et al. 2023), the value at which this turnover occurs for our model is very close to the resolution limit of the Millennium Simulation. In order to quantify the impact of resolution, we have run our GAEA2023 model on the MillenniumII Simulation (Boylan-Kolchin et al. 2009). As mentioned above, this simulation adopts the same cosmological model of the Millennium but has 125 times better mass resolution, allowing to resolve galaxies down to stellar masses ∼108 M⊙.

Figure 11 shows the predicted galaxy stellar mass functions for all (black and grey lines for model data and blue symbols for observations) and quiescent (dark red and red lines for models and red symbols for observational data) galaxies. Solid lines correspond to the Millennium Simulation while thick dashed lines show the corresponding predictions for the Millennium II. All model predictions assume an observational uncertainty on galaxy stellar masses of 0.25 dex. Symbols with error bars show different available observational estimates, as indicated in the legend. Our GAEA2023 model exhibits a very good convergence over the entire galaxy mass range considered for the Millennium Simulation when considering all model galaxies. The low-mass end of the quiescent galaxy stellar mass function is somewhat more affected by resolution. However, a clear turn-over is still evident at galaxy stellar masses between ∼109 and ∼1010 M⊙. Below this mass limit, our model appears to over-predict the number densities measured for passive galaxies, and even more so when considering predictions from the Millennium II simulation. As noted earlier, this mass regime can be significantly affected by a different selection of quiescent galaxies, with a colour-colour cut typically selecting less galaxies than a cut based on sSFR. We defer a more detailed comparison to a future study, where we will also push our model predictions to earlier cosmic epochs. The last two panels of Fig. 11 show that our new model under-predicts the number densities of the most massive quiescent galaxies with respect to observational estimates by Weaver et al. (2023). We look into this in more detail below.

thumbnail Fig. 11.

Galaxy stellar mass function for all galaxies (grey and lack) and for quiescent galaxies selected on the basis of their sSFR (red and dark red) as predicted from the GAEA2023 model run on the Millennium Simulation (solid lines) and on the Millennium II (thick dashed lines). Symbols with error bars show different observational predictions as indicated in the legend. All model predictions shown in this figure assume a stellar mass uncertainty of 0.25 dex. The redshift bins used for the model have been optimized to match those employed by Weaver et al. (2023), that are slightly offset from that adopted in Santini et al. (2022) as indicated in the legend.

6.3. The cosmic stellar mass density evolution

The number density of massive passive galaxies at early cosmic epochs has gathered significant attention in recent years. In fact, the assembly of massive galaxies in a relatively short cosmic time (∼1.5 − 2 billion years when focussing on z > 3) places strong constraints on galaxy formation models. The topic has been subject to an intense debate over the years, and the debate has been revived by recent studies based on JWST observations, showing an excess of massive quiescent galaxies at z > 3 with respect to previous estimates (e.g. Carnall et al. 2023; Valentino et al. 2023). Predictions from our GAEA2023 model are shown in Fig. 12, together with a collection of recent pre- and post-JWST estimates. We show model predictions corresponding to different stellar mass cuts, as indicated in the legend, and assuming an uncertainty on stellar mass of 0.25 dex. The thick solid line is the one closest to the cut adopted in all observational studies considered. While our model is in quite good agreement with observational estimates based on pre-JWST data, the figure confirms that it falls short of massive galaxies at z > 3 when considering the most recent observations. The grey lines show model predictions obtained considering 125 independent sub-volumes of 100 Mpc h−1 on a side. These lines show that the expected cosmic variance is large, which is not surprising considered that these are extreme and rare galaxies. In fact, when taking the expected variance into account, our model predictions are consistent with estimates by Valentino et al. (2023), while still significantly below observational estimates by Carnall et al. (2023). As for most of the analysis presented in previous sections, we have selected quiescent galaxies considering a time-dependent sSFR cut. However, results do not vary significantly when considering a flat sSFR cut corresponding to 10−10 yr−1.

thumbnail Fig. 12.

Cumulative number density of massive quiescent galaxies as a function of redshift. Lines show model predictions corresponding to stellar mass cuts of 1010 (dashed), 1010.5 (solid – this is the closest to the observational measurements), and 1011 M⊙ (dotted). Symbols with error bars show different observational estimates, as indicated in the legend. Measurements based on JWST are indicated in green (Carnall et al. 2023) and magenta (Valentino et al. 2023). The grey lines show the cumulative number density evolution estimated in 125 independent sub-volumes of 100 Mpc h−1 on a side. All model predictions assume an uncertainty on galaxy stellar masses of 0.25 dex.

The comparison shown in Fig. 12 should be considered with caution: the data shown correspond to different selections of quiescent galaxies, and rather different volumes. Estimates by Carnall et al. (2023) are based on a total effective area of ∼30 arcmin2, while Valentino et al. (2023) base their estimates on 11 JWST fields with publicly available observations covering a total effective sky area of ∼145 arcmin2. The latter roughly correspond, at 3 < z < 5, to a total volume of 3.8 × 106 Mpc3, approximately equal to the individual sub-volumes considered to show the expected cosmic variance from our model. Valentino et al. (2023) also find a significant field-to-field variation in the estimated number densities, up to a factor 2−3. The presence of over-densities in some of the areas studied could further decrease the tension between our model predictions and the JWST estimates considered in Fig. 12. As mentioned, the data shown in Fig. 12 are based on different colour-colour selections and these could have a non negligible impact at these redshift. Uncertainties in photometric redshifts could also have an important impact on the estimated number densities. Finally, it is also important to bear in mind that the Millennium Simulation considered in this work has only 27 and 20 snapshots at z > 3 and z > 5, respectively. It remains to be demonstrated that our model results are stable on dark matter merger trees constructed on such a low number of snapshots. We defer this to a future work.

In order to put our model predictions in a context, we compare them to those from independent recently published theoretical models in Fig. 13. Specifically, we compare our model predictions with results from the latest version of the SHARK semi-analytic model (in blue, Lagos et al. 2024), from TNG100 and TNG300 (light and dark green respectively, Nelson et al. 2019), two different versions of Simba (magenta and purple, Davé et al. 2019; Hough et al. 2023), and EAGLE (brown, McAlpine et al. 2016 and references therein). We show only number densities corresponding to quiescent galaxies more massive than 1010.5 M⊙. We have considered, in all cases, two flat cuts in sSFR (top and bottom panels, as indicated in the legend). Predictions shown in this figure do not assume any uncertainty in the physical properties used (i.e., intrinsic values of SFR and stellar mass are used). The models considered cover a range of number densities, at all redshifts. As discussed above, this is partially due to different simulated volumes. However, based on results shown in Fig. 12, volume differences cannot be the main drivers of the differences between the models shown in Fig. 13: TNG100 and TNG300, based on boxes of 100 and 300 Mpc respectively, predict number density evolutions that are much closer to each other than to predictions from any other simulation considered. The same is true for Simba and Simba-c (corresponding to boxes of size ∼150 and ∼74 Mpc on a side, respectively). Predictions from EAGLE correspond to a simulated volume comparable to that of TNG100 and Simba, while predictions from SHARK are based on a box of ∼300 Mpc, comparable to that of TNG300. While this is smaller than the volume of the Millennium Simulation, it is larger than the sub-boxes considered in Fig. 12. In summary, we believe that the differences shown in Fig. 13 are largely driven by a different treatment of the baryonic processes.

thumbnail Fig. 13.

Cumulative number density of massive quiescent galaxies as a function of redshift. Predictions from GAEA2023 (shown by the solid black line in each panel) are compared with predictions from independent recently published models. The latter are shown as coloured lines, as indicated in the legend in the top panel. All model predictions correspond to galaxies more massive than 1010.5 M⊙ and do not include any uncertainty on the properties considered (intrinsic predictions are shown). Top and bottom panels correspond to two different (flat) cuts in sSFR to select quiescent galaxies, as indicated in the bottom left of each panel.

At z = 0, the EAGLE simulation predicts the lowest number density of quiescent galaxies. GAEA2023 and the TNG simulations predict the largest values of number densities – these are larger than those predicted by EAGLE by about one order of magnitude. The number densities of massive quiescent galaxies decrease rapidly with increasing redshift, for all models. There are no massive quiescent galaxies above z ∼ 3.7 in TNG100, and above z ∼ 4 in TNG300. The decrease is smoother for the EAGLE simulation that however also predicts number densities of massive quiescent galaxies of only ∼10−6 Mpc−3 at z = 3 when considering a sSFR cut < 10−10 yr−1. Among all models considered, GAEA2023 predicts the largest number densities at z > 3.

7. Summary and conclusions

In this paper, we present an updated version of our GAlaxy Evolution and Assembly (GAEA) theoretical model of galaxy formation (De Lucia et al. 2014; Hirschmann et al. 2016). In particular, the latest version of the model (GAEA2023) now combines for the first time (i) an updated treatment of AGN feedback that includes an improved modelling of cold gas accretion on super-massive black holes and an explicit implementation of AGN-driven outflows (Fontanot et al. 2020); and (ii) an improved modelling of both cold and hot gas stripping from satellite galaxies (Xie et al. 2020). We have shown that our latest model version predicts specific star formation rate (sSFR) distributions, for galaxies of different stellar mass, that are in remarkable agreement with observational measurements in the local Universe.

Comparing with previous versions of our model, we can characterize the physical drivers of the improved agreement with data: the updated treatment of satellite galaxies, featuring a non-instantaneous gas stripping, leads to relatively long quenching timescales for low-mass galaxies and a sSFR distribution with a single peak that coincides with the value observed in the local Universe. The inclusion of quasar winds leads to a rapid suppression of the star formation, that manifests itself as a clear bimodality in the sSFR distributions for galaxies with stellar masses between 109.5 to 1011 M⊙, and no excess of star formation for more massive galaxies. We find that the improved treatment of the density threshold for star formation, resulting from our explicit modelling of the partition between atomic and molecular hydrogen (Xie et al. 2017, 2020), does play an important role for the most massive galaxies by avoiding an extremely bursty and unrealistic behaviour for their star formation histories.

Predictions from our updated model version are also in quite nice agreement with observational measurements of the quenched fractions up to z  ∼  3 − 4. This again improves significantly over our previous models that were predicting either too-low fractions of quiescent galaxies at z > 1 (Xie et al. 2020, see also De Lucia et al. 2019), or an important excess of quenched galaxies below stellar masses ∼1010 M⊙ (Fontanot et al. 2020). We show that all models discussed in this work predict a turn-over of the number densities of quiescent galaxies at low-masses, with the exact location of this turn-over depending on the model considered and being closest to observational estimates for our latest model version. This is a robust prediction of our models, and it is not affected by limited resolution. Deeper observations, that can attempt to constrain the location of this turnover at intermediate to high redshift, can therefore provide important constraints on the physical processes regulating the number densities of low-mass quiescent galaxies, and in particular on the relative importance of stellar feedback and environmental processes in suppressing star formation in these systems.

As discussed in recent work (e.g. Kukstas et al. 2023), state-of-the-art theoretical models struggle to reproduce the observed quenching of galaxies in environments that are satellites dominated, that is, in galaxy clusters. Our updated model is in remarkable agreement with observational estimates at z ∼ 1, when the (naturally large) expected halo-to-halo variation is taken into account. However, a comparison with available observational measurements of the quiescent galaxy mass function suggests that the model still slightly over-predicts the fractions of quenched galaxies at the lowest stellar masses sampled, and under-predicts the quenched number densities around the knee. That is, the overall shape of the galaxy stellar mass function for quiescent and star-forming galaxies is very similar at intermediate galaxy stellar masses, in contrast with available data. This regime is sensitive to both stellar and AGN feedback, so accurate measurements for larger samples of clusters might provide very useful constraints on these physical processes.

The bulk of the analysis just summarized is based on an evolving cut in sSFR to select quiescent galaxies, while most observational estimates at z > 0 are based on some colour-colour selection. Using a sub-volume of the simulation, we show that a UVJ colour-colour selection would not have a significant impact on the discussed agreement between model predictions and observational estimates up to z ∼ 3. At larger redshift, a classical UVJ diagram fails to identify increasingly larger fractions of galaxies that are intrinsically quiescent (i.e., have very low levels of sSFR) in GAEA2023. We also demonstrate that a different selection might impact significantly the number densities of quiescent galaxies at stellar masses < 1010 M⊙. In future studies, we plan to carry out a more detailed comparison with observational data by mimicking the adopted colour-colour selections and extending the analysis to higher redshifts than considered in the present work.

Our latest model predicts number densities of massive quiescent galaxies at z > 3 that are larger than those predicted by a number of recently published state-of-the-art models. Yet, our model predictions appear to be below JWST observational measurements (Carnall et al. 2023; Valentino et al. 2023). We show that the expected cosmic variance is very large, and it can easily accomodate at least some of the recently published measurements. More careful comparisons, based on the same criteria to select quiescent galaxies, are needed. Caution is also needed when interpreting these results due to the fact that the simulations used in this work do not have a fine sampling of the dark matter merger trees at these early cosmic epochs. In future work, we will re-address this issue by taking advantage of a larger simulation with much finer time sampling at high redshift.

1

Details about the model, and access to a selection of data products, can be found at: https://sites.google.com/inaf.it/gaea/

5

This is the mass within a radius corresponding to an overdensity 200 times the critical density of the Universe.

Acknowledgments

We acknowledge the use of INAF-OATs computational resources within the framework of the CHIPP project (Taffoni et al. 2020). We thank Claudia Lagos, Romeel Davé, and Jakub Szpila for making results from SHARK v2.0, Simba, and Simba-c available in electronic format. We acknowledge stimulating discussions with Adam Muzzin, Veronica Strazzullo, and John Weaver. M.H. acknowledges funding from the Swiss National Science Foundation (SNSF) via a PRIMA grant PR00P2-193577 “From cosmic dawn to high noon: the role of BHs for young galaxies”. We acknowledge the anonymous referee and Veronica Strazzullo, for careful reading of our manuscript and many useful suggestions that helped improving the clarity of our presentation.

References

  1. Andreani, P., Miyamoto, Y., Kaneko, H., et al. 2020, A&A, 643, L11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Bahé, Y. M., Barnes, D. J., Dalla Vecchia, C., et al. 2017, MNRAS, 470, 4186 [Google Scholar]
  3. Baldry, I. K., Glazebrook, K., Brinkmann, J., et al. 2004, ApJ, 600, 681 [Google Scholar]
  4. Balogh, M. L., Gilbank, D. G., Muzzin, A., et al. 2017, MNRAS, 470, 4168 [NASA ADS] [CrossRef] [Google Scholar]
  5. Balogh, M. L., van der Burg, R. F. J., Muzzin, A., et al. 2021, MNRAS, 500, 358 [Google Scholar]
  6. Benson, A. J., Bower, R. G., Frenk, C. S., et al. 2003, ApJ, 599, 38 [NASA ADS] [CrossRef] [Google Scholar]
  7. Biviano, A., van der Burg, R. F. J., Balogh, M. L., et al. 2021, A&A, 650, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Blanton, M. R., Hogg, D. W., Bahcall, N. A., et al. 2003, ApJ, 594, 186 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bower, R. G., Benson, A. J., Malbon, R., et al. 2006, MNRAS, 370, 645 [Google Scholar]
  10. Boylan-Kolchin, M., Springel, V., White, S. D. M., Jenkins, A., & Lemson, G. 2009, MNRAS, 398, 1150 [Google Scholar]
  11. Brammer, G. B., Whitaker, K. E., van Dokkum, P. G., et al. 2009, ApJ, 706, L173 [Google Scholar]
  12. Bravo, M., Robotham, A. S. G., Lagos, C. D. P., et al. 2023, MNRAS, 522, 4481 [NASA ADS] [CrossRef] [Google Scholar]
  13. Brinchmann, J., Charlot, S., White, S. D. M., et al. 2004, MNRAS, 351, 1151 [Google Scholar]
  14. Carnall, A. C., McLeod, D. J., McLure, R. J., et al. 2023, MNRAS, 520, 3974 [NASA ADS] [CrossRef] [Google Scholar]
  15. Conroy, C., Gunn, J. E., & White, M. 2009, ApJ, 699, 486 [Google Scholar]
  16. Croton, D. J., Springel, V., White, S. D. M., et al. 2006, MNRAS, 365, 11 [Google Scholar]
  17. Davé, R., Anglés-Alcázar, D., Narayanan, D., et al. 2019, MNRAS, 486, 2827 [Google Scholar]
  18. Davies, L. J. M., Robotham, A. S. G., Lagos, C. D. P., et al. 2019, MNRAS, 483, 5444 [Google Scholar]
  19. De Lucia, G., & Blaizot, J. 2007, MNRAS, 375, 2 [Google Scholar]
  20. De Lucia, G., Kauffmann, G., & White, S. D. M. 2004, MNRAS, 349, 1101 [Google Scholar]
  21. De Lucia, G., Tornatore, L., Frenk, C. S., et al. 2014, MNRAS, 445, 970 [Google Scholar]
  22. De Lucia, G., Hirschmann, M., & Fontanot, F. 2019, MNRAS, 482, 5041 [Google Scholar]
  23. De Lucia, G., Xie, L., Fontanot, F., & Hirschmann, M. 2020, MNRAS, 498, 3215 [NASA ADS] [CrossRef] [Google Scholar]
  24. Donnari, M., Pillepich, A., Nelson, D., et al. 2021, MNRAS, 506, 4760 [NASA ADS] [CrossRef] [Google Scholar]
  25. Dubois, Y., Gavazzi, R., Peirani, S., & Silk, J. 2013, MNRAS, 433, 3297 [CrossRef] [Google Scholar]
  26. Fontanot, F., De Lucia, G., Monaco, P., Somerville, R. S., & Santini, P. 2009, MNRAS, 397, 1776 [NASA ADS] [CrossRef] [Google Scholar]
  27. Fontanot, F., De Lucia, G., Hirschmann, M., et al. 2017a, MNRAS, 464, 3812 [NASA ADS] [CrossRef] [Google Scholar]
  28. Fontanot, F., Hirschmann, M., & De Lucia, G. 2017b, ApJ, 842, L14 [Google Scholar]
  29. Fontanot, F., De Lucia, G., Hirschmann, M., et al. 2020, MNRAS, 496, 3943 [CrossRef] [Google Scholar]
  30. Fontanot, F., Calabrò, A., Talia, M., et al. 2021, MNRAS, 504, 4481 [CrossRef] [Google Scholar]
  31. Fontanot, F., La Barbera, F., De Lucia, G., et al. 2024, A&A, 686, A302 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Fossati, M., Wilman, D. J., Mendel, J. T., et al. 2017, ApJ, 835, 153 [Google Scholar]
  33. Franx, M., van Dokkum, P. G., Förster Schreiber, N. M., et al. 2008, ApJ, 688, 770 [NASA ADS] [CrossRef] [Google Scholar]
  34. Guo, Q., White, S., Angulo, R. E., et al. 2013, MNRAS, 428, 1351 [NASA ADS] [CrossRef] [Google Scholar]
  35. Haynes, M. P., Giovanelli, R., Martin, A. M., et al. 2011, AJ, 142, 170 [Google Scholar]
  36. Hirschmann, M., De Lucia, G., Wilman, D., et al. 2014, MNRAS, 444, 2938 [NASA ADS] [CrossRef] [Google Scholar]
  37. Hirschmann, M., De Lucia, G., & Fontanot, F. 2016, MNRAS, 461, 1760 [Google Scholar]
  38. Hough, R. T., Rennehan, D., Kobayashi, C., et al. 2023, MNRAS, 525, 1061 [NASA ADS] [CrossRef] [Google Scholar]
  39. Kauffmann, G., Colberg, J. M., Diaferio, A., & White, S. D. M. 1999, MNRAS, 303, 188 [NASA ADS] [CrossRef] [Google Scholar]
  40. Kauffmann, G., Heckman, T. M., White, S. D. M., et al. 2003, MNRAS, 341, 33 [Google Scholar]
  41. Keres, D., Yun, M. S., & Young, J. S. 2003, ApJ, 582, 659 [NASA ADS] [CrossRef] [Google Scholar]
  42. Kukstas, E., Balogh, M. L., McCarthy, I. G., et al. 2023, MNRAS, 518, 4782 [Google Scholar]
  43. Lagos, C. D. P., Cora, S. A., & Padilla, N. D. 2008, MNRAS, 388, 587 [Google Scholar]
  44. Lagos, C. D. P., Bravo, M., Tobar, R., et al. 2024, MNRAS, 531, 3551 [CrossRef] [Google Scholar]
  45. Lustig, P., Strazzullo, V., Remus, R.-S., et al. 2023, MNRAS, 518, 5953 [Google Scholar]
  46. McAlpine, S., Helly, J. C., Schaller, M., et al. 2016, Astron. Comput., 15, 72 [NASA ADS] [CrossRef] [Google Scholar]
  47. Mitchell, P. D., Lacey, C. G., Baugh, C. M., & Cole, S. 2013, MNRAS, 435, 87 [NASA ADS] [CrossRef] [Google Scholar]
  48. Muzzin, A., Marchesini, D., Stefanon, M., et al. 2013, ApJ, 777, 18 [NASA ADS] [CrossRef] [Google Scholar]
  49. Nelson, D., Springel, V., Pillepich, A., et al. 2019, Comput. Astrophys. Cosmol., 6, 2 [Google Scholar]
  50. Puchwein, E., & Springel, V. 2013, MNRAS, 428, 2966 [NASA ADS] [CrossRef] [Google Scholar]
  51. Salim, S., Rich, R. M., Charlot, S., et al. 2007, ApJS, 173, 267 [NASA ADS] [CrossRef] [Google Scholar]
  52. Santini, P., Castellano, M., Fontana, A., et al. 2022, ApJ, 940, 135 [NASA ADS] [CrossRef] [Google Scholar]
  53. Shen, X., Hopkins, P. F., Faucher-Giguère, C.-A., et al. 2020, MNRAS, 495, 3252 [Google Scholar]
  54. Sijacki, D., Springel, V., Di Matteo, T., & Hernquist, L. 2007, MNRAS, 380, 877 [NASA ADS] [CrossRef] [Google Scholar]
  55. Springel, V., White, S. D. M., Jenkins, A., et al. 2005, Nature, 435, 629 [Google Scholar]
  56. Taffoni, G., Becciani, U., Garilli, B., et al. 2020, in Astronomical Data Analysis Software and Systems XXIX, eds. R. Pizzo, E. R. Deul, J. D. Mol, J. de Plaa, & H. Verkouter, ASP Conf. Ser., 527, 307 [NASA ADS] [Google Scholar]
  57. Valentino, F., Brammer, G., Gould, K. M. L., et al. 2023, ApJ, 947, 20 [NASA ADS] [CrossRef] [Google Scholar]
  58. van der Burg, R. F. J., Rudnick, G., Balogh, M. L., et al. 2020, A&A, 638, A112 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Wang, J., De Lucia, G., Kitzbichler, M. G., & White, S. D. M. 2008, MNRAS, 384, 1301 [CrossRef] [Google Scholar]
  60. Wang, E., Wang, H., Mo, H., et al. 2018, ApJ, 864, 51 [NASA ADS] [CrossRef] [Google Scholar]
  61. Weaver, J. R., Kauffmann, O. B., Ilbert, O., et al. 2022, ApJS, 258, 11 [NASA ADS] [CrossRef] [Google Scholar]
  62. Weaver, J. R., Davidzon, I., Toft, S., et al. 2023, A&A, 677, A184 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  63. Weinmann, S. M., van den Bosch, F. C., Yang, X., et al. 2006, MNRAS, 372, 1161 [CrossRef] [Google Scholar]
  64. Wetzel, A. R., Tinker, J. L., & Conroy, C. 2012, MNRAS, 424, 232 [NASA ADS] [CrossRef] [Google Scholar]
  65. Whitaker, K. E., Labbé, I., van Dokkum, P. G., et al. 2011, ApJ, 735, 86 [Google Scholar]
  66. Williams, R. J., Quadri, R. F., Franx, M., van Dokkum, P., & Labbé, I. 2009, ApJ, 691, 1879 [NASA ADS] [CrossRef] [Google Scholar]
  67. Xie, L., De Lucia, G., Hirschmann, M., Fontanot, F., & Zoldan, A. 2017, MNRAS, 469, 968 [Google Scholar]
  68. Xie, L., De Lucia, G., Hirschmann, M., & Fontanot, F. 2020, MNRAS, 498, 4327 [NASA ADS] [CrossRef] [Google Scholar]
  69. Xie, L., De Lucia, G., Fontanot, F., et al. 2024, ApJ, 966, L2 [NASA ADS] [CrossRef] [Google Scholar]
  70. Yang, X., Mo, H. J., van den Bosch, F. C., et al. 2007, ApJ, 671, 153 [Google Scholar]

Appendix A: Model calibration

As mentioned in Section 2, the primary observables that have been considered to recalibrate our model are: the observed galaxy stellar mass function and its evolution up to z ∼ 3, the HI and H2 mass functions in the local Universe, the evolution of the AGN luminosity function up to z ∼ 4. For completeness, we show in this Appendix the main calibration plots, including all three models that have been presented in this paper.

Fig. A.1 shows the predicted galaxy stellar mass function at different redshifts, from all three models considered in this paper, compared with a compilation of observational data as in Hirschmann et al. (2016). Fig. A.2 shows the HI galaxy mass function in the left panel, compared with observational estimates by Haynes et al. (2011), and the H2 galaxy mass function, compared with observational estimates by Keres et al. (2003) and Andreani et al. (2020). Finally, Fig. A.3 shows the predicted AGN luminosity function from the GAEA2023 and GAEA2020F models, compared with observational estimates by Shen et al. (2020, we consider an uncertainty of +/-0.3 dex).

thumbnail Fig. A.1.

Galaxy stellar mass functions at different redshifts. Solid lines show predictions all three models considered in this paper, while symbols show the same compilation of observational data considered in Hirschmann et al. (2016).

thumbnail Fig. A.2.

HI (left panel) and H2 (right panel) mass function at z = 0, compared with observational estimates. The solid black line corresponds to our GAEA2023 model, while the dark orange line shows model predictions from Xie et al. (2020).

thumbnail Fig. A.3.

AGN bolometric luminosity function at different redshifts. Solid lines show predictions from the GAEA2023 and GAEA2020F models, while the solid lines correspond to observational estimates by Shen et al. (2020) assuming an uncertainty of +/-0.3 dex.

All Tables

Table 1.

Values of the relevant parameters used in Xie et al. (2020) and Fontanot et al. (2020), and the corresponding values adopted in the model considered in this study.

In the text

All Figures

thumbnail Fig. 1.

Specific star formation rate distribution as predicted by the models indicated in the legend, and compared to observational estimates based on SDSS DR8 (grey shaded histograms – see text for details). Each panel corresponds to a different bin in galaxy stellar mass, as indicated in the top-right legend. The histograms corresponding to the GAEA2020X and GAEA2020F models have been slightly shifted (by 0.02 dex) to the left and right with respect to the GAEA2023 model for clarity.
In the text
thumbnail Fig. 2.

Passive fractions for all galaxies (left panel), central galaxies (middle panel), and satellites (right panels). Symbols correspond to observational measurements. In particular: squares correspond to estimates from Davies et al. (2019) based on the GAMA survey, and adopt a log(sSFR) threshold of −10.5 yr−1 to identify passive fractions. Filled and empty circles correspond to measurements based on SDSS DR8 and are obtained assuming a threshold of −11 yr−1 and −10.6 yr−1, respectively. Solid lines correspond to predictions from different models, as indicated in the legend, assuming a log(sSFR) threshold of −11.0 yr−1 to select passive galaxies. The grey regions show the minimum-maximum fractions obtained when considering predictions from the GAEA2023 model from 125 independent boxes, each of 100 Mpc h−1 on a side. The dashed line shows predictions from the GAEA2020X model corresponding to a log(sSFR) threshold of −10.6 yr−1, as adopted in Xie et al. (2020).
In the text
thumbnail Fig. 3.

Fraction of passive centrals (solid lines) and satellites (dashed lines) as predicted by the models considered in this study. Symbols with error bars (squares are for centrals and triangles for satellites) show observational estimates by Fossati et al. (2017), based on five CANDELS/3D-HST fields.
In the text
thumbnail Fig. 4.

Passive fraction at different cosmic epochs. Lines correspond to predictions from the models considered in this study. The shaded area and red symbols with error bars show observational estimates from Muzzin et al. (2013). Specifically, the shaded areas show the 1/Vmax measurements, while the open symbols with error bars show the maximum-likelihood measurements (see bottom panels of Fig. 6 in Muzzin et al. 2013). The blue symbols with error bars correspond to the observational estimates by Weaver et al. (2023). A sSFR cut has been employed to select model passive galaxies (see text for details), while the observational samples are based on different colour-colour selections, a SED fitting estimate of the galaxy stellar mass, and photometric redshifts estimates.
In the text
thumbnail Fig. 5.

Galaxy stellar mass function at different cosmic epochs. Symbols with error bars correspond to the observational estimates published in Muzzin et al. (2013, in this case we are showing only the 1/Vmax estimates published in this study) and Weaver et al. (2023). Solid lines correspond to predictions from the models considered in this study, as indicated in the legend. The dashed line in each panel shows, for the GAEA2023 model only, the corresponding predictions obtained assuming that the stellar masses have uncertainties that are Gaussian distributed with a width of 0.25 dex.
In the text
thumbnail Fig. 6.

As in Fig. 5, but for quiescent galaxies only.
In the text
thumbnail Fig. 7.

Galaxy stellar mass function for the star-forming (top panel) and quiescent galaxies (bottom panel) in clusters at 1 < z < 1.4. Solid lines show predictions from the models considered in this study, as indicated in the legend. Symbols with error bars show corresponding observational estimates published in van der Burg et al. (2020). The shaded regions show the area enclosed by the 10th and 90th percentiles of the distributions of number densities obtained, for each galaxy stellar mass bin, when considering 25 samples of simulated haloes that match the mass and redshift distributions of the GOGREEN sample (see text for details). All model lines shown in this figure assume an uncertainty in stellar mass of 0.25 dex.
In the text
thumbnail Fig. 8.

Passive fraction in the cluster environment at z ∼ 1. Lines show predictions from the models considered in this study, while symbols with error bars show the observational estimates from van der Burg et al. (2020). The dashed lines show the area enclosed by the 20th and 80th percentiles of the quiescent fractions obtained, for each galaxy stellar mass bin, when considering 25 samples of simulated haloes that match the mass and redshift distributions of the GOGREEN sample (see text for details). We show the dashed lines only for the GAEA2023 model for clarity. All model lines shown in this figure assume an uncertainty in stellar mass of 0.25 dex.
In the text
thumbnail Fig. 9.

U − V versus V − J diagram for model galaxies at different cosmic epochs. Grey shaded regions show the distributions of all model galaxies more massive than ∼109 M⊙ at each redshift. Contours mark the regions enclosing 90, 60, and 30% of the galaxies. The red contours show the corresponding contour levels for galaxies that are classified as quiescent according to their sSFR. Dashed lines indicate the regions typically considered for a photometric selection of passive galaxies, in different studies. The fractions given in each panel indicate the completeness and contamination of the passive and active samples selected using the UVJ diagram (see text for details).
In the text
thumbnail Fig. 10.

As in Fig. 9 but only considering model galaxies more massive than ∼1010 M⊙ and for the three highest redshifts considered.
In the text
thumbnail Fig. 11.

Galaxy stellar mass function for all galaxies (grey and lack) and for quiescent galaxies selected on the basis of their sSFR (red and dark red) as predicted from the GAEA2023 model run on the Millennium Simulation (solid lines) and on the Millennium II (thick dashed lines). Symbols with error bars show different observational predictions as indicated in the legend. All model predictions shown in this figure assume a stellar mass uncertainty of 0.25 dex. The redshift bins used for the model have been optimized to match those employed by Weaver et al. (2023), that are slightly offset from that adopted in Santini et al. (2022) as indicated in the legend.
In the text
thumbnail Fig. 12.

Cumulative number density of massive quiescent galaxies as a function of redshift. Lines show model predictions corresponding to stellar mass cuts of 1010 (dashed), 1010.5 (solid – this is the closest to the observational measurements), and 1011 M⊙ (dotted). Symbols with error bars show different observational estimates, as indicated in the legend. Measurements based on JWST are indicated in green (Carnall et al. 2023) and magenta (Valentino et al. 2023). The grey lines show the cumulative number density evolution estimated in 125 independent sub-volumes of 100 Mpc h−1 on a side. All model predictions assume an uncertainty on galaxy stellar masses of 0.25 dex.
In the text
thumbnail Fig. 13.

Cumulative number density of massive quiescent galaxies as a function of redshift. Predictions from GAEA2023 (shown by the solid black line in each panel) are compared with predictions from independent recently published models. The latter are shown as coloured lines, as indicated in the legend in the top panel. All model predictions correspond to galaxies more massive than 1010.5 M⊙ and do not include any uncertainty on the properties considered (intrinsic predictions are shown). Top and bottom panels correspond to two different (flat) cuts in sSFR to select quiescent galaxies, as indicated in the bottom left of each panel.
In the text
thumbnail Fig. A.1.

Galaxy stellar mass functions at different redshifts. Solid lines show predictions all three models considered in this paper, while symbols show the same compilation of observational data considered in Hirschmann et al. (2016).
In the text
thumbnail Fig. A.2.

HI (left panel) and H2 (right panel) mass function at z = 0, compared with observational estimates. The solid black line corresponds to our GAEA2023 model, while the dark orange line shows model predictions from Xie et al. (2020).
In the text
thumbnail Fig. A.3.

AGN bolometric luminosity function at different redshifts. Solid lines show predictions from the GAEA2023 and GAEA2020F models, while the solid lines correspond to observational estimates by Shen et al. (2020) assuming an uncertainty of +/-0.3 dex.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.