# Megatron: how the first stars create an iron metallicity plateau in the smallest dwarf galaxies

Martin Rey, Harley Katz, Corentin Cadiou, Mahsa Sanati, Oscar Agertz, Jeremy Blaizot, Alex Cameron, Nicholas Choustikov, Julien Devriendt, Uliana Hauk, et al.

► To cite this version:

Martin Rey, Harley Katz, Corentin Cadiou, Mahsa Sanati, Oscar Agertz, et al.. Megatron: how the first stars create an iron metallicity plateau in the smallest dwarf galaxies. 2025. hal-05352224

HAL Id: hal-05352224

https://hal.science/hal-05352224v1

Preprint submitted on 6 Nov 2025

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
MEGATRON: HOW THE FIRST STARS CREATE AN IRON METALLICITY PLATEAU
IN THE SMALLEST DWARF GALAXIES

MARTIN P. REY¹*, HARLEY KATZ²,³, CORENTIN CADIOU⁴,

MAHSA SANATI⁵, OSCAR AGERTZ⁶, JEREMY BLAIZOT⁷, ALEX J. CAMERON⁵, NICHOLAS CHOUSTIKOV⁵, JULIEN DEVRIENDT⁵, ULIANA HAUK², ALEXANDER P. JI²,³, GARETH C. JONES⁸,⁹, TAYSUN KIMM¹⁰, ISAAC LASETER¹¹, SERGIO

MARTIN-ALVAREZ¹², KOSEI MATSUMOTO¹³, AUTUMN PEARCE², YVES REVAZ¹⁴, FRANCISCO RODRÍGUEZ MONTERO²,³, JOKI

ROSDAHL⁷, AAYUSH SAXENA⁵, ADRIANNE SLYZ⁵, RICHARD STISKALEK⁵, ANATOLE STORCK⁵, OSCAR VEENEMA⁵, WONJAE

¹University of Bath, Department of Physics, Claverton Down, Bath, BA2 7AY, UK

²Department of Astronomy & Astrophysics, University of Chicago, 5640 S Ellis Avenue, Chicago, IL 60637, USA

³Kavli Institute for Cosmological Physics, University of Chicago, Chicago IL 60637, USA

⁴Institut d'Astrophysique de Paris, Sorbonne Universités, CNRS, UMR 7095, 98 bis bd Arago, 75014 Paris, France

⁵Sub-department of Astrophysics, University of Oxford, Keble Road, Oxford OX1 3RH, United Kingdom

⁶Lund Observatory, Division of Astrophysics, Department of Physics, Lund University, Box 43, SE-221 00 Lund, Sweden

⁷Université Claude Bernard Lyon 1, CRAL UMR5574, ENS de Lyon, CNRS, Villeurbanne, F-69622, France

⁸Kavli Institute for Cosmology, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK

⁹Cavendish Laboratory, University of Cambridge, 19 JJ Thomson Avenue, Cambridge CB3 0HE, UK

¹⁰Department of Astronomy, Yonsei University, 50 Yonsei-ro, Seodaemun-gu, Seoul 03722, Republic of Korea

¹¹Department of Astronomy, University of Wisconsin-Madison, Madison, WI 53706, USA

¹²Kavli Institute for Particle Astrophysics & Cosmology (KIPAC), Stanford University, Stanford, CA 94305, USA

¹³Sterrenkundig Observatorium Department of Physics and Astronomy Universiteit Gent, Krijgslaan 281 S9, B-9000 Gent, Belgium and

¹⁴Institute of Physics, Laboratoire d'Astrophysique, École Polytechnique Fédérale de Lausanne (EPFL), CH-1015 Lausanne, Switzerland

Version October 8, 2025

ABSTRACT

We study the stellar mass–iron metallicity relation of dwarf galaxies in the new high-resolution MEGATRON cosmological radiation-hydrodynamics simulations. These simulations model galaxy formation up to $z \approx 8$ in a region that will collapse into a Milky-Way-like galaxy at $z = 0$, while self-consistently tracking Population III and II (Pop. III, Pop. II) star formation, feedback and chemical enrichment. MEGATRON dwarf galaxies are in excellent agreement with the observed stellar mass–metallicity relation at $z = 0$, including an over-abundance of dwarfs along a flat plateau in metallicity ($\langle[\text{Fe/H}]\rangle \approx -2.5$) at low stellar masses ($M_* \le 10^5 M_\odot$). We tie this feature to the chemical enrichment of dwarf galaxies by Pop. III pair-instability supernova (PISN) explosions. The strong Lyman-Werner background (LW) from the protogalaxy ensures that PISNs occur in haloes massive enough ($\approx 10^7 M_\odot$) to retain their ejecta. We also predict a tail of $\approx 20\%$ of iron-deficient ($\langle[\text{Fe/H}]\rangle \le -3$) dwarf galaxies. We show that both plateau and tail (i) are robust to large variations in Pop. II feedback assumptions, and (ii) survive in bound satellites surrounding the central galaxy at $z = 0$.

Subject headings: galaxies: dwarf, abundances, formation – stars: Population III – methods: numerical

1. INTRODUCTION

Characterizing how the first generation of stars enriched the pristine Universe with chemical elements remains a key question of modern astrophysics. The advent of the *James Webb Space Telescope* (JWST) now allows us to observe this process in-situ. With its new infrared spectroscopic capabilities, we can observe emission lines originating in the interstellar medium (ISM) of high-redshift galaxies. This allows us to, for example, constrain the relationship between stellar mass and gas-phase oxygen abundance at $z \ge 6$ (e.g. Curti et al. 2023; Heintz et al. 2023; Curti et al. 2024; Laseter et al. 2024; Sanders et al. 2024; Li et al. 2025; Scholte et al. 2025), as well as the relative abundances between individual chemical elements (e.g. Cameron et al. 2023; Isobe et al. 2023; Stiavelli et al. 2023; Hsiao et al. 2024; Schaerer et al. 2024; Stanton et al. 2025; Topping et al. 2024, 2025; Nakane et al. 2025).

The same star formation episodes that spawn the short-lived ($\approx 10$ Myr) massive stars which power the

observed emission lines at $z \ge 6$ also create long-lived ($\ge 10$ Gyr) low-mass stars that survive to the present day. The chemistry of these metal-poor stars keeps a record of the gas conditions in which they formed, providing us with a complementary view of the chemical enrichment history of the Universe at high redshift. We have long been able to identify and characterize these remnant metal-poor stars in the Milky Way and nearby satellite galaxies (e.g. Tolstoy, Hill & Tosi 2009; Frebel & Norris 2015; Bonifacio et al. 2025 for reviews). In turn, the chemistry of $z = 0$ metal-poor stars has allowed us to, for example, constrain the mass function of primordial Population III (Pop. III) stars (e.g. Hartwig et al. 2015; Ishigaki et al. 2018; Jiang et al. 2024) or the yields and explosion properties of metal-poor massive stars (e.g. Kobayashi et al. 2014; Koutsouridou et al. 2023).

However, the connection between the observed spectroscopic properties of $z > 6$ galaxies and the local metal-poor stars is challenging. Not only because of the $\ge 13$ billion years gap in cosmic evolution between the two regimes, but also because the elements accessible to galaxy and stellar spectroscopy tend to differ

<footnote>*E-mail: mpr47@bath.ac.uk</footnote>
(e.g. oxygen in the former, iron in the latter). Furthermore, metal-poor stars can only be resolved in the immediate vicinity of the Milky Way, providing us with an important, but biased, view of low-metallicity chemical enrichment in one specific environment. How to generalize Milky-Way data to the high-redshift Universe and vice-versa remains an open question.

In this paper, we introduce the MEGATRON project, which aims to clarify this connection between high-redshift emission lines and low-redshift metal-poor stars.

This paper initiates this effort by focusing on the properties of ultra-faint dwarf galaxies (UFDs; $M_V \ge -6$; $M_* \le 10^5 M_\odot$ at $z=0$). These galaxies populate the faintest end of the galaxy luminosity function and are the most iron-poor galaxies known at $z=0$ ($\langle[\text{Fe/H}]\rangle \le -2$; see Simon 2019 for a review). Their star formation histories almost invariably favour a high-redshift formation, with limited or non-existent star formation after $z \le 4$ (e.g. Okamoto et al. 2012; Brown et al. 2014; Weisz et al. 2014; Savino et al. 2025; Durbin et al. 2025). This truncation in star formation is usually attributed to the low dynamical masses of UFDs. Once the cosmic ultraviolet (UV) background following cosmic reionization has heated the intergalactic medium (IGM), their gravitational potential wells are too shallow to accrete fresh gas from the now-warm IGM (Efstathiou 1992; Shapiro, Giroux & Babul 1994; Gnedin 2000; Noh & McQuinn 2014), terminating star formation activity (e.g. Bullock, Kravtsov & Weinberg 2000; Benson et al. 2002; Somerville 2002). As a result, UFDs are ‘fossil relics’ of the high-redshift Universe (e.g. discussions in Ricotti & Gnedin 2005; Bovill & Ricotti 2009; Salvadori & Ferrara 2009; Bland-Hawthorn, Sutherland & Webster 2015), with their lack of star formation after $z \approx 3^1$ making them ideal to isolate and constrain metal-free and metal-poor stellar evolution with observed abundances (e.g. Magg et al. 2018; Rossi, Salvadori & Skúladóttir 2021; Rossi et al. 2023, 2025; Chiti et al. 2025).

In particular, the iron metallicity of UFDs is a sensitive probe of high-redshift star formation and feedback conditions. Observed UFDs around the Milky Way show a break in slope in their mass-metallicity relation<fnref n="2" />. All UFDs with $M_* \le 10^5 M_\odot$ cluster around $\langle[\text{Fe/H}]\rangle \approx -2.5$, in a flat ‘plateau’ as stellar masses decrease to $M_* \le 10^2 M_\odot$ (e.g. Simon 2019; Fu et al. 2023). This plateau has been a challenge for recent models of UFDs, which struggle to reproduce the comparatively-high iron content of UFDs given their very low stellar masses (e.g. Munshi et al. 2019; Wheeler et al. 2019; Applebaum et al. 2021; Sanati et al. 2023; Go et al. 2025). The origin of this plateau is debated, with studies showing that it could stem from modifications of the Pop. II IMF (Prgomet et al. 2022), Pop. III physics (e.g. Jeon, Besla & Bromm 2017), a threshold metallicity in the intergalactic medium set by external enrichment (e.g. Jeon, Besla & Bromm 2017; Ahvazi et al. 2024, though see Wheeler et al. 2025), or the details of how galactic

outflows transport iron out of the ISM of dwarf galaxies (e.g. Agertz et al. 2020; Rey et al. 2025).

In this paper, we show that a plateau in iron metallicity is naturally explained from internal chemical enrichment following the birth of high-mass Pop. III stars. To demonstrate this, we leverage the population of faint dwarf galaxies formed in the new suite of high-redshift MEGATRON simulations (Katz et al. 2025).

The novelty of the MEGATRON simulations is multifold:

* We solve for both radiative transfer and non-equilibrium chemistry of $\ge 80$ primordial species, molecules, and metal ions. This allows us to make ab-initio, robust predictions of emission lines and absorption lines (Katz et al. 2025; Cadiou et al. 2025), in turn allowing a direct quantification of how these observables respond to changes in feedback modelling (e.g. N. Choustikov et al. in prep).

* We model small enough spatial and mass scales to resolve most Pop. III star formation sites in the smallest first galaxies (dark matter particle mass $m_{DM} = 2.5 \times 10^4 M_\odot$, spatial resolution $\Delta x \approx 3$ pc), allowing us to quantify their birth properties and observables (A. Storck et al. in prep).

* The simulation volume is chosen to be a protosolar Milky-Way-like environment, with the final progenitor reaching a halo mass $\approx 10^{12} M_\odot$ at $z=0$. These initial conditions form a large gas and stellar disc ($M_* \approx 5 \times 10^{10} M_\odot$) when evolved to $z=0$ with well-tested low-redshift models (Rey et al. 2023; Joshi et al. 2025). Leveraging this last aspect to bridge high- and low-redshift observables of early chemical enrichment is the focus of this paper.

We present the MEGATRON simulations and our galaxy formation model in Section 2. We then showcase the mass-metallicity relation of the population of simulated MEGATRON dwarf galaxies ($\ge 500$) in Section 3, highlighting the emergence of a plateau around $\langle[\text{Fe/H}]\rangle \approx -2.5$ at low stellar masses ($M_* \le 10^5 M_\odot$) in excellent agreement with local Universe observations. We attribute the origin of this plateau to Pop. III physics in Section 4 and justify the connection between high-redshift simulations and low-redshift data in Section 5. We conclude in Section 6.

## 2. THE MEGATRON HIGH-REDSHIFT SUITE

We use the high-redshift MEGATRON simulations that are evolved to $z \approx 8$ (Katz et al. 2025). This suite consists of four simulations of the same initial conditions (ICs), each with a different Pop. II star formation and feedback model (but identical Pop. III assumptions). We describe the ICs and their link to the Milky Way in Section 2.1, and summarize the MEGATRON galaxy formation model in Section 2.2 (see Katz et al. 2024, 2025 for a more in-depth description).

