Unveiling the properties and origin of massive quenched galaxies at $z\ge2$ in the COLIBRE hydrodynamical simulations

TL;DR

This work uses the COLIBRE cosmological hydrodynamical simulations to study massive quenched galaxies (MQGs) at z≥2, addressing their abundance, formation and quenching timescales, and ISM properties. By implementing explicit cold gas physics, dust, and two AGN feedback modes (thermal and jet-inclusive hybrid), the authors show MQGs form rapidly in overdense environments and are quenched primarily by AGN feedback, with BH growth fueling powerful outflows that deplete dust and H2 reservoirs. The study finds COLIBRE predictions align with JWST measurements once observational uncertainties are convolved, and it highlights the roles of environment and AGN feedback in driving MQG evolution, including a sizable fraction that survive to z=0 and substantial rejuvenation episodes. These results offer robust, testable predictions for MQG properties and provide a pathway to reconcile observations with simulations, guiding future JWST/ALMA studies.

Abstract

JWST has uncovered a substantial population of Massive ($M_{\star} \gtrsim 10^{10 }\mathrm{M_{\odot}}$), Quenched Galaxies (MQGs) in the early Universe ($z \gtrsim 2$), whose properties challenge current galaxy formation models. In this paper, we examine this population of MQGs within the new COLIBRE cosmological hydrodynamical simulations. We report number densities and stellar mass functions in broad agreement with the latest observations. The predicted quenching and formation timescales are qualitatively consistent with observational inferences. Leveraging the state-of-the-art physics in COLIBRE, the model predicts that MQGs have dust and $\rm H_{2}$ fractions more than $1$ dex lower than their massive star-forming counterparts; while their sizes and kinematics remain broadly similar. We further explore the processes driving galaxies to become massive and quenched in COLIBRE, identifying active galactic nucleus (AGN) feedback as the primary quenching mechanism. Compared to star-forming galaxies of similar mass, MQGs host more massive black holes (BHs) and exhibit higher star formation efficiencies. These differences arise from their environments, particularly at local ($\rm 0.3\,cMpc$) to intermediate scales ($\rm 1.0\,cMpc$) before quenching, where overdense regions are associated with enhanced gas inflows, higher BH accretion and, hence, feedback power. We find that about $55\%$ of MQGs survive as the main progenitors of $z=0$ galaxies when they are selected at $z=3$, although up to $55\%$ experience rejuvenation episodes. Our results provide robust predictions for MQGs, show that tensions with observations are reduced when an effective observational uncertainty is forward-modelled, and clarify the mechanisms behind their origin.

Figures (21)

• Figure 1: Top panel: normalised SFHs ($\rm SFR/SFR_{max}$ versus lookback time, with redshift indicated by the top $x$-axis) of two example galaxies identified to be massive and quenched at $z=3$ in L200m6. The horizontal solid line corresponds to the peak SFR value, $\rm SFR_{max}$, while the dotted line shows the $\rm 0.1\times SFR_{max}$ threshold used to define the quenching timescale, $t_{\rm q}$. Arrows indicate $t_{\rm q}$ for each galaxy, defined as the latest time the galaxy's SFR drops below this threshold; as well as the SFR peak time, $t_{\rm SFR_{max}}$. The green curve represents a galaxy that, once it has crossed the threshold indicated by the dotted line, never rises above it again. The orange curve represents a galaxy that experiences several star formation episodes above the threshold after its SFR peak. Bottom panel: normalised stellar mass assembly histories ($M_{\star}/M_{\star, t=0}$ versus lookback time, with redshift indicated by the top $x$-axis) for the same galaxies. The horizontal solid line shows the maximum stellar mass value, $M_{\star,t=0}$; while the dotted line shows the $0.5\times M_{\star,t=0}$ threshold used to define the formation timescale, $t_{\rm 50}$. Arrows indicate $t_{\rm 50}$.
• Figure 2: Comoving number density of MQGs, defined by $M_{\star}>10^{10}$ and $\mathrm{sSFR}<0.2/t_{\rm age}$, as a function of redshift. Lines represent predictions from various COLIBRE boxes (spanning different volumes, resolutions, and AGN feedback models) in Table \ref{['tab:runs']}, with the selection at $2 \le z \le 8$. These are compared with the latest JWST observational estimates: black (smaller area) and dark red (larger area) points show spectroscopic measurements, while grey (smaller area) and red (larger area) points correspond to photometric measurements. Observational data are taken from carnall23valentino23alberts24nanayakkara25weibel25baker25zhang25stevenson25. Top panel: predicted simulation values. Bottom panel: $M_{\star}$ and SFR values in the simulations are convolved independently with a Gaussian distribution (mean 0, standard deviation $0.3$ dex) representing a reasonable error budget for these quantities.
• Figure 3: Comoving number density of MQGs, defined by $M_{\star}>10^{10}$ and $\mathrm{sSFR}<0.2/t_{\rm age}$, as a function of redshift. Light-blue line shows predictions from our fiducial L200m6 COLIBRE simulation, with the selection at $2 \le z \le 8$, compared to results from other galaxy formation and evolution models in the literature ( GAEA defines quenching using a fixed cut of $\rm sSFR = 10^{-10}yr^{-1}$). For all simulation results (except GAEA), $M_{\star}$ and SFR values in the simulations are convolved independently with a Gaussian distribution (mean 0, standard deviation $0.3$ dex) representing a reasonable error budget for these quantities. These are compared with the latest JWST observational estimates: black (smaller area) and dark red (larger area) points show spectroscopic measurements, while grey (smaller area) and red (larger area) points correspond to photometric measurements. Observational data are taken from carnall23valentino23alberts24nanayakkara25weibel25baker25zhang25stevenson25.
• Figure 4: SMF of MQGs, defined by $M_{\star}>10^{10}$ and $\mathrm{sSFR}<0.2/t_{\rm age}$, from the fiducial L200m6 simulation. Left panel: predicted simulation values at different redshifts $2 \le z \le 7$ (colour-coded as labelled in the left panel), compared with observational estimates from baker25 and stevenson25 within the same redshift range, as labelled. Right panel: the SMF at $z=3$ from COLIBRE in light blue (fiducial L200m6) compared with results from other galaxy formation and evolution models, where $M_{\star}$ and SFR values in the simulations are convolved independently with a Gaussian distribution (mean 0, standard deviation $0.3$ dex) representing a reasonable error budget for these quantities ( GAEA defines quenching using a fixed cut of $\rm sSFR = 10^{-10}\,yr^{-1}$ and its results are not convolved with potential observational errors). The horizontal dotted line denotes the threshold of 10 galaxies, below which the statistics become unreliable.
• Figure 5: SFH properties of MQGs as a function of $M_{\star}$. Thick solid lines show the median predictions from L200m6 at different redshifts $2 \le z \le 5$ (as labelled in the bottom panel), with shaded regions indicating the 16th and 84th percentile range. At $z=5$, MQGs are too few to split in stellar mass bins (9 galaxies). Thin solid lines show the corresponding medians from L400m7. These are compared to recent JWST spectroscopic measurements carnall24nanayakkara25baker25b, shown as symbols colour-coded by observed galaxy redshift, with the corresponding colorbar at the top. Top panel: formation times ($t_{50}$), defined in § \ref{['ssec:prop-sfh']}. Middle panel: quenching times ($t_{\rm q}$), defined in § \ref{['ssec:prop-sfh']}carnall24. Bottom panel: SFR peak ($\mathrm{SFR}_{\rm max}$) carnall24.