alphaFeedback

WorksPhysics and Astronomy

How Quickly Did the Giants Form? Constraining the Formation Timescales of Massive Quiescent Galaxies at Cosmic Noon Using Chemical Evolution Modelling

AI contributions
  • Gemini · multiple versions
  • Claude by Anthropic
View details

Abstract

Massive quiescent galaxies (MQGs; M⋆∼1011 M⊙M_\star \sim 10^{11} \, \mathrm{M_\odot}) at cosmic noon (z∼2z \sim 2) host αα-enhanced stellar populations ([Mg/Fe]∼0.4−0.5\mathrm{[Mg/Fe]} \sim 0.4 - 0.5) at near-solar [Fe/H], implying rapid formation histories that spectral energy distribution (SED) fitting alone constrains poorly. To determine the physical conditions and formation timescales that generate these abundances, we employ a one-zone chemical evolution model with nucleosynthetic yields calibrated against resolved stars in the Sculptor dwarf spheroidal. Comparing our model to a compilation of observed abundances, we find MQGs are consistent with rapid formation timescales of ΔtSF∼0.1−0.4 GyrΔt_\mathrm{SF} \sim 0.1 - 0.4 \, \mathrm{Gyr}. The progenitors must sustain high star-formation efficiencies (approximately the molecular gas depletion time; τ⋆∼0.05 Gyrτ_\star \sim 0.05 \, \mathrm{Gyr}) and retain most of their metal-enriched ejecta (η∼0η\sim 0), consistent with high-redshift dusty star-forming galaxies (DSFGs). However, our model favours moderately active progenitors with typical star-formation rates (SFRs) of 200−400 M⊙ yr−1200-400 \, \mathrm{M_\odot \, yr^{-1}}, matching recently identified `NIRCam-dark' galaxies more closely than the broader DSFG population. Sustained starbursts (SFR≳1000 M⊙ yr−1\mathrm{SFR} \gtrsim 1000 \, \mathrm{M_\odot \, yr}^{-1}) are disfavoured, since they would form most stellar mass before significant SN Ia enrichment and underproduce the observed [Fe/H]. The same framework also reproduces the central abundances of local early-type galaxies over longer timescales (ΔtSF∼0.6−1 GyrΔt_\mathrm{SF} \sim 0.6 - 1 \, \mathrm{Gyr}). Chemical abundances can therefore provide powerful constraints on the formation histories of MQGs, and incorporating them into SED-fitting routines will ensure inferred histories remain consistent with the observed chemistry.

Introduction

The chemical composition of a galaxy’s stars preserves a fossil record of its star formation and enrichment history (Tinsley 1980; Maiolino & Mannucci 2019). By interpreting the elemental abundance distributions of stellar populations, we can begin to reconstruct when, how rapidly, and from what sources a galaxy formed and assembled its mass.

Of the abundance ratios accessible in integrated-light spectroscopy, [α\alpha/Fe]1 is the most informative tracer of how rapidly a stellar population formed. α\alpha-elements are produced almost exclusively by core-collapse supernovae (CCSNe) on timescales of 1−10 Myr1 - 10\ \mathrm{Myr}, while Fe production is dominated by Type Ia supernovae (SNe Ia), which contribute on delayed timescales of 0.1−1.0 Gyr0.1 - 1.0\ \mathrm{Gyr} (Matteucci & Greggio 1986). As a result, stellar populations whose star formation completed before significant SN Ia enrichment could preserve enhanced [α\alpha/Fe], while more extended formation histories drive [α\alpha/Fe] toward solar values (e.g. Tinsley 1979). [α\alpha/Fe] hence acts as a chronometer of star formation duration, complementing [Fe/H] as a tracer of total enrichment (metallicity).

Massive quiescent galaxies (MQGs) observed near cosmic noon (z≳2z \gtrsim2) are a powerful target for this chemical analysis. Having ceased significant star formation, their stellar populations are largely undisturbed by ongoing enrichment leaving a clean record of the conditions under which they formed. Observed close to their quenching epoch, they are also caught before significant ex-situ accretion expected from dry-minor merging events (e.g. Naab et al. 2009), isolating the in-situ stellar population assembled during the galaxy’s formation. MQGs also dominate the high-mass end of the galaxy population, where the physical processes that regulate growth through the gas supply – feedback, accretion and quenching – operate most strongly and where current galaxy formation models remain most uncertain (Naab & Ostriker 2017). Understanding how MQGs formed the bulk of their stellar mass is therefore central to a broader picture of how the most massive structures in the Universe came to be.

Most of what is known about the chemistry of massive galaxies has been established not from MQGs directly, but from their assumed descendants in the local universe, early-type galaxies (ETGs; [171] (Patel et al. 2013)), whose proximity allows the spectroscopic detail needed to derive stellar abundances. Archaeological studies find the properties of their stellar populations correlate strongly with the stellar mass (M⋆M_\star) and velocity dispersion (σ\sigma), where the most massive, highest-σ\sigma galaxies host old stellar populations with systematically high metallicities and enhanced α\alpha-element abundances (typically [α/Fe]∼0.2[\alpha/\mathrm{Fe}] \sim0.2; Thomas et al. 2005; Graves et al. 2009a,b; Smith et al. 2009; Conroy et al. 2014; Greene et al. 2019). Even low-mass ETGs at z∼0z \sim0 show super-solar [Mg/Fe], with the enhancement growing systematically with velocity dispersion (Conroy et al. 2014). Three explanations were originally proposed by Faber et al. (1992) and Worthey et al. (1992): (i) variation in the timescales over which star formation proceeded, (ii) a non-universal initial mass function (IMF), and (iii) differential entrainment or loss of Mg and Fe in winds. However, the latter explanation is disfavoured due to the absence of plausible physical mechanisms for ejecting certain species (Matteucci 1994). Instead, variations in SF timescales are the favoured interpretation of [Mg/Fe] variations in both galactic and extragalactic archaeology (e.g. [166] (Matteucci 1994; Thomas et al. 1999; Thomas et al. 2005)); though, variation of the IMF in early-type galaxies remains a subject of ongoing debate (see Bastian et al. 2010; Smith 2020, for a review).

Whether the chemical patterns observed in local ETGs were already in place at the formation epoch of their MQG progenitors can only be addressed by direct observation at high redshift. Extending these abundance measurements to the quiescent population at z≳2z \gtrsim2, close to when they actually formed and quenched, has only recently become feasible with near-infrared spectroscopy from instruments such as JWST/NIRSpec and Keck/MOSFIRE. If MQGs evolved purely passively into local ETGs, their chemistry should mirror that of their descendants. Instead, observations of z∼2z \sim2 MQGs reveal more extreme α\alpha-enhancement, with [Mg/Fe] often reaching 0.4 − 0.5 dex (Kriek et al. 2016; Jafariyazani et al. 2020; Beverage et al. 2024; Carnall et al. 2024; Beverage et al. 2025; Jafariyazani et al. 2025; McConachie et al. 2026) compared to the typical ∼0.2\sim0.2 dex found in local ETGs. Taking [α/Fe][\alpha/\mathrm{Fe}] to be a measure of the duration of star formation, the enhanced [Mg/Fe] implies that these galaxies formed the bulk of their stellar mass on remarkably short timescales (≲200\lesssim200 Myr) before quenching (Kriek et al. 2016), a scenario reinforced by abundance patterns across the broader nucleosynthetic groups (Beverage et al. 2025). It remains an open question, however, whether this offset in abundances reflects the typical formation history of all local ETG progenitors, or only of an early, rapidly-forming subset that we preferentially observe at z≳2z \gtrsim2.

The case for rapid formation is supported by recent findings from JWST which reveal that the number density of MQGs in the early universe is substantially higher than predicted prior to launch (Valentino et al. 2023; Carnall et al. 2024; Nanayakkara et al. 2025; Merlin et al. 2025; Stevenson et al. 2026; Zhang et al. 2026), with individual galaxies in place as early as ∼500\sim500 Myr after the Big Bang according to their star-formation histories (SFHs) (Carnall et al. 2023; de Graaff et al. 2025; McConachie et al. 2026; Weibel et al. 2025). However, reproducing such galaxies requires extreme formation scenarios with baryon conversion efficiencies near the theoretical limit (Glazebrook et al. 2024; McConachie et al. 2026). Furthermore, existing cosmological and semi-analytic galaxy formation models are unable to reproduce the high number density of high-mass MQGs at z≳3z \gtrsim3 ([83] (Stevenson et al. 2026)), exposing a gap in our understanding of how the most massive galaxies form and quench so rapidly (Yang et al. 2026).

One route to closing the gap is to directly identify the progenitor population observationally, since the high number densities, SFHs and host halo masses of MQGs should impose strict criteria on their progenitors. For example, the formation times derived from the SFHs of several of the most extreme galaxies imply the bulk of their stellar mass formed prior to z∼10z \sim10 (e.g. de Graaff et al. 2025; Glazebrook et al. 2024; McConachie et al. 2026), while other population-level studies favour formation redshifts between z∼5−8z \sim5 - 8 (e.g. Park et al. 2024; Nanayakkara et al. 2025). Together, these constraints point towards massive, dust-obscured star-forming galaxies bright in the sub-millimetre as the primary candidates for the progenitors of MQGs (Casey et al. 2014; Wilkinson et al. 2017; Hodge & da Cunha 2020; Valentino et al. 2020; Araya-Araya et al. 2025), since they combine the high star formation rates (SFR ≳100\gtrsim100 M⊙_\odot yr−1^{-1}), high stellar masses (M⋆∼1010−1011M_\star\sim10^{10} - 10^{11} M⊙_\odot) and compact morphologies (Re,circ∼1R_{\mathrm{e,circ} \sim1} kpc) needed to build up a MQG within the ∼100−500\sim100 - 500 Myr implied by their SFHs (Toft et al. 2014). More recently, so-called ‘NIRCam-dark’ galaxies – extremely faint or undetected galaxies in deep NIRCam imaging – have been proposed as promising candidates, with redshifts and inferred properties that link them naturally to both MQGs and the JWST-identified ‘Blue Monster’ population at earlier epochs (Ziparo et al. 2023; Bing et al. 2026; Sun et al. 2026). However, these links rely on the inferred properties of MQGs themselves, which remain uncertain.

Our understanding of when and how rapidly MQGs formed their stars, and consequently what their progenitors must have been, rests on the accuracy of the underlying modelling of the spectral energy distribution (SED), which is known to break down for galaxies with the abundance patterns observed in MQGs. For example, standard SED-fitting routines use stellar population libraries that assume solar-scaled abundance patterns which can introduce systematic biases in the derived metallicity when applied to extremely α\alpha-enhanced populations (Beverage et al. 2025). As a result, the derived formation times and ages can vary significantly depending on the assumed abundance pattern (e.g. McConachie et al. 2026), which is compounded due to the degeneracy between stellar age and metallicity (Vazdekis et al. 2015). A diagnostic that constrains formation timescales directly from the observed chemistry, independent of SED-derived SFHs, is therefore needed to test the rapid-formation scenario.

Galactic chemical evolution (GCE) modelling provides such a diagnostic by using the observed abundance ratios as direct constraints on the duration of star formation (for a review see Matteucci 2021). By comparing the relative abundances of α\alpha-elements and Fe to the timescales of their nucleosynthetic progenitors while also accounting for overall metallicity evolution, this approach allows for the formation timescales to be derived independently of the systematics inherent to SED modelling. Beyond the formation timescale, GCE models also constrain the star-formation efficiency (SFE) timescale (τ⋆\tau_\star) and mass-loading factor (η\eta) required to reproduce the observed abundances, parameters that can be measured independently in candidate progenitor populations such as DSFGs. With JWST now producing multi-element abundances for MQGs at z≳2z \gtrsim2, GCE models are well-placed to make use of this data.

In this work, we use the Versatile Integrator for Chemical Evolution (VICE2^2; Johnson & Weinberg 2020) GCE code to explore trends in abundance ratios [X/Fe] as a function of the star-formation timescale (ΔtSF\Delta t_{\mathrm{SF}}) and to investigate the parameters (τ⋆,η\tau_\star,\eta) and formation timescales required to reproduce the abundances of observed MQGs. Our analysis is guided by the following questions:

  • How well can a one-zone GCE model reproduce the abundance patterns of MQGs at cosmic noon?

  • What implications do our results have for the formation timescales and properties of MQGs, including their likely progenitors?

  • How sensitive are our conclusions to systematic uncertainties in the underlying nucleosynthetic yields?

We describe the model framework and global parameter calibration using the nearby dwarf galaxy Sculptor in § 2 and apply it to MQGs in § 3. Results are presented in § 4 and § 5 and we discuss implications of MQG formation in § 6. Throughout this paper we adopt a Kroupa (2001) IMF, solar abundances from Asplund et al. (2009) and a flat Lambda cold dark matter (Λ\LambdaCDM) cosmology with Ωm=0.3\Omega_{\mathrm{m}} = 0.3 and H0=70 km s−1 Mpc−1H_0 = 70\ \mathrm{km\ s^{-1}\ Mpc^{-1}}.

Chemical Evolution Model and Calibration

Our chemical evolution modelling uses the publicly available VICE GCE code which is described in great detail in Johnson & Weinberg (2020). In § 2.1, we provide a short overview of the one-zone implementation adopted in this work and in § 2.2, we describe the determination of empirical yields using the dwarf galaxy Sculptor.

We utilise VICE for our analysis due to its programmatic flexibility. The code allows for arbitrary star-formation histories, both theoretical yield tables and user-defined nucleosynthetic yields, and custom IMFs with self-consistent continuous recycling. We leverage these capabilities to investigate abundance variations in MQGs at cosmic noon.

Chemical Evolution Model Framework

We employ a single-zone model to describe the chemical evolution of a galaxy treating it as a single, well-mixed reservoir that accretes and ejects gas, forms stars, and is enriched by their subsequent return of metals to the ISM.

The evolution of the system is governed by a set of ordinary differential equations, which we integrate numerically using VICE. We choose to define a SFH parametrised as

M˙⋆(t)=A⋅(1−e−t/τ1)⋅e−t/τ2,(1)\dot{M}_{\star}(t) = A \cdot\left(1 - e^{-t/\tau_1}\right) \cdot e^{-t/\tau_2}, \tag*{(1)}

where M˙⋆(t)\dot{M}_{\star}(t) is the SFR, AA is a normalisation factor we fit for a given stellar mass M⋆M_{\star}, and τ1\tau_1 and τ2\tau_2 are free parameters that govern the timescales over which the SFR rises and falls respectively. This is the rise-and-fall parametrisation adopted in Johnson et al. (2021) and used by Gountanis et al. (2025) to model the chemical evolution of elliptical galaxies at redshifts z∼0z \sim0 and z∼0.7z \sim0.7, making it useful for our application to MQGs.

Enrichment of the stellar population is driven by the gas supply, which provides the fuel for star formation and the reservoir from which stars obtain their chemical composition. The gas mass evolves according to the balance between sources (inflows M˙in\dot{M}_{\mathrm{in}}, stellar recycling M˙returned\dot{M}_{\mathrm{returned}}) and sinks (star formation M˙⋆\dot{M}_{\star}, outflows M˙out\dot{M}_{\mathrm{out}}):

M˙g=M˙in−M˙⋆−M˙out+M˙returned.(2)\dot{M}_{g} = \dot{M}_{\mathrm{in}} - \dot{M}_{\star} - \dot{M}_{\mathrm{out}} + \dot{M}_{\mathrm{returned}}. \tag*{(2)}

Since accretion histories cannot be directly derived from observations, it is more useful to describe Eqn. (2) using observable parameters (e.g. Weinberg et al. (2017)), namely the SFE timescale (τ⋆=Mg/M˙⋆\tau_{\star} = M_g/\dot{M}_{\star}) and mass-loading factor (η=M˙out/M˙⋆\eta= \dot{M}_{\mathrm{out}}/\dot{M}_{\star}). Implicitly, our model assumes instantaneous cooling—the gas available for star formation is equal to the gas present in the zone—and hence the SFE timescale effectively describes the molecular gas depletion timescale τ⋆≈τmol\tau_{\star} \approx\tau_{\mathrm{mol}}. This assumption appears reasonable given a typical cooling time of <105 yr<10^5\ \mathrm{yr} for dense molecular clouds (n(H2)∼103−104 cm−3n(\mathrm{H}_2) \sim10^3-10^4\ \mathrm{cm}^{-3}; Goldsmith & Langer (1978)), which is much shorter than the depletion times considered in this work. Mass return from stars is modelled using VICE’s time-dependent continuous recycling prescription (Eqn. (2) of Johnson & Weinberg (2020)) rather than assuming instantaneous. With these substitutions, Eqn. (2) can be solved for the gas accretion rate M˙in\dot{M}_{\mathrm{in}} in terms of observable quantities providing a complete description of the gas supply fuelling star formation.

Chemical enrichment is similarly governed by a balance of metal sources and sinks. Source terms include the IMF-integrated nucleosynthetic yields from core-collapse supernovae (CCSNe), Type Ia supernovae (SNe Ia) and asymptotic giant branch (AGB) stars, along with metals carried by accreted gas and returned via stellar recycling. Conversely, sink terms account for metals removed in outflows and those locked into newly formed stars. All yields in VICE are IMF-integrated, fractional net yields (Johnson & Weinberg (2020)) specifying the net amount of that element gained by the ISM from a single star or stellar population yXΛy_X^{\Lambda}, where Λ\Lambda is the enrichment channel (e.g. CCSNe) and XX is the element. In principle, these yields can be adopted from standard theoretical yield tables, however in practice, these values are highly uncertain and struggle to reproduce observed abundances in the Milky Way (e.g. Andrews et al. (2017)).

The framework above defines two distinct classes of parameter that control chemical enrichment. We refer to these as global parameters (those assumed universal across galactic environments, including the IMF-averaged nucleosynthetic yields yXΛy_X^{\Lambda} and the minimum SNe Ia delay time tdelayt_{\mathrm{delay}}) and local parameters (those that are specific to a given galaxy, including the gas depletion timescale τ⋆\tau_{\star}, the mass-loading factor η\eta and the SFH timescales τ1,τ2\tau_1,\tau_2). The assumption of global yields requires a universal IMF and metallicity independence, whose implications we revisit in § 6.2. In § 2.2 we calibrate the global parameters empirically using the dwarf spheroidal (dSph) galaxy Sculptor, exploiting the fact that yields fixed against one well-characterised system can be applied to others. The local parameters τ⋆\tau_{\star} and η\eta are then modified between galaxies for application to MQGs.

Yield Calibration with Sculptor

Calibrating the global parameters defined in § 2.1 requires a galaxy whose chemistry can be modelled cleanly with the one-zone approximation. Although the Milky Way offers extensive data sets for such calibration, its complex assembly and star formation history make it a poor benchmark for a one-zone model. In contrast, Sculptor is a single, well-mixed galaxy whose stellar population has been well-studied (e.g. Kirby et al. 2009; Hill et al. 2019; Henderson et al. 2025) and is more suited to the assumptions made by single-zone models (e.g. Escala et al. 2018; de los Reyes et al. 2022). Moreover, Sculptor’s SFH (e.g. Kirby et al. 2011; Vincenzo et al. 2014; Weisz et al. 2014; Bettinelli et al. 2019; de los Reyes et al. 2022) appears analogous to those assumed for MQGs, whereby star formation proceeds in a rapid burst at very-early times (e.g. Gountanis et al. 2025).

We fit the CCSN yields of Mg, Ca, Ti and Fe against the observed abundances of the Sculptor red-giant branch (RGB) stars. These four elements span the dominant nucleosynthetic channels relevant to our framework: the α\alpha-elements, Mg, Ca and Ti, are produced predominantly by CCSNe with minor contributions from SNe Ia for Ca and Ti. In contrast, Fe receives substantial contributions from both channels, with SNe Ia dominating its yield in stellar populations old enough for the SN Ia delay-time distribution to have significant effect (Nomoto et al. 2013; Kobayashi et al. 2020). Fitting the three abundance ratios [Mg/Fe], [Ca/Fe] and [Ti/Fe] simultaneously constrains yFeCCy_{\mathrm{Fe}}^{\mathrm{CC}} better than any single ratio could, since the height of the early-time ‘CCSN-plateau’ ([α/Fe]cc[\alpha/\mathrm{Fe}]_{\mathrm{cc}}; the plateau in [α/Fe][\alpha/\mathrm{Fe}] at low [Fe/H] prior to the onset of SNe Ia) is set by the ratio yαCC/yFeCCy_{\alpha}^{\mathrm{CC}}/y_{\mathrm{Fe}}^{\mathrm{CC}} for each α\alpha-element.

The remaining yields enter the model as fixed inputs. AGB yields are strongly mass- and metallicity-dependent and cannot be appropriately described by a single population-averaged yield. For this reason, we select to model elements that are not enriched by AGB stars. For SN Ia, we fix the Fe yield yFeIay_{\mathrm{Fe}}^{\mathrm{Ia}} to the empirical normalisation of the SN Ia delay-time distribution (DTD) rather than to a single theoretical model. Following Weinberg et al. (2017) (and similarly Johnson & Weinberg (2020)), the IMF-integrated, time-averaged yield of Fe from SNe Ia is given by

yFeIa=MFe∫0∞RIa(t) dt(3)y_{\mathrm{Fe}}^{\mathrm{Ia}} = M_{\mathrm{Fe}} \int_{0}^{\infty} R_{\mathrm{Ia}}(t)\,\mathrm{d}t \tag*{(3)}

where MFeM_{\mathrm{Fe}} is the average mass of Fe produced in a single SN Ia event. The integral over the DTD then represents the number of SN Ia events per stellar mass formed, ∫0∞RIa(t) dt=NIa/M⋆\int_{0}^{\infty} R_{\mathrm{Ia}}(t)\,\mathrm{d}t = N_{\mathrm{Ia}}/M_{\star}. Using the constraint NIa/M⋆=(2±1)×10−3 M⊙−1N_{\mathrm{Ia}}/M_{\star} = (2 \pm1) \times10^{-3}\,M_{\odot}^{-1} found by Maoz et al. (2012), one can solve for yFeIay_{\mathrm{Fe}}^{\mathrm{Ia}} given MFeM_{\mathrm{Fe}}. Adopting an average iron mass of MFe=0.78 M⊙M_{\mathrm{Fe}} = 0.78\,M_{\odot} from the W70 explosion model of Iwamoto et al. (1999), we fix the SN Ia iron yield to yFeIa=0.0017y_{\mathrm{Fe}}^{\mathrm{Ia}} = 0.0017, consistent with the value used by Weinberg et al. (2017) and Johnson & Weinberg (2020). For the remaining elements, we take the SN Ia yields directly from the yield table produced in Iwamoto et al. (1999), with the exception of Mg, for which we assume negligible SN Ia enrichment setting yMgIa=0y_{\mathrm{Mg}}^{\mathrm{Ia}} = 0.

The DTD itself is taken to be a power law of the form RIa∝t−1.1R_{\mathrm{Ia}} \propto t^{-1.1}, matching the slope inferred by Maoz et al. (2012). Related to the DTD, the minimum delay time tdelayt_{\mathrm{delay}} then sets the interval between the formation of the first stars and the onset of Type Ia supernovae. The delay time is inherently sensitive to the lifetimes of SN Ia progenitors which remain a subject of debate. However, both observational and theoretical constraints suggest prompt SNe Ia are plausible as early as ∼50\sim50 Myr (e.g. (Mannucci et al. 2006; Aubourg et al. 2008; Totani et al. 2008; Mennekens et al. 2010; Wang & Han 2012; Maoz et al. 2014; Castrillo et al. 2021; Chen et al. 2021)). To this end, we fix the minimum delay time to tdelay=50t_{\mathrm{delay}} = 50 Myr, motivated by literature findings and tests that show our posterior chains run up the lower wall when a lower bound on the prior is set. Like the yields, tdelayt_{\mathrm{delay}} is assumed to be a global parameter and we apply the same delay time to both Sculptor and MQGs.

The remaining parameters governing the enrichment of the ISM are the SFE timescale τ⋆\tau_{\star} and the mass-loading factor η\eta, which are local properties. Although their best-fit values apply only to Sculptor and are not applied in the MQG analysis, both must be fit simultaneously with the global parameters since they regulate the gas reservoir. We adopt uniform priors for both, in log-space for τ⋆\tau_{\star} and linear-space for η\eta, with bounds given in Table 1.

ParameterDescriptionPriorDefault NKT13 Net YieldsBest-fitScaling
LocalAA Normalisation of the SFH (M⊙ yr−1M_\odot\,\mathrm{yr}^{-1})log⁡10U(10−4,10−1)\log_{10}\mathcal{U}(10^{-4},10^{-1})—7.17−1.32+1.45×10−37.17^{+1.45}_{-1.32}\times10^{-3}—
Localτ⋆\tau_\star Star formation efficiency timescale (Gyr)log⁡10U(1,20)\log_{10}\mathcal{U}(1,20)—1.42−0.24+0.291.42^{+0.29}_{-0.24}—
Localη\eta Mass loading factor (outflow efficiency)U(5,60)\mathcal{U}(5,60)—38.62−2.33+2.4538.62^{+2.45}_{-2.33}—
Localτ1\tau_1 SFH rise timescale (Gyr)log⁡10U(0.1,1)\log_{10}\mathcal{U}(0.1,1)—0.02−0.01+0.020.02^{+0.02}_{-0.01}—
Localτ2\tau_2 SFH fall timescale (Gyr)log⁡10U(0.1,1)\log_{10}\mathcal{U}(0.1,1)—0.37−0.05+0.050.37^{+0.05}_{-0.05}—
GlobalyFeCCy_{\mathrm{Fe}}^{\mathrm{CC}} Net Fe yield from CCSNelog⁡10U(0.1 ydef,10 ydef)\log_{10}\mathcal{U}(0.1\,y_{\mathrm{def}},10\,y_{\mathrm{def}})5.56×10−45.56\times10^{-4}6.64−0.61+0.62×10−46.64^{+0.62}_{-0.61}\times10^{-4}1.19
GlobalyMgCCy_{\mathrm{Mg}}^{\mathrm{CC}} Net Mg yield from CCSNelog⁡10U(0.1 ydef,10 ydef)\log_{10}\mathcal{U}(0.1\,y_{\mathrm{def}},10\,y_{\mathrm{def}})1.51×10−31.51\times10^{-3}1.33−0.09+0.09×10−31.33^{+0.09}_{-0.09}\times10^{-3}0.88
GlobalyCaCCy_{\mathrm{Ca}}^{\mathrm{CC}} Net Ca yield from CCSNelog⁡10U(0.1 ydef,10 ydef)\log_{10}\mathcal{U}(0.1\,y_{\mathrm{def}},10\,y_{\mathrm{def}})7.67×10−57.67\times10^{-5}6.94−0.57+0.60×10−56.94^{+0.60}_{-0.57}\times10^{-5}0.90
GlobalyTiCCy_{\mathrm{Ti}}^{\mathrm{CC}} Net Ti yield from CCSNelog⁡10U(0.1 ydef,10 ydef)\log_{10}\mathcal{U}(0.1\,y_{\mathrm{def}},10\,y_{\mathrm{def}})1.30×10−61.30\times10^{-6}2.83−0.28+0.30×10−62.83^{+0.30}_{-0.28}\times10^{-6}2.18

Table 1. Best-fit parameters for the Sculptor dSph derived from MCMC modelling. Priors are uniform in log-space for all parameters except η\eta, with bounds given in linear units. Tabulated CCSN yields are taken from the Nomoto et al. (2013, NKT13) grid at log⁡(Z/Z⊙)=−0.24\log(Z/Z_\odot)=-0.24 and provided as a comparison to the fitted yields, with the scaling shown in the last column. Best-fit values indicate the median of the posterior distributions with uncertainties representing the 16th and 84th percentiles.

Since our model’s gas component physically corresponds to the gas available for star formation, the prior on τ⋆\tau_{\star} is motivated by analyses of H2_2 (molecular gas) depletion times (τmol\tau_{\mathrm{mol}}) in nearby star-forming galaxies and in low-metallicity dwarf galaxies, which represent the closest present-day analogues to dSph progenitors at their formation epochs. In MW-like spiral galaxies, the molecular gas depletion time is approximately τmol∼2\tau_{\mathrm{mol}} \sim2 Gyr (Bigiel et al. 2008) whereas in star-forming dwarfs, τmol\tau_{\mathrm{mol}} is more variable and on average shorter, spanning ∼0.1−6\sim0.1 - 6 Gyr (Bolatto et al. 2011; Grossi et al. 2016) and reaching as low as ∼100\sim100 Myr in starbursting galaxies (Hunt et al. 2015; Amorín et al. 2016). In contrast, depletion times for HI in dwarfs are significantly longer, τHI∼10−100\tau_{\mathrm{HI}} \sim10 - 100 Gyr (e.g. (Bigiel et al. 2010; Roychowdhury et al. 2014)). Since we do not explicitly distinguish H2_2 and HI in our models, our prior spans the range of observed τmol,τHI\tau_{\mathrm{mol}}, \tau_{\mathrm{HI}} values allowing the data to constrain the physical representation of the gas reservoir.

The prior on η\eta is motivated by GCE analyses of dSphs, which require mass-loading factors higher than in typical star-forming galaxies in order to support the efficient removal of gas necessary to reproduce the low-metallicity populations hosted by these galaxies (Lanfranchi & Matteucci 2004; Vincenzo et al. 2014). While direct observations of outflows in dwarfs are limited, available measurements for M⋆∼106−7 M⊙M_{\star} \sim10^{6-7}\,M_{\odot} galaxies suggest η\eta is of the order ∼101−102\sim10^{1} - 10^{2} (e.g. [158]). Similar predictions for η\eta are made using the observed mass-metallicity relation (Peeples & Shankar 2011; Lu et al. 2015), where η\eta is found to increase for lower-mass galaxies. We therefore adopt a prior range of 5 – 60 for Sculptor.

For our observational reference, we adopt elemental abundances from the medium- and high-resolution spectroscopic catalogues of RGB stars from Henderson et al. (2025) and Hill et al. (2019), observed using Keck/DEIMOS and VLT/FLAMES respectively. Both surveys report [Fe/H] and [Mg/Fe] and we additionally include [Ca/Fe] and [Ti/Fe] from Hill et al. (2019) to broaden the CCSN-yield constraints, since these elements are not reported in the Henderson et al. (2025) catalogue. To ensure data quality, we restrict the sample to stars with measurement uncertainties σ<0.3\sigma< 0.3 dex across all considered abundance ratios and exclude outliers with [Mg/Fe] >0.75> 0.75. The latter cut is a deliberately liberal upper bound given that the Milky Way CCSN plateau is generally taken to be [α/Fe]cc∼0.45[\alpha/\mathrm{Fe}]_{\mathrm{cc}} \sim0.45, however, this is chosen to remove likely outliers without aggressively cutting high-α\alpha stars. We discuss the implications of this in § 6.2.

We infer the posterior distributions for the 9 free model parameters (τ⋆\tau_{\star}, η\eta, τ1\tau_1, τ2\tau_2, yFeCCy_{\mathrm{Fe}}^{\mathrm{CC}}, yMgCCy_{\mathrm{Mg}}^{\mathrm{CC}}, yCaCCy_{\mathrm{Ca}}^{\mathrm{CC}}, yTiCCy_{\mathrm{Ti}}^{\mathrm{CC}}) using an MCMC sampling routine, following Kirby et al. (2011) and de los Reyes et al. (2022). We adopt the medians of the resulting posteriors as our parameter estimates, with uncertainties derived from the 16th and 84th percentiles. Each observed star ii is treated as a probability density in a multi-dimensional chemical space defined by ϵj,i\epsilon_{j,i}, where jj denotes the abundance ratio (e.g. [Mg/Fe]) and σj,i\sigma_{j,i} is the associated measurement uncertainty. The likelihood LiL_i that a single star originates from a specific model is calculated by integrating the star formation rate along the predicted chemical evolution track ϵj(t)\epsilon_j(t) from VICE such that

Li=∫0tend∏j12πσi,jexp⁡(−(ϵi,j−ϵj(t))22(σi,j)2)M˙⋆(t)M⋆,final dt(4)L_i = \int_{0}^{t_{\mathrm{end}}} \prod_j \frac{1}{\sqrt{2\pi}\sigma_{i,j}} \exp\left(\frac{-(\epsilon_{i,j} - \epsilon_j(t))^2}{2(\sigma_{i,j})^2}\right) \frac{\dot{M}_{\star}(t)}{M_{\star,\mathrm{final}}}\,\mathrm{d}t \tag*{(4)}

This effectively weights the probability that a star forms at time tt, given by the mass of stars formed in each time step M˙⋆/M⋆,final\dot{M}_{\star}/M_{\star,\mathrm{final}} with an abundance ratio ϵj(t)\epsilon_j(t). The maximum runtime of the model set by tendt_{\mathrm{end}} is nominally 3 Gyr, motivated by other models which suggest Sculptor forms in a rapid burst (e.g. Vincenzo et al. 2014). However in practice, this acts as an upper limit and the integration terminates earlier, at the time step when the cumulative formed stellar mass reaches 99% to avoid artefacts from near-zero gas mass.

The total likelihood LL is then given by the product of the individual likelihoods for NN stars, combined with Gaussian penalty terms Notes.

included for the total stellar (M⋆,obsM_{\star,\mathrm{obs}}) and gas masses (Mg,obsM_{\mathrm{g,obs}}):

L=∏i=1N[Li×12πδM⋆,obsexp⁡(−(M⋆,obs−M⋆,VICE)22(δM⋆,obs)2)×12πδMg,obsexp⁡(−Mg,VICE22(δMg,obs)2)].(5)L=\prod_{i=1}^{N}\left[L_i\times\frac{1}{\sqrt{2\pi}\delta M_{\star,\mathrm{obs}}}\exp\left(-\frac{(M_{\star,\mathrm{obs}}-M_{\star,\mathrm{VICE}})^2}{2(\delta M_{\star,\mathrm{obs}})^2}\right)\times\frac{1}{\sqrt{2\pi}\delta M_{\mathrm{g,obs}}}\exp\left(-\frac{M_{\mathrm{g,VICE}}^2}{2(\delta M_{\mathrm{g,obs}})^2}\right)\right]. \tag*{(5)}

Unlike Kirby et al. (2011) and de los Reyes et al. (2022), we do not include an exponent weighting for the mass penalties since we observe the likelihood contributions from stellar abundances and these penalties to be on the same order of magnitude. We adopt the stellar mass from McConnachie (2012) as the observational reference and assume the gas mass at the end of formation is zero. Additionally, we assume arbitrary observational uncertainties of δM⋆,obs=0.1M⋆\delta M_{\star,\mathrm{obs}}=0.1M_\star and δMg,obs=0.01M⋆\delta M_{\mathrm{g,obs}}=0.01M_\star for the stellar and gas masses respectively.

We explore the parameter space using emcee (Foreman-Mackey et al. (2013)) with an ensemble of 64 walkers evolved for 15000 steps. To accelerate convergence and avoid chains getting stuck in local maxima, we employ a resampling strategy following the ‘burn-in’ phase. After 5000 iterations, walkers falling in the bottom 25th percentile of likelihoods are re-initialised to the position of the highest-likelihood walker. The walkers are then evolved for a further 1000 steps to allow those that have been re-initialised to diffuse in the posterior space. We discard these first 6000 steps as burn-in and use the subsequent 9000 steps for posterior sampling. Convergence is verified by monitoring the evolution of autocorrelation time across the production chains. We adopt log-uniform priors on all free parameters, except η\eta, with bounds motivated by previous chemical evolution models of Sculptor (e.g. τ⋆\tau_\star) or theoretical/observational constraints (e.g. η\eta). For the four CCSN IMF-averaged yields, the bounds span ±1\pm1 dex around the Nomoto et al. (2013, NKT13) theoretical values for each element. A summary of the priors, the default VICE settings and our resulting best-fit values are provided in Table 1.

Calibration Results and Validation

In Figure 1 we show the best-fit VICE model compared to the observational data, in terms of the abundance and metallicity distributions (MDFs), as well as the abundance trends themselves. The model successfully reproduces the observed abundance trends across the three fitted α\alpha-elements with the chemical evolution tracks capturing the ‘knee’ in each of the [α\alpha/Fe]-[Fe/H] planes, and predicted MDFs that closely match the observed distributions. Together, these confirm that the fit has effectively constrained the interplay between the CCSN yields, τ⋆\tau_\star, and η\eta. This strong agreement provides confidence that our empirically-derived yields are physically reasonable. As detailed in Table 1, the best-fit CCSN yields for Mg and Ca are scaled down by a factor of ∼0.9\sim0.9 relative to the NKT13 baseline, while yFeCCy_{\mathrm{Fe}}^{\mathrm{CC}} is scaled up by ∼1.2\sim1.2. Conversely, the Ti yield is scaled up by a factor of ∼2.2\sim2.2, which is expected given theoretical yields tend to underproduce Ti (e.g. Timmes et al. 1995; Romano et al. 2010; Kobayashi et al. 2020).

Best-fit chemical evolution model for the Sculptor dSph

Figure 1. Best-fit chemical evolution model for the Sculptor dSph. The main panels display the observed abundances for three datasets: Keck-II/DEIMOS medium-resolution data from Henderson et al. (2025) (blue), VLT/FLAMES high-resolution data from Hill et al. (2019) (orange) and Keck-II/DEIMOS medium-resolution data from Kirby et al. (2010) (aqua). Overlaid in red are the evolutionary tracks from the best-fit VICE model. The top and righthand marginal panels show the corresponding MDFs for each dataset in their respective colours, alongside the predicted MDF from VICE. The pink histograms (labelled ‘Total’) in the [Fe/H] and [Mg/Fe] panels represent the combined abundance sample used to constrain the model for those particular abundances (i.e. Henderson et al. (2025) and Hill et al. (2019)).

For application to the multi-element abundances reported for MQGs at cosmic noon, we also extend our calibrated yield set to Si, another α\alpha-element with contributions from both CCSNe and SNe Ia. We note ySiCCy_{\mathrm{Si}}^{\mathrm{CC}} by visual inspection to match observed abundances in the [Si/Fe]-[Fe/H] space from Kirby et al. (2010) and Hill et al. (2019). Because the primary [Fe/H] measurements used in our likelihood come from a combined sample of observations from Hill et al. (2019) and Henderson et al. (2025), we treat [Si/Fe] as an external constraint to avoid mixing additional samples in a single MCMC routine. We find ySiCC=0.5⋅ySi,NKT13CC=7.3×10−4y_{\mathrm{Si}}^{\mathrm{CC}}=0.5\cdot y_{\mathrm{Si,NKT13}}^{\mathrm{CC}}=7.3\times10^{-4} best reproduces the observed [Si/Fe] as shown in Figure 1. This corresponds to a ∼0.26\sim0.26 dex reduction in the mass-weighted mean stellar [Si/Fe] when compared to the same model using the NKT13 yield, which is consistent with the ∼0.2\sim0.2 dex overproduction of [Si/Fe] reported in models applying theoretical yields to MW stars (e.g. Kobayashi et al. 2011, 2020).

In Figure 2 we compare our derived cumulative SFH to those from the literature, which include from other GCE models (Kirby et al. 2011; Vincenzo et al. 2016; de los Reyes et al. 2022) and those derived from photometric observations ([164] (Bettinelli et al. 2019)). We note that GCE models are agnostic to cosmic time and hence the start time of our model is not constrained. For a meaningful physical comparison with other models, we instead align our model to the formation epoch assumed by de los Reyes et al. (2022). Generally SFHs derived from observations predict longer formation timescales for Sculptor, however as noted by de los Reyes et al. (2022), the photometry used for these measurements tend to not be deep enough to observe stars below the main-sequence turn-off. When observations sufficiently probe stars below the turn-off, they find that Sculptor formed the majority of its stars in a short burst (∼2.2\sim2.2 Gyr at FWHM; Bettinelli et al. (2019)) with the majority forming over 10 Gyr ago (Weisz et al. 2014).

Cumulative SFH as a function of lookback time for the Sculptor dSph

Figure 2. Cumulative SFH as a function of lookback time for the Sculptor dSph. The predicted SFH from VICE is shown in red, with the shaded region representing the 1σ1\sigma uncertainty derived from the SFH parameters τ1=0.02−0.01+0.02\tau_1 = 0.02^{+0.02}_{-0.01} Gyr and τ2=0.37−0.05+0.05\tau_2 = 0.37^{+0.05}_{-0.05} Gyr. SFHs from chemical evolution models are shown in light grey: Kirby et al. (2011) (dotted), Vincenzo et al. (2016) (dashed), and de los Reyes et al. (2022) (solid). SFHs from observational studies are shown in dark grey: de Boer et al. [164] (dotted) and Bettinelli et al. (2019) (solid).

The formation timescale we recover is consistent with previous GCE models of Sculptor. We find 99% of the stellar mass formed in 1.72 Gyr (τ1=0.02−0.01+0.02\tau_1 = 0.02^{+0.02}_{-0.01} Gyr, τ2=0.37−0.05+0.05\tau_2 = 0.37^{+0.05}_{-0.05} Gyr) after the onset of star formation, in rough agreement with the 2.16 Gyr reported by Vincenzo et al. (2016) using the same timescale. Other GCE models report shorter bursts of ∼1\sim1 Gyr (Kirby et al. 2011; de los Reyes et al. 2022) which may be attributed to a combination of methodological differences: our model uses a SFH parametrisation rather than a prescribed gas-inflow history; we adopt Henderson et al. (2025) as our primary abundance sample, which extends to lower metallicities than Hill et al. (2019) and reports systematically higher [Mg/Fe] and lower [Fe/H] than Kirby et al. (2010); and our outflow prescription, while broadly consistent with the wind models presented in Kirby et al. (2011) and subsequently de los Reyes et al. (2022), closely follows the form in Vincenzo et al. (2014). We test whether adopting observational data from Kirby et al. (2010) has a significant impact on the resulting SFH and parameters but find no significant difference. Otherwise, we do not attempt to isolate which methodological difference contributes the most to the timescale offset, but identify them as the principal differences between the models.

In addition to the abundances and SFH parameters, we also derive estimates of the SFE timescale and mass-loading. While τ⋆\tau_\star does not explicitly define a depletion timescale in terms of H2_2 or HI, our bestfit value (τ⋆=1.42−0.24+0.29\tau_\star= 1.42^{+0.29}_{-0.24} Gyr) is physically consistent with observations of τmol\tau_{\rm mol} in local star-forming dwarfs (τmol∼0.1−6\tau_{\rm mol} \sim0.1 - 6 Gyr; Grossi et al. 2016) and is at least an order of magnitude shorter than typical HI depletion times (e.g. Bigiel et al. 2010; Roychowdhury et al. 2014). Given this agreement, we interpret τ⋆\tau_\star to physically represent the depletion timescale of star-forming molecular gas and adopt this interpretation in our application to MQGs. This interpretation of τ⋆\tau_\star supports the assessment that MgM_{\rm g} in our models only represents the star-forming gas, as expected. Furthermore, the mass-loading factor in our models is well-constrained (η=38.62−2.33+2.45\eta= 38.62^{+2.45}_{-2.33}) and in agreement with literature expectations (η∼101−102\eta\sim10^1 - 10^2; [158] (Lu et al. 2015)).

In summary, while our primary focus is not to constrain Sculptor’s formation and assembly but to use it as a tool to empirically calibrate the global yields that will be applied to MQGs in § 3, the results of our fitting support the validity of our calibration. The inferred SF timescale agrees with independent constraints derived from deep CMD-based studies of Sculptor’s stellar population (e.g. Bettinelli et al. 2019). Following a scaling of ySiCCy_{\rm Si}^{\rm CC} to match the observed plateau, the metallicity evolution of our model predicts the observed [Si/Fe]-[Fe/H] trend for an element not included in the fit, indicating that the recovery of model parameters beyond the yields (τ⋆\tau_\star, η\eta) is successful. Lastly, our derived τ⋆\tau_\star and η\eta values are physically consistent with observations of τmol\tau_{\rm mol} and mass-loading in low-mass galaxies. Together, these results justify adopting the calibrated CCSN yields (including the scaled Si yield), SN Ia yields and minimum delay time as the global parameters in our MQG analysis.

Application to Massive Quiescent Galaxies

Having calibrated empirically-motivated yields using the Sculptor dSph, we now consider the main aim of this work: modelling the chemical evolution of massive quiescent galaxies in the early universe. While we retain the global parameters (yields and SNe Ia delay time tdelay=50t_{\rm delay} = 50 Myr) from fitting Sculptor, the local properties governing ISM enrichment (SFE timescale τ⋆\tau_\star, mass-loading factor η\eta and SFHs) take values appropriate to MQGs and must be reset accordingly. Here, we establish our observational reference sample of MQGs in § 3.1 and then describe our treatment of the local parameters in § 3.2.

Observed Sample of MQGs

To place our models in context, we compare our predicted abundances to a compilation of observational abundance measurements from the literature. Until the launch of JWST, only a handful of [α/Fe] abundances were available for MQGs at z>1z > 1, with observations drawn from either lensing studies or deep near-infrared observations (e.g. Kriek et al. 2016; Jafariyazani et al. 2020; Carnall et al. 2024; Jafariyazani et al. 2025; McConachie et al. 2026). Recent observations with the Heavy Metal and JWST-SUSPENSE surveys have greatly expanded upon this sample providing abundance measurements for several MQGs, with multi-element abundances available from JWST-SUSPENSE in particular (Kriek et al. 2024; Beverage et al. 2024, 2025). In § 5 we primarily compare our model to these surveys and therefore summarise their selection criteria below.

The Heavy Metal survey (Kriek et al. 2024) is an ultra-deep, rest-frame optical spectroscopic survey that targeted MQGs in the early universe across two redshift ranges (1.30<z<1.501.30 < z < 1.50, 1.92<z<2.281.92 < z < 2.28) in the common area covering the UltraVISTA (McCracken et al. 2012), COSMOS (Scoville et al. 2007) and COSMOS-DASH (Momcheva et al. 2017; Mowla et al. 2019) surveys. Targets in this sample were first selected as quiescent based on their UVJ colours from the UltraVISTA survey (Muzzin et al. 2013) and observed using Keck/LRIS and Keck/MOSFIRE. The second selection criterion requires the 1.30<z<1.501.30 < z < 1.50 MQGs to be brighter than J=21.6J = 21.6 and the 1.92<z<2.281.92 < z < 2.28 MQGs to be brighter than H=21.8H = 21.8 to meet S/N requirements (Kriek et al. 2024). With this selection, the Heavy Metal survey consists of 21 massive (log⁡M/M⊙∼10.60−11.35\log M/M_\odot\sim10.60 - 11.35) quiescent galaxies across 1.3<z<2.11.3 < z < 2.1, with chemical abundances available for 19 of these (Kriek et al. 2024; Beverage et al. 2024).

Complementing this, we also draw data from the JWST-SUSPENSE program (Slob et al. 2024; Beverage et al. 2025). SUSPENSE is an ultra-deep, rest-frame optical spectroscopic survey using NIRSpec/MSA on JWST. Targets were identified from the same overlapping field as Heavy Metal, with similar selection criteria applied; UVJ selection for quiescence and a minimum brightness limit of H<22.6H < 22.6 (Slob et al. 2024). In total, 20 quiescent galaxy candidates were observed with multi-element abundances available for 15 massive (log⁡M/M⊙∼10.2−11.2\log M/M_\odot\sim10.2 - 11.2) quiescent galaxies at z∼1−3z \sim1 - 3 (Slob et al. 2024; Beverage et al. 2025). We combine these two samples with other available abundance measurements from single galaxy studies and lensing work to form a z∼2z \sim2 reference sample of ∼35\sim35 galaxies with available [Mg/Fe], [Fe/H] measurements (Kriek et al. 2016; Jafariyazani et al. 2020; Carnall et al. 2024; Jafariyazani et al. 2025; McConachie et al. 2026). We also adopt abundance ratios for other elements beyond [α/Fe] from the analysis of SUSPENSE galaxies by Beverage et al. 2025.

While we focus our modelling efforts on the early MQGs, the recent tension in their observed abundances arises due to the discrepancy between their [α/Fe]-enhanced stellar population and the abundances of their assumed descendants (under the two-phase formation scenario) in the local universe (e.g. Oser et al. 2010; Conroy et al. 2014; Jafariyazani et al. 2025; Beverage et al. 2025). Given these galaxies at cosmic noon are expected to evolve passively, the chemistry of their stellar population should not change significantly over cosmic time. To quantify this discrepancy and place these galaxies in a broader context, we also compare our modelling results to the local ETG population. We benchmark the z∼0z \sim0 population using two catalogues of local ellipticals drawn from the SAURON and MASSIVE galaxy surveys ([167, 170] (Bacon et al. 2001; Kuntschner et al. 2010; Ma et al. 2014)). The SAURON galaxies are an extensively studied sample selected to be a representative sample of E, S0 and Sa galaxies in the local universe (details within [167]). The stellar population analysis of 48 early-type SAURON galaxies by Kuntschner et al. (2010) covers a range of ages and metallicities and is used by [170] who observe 34 of these galaxies with the Keck/LRIS instrument. Specifically, we adopt the catalogue of ETGs from [170] as they preference older galaxies with metallicities 0<Z<0.200 < Z < 0.20 and α-enhancement 0.15<[α/Fe]<0.250.15 < [\alpha/\mathrm{Fe}] < 0.25. While these galaxies are partially selected for their chemistry, their abundance patterns reflect the spread of the broader population of nearby ETGs when compared to stacks of high-quality SDSS data (Conroy et al. 2014). Additionally, we draw a second local reference sample (z<0.025z < 0.025) from the MASSIVE survey (selection details within Ma et al. 2014); a sample of ∼100\sim100 of the most massive ETGs (M⋆≥1011.5M⊙M_\star\ge10^{11.5} M_\odot) within 108 Mpc. Chemical abundances are available from Gu et al. (2022) who observe 41 of the MASSIVE galaxies using deep spectroscopy on the Magellan/LDSS3. These 41 galaxies were selected based only on their declination, greatly improving the S/N and wavelength coverage on previously available spectroscopy (Gu et al. 2022). sample only the inner regions of these galaxies (within <Re/8< R_{\mathrm{e}}/8), providing an excellent comparison to the cosmic noon sample of MQGs, as the centres of local ETGs are their assumed descendants according to the two-phase evolution model (e.g. Oser et al. 2010).

The aforementioned cosmic noon sample of MQGs and local sample of ETGs constitute our compilation of literature abundances. All works measure abundance ratios using full-spectrum fitting of rest-frame optical stellar absorption features. For specific details on the individual samples and fitting, we refer the reader to the original studies.

MQG Model Framework

While the extended, inefficient star formation typical of dSphs is generally not applicable to these massive, rapidly-forming quiescent galaxies, Sculptor’s assembly serves as useful qualitative analogue. The formation histories of both Sculptor and MQGs can be characterised by early, relatively rapid stellar assembly followed by rapid quenching. For this reason, we adopt the same continuous rise-and-fall parametrisation for the SFH (e.g. Johnson et al. 2021; Gountanis et al. 2025). However, rather than adopting fixed values for the rise or fall time, we vary these parameters as needed to investigate the connection between SF timescale and the expected MQG abundances.

More specifically, SFHs with shorter durations are often invoked to explain the deficiency in Type Ia products, both locally and in the distant universe (e.g. Worthey et al. 1992; Conroy et al. 2014; Beverage et al. 2025). To this end, we use a forward-modelling approach to investigate the SF timescales necessary to reproduce the observed abundances of MQGs at cosmic noon.

We explore stellar chemical abundance ratios for a fiducial MQG (M⋆=1011M⊙M_{\star} = 10^{11} M_{\odot}; formed mass) as a function of 6 different star-formation timescales: 0.1, 0.2, 0.4, 0.6, 0.8 and 1.0 Gyr. We parametrise this timescale by ΔtSF=t90−t10\Delta t_{\mathrm{SF}} = t_{90} - t_{10}, defined as the time interval between the formation of 10% (t10t_{10}) and 90% (t90t_{90}) of the galaxy’s formed stellar mass. In this way, ΔtSF\Delta t_{\mathrm{SF}} effectively describes the width of star formation and hence the available time for SN Ia products to be incorporated into newly formed stars before quenching. Since the rise-and-fall model (Eqn. (1)) is not directly defined by ΔtSF\Delta t_{\mathrm{SF}}, we must map this SF duration to the parameters τ1\tau_{1} (rise timescale) and τ2\tau_{2} (decay timescale). This mapping is degenerate: there exists a range of τ1,τ2\tau_{1}, \tau_{2} parameter combinations that could satisfy a fixed ΔtSF\Delta t_{\mathrm{SF}}. To resolve this, we perform a sensitivity analysis by sampling across a grid of τ1\tau_{1} values for each fixed formation duration and solve for the corresponding τ2\tau_{2}. We find that variations in the shape of the underlying SFH produce minimal variation in the [Mg/Fe]-[Fe/H] plane compared to changes in the total duration of star formation. Given this insensitivity, we choose τ1=1/2⋅ΔtSF\tau_{1} = 1/2 \cdot\Delta t_{\mathrm{SF}} as a reasonable approximation and solve for the corresponding τ2\tau_{2}. The resulting SFHs for each of the SF timescales are shown in Figure 3.

Adopted star-formation histories for fiducial massive quiescent galaxy models

Figure 3. Adopted SFHs for our fiducial MQG models. The tracks are parametrised by the star-formation timescale ΔtSF=t90−t10\Delta t_{\mathrm{SF}} = t_{90} - t_{10}. The colour gradient from yellow to maroon corresponds to progressively longer durations of SF, ranging from 0.1–1.0 Gyr.

We must also choose appropriate values for the local parameters: the SFE timescale and mass-loading factor. For Sculptor, we found that τ⋆=1.42−0.24+0.29\tau_{\star} = 1.42^{+0.29}_{-0.24} and η=38.62−2.33+2.45\eta= 38.62^{+2.45}_{-2.33} best described the observed abundances, consistent with observations of star-forming dwarfs (§ 2.3). As established in our modelling, τ⋆\tau_{\star} physically corresponds to the depletion timescale for star-forming gas, i.e. molecular gas. For typical star-forming spiral galaxies, the molecular gas depletion timescale is observed to be ∼1\sim1–22 Gyr (e.g. Bigiel et al. 2008; Genzel et al. 2010) with scaling relations finding shorter depletion times at higher redshift and at greater offsets above the star-forming main sequence (SFMS; e.g. Saintonge et al. 2013; Tacconi et al. 2013, 2018; Saintonge & Catinella 2022). However, the progenitors of MQGs are likely not typical star-forming galaxies, but rather experience significant bursts of star formation (e.g. ΔtSF∼0.1\Delta t_{\mathrm{SF}} \sim0.1–1.01.0 Gyr) driving them further above the SFMS. In cases of starburst galaxies – for example local ultra-luminous infrared galaxies (ULIRGs), distant submillimetre galaxies and simulated starbursting galaxies – molecular gas depletion timescales of ∼10\sim10–100100 Myr have been identified suggesting extremely rapid consumption of gas (e.g. Greve et al. 2005; Pereira-Santaella et al. 2021; Cenci et al. 2024). Our preliminary grid search indicates that to reproduce the elevated abundances of cosmic noon MQGs, a depletion timescale as short as τ⋆=0.05\tau_{\star} = 0.05 Gyr is required, and is in good agreement with the observed depletion times in starbursts. We explore the sensitivity to this parameter further in § 4.2 but adopt τ⋆=0.05\tau_{\star} = 0.05 Gyr as the fiducial value for our MQG models.

We adopt a similar grid-search approach for η\eta. Our analysis favours models with negligible mass loading, consistent with the deep gravitational potential wells of massive galaxies, which inhibit gas loss from the halo. The expected decrease of η\eta with increasing M⋆M_{\star} is supported by inferences from the stellar mass-metallicity relation (Lu et al. 2015). We therefore adopt a minimal mass-loading factor of η=0.01\eta= 0.01 as a fiducial choice, effectively treating these systems as ‘closed’ boxes during their formation. VICE also permits the user to specify an outflow ‘enhancement’ factor ξ\xi, which defines the ratio of the metallicity in the outflow to the ISM. This enhancement factor and the mass-loading factor have a degenerate effect on the stellar metallicity and hence we assume the outflow metallicity is the same as the ISM metallicity (ξ=1\xi= 1) at any given time step.

Finally, we compute our models across a finite range in cosmic time relevant to the high-zz universe. Nominally, the models are allowed a maximum runtime of 3.2 Gyr, corresponding to the age of the universe at z∼2z \sim2. However to ensure numerical stability, we impose a termination condition where the simulation ends once 97.5% of the stellar mass has formed. This prevents the ISM abundances from diverging to unphysically high values as the gas reservoir is exhausted. We use a Kroupa (2001) IMF sampled across the range 0.08<M⋆/M⊙<1000.08 < M_{\star}/M_{\odot} < 100, consistent with the Sculptor calibration in Section 2.2. We also assume a post-main sequence lifetime equivalent to 10% of the main sequence lifetime. Lastly, given the short star-formation timescales modelled, the assumption that gas is returned to the ISM instantaneously is likely invalid. Therefore, our models use time-dependent recycling, explicitly tracking the delayed return of gas mass from dying stars.

The Impact of Model Parameters on Stellar Abundances

In this section, we investigate how the predicted stellar abundances depend on model parameters, with a particular focus on the star-formation timescale ΔtSF\Delta t_{\mathrm{SF}}. We present our results when varying ΔtSF\Delta t_{\mathrm{SF}} in § 4.1, and when varying the two local parameters controlling ISM enrichment, the SFE timescale (τ⋆\tau_\star) and the mass-loading factor (η\eta), in § 4.2. For completeness, we also explore sensitivity to assumptions about the minimum SNe Ia delay time tdelayt_{\mathrm{delay}}.

Star-Formation Timescales

Variations in ΔtSF\Delta t_{\mathrm{SF}} shape the chemical evolution of MQGs through the relative timing of CCSN and SN Ia enrichment. CCSN operate almost instantaneously, while SN Ia contributions are delayed by tdelayt_{\mathrm{delay}}, so the SN Ia contribution to the stellar population depends on how much star formation continues past their onset.

In Figure 4, we show the fraction of Fe produced by SNe Ia as a function of time, in both the ISM (dashed) and stellar population (solid), for varying ΔtSF\Delta t_{\mathrm{SF}}. In all models, the early build-up of stellar mass produces a population dominated by CCSN enrichment. By varying the amount of SF occurring after tdelayt_{\mathrm{delay}} (i.e. increasing ΔtSF\Delta t_{\mathrm{SF}}), there is more available time for SN Ia products to be incorporated into newly formed stars resulting in greater contributions from SNe Ia. While the ISM abundances do receive greater contributions from SNe Ia, the amount that can be incorporated into the stellar population is limited by the SFR at any given time. This is particularly true for

Fraction of Fe produced by SNe Ia in the ISM and stellar population

Figure 4. Fraction of Fe produced by SNe Ia in the ISM (dashed) and in the stellar population (solid). The models are coloured by their associated SF timescale ΔtSF\Delta t_{\mathrm{SF}} as indicated by the colour bar.

short SF timescales, where the SFR is already beginning to decrease by the time SN Ia onset occurs. The resulting CCSN-dominated stellar enrichment is the primary consequence of varying ΔtSF\Delta t_{\mathrm{SF}} and is the likely contributor to the observed enhancement of α\alpha-elements (produced primarily by CCSNe) in MQG spectra.

The effects of varying ΔtSF\Delta t_{\mathrm{SF}} can clearly be observed in the abundances of Mg and Fe. Figure 5 shows the evolution of the ISM abundances ([Mg/Fe], [Mg/H] and [Fe/H]) as a function of time, alongside the resulting mass-weighted stellar abundances for each ΔtSF\Delta t_{\mathrm{SF}}. The reported stellar abundances are the medians of the final stellar MDFs reported by VICE, with error bars showing the 16th and 84th percentiles. As expected, shorter ΔtSF\Delta t_{\mathrm{SF}} produces less SN Ia enrichment, resulting in a stellar population enhanced in Mg relative to Fe. In the most extreme case (ΔtSF=0.1\Delta t_{\mathrm{SF}} = 0.1 Gyr), the star formation ceases just after the SNe Ia channel becomes active (tdelay=0.05t_{\mathrm{delay}} = 0.05 Gyr) and the stellar [Mg/Fe] remains near the [α\alpha/Fe] plateau. As ΔtSF\Delta t_{\mathrm{SF}} increases, the stellar [Mg/H] and [Fe/H] rise toward equilibrium values, with Mg reaching equilibrium more quickly than Fe, reflecting the prompt timescale of CCSN enrichment.

Chemical evolution of the ISM and resulting stellar populations

Figure 5. Chemical evolution of the ISM (dashed lines) and resulting stellar populations (diamonds) as a function of time, for varying star-formation timescales ΔtSF\Delta t_{\mathrm{SF}}. Panels show (a) [Mg/Fe], (b) [Mg/H] and (c) [Fe/H]. Both ISM tracks and stellar abundances are colour-coded by ΔtSF\Delta t_{\mathrm{SF}}, as indicated by the colour bar. ISM tracks are truncated at the epoch when 97.5% of the stellar mass has formed, to avoid the numerical artefacts discussed in the text. Stellar abundances are the medians of the final stellar MDFs, with error bars showing the 16th and 84th percentiles.

Figure 5 also shows a distinct discrepancy between the final ISM and median stellar abundances across all ΔtSF\Delta t_{\mathrm{SF}}. This discrepancy arises because the ISM continues to be enriched by dying stars, even after star formation has slowed substantially. The continuing evolution of the ISM abundances is partly unphysical and driven by the depletion of star-forming gas reservoir; as the mass of hydrogen available for star formation is consumed, the metals returned by evolved stars are concentrated into an increasingly small gas mass, driving high elemental mass fractions. This results in a dramatic boost in the measured [X/H] ratios, which is effectively a numerical artefact distinct from true chemical enrichment. This effect is a consequence of the treatment of gas in the model, since the gas reservoir does not explicitly distinguish between molecular and atomic phases (H2_2, HI). While the star-forming gas mass approaches zero at late times, the total gas mass in a real galaxy is likely non-negligible, with the returned metals distributed across both phases. To aid in visual interpretation and mitigate this effect, the ISM tracks in this figure are truncated at the epoch when 97.5% of the stellar mass has formed. The stellar abundances are also reported at this epoch and shift by ≲0.02\lesssim0.02 dex as a result; elsewhere we adopt the abundances at the 99.9% epoch.

SFE Timescales, Mass-Loading Factors and Delay Times

Although the effects of varying τ⋆\tau_\star and η\eta are well-established (see Andrews et al. 2017), the intense regime of star formation the progenitors of MQGs likely operate in requires us to explore extreme ranges of these parameters. Furthermore, since we deal with integrated-light spectra, we must investigate the impact on the total stellar population. We are primarily interested in whether τ⋆\tau_\star or η\eta can produce stellar abundance variations comparable to those induced by ΔtSF\Delta t_{\mathrm{SF}}. If τ⋆\tau_\star or η\eta produce larger variations than ΔtSF\Delta t_{\mathrm{SF}}, the physical realisations of these parameters may instead be the dominant driver of the α\alpha-enhanced MQG population. Furthermore, while we choose tdelayt_{\mathrm{delay}} to be a global parameter, it has direct impact on the interpretation of ΔtSF\Delta t_{\mathrm{SF}} and hence we investigate it here.

Figure 6 shows how the mass-weighted stellar abundances respond to variations in τ⋆\tau_\star, η\eta and tdelayt_{\mathrm{delay}}, at four of the SF timescales (ΔtSF=0.1,0.2,0.4,0.6\Delta t_{\mathrm{SF}} = 0.1, 0.2, 0.4, 0.6 Gyr) explored in § 4.1, with arrows indicating the direction in which the abundance shifts as each parameter is increased, averaged over ΔtSF\Delta t_{\mathrm{SF}}. Three distinct results are apparent. First, τ⋆\tau_\star and η\eta shift stellar abundances along different axes in [Mg/Fe]–[Fe/H] space: shorter τ⋆\tau_\star moves abundances along the canonical [Fe/H] track to higher [Fe/H], pushing the SN Ia transition (the ‘knee’) to higher metallicity, while lower η\eta shifts abundances orthogonally to this track, increasing both [Mg/Fe] and [Fe/H]. Second, in [Mg/Fe] the magnitude of these shifts are generally smaller than the variation induced by ΔtSF\Delta t_{\mathrm{SF}}, although in [Fe/H] the τ⋆\tau_\star response is much larger. Third, varying τ⋆\tau_\star primarily shifts the stellar abundances of shorter SF timescales (ΔtSF≲0.2\Delta t_{\mathrm{SF}} \lesssim0.2 Gyr), while η\eta dominates on longer timescales (ΔtSF≥0.4\Delta t_{\mathrm{SF}} \ge0.4 Gyr). At ΔtSF=0.1\Delta t_{\mathrm{SF}} = 0.1 Gyr,

Sensitivity of the mass-weighted stellar abundances

Figure 6. Sensitivity of the mass-weighted stellar abundances ([Mg/Fe], [Fe/H]) to variations in the model parameters. Panels (a), (b) and (c) show the changes in abundance ratios in response to variations in the star-formation efficiency timescale τ⋆\tau_\star, mass-loading factor η\eta and SNe Ia delay time tdelayt_{\mathrm{delay}} respectively. In panels (a) and (b), tdelayt_{\mathrm{delay}} is fixed while both τ⋆\tau_\star and η\eta vary across a 2D grid; the same points are coloured by τ⋆\tau_\star in panel (a) and by η\eta in panel (b) to isolate each parameter’s effect. In panel (c), τ⋆\tau_\star and η\eta are fixed at the values shown in the lower left while tdelayt_{\mathrm{delay}} is varied. The points are coloured by the value of the parameter being varied, as indicated by the colour bar in each panel. The coloured arrows indicate the directional shifts in the predicted abundances as the respective parameters (τ⋆\tau_\star, η\eta, tdelayt_{\mathrm{delay}}) are increased.

shortening τ⋆\tau_\star from 1.0 Gyr to 0.05 Gyr raises [Fe/H] by up to ∼1.1\sim1.1 dex, whereas for ΔtSF=0.6\Delta t_{\mathrm{SF}} = 0.6 Gyr, decreasing η\eta from 1.0 to 0.01 produces a ∼0.3\sim0.3 dex increase in [Fe/H] and a 0.2 dex increase in [Mg/Fe].

The response in the chemical abundances to changes in τ⋆\tau_\star and η\eta can be explained by their impact on the gas reservoir. The SFE timescale relates the star-forming gas to the star formation rate via Mg=τ⋆M˙⋆M_g = \tau_\star\dot{M}_\star. At a given SFR, a short τ⋆\tau_\star results in a smaller gas reservoir, meaning that for a fixed yield, the resulting elemental mass fraction is higher, raising both [Mg/H] and [Fe/H]. This effect is most pronounced at short ΔtSF\Delta t_{\mathrm{SF}}, since this is where the gas reservoir is less metal-rich. Conversely, η\eta defines the magnitude of outflows from the gas reservoir, without differential entrainment (ξ=1\xi= 1), meaning that the ISM abundances are preserved in the outflow. Increasing η\eta therefore drives more metals out of the galaxy lowering the mass fractions of all species ([Mg/H], [Fe/H]). Furthermore, since the gas mass MgM_g is fixed as per τ⋆\tau_\star, higher gas inflows are required to sustain the gas reservoir driving in more pristine gas. Together, this lowers the metallicity of the star-forming gas resulting in a stellar population with lower [Fe/H]. This effect is enhanced for longer formation timescales since the impact of η\eta accumulates over the SFH. For short ΔtSF\Delta t_{\mathrm{SF}}, little Fe has been produced and so little has been lost, while for longer formation timescales, the cumulative effect of high η\eta grows, widening the [Fe/H] spread.

The change in relative abundances when modifying η\eta is more subtle. Since each element is equally affected, for short formation timescales (ΔtSF≤0.2 Gyr\Delta t_{\mathrm{SF}} \leq0.2\,\mathrm{Gyr}) where Mg and Fe abundances are dominated by CCSN enrichment, [Mg/Fe] remains fixed. As SNe Ia begin to preferentially enrich Fe, the faster turnover of gas driven by higher η\eta, means that the stellar abundances more closely track recent ISM enrichment. Consequently, the system reaches a lower [Mg/Fe] equilibrium value. Despite these impacts, the variation in the [Mg/Fe] abundance as a consequence of modifying these parameters is smaller than the variation introduced by varying ΔtSF\Delta t_{\mathrm{SF}} on short formation timescales. This ensures that while our choice of τ⋆\tau_\star, η\eta may introduce some scatter, the dominant contribution to variations in this space is robustly driven by differences in star-formation durations.

The same is not true for variations in tdelayt_{\mathrm{delay}}, as shown in panel (c) of Figure 6. As the time delay between the beginning of star formation and the onset of SNe Ia decreases, at fixed ΔtSF\Delta t_{\mathrm{SF}} there is more time for the products of this enrichment channel to pollute the gas in the system effectively shifting the stellar abundance further along the canonical [Mg/Fe] track. This results in a large degree of scatter along the direction of the [α\alpha/Fe] track when tdelay<0.1 Gyrt_{\mathrm{delay}} < 0.1\,\mathrm{Gyr}, which exceeds the variation in [Mg/Fe] across ΔtSF\Delta t_{\mathrm{SF}}, complicating the interpretation of ΔtSF\Delta t_{\mathrm{SF}} from [Mg/Fe] alone. This motivates the use of a SN Ia-independent measure of the SF duration, which may come from elements whose production is dominated by AGB stars (e.g. neutron-capture elements). Although not investigated in this work, we note that similar timescale arguments are relevant for the incorporation of AGB enrichment products, which are only expected to contribute significantly to populations with formation timescales longer than ∼0.3 Gyr\sim0.3\,\mathrm{Gyr} (e.g. Cristallo et al. 2011; Karakas & Lugaro 2016).

Comparison to MQGs at Cosmic Noon

In Figure 7, we compare our modelled [Mg/Fe]-[Fe/H] against the compilation of observed abundances described in § 3.1, with the mass-weighted median abundances (including [Mg/Fe] and other elements) predicted by our models tabulated in Table 2. This dataset spans two distinct epochs and galaxy types. The first consists of individual MQGs (Kriek et al. 2016; Carnall et al. 2024; McConachie et al. 2026) and population surveys (Beverage et al. 2024, 2025) observed near cosmic noon (z∼2z \sim2). The second draws from local ETG observations, including both integrated measurements (Conroy et al. 2014) and central abundances measured within R<Re/8R < R_e/8, where the sample has a mean effective radius ⟨Re⟩=7.2±2.8 kpc\langle R_e \rangle= 7.2 \pm2.8\,\mathrm{kpc} along the semi-major axis (Gu et al. 2022). Following the two-phase formation scenario (e.g. Oser et al. 2010), the assembled MQGs observed at cosmic noon are assumed to evolve primarily by dry minor merging (e.g. Bezanson et al. 2009; Naab et al. 2009) to become the cores of local ETGs. This evolutionary link implies that their stellar populations should share similar properties, motivating a comparison of their chemical abundances.

Figure 7. Predicted [Mg/Fe] vs [Fe/H] for fiducial MQG models as a function of SF timescale. The diamonds show the mass-weighted median stellar abundances reported by the stellar MDFs, coloured by the duration of their star formation (ΔtSF=t90−t10\Delta t_{\mathrm{SF}} = t_{90} - t_{10}) as indicated by the colour bar. The coloured arrows indicate the directional shifts in the predicted abundances as the respective parameters (tdelayt_{\mathrm{delay}}, τ⋆\tau_\star, η\eta) are increased. The model is overlaid on a sample of MQGs observed at cosmic noon (light grey symbols), comprising both individual measurements: ZF-UDS-7329 (pentagon; (Carnall et al. 2024)), Eridu (hexagon; (McConachie et al. 2026)), COSMOS 11494 (diamond; (Kriek et al. 2016)), as well as population surveys from JWST-SUSPENSE (squares; (Beverage et al. 2025)), Heavy Metal (triangles; (Beverage et al. 2024)), and lensed MQGs (octagons; (Jafariyazani et al. 2020, 2025)). For comparison, the dark grey symbols show measurements of local ETGs, including both integrated light (circles; (Conroy et al. 2014)) and central (R<Re/8R < R_{\mathrm{e}}/8) abundances (upside-down triangles; (Gu et al. 2022)).

ΔtSF\Delta t_{\mathrm{SF}}[M/H]ᵃ[Fe/H][Mg/Fe][Si/Fe][Ca/Fe][Ti/Fe]
0.1-0.23-0.500.540.280.290.21
0.2-0.04-0.280.500.270.270.20
0.40.10-0.050.390.170.190.10
0.60.150.060.320.120.140.04
0.80.180.120.270.080.110.00
1.00.200.160.240.060.09-0.02

Table 2. Mass-weighted median chemical abundances of MQGs from VICE as a function of SF formation timescale (ΔtSF=t90−t10\Delta t_{\mathrm{SF}} = t_{90} - t_{10}).

Figure 7 shows that our model reproduces the observed [Mg/Fe]-[Fe/H] distribution of cosmic noon MQGs across a range of SF timescales. A formation duration of ΔtSF∼0.2−0.4 Gyr\Delta t_{\mathrm{SF}} \sim0.2 - 0.4\,\mathrm{Gyr} best matches the average observed abundance (in both [Mg/Fe] and [Fe/H]), with the most α\alpha-enhanced MQGs requiring rapid timescales around ∼0.1 Gyr\sim0.1\,\mathrm{Gyr}, suggesting variation in the formation across the observed population. As discussed in § 4.2 and indicated by the parameter arrows, variations in τ⋆\tau_\star tend to shift our predictions along the tracks mapped out by changes in ΔtSF\Delta t_{\mathrm{SF}} but with smaller magnitude, and so different fiducial values do not change the qualitative preference for short-formation timescales.

Several MQGs at cosmic noon, particularly in the survey samples, host stellar populations with lower [Mg/Fe] and [Fe/H] than our fiducial models predict. This offset is consistent with intrinsic galaxy-to-galaxy variation in the local parameters across the MQG progenitor population: increasing η\eta at fixed ΔtSF\Delta t_{\mathrm{SF}} moves the predicted tracks toward lower [Fe/H] and [Mg/Fe], as indicated by the parameter arrows in Figure 7 and the trends shown in Figure 6. Despite the diversity of outflow magnitudes experienced by these galaxies, the formation timescales implied by higher η\eta are still consistent with ΔtSF∼0.1−0.4 Gyr\Delta t_{\mathrm{SF}} \sim0.1 - 0.4\,\mathrm{Gyr}. While the [Fe/H] axis depends on CCSN and SN Ia Fe yields and hence carries uncertainty, the [Mg/Fe] axis is controlled by the ratio yMgCC/yFeCCy_{\mathrm{Mg}}^{\mathrm{CC}}/y_{\mathrm{Fe}}^{\mathrm{CC}}, which is well-constrained by our prior calibration with Sculptor. Even allowing for diversity in η\eta across the progenitor population, the formation timescales implied by the [Mg/Fe] axis remain robust.

These rapid formation timescales are consistent with multiple lines of evidence in the literature. The high [Mg/Fe] enhancement observed in cosmic noon MQGs has previously been argued to require formation durations of 100−500 Myr100 - 500\,\mathrm{Myr} from a consideration of SN Ia enrichment timescales alone (Kriek et al. 2016; Beverage et al. 2025). Recent SFH-based inferences from SED fitting of other MQG samples also return formation durations in a similar range, though they generally predict shorter timescales (e.g. Hamadouche et al. 2026; Leung et al. 2026). Now using a full treatment of the chemical evolution, our analysis confirms this rapid-formation picture beyond timing arguments and SFH inferences, supporting the view that MQGs at cosmic noon can plausibly form the bulk of their stellar mass on timescales of order a few hundred Myr.

While our model captures the average behaviour of the population, with η\eta accounting for intrinsic variation between galaxies, reproducing the most [Mg/Fe]-enhanced galaxies at their observed [Fe/H] (e.g. [Mg/Fe] ∼0.6\sim0.6 at [Fe/H] ∼−0.4\sim-0.4) pushes the fiducial framework close to its physical limits. As indicated by the parameter arrows in Figure 7, reaching this extreme region of chemical space requires not only the shortest formation timescales (ΔtSF∼0.1 Gyr\Delta t_{\mathrm{SF}} \sim0.1\,\mathrm{Gyr}), but also driving τ⋆\tau_\star and η\eta to their lower extremes, with η→0\eta\rightarrow0 corresponding to the closed-box limit. Further tuning of the yields offers one route to better match these galaxies, however, the yields implemented here are motivated by the chemistry of the Sculptor dSph, which shows similar [Mg/Fe] enhancements as MQGs in metal-poor stars. The yields adopted here are therefore already favourable to producing [Mg/Fe]-enhanced populations compared to the MW, which is typically assumed to have [Mg/Fe]CC=0.45_{\mathrm{CC}} = 0.45 (e.g. Weinberg et al. 2024). Another alternative is a steeper IMF, which would boost the relative

Notes. ᵃ The metallicity [M/H] = log⁡(Z/Z⊙)\log(Z/Z_{\odot}) is the total metallicity by mass of both gas and stars tracked by the model. We use Z⊙=0.014Z_{\odot} = 0.014 from Asplund et al. (2009).

CCSN (and hypernova) contribution and increase the yield of Mg relative to Fe, without requiring extreme values of τ⋆\tau_{\star} or η\eta. Figure 7 also shows that the intrinsic evolution of MQGs formed over longer formation timescales may naturally reproduce the observed abundances of the cores of local ETGs. As ΔtSF\Delta t_{\mathrm{SF}} increases, our model approaches an equilibrium value in [Mg/Fe]-[Fe/H] space that is close to the value measured for the centres of ETGs. Two interpretations follow. First, the observed cosmic noon MQG population may not be the direct progenitors of ETG cores, but rather the earliest-forming, most extreme members of the broader quiescent population; the average ETG core would then be a quiescent galaxy that formed over a longer timescale. Alternatively, if the currently selected and observed population of MQGs at cosmic noon are the progenitors of local ETG cores, their stellar populations must undergo additional chemical evolution. This evolution could occur through rejuvenated or hidden ongoing star formation or dry major merging with galaxies formed on longer timescales, both of which are seen in simulated galaxies (e.g. Lagos et al. [82]; Remus & Kimmig 2025; Chandro-Gómez et al. 2026). These interpretations are consistent with our model and the observed chemistry and both may be partly responsible; however, dry minor merging, which adds metal-poor stars from low-mass galaxies with lower [Mg/H] and [Fe/H] (Kirby et al. 2013), is unlikely to significantly alter the chemistry in the required direction, as discussed in Beverage et al. (2024).

Figure 8 compares our model predictions to the observed multi-element abundances from the JWST-SUSPENSE survey (Slob et al. 2024; Beverage et al. 2025) for Mg, Si, Ca and Ti as a function of atomic number. Consistent with Figure 7, we find a SF duration of ∼0.4\sim0.4 Gyr reproduces the median [Mg/Fe] of cosmic noon MQGs observed by SUSPENSE, with a slight preference towards longer ΔtSF\Delta t_{\mathrm{SF}}, and this agreement extends to [Si/Fe] and [Ca/Fe]. The agreement in [Si/Fe] between our model and observations is particularly interesting since we did not rigorously calibrate the Si yield, beyond fitting by visual inspection. The broad agreement in formation timescales for these α\alpha-elements suggests that their yields are mostly metallicity-independent and applicable across different star-forming environments.

Model abundance patterns for Mg, Si, Ca and Ti as a function of atomic number

Figure 8. Model abundance patterns ([X/Fe] for Mg, Si, Ca and Ti) for fiducial MQG models as a function of atomic number. The diamonds show the mass-weighted median stellar abundances reported by the stellar MDFs, coloured by the duration of their star formation (ΔtSF=t90−t10\Delta t_{\mathrm{SF}} = t_{90} - t_{10}) as indicated by the colour bar. Overlaid as light grey squares are the observed median abundance patterns of 1≤z≤31 \le z \le3 MQGs from the JWST-SUSPENSE survey ((Slob et al. 2024; Beverage et al. 2025)), with individual galaxy abundances shown in the background as smaller grey squares.

Despite this broad agreement, Sculptor appears to be a poor model for Ti production in MQGs, with our calibrated yield underproducing the observed [Ti/Fe] by 0.25 dex at ΔtSF=0.4\Delta t_{\mathrm{SF}} = 0.4 Gyr. This offset could be due to several effects which may compound. First, Ti is susceptible to departures from local thermodynamic equilibrium (LTE), with [Ti/Fe] abundances differing by up to ∼0.1\sim0.1 dex between LTE and non-LTE treatments (Bergemann 2011), though it is unclear how this translates from single stars to integrated stellar populations across different regimes of galaxies. In the case of ionisation states, we adopt [Ti/Fe] abundances from Hill et al. (2019) as these are considered more robust than those derived from Ti I (Bergemann 2011; Mallinson et al. [89]), whereas Beverage et al. (2025) derive [Ti/Fe] from full-spectrum fitting in a wavelength range where a blended Fe I + Ti I feature at 4540 Å is the primary Ti diagnostic. However, there is no major difference in the MDFs of [Ti/Fe] and [TiII/Fe] abundances from Hill et al. (2019) that could significantly change the calibrated yields from Sculptor, although we cannot apply this test to MQGs. Second, stars in Sculptor may be deficient in [Ti/Fe]. The core-collapse plateau [Ti/Fe]cc_{\mathrm{cc}} in Sculptor lies ∼0.15\sim0.15 dex below the plateau in the Milky Way (Hill et al. (2019)). Such a difference would lower our calibrated yTiccy_{\mathrm{Ti}}^{\mathrm{cc}} and result in even our shortest SF timescale, which reflects enrichment from only CCSNe, underproducing [Ti/Fe] as shown in Figure 8. Lastly, observational uncertainties likely contribute to these discrepancies; as noted by Beverage et al. (2025), individual measurements for Ti are subject to large errors of ∼0.15−0.2\sim0.15 - 0.2 dex and with a small number of galaxies, this leads to a larger uncertainty of 0.08 dex on the median [Ti/Fe] compared to that of other elements. Altogether, these effects are comparable to the 0.25 dex offset, suggesting that the disagreement with our model reflects systematics in abundance measurements of Ti and Sculptor’s lower [Ti/Fe]cc_{\mathrm{cc}} rather than a difference in ΔtSF\Delta t_{\mathrm{SF}}.

Next we consider the absolute abundances of [Mg/H] and [Fe/H] as predicted by our model, presented in Figure 9. While the [Mg/Fe]-[Fe/H] plane tests the relative timing of CCSN and SN Ia enrichment, [Mg/H]-[Fe/H] adds a constraint on the absolute Mg and Fe predicted by our models. Here we show the usual mass-weighted median stellar abundances from the model MDFs, overlaid on the abundances of z∼2z \sim2 MQGs observed by Heavy Metal and SUSPENSE and the centres of ETGs from Gu et al. (2022). Our fiducial models show consistent behaviour with the observational data, reproducing the trends in the observed population across the parameter space. However at fixed [Fe/H], we tend to slightly overproduce [Mg/H] for several galaxies in the cosmic noon population.

Figure 9. Predicted [Mg/H] vs [Fe/H] for fiducial MQG models, with colours, symbols and arrows following the formatting of Figure 7.

This offset could reflect intrinsic variations in the local parameters of these galaxies, namely η\eta. This explanation is feasible for MQGs with [Fe/H] ∼−0.25\sim-0.25, where outflows can drive out significant amounts of SNe Ia products, however, galaxies with low [Mg/H] at [Fe/H] ∼−0.5\sim-0.5 are difficult to reconcile under this interpretation. On such short formation timescales, Mg and Fe are produced almost exclusively by CCSNe in roughly fixed proportion, so changes in τ⋆\tau_\star or η\eta shift [Mg/H] and [Fe/H] together rather than separately. The observed excess of [Mg/H] at fixed [Fe/H], for ΔtSF∼0.1\Delta t_{\mathrm{SF}} \sim0.1 Gyr, hence cannot be reconciled by varying the local parameters alone. Alternatively, a more bottom-heavy IMF could produce less Mg at fixed Fe by reducing the number of high-mass CCSN progenitors. A shorter tdelayt_{\mathrm{delay}} could also help by allowing SN Ia Fe enrichment to begin earlier, though it is not clear this is physically plausible and it would require modifications to τ⋆\tau_\star and η\eta. The overall offset is also partly attributable to uncertainties in the measurement of α\alpha abundances in MQGs. As with other abundance comparisons in this work, the discrepancy most likely reflects a combination of these factors, with the dominant contributions coming from uncertainties in the CCSN Mg and Fe yields and in the abundance measurements themselves.

Recently, the opposite was identified by Gountanis et al. (2025) who found it difficult to reproduce highly-enhanced [Mg/H] at near-solar [Fe/H] for z∼0z \sim0 and z∼0.7z \sim0.7 ellipticals. In particular, using Milky Way-motivated, empirical nucleosynthetic yields, they could not reproduce the high [Mg/Fe] and [Mg/H] of ETGs with any combination of τ⋆\tau_\star or η\eta. To resolve this, they increased the CCSNe Mg yield to achieve a plateau of [Mg/Fe]CCSNe=0.55_{\mathrm{CCSNe}} = 0.55 rather than the standard 0.45. In our models, because our framework adopts yields physically motivated by the chemistry of Sculptor, which is enhanced in [Mg/Fe] relative to the MW as previously discussed, our models naturally capture this high [Mg/H] space. Hence, while we approach the problem from a different angle, we arrive at a similar conclusion; to reproduce the abundance trends seen in [Mg/Fe] and [Mg/H] at near-solar [Fe/H] for MQGs, the core-collapse contribution of Mg relative to Fe must be greater than what is generally observed/assumed in the Milky Way. Furthermore, the SF must be highly efficient (short τ⋆\tau_\star) with significant metal retention (η∼0\eta\sim0). This is also consistent with Gountanis et al. (2025) finding that high-σ\sigma ETGs – which formed earlier and over shorter timescales – require short SFE timescales and negligible mass-loading to match observed abundances. We discuss the physical viability of this in § 6, however it is worth noting that other changes, for example modifying the IMF or self-terminating bursts (both discussed by Gountanis et al. (2025)), may provide better agreement with the observed abundances when adopting standard MW-motivated yields.

Finally, beyond these abundance trends for individual elements, we find that our modelled metallicities ([M/H] = log⁡(Z/Z⊙)\log(Z/Z_\odot); see Table 2) are broadly consistent with observed samples near or across these epochs (e.g. (Beverage et al. 2025; Nanayakkara et al. 2025; Hamadouche et al. 2026; Leung et al. 2026)). However, we caution that direct comparisons to existing literature metallicities are complicated. Many current metallicity estimates for high-redshift MQGs are derived using stellar population models that assume solar-scaled abundance patterns. However, because these galaxies are α\alpha-enhanced, using solar-scaled templates can introduce systematic biases into the derived total metallicities (e.g. (Beverage et al. 2025)). Therefore, while the consistency between our models and available metallicity measurements is encouraging, the observed metallicities may not be representative of the true chemistry of these galaxies.

Discussion

Implications for the Formation Timescales and Properties of MQGs

Our analysis demonstrates that it is possible to reproduce the chemical abundances of z∼2z \sim2 massive quiescent galaxies by assuming they form the majority of stars in rapid, intense star-formation episodes. By limiting the time available for star formation following the onset of SNe Ia, the incorporation of Fe-rich ejecta into the ISM is slowed, producing the observed [Mg/Fe]-enhanced populations. On average, we find a timescale of ΔtSF∼0.1−0.4\Delta t_{\mathrm{SF}} \sim0.1 - 0.4 Gyr successfully produces mass-weighted chemical abundance ratios ([Mg/Fe], [Si/Fe], [Ca/Fe]) consistent with current observational samples at cosmic noon. For the most extremely α\alpha-enhanced galaxies, shorter formation timescales are feasible (e.g. ΔtSF∼0.1\Delta t_{\mathrm{SF}} \sim0.1 Gyr); however, arbitrarily short timescales are precluded given the non-negligible Fe enrichment already present in these galaxies. Broadly, this is consistent with the estimates given by Kriek et al. (2016) and Beverage et al. (2025), who use SNe Ia enrichment timescales to constrain formation durations, as well as those derived from inferred SFHs (e.g. (Glazebrook et al. 2017; Hamadouche et al. 2026; Leung et al. 2026)).

Simulations also paint a similar picture; in their comparison of the SFHs of MQGs from several existing simulations, Lagos et al. [83] found t90−t50∼100−300t_{90} - t_{50} \sim100 - 300 Myr. Although t90−t50t_{90} - t_{50} samples only half the burst, an approximately symmetric single-burst SFH implies ΔtSF≈2(t90−t50)∼200−600\Delta t_{\mathrm{SF}} \approx2(t_{90} - t_{50}) \sim200 - 600 Myr, consistent with our derived range. While we cannot entirely rule out the very short bursts of SF suggested by some parametrisations of the SFH (e.g. 50 Myr; Valcinom et al. (2020)), we do not require the implied SFRs of M˙≳1000 M⊙ yr−1\dot{M} \gtrsim1000\,M_\odot\,\mathrm{yr}^{-1} to reproduce the observed chemical abundances.

We note this timescale relies on several assumptions in our modelling, including the choice of SN Ia delay time. Figure 6 shows that when tdelayt_{\mathrm{delay}} becomes sufficiently short (≲0.1\lesssim0.1 Gyr), the shift in [Mg/Fe] can be on the order of 0.1 dex for the shortest ΔtSF\Delta t_{\mathrm{SF}}. Our adoption of tdelay=0.05t_{\mathrm{delay}} = 0.05 Gyr is motivated by theoretical and observational expectations that SNe Ia progenitor lifetimes can be as short as 50 Myr (e.g. Maoz et al. 2014; Castrillo et al. 2021), and by our modelling of the Sculptor dSph, which requires a prompt enrichment channel to reproduce its observed abundances. This is broadly consistent with los Reyes et al. (2022) and Tang et al. (2023) who find delay times of 100 Myr for Sculptor, though los Reyes et al. (2022) also show that shorter delay times produce worse fits to the observed abundances. Ultimately, the minimum delay time remains a subject of strong debate, with uncertainties in the progenitor channels allowing for a range of values (see Ruiter & Seitenzahl (2025), and references therein); however, an onset of SNe Ia as early as ∼0.05\sim0.05 Gyr remains physically plausible. These systematic uncertainties highlight a need for a SNIa-independent clock, which we suggest could come from neutron-capture elements, in agreement with Conroy et al. (2013), and we leave this for future work.

The requirement for such rapid formation timescales (ΔtSF∼0.1−0.4\Delta t_{\mathrm{SF}} \sim0.1 - 0.4 Gyr) raises the question of whether massive galaxies can physically assemble and quench so quickly. Moreover, the inferred SFHs of several MQGs observed using JWST continue to suggest that their progenitors had assembled a bulk of their mass (M⋆∼1010−1011 M⊙M_\star\sim10^{10} - 10^{11}\,M_\odot) beyond z≳10z \gtrsim10, when the universe was less than 500 Myr old (e.g. Glazebrook et al. 2024; Carnall et al. 2024; de Graaff et al. 2025; McConachie et al. 2026). While the physical mechanisms for how this may occur are not understood, especially given dark matter halos of a compatible size had not yet formed by this time (e.g. Behroozi & Silk 2018; Boylan-Kolchin 2023; Glazebrook et al. 2024), the chemical abundances are con-

ends_mid=1 sistent with SF timescales that would allow these galaxies to form prior to z∼10z \sim10. However, galaxies with such masses have yet to be observed and require SFRs significantly higher (100−500×100 - 500\times) than other spectroscopically-confirmed galaxies (M⋆∼107−108M⊙M_\star\sim10^7 - 10^8 M_\odot) at a similar redshift (e.g. Curtis-Lake et al. 2023). Furthermore, modelling of the progenitor SED for the galaxy ‘Eridu’ by McConachie et al. (2026) revealed no existing JWST photometric candidates in either GOODS field with comparable brightness or properties.

While the implied SFRs are inconsistent with current observations of galaxies at z>10z > 10, our chemical evolution modelling has no constraining power on the formation ages of these galaxies. Consequently, the tension with existing observations is largely dependent on the accuracy of the SFHs and absolute ages derived from SED fitting, which have known degeneracies with metallicity (Worthey 1994; Conroy 2013) that can impact the formation ages of MQGs (e.g. de Graaff et al. 2025). If the formation redshifts are overestimated and these galaxies formed at slightly later epochs (z∼5−8z \sim5 - 8), this tension is significantly alleviated. Alternatively, these galaxies may represent the extremes of the MQG population and do not actually reflect the typical formation times of the bulk of these galaxies. Indeed, in the analyses of the quiescent population by Beverage et al. (2024); Park et al. (2024); Beverage et al. (2025); Nanayakkara et al. (2025), the majority of MQGs exhibit typical formation redshifts of z<10z < 10 with only a few having formed the bulk of their stars earlier.

The later formation epoch points to dusty star-forming galaxies (DSFGs) or submillimetre galaxies (SMGs; Casey et al. 2014; Hodge & da Cunha 2020) as the potential progenitors of MQGs at cosmic noon, as has been proposed previously (e.g. Toft et al. 2014; Forrest et al. 2020; Valentino et al. 2020; Araya-Araya et al. 2026). Such galaxies have sufficiently high star-formation rates ($10^2 - 10^3 M_\odot,\mathrm{yr}^{-1}$) – also predicted by simulations of MQG progenitors (Lagos et al. [83]) – while also matching the comoving number densities (Valentino et al. 2020) and stellar masses (e.g. Zavala et al. 2026) of MQGs around z∼2−4z \sim2 - 4. Very recently, so-called ‘NIRCam-dark’ DSFGs have been shown by Sun et al. (2026) to exactly match the stellar mass and SFRs (80−250M⊙ yr−180 - 250 M_\odot\,\mathrm{yr}^{-1}) of the assumed progenitors of z∼4z \sim4 MQGs, while also aligning well as possible descendants of the JWST ‘Blue Monster’ population at $z > 10 (e.g. Ziparo et al. 2023). Our chemical evolution modelling demonstrates that the SFHs required to reproduce the bulk of the observed chemical abundances (ΔtSF=0.2−0.4 Gyr\Delta t_{\mathrm{SF}} = 0.2 - 0.4\,\mathrm{Gyr}; see Figure 3) imply average SFRs (200−400M⊙ yr−1200 - 400 M_\odot\,\mathrm{yr}^{-1}) that are consistent with the measured SFRs of these NIRCam-dark galaxies (Bing et al. 2026; Sun et al. 2026). However, the rapid accumulation of mass in our model challenges the connection to z>10z > 10 ‘Blue Monsters’, as an evolutionary link to these earlier structures would imply a more extended star-formation history than our chemical constraints suggest. Yet, Blue Monsters have low stellar masses (M⋆∼108M⊙M_\star\sim10^8 M_\odot) contributing few stars to the final MQG stellar budget (M⋆∼1011M⊙M_\star\sim10^{11} M_\odot), and hence a significant burst at later epochs could still produce the observed chemistry of MQGs.

Our modelling also demonstrates that MQG formation timescales ΔtSF≲0.1 Gyr\Delta t_{\mathrm{SF}} \lesssim0.1\,\mathrm{Gyr} are not favoured by the measured chemical abundances around cosmic noon. Almost all MQGs in our collated sample show signs of enrichment from SNe Ia, which would be absent if the SF-burst is rapid enough and forms the majority of stars within the galaxy. Instead, our modelling implies that hyper-luminous DSFGs with sustained SFRs significantly exceeding 1000M⊙ yr−11000 M_\odot\,\mathrm{yr}^{-1} (e.g. Riechers et al. 2013; Toft et al. 2014; Marrone et al. 2018) are effectively ruled out as the primary progenitors for the bulk of the population, except perhaps for the most extremely [α/Fe][\alpha/\mathrm{Fe}]-enhanced galaxies. To build a typical stellar mass of ∼1011M⊙\sim10^{11} M_\odot at such extreme SFRs would require less than 100 Myr resulting in a near-purely α\alpha-enhanced population that lacks the observed enrichment in Fe from SNe Ia. However, we note that shorter formation timescales are plausible if the minimum SN Ia delay time is shorter than our adopted value, or the true CCSNe Fe yield, yFeCCy_{\mathrm{Fe}}^{\mathrm{CC}}, is higher. While our model supports the general DSFG-to-MQG connection, it suggests the average MQG forms from the moderately active, dusty population rather than the extreme hyper-starbursts.

Our finding that extremely short SF timescales are disfavoured is in tension with the timescales inferred from recent population studies of MQGs (e.g. Nanayakkara et al. 2025; Leung et al. 2026). While their SED-derived SFHs recover individual durations consistent with ours, the population generally trends towards shorter formation periods. However, we caution that the SFHs derived purely from SED modelling are fit without specific chemical constraints, which can lead to misalignment between the inferred formation history and the observed chemical abundances. Leung et al. (2026) highlight this for the galaxy PRIMER EXCELS-55410 (ZF-UDS-7329), which has an inferred SF timescale of ∼300\sim300 Myr despite showing signs of moderate enrichment from SNe Ia. While this scenario is plausible under our model, Leung et al. (2026) find other galaxies with near-solar [Fe/H] with even shorter formation timescales. Neglecting the uncertainty in nucleosynthetic yields, achieving this significant [Fe/H] enrichment on such rapid timescales would require an extremely prompt SN Ia component (<50< 50 Myr), which is not supported observationally (Maoz et al. 2014, and references therein). This points to the need for SED modelling to account for chemical constraints as a means of inferring the formation of such massive galaxies.

Another useful physical comparison from our modelling results involves the local properties governing the enrichment of the ISM – the SFE timescale (approximately τmol\tau_{\mathrm{mol}}) and the mass-loading factor. In general, we find that very short SFE timescales (τ⋆∼0.05 Gyr\tau_\star\sim0.05\,\mathrm{Gyr}) are necessary to reproduce the observed abundance patterns at cosmic noon. Such depletion times are commonly observed in galaxies experiencing significant star formation, both locally in ULIRGs (e.g. Solomon et al. 1997; Combes et al. 2013) and at early times in DS-FGs (e.g. Aravena et al. 2016; Reuter et al. 2020; Ikarashi et al. 2022; Prajapati et al. 2026), implying that the extreme efficiencies required by our model are plausible, and that comparisons between τ⋆\tau_\star and τmol\tau_{\mathrm{mol}} are appropriate. We note that the depletion time derived for the individual NIRCam-dark galaxy AC-2168 is slightly longer (τmol∼150\tau_{\mathrm{mol}} \sim150 Myr; Bing et al. 2026) than our fiducial 50 Myr, however, given the scatter in the abundance planes of the MQG population, which we suspect to be driven by variations in gas fractions and star-formation efficiencies, this factor of ∼3\sim3 difference likely reflects intrinsic diversity rather than inconsistency with the model.

Our models allow for negligible mass loading (η∼0\eta\sim0), implying MQG progenitors evolve as effectively closed boxes retaining the bulk of their metal-enriched ejecta. Such strong retention is necessary to reproduce the near-solar [Fe/H] at high [Mg/Fe] (see Figure 6), however η∼1\eta\sim1 is suitable for MQGs at lower [Fe/H]. For the most α\alpha-enhanced galaxies, η∼0\eta\sim0 is in tension with observational constraints from DSFGs at z>4z > 4, where typical mass-loading factors are near unity (Ginolfi et al. 2020; Spilker et al. 2020; Herrera-Camus et al. 2021). However, it is important to distinguish between outflows leaving the galactic disk and those escaping the galaxy’s halo.

Observed mass-loading factors are typically derived from outflows on the disk-scale rather than those leaving the virial radius, with many studies finding most of the outflowing gas does not reach escape velocity and is actually retained by the galaxy (e.g. Davies et al. 2019; Fluetsch et al. 2019; Rodríguez Del Pino et al. [120]). In the implementation in VICE, η\eta effectively describes the loss of metals from the box, rather than the movement of gas within the halo of a galaxy. Consequently, our requirement for η∼0\eta\sim0 does not rule out the existence of strong stellar feedback driving galactic winds. Instead, it aligns with a ‘galactic fountain’ scenario where gas ejected from the disk rains back down on the galaxy and is re-accreted into the ISM on timescales rapid enough to fit within the SFHs of our fiducial MQGs (e.g. 100 – 300 Myr; Oppenheimer & Davé 2008; Anglés-Alcázar et al. 2017). Therefore, our modelling suggests that while the progenitors of MQGs retain nearly all of their metals, they may still have intense disk-scale outflows driven by stellar feedback.

Overall, our chemical evolution modelling implies a coherent picture for the formation of MQGs, summarised in Figure 10: such galaxies likely formed over short timescales (ΔtSF=0.1−0.4\Delta t_{\mathrm{SF}} = 0.1 - 0.4 Gyr), sustained by rapid gas consumption (τ⋆∼0.05\tau_\star\sim0.05 Gyr) and effectively ‘closed-box’ metal retention (η∼0\eta\sim0). This picture of formation aligns with the moderately active DSFG population, and in particular with the recently observed NIRCam-dark galaxies identified as candidate MQG progenitors (Bing et al. 2026; Sun et al. 2026), while disfavouring sustained, hyper-luminous starbursts as primary progenitors of the typical MQG.

Illustrative summary of the primary parameters governing chemical enrichment

Figure 10. Illustrative summary of the primary parameters (ΔtSF\Delta t_{\mathrm{SF}}, τ⋆\tau_\star, η\eta) governing the chemical enrichment of MQGs in our model. Each parameter is defined and its physical role described in a text box (ΔtSF\Delta t_{\mathrm{SF}} top, τ⋆\tau_\star left, η\eta right), while the central graphic illustrates their physical effect on a model galaxy. The side panels show the effect of increasing τ⋆\tau_\star (left) and η\eta (right) on the median stellar abundances of a model galaxy, relative to the observed MQGs (grey contour enclosing 83% of the sample). Each curve traces a sequence of independent models spanning a range of ΔtSF\Delta t_{\mathrm{SF}} (colour bar) for the fiducial value of each parameter (τ⋆=0.05\tau_\star= 0.05 Gyr, η=0.01\eta= 0.01) and a secondary representative value (τ⋆=1.0\tau_\star= 1.0 Gyr, η=2.0\eta= 2.0); the curves are the locations of model medians and are not evolutionary tracks of a single galaxy. ΔtSF\Delta t_{\mathrm{SF}} defines the formation timescale in terms of t90−t10t_{90} - t_{10}, controlling the balance of CCSN enrichment (large blue supernovae) vs. SN Ia enrichment (small orange supernovae). For the fiducial τ⋆\tau_\star and η\eta, our modelling shows that ΔtSF∼0.1−0.4\Delta t_{\mathrm{SF}} \sim0.1 - 0.4 Gyr is compatible with the bulk of the MQG population. In terms of the physical parameters τ⋆\tau_\star and η\eta, shorter SFE timescales τ⋆\tau_\star (approximately τmol\tau_{\mathrm{mol}}) result in smaller gas reservoirs and a quicker turnover of gas, leading to more rapid enrichment, particularly in [Fe/H]; rapid depletion times of ∼0.05\sim0.05 Gyr are necessary to produce the metal-rich stellar populations hosted by MQGs at z∼2z \sim2. Increasing η\eta allows a larger fraction of metals to escape the box (halo), resulting in a significant decrease in [Fe/H] and a moderate decrease in [Mg/Fe] at fixed ΔtSF\Delta t_{\mathrm{SF}}. To reach the metallicities observed in MQGs at cosmic noon, a high degree of metal retention is required (η∼0.01\eta\sim0.01) with moderate variations up to values of order unity able to reproduce the scatter in the observed population.

Systematics in Nucleosynthetic Yields and the High- Plateau

Constraining the star-formation histories of galaxies from chemical abundances is fundamentally limited by the systematic uncertainties inherent to theoretical nucleosynthetic yields (Kobayashi et al. 2020). Calculating these yields from parametrised stellar evolution models is computationally challenging and relies on physical assumptions regarding the treatment of convection, mass loss and explosion mechanisms. Such assumptions and unconstrained parameters drive the significant uncertainties in stellar yield calculations, to the extent that different yield sets can predict vastly different abundance ratios for the same stellar progenitor (e.g. Romano et al. 2010; Karakas & Lattanzio 2014). Consequently, the choice of yield table introduces a systematic error that may limit the interpretation of inferred SFHs.

In the context of our GCE modelling, the adopted nucleosynthetic yields and IMF control the relative importance of different enrichment channels on varying timescales, as well as the individual enrichment of each element. Typically in GCE models, yields are chosen from supernovae tables such as the ones adopted in this work for CCSNe (NKT13) and SNe Ia (Iwamoto et al. 1999). Rather than rely directly on these tables, we motivate our adopted yields for the CCSN and SN Ia contributions to Mg, Ca, Si, Ti and Fe based on our chemical analysis of the dSph galaxy Sculptor. This calibration using a resolved stellar population allows us to empirically anchor the population-averaged yields of these key elements, reducing the reliance on theoretically calculated quantities. However, the adopted yields are IMF-averaged constant values, meaning that mass-dependent enrichment channels (e.g. hypernovae), metallicity dependencies and rotational mixing are neglected, adding some uncertainty. The latter two effects are expected to have minimal impact on the yields of the elements considered here, whereas excluding mass-dependent yields may have moderate consequences if the IMF were to vary between galaxies, a possibility raised by some studies of massive ETGs (see Smith 2020).

While Sculptor offers an ideal testbed for developing our model and constraining empirically-motivated yields, this choice comes with its own caveats. Despite the theoretical assumption that the nucleosynthetic yields between galaxies should be constant, there are notable differences in the abundances of stars in the MW and Sculptor on the [α/Fe]cc[\alpha/\mathrm{Fe}]_{\mathrm{cc}} plateau. In the MW, the typical [Mg/Fe][\mathrm{Mg}/\mathrm{Fe}] of stars prior to SNe Ia enrichment is [Mg/Fe]cc∼0.45[\mathrm{Mg}/\mathrm{Fe}]_{\mathrm{cc}} \sim0.45 (e.g. Weinberg et al. 2024), whereas in Sculptor we identify a higher plateau around [Mg/Fe]cc=0.6[\mathrm{Mg}/\mathrm{Fe}]_{\mathrm{cc}} = 0.6 from our modelling. This higher plateau is not necessarily limited to Sculptor, but is also observed among other Local Group dSphs (Kirby et al. 2011). Though from the VLT/FLAMES data which was used as our observational reference, Hill et al. (2019) found the observed [Mg/Fe][\mathrm{Mg}/\mathrm{Fe}] to be consistent with the MW, with scatter around the plateau ranging from [Mg/Fe]∼0.3−0.75[\mathrm{Mg}/\mathrm{Fe}] \sim0.3 - 0.75. There is less scatter observed in the Keck-II/DEIMOS data from Henderson et al. (2025), however with a larger sample of stars they find many stars with [Mg/Fe]>0.4[\mathrm{Mg}/\mathrm{Fe}] > 0.4 with some as high as [Mg/Fe]∼1.0[\mathrm{Mg}/\mathrm{Fe}] \sim1.0 at low-metallicity. To avoid over-fitting to the stars at low-metallicity, we have employed a high-[Mg/Fe] cut at [Mg/Fe]<0.75[\mathrm{Mg}/\mathrm{Fe}] < 0.75, however we still obtain a high plateau. Despite this, our model trends in the [Mg/Fe[\mathrm{Mg}/\mathrm{Fe}]-[Fe/H][\mathrm{Fe}/\mathrm{H}] space are consistent with other GCE models of Sculptor and show an excellent fit to the observed data (Kirby et al. 2011; de los Reyes et al. 2022).

Although the physical drivers of the elevated [Mg/Fe]cc[\mathrm{Mg}/\mathrm{Fe}]_{\mathrm{cc}} in local dSphs are not fully understood, Kirby et al. (2011) proposed a coupling between early-time gas masses and IMF sampling of CCSN progenitors. In this scenario, the initial gas reservoir is sufficiently small that the first generation of supernovae, particularly the most massive, short-lived progenitors that produce relatively more α\alpha-elements per Fe (NKT13), rapidly enriches the ISM before lower-mass CCSNe contribute. Consequently, these very early SNe imprint a very high [Mg/Fe][\mathrm{Mg}/\mathrm{Fe}] on the ISM from which new stars may form, resulting in the elevated [α/Fe][\alpha/\mathrm{Fe}] plateau. This scenario is consistent with the recent multicomponent decomposition model of Yang et al. (2025), who find that the early phase of Sculptor’s evolution was dominated by yields from extremely massive stars.

Although our one-zone framework assumes instantaneous spatial mixing across the galaxy, the elevated [Mg/Fe]cc[\mathrm{Mg}/\mathrm{Fe}]_{\mathrm{cc}} we recover in both dSphs and MQGs is consistent with what an inhomogeneous, clumpy star-formation mode would naturally produce. In the MW, larger gas reservoirs, well-sampled IMFs and gradual star formation would average each CCSN population toward the IMF-averaged yield whereas in dSphs, the smaller gas reservoir (and in the case of Sculptor, relatively-rapid growth in stellar mass) leads to globally under-sampled IMFs, with the most massive progenitors dominating early enrichment, as mentioned earlier.

A similar IMF-undersampling regime may operate in MQG progenitors, not through globally small gas reservoirs, but through localised small reservoirs within individual star-forming clumps. If high-zz star formation proceeds in a fragmented, clumpy mode (as suggested by observations; [168] (Guo et al. 2015; Harikane et al. 2025; Fujimoto et al. 2025; Hainline et al. 2026; Mawatari et al. 2026)), enrichment within each clump would be dominated by the most prompt, most massive progenitors (e.g. hypernovae) producing α\alpha-rich ejecta before significant dilution. The high average SFRs (200 – 800 M⊙M_\odot yr−1^{-1} across our full range of ΔtSF\Delta t_{\mathrm{SF}}) and short SFE timescales we infer (τ⋆∼0.05\tau_\star\sim0.05 Gyr) are consistent with the rapid local gas consumption and prompt formation of stars from the resulting enriched material in this picture. Simulations attempting to reconstruct the MW’s bimodality in [α/Fe][\alpha/\mathrm{Fe}] (Clarke et al. 2019) also support this, finding that clumpy SF with high SFR density naturally produces high-α\alpha abundances. While clumpy SF and inhomogeneous mixing on the galactic scale contradict the instantaneous mixing approximation, our model captures the population-averaged outcome of these sub-zone processes through its empirical Sculptor calibration.

The consequence of adopting yields from Sculptor given the offset in the [α/Fe][\alpha/\mathrm{Fe}]-plateau is that our models are more favourable to producing the elevated [Mg/Fe][\mathrm{Mg}/\mathrm{Fe}] abundances observed in MQGs around cosmic noon. Yet without modifying the plateau in some way, it is impossible to reproduce the observed abundances of the most extremely α\alpha-enhanced ([Mg/Fe]∼0.5−0.6[\mathrm{Mg}/\mathrm{Fe}] \sim0.5 - 0.6) MQGs. Similarly for early-type galaxies at z∼0.7z \sim0.7, Gountanis et al. (2025) find that MW-motivated yields alone are unable to explain the α/Fe\alpha/\mathrm{Fe} abundances, demonstrating that a plateau around [α/Fe]cc∼0.55[\alpha/\mathrm{Fe}]_{\mathrm{cc}} \sim0.55 provides better agreement while still slightly under-predicting the observed abundances. Despite the offset between the MW and dSph plateaus, it remains possible that the true CCSN plateau in the Milky Way is higher than typically assumed. The analysis of metal-poor stars by Conroy et al. (2022) suggests that the intrinsic plateau could be around 0.6 dex. Similarly, Sit et al. (2025) argue that even stars on the observed MW high-α\alpha plateau (0.25≤[Mg/Fe]≤0.350.25 \leq[\mathrm{Mg}/\mathrm{Fe}] \leq0.35) already reflect varying SN Ia contributions rather than pure CCSN enrichment, proposing that the true CCSN plateau lies above this observed range, potentially as high as [Mg/Fe]cc≈0.6[\mathrm{Mg}/\mathrm{Fe}]_{\mathrm{cc}} \approx0.6. These findings are consistent with the plateau found in this work ([Mg/Fe]cc≈0.6[\mathrm{Mg}/\mathrm{Fe}]_{\mathrm{cc}} \approx0.6), lending support to the adopted nucleosynthetic yields and their application to MQGs.

Summary

In this work, we present a one-zone chemical evolution model for the formation of a massive quiescent galaxy (M⋆=1011 M⊙M_\star= 10^{11}\ M_\odot) and compare the predicted chemical abundances to a compilation of observed abundances around cosmic noon (z∼2z \sim2). The adopted nucleosynthetic yields (Mg, Si, Ca, Ti, Fe) are calibrated against a chemical analysis of stars in the Sculptor dSph. We use this model to investigate the star-formation timescales and ISM properties of the progenitors of MQGs, which are sensitive to the chemical enrichment. Our main findings are summarised as follows:

  • The duration of star-formation (ΔtSF\Delta t_{\mathrm{SF}}) directly regulates the resulting stellar [Mg/Fe] abundances, with shorter timescales neces-

sary to limit the enrichment of Fe from SNe Ia after the characteristic delay time tdelayt_{\mathrm{delay}}, consistent with previous studies.

  • We find that rapid molecular gas depletion timescales (SFE timescale; τ⋆∼0.05\tau_\star\sim0.05 Gyr) and negligible mass-loading (η∼0\eta\sim0) are required to reproduce the high [Fe/H] and [Mg/Fe] abundances observed in MQGs at cosmic noon. However, the observed scatter in the [Mg/Fe]-[Fe/H] space may reflect intrinsic variations in the efficiency and strength of outflows across the population, with the shortest τ⋆\tau_\star and lowest η\eta found here representing the extremely α\alpha-enhanced subset rather than the typical MQG.

  • A star-formation timescale of ΔtSF∼0.1−0.4\Delta t_{\mathrm{SF}} \sim0.1 - 0.4 Gyr provides the best agreement in the [Mg/Fe]-[Fe/H] plane for the bulk of MQGs observed at cosmic noon, with the most [α\alpha/Fe]-enhanced galaxies consistent with the shortest timescales (ΔtSF∼0.1\Delta t_{\mathrm{SF}} \sim0.1 Gyr). Additionally, the observed [Si/Fe] and [Ca/Fe] abundances measured from the JWST-SUSPENSE survey (Slob et al. 2024; Beverage et al. 2025) are also consistent with ΔtSF=0.4\Delta t_{\mathrm{SF}} = 0.4 Gyr. However, we are unable to reproduce the abundances of [Ti/Fe] despite motivating its yield using Sculptor, possibly due to physical differences in Ti production between galaxies, or more likely, uncertainty in measuring [Ti/Fe].

  • Our model is capable of reproducing the observed abundances of the centres of local ETGs (R<Re/8R < R_e/8) (Gu et al. 2022), naturally explaining them as quiescent galaxies that formed over longer timescales (ΔtSF∼1.0\Delta t_{\mathrm{SF}} \sim1.0 Gyr). The same chemistry also follows naturally from dry merging of early-forming MQGs with slower-forming, less α\alpha-enhanced galaxies, which would dilute the original extreme α\alpha-enhancement.

  • Our model framework captures the high [Mg/H] space ([Mg/H] ∼0.2−0.4\sim0.2 - 0.4 at [Fe/H] ∼0\sim0) without direct modifications of the IMF, in contrast to previous models in the literature that struggle to reproduce absolute [Mg/H] abundances with Milky Way-motivated yields. Although our adopted yields are calibrated to the Sculptor dSph, which is enhanced in [Mg/Fe] relative to the MW, recent analyses of MW stars suggest the true CCSN [Mg/Fe] plateau may be higher than the canonical value, alleviating this tension.

Collectively, these findings demonstrate that the formation histories of massive quiescent galaxies observed at cosmic noon are rapid. Viewed through the lens of stellar abundances, our analysis indicates that these galaxies formed over very short timescales, with the most extreme cases forming the bulk of their stars in ∼100\sim100 Myr. This rapid assembly has important implications for the progenitors of these galaxies. They likely experienced highly-efficient star formation, with molecular gas depletion timescales comparable to those of local ULIRGs and high-redshift DSFGs, while simultaneously retaining a significant fraction of their metals. Even in the presence of strong stellar feedback, these metals must be rapidly recycled back into the star-forming ISM to reproduce the observed abundances. Furthermore, while extreme hyper-starbursts are observed among some DSFGs, our chemical constraints suggest the progenitors of typical MQGs could not have sustained such extreme SFRs without underproducing the observed Fe enrichment, unless the SN Ia delay time is unexpectedly short (tdelay<50t_{\mathrm{delay}} < 50 Myr) or these bursts are short-lived. Instead, our preferred galaxy progenitors with typical star-formation rates of 200−400 M⊙ yr−1200 - 400\,M_\odot\,\mathrm{yr}^{-1} for the bulk of the population (ΔtSF=0.2−0.4\Delta t_{\mathrm{SF}} = 0.2 - 0.4 Gyr), rising to 800 M⊙ yr−1800\,M_\odot\,\mathrm{yr}^{-1} for the most α\alpha-enhanced galaxies at ΔtSF=0.1\Delta t_{\mathrm{SF}} = 0.1 Gyr. This evolutionary picture makes the typical MQG progenitor physically consistent with recent JWST observations of ‘NIRCam-dark’ galaxies, which have SFRs and depletion timescales that align more closely with our model framework (Bing et al. 2026; Sun et al. 2026).

Finally, our chemical evolution modelling addresses just one piece of the broader puzzle regarding MQG formation and evolution from a chemical perspective. An open question remains regarding the chemical abundance discrepancy between MQGs observed at cosmic noon and their assumed descendants, local ETG cores. As previously highlighted, and discussed by Jafariyazani et al. (2025) and Beverage et al. (2025), passive evolution alone cannot account for this difference, though mechanisms such as obscured star formation or subsequent major merging may alleviate this tension. Our model supports both scenarios: the abundances of local ETG centres are reproduced when MQGs form over ∼1\sim1 Gyr timescales, which dry mergers and continued enriched star formation would equally produce. We note, however, that obscured star formation must proceed from the already-enriched reservoir, since pristine inflows would lower [Fe/H] below the observed values.

We plan to address these broader evolutionary pathways in future work using hydrodynamical simulations to explore the chemistry of MQGs at cosmic noon and trace their subsequent evolution to z=0z = 0.

Acknowledgements

We acknowledge the traditional owners of the land on which the Australian National University stands, the Ngunnawal and Ngambri people, who have surveyed the sky for thousands of years. We pay our respects to elders past, present, and those emerging. PO acknowledges support from the Australian Government Research Training Program (RTP) Scholarship. SB acknowledges support from the Australian Research Council under grant number DE240100150. MKN acknowledges funding from the ARC for project DP260102813. Refinements to parts of the text and code were developed through highly-supervised engagement with AI large language models (Gemini by Google and Claude by Anthropic).

Data availability

The chemical evolution models presented in this work are computed using the publicly available VICE code (Johnson & Weinberg 2020). The observational abundance data used for comparison are compiled from previously published, publicly-available sources, as referenced throughout this work. Post-processed data underlying results of this work can be requested from the lead author.

References

  1. [82]Lagos C. d. P., et al., 2018, MNRAS, 473, 4956
  2. [83]Lagos C. d. P., et al., 2025, MNRAS, 536, 2324
  3. [89]Mallinson J. W. E., Lind K., Amarsi A. M., Barklem P. S., Grumer J., Belyaev A. K., Youakim K., 2022, A&A, 668, A103
  4. [120]Rodríguez Del Pino B., et al., 2026, A&A, 711, A66
  5. [158]Xu X., et al., 2022, ApJ, 933, 222
  6. [164]de Boer T. J. L., et al., 2012, A&A, 539, A103
  7. [166]de La Rosa I. G., La Barbera F., Ferreras I., de Carvalho R. R., 2011, MNRAS, 418, L74
  8. [167]de Zeeuw P. T., et al., 2002, MNRAS, 329, 513
  9. [168]de la Vega A., Mobasher B., Sattari Z., Chartab N., Manesh F., Sharei N., 2026, ApJ, 1005, 169
  10. [170]van Dokkum P. G., Conroy C., 2012, ApJ, 760, 70
  11. [171]van Dokkum P. G., et al., 2010, ApJ, 709, 1018

Paper details

vlnc:W_q3pfrxqv
v1 · latestHistory (1)

published · Fresh — added to the Hub within the last 30 days

added ·

Usage

views · 0

downloads · 0

A view counts once per browser per day, a download each time the PDF is requested.

CC BY 4.0

Original source

  • Physics and Astronomy

Gemini (multiple versions), Claude by Anthropic · View statement →

Looking for related papers…
Contents