### 2.1. *Proto-Milky Way initial conditions*

The MEGATRON suite is based on the ICs of the VINTERGATAN-GM suite of hydrodynamical simulations (Rey et al. 2023; Joshi et al. 2025). These initial conditions are selected to sample a range of Milky-Way-mass

<footnote n="1">¹ Reionization quenching is not instantaneous, and residual star formation can persist after reionization until galactic outflows have expelled self-shielded gas within haloes (e.g. Susa & Umemura 2004; Oñorbe et al. 2015; Rey et al. 2020)</footnote>

<footnote n="2">² Hereafter, the mass-metallicity relation refers to the relationship between $M_*$ and stellar-based metallicities $\langle[\text{Fe/H}]\rangle$.</footnote>
haloes ($M_{200}(z=0) \approx 10^{12} M_{\odot}$), where $M_{200}$ defines the mass enclosed within the radius $r_{200}$ where the density equals 200 times the critical density of the Universe.

Out of the library of VINTERGATAN-GM ICs, MEGATRON uses an early-forming Milky-Way-mass halo (introduced as ‘Halo 599’ in Rey & Starkenburg 2022). The mass accretion history is then modified using the ‘geNETIC modification’ technique (Roth, Pontzen & Peiris 2016; Rey & Pontzen 2018; Stopyra et al. 2021). First, we genetically modify the ICs to zero the velocity of the Lagrangian patch with respect to the grid. This enhances the precision of numerical integration with a grid code by reducing advection errors (Pontzen et al. 2021). Then, we engineer two genetically-modified variants to have more or less high-redshift substructures using quadratic modifications (Rey & Pontzen 2018), but always maintaining a very similar mass growth between $z=2$ and $z=0$ (see Katz et al. 2025, fig. 1).

To benchmark these ICs, the reference, unmodified formation history is evolved with two galaxy formation models that are well tested at $z=0$: the VINTERGATAN model (Agertz et al. 2021) and the ILLISTRISTNG model (Pillepich et al. 2018). With both models, the MEGATRON ICs lead a large disc galaxy at $z=0$ ($M_{\star}=3-7 \times 10^{10} M_{\odot}$) hosted in a Milky-Way-mass halo (Joshi et al. 2025, fig.2, left column for a visual). These low-redshift runs serve two purposes: (i) they provide us with well-defined benchmarks on which to compare our high-redshift galaxy formation efficiencies for which data is currently sparse (see discussion in Katz et al. 2025); and (ii) they strongly establish that MEGATRON is studying high-redshift galaxy formation in an environment that represents a reasonable realization of the proto-Milky Way.

## 2.2. The high-redshift galaxy formation model

Each MEGATRON simulation is a zoomed, cosmological simulation with $m_{\text{DM}} = 2.5 \times 10^4 M_{\odot}$. The spatial resolution evolves with redshift, from $\Delta x = 2.5$ pc at $z = 25$ to $\Delta x \approx 5$ pc at $z = 8$ (see Katz et al. 2025 for more details).

The simulations are evolved with the adaptive mesh refinement RAMSES-RTZ code (Teyssier 2002; Rosdahl et al. 2013; Katz 2022). We solve for gravity, hydrodynamics, radiative transfer, non-equilibrium primordial chemistry (including H$_2$), and non-equilibrium metal chemistry of $\ge 70$ heavy elements and molecules. Non-equilibrium metal chemistry and radiative transfer are coupled through the PRISM model (Katz et al. 2022; see Katz et al. 2024 for updates).

MEGATRON employs a detailed galaxy formation model including the treatment of Pop. III star formation, supernova feedback, radiative feedback, and metal enrichment from individual stars. Briefly, Pop. III star formation is modelled as in Kimm et al. (2017); Katz et al. (2023), spawning individual Pop. III stars following a lognormal IMF with a characteristic mass $100 M_{\odot}$ (Wise et al. 2012). Pop. III stars with masses $10 < m_{\star} \le 20 M_{\odot}$ explode as core-collapse supernovae (CCSNe), those with masses $20 < m_{\star} \le 40 M_{\odot}$ explode as hypernovae (HNe), and those with masses $140 \le m_{\star} \le 300 M_{\odot}$ explode as PISNe. Each of these stellar types has a different explosion energy, namely $10^{51}$, $\approx 2 \times 10^{52}$, $\approx 5 \times 10^{52}$ erg for CCSNe, HNe and PISNe, respectively (see Kimm et al.

2017, section 2.3.3 for the exact scaling with progenitor mass or core mass). Element-by-element Pop. III yields are injected according to Nomoto, Kobayashi & Tominaga (2013). Radiative feedback from Pop. III stars follow the models of Schaerer (2002). Pop. III stars whose progenitor masses are outside the mass range mentioned are assumed to collapse directly as a black hole (no explosion, no injection of chemical elements). More details of the implementation can be found in A. Storck et al. in prep, with their Fig. 1 providing a visualization of the IMF and stellar evolution assumptions.

Once gas is above $Z \ge 10^{-6} Z_{\odot}$ where $Z_{\odot} = 0.02$, we transition to a Pop. II star formation mode. The star formation algorithm remains the same. We use a multi-freefall approach based on a local gravo-thermoturbulent condition (Kimm et al. 2017), and stellar particles represent a simple stellar population rather than individual stars. We assume a fiducial Pop. II IMF from Kroupa (2001) and a minimum stellar particle mass of $500 M_{\odot}$. Pop. II populations eject stellar winds, undergo supernova explosions, and inject radiative feedback and metal enrichment following the approach in Agertz et al. (2020). Briefly, supernova feedback is injected as energy ($10^{51}$ erg) if the cooling radius is resolved, and as momentum otherwise ($3 \times 10^5 M_{\odot} \text{ km s}^{-1}$). Radiative feedback is injected following the BPASS v2.2.1 binary stellar evolution model (Stanway, Eldridge & Becker 2016). CCSNe inject chemical elements following the yields of Limongi & Chieffi (2018) assuming no rotation. We choose this Pop. II feedback budget as it is well tested at $z=0$ in both Milky-Way galaxies (e.g. Agertz et al. 2021) and dwarf galaxies (e.g. Rey et al. 2025).

More details about the implementation of MEGATRON can be found in Katz et al. 2024, 2025.

## 2.3. The suite of four simulations

The MEGATRON high-redshift suite consists of 4 simulations of the same initial conditions varying assumptions in Pop. II star formation and feedback modelling ('Efficient SF', 'Bursty SF', 'Varying IMF' and 'HN, $\epsilon_{eff}$').

Briefly, the first two models differ in the strength of Pop. II supernova feedback ($10^{51}$ erg for 'Efficient SF' vs. $5 \times 10^{51}$ erg for 'Bursty SF', respectively), leading to efficient gas conversion into stars in the first case and burstier star formation histories in the second case. The 'Varying IMF' implements a (primarily) density-dependent Pop. II IMF from Marks et al. (2012), which becomes more top-heavy in higher-density star formation environments and converges to our fiducial Kroupa (2001) otherwise. In addition to a varying-IMF, this run assumes that low-metallicity Pop. II stars can undergo HN explosions, with a varying energy explosion depending on mass (Nomoto et al. 2006) and a fraction of HNe that strongly depends on metallicity (Kobayashi et al. 2006). The last model changes parameters of the star formation algorithm to ensure that star formation occurs is more efficient (the star formation efficiency $\epsilon_{ff}$ is fixed at 100% rather than varying with Mach number and turbulence) and more clustered (stellar particle mass of $2,000 M_{\odot}$ instead of our fiducial $500 M_{\odot}$). This last model also implements the same HN feedback for Pop. II low-metallicity stars as the 'Varying IMF'.

We refer the reader to Katz et al. (2025) for the exact
<figure><img src="image_2.png" /><figcaption>FIG. 1.— MEGATRON cosmological radiation-hydrodynamics simulations (left panels, gas density, gas temperature, gas iron metallicity, LW background at $z \approx 8$) model galaxy formation of a proto-galaxy that will turn into a Milky-Way-like galaxy at $z = 0$ (right panels). High-redshift dwarf galaxies in MEGATRON are mapped onto a dark matter only simulation (right, top panel, orange circles highlighting dwarf galaxies with $M_{\star} \le 10^5 M_{\odot}$) to predict their dynamical future around the Milky-Way-like galaxy at $z = 0$ (e.g. first infall, disrupted in the stellar halo, bound in orbiting satellites, bottom right panel, see Section 5 for details of the particle tagging technique).</figcaption></figure>

parameter choices of each model and to Choustikov et al. (2025) for how they affect the ISM structure at high redshift. In the context of this paper, comparing between these four simulations simply allows us to verify the robustness of our findings when making large variations to Pop. II feedback implementation.

## 2.4. Halo finding and dynamical future of high-redshift dwarfs

The MEGATRON simulations are only evolved to $z \approx 8$, but we wish to connect the properties of their simulated dwarfs to local data around the Milky Way at $z = 0$.

To perform this connection, we leverage the dark-matter-only counterpart of the MEGATRON ICs, which we evolve to $z = 0$ with the RAMSES code (Teyssier 2002) under gravity and the same $m_{DM}$. Even though we have access to hydrodynamical runs of these initial conditions at $z = 0$ (see Section 2.1), these runs are evolved with 8x coarser dark matter resolution than the fiducial MEGATRON simulations ($m_{DM} \approx 2 \times 10^5 M_{\odot}$ rather than $m_{DM} = 2.5 \times 10^4 M_{\odot}$). As a result, the dark-matter-only run with matching resolution is better suited to track the fate of faint objects.

In all simulations, we identify dark matter haloes and subhaloes with the ROCKSTAR halo finder (Behroozi, Wechsler & Wu 2013). We then match galaxies in the hydrodynamical MEGATRON simulations to their dark-matter-only counterparts at $z \approx 8$ (Figure 1, upper right panel). The IDs of dark matter particles are generated self-consistently across all MEGATRON simulations, allowing us to cross-match the haloes with the most particles in common between simulations. We discard haloes if $M_{200}$ differ by more than a factor two between the hydrodynamical and dark-matter-only runs.

Once a halo is matched between the hydrodynamical and dark-matter-only runs, we use the particle-tagging approach from Rey & Starkenburg (2022) to estimate its dynamical fate at $z = 0$. We identify the 5% most bound dark matter particles of the host halo at $z \approx 8$ and ‘tag’ them as dynamical proxies for the future evolution of the stars (Rey & Starkenburg 2022; see also Bullock & Johnston 2005; De Lucia & Helmi 2008; Cooper et al. 2010, 2017 for similar approaches).

Tagged particles are then tracked to $z = 0$, allowing us to predict whether the dwarf galaxy is disrupted in the stellar halo of the Milky Way, is on first infall, or is bound in orbiting satellites (Figure 1, lower right panel). The dwarf is labelled as a surviving $z = 0$ structure if at least 50% of the tagged particles belong to a bound structure identified by the ROCKSTAR halo finder. We verified that using 2% of the most bound particles does not qualitatively change our results (see also discussion in Le Bret et al. 2017; Rey & Starkenburg 2022).

Another technicality in the comparison between high and low redshift is that stellar evolution will induce mass-loss over the Hubble time, reducing stellar masses between $z \approx 8$ and $z = 0$. To estimate the importance of this effect, we compute the projected stellar mass at $z = 0$ for each dwarf galaxy by counting which stars would still be on the main sequence at $z = 0$ given our stellar evolution model (mass-loss rate for low-mass stars and explodability range for massive stars; see Katz et al. 2025). Using the projected stellar mass at $z = 0$ or the simulated one at $z \approx 8$ does not impact our conclusions (see Section 5.1).

## 3. THE STELLAR MASS-METALLICITY RELATION OF DWARF GALAXIES
<figure><img src="image_3.png" /><figcaption>FIG. 2.— The stellar mass-mass metallicity relation of dwarf galaxies, observed at $z = 0$ (grey points, Pace 2025), observed in $z = 3$ star-forming galaxies (Stanton et al. 2024), and simulated MEGATRON galaxies at $z \approx 8$ (color points). The left panel shows all individual simulated dwarfs, while the right panel shows the mode in bins of $M_\star$. The agreement in slope, normalization, and scatter over $10^5 \le M_\star \le 10^8 M_\odot$ is striking. Only one simulation ('Efficient SF', purple) overshoots the data, due to inefficient galactic outflows in higher-$M_\star$ galaxies. All MEGATRON simulations plateau in iron metallicities at low masses ($M_\star \le 10^5 M_\odot$), with an extended tail of iron-poor galaxies (Figure 3).</figcaption></figure>

Figure 2 shows the average iron content in the stars of MEGATRON galaxies at $z \approx 8$ as a function of their stellar masses. We select only central galaxies (i.e. excluding galaxies that are already satellites at $z \approx 8$), and galaxies that contain at least 1 Pop. II star particle (excluding purely Pop. III objects which are still forming at $z \approx 8$; A. Storck et al. in prep, fig.2). We compute $\langle [\text{Fe/H}] \rangle$ as the mass-weighted average over stellar particles (see Escala et al. 2018, eq. 3 and 4) and verified that using the median does not impact our results.

To obtain $M_\star$, we sum the mass of all Pop. II stellar particles within $r_{200}$, removing the mass contributed by non-luminous stellar evolution remnants (e.g. black holes). $\langle [\text{Fe/H}] \rangle$ is computed by averaging the metallicity of all star particles within $r_{200}$ weighted by the current mass of each star particle. We use solar abundance patterns from Asplund et al. (2009) to calculate $\langle [\text{Fe/H}] \rangle$.

The left panel shows all individual simulated dwarf galaxies per simulation, while the right panel shows the mode of the $\langle [\text{Fe/H}] \rangle$ distribution in bins of $M_\star$ for each simulation (where individual simulations have been offset from the centre of the bin for clarity). We use the mode, rather than the median, to highlight where most dwarfs cluster. (We show the full shape of the distribution at the faint end in Figure 3). The grey points show the observed stellar mass-metallicity in local dwarf galaxies taken from the database of Pace (2025) (see Appendix A for individual citations). The black line shows the best-fit to the stellar mass-mass metallicity from $z = 3$ star-forming galaxies from Cullen et al. (2021); Stanton et al. (2024).

The agreement between the simulated and observed

stellar mass-mass relation is striking. For all simulations, the population of simulated dwarf galaxies overlap with the $z = 0$ data, particularly in the range $10^4 \le M_\star \le 10^6 M_\odot$ where they match the normalization and slope of the observed data. This is despite the fact that the simulated dwarf galaxies are formed at $z \approx 8$, while observed dwarf galaxies are at $z = 0$ or $z = 3$. This confirms the view that the stellar mass-metallicity relation of small dwarf galaxies is already set at high redshift and does not significantly evolve thereafter (see Section 5 for further discussion).

Towards the high-mass end ($M_\star \ge 10^6 M_\odot$), the 'Efficient SF' simulation (purple) overshoots the relation. The 'Efficient SF' simulation drives the weakest galactic outflows (Katz et al. 2025), allowing a higher iron retention in the ISM and iron over-abundance in the stars. The other three simulations, with much more efficient outflows, are in close agreement to the observed data. These results re-affirm that the slope and normalization of the stellar-mass-metallicity relation of dwarf galaxies is particularly sensitive to the strength of galactic outflows (see also Agertz et al. 2020; Rey et al. 2025), for the first time demonstrated over a large population of simulated dwarf galaxies ($N \ge 200$ per simulation) rather than a handful of objects.

At the low-mass end of the relation, the observed data in the Local Volume show a break in slope, turning into a 'plateau' of $\langle [\text{Fe/H}] \rangle \approx -2.5$ for $M_\star \le 10^5 M_\odot$ and a dearth of iron-poor galaxies at low-$M_\star$. All simulations show a similar break in slope in the mode, with the majority of simulated dwarf galaxies clustering around
<figure><img src="image_4.png" /><figcaption>FIG. 3.— Same as Figure 2 but zooming on the low-mass end of the stellar mass-metallicity relation. In all simulations, the distribution of faint dwarf galaxies is strongly peaked around $\langle[\text{Fe/H}]\rangle \approx -2.5$ over a broad range of $M_*$, reproducing the observed ‘plateau’ of local dwarf galaxies (grey histogram and data points), and a long tail of iron-deficient galaxies. Both features emerge from our choices of Pop. III star formation and stellar evolution (Figure 4 and Figure 5).</figcaption></figure>

the same plateau as the data. However, the simulations also show a more continuous and scattered distribution of iron metallicities and a tail extending below the plateau. To visualize this distribution better, Figure 3 shows the stellar mass-metallicity zoomed in on the low-mass end of the relation with aggregated histograms (right panel). Note that the lack of very-low-mass dwarf galaxies in the ‘HN, $\epsilon_{eff}$’ simulation is due to the increased stellar particle mass in this run only allowing us to start resolving Pop. II galaxies with $M_* \ge 10^{3.5} M_\odot$.

From Figure 3, we quantitatively confirm that 78% of simulated faint dwarf galaxies cluster between $-3 \le \langle[\text{Fe/H}]\rangle \approx -2.0$ (hereafter referred to as the ‘iron plateau’), similar to the observed ultra-faint dwarfs in the Local Group (grey points and histogram). Unlike the data, however, 22% of the dwarf population populates a long tail of iron-poor galaxies below $\langle[\text{Fe/H}]\rangle \le -3.0$ with a handful of objects having higher-than-plateau iron metallicities. We discuss the observability of this tail around the Milky Way and its implications in Section 6.

All simulations showcase similar distributions of dwarf galaxies, despite large variations in Pop. II assumptions within them. This hints that the iron distribution at the faintest end is rather set by Pop. III modelling which stays the same across all four simulations. We explore this further in the next section.

# 4. THE ORIGIN OF THE IRON PLATEAU

## 4.1. *The importance of Pop. III microphysics*

To understand the origin of the distribution of $\langle[\text{Fe/H}]\rangle$, we break the sample of faint dwarf galaxies into three populations: above, on, and below the plateau. We define the boundaries of the plateau as $-3.0 \le \langle[\text{Fe/H}]\rangle \le -2.0$ (horizontal lines; Figure 3). The lower boundary is somewhat arbitrary, since the simulations show a continuous tail of iron-poor dwarfs,

<figure><img src="image_5.png" /><figcaption>FIG. 4.— The number of Pop. III explosions within MEGATRON dwarf galaxies broken down by their position on, above, or below the plateau (red, brown, and gold, respectively). Almost all faint dwarf galaxies undergo exactly one Pop. III explosion. A small fraction (2% of the total population) of dwarf galaxies are externally enriched by diluted ejecta from nearby objects and have had no Pop. III explosions. Those galaxies are exclusively found below the plateau but cannot explain the whole tail of iron-deficient dwarf galaxies.</figcaption></figure>

but we check that all results are qualitatively unchanged if we instead define the lower limit of the plateau as $\langle[\text{Fe/H}]\rangle \ge -2.7$ or $\langle[\text{Fe/H}]\rangle \ge -3.2$.

Figure 4 then shows the number of exploding Pop. III stars in each subpopulation of dwarf galaxies stacked across all four simulations. Irrespective of whether they are on, above, or below the plateau, the majority of dwarf galaxies undergo precisely one Pop. III explosion. As we will see in Section 4.3, chemical elements are being retained in the ISM of our Pop. III formation sites, despite the large energies associated with Pop. III explosions. This single explosion is thus enough to enrich the minihalo of the nascent galaxy above our critical metallicity threshold ($Z \ge 10^{-6}Z_\odot$) and transition to Pop. II star formation.

A small fraction (6%) of dwarf galaxies have transitioned to Pop. II star formation with multiple Pop. III explosions. On close inspection, we find that this arises roughly half of the time when Pop. III stars are rapidly spawned in a small cluster with an age spread $\le 1$ Myr, leading to nearly simultaneous explosions and enrichment. The other half of the time, the age spreads are much larger ($\ge 5$ Myr) which we hypothesize can occur either when all metals are vacated from the minihalo after the first Pop. III explosion (requiring a second explosion before transitioning to Pop. II star formation), or when Pop. III form in two independent minihaloes that merge together into the faint dwarf we are studying at $z=8$.

An even smaller fraction (2%) of the dwarf galaxy population have transitioned to Pop. II star formation without any Pop. III explosions. These exclusively sit
<figure><img src="image_6.png" /><figcaption>FIG. 5.— Same as Figure 4, but breaking down each sub-population of dwarfs by the exact progenitor of Pop. III explosions. MEGATRON galaxies experience almost no Pop. III CCSNe and HNe explosions due to the choice of a top-heavy Pop. III IMF (first and second panels). But dwarf galaxies on the metallicity plateau are mostly enriched by high-mass PISNe (fourth panel), while galaxies below the plateau have a significant contribution from low-mass PISNe (third panel). This split is naturally explained by the vastly different iron yields in explosions of low-mass and high-mass PISNe (see text for discussions).</figcaption></figure>

below the plateau of the mass-metallicity relation (red line) and have likely been enriched externally by diluted ejecta from nearby star formation. However, the small numbers of such externally-enriched dwarfs cannot explain the whole tail of iron-deficient dwarfs noted in Figure 3. Furthermore, we find no difference in the spatial distribution or chemical abundances of externally-enriched dwarfs compared to the rest of the population (Appendix B). Given their low occurrence rate which will be further reduced by tidal disruption and stripping by $z=0$, this supports the picture that the chemical content of the Milky-Way's ultra-faint dwarfs is dominated by internal Pop. III or Pop. II events, rather than external contamination (see also discussion in Griffen et al. 2018; Wheeler et al. 2025).

Next, we break the populations by each type of Pop. III explosions in Figure 5. Beyond CCSNe (first panel), HNe (second panel), we divide PISNe into low-mass progenitors ($140 \le m_* \le 160 M_{\odot}$; third panel) and high-mass progenitors ($160 < m_* \le 300 M_{\odot}$; fourth panel). The precise value for this split is motivated by a strong evolution in carbon and iron production below and above this value (Nomoto, Kobayashi & Tominaga 2013), although the precise value is somewhat arbitrary since we sample a continuum of progenitor masses and (linearly) interpolate between progenitor mass bins of the yield tables.

First, almost no MEGATRON dwarf galaxies undergo Pop. III CCSNe and HNe (first and second panels). This is due to our choice of a top-heavy Pop. III IMF, favouring the formation of more massive progenitors (see A. Storck et al. in prep, fig. 1). As a result, the vast majority of explosions experienced by our dwarf galaxies are PISNe.

Dwarf galaxies that sit on the plateau of the mass-metallicity relation (brown line) almost exclusively host high-mass PISN explosions (fourth panel). In contrast, dwarf galaxies below the plateau (gold curve) show a mix of explosions and have a significant contribution from low-mass PISNe. This split is naturally explained when considering the iron yield of each of these progenitors. According to Nomoto, Kobayashi & Tominaga (2013), an explosion with progenitor mass $m_* = 140 M_{\odot}$ yields $\approx 0.5 M_{\odot}$ of iron, while a progenitor with $m_* = 180 M_{\odot}$

yields $\approx 5 M_{\odot}$ of iron, with little evolution as progenitor mass keeps increasing. Thanks to their much more extended progenitor mass range ($160 < m_* \le 300 M_{\odot}$), high-mass PISNe dominate the number of explosions once integrated over the IMF. In all cases, the explosion yields enough carbon, oxygen and other elements for the total metallicity of the ISM to transition to Pop. II star formation, but the difference in iron leads to the position on, or below, the plateau.

These results clearly show that the low-mass end of the mass-metallicity is highly sensitive to the exact assumptions and details of Pop. III physics (see also Jeon et al. 2021; Sanati et al. 2023). In particular, we find that, everything else being fixed, the number of dwarf galaxies on the plateau is tied to high-mass Pop. III PISNe, while the number below the plateau (or lack thereof) are tied to low-mass Pop. III PISNe. This is highly impactful, as one could leverage this knowledge and the distribution of observed iron metallicities into an effective constraint on the Pop. III IMF and PISN iron yields. We leave this to future work and show next how, beyond iron, this dichotomy of PISN progenitors is reflected on the chemical abundances of our simulated dwarfs.

## 4.2. Chemical abundances beyond the plateau

Figure 6 shows the average carbon, nitrogen and magnesium content of all MEGATRON dwarf galaxies (left, middle, right panels, respectively). The contours show the 2D distribution of simulated galaxies in average abundance versus $\langle[\text{Fe/H}]\rangle$. This is contrasted with the abundances of individual stars in ultra-faint dwarf galaxies (see Appendix A for exact citations).

Carbon and magnesium show broad agreement between the simulated and observed data, providing encouraging support for the PISN enrichment scenario we highlighted above. Nitrogen, however, is under-produced in our simulated dwarf galaxies compared to observed stars. Most stars in UFDs only have upper limits in nitrogen content (triangles), but the few with detections (diamonds) are several order of magnitude above the simulated data. PISNe are expected to have strong odd-even anti-correlation in chemical yields (e.g. Heger & Woosley 2002), with the lack of nitrogen yield being a defining characteristic (e.g. Salvadori et al. 2019). This is directly reflected in our simulated data for which PISN
<figure><img src="image_7.png" /><figcaption>FIG. 6.— Individual abundances in simulated MEGATRON dwarfs (red) compared to local data (gray). There is reasonable agreement for carbon and magnesium, although simulated data favour lower [C/Fe] and [Mg/Fe] abundances at [Fe/H] ≈ −2.5. This reflects our choice of yield models for high-mass PISNe which dominate enrichment close to plateau metallicities (see text for discussion). Nitrogen is vastly under-produced in our simulated dwarfs since PISNe are the sole enrichment channel in MEGATRON. Additional channels (stellar winds, rotation) could help alleviate this tension (Section 4.2).</figcaption></figure>

explosions dominate the enrichment budget.

We do not see this discrepancy in nitrogen content as invalidating the claim that PISNe drive the iron plateau. MEGATRON lacks many additional channels of nitrogen production, for example pre-explosion winds which could significantly boost nitrogen production in such massive low-metallicity stars (e.g. Hirschi 2007; Vink 2023; Charbonnel et al. 2023). We also assume nucleosynthetic yields from non-rotating stellar models, while stellar rotation also significantly promotes internal mixing and the nitrogen yield (e.g. Meynet & Maeder 2002; Chieffi & Limongi 2013; Limongi & Chieffi 2018). Lastly, we stress that high-redshift galaxies with surprisingly high nitrogen contents are currently being revealed by JWST (Cameron et al. 2023; Topping et al. 2024; Schaerer 2002), highlighting gaps in our understanding of nitrogen enrichment in metal-poor environments. As a result, we see the discrepancy in nitrogen content as a strong motivation for future studies to explore the impact of additional stellar evolution channels for nitrogen production.

The signatures of PISNe are also clearly visible in the distribution of carbon and magnesium abundances. Simulated dwarf galaxies that have plateau-like metallicities ($\langle[\text{Fe/H}]\rangle \approx -2.5$) have a solar to sub-solar carbon content ($\langle[\text{C/Fe}]\rangle < 0$). In contrast, the tail of iron-poor galaxies have strong carbon enhancements. Again, these differences reflects the dichotomy in their respective PISN progenitors. High-mass PISNe that dominate the plateau population (Figure 5) yield sub-solar [C/Fe] $\approx -0.5$ for a $m_\star = 180\,\text{M}_\odot$ PISN progenitor (Nomoto, Kobayashi & Tominaga 2013). Similarly, low-mass PISN progenitors that contribute to the tail yield super-solar [C/Fe] (e.g. [C/Fe] $\approx +2.5$ for a $m_\star = 140\,\text{M}_\odot$ PISN progenitor). Patterns in [Mg/Fe] (right panel) closely follow those in [C/Fe].

An encouraging success is that the shape of the distribution of carbon abundances matches well that of the observed data. However, the simulated data seems systematically shifted down, by $\approx 0.5$ dex. Such shifts are greater than typical uncertainties in observational data

(e.g. non-local thermodynamic equilibrium corrections can shift observed abundances by up to $\approx 0.3$ dex; discussion in Koutsouridou, Skúladóttir & Salvadori 2025) but well within the uncertainties in nucleosynthetic yields from PISNe. For example, differences in the treatment of convection and mixing in 1D stellar evolution codes can change the carbon yields of PISNe by a factor of a few (e.g. Heger & Woosley 2002; Nomoto, Kobayashi & Tominaga 2013; Takahashi, Yoshida & Umeda 2018), and uncertainties in nuclear reactions rate directly impacts predictions of the iron yields, the explosion energy and the window of explodability of PISNe (e.g. Farmer et al. 2019; Kawashimo et al. 2024).

Overall, we conclude that, given the chemical information available in MEGATRON, the PISN scenario is in reasonable agreement with observed data. However, we note that obtaining conclusive evidence for the full chemical signature of PISNe is likely to require a detailed study of many more chemical elements (e.g. Salvadori et al. 2019; Koutsouridou, Salvadori & Skúladóttir 2024) as the odd-even effects can be rapidly masked by additional enrichment channels (e.g. CCSNe; Ji, Frebel & Bromm 2015). Unfortunately, the choice of chemical elements tracked by our simulations is motivated by the wish to model gas cooling, limiting us to C, O, Ne, Mg, Si, S, Ca, Fe with even atomic numbers and only a single odd element (N). Future MEGATRON simulations will track additional elements relevant for stellar spectroscopy and abundance constraints (e.g. Na, Al, Mn, Cu, Zn) to enable such studies in the future.

So far, we have established that Pop. III PISNe are the key driver of the iron plateau and chemical abundances of UFDs in MEGATRON. However, given the large explosion energies, metals from PISN explosions are easily expelled from minihaloes (e.g. Bromm, Yoshida & Hernquist 2003; Greif et al. 2007) and can be ineffective in polluting metal-poor stars (e.g. Cooke & Madau 2014). In the next section, we explore what enables MEGATRON Pop. III star formation sites to retain the signatures of PISN enrichment.
### 4.3. The coupling between the plateau and Pop. III star formation sites

In this section, we connect our findings to the properties of Pop. III star formation sites in MEGATRON. Figure 7 shows the distribution of dark matter halo masses hosting newborn Pop. III stars<fnref n="3" /> (top panel) and the average iron metallicity when diluting typical Pop. III iron yields in the gas reservoir of a halo of this mass (bottom panel).

Focussing first on the top panel, we find that Pop. III stars in MEGATRON predominantly form in haloes with $M_{200} \ge 5 \times 10^6 M_{\odot}$, with a median of $M_{200} \approx 1-3 \times 10^7 M_{\odot}$ (vertical coloured lines). In A. Storck et al. (in prep), we explore the physical mechanisms setting this mass-scale in detail. Briefly, we find that a strong LW radiation field sourced by the central galaxy (Figure 1) suppresses $H_2$ formation and cooling in small minihaloes and elevates the mass-scale of Pop. III star formation (see also Ahn et al. 2009; Xu, Wise & Norman 2013; Wise et al. 2014; Xu et al. 2016; Zier et al. 2025 for similar results).

Such $10^7 M_{\odot}$ haloes can retain metals from much larger explosion energies ($\ge 10^{53}$ erg) due to the deeper gravitational potential wells and increased cooling losses (e.g. Kitayama & Yoshida 2005; Whalen et al. 2008; Cooke & Madau 2014). A single PISN is insufficient to evacuate metals, which are then readily re-incorporated in the next generation of (Pop. II) stars. This aligns with our findings in Section 4.1 and 4.2 that most MEGATRON dwarf galaxies undergo exactly one Pop. III explosion, retain the metals from this enrichment event, and form Pop. II stars with the chemical signature of PISNe.

We further connect this point to the mass-metallicity plateau in the bottom panel of Figure 7 in which we show the iron metallicity obtained when diluting three typical iron yields into the primordial gas of a dark matter halo. Here, we convert a virial mass $M_{200}$ to a total hydrogen mass assuming a 100% baryon fraction with $\Omega_b = 0.04916$ and a primordial mass fraction of hydrogen $Y_H = 0.76$ (e.g. Planck Collaboration et al. 2020). We then uniformly dilute the iron yield across this gas reservoir assuming that the explosion shock is contained within the virial radius (e.g. Whalen et al. 2008; Cooke & Madau 2014). We use $Y_{\text{Fe}} = 5 M_{\odot}, 0.5 M_{\odot}, 0.07 M_{\odot}$ to represent a high-mass Pop. III PISN progenitor, a low-mass Pop. III PISN progenitor, and a Pop. III or Pop. II CCSN respectively.

Diluting a fixed yield in a halo of increasing mass naturally leads to lower iron metallicities, as expected. But more importantly, diluting the yield of a high-mass PISN (brown line) in a $M_{200} \approx 10^7 M_{\odot}$ minihalo naturally leads to $\langle [\text{Fe/H}] \rangle \approx -2.5$, the average position of the plateau in the observed stellar mass-metallicity relation. In contrast, diluting the yield of a low-mass PISN (gold line) in the same halo leads to $\langle [\text{Fe/H}] \rangle \approx -3.5$, while a Pop. II CCSN (red line) leads to $\langle [\text{Fe/H}] \rangle \approx -4.5$. Since the bulk of Pop. III star formation sites in MEGATRON have $M_{200} \approx 10^7 M_{\odot}$ haloes (top panel) and undergo PISN explosions (Section 4.1), this naturally explains the emergence of the plateau and the tail of iron-poor galaxies in

<figure><img src="image_8.png" /><figcaption>FIG. 7.— Top: Distribution of dark matter halo masses hosting a Pop. III star formation event. Bottom: iron metallicity resulting from the dilution of an iron yield in the gas reservoir of a halo as a function of mass. In MEGATRON, Pop. III stars form in haloes with $M_{200} \ge 5 \times 10^6 M_{\odot}$ (medians marked at the top). This mass allows the retention of metals from high-energy PISN explosions (bottom shading, see Cooke & Madau (2014), fig. 3 and text for discussion). Diluting the iron yield of a high-mass PISN (brown) in a $10^7 M_{\odot}$ halo naturally leads to the emergence of the iron plateau at $\langle [\text{Fe/H}] \rangle \approx -2.5$ (horizontal dashed line). Low-mass PISNe (gold) and CCSNe (red) lead to lower iron metallicities, populating the tail of iron-poor galaxies below the plateau.</figcaption></figure>

### MEGATRON.

We stress that this is an approximate and heavily simplified calculation bypassing many of the complexities of metal mixing and retention in realistic haloes. In particular, we have assumed that all metals are retained within the virial radius and uniformly mixed within the gas reservoir. This is unlikely to be the case in reality (e.g. Ritter et al. 2015). Such additional effects are likely to explain the non-zero width of the distribution, but our first-order calculation readily explains how the plateau arises from the coupling between Pop. III yields and their halo formation sites.

According to this same argument, PISNe explosions in higher-mass haloes are not the only solution to produce plateau-like metallicities. In particular, a single CCSN event (red) in a much smaller minihalo ($M_{200} \approx 10^5 M_{\odot}$) could also yield $\langle [\text{Fe/H}] \rangle \approx -2.5$, especially since metal retention is also expected to be efficient given the lower explosion energy ($10^{51}$ erg; Cooke & Madau 2014, fig. 3). This highlights the degeneracies between the Pop. III IMF and yields (controlling iron production) and the halo mass-scale for Pop. III star formation (controlling iron retention and dilution). Furthermore, these effects are coupled together – as just one example, a more top-heavy Pop. III IMF raises the effective iron yield, but

<footnote n="3">³ Our output cadence ($\approx 5$ Myr) ensures that halo masses evolve minimally even if the Pop. III star is born in-between outputs</footnote>
also the amount of LW background and in turn the mass-scale for Pop. III star formation in nearby haloes (see e.g. Brauer et al. 2025b).

The emergence of the plateau in MEGATRON is thus neither fine-tuned nor obvious, rather resulting from the choices of Pop. III implementation we have made and the specific simulated environment of a Milky-Way-mass protogalaxy. Our results highlight the importance of accounting for these coupled effects self-consistently to predict the distribution of UFD metallicities and leverage local data to constrain Pop. III microphysics. We continue to explore these prospects in the next section, quantifying how much of the features we have highlighted at $z \approx 8$ would remain observable at $z = 0$ in local data.

### 5. SHOULD WE EXPECT AGREEMENT BETWEEN LOCAL DATA AND HIGH-REDSHIFT SIMULATIONS?

By comparing simulated dwarf galaxies at $z \approx 8$ with observed dwarf galaxies at $z = 0$, we have jumped 12 Gyr of cosmic evolution. In this Section, we leverage a defining feature of MEGATRON simulations – that they follow the formation of a Milky-Way-like galaxy – to show that bridging between these two epochs is quantitatively justified.

#### 5.1. Linking high- and low-redshift star formation histories and stellar masses

Stellar masses and star formation at the faint end of the stellar mass function of dwarf galaxies is expected to be primarily shaped by reionization quenching. In fact, all ultra-faint dwarfs observed at $z = 0$ are quenched and have formed most of their stars $\ge 10$ Gyr ago (e.g. Savino et al. 2025; Durbin et al. 2025). By contrast, most dwarfs plotted in Fig 2 are actively forming stars since none of the MEGATRON simulations are fully reionized by $z \approx 8$ (see Appendix C for more details). Catching these dwarfs before their ‘fossilization’ by reionization quenching could mean that their properties significantly evolve after $z \approx 8$ and they are not progenitors of local UFDs.

However, reionization quenching is strongly modulated by the specific spatial location of a given dwarf. Dwarfs closer to photon-producing regions will be reached first by ionization fronts and see their halo gas reservoir climb in temperature (e.g. Aubert et al. 2018; Katz et al. 2020). We thus expect a subpopulation of dwarfs in MEGATRON to already be reionization-quenched by $z \approx 8$, particularly those closer to the centre of the simulated volume where most star formation and photon production occurs (see temperature map in Figure 1).

We extract this subpopulation by selecting central faint dwarf galaxies that have formed no new stars for the last 100 Myrs (since $z \approx 10$). Figure 8 shows their mass-metallicity relation. (Restricting to dwarfs quenched for 200 Myrs – since $z \approx 11.5$ – shows the same trend, but prevents statistical tests by decreasing the sample size to too few objects).

As expected, the subpopulation of already-quenched dwarfs is much smaller than in Figure 2, but the overabundance of dwarfs along the plateau remains (see histogram on the right-hand panel). The sample sizes decrease from $N \ge 250$ per simulation to 3, 28, 17, and

51 quenched dwarfs for the ‘Efficient SF’, ‘Bursty SF’, ‘Varying IMF’ and ‘HN, $\epsilon_{eff}$’, respectively. The lowest number of quenched dwarfs in the ‘Efficient SF’ simulation correlates with the fact that this simulation is the least reionized at $z \approx 8$ (Figure 11). Combining all simulations together, we run a Kolmogorov-Smirnov test to quantitatively verify that the distribution of $\langle[\text{Fe/H}]\rangle$ in the quenched sample is consistent with the distribution of $\langle[\text{Fe/H}]\rangle$ in the overall sample. We find a p-value $\ge 0.1$, indicating that the two samples are statistically comparable at the $\approx 1\sigma$ level.

Importantly, Figure 8 also shows dwarfs over a range of masses, up to $M_* \approx 10^5 M_\odot$. MEGATRON dwarfs have thus already undergone significant evolution and several star formation episodes. We confirm that the median distance of the quenched sample is at most 25% away from the center of mass of the simulated volume. These dwarfs are also unlikely to restart star formation quickly, since their gas reservoir has a volume average temperature of $T = 1.1 \times 10^4$ K (20% increased compared to the median of the overall faint dwarf population) and a volume average H I fraction $\le 0.5$ (30% decreased).

We conclude that reionization feedback is already acting on the subpopulation of dwarfs closer to the centre, where the protogalaxy is collapsing and most ionizing photons are produced. These reionization-quenched dwarfs show no difference in $\langle[\text{Fe/H}]\rangle$ compared to the rest of the population, as expected since the physics setting their iron content is primarily internal (Section 4).

We further match these dwarfs to their counterpart haloes in the dark matter only simulation and find that their median peak halo mass over the whole Hubble time is $8 \times 10^8 M_\odot$. This is well below the threshold to reaccrete gas at later times once the UV background is homogeneous (see e.g. discussions in Rey et al. 2020; Benitez-Llambay & Frenk 2020). The dwarfs plotted in Figure 8 are thus already ‘fossilized’ at $z \approx 8$, and would stay quenched at all times like observed UFDs at $z = 0$.

Another bias in the comparison between high and low redshift is that stellar evolution will induce mass-loss over the Hubble time. In Figure 8, dots show the stellar mass at $z \approx 8$ while diamonds show the projected stellar mass at $z = 0$ for quenched dwarf galaxies (calculated as explained in Section 2.4). Accounting for this mass loss shifts stellar masses by $0.2 - 0.3$ dex, but makes little difference to the shape of the low-mass end of the mass-metallicity relation.

#### 5.2. Future environmental post-processing by a Milky Way's tidal field

The second major physical process at play between $z \approx 8$ and $z = 0$ is environmental processing of the dwarf galaxy population by the tidal field of the central galaxy. Figure 9 shows the $z \approx 8$ mass-metallicity relation of dwarf galaxies that we identify as survivors at $z = 0$ according to the procedure described in Section 2.4.

As when restricting to quenched dwarfs, the sample size significantly decreases. Overall, 25% of the faint dwarf galaxies at $z \approx 8$ survive in a bound remnant at $z = 0$, with the rest being disrupted by the central object and instead populating the diffuse halo (Figure 1). There is little variance in this reduction between individual MEGATRON simulations, with numbers going from $N \ge 250$ to (66, 76, 92, 70) survivors for the ‘Efficient
<figure><img src="image_9.png" /><figcaption>FIG. 8.— Same as Figure 2 but for the subpopulation of dwarf galaxies already quenched by reionization at $z \approx 8$. Dots show their stellar mass at $z \approx 8$, while diamonds show the projected stellar mass at $z = 0$ after accounting for a Hubble time of stellar mass loss. The faint end of the mass-metallicity relation is unchanged by either effects.</figcaption></figure>

SF', 'Bursty SF', 'Varying IMF' and 'HN, $\epsilon_{eff}$', respectively. This is to be expected since the dark-matter-only history of the central object is the same in all cases.

The shape of the distribution of $\langle[\text{Fe/H}]\rangle$ is unchanged by considering only surviving dwarfs. As for the quenched sample, there is no statistical difference between this subpopulation and the parent sample (p-value = 0.55, i.e. the two populations are $1\sigma$ compatible). Together with Section 5.1, this firmly establishes that the comparison between $z \approx 8$ simulated dwarfs and $z = 0$ observed dwarfs is quantitatively justified. Features in the mass-metallicity relation identified and explained at $z \approx 8$ are not transient at high redshift, but rather robust predictions that would survive to $z = 0$.

Our analysis highlights the power of using highly-resolved and detailed high-redshift simulations to inform Galactic archaeology studies aiming to constrain how early chemical enrichment proceeded in the proto-Milky Way. Particle-tagging techniques based on dark-matter-only simulations are a useful and cheap complementing tool to the expensive high-redshift hydrodynamical simulations like MEGATRON. But they are also limited in their accuracy (e.g. discussions in Bailin et al. 2014; Le Bret et al. 2017; Cooper et al. 2017). In the future, we plan on strengthening our low-redshift forecasts by mapping onto another hydrodynamical simulation of the same initial conditions evolved to $z = 0$ with the VINTERGATAN model. This will allow us to account for the potential impact of disc formation on the survival of faint dwarfs (e.g. D'Onghia et al. 2010).

## 6. DISCUSSION AND CONCLUSION

We have presented new results from the MEGATRON simulation suite showing how Pop. III stars drive the distribution of iron metallicities in ultra-faint dwarf galaxies. The MEGATRON offers a new theoretical framework to connect high-redshift and low-redshift chemical enrichment data by (i) offering high resolution in the ISM

<figure><img src="image_10.png" /><figcaption>FIG. 9.— Same as Figure 2 but for the subpopulation of dwarf galaxies that survive as $z = 0$ bound satellites. The plateau at $\langle[\text{Fe/H}]\rangle \approx -2.5$ and the tail of iron-deficient dwarfs both survive in bound satellites around the host Milky Way at $z = 0$.</figcaption></figure>

of small dwarfs ($\approx 5$pc) to resolve Pop. III and Pop. II star formation sites; (ii) an extensive physical model (hydrodynamics, radiative transfer, non-equilibrium primordial and metal chemistry) to self-consistently capture the build-up of radiation backgrounds from the first stars to hydrogen reionization; and (iii) simulating the environment of a protogalaxy that will collapse into a Milky-Way-mass object to allow connections with local Galactic archaeology data.

Out of four MEGATRON simulations varying Pop. II star formation and feedback models, dwarf galaxies in three simulations provide an excellent match to the observed mass-metallicity relation of UFDs at $z = 0$, capturing the slope and normalization of the relation over several decades in stellar mass (Figure 2). This agreement is non-trivial, as shown by the fourth simulation ('Efficient SF') which over-shoots in iron content for $M_* \ge 10^6 \text{ M}_\odot$ due to weaker galactic outflows. These results reaffirm that the mass-metallicity relation is a sensitive probe of the star formation and feedback assumptions in dwarf galaxy formation models (see also Agertz et al. 2020; Prgomet et al. 2022; Sanati et al. 2023; Go et al. 2025; Rey et al. 2025).

At the low-mass end of the relation ($M_* \le 10^5 \text{ M}_\odot$), all MEGATRON simulations showcase an over-abundance of galaxies with $\langle[\text{Fe/H}]\rangle \approx -2.5$, similar to the observed plateau at $z = 0$ in the local data, and a tail of iron-deficient dwarf galaxies ($\langle[\text{Fe/H}]\rangle \le -3.0$) that is yet to be observed in local data. Both features are robust to large variations in Pop. II star formation and feedback variations and survive as observable features at $z = 0$ (Figure 8, 9).

We causally connect the iron plateau to internal enrichment from a single explosion from a high-mass Pop. III PISNe ($\ge 180 \text{ M}_\odot$, Figure 4, 5). The over-dense MEGATRON environment drives a strong LW background, ensuring that such explosions occur in haloes ($M_{200} \approx 10^7 \text{ M}_\odot$) where metals are retained efficiently and incorporated into subsequent generations of stars (Figure 7).
This scenario naturally explains the normalization of the plateau and provides qualitatively consistent carbon and magnesium abundances, although offset by ~0.5 dex relative to observations (Figure 6). The tail of iron-poor galaxies is instead created by low-mass Pop. III PISNe explosions (≈ 170 $M_\odot$) and to a lesser extent by external pollution.

Our results support the interpretation that local UFDs are internally enriched at high redshift, with their chemical properties closely related to the microphysics of early chemical enrichment by the first stars (see also Jeon, Besla & Bromm 2017; Sanati et al. 2023; Go et al. 2025). As an example, following the MEGATRON results, a complete census of iron metallicities in UFDs could be leveraged to constrain the Pop. III IMF. In particular, the relative frequency of high-mass to low-mass PISNe would directly set the ratio between dwarf galaxies on the plateau and in the tail. A lack of such iron-poor tail might instead push the IMF towards extremely top-heavy values (i.e. only high-mass PISNe), or instead point to another origin to the UFD iron plateau (e.g. CCSNe in a lower-mass minihalo, see discussion in Section 4.3).

Current local data is yet to observe a UFD with galaxy-averaged $\langle[\text{Fe/H}]\rangle \le -3.0$. If the underlying distribution of $\langle[\text{Fe/H}]\rangle$ in UFDs is similar to that predicted by MEGATRON, we expect $\approx 20\%$ of UFDs to have such low iron contents. The probability of not observing any dwarfs with $\langle[\text{Fe/H}]\rangle \ge -3.0$ out of 38 independent samples with is $(1-0.2)^3 \approx 0.02\%$, i.e. disfavouring the existence of an iron-poor tail. However, the completeness of iron metallicity measurements in UFDs is sparse, both in surveying different objects and in the number of stars measured per galaxy. For example, Pace et al. (2025) recently reported the first characterization of Pictor II, finding $\langle[\text{Fe/H}]\rangle \approx -2.99$ pushing down towards the tail identified in MEGATRON. Furthermore, in the same galaxy, Chiti et al. (2025) identified the first ultra iron-poor star ($\langle[\text{Fe/H}]\rangle \le -4.0$, not included in the average above), proving that (i) such low iron contents can be expected in UFDs; and (ii) that galaxy averages derived from low numbers of iron measurements could be biased high. To further illustrate this point, Appendix D shows the composite metallicity distribution functions of MEGATRON dwarfs. Even for galaxies with an average iron content close to the plateau, individual stars can populate a long tail with $\langle[\text{Fe/H}]\rangle \le -3.0$ due to the combination of low-mass Pop. III PISN and Pop. II CCSN enrichment. The MEGATRON tail of low $\langle[\text{Fe/H}]\rangle$ is just an extension in galactic average of a tail of iron-poor stars in all dwarfs that might currently be hidden by the limited number of stars observed per galaxy.

Our results thus provide strong theoretical incentives to complete the census of UFD iron metallicities. Major efforts are already underway in this direction and will likely yield fruits in the coming years. The now-online Vera C. Rubin Telescope is expected to propose many new photometric candidates (e.g. Mutlu-Pakdil et al. 2021) and better-multiplexed spectroscopic surveys will greatly enhance the number of chemical abundances measured in dwarf galaxies (e.g. Skúladóttir et al. 2023). Simulation efforts like MEGATRON will be key to interpret these data and leverage them to constrain the physics of

early chemical enrichment.

In particular, the ability of MEGATRON to self-consistently model both Pop. III and Pop. II star formation and feedback in a cosmological-radiation hydrodynamics framework is key to our results. The emergence of the iron plateau is not a fine-tuned result, but rather stems from the combination of (i) a strong LW background forcing Pop. III stars to form in higher-mass haloes that can retain metals; and (ii) a top-heavy IMF favouring PISNe explosions to obtain the right iron content and a reasonable chemistry for faint dwarfs.

Capturing the first aspect is only possible with cosmological-radiation hydrodynamics. Such setups have become increasingly common in high-redshift simulations (e.g. Xu, Wise & Norman 2013; Wise et al. 2014; Kimm et al. 2017; Brauer et al. 2025a; Kang et al. 2025; Zier et al. 2025), although the requirement to resolve both Pop. III and Pop. II star formation sites over a large volume still poses computational challenges. Such simulations have repeatedly highlighted the inhomogeneity and environmental dependence of the LW background (e.g. Xu et al. 2016; Zier et al. 2025). In the case of MEGATRON, the LW background is always elevated due to the specific cosmic environment we have picked – an overdense region collapsing early to form a disc galaxy at $z=0$ and which, as a result, contains more stars for its halo mass than the cosmic average at $z=8$ (Katz et al. 2025). As a result, irrespective of the four Pop. II feedback models we have explored, Pop. III stars always form in higher-mass haloes. This makes it clear that accounting, at least partially, for the special environment of the proto Milky Way will be key for accurate predictions of early chemical enrichment observables in Galactic data.

Another important choice driving our results lies in a top-heavy Pop. III IMF that favours PISN explosions. This is a classical assumption motivated by early simulations pointing out that the lack of metal cooling in primordial gas leads to gas clumps with characteristic masses $\ge 100 M_\odot$ close to PISN progenitors (e.g. Abel et al. 1998; Bromm et al. 2001; Nakamura & Umemura 2001; Abel, Bryan & Norman 2002; Bromm, Coppi & Larson 2002). More recent models, however, show that hydrodynamical instabilities and turbulence during minihalo collapse likely drives additional gas fragmentation, potentially lowering the characteristic mass of primordial star forming clouds to $\le 100 M_\odot$ (e.g. Turk, Abel & O'Shea 2009; Greif et al. 2012; Sugimura et al. 2020; Wollenberg et al. 2020). This in turn would lead to much less top-heavy IMF predictions than we have assumed here (e.g. Susa, Hasegawa & Tominaga 2014; Stacy, Bromm & Lee 2016; Prole et al. 2022). However, turbulence-amplified magnetic fields (e.g. Sharda, Federrath & Krumholz 2020) and the radiative input of new-born Pop. III stars (e.g. Hirano et al. 2014, 2015, though see also Jaura et al. 2022) can prevent fragmentation in the primordial minihalo and raise the characteristic mass again. Overall, the uncertainties surrounding predictions of the Pop. III IMF are significant (see Klessen & Glover 2023 for a more in-depth review), making our choice justified, but one among many.

Despite these uncertainties, our results show that the Pop. III IMF leaves a clear imprint on the iron metallicity distribution of UFDs. Strengthening the connection between high-redshift simulations and local archeologi-
cal data as we are doing here is thus full of promises but will require a more systematic exploration of parameter space. We plan on addressing this in future work re-simulating our volume with a different Pop. III IMF and studying the abundances and spatial distribution of all metal-poor stars, rather than just those inside dwarf galaxies. This will be highly complementary to galactic chemical enrichment models that aim to constrain the Pop. III IMF from stellar abundances (e.g. de Bennasuti et al. 2017; Ishigaki et al. 2018; Rossi, Salvadori & Skúladóttir 2021), for example allowing to calibrate some of their free parameters (e.g. metal retention factors) from first principles. The end goal of this exercise will be to jointly leverage high-redshift emission line data and local chemical enrichment data to constrain the physics of early chemical enrichment. This is a timely and exciting endeavour since both JWST and Galactic spectroscopic surveys are maturing at the same time.

## ACKNOWLEDGEMENTS

AJC acknowledges funding from the 'First Galaxies' Advanced Grant from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (Grant agreement No.789056). AFS acknowledges support from the Science and Technology Facilities Council (STFC) for a PhD studentship. FRM is supported by the Kavli Institute for Cosmological physics at the University of Chicago through an endowment from the Kavli Foundation and its founder Fred Kavli. GCJ acknowledges support by the Science and Technology Facilities Council (STFC), by the ERC through Advanced Grant 695671 'QUENCH', and by the UKRI Frontier Research grant RISEandFALL. HK acknowledges support from FACCTS. The material in this manuscript is based upon work supported by NASA under award No. 80NSSC25K7009. KM acknowledges the Flemish Fund for Scientific Research (FWO-Vlaanderen), Grant number 1169822N. M.S. acknowledges the support from the Swiss National Science Foundation under Grant No. P500PT.214488. NC acknowledges support from the Science and Technology Facilities Council (STFC) for a PhD studentship. OA acknowledges support from the Knut and Alice Wallenberg Foundation, the Swedish Research Council (grant 2019-04659), the Swedish National Space Agency (SNSA Dnr 2023-00164), and the LMK foundation. TK

is supported by the National Research Foundation of Korea (RS-2022-NR070872 and RS-2025-00516961) and by the Yonsei Fellowship, funded by Lee Youn Jae.

This work made extensive use of the dp265, dp016, dp373, and dp379 projects on the DiRAC ecosystem. This work was performed using the DiRAC Data Intensive service at Leicester, operated by the University of Leicester IT Services, which forms part of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The equipment was funded by BEIS capital funding via STFC capital grants ST/K000373/1 and ST/R002363/1 and STFC DiRAC Operations grant ST/R001014/1. This work used the DiRAC@Durham facility managed by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The equipment was funded by BEIS capital funding via STFC capital grants ST/P002293/1, ST/R002371/1 and ST/S002502/1, Durham University and STFC operations grant ST/R000832/1. This work was performed using resources provided by the Cambridge Service for Data Driven Discovery (CSD3) operated by the University of Cambridge Research Computing Service (www.csd3.cam.ac.uk), provided by Dell EMC and Intel using Tier-2 funding from the Engineering and Physical Sciences Research Council (capital grant EP/T022159/1), and DiRAC funding from the Science and Technology Facilities Council (www.dirac.ac.uk). DiRAC is part of the National e-Infrastructure. This work has made use of the Infinity Cluster hosted by Institut d'Astrophysique de Paris. We thank Stephane Rouberol for running smoothly this cluster for us. The authors thank Jonathan Patterson for smoothly running the Glamdring Cluster hosted by the University of Oxford, where part of the data processing was performed. The authors also acknowledge financial support from Oriel College's Research Fund.

We thank the developers and maintainers of PYNBODY (Pontzen et al. 2013), TANGOS (Pontzen & Tremmel 2018), YT (Turk et al. 2011), NUMPY (van der Walt, Colbert & Varoquaux 2011; Harris et al. 2020), SciPY (Virtanen et al. 2020), JUPYTER (Ragan-Kelley et al. 2014), MATPLOTLIB (Hunter 2007), the Astrophysics Data Service and the arXiv preprint repository for providing open-source softwares that were used extensively in this work.

## REFERENCES

Abel T., Anninos P., Norman M. L., Zhang Y., 1998, ApJ, 508, 518

Abel T., Bryan G. L., Norman M. L., 2002, Science, 295, 93

Agertz O. et al., 2020, MNRAS, 491, 1656

Agertz O. et al., 2021, MNRAS, 503, 5826

Ahn K., Shapiro P. R., Iliev I. T., Mellema G., Pen U.-L., 2009, ApJ, 695, 1430

Ahvazi N., Benson A., Sales L. V., Nadler E. O., Weerasooriya S., Du X., Bovill M. S., 2024, MNRAS, 529, 3387

Andersson E. P., Rey M. P., Pontzen A., Cadiou C., Agertz O., Read J. I., Martin N. F., 2025, ApJ, 978, 129

Applebaum E., Brooks A. M., Christensen C. R., Munshi F., Quinn T. R., Shen S., Tremmel M., 2021, ApJ, 906, 96

Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ARA&A, 47, 481

Aubert D. et al., 2018, ApJL, 856, L22

Bailin J., Bell E. F., Valluri M., Stinson G. S., Debattista V. P., Couchman H. M. P., Wadsley J., 2014, ApJ, 783, 95

Behroozi P. S., Wechsler R. H., Wu H.-Y., 2013, ApJ, 762, 109

Bellazzini M., Gennari N., Ferraro F. R., 2005, MNRAS, 360, 185

Belokurov V. et al., 2007, ApJ, 654, 897

Benitez-Llambay A., Frenk C., 2020, MNRAS, 498, 4887

Benson A. J., Frenk C. S., Lacey C. G., Baugh C. M., Cole S., 2002, MNRAS, 333, 177

Bhardwaj A., Rejkuba M., Ngeow C.-C., Marconi M., Ripepi V., Samantaray A. S., Singh H. P., 2024, AJ, 167, 247

Bland-Hawthorn J., Sutherland R., Webster D., 2015, ApJ, 807, 154

Boettcher E. et al., 2013, AJ, 146, 94

Bonifacio P., Caffau E., François P., Spite M., 2025, ARA&A, 33, 2

Bovill M. S., Ricotti M., 2009, ApJ, 693, 1859

Brauer K. et al., 2025a, ApJ, 980, 41

Brauer K. et al., 2025b, Aeos: The Impact of Pop III Initial Mass Function and Star-by-Star Models in Galaxy Simulations

Bromm V., Coppi P. S., Larson R. B., 2002, ApJ, 564, 23

Bromm V., Ferrara A., Coppi P. S., Larson R. B., 2001, MNRAS, 328, 969

Bromm V., Yoshida N., Hernquist L., 2003, ApJ, 596, L135

Brown T. M. et al., 2014, ApJ, 796, 91
Bruce J., Li T. S., Pace A. B., Heiger M., Song Y.-Y., Simon J. D., 2023, ApJ, 950, 167

Bullock J. S., Johnston K. V., 2005, ApJ, 635, 931

Bullock J. S., Kravtsov A. V., Weinberg D. H., 2000, ApJ, 539, 517

Cadiou C., Dubois Y., Pichon C., 2019, A&A, 621, A96

Cadiou C., Katz H., Rey M. P., et al., 2025, see arXiv Oct 1st

Cameron A. J., Katz H., Rey M. P., Saxena A., 2023, MNRAS, 523, 3516

Cantu S. A. et al., 2021, ApJ, 916, 81

Carlin J. L., Grillmair C. J., Muñoz R. R., Nidever D. L., Majewski S. R., 2009, ApJL, 702, L9

Carlin J. L., Sand D. J., 2018, ApJ, 865, 7

Carlin J. L. et al., 2017, AJ, 154, 267

Casey Q. O. et al., 2025, ApJ, 984, 148

Cerny W. et al., 2025, ApJ, 979, 164

Cerny W. et al., 2023a, ApJL, 953, L21

Cerny W. et al., 2023b, ApJ, 953, 1

Cerny W. et al., 2021a, ApJ, 910, 18

Cerny W. et al., 2021b, ApJL, 920, L44

Cerny W. et al., 2023c, ApJ, 942, 111

Charbonnel C., Schaeerer D., Prantzos N., Ramírez-Galeano L., Fragos T., Kuruvandothi A., Marques-Chaves R., Gieles M., 2023, A&A, 673, L7

Chieffi A., Limongi M., 2013, ApJ, 764, 21

Chiti A. et al., 2023, AJ, 165, 55

Chiti A. et al., 2021, Nat Ast, 5, 392

Chiti A. et al., 2025, A second-generation star in a relic dwarf galaxy

Chiti A., Simon J. D., Frebel A., Pace A. B., Ji A. P., Li T. S., 2022, ApJ, 939, 41

Choustikov N., Rey M. P., Cadiou C., et al., 2025, Meg P1

Cicuéndez L. et al., 2018, A&A, 609, A53

Cooke R. J., Madau P., 2014, ApJ, 791, 116

Cooper A. P., Cole S., Frenk C. S., Le Bret T., Pontzen A., 2017, MNRAS, 469, 1691

Cooper A. P. et al., 2010, MNRAS, 406, 744

Correnti M., Bellazzini M., Ferraro F. R., 2009, MNRAS, 397, L26

Crnojević D., Sand D. J., Zaritsky D., Spekkens K., Willman B., Hargis J. R., 2016, ApJ, 824, L14

Cullen F. et al., 2021, MNRAS, 505, 903

Curti M. et al., 2023, MNRAS, 518, 425

Curti M. et al., 2024, A&A, 684, A75

Dall'Ora M. et al., 2006, ApJ, 653, L109

Dall'Ora M. et al., 2012, ApJ, 752, 42

de Bennassuti M., Salvadori S., Schneider R., Valiante R., Omukai K., 2017, MNRAS, 465, 926

De Lucia G., Helmi A., 2008, MNRAS, 391, 14

D'Onghia E., Springel V., Hernquist L., Keres D., 2010, ApJ, 709, 1138

Drlica-Wagner A. et al., 2015, ApJ, 813, 109

Durbin M. J. et al., 2025, The HST Legacy Archival Uniform Reduction of Local Group Imaging (LAURELIN). I. Photometry and Star Formation Histories for 36 Ultra-faint Dwarf Galaxies

Efstathiou G., 1992, MNRAS, 256, 43P

Escala I. et al., 2018, MNRAS, 474, 2194

Farmer R., Renzo M., de Mink S. E., Marchant P., Justham S., 2019, ApJ, 887, 53

Frebel A., Norris J. E., 2015, ARA&A, 53, 631

Fritz T. K., Carrera R., Battaglia G., Taibi S., 2019, A&A, 623, A129

Fu S. W. et al., 2023, ApJ, 958, 167

Garofalo A., Clementini G., Cusano F., Muraveva T., Monti L., 2025, A&A, 695, A88

Garofalo A. et al., 2013, ApJ, 767, 62

Gnedin N. Y., 2000, ApJ, 542, 535

Go M. et al., 2025, ApJ, 986, 214

Greco C. et al., 2008, ApJ, 675, L73

Greif T. H., Bromm V., Clark P. C., Glover S. C. O., Smith R. J., Klessen R. S., Yoshida N., Springel V., 2012, MNRAS, 424, 399

Greif T. H., Johnson J. L., Bromm V., Klessen R. S., 2007, ApJ, 670, 1

Griffen B. F., Dooley G. A., Ji A. P., O'Shea B. W., Gómez F. A., Frebel A., 2018, MNRAS, 474, 443

Hansen T. T., Simon J. D., Li T. S., Sharkey D., Ji A. P., Thompson I. B., Reggiani H. M., Galarza J. Y., 2024, ApJ, 968, 21

Harris C. R. et al., 2020, Nature, 585, 357

Hartwig T., Bromm V., Klessen R. S., Glover S. C. O., 2015, MNRAS, 447, 3892

Heger A., Woosley S. E., 2002, ApJ, 567, 532

Heiger M. E. et al., 2024, ApJ, 961, 234

Heintz K. E. et al., 2023, Nat Ast, 7, 1517

Hirano S., Hosokawa T., Yoshida N., Omukai K., Yorke H. W., 2015, MNRAS, 448, 568

Hirano S., Hosokawa T., Yoshida N., Umeda H., Omukai K., Chiaki G., Yorke H. W., 2014, ApJ, 781, 60

Hirschi R., 2007, A&A, 461, 571

Homma D. et al., 2019, PASJ, 71, 94

Homma D. et al., 2024, PASJ, 76, 733

Homma D. et al., 2018, PASJ, 70, S18

Hsiao T. Y.-Y. et al., 2024, First direct carbon abundance measured at $z>10$§ in the lensed galaxy MACS0647$-$JD

Hunter J. D., 2007, CiSE, 9, 90

Ishigaki M. N., Tominaga N., Kobayashi C., Nomoto K., 2018, ApJ, 857, 46

Isobe Y. et al., 2023, ApJ, 959, 100

Jaura O., Glover S. C. O., Wollenberg K. M. J., Klessen R. S., Geen S., Haemmerlé L., 2022, MNRAS, 512, 116

Jenkins S. A., Li T. S., Pace A. B., Ji A. P., Koposov S. E., Mutlu-Pakdil B., 2021, ApJ, 920, 92

Jeon M., Besla G., Bromm V., 2017, ApJ, 848, 85

Jeon M., Bromm V., Besla G., Yoon J., Choi Y., 2021, MNRAS, 502, 1

Ji A. P., Frebel A., Bromm V., 2015, MNRAS, 454, 659

Ji A. P. et al., 2021, ApJ, 921, 32

Jiang R., Zhao G., Li H., Xing Q., 2024, ApJ, 976, 68

Joshi G. D., Pontzen A., Agertz O., Read J., Rey M. P., 2025, eprint arXiv:2507.05401, arXiv:2507.05401

Kang C., Kimm T., Han D., Katz H., Devriendt J., Slyz A., Teysier R., 2025, A&A, 693, A149

Karczmarek P. et al., 2015, AJ, 150, 90

Katz H., 2022, MNRAS, 512, 348

Katz H., Kimm T., Ellis R. S., Devriendt J., Slyz A., 2023, MNRAS, 524, 351

Katz H. et al., 2022, arXiv e-prints, 2211.04626

Katz H. et al., 2020, MNRAS, 494, 2200

Katz H., Rey M. P., Cadiou C., et al., 2025, see arXiv Oct 1st

Katz H., Rey M. P., Cadiou C., Kimm T., Agertz O., 2024, arXiv e-prints, 2411.07282

Kawashimo H., Sawada R., Suwa Y., Moriya T. J., Tanikawa A., Tominaga N., 2024, MNRAS, 531, 2786

Kim D., Jerjen H., 2015, ApJ, 808, L39

Kim D. et al., 2016, ApJ, 833, 16

Kimm T., Katz H., Haehnelt M., Rosdahl J., Devriendt J., Slyz A., 2017, MNRAS, 466, 4826

Kirby E. N., Cohen J. G., Guhathakurta P., Cheng L., Bullock J. S., Gallazzi A., 2013, ApJ, 779, 102

Kirby E. N., Cohen J. G., Simon J. D., Guhathakurta P., 2015, ApJ, 814, L7

Kirby E. N., Rizzi L., Held E. V., Cohen J. G., Cole A. A., Manning E. M., Skillman E. D., Weisz D. R., 2017, ApJ, 834, 9

Kitayama T., Yoshida N., 2005, ApJ, 630, 675

Klessen R. S., Glover S. C. O., 2023, ARA&A, 61, 65

Kobayashi C., Ishigaki M. N., Tominaga N., Nomoto K., 2014, ApJL, 785, L5

Kobayashi C., Umeda H., Nomoto K., Tominaga N., Ohkubo T., 2006, ApJ, 653, 1145

Koposov S. E., Belokurov V., Torrealba G., Evans N. W., 2015a, ApJ, 805, 130

Koposov S. E. et al., 2015b, ApJ, 811, 62

Koposov S. E. et al., 2011, ApJ, 736, 146

Koposov S. E. et al., 2018, MNRAS, 479, 5343

Koutsouridou I., Salvadori S., Skúladóttir Á., 2024, ApJL, 962, L26

Koutsouridou I., Salvadori S., Skúladóttir Á., Rossi M., Vanni I., Pagnini G., 2023, MNRAS, 525, 190

Koutsouridou I., Skúladóttir Á., Salvadori S., 2025, A&A, 699, 32

Kroupa P., 2001, MNRAS, 322, 231

Kuehn C. et al., 2008, ApJL, 674, L81

Kulkarni G., Keating L. C., Haehnelt M. G., Bosman S. E. I., Puchwein E., Chardin J., Aubert D., 2019, MNRAS, 485, L24

Laseret I. H. et al., 2024, A&A, 681, A70

Le Bret T., Pontzen A., Cooper A. P., Frenk C., Zolotov A., Brooks A. M., Governato F., Parry O. H., 2017, MNRAS, 468, 3212

Lee M. G., Yuk I.-S., Park H. S., Harris J., Zaritsky D., 2009, ApJ, 703, 692

Li T. S. et al., 2017, ApJ, 838, 8

Li T. S. et al., 2018, ApJ, 857, 145

Li Z. et al., 2025, Insights on Metal Enrichment and Environmental Effect at $z\approx5-7$§ with JWST ASPIRE/EIGER and Chemical Evolution Model

Limongi M., Chieffi A., 2018, ApJS, 237, 13

Longeard N. et al., 2018, MNRAS, 480, 2609

Madau P., Giallongo E., Grazian A., Haardt F., 2024, ApJ, 971, 75

Magg M., Hartwig T., Agarwal B., Frebel A., Glover S. C. O., Griffen B. F., Klessen R. S., 2018, MNRAS, 473, 5308

Marks M., Kroupa P., Dabringhausen J., Pawlowski M. S., 2012, MNRAS, 422, 2246
Martínez-Vázquez C. E. et al., 2021, AJ, 162, 253

Martínez-Vázquez C. E. et al., 2015, MNRAS, 454, 1509

Martínez-Vázquez C. E. et al., 2019, MNRAS, 490, 2183

Mateo M., Olszewski E. W., Walker M. G., 2008, ApJ, 675, 201

McConnachie A. W., 2012, AJ, 144, 4

Medina G. E. et al., 2018, ApJ, 855, 43

Meynet G., Maeder A., 2002, A&A, 381, L25

Moskowitz A. G., Walker M. G., 2020, ApJ, 892, 27

Muñoz R. R., Côté P., Santana F. A., Geha M., Simon J. D., Oyarzún G. A., Stetson P. B., Djorgovski S. G., 2018, ApJ, 860, 66

Munshi F., Brooks A. M., Christensen C., Applebaum E., Holley-Bockelmann K., Quinn T. R., Wadsley J., 2019, ApJ, 874, 40

Musella I. et al., 2009, ApJ, 695, L83

Mutlu-Pakdil B. et al., 2018, ApJ, 863, 25

Mutlu-Pakdil B. et al., 2021, ApJ, 918, 88

Mutlu-Pakdil B. et al., 2020, ApJ, 902, 106

Nakamura F., Umemura M., 2001, ApJ, 548, 19

Nakane M. et al., 2025, Fe Abundances of Early Galaxies at $\S z=9-12$ Derived with Deep JWST Spectra

Noh Y., McQuinn M., 2014, MNRAS, 444, 503

Nomoto K., Kobayashi K., Tominaga N., 2013, ARA&A, 51, 457

Nomoto K., Tominaga N., Umeda H., Kobayashi C., Maeda K., 2006, NuPhA, 777, 424

Oakes E. K., Hoyt T. J., Freedman W. L., Madore B. F., Tran Q. H., Cerny W., Beaton R. L., Seibert M., 2022, ApJ, 929, 116

Okamoto S., Arimoto N., Yamada Y., Onodera M., 2012, ApJ, 744, 96

Oñorbe J., Boylan-Kolchin M., Bullock J. S., Hopkins P. F., Kereš D., Faucher-Giguère C.-A., Quataert E., Murray N., 2015, MNRAS, 454, 2092

Paardekooper J.-P., Khochfar S., Dalla Vecchia C., 2015, MNRAS, 451, 2544

Paardekooper J.-P., Pelupessy F. I., Altay G., Kruip C. J. H., 2011, A&A, 530, A87

Pace A. B., 2025, OpJA, 8

Pace A. B. et al., 2020, MNRAS, 495, 3022

Pace A. B. et al., 2025, OpJA, 8, 112

Pillepich A. et al., 2018, MNRAS, 473, 4077

Planck Collaboration et al., 2020, A&A, 641, A6

Pontzen A., Rey M. P., Cadiou C., Agertz O., Teissier R., Read J., Orkney M. D. A., 2021, MNRAS, 501, 1755

Pontzen A., Roškar R., Stinson G., Woods R., 2013, Astrophysics Source Code Library, ascl:1305.002

Pontzen A., Tremmel M., 2018, ApJS, 237, 23

Prgomet M., Rey M. P., Andersson E. P., Segovia Otero A., Agertz O., Renaud F., Pontzen A., Read J. I., 2022, MNRAS, 513, 2326

Prole L. R., Clark P. C., Klessen R. S., Glover S. C. O., 2022, MNRAS, 510, 4019

Ragan-Kelley M., Perez F., Granger B., Kluyver T., Ivanov P., Frederic J., Bussonnier M., 2014, Am. Geophys. Un., 2014, H44D

Rey M. P. et al., 2023, MNRAS, 521, 995

Rey M. P., Pontzen A., 2018, MNRAS, 474, 45

Rey M. P., Pontzen A., Agertz O., Orkney M. D. A., Read J. I., Rosdahl J., 2020, MNRAS, 497, 1508

Rey M. P., Starkenburg T. K., 2022, MNRAS, 510, 4208

Rey M. P. et al., 2025, MNRAS, 541, 1195

Richstein H. et al., 2024, ApJ, 967, 72

Richstein H. et al., 2022, ApJ, 933, 217

Ricotti M., Gnedin N. Y., 2005, ApJ, 629, 259

Ritter J. S., Sluder A., Safranek-Shrader C., Milosavljević M., Bromm V., 2015, MNRAS, 451, 1190

Rosdahl J., Blaizot J., Aubert D., Stranex T., Teissier R., 2013, MNRAS, 436, 2188

Rosdahl J. et al., 2018, MNRAS, 479, 994

Rossi M., Salvadori S., Skúladóttir Á., 2021, MNRAS, 503, 6026

Rossi M., Salvadori S., Skúladóttir Á., Vanni I., 2023, MNRAS, 522, L1

Rossi M., Salvadori S., Skúladóttir Á., Vanni I., Koutsouridou I., 2025, ApJ, 987, 121

Roth N., Pontzen A., Peiris H. V., 2016, MNRAS, 455, 974

Salvadori S., Bonifacio P., Caffau E., Korotin S., Andreevsky S., Spite M., Skuladottir Á., 2019, MNRAS, 487, 4261

Salvadori S., Ferrara A., 2009, MNRAS, 395, L6

Sanati M., Jeanquartier F., Revaz Y., Jablonka P., 2023, A&A, 669, A94

Sand D. J., Strader J., Willman B., Zaritsky D., McLeod B., Caldwell N., Seth A., Olszewski E., 2012, ApJ, 756, 79

Sanders R. L., Shapley A. E., Topping M. W., Reddy N. A., Brammer G. B., 2024, ApJ, 962, 24

Sandford N. R. et al., 2025, Chemodynamics of BoötesI with $\$$\{5\}$\$§: Revised Velocity Gradient, Dark Matter Density, and Galactic Chemical Evolution Constraints

Savino A. et al., 2025, ApJ, 979, 205

Schaerer D., 2002, A&A, 382, 28

Schaerer D., Marques-Chaves R., Xiao M., Korber D., 2024, A&A, 687, L11

Scholte D. et al., 2025, MNRAS, 540, 1800

Shapiro P. R., Giroux M. L., Babul A., 1994, ApJ, 427, 25

Sharda P., Federrath C., Krumholz M. R., 2020, MNRAS, 497, 336

Simon J. D., 2019, ARA&A, 57, 375

Simon J. D. et al., 2015, ApJ, 808, 95

Simon J. D., Geha M., 2007, ApJ, 670, 313

Simon J. D. et al., 2011, ApJ, 733, 46

Simon J. D. et al., 2017, ApJ, 838, 11

Simon J. D. et al., 2020, ApJ, 892, 137

Skúladóttir Á. et al., 2023, The Messenger, 190, 19

Smith S. E. T. et al., 2023, AJ, 166, 76

Somerville R. S., 2002, ApJ, 572, L23

Spencer M. E., Mateo M., Olszewski E. W., Walker M. G., McConnachie A. W., Kirby E. N., 2018, AJ, 156, 257

Spencer M. E., Mateo M., Walker M. G., Olszewski E. W., 2017, ApJ, 836, 202

Stacy A., Bromm V., Lee A. T., 2016, MNRAS, 462, 1307

Stanton T. M. et al., 2025, MNRAS, 537, 1735

Stanton T. M. et al., 2024, MNRAS, 532, 3102

Stanway E. R., Eldridge J. J., Becker G. D., 2016, MNRAS, 456, 485

Stetson P. B., Fiorentino G., Bono G., Bernard E. J., Monelli M., Irannicola G., Gallart C., Ferraro I., 2014, PASP, 126, 616

Stiavelli M., Morishita T., Chiaberge M., Grillo C., Rosati P., Schuldt S., Trenti M., Treu T., 2023, ApJL, 957, L18

Storypa S., Pontzen A., Peiris H., Roth N., Rey M. P., 2021, ApJS, 252, 28

Sugimura J., Matsumoto T., Hosokawa T., Hirano S., Omukai K., 2020, ApJL, 892, L14

Susa H., Hasegawa K., Tominaga N., 2014, ApJ, 792, 32

Susa H., Umemura M., 2004, ApJ, 600, 1

Takahashi K., Yoshida T., Umeda H., 2018, ApJ, 857, 111

Tan C. Y. et al., 2025, ApJ, 979, 176

Teissier R., 2002, A&A, 385, 337

Tolstoy E., Hill V., Tosi M., 2009, ARA&A, 47, 371

Topping M. W. et al., 2025, ApJ, 980, 225

Topping M. W. et al., 2024, MNRAS, 529, 3301

Torrealba G. et al., 2018, MNRAS, 475, 5085

Torrealba G., Koposov S. E., Belokurov V., Irwin M., 2016a, MNRAS, 459, 2370

Torrealba G. et al., 2016b, MNRAS, 463, 712

Trebitsch M., Blaizot J., Rosdahl J., Devriendt J., Slyz A., 2017, MNRAS, 470, 224

Turk M. J., Abel T., O'Shea B., 2009, Science, 325, 601

Turk M. J., Smith B. D., Oishi J. S., Skory S., Skillman S. W., Abel T., Norman M. L., 2011, ApJS, 192, 9

van der Walt S., Colbert S. C., Varoquaux G., 2011, Comput. Sci. Eng., 13, 22

Vink J. S., 2023, A&A, 679, L9

Virtanen P. et al., 2020, Nat Methods, 17, 261

Vivas A. K., Martínez-Vázquez C., Walker A. R., 2020, ApJS, 247, 35

Vivas A. K., Martínez-Vázquez C. E., Walker A. R., Belokurov V., Li T. S., Erkal D., 2022, ApJ, 926, 78

Vivas A. K. et al., 2016, AJ, 151, 118

Walker M. G., Mateo M., Olszewski E. W., Bailey J. I., Koposov S. E., Belokurov V., Evans N. W., 2015, ApJ, 808, 108

Walker M. G., Mateo M., Olszewski E. W., Sen B., Woodroofe M., 2009, AJ, 137, 3109

Walker M. G., Olszewski E. W., Mateo M., 2015, MNRAS, 448, 2717

Walsh S. M., Willman B., Sand D., Harris J., Seth A., Zaritsky D., Jerjen H., 2008, ApJ, 688, 245

Wang M. Y. et al., 2019, ApJ, 881, 118

Weisz D. R., Dolphin A. E., Skillman E. D., Holtzman J., Gilbert K. M., Dalcanton J. J., Williams B. F., 2014, ApJ, 789, 148

Whalen D., van Veelen B., O'Shea B. W., Norman M. L., 2008, ApJ, 682, 49

Wheeler C. et al., 2019, MNRAS, 490, 4447

Wheeler V., Kravtsov A., Chiti A., Katz H., Semenov V. A., 2025, eprint arXiv:2507.03182, arXiv:2507.03182

Willman B., Geha M., Strader J., Strigari L. E., Simon J. D., Kirby E., Ho N., Warres A., 2011, AJ, 142, 128

Willman B. et al., 2006, arXiv e-prints, astro

Wise J. H., Cen R., 2009, ApJ, 693, 984

Wise J. H., Demchenko V. G., Halicek M. T., Norman M. L., Turk M. J., Abel T., Smith B. D., 2014, MNRAS, 442, 2560

Wise J. H., Turk M. J., Norman M. L., Abel T., 2012, ApJ, 745, 50

Wollenberg K. M. J., Glover S. C. O., Clark P. C., Klessen R. S., 2020, MNRAS, 494, 1871
# APPENDIX

## A. DATA COMPILATION

Dwarf galaxy data presented here is drawn from the Local Volume Database (Pace 2025) which compiles UFD data (including the LMC system) from Bellazzini, Gennari & Ferraro (2005); Belokurov et al. (2007); Bhardwaj et al. (2024); Boettcher et al. (2013); Bruce et al. (2023); Cantu et al. (2021); Carlin et al. (2009, 2017); Carlin & Sand (2018); Casey et al. (2025); Cerny et al. (2021a,b, 2023b,a,c, 2025); Chiti et al. (2021, 2023, 2022); Cicuéndez et al. (2018); Correnti, Bellazzini & Ferraro (2009); Crnojević et al. (2016); Dall'Ora et al. (2006, 2012); Drlica-Wagner et al. (2015); Fritz et al. (2019); Garofalo et al. (2013, 2025); Greco et al. (2008); Hansen et al. (2024); Heiger et al. (2024); Homma et al. (2018, 2019, 2024); Jenkins et al. (2021); Ji et al. (2021); Karczmarek et al. (2015); Kim & Jerjen (2015); Kim et al. (2016); Kirby et al. (2013, 2015, 2017); Koposov et al. (2011, 2015a,b, 2018); Kuehn et al. (2008); Lee et al. (2009); Li et al. (2017, 2018); Longeard et al. (2018); Martínez-Vázquez et al. (2015, 2019, 2021); Mateo, Olszewski & Walker (2008); McConnachie (2012); Medina et al. (2018); Moskowitz & Walker (2020); Muñoz et al. (2018); Musella et al. (2009); Mutlu-Pakdil et al. (2018, 2020); Oakes et al. (2022); Pace et al. (2020); Richstein et al. (2022, 2024); Sand et al. (2012); Simon & Geha (2007); Simon et al. (2011, 2015, 2017); Simon (2019); Simon et al. (2020); Smith et al. (2023); Spencer et al. (2017, 2018); Stetson et al. (2014); Tan et al. (2025); Torrealba et al. (2016b,a, 2018); Vivas et al. (2016); Vivas, Martínez-Vázquez & Walker (2020); Vivas et al. (2022); Walker et al. (2009); Walker, Olszewski & Mateo (2015); Walker et al. (2015); Walsh et al. (2008); Wang et al. (2019); Willman et al. (2006, 2011).

## B. EXTERNALLY-ENRICHED DWARFS

Figure 4 shows that a small fraction of faint dwarfs have transitioned to Pop. II star formation without ever hosting a Pop. III star.

Figure 10 shows their spatial distribution across the simulated volume (left panel) and the carbon abundances of the oldest Pop. II stellar population in these externally-enriched dwarfs (i.e. their first star formed since no Pop. III are present; right panel). Black contours show the rest of the MEGATRON population. Our externally-enriched dwarfs are not particularly clustered in space, tracking the broader large-scale structure. The carbon abundances of their oldest Pop. II stars (red dots) are also representative of the rest of the dwarf population (black contour), with most stars having $-0.5 \le [C/Fe] \approx 0.5$. With our assumptions, such ratios are typical for yields of low-mass Pop. III PISNe ($[C/Fe] \approx -0.5$) and low-metallicity Pop. II CCSNe ($[C/Fe] \approx 0.5$). We conclude that, occasionally, a star-forming galaxy can enrich a neighbouring minihalo to the critical metallicity for Pop. II star formation before it forms its own Pop. III star. Determining whether the ejecta is primarily from Pop. II or Pop. III stars will require more detailed tracking of the metal flows across our simulated volume, which we plan to tackle in future work using our Lagrangian particle tracers (Cadiou, Dubois & Pichon 2019).

Using the same methodology as in Section 5.2 and across all four simulations, we quantify that only 3 externally-enriched dwarfs out of 22 will survive at $z=0$ after processing by the tidal field of the Milky-Way-like host. This fraction (13%) is small but not negligible, suggesting that several observed UFDs around the Milky Way could have been externally-enriched at high redshift. However, since their chemical abundances are indistinguishable from the rest of the population, the observational prospects to identify them would be challenging.

## C. REIONIZATION HISTORY OF MEGATRON SIMULATIONS

We show the reionization history of the four MEGATRON simulations in Figure 11. We compute the ionized fraction of the gas within the high-resolution region of the simulation and show the volume-weighted ionized fractions as a function of redshift. All simulations are more reionized than the cosmic average inferred from data (e.g. Kulkarni et al. 2019; Madau et al. 2024), as expected from sampling the overdense environment that will collapse into a Milky-Way-like galaxy.

We observe significant differences between models, with the 'Efficient SF' simulation being the least reionized at a given redshift, and the 'HN, $\epsilon_{eff}$' simulation being completely reionized by $z \approx 8$. These differences are not driven by the number of ionizing photons produced, but by feedback-modulated escape fractions. For example, the 'Efficient SF' simulation produces significantly more ionizing photons due to its higher star formation rates (Katz et al. 2025, fig. 6). But the higher-density ISM drives a much lower escape fraction of these ionizing photons (see N. Choutikov et al. in prep, fig. 10), in turn reducing the ability to reionize the surrounding IGM. These findings align with other numerical simulations predicting that the escape of ionizing photons is feedback regulated and strongly dependent on the ISM structure, with mechanical feedback from SNe playing a larger role than radiative feedback (e.g. Wise & Cen 2009; Paardekooper et al. 2011; Paardekooper, Khochfar & Dalla Vecchia 2015; Trebitsch et al. 2017; Kimm et al. 2017; Rosdahl et al. 2018). Differences in reionization history explains why the fraction of dwarf galaxies that are reionization-quenched varies across simulations in Section 5.1 and is highest in the 'HN, $\epsilon_{eff}$' simulation.
<figure><img src="image_11.png" /></figure>

FIG. 10.— The spatial (left) and chemical (right) distribution of externally-enriched dwarfs that contain no Pop. III stars within them. These dwarfs are compatible with the rest of the population in both spaces, reflecting the stochasticity of external enrichment by either Pop. III- or Pop. II-driven outflows.

<figure><img src="image_12.png" /></figure>

FIG. 11.— The volume-averaged ionized hydrogen fraction as a function of redshift in the four MEGATRON simulations. All four simulations are more reionized than the cosmic average (e.g. Kulkarni et al. 2019; Madau et al. 2024), with the differences between models being driven by feedback-modulated escape fractions of ionizing photons.

## D. IRON DISTRIBUTION FUNCTIONS

Figure 3 and Figure 6 shows the galaxy-averaged abundances of faint dwarfs. However, each galaxy contains a metal distribution function (MDF) from which this average is derived. The shape of this distribution and the number of stars for which we have spectroscopic abundance measurements can introduce biases in the measurement of $\langle[\text{Fe/H}]\rangle$.
<figure><img src="image_13.png" /><figcaption>FIG. 12.— The iron distribution functions of individual stars within faint dwarf galaxies on, above and below the plateau. As expected, dwarfs on the plateau have a peak around [Fe/H] ≈ −2.5, but also show a tail of iron-poor stars down to [Fe/H] ≤ −4.0. Similarly, dwarfs with an average metallicity below the plateau can still have individual stars with higher iron contents.</figcaption></figure>

(see e.g. Andersson et al. 2025 and Sandford et al. 2025 for theoretical and observational perspectives, respectively).

Figure 12 exemplifies this point, showing the composite iron (left) and carbon (right) distribution functions of individual stars within dwarfs with $M_* \le 10^5 M_\odot$ that are on, above, and below the iron plateau. The shape of the distribution of individual stars across all galaxies clearly follows the shape of the distribution of galaxy-averaged $\langle[\text{Fe/H}]\rangle$ in Figure 3, as expected given the large numbers of both stars and galaxies involved. Also expected, dwarfs on the plateau have a peak around $[\text{Fe/H}] \approx -2.5$. But Figure 12 also shows that plateau dwarfs have a tail of iron-poor stars down to $[\text{Fe/H}] \le -4.0$ similar to the object recently discovered by Chiti et al. (2025). Similarly, dwarfs with an average metallicity below the plateau can still have individual stars with plateau-like iron contents. The same picture holds for $[\text{C/Fe}]$ (right panel), where significant overlap exists between the distributions of dwarfs on, above and below the plateau.

From this, it is clear that a complete census of MDFs and a careful unpicking of their distribution over the full population of dwarfs will be key to accurately measure their galaxy-averaged $\langle[\text{Fe/H}]\rangle$ and capture the exact relationship between Pop. III enrichment and UFD chemistry. Small-sample statistics and individual detections could otherwise bias the interpretation.