# On the dearth of C-enhanced metal-poor stars in the galactic bulge

G. Pagnini, S. Salvadori, M. Rossi, D. Aguado, I. Koutsouridou, Á.

Skúladóttir

► To cite this version:

G. Pagnini, S. Salvadori, M. Rossi, D. Aguado, I. Koutsouridou, et al. On the dearth of C-enhanced metal-poor stars in the galactic bulge. Monthly Notices of the Royal Astronomical Society, 2023, 521, pp.5699-5711. 10.1093/mnras/stad912. insu-04848732

HAL Id: insu-04848732

https://insu.hal.science/insu-04848732v1

Submitted on 19 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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
# On the dearth of C-enhanced metal-poor stars in the galactic bulge

G. Pagnini,<fnref n="1" />* S. Salvadori,<fnref n="2" />,<fnref n="3" /> M. Rossi,<fnref n="2" />,<fnref n="3" /> D. Aguado,<fnref n="2" />,<fnref n="3" />,<fnref n="4" /> I. Koutsouridou,<fnref n="2" />,<fnref n="3" /> and Á. Skúladóttir<fnref n="2" />,<fnref n="3" />

Accepted 2023 March 22. Received 2023 March 20; in original form 2022 October 7

## ABSTRACT

The chemical fingerprints of the first stars are retained within the photospheres of ancient unevolved metal-poor stars. A significant fraction of these stellar fossils is represented by stars known as Carbon-Enhanced Metal-Poor (CEMP), $[C/Fe] > +0.7$ and $[Fe/H] < -2$, which are likely imprinted by low-energy primordial supernovae. These CEMP stars are largely observed in the Galactic halo and ultrafaint dwarf galaxies, with values reaching $[C/Fe] = +4.5$. The Galactic bulge is predicted to host the oldest stars, but it shows a striking dearth of CEMP stars with $[C/Fe] \gtrsim +2.0$. Here, we explore the possible reasons for this anomaly by performing a statistical analysis of the observations of metal-poor stars in combination with the predictions of Lambda cold dark matter models. We suggest that the dearth of CEMP stars with high $[C/Fe]$ is not due to the low statistics of observed metal-poor stars but is the result of the different formation process of the bulge. N-body simulations show that the first star-forming haloes which end up in the bulge are characterized by the highest star formation rates. These rates enable the formation of rare massive first stars exploding as pair-instability supernovae (PISNe), which wash out the signature of primordial faint supernovae. We demonstrate that the mean $[C/Fe]$ of first stars polluted environments decreases with the increasing contribution of PISNe. We conclude that the dearth of CEMP stars in the Galactic bulge indirectly probes the existence of elusive PISNe, and propose a novel method which exploits this lack to constrain the mass distribution of the first stars.

**Key words:** stars: carbon – stars: Population III – Galaxy: bulge – galaxies: formation – galaxies: high-redshift.

## 1 INTRODUCTION

Today we are surrounded by stars – the Milky Way galaxy alone contains hundreds of billions of stars – but there was a time, billions of years ago, when stars were absent and the Universe was extremely simple. At that time the Universe was mostly neutral and mainly composed by hydrogen and helium produced during the big bang Nucleosynthesis. The first stars played a fundamental role in the evolution of the Universe, as they were responsible for the transition from this very simple and early stage to the more complex one visible today. In fact, the first stars were the sources of the first chemical elements heavier than lithium and of the first hydrogen-ionizing photons. Hence, they initiated the extended processes of reionization and metal-enrichment.

Within the standard Lambda cold dark matter (ΛCDM) model of structure formation, the first stars (referred to as Population III or Pop III stars) are predicted to have formed within the first few hundred million years after the big bang, corresponding to redshifts of $z \sim 15–20$. The primordial birth environment of Pop III stars, characterized by lack of heavy elements and dust, may have resulted in a higher Pop III characteristic stellar mass with respect to present-day stars (e.g. Tan & McKee 2004; Hosokawa et al. 2011) although sub-solar mass stars might also have been able to form (e.g. Greif

et al. 2011; Stacy et al. 2016). Since the primordial star formation process is still very poorly understood, we can state that the initial mass function (IMF) of Pop III stars is almost completely unknown. Despite long searches (Beers & Christlieb 2005; Caffau et al. 2013), zero-metallicity stars have not yet been observed, confirming the hypothesis of a primordial IMF biased towards more massive stars than the present-day IMF (e.g. Salvadori et al. 2007; Magg et al. 2019; Rossi et al. 2021). Massive Pop III stars explode as supernovae (SNe) polluting the surrounding medium with their chemical products, whose yields depend upon the mass of the progenitor star along with the explosion energy (e.g. Nomoto et al. 2004; Heger & Woosley 2010). Hence, even if we cannot directly observe short-lived zero-metallicity stars, we can still catch their long-lived descendants (Keller et al. 2007; Starkenburg et al. 2017b). Low-mass (Pop II) stars formed in environments enriched by the chemical products of the first stars to the critical metallicity value, $Z > Z_{cr}$, at which normal star formation is expected to proceed (Bromm et al. 2001; Schneider et al. 2002). In this context Stellar Archaeology operates: searching for the chemical signatures of the first stellar generations in the photospheres of old (> 12 Gyr) and metal-poor stars that dwell in our Galaxy and its ancient dwarf satellites.

Which components of the Milky Way should be examined to find the oldest living stars? The first surveys looking for first star descendants have attempted to select metal-poor stars in the Galactic stellar halo (e.g. Beers et al. 1992; Christlieb 2003), which is expected to be the most metal-poor component. Typically the iron-

<footnote n="1">¹GEPI, Observatoire de Paris, PSL Research University, CNRS, Place Jules Janssen, F-92195 Meudon, France</footnote>

<footnote n="2">²Dipartimento di Fisica e Astronomia, Università degli Studi di Firenze, via G. Sansone 1, I-50019 Sesto Fiorentino, Italy</footnote>

<footnote n="3">³INAF/Osservatorio Astrofisico di Arcetri, Largo E. Fermi 5, I-50125 Firenze, Italy</footnote>

<footnote n="4">⁴Instituto de Astrofísica de Canarias, Vía Láctea, E-38205 La Laguna, Tenerife, Spain</footnote>

<footnote>* E-mail: giulia.pagnini@obspm.fr</footnote>
abundance<fnref n="1" />, [Fe/H], is measured in these surveys and used as a metallicity indicator. The most metal-poor stars are then selected for high-resolution follow-up, often revealing chemically peculiar stars with strong enhancements or deficiencies of particular elements (e.g. Christlieb 2003; Caffau et al. 2011; Keller et al. 2014; Bonifacio et al. 2015; Caffau et al. 2016; Aguado et al. 2018; François et al. 2018; González Hernández et al. 2020).

The best-known type of chemically interesting stellar object at $[Fe/H] < -2$ is the carbon-enhanced metal-poor (CEMP) class, which has $[C/Fe] > +0.7$ (e.g. Beers & Christlieb 2005; Aoki et al. 2007). This class can be divided into two main populations: (i) carbon-rich stars that also exhibit an excess in heavy elements formed by the slow neutron-capture processes, having $[Ba/Fe] > 1$ and named CEMP-s stars, and (ii) carbon-rich stars with no excess of the heavy elements, having $[Ba/Fe] < 0$ and known as CEMP-no stars. The CEMP-s stars are commonly assumed to be chemically enriched by mass transfer from a binary companion star that has gone through the asymptotic giant branch (AGB) phase (Abate et al. 2015), and these objects are preferentially found in binary systems (e.g. Suda et al. 2004; Lucatello et al. 2005; Starkenburg et al. 2014; Hansen et al. 2016b). On the other hand, CEMP-no stars are not primarily found in binary systems (Lucatello et al. 2005; Norris et al. 2012; Starkenburg et al. 2014; Hansen et al. 2016a) and even those that have a binary companion (e.g. Arentsen et al. 2019) show high values of $^{12}$C/$^{13}$C, which implies that the surface composition has not been altered by mass transfer (see Aguado et al. 2022, 2023b). Hence the C-excess in CEMP-no stars is expected to be representative of the interstellar medium (ISM) out of which they formed, which was likely primarily polluted by the first stellar generations (e.g. Salvadori et al. 2015; De Bennassuti et al. 2017).

The observed chemical abundance patterns of the most Fe-poor, $[Fe/H] < -4$, CEMP-no stars are indeed consistent with the yields of Pop III stars exploding as a faint supernovae and experiencing mixing and fallback (e.g. Iwamoto et al. 2005; Marassi et al. 2014). Because of their low explosion energy, only the outer layers of the Pop III progenitor star, which are rich in C and other light elements, can be expelled by faint SNe. On the other hand, the inner part, which is rich in Fe-peak elements, falls back into the center forming a neutron star or a black hole (e.g. Heger & Woosley 2010). The increased fraction of CEMP-no stars with decreasing $[Fe/H]$ further supports such a link with Pop III star pollution (De Bennassuti et al. 2017). Furthermore, the high values of $^{12}$C/$^{13}$C recently observed in CEMP-no stars are consistent with an imprint from low-energy Pop III SNe (Aguado et al. 2023b) and rule out the so-called 'spinstars' (e.g. Meynet et al. 2006) as main pollutants. CEMP-no stars have been found in a significant fraction in the Galactic halo (Yong et al. 2012; Carollo et al. 2014; Placco et al. 2014; Lee et al. 2017; Yoon et al. 2018; Lee et al. 2019) and in the faintest satellites of the Milky Way, the so-called ultrafaint dwarf galaxies (UFDs; Norris et al. 2010; Lai et al. 2011; Gilmore et al. 2013; Spite et al. 2018), which are the oldest galaxies in the Local group (e.g. Simon et al. 2010; Gallart et al. 2021). On the contrary, CEMP-no stars seem to be quite rare in the more luminous classical dwarf spheroidal (dSph) galaxies (e.g. Skúladóttir et al. 2015; Yoon et al. 2020), which show more complex and longer star formation histories with respect to UFDs (see Salvadori et al. 2015, for a global view). Ultimately, these

observational results confirm that the descendants of the first stars are preferentially found in purely ancient environments.

Relying on the $\Lambda$CDM model and hierarchical clustering, White & Springel (2000) first predicted the oldest stars to be in the inner part of the Milky Way, i.e. the Galactic bulge. This idea was confirmed in the following years through different numerical simulations of the Milky Way (Diemand et al. 2005; Tumlinson 2009; Salvadori et al. 2010; Starkenburg et al. 2017a). More recent simulations have shown this to also hold for different galaxies, including the dwarf satellites of Lyman Break galaxies at redshift $z \approx 6$ (Gelli et al. 2020). Unfortunately, the Galactic bulge is a dusty and overcrowded region, predominantly populated by metal-rich stars, so metal-poor objects are difficult to find (Zoccali et al. 2008; Ness et al. 2013; Howes et al. 2016). The EMBLA Survey (Extremely Metal-poor bulge stars with AAOmega), has been the first in attempting to discover candidate metal-poor stars in this inner region (Howes et al. 2014, 2016). Although it successfully identified $\approx 30$ stars with $[Fe/H] < -2$, from its observations it is clear that there is a dearth of CEMP-no stars in this region; in fact, only one CEMP-no star was found, having $[Fe/H] = -3.48$ and $[C/Fe] = +0.98$ (Howes et al. 2015). A reason for a lack of CEMP stars in the EMBLA sample could be a selection bias introduced by the SkyMapper photometric selection since photometric colours may have been affected by the strong CH absorption of CEMP stars with $[C/Fe] > +2$, placing them outside the selected region in the colour-colour diagram. Therefore, the EMBLA Survey could have missed stars with extremely large C enhancements but stars with mild C enhancements should still have been found.

More recently, the Pristine Inner Galaxy Survey (PIGS; Arentsen et al. 2020) targeted the Galactic bulge with low/intermediate resolution spectroscopy ($R \approx 1300$ at 3700–5500 Å, and $R \approx 11000$ at 8400–8800 Å), collecting 1900 stars with $[Fe/H] < -2.0$, which is currently the largest sample of confirmed very metal-poor stars in the inner Galaxy. Since s-process abundance measurements for this sample are unavailable, a different CEMP classification was applied based on the absolute carbon abundance<fnref n="2" />, A(C), and the $[Fe/H]$ of the stars (Yoon et al. 2016). According to Bonifacio et al. (2015) (see also Spite et al. 2013) CEMP stars can be divided in two main groups: (i) the high-carbon band, A(C) $> 7.4$, largely containing CEMP-s stars that typically have higher iron-abundance; (ii) the low-carbon band, A(C) $< 7.4$, predominantly containing CEMP-no stars at lower $[Fe/H]$. Using this classification while correcting for evolutionary effects in Red Giant Branch (RGB) stars (e.g. Placco et al. 2014), the PIGS survey identified 24 new CEMP-no candidate stars in the Galactic bulge. Still, the overall fraction of CEMP-no stars obtained by PIGS is only $\lesssim 6$ per cent at $[Fe/H] < -2$, i.e. much lower than what is found in the Galactic halo ($\approx 20$ per cent, see Arentsen et al. 2021). Furthermore, $\simeq 96$ per cent of their CEMP-no candidates have only a moderate C-enhancement, $+0.7 < [C/Fe] < +1.2$, with only one CEMP-no star at $[Fe/H] \simeq -3.5$ that shows a high $[C/Fe] \simeq +2.2$, while similar values are more frequently observed both in the Galactic halo and in UFDs (e.g. fig. 1 from Salvadori et al. 2015). As explained in Arentsen et al. (2021) and Arentsen et al. (2022), the selection of metal-poor stars in PIGS could be significantly affected by biases in the photometry against CEMP stars with strong CH lines, i.e. those with very high carbon abundances and/or cooler temperatures. However, stars with $[Fe/H] < -3.0$ are not expected to be substantially biased, except for

<footnote n="1">¹Throughout this paper we will be using the notation $[A/B] \equiv \log_{10}(N_A/N_B)_*$ – $\log_{10}(N_A/N_B)_\odot$, where $N_A$ and $N_B$ refer to the numbers of atoms of elements A and B, respectively.</footnote>

<footnote n="2">²Elemental abundances can also be referred to as an 'absolute' scale, relative to the number of hydrogen atoms, defined as $A(C) \equiv \log_{10}(N_C/N_H)_* + 12$.</footnote>
the coolest ($T_{eff} < 4750$ K) and most carbon-rich ($[C/Fe] > +2.0$) stars. For stars with $[Fe/H] < -2.0$ and $<-2.5$, slightly warmer and less carbon-rich stars are also affected, but still most of the CEMP-no stars are expected to stay within their selection. Therefore they expect that the CEMP-no fraction is less impacted by a photometric bias compared to that of CEMP-s stars.

These observational results raise several questions: Why is there an apparent dearth of CEMP-no stars in the Galactic bulge? And why are CEMP-no stars with high $[C/Fe]$ so rare in this ancient environment? Is it entirely due to a bias in the selection of metal-poor stars or is there a process that reduces the CEMP-no fraction in the Galactic bulge? The aim of this paper is to answer these questions, and it is structured as follows. In Section 2, we will carry out a preliminary analysis on the observational data of metal-poor stars in the different regions of the Local Group in order to understand whether the dearth of C-enhanced stars in the Galactic bulge is due to a statistical effect. In Section 3, we will illustrate the cosmological model used to follow the evolution of the Galaxy that combines a N-body simulation, following the hierarchical assembly of a Milky Way (MW)-like galaxy, with a semi-analytical model that follows the evolution of baryons. The results obtained from this model and the additional analytical calculations performed will be described in Section 4. Finally, in Section 5 we will draw the conclusions of our analysis with the associated implications, and we will list the prospects for future works.

## 2 OBSERVATIONAL DATA ANALYSIS

To understand why CEMP-no stars are less common in the Galactic bulge with respect to other ancient environments of the Local Group, we should first analyse if this dearth can be just a statistical effect. Indeed, since CEMP-no stars are more common at decreasing $[Fe/H]$, the probability to discover them in environments dominated by metal-rich stars, such as the Galactic bulge, can be intrinsically very low (Salvadori et al. 2015).

### 2.1 The fraction of CEMP stars

Henceforth, we will simply speak of CEMP stars when referring to the subclass of CEMP-no stars, i.e to stars with $[Fe/H] < -2$, $[C/Fe] > +0.7$ within the low-carbon band, A(C) < 7.4. The upper panel of Fig. 1 shows the measured $[C/Fe]$ versus $[Fe/H]$ for a sample of 984 halo stars with carbon measurements from the JINAbase<fnref n="3" /> (Abohalima & Frebel 2018). The abundance and stellar parameter data collected there are based on high-resolution spectroscopic studies ($R = \lambda/\Delta\lambda \gtrsim 15\,000$, with the majority having $R = 30\,000-40\,000$) found in the literature. In the same plot, we also show the measurements for the sample of bulge stars (Arentsen, private communication) as selected in Arentsen et al. (2021) in their fig. 9. In detail they select measurements from a FERRE analysis with the following cuts on gravity, temperature and $[C/Fe]$ uncertainty: $\log g < 3.5$, $4600\,K < T_{eff} < 5500\,K$ and $\epsilon_{[C/Fe]} < 0.5$. As stated in Section 1, for bulge stars there are not available measurements of the s-process element Ba. Thus, the A(C) classification is used to distinguish among CEMP-no and CEMP-s candidates. Although many halo stars have Ba measurements, here we adopt the same definition to separate CEMP-no and CEMP-s stars to have a self-consistent classification. As explained in detail in Yoon et al. (2016), despite the differences in the two criteria, the clear distinction in A(C)

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Top panel: Measured $[C/Fe]$ versus $[Fe/H]$ for stars in the Galactic halo (black and grey points; JINAbase) and the bulge PIGS sample (orange and red points; Arentsen et al. 2021). Bottom panel: Fraction of CEMP-no stars in the Galactic halo in different $[Fe/H]$ ranges: grey includes all CEMP-no stars, while cyan is only CEMP-no stars with carbon-excess, $[C/Fe] > +2.0$. Error bars are Poissonian errors derived from Gehrels (1986). The number of stars in each bin, $N_*([Fe/H])$, is listed on top.</figcaption></figure>

between CEMP-s and the CEMP-no stars in their halo sample appears to be as successful, and likely more astrophysically fundamental, for the separation of these sub-classes as the previously recommended criterion based on $[Ba/Fe]$ abundance ratios. Fig. 1 (upper panel) clearly shows how the $[C/Fe]$ values in Galactic halo stars are higher at $[Fe/H] < -3$ compared to stars in the bulge, and how in general CEMP-no stars are more frequent in the halo, although we cannot exclude that this can be due to a statistical effect given the rarity of extremely metal-poor stars in the Galactic bulge (see also Section 2.3). To quantify this, we compute the fraction of CEMP-no stars, defined as the ratio of CEMP-no stars over the total number of stars at a given $[Fe/H]$:

$$F_{\text{CEMP}}([\text{Fe/H}]) = \frac{N_{\text{CEMP}}([\text{Fe/H}])}{N_*([\text{Fe/H}])} \quad (1)$$

The Galactic halo is currently the best observed ancient environment of the Local Group and for which we have the largest statistics regarding the chemical properties of the metal-poor stellar population (e.g. Yong et al. 2012; Bonifacio et al. 2021). Therefore, we compute $F_{\text{CEMP}}^{\text{halo}}$ for the sample of halo stars in Fig. 1. The results are shown as grey histogram in the bottom panel of Fig. 1

<footnote n="3">³https://jinabase.pythonanywhere.com/</footnote>
<table><caption>**Table 1.** Luminosities and average stellar metallicities of the Local Group environments taken into account in this work.</caption><thead><tr><th></th><th>Least luminous<br>UFDs</th><th>Bootes I</th><th>Sculptor</th><th>Fornax</th><th>Galactic<br>bulge</th></tr></thead><tbody><tr><td>L<sub>*</sub></td><td>&le;10<sup>4</sup> L<sub>&odot;</sub></td><td>10<sup>4.5</sup> L<sub>&odot;</sub></td><td>10<sup>6.34</sup> L<sub>&odot;</sub></td><td>10<sup>7.25</sup> L<sub>&odot;</sub></td><td>10<sup>10</sup> L<sub>&odot;</sub></td></tr><tr><td>&lang;[Fe/H]&rang;</td><td>&le;-2.2</td><td>-2.1</td><td>-1.8</td><td>-1.0</td><td>0.0</td></tr></tbody></table>

where we have also computed $F_{\text{CEMP}}^{\text{halo}}$ for CEMP stars with high carbon-excess, $[C/Fe] > +2.0$ (cyan) (cyan histogram). As we can see, $F_{\text{CEMP}}$ increases as $[Fe/H]$ decreases and reaches the value of ≈1 for $[Fe/H] ≤ -5$. Because of the large number of CEMP-no stars at $[Fe/H] ≤ -4$ discovered during the last 10 yr (see e.g. Norris et al. 2012; Yong et al. 2012; Keller et al. 2014; Bonifacio et al. 2015; Aguado et al. 2017; Da Costa et al. 2019; Nordlander et al. 2019; Li et al. 2022, and references therein), we can now derive $F_{\text{CEMP}}$ at different $[Fe/H]$ values, equally spaced every 0.5 dex. Still, we should remind that these values suffer various systematic uncertainties. On one hand, as suggested in Arentsen et al. (2022), high-resolution spectroscopic compilations of very metal-poor stars ($[Fe/H] < -2.0$) might be biased towards (very) carbon-rich objects at the 'high' metallicity end ($[Fe/H] > -3.0$) due to their follow-up strategies. This leads to an overestimate of the CEMP fraction in this metallicity range. Furthermore, $F_{\text{CEMP}}^{\text{halo}}$ might increase when correcting the carbon measurements to account for the internal depletion of carbon occurring in RGB stars (Placco et al. 2014). On the other hand, $F_{\text{CEMP}}^{\text{halo}}$ might decrease when accounting for non-local thermodynamic equilibrium (non-LTE) effects, which can drastically lower the measured $[C/Fe]$ value (Amarsi et al. 2019). Finally, the selection function of the data might suffer some bias in either direction. Since the main two effects compensate and non-LTE corrections are only available for a handful of stars (Amarsi et al. 2019), we will consider the derived $F_{\text{CEMP}}^{\text{halo}}$ as the reference fraction of Galactic halo stars at different iron abundances. It is also worth mentioning some work suggesting that the frequencies of CEMP varies with location in the Galactic halo (Yoon et al. 2018; Lee et al. 2019). Despite the mentioned caveats, as an initial hypothesis, we will assume that the fraction of CEMP-no stars is the same in all environments and equal to that of the Galactic halo, $F_{\text{CEMP}} \equiv F_{\text{CEMP}}^{\text{halo}}$.

## 2.2 Ancient environments in the Local Group

To understand whether CEMP-no stars might simply be hidden in the Galactic bulge, we need to derive the probability to catch CEMP-no stars while blindly observing the bulge and then compare to other ancient environments in the Local Group (Salvadori et al. 2015). To this end we analyzed a number of environments with increasing stellar mass (luminosity) and correspondingly an increasing average stellar metallicity. In particular we considered (see Table 1 for details): the least luminous UFDs, the most luminous UFD *Bootes I*, the dSph galaxy *Sculptor*, the dSph galaxy *Fornax*, and the Galactic bulge.

The carbon measurements in the least luminous UFDs are from high-resolution spectroscopic studies (Segue 1: Norris et al. (2010); Frebel et al. (2014). Pisces II: Spite, M. et al. (2018). Ursa Major II and Coma Berenice: Frebel et al. (2009). Leo IV: Simon et al. (2010)). Only few stars have been observed in each UFD because these galaxies are faint and distant so that only few RGB stars are available for spectroscopic observations. For this reason, we consider stars belonging to the faintest UFDs all together. Stars in *Boötes I*, the most luminous UFD, are studied with low-resolution spectroscopy (Norris et al. 2010; Lai et al. 2011) with the only exception of seven

stars (Gilmore et al. 2013). In the *Sculptor* dSph galaxy, many carbon measurements are available from both low- (Kirby et al. 2015) and high-resolution spectroscopic studies (Frebel et al. 2010; Tafelmeyer et al. 2010; Starkenburg et al. 2014; Jablonka et al. 2015; Simon et al. 2015; Skúladóttir et al. 2015). Finally, for the *Fornax* dSph we include the low-resolution measurements from Kirby et al. (2015).

For the Galactic bulge, we use data from the JINAbase and those provided by the high-resolution, large multi-object APOGEE spectroscopic survey (Data Release 16, Jönsson et al. 2020).

## 2.3 Are CEMP stars hidden in the Galactic bulge?

For a 'blind' survey that does not pre-select the most metal-poor stars, the ability to find CEMP stars is naturally limited by the number of stars that exist at each $[Fe/H]$, i.e. the metallicity distribution function (MDF; Salvadori et al. 2015). In Fig. 2 (left), we show the normalized MDFs for the different environments studied. Note that the total number of stars observed in these environments strongly increases with galaxy luminosity: only 19 stars constitute the MDF in the faintest UFDs, while for the Galactic bulge we have >17 000 stars.

We then ask ourselves: if we assume that the fraction of CEMP stars is 'Universal' and equal to that of the Galactic halo, $F_{\text{CEMP}}^{\text{halo}}$, what is the joint probability, $P_{\text{obs}}$, to observe a star that has a given $[Fe/H]$ value, which is also carbon-enhanced? This is done by combining two independent functions: the halo fraction of CEMP stars and the normalized MDF derived for each environment:

$$P_{\text{obs}}([\text{Fe/H}]) = F_{\text{CEMP}}^{\text{halo}}([\text{Fe/H}]) \times \frac{N_{*}([\text{Fe/H}])}{N_{\text{tot}}}. \quad (2)$$

In the right-hand panel of Fig. 2, we show $P_{\text{obs}}$ for the different environments studied, and we note a clear trend: the brighter is the galaxy, and thus more shifted towards higher $[Fe/H]$ values, the lower is the overall probability to observe a CEMP star. In other words, in the most massive galaxies, stars in the low-$[Fe/H]$ tail become almost 'invisible' with respect to the metal-rich stars since they are much rarer. This is especially evident in the bulge since the peak of the resulting probability ($P_{\text{obs}} < 0.01$ per cent) is two orders of magnitude lower than the values obtained in UFDs ($P_{\text{obs}} \sim 3$ per cent), where more CEMP stars have actually been discovered. This peak of $P_{\text{obs}}$ in the bulge, at $[Fe/H] \simeq -3.5$, is almost flat for $-3.5 < [Fe/H] < -2.0$, and drops significantly at $[Fe/H] < -3.5$.

In conclusion, the scarcity of stars at low-$[Fe/H]$, makes it very difficult to search for CEMP stars in luminous environments, such as the Galactic bulge, that are dominated by metal-rich stars. The dearth of CEMP stars in the Galactic bulge can therefore partially be ascribed to this reason. However, neither EMBLA or PIGS are 'blind' surveys since they specifically search for the most metal-poor stars in the Galactic bulge. PIGS, in particular, obtained a sample of 1900 stars with $[Fe/H] < -2$. Thus, if we hypothesize that the CEMP fraction within the bulge is equal to the one we obtained for the Galactic halo, namely $F_{\text{CEMP}}^{\text{bulge}}([Fe/H] < -2) = F_{\text{CEMP}}^{\text{halo}}([Fe/H] < -2) \approx 24$ per cent, then we can estimate that $\sim 450$ CEMP stars
<figure><img src="image_3.png" /><figcaption>**Figure 2.** The MDF (left-hand panel) and the probability of observing a CEMP star (right-hand panel) for the different environments displayed in order of increasing luminosity (mass): ultrafaint galaxies (orange), Boötes (red), Sculptor (pink), Fornax (green), and the Galactic bulge (blue).</figcaption></figure>

should have been observed by PIGS. Instead, <62 CEMP(-no) candidates have been identified so far (Arentsen et al. 2021). In agreement with the PIGS papers, we thus conclude that the dearth of CEMP stars in the Galactic bulge cannot be only a consequence of a statistical effect in this metal-rich environment. It is also interesting that we derive a peak probability of observing CEMP stars in the Galactic bulge at [Fe/H] ≃ −3.5, i.e. at the metallicity of the only CEMP star with high [C/Fe] observed by Arentsen et al. (2021).

## 3 COSMOLOGICAL MODEL DESCRIPTION

In this section, we will briefly recap the cosmological model used to identify the first star-forming haloes whose stars are currently located in the Galactic bulge and to investigate their properties. The model combines a N-body simulation that follows the hierarchical assembly of a Milky Way (MW) – like galaxy (used also in Salvadori et al. 2010; Pacucci et al. 2017) with a semi-analytical model (GAMETE; Salvadori et al. 2007, 2015) that follows the evolution of baryons, from the formation of the first stars to the present day. The semi-analytical model allows us to follow the star formation and metal enrichment history of the Milky Way from early times ($z \sim 20$) until now ($z = 0$) and thus to link the chemical abundances of present-day stars with the properties of the first stellar generations.

(i) **N-body simulation:** The N-body simulation used to study the hierarchical formation of the MW has a low-resolution region corresponding to a sphere of radius $10 h^{-1}$ Mpc (see Scannapieco et al. 2006). The region of the highest resolution is a sphere of radius $1 h^{-1}$ Mpc, i.e. four times the virial radius of the MW at $z = 0$ ($r_{vir} = 239$ kpc). A low-resolution simulation including gas physics and star formation has been used to confirm that the initial conditions will lead to a disc galaxy like the Milky Way. The system consists on about $10^6$ dark matter (DM) particles within $r_{vir}$ with masses of $7.8 \times 10^5$ M$_\odot$; its virial mass and radius are respectively $M_{vir} = 7.7 \times 10^{11}$ M$_\odot$ and $r_{vir} = 239$ kpc, roughly consistent with the observational estimates of the MW, $M_{vir} \approx 10^{12}$ M$_\odot$ (e.g. McMillan 2011). The softening length is 540 kpc. The simulation data is output every 22 Myr between $z = 8-17$ and every 110 Myr for $z < 8$. At each output, the virialized DM haloes have been identified using a *friend-of-friend* group finder with a linking parameter $b = 0.15$, and

a threshold number of particles constituting virialized haloes equal to 50. The N-body simulation enables to reconstruct a hierarchical history of the MW that proceeds through the consecutive merging of DM haloes maintaining the information about the spatial distribution of DM particles belonging to them.

(ii) **Star formation:** At the initial redshift of the simulation, $z = 20$, gas in DM haloes is assumed to have a primordial composition, and only objects with virial temperatures $T_{vir} \ge 10^4$ K form stars. This choice, which is equivalent to assuming that the star formation activity is rapidly quenched in mini-haloes, is dictated by the limited DM resolution of the N-body simulation (Graziani et al. 2015). At each time-step, stars are assumed to form in a single burst, proportional to the available cold gas mass, $M_{gas}$. The constant of proportionality is a redshift-dependent star-formation efficiency $f_*(z) = \epsilon_* \Delta t(z)/t_{ff}$, where $t_{ff}$ is the free-fall time, $\Delta t(z)$ the N-body time-step and $\epsilon_*$ a free parameter of the model, physically corresponding to a 'local' star formation efficiency.

(iii) **Pop III stars:** The model simply assumes that the Pop III IMF is a delta function centered either on the average mass value of the pair-instability SNe (PISNe; Salvadori et al. 2010) or faint SNe (Pacucci et al. 2017), i.e:

$$ \Phi(m) = \frac{dN}{dm} \propto \delta(m_{\text{Pop III}}). \qquad (3) $$

To give an idea of the amount of chemical elements released by these different Pop III stars, an 'average' PISN ($m_{\text{Pop III}} = 200$ M$_\odot$) yields of metals, iron, and carbon, respectively, is equal to: $Y = M_Z/M_{\text{Pop III}} = 0.45$, $Y_{\text{Fe}} = 0.022$, and $Y_C = 0.02$ (Heger & Woosley 2002). Conversely, an average faint SN ($m_* = 25$ M$_\odot$) gives $Y \sim 0.1$, $Y_{\text{Fe}} = 4 \times 10^{-7}$, and $Y_C = 9.98 \times 10^{-3}$ (Iwamoto et al. 2005; Marassi et al. 2014).

(iv) **Pop III-to-Pop II transition:** Following the critical metallicity scenario (e.g. Bromm et al. 2001; Schneider et al. 2002), we assume that the IMF of newly formed stars depends upon the initial metallicity of the star-forming clouds. Therefore, a star forming halo with a gas metallicity $Z \le Z_{cr}$ will form Pop III stars, otherwise, if $Z > Z_{cr}$, it will form Pop II/I stars. Exploiting the data driven constraints of De Bennassuti et al. (2017), we set $Z_{cr} = 10^{-4.5} Z_\odot$. Normal Pop II/I stars are assumed to have masses in the range [0.1, 100] M$_\odot$ and to form according to a standard Larson IMF (Larson 1998), which
peaks at the characteristic mass $m_{ch} = 0.35M_{\odot}$ and rapidly declines with a Salpeter-like shape towards larger masses (Salpeter 1955).

(v) **Metal enrichment:** The newly formed stars are assumed to instantaneously evolve at the following snapshot of the simulation (Instantaneous Recycling Approximation, IRA) since the time elapsed between two neighbouring steps is larger than the lifetime of the least massive stars evolving as SNe. For Pop III stars exploding as PISN we assume the yields of Heger & Woosley (2002), while for primordial faint SN we assume those of Iwamoto et al. (2005). For Pop II/I stars we adopt the yields from Woosley & Weaver (1995). To not overestimate the contribution of carbon due to AGB stars that also produce slow neutron-capture process elements and, even if not in binary systems, can lead to the formation of 'moderate' CEMP-s stars (Rossi et al. 2023), we only considered the chemical products of Pop II/I that evolve as core-collapse SNe in short time-scales (3–30 Myr). After SN explosions the newly produced/injected metals are assumed to be instantaneously and homogeneously mixed within the inter stellar medium (ISM) eventually reaching the inter galactic medium (IGM) via supernova driven outflows.

(vi) **Gas and metals dispersal:** Supernovae may release a high amount of energy, which may overcome the binding energy of the hosting halo leading to a partial gas and metals removal from the galaxy. The mass of gas ejected into the IGM, $M_{\text{gas}}^{\text{ej}}$, depends on the balance between the escape velocity of the halo, $v_{\text{esc}}$, and the kinetic energy released during the explosion, namely:

$$M_{\text{gas}}^{\text{ej}} = (2E_{\text{SN}})/v_{\text{esc}}^2, \qquad (4)$$

where $E_{\text{SN}} = \epsilon_w N_{\text{SN}} \langle E_{\text{SN}} \rangle$, with $N_{\text{SN}}$ number of SN explosions and $\langle E_{\text{SN}} \rangle$ the average explosion energy. For a typical PISN $\langle E_{200M_{\odot}}^{\text{PISN}} \rangle \sim 2.7 \times 10^{52}$ erg (Heger & Woosley 2002), while a typical faint SN provides $\langle E_{25M_{\odot}}^{\text{faint}} \rangle \sim 0.7 \times 10^{51}$ erg (Iwamoto et al. 2005; Marassi et al. 2014), which is lower than the energy released by a normal $25 M_{\odot}$ core-collapse SN, $E_{25M_{\odot}}^{\text{cc}} \approx 10^{51}$ erg. The quantity $\epsilon_w$, representing the second free parameter of the model, is the wind efficiency, that is the fraction of the explosion energy converted into kinetic form.

(vii) **Chemical evolution:** Due to mechanical feedback, the mass of gas and metals in a halo can decrease substantially. At each simulation time-step the gas mass reservoir in each halo, $M_{\text{gas}}$, is updated with respect to the initial gas mass, $M_{\text{gas}}^{\text{in}}$, to account for the mass of stars locked into newly formed stars, $M_{*}$, and the gas mass ejected out of the halo, $M_{\text{gas}}^{\text{ej}}$:

$$M_{\text{gas}} = M_{\text{gas}}^{\text{in}} - (1+R)M_{*} - M_{\text{gas}}^{\text{ej}}, \qquad (5)$$

where R is the returned fraction, which is equal to 1 only for PISN and lower than unity otherwise. Similarly, the mass of metals, $M_{Z}$, in each hosting halo is updated as follows:

$$M_{Z} = M_{Z}^{\text{in}} + Y M_{*} - Z_{\text{ISM}}^{\text{in}} M_{*} - Z_{\text{ISM}}^{\text{in}} M_{\text{gas}}^{\text{ej}}, \qquad (6)$$

where $M_{Z}^{\text{in}}$ is the initial mass of metals, $Y_{Z}$ the metal yield, and $Z_{\text{ISM}}^{\text{in}}$ the initial metallicity of the ISM.

(viii) **Model calibration:** To set the best values of the two model free parameters, $\epsilon_{*}$ and $\epsilon_{w}$, the observed properties of the MW have been used as a benchmark. In particular, the results of the simulations at redshift $z=0$ have been compared with the gas/stellar mass and metallicity of the MW, the baryon-to-dark matter ratio, and the metallicity of high-velocity clouds (see Salvadori et al. 2007, 2010, for more details).

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Present-day mass fraction of CEMP stars with [C/Fe] > +2, relative to all stars at [Fe/H] < −2, contained in an annulus of radial width within 1 kpc (see colour map) plotted in the $\rho-\zeta$ cylindrical coordinate plane of our simulated MW-analogue galaxy. The results have been obtained from our semi-analytical model by assuming that all Pop III stars have $m_{\text{Pop III}} = 25M_{\odot}$ and explode as faint SNe.</figcaption></figure>

# 4 SIMULATION ANALYSIS AND RESULTS

## 4.1 Extreme cases of Pop III enrichment

To understand how much the properties of Pop III stars can affect the present-day distribution of CEMP stars, we first analyse the $z=0$ outputs of the simulation for two extreme cases: (i) when all Pop III stars explode as an average PISN of $200 M_{\odot}$, or (ii) when they all explode as an average faint SN of $25 M_{\odot}$ (see Section 3). In the first case, we find that CEMP stars are never produced by the model, the carbon-to-iron value of present-day stars being [C/Fe] $\lesssim 0.2$. This is perfectly consistent with the low [C/Fe] value of an ISM imprinted by massive $> 160 M_{\odot}$ Pop III stars evolving as PISN (e.g. Salvadori et al. 2019, left-hand panel of fig. 2).

Conversely, CEMP stars are formed when all Pop III stars are assumed to explode as $25 M_{\odot}$ faint SNe and they can have huge C-excess, $[\text{C/Fe}] \ge +4$, at the lowest $[\text{Fe/H}] < -4$. This result is in line with the idea that the most iron-poor CEMP(-no) stars are indeed the descendants of this kind of zero-metallicity SNe (Section 1). Under this assumption, we can also compute the present-day distribution of CEMP stars which predicts that at $[\text{Fe/H}] < -2$ the innermost region hosts the largest mass fraction of CEMP stars. This remains true even if we only consider CEMP stars with $[\text{C/Fe}] > +2.0$, see Fig. 3.

Although we cannot completely rule out a statistical effect (see Section 2.3), these results are at odds with observations, which actually suggest a lower fraction of CEMP stars in the bulge than in the outer regions (Section 2.1). Therefore it is clear that to model the early chemical enrichment in the first star-forming haloes it is necessary to consider a more realistic Pop III IMF, which accounts for the contribution of both faint SNe and PISNe.

## 4.2 The first star-forming haloes

To identify the first star-forming haloes and see where the oldest stars are currently located, we reconstructed the hierarchical tree in the MW-analogue by combining the positions of the DM particles within which star formation occurs, and their belonging to a specific halo at different redshifts. A halo is identified as the main progenitor of another at the following time-step, if the latter inherits at least 90 per cent of its DM particles. If a halo has two or more progenitors then we can assume that a merging process has occurred.
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Present-day positions, $z = 0$, in the $\rho-\zeta$ cylindrical coordinate plane of Pop II stars formed in a specific halo of our simulation (id number and colour in label), at the specified $z$ which decreases from the top to the bottom panel. Identified Galactic regions are: the inner halo, 7 kpc $< r \le 20$ kpc (grey circles); and the bulge, $r \le 2.5$ kpc (green circle).</figcaption></figure>

Fig. 4 shows the present-day, $z = 0$, positions of Pop II stars formed in the first star-forming haloes, at $z \ge 13$. As we can see from the top panel of Fig. 4, the oldest Pop II stars are currently located in the innermost Galactic region, i.e. the bulge. The position

of DM particles belonging to haloes that begin forming stars at lower redshifts (lower panels of Fig. 4) is gradually spread towards the outermost regions, filling first the inner and then the outer Galactic halo. These results confirm the idea that the most ancient stellar populations, some of which should have been formed from the ashes of the first stars, must dwell in the bulge (see Section 1).

Fig. 5 shows the present-day, $z = 0$, spatial distribution of DM particles, colour coded with their SFR during the first $\approx 800$ Myr, i.e. averaged over the cosmic time between $z = 15$ and $z = 8$, which corresponds to the lowest redshift for which we have Pop III star formation. Thus we assigned to each DM particle the mean SFR of the halo it belongs to, so that:

$$ \langle \text{SFR} \rangle = \frac{\sum_{i=1}^{N} M_{\text{si}}^{\text{tot}}}{t(z=8) - t(z=15)}, \qquad (7) $$

where $M_{\text{si}}^{\text{tot}}$ is the total stellar mass formed in the halo at the $i$th time-step, and $N$ is the number of time-steps between $z = 15$ and $z = 8$. Then we computed the final $\langle \text{SFR} \rangle$ by averaging among the results of all DM particles in the considered pixel. This results in a wide range of average star formation rates (Fig. 5), $\langle \text{SFR} \rangle \approx (10^{-4}-1.4) M_{\odot}/\text{yr}$.

The spatially resolved region with the highest mean SFR at $z > 8$ is the innermost one, i.e. the Galactic bulge. This result, which is in agreement with previous theoretical (e.g. Cescutti & Matteucci 2011) and observational (e.g. Lucertini et al. 2022) findings, suggests that the progenitors of the Galactic bulge experience, on average, a more intense star formation at early times with respect to the progenitors of the Galactic halo. A higher star formation rate in a metal-poor environment might have strong consequences: in particular, it could allow to produce many rare stellar populations, most likely including the massive progenitors of PISNe (Rossi et al. 2021).

## 4.3 Pop III enrichment in the first star-forming haloes

We will now analytically model the Pop III star enrichment within the first star-forming haloes of the cosmological simulation by using more realistic Pop III IMFs, and by accounting for their incomplete sampling. To make our analytical calculations, we selected those pristine and star-forming haloes for which more than 80 per cent of particles are predicted to dwell in the Galactic bulge, i.e. with $r = \sqrt{\rho^2 + \zeta^2} \le 2.5$ kpc, at the present time (8 haloes). Their halo mass and redshift of formation range between $m_h = (5.5 \times 10^7 - 1.7 \times 10^8) M_{\odot}$, and $11.4 < z < 16.3$. The predicted SFR for these first star-forming haloes is $SFR < 1.8 \times 10^{-2} M_{\odot}\text{yr}^{-1}$, which in all cases is less than the threshold value for a fully populated Pop III IMF, $SFR_{\min} \sim 10^{-1} M_{\odot}\text{yr}^{-1}$ (Rossi et al. 2021). Thus, it is fundamental to include incomplete sampling of the adopted Pop III IMFs.

### 4.3.1 Realistic IMFs of Pop III stars

In addition to the delta functions considered so far (Section 4.1), we will assume that Pop III stars form according to a Larson IMF:

$$ \Phi(m) = \frac{dN}{dm} \propto m^{-2.35} \exp(-m_{ch}/m), \qquad (8) $$

biased towards more massive stars as it is expected for the first stellar generations. In particular, following the latest data-driven results from stellar archaeology (Rossi et al. 2021), we assume a minimum mass of Pop III stars equal to $m_{\min} = 0.8 M_{\odot}$, a maximum mass $m_{\max} = 1000 M_{\odot}$, and we explore different characteristic masses, $m_{ch} \in \{1, 10, 100 M_{\odot}\}$. This allows us to account for the contribution to chemical enrichment of Pop III stars exploding as both faint SNe and PISNe and to vary their relative proportions. All the Pop III IMFs considered in this study are shown in Fig. 6.
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Present-day spatial distribution of DM particles in the cylindrical coordinate plane of the simulated galaxy, colour coded with their SFR averaged over the cosmic time between $z = 15$ and $z = 8$.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 6.** The Pop III IMFs used in this work: delta functions centered on the average mass values of faint SNe ($m_{\text{Pop III}} = 25\ M_{\odot}$; orange) and PISNe ($m_{\text{Pop III}} = 200\ M_{\odot}$; red), Larson IMFs with $m_{\text{min}} = 0.8\ M_{\odot}$, $m_{\text{max}} = 1000\ M_{\odot}$, and $m_{\text{ch}} \in \{1, 10, 100\ M_{\odot}\}$ are, respectively, in cyan, green, and magenta.</figcaption></figure>

For each of our selected first star-forming halo, we derived the *effective* Pop III IMF associated to each burst of primordial star-formation to account for the stochastic and incomplete IMF sampling. This was done using the Monte Carlo procedure developed by Rossi et al. (2021). Starting from the total mass of Pop III stars formed, we

generate a random number of stars, which are distributed according to the assumed stellar IMF. Due to the stochastic nature of this sampling, every time that stars are formed, the effective stellar mass distribution is differently populated, especially at the higher masses (Fig. 7). Following Rossi et al. (2021) we assumed a time-scale for star formation equal to $\Delta t = 1$ Myr and then we computed the total mass of Pop III stars formed in the burst accordingly, $M_{\text{Pop III, burst}} = \text{SFR} \times \Delta t$. Then we applied a statistical approach by averaging among the results of 50 incomplete random sampling of the IMF and quantifying the scatter among them.

Fig. 7 shows, for different $m_{\text{ch}}$, the comparison between the theoretical Pop III IMFs and effective ones, which have been obtained in one run of the random sampling. The mass range covered by faint SNe, (8–40) $M_{\odot}$; PISNe, (140–260) $M_{\odot}$; and other stars that do not end their lives as supernovae are identified using different colours. Among these, stars with $m < 8\ M_{\odot}$ lose their external envelope during the AGB phase, and stars in the two mass ranges, (40–140) $M_{\odot}$ and (260–1000) $M_{\odot}$ are predicted to directly end their life collapsing in a black hole (Heger & Woosley 2002). Note that the average mass fraction of stars exploding as faint SNe (and PISNe) strongly varies, not only with $m_{\text{ch}}$ but also as a consequence of the incomplete IMF sampling. In the mass range of faint SNe, the Pop III IMFs are almost completely sampled except in the extreme case $m_{\text{ch}} = 100\ M_{\odot}$. On the contrary, in the typical mass range of PISNe the Pop III IMF is only partially populated. As expected, however, the mass range of PISNe becomes more densely populated as the characteristic mass increases.

The mass of metals, iron, and carbon produced by the total number of Pop III stars formed was computed by summing up the contribution of Pop III stars with different masses, i.e. $M_{\text{X}} = \sum_{i}^{N} Y_{\text{X}}(m_{\text{Pop III, i}})m_{\text{Pop III, i}}$, where $Y_{\text{X}}(m_{\text{Pop III, i}})$ is the yield of the
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Theoretical (solid lines) and effective (coloured histogram) Larson IMFs of Pop III stars with increasing characteristic masses from top to bottom ($m_{\text{ch}} = 1, 10, 100 M_{\odot}$) obtained from one run of the random sampling. In each panel, the mass range covered by faint SNe, (8–40) $M_{\odot}$ (orange); PISNe, (140–260) $M_{\odot}$ (red); and other stars that do not end their lives as supernovae (grey), are specified. The total mass fractions are also listed for each mass range. The total mass of the burst of Pop III is $7.7 \times 10^3 M_{\odot}$.</figcaption></figure>

element X produced by Pop III stars in the *i*-th mass bins. For faint SNe, we followed De Bennassuti et al. (2017) and assumed that the yield corresponding to 25 $M_{\odot}$ is simply re-scaled to the mass of the Pop III star exploding as faint SN, i.e. $Y_X(m_{\text{Pop III}}) = Y_X(25 M_{\odot}) \times (m_{\text{Pop III}}/25 M_{\odot})$. This scaling is quite a good approximation as shown by Marassi et al. (2014) and it has been used in different works investigating stellar nucleosynthesis (e.g. Nomoto et al. 2013; De Bennassuti et al. 2017). For PISNe we exploited the yields and SN explosion energies provided by Heger & Woosley (2002). Similarly, we computed the total amount of energy released by SNe by accounting for the number of Pop III stars effectively formed in different mass bins, $E_{\text{SN}} = \epsilon_w \sum_{i}^{N} N_{\text{SN}}(m_{\text{Pop III}}) \times E_{\text{SN}}(m_{\text{Pop III}})$, as done in Rossi et al. (2021). Finally, by exploiting equation (5) we computed the ISM metallicity along with [Fe/H] and [C/Fe] after the Pop III star enrichment.

### 4.3.2 Average and maximum [C/Fe] of bulge stars

In Fig. 8, we show the [Fe/H] and [C/Fe] of the ISM in the 8 first star-forming haloes currently dwelling in the Galactic bulge as obtained with our incomplete IMF sampling procedure (Section 4.3.1) and as a function of the five different Pop III IMFs assumed (see Fig. 6). The abundance ratios for the different haloes are colour coded according to the halo mass and, for comparison, we also show the results obtained considering the fully sampled IMFs as solid grey lines. For each $m_{\text{ch}}$, the average abundances over all haloes are shown as empty black symbols together with their errors. The observed ranges of [C/Fe] within the Galactic halo and bulge (from Fig. 1) are highlighted as well in Fig. 8 (respectively, as grey and orange shaded areas).

First, we note that for all IMFs, the total gas metallicity after the Pop III star enrichment is $Z_{\text{ISM}} > Z_{\text{cr}} = 10^{-4.5}Z_{\odot}$, which implies that normal low-mass long-lived stars will be able to form in such environments. Furthermore, when we only account for the contribution of faint SNe, we get a very low value of [Fe/H] $\simeq -5.6$ and an extremely high [C/Fe] $\simeq +4$, which is close to the maximum value observed in the Milky Way halo (Keller et al. 2014). As soon as the chemical contribution of PISNe is also considered, namely for $m_{\text{ch}} = 1 M_{\odot}$, the [C/Fe] value drops dramatically by at least 2 orders of magnitude. As we can see at fixed $m_{\text{ch}}$, as the halo mass increases, the abundances approach the fully sampled case, that is, [Fe/H] increases and [C/Fe] decreases. This result can easily be explained as a consequence of the Pop III IMF sampling: in the lower mass haloes, which have lower star formation rates, on average only one PISN explodes, thus only partially lowering the [C/Fe] value obtained in the case of faint SNe only. Still, even a single PISN is able to inject into the ISM 50 per cent of its total stellar mass in form of heavy elements, thus strongly affecting the final [C/Fe] value of the ISM. In the most massive halo more PISNe can actually form and explode, thus further lowering the expected [C/Fe] value. For increasing characteristic mass, the mean [C/Fe] values become even lower, since PISN production starts to dominate the chemical enrichment. When only PISNe explode, no CEMP stars are able to form, regardless of the halo mass.

In conclusion, the measured [C/Fe] values in the Galactic bulge suggest that massive Pop III stars exploding as PISNe have likely formed in this environment, washing out the high [C/Fe] signature left by low-energy primordial SNe. Furthermore, by comparing our predicted [C/Fe] with the values observed in the bulge (Fig. 8), we can exclude the two extreme IMF cases, where Pop III stars explode
halo, then biased surveys specifically searching for the most metal-poor stars in the Galactic bulge – such as the PIGS survey (Arentsen et al. 2021) – should have observed many more CEMP stars (~450) than what they actually found (<62, Section 2). This discrepancy could be due to selection biases in the photometric selection, but could there be another mechanism that reduces the [C/Fe] values in the bulge?

(ii) We then asked whether the lower [C/Fe] values of CEMP stars in the Galactic bulge reflects the different formation and evolution of this ancient environment. To address this question, we focused on the predictions derived from the ΛCDM cosmological model, through the use of a *N*-body simulation that follows the hierarchical formation of a MW-like galaxy combined with the semi-analytical model GAMETE (e.g. Salvadori et al. 2007, 2010, 2015; Graziani et al. 2015; Pacucci et al. 2017), which is required to follow the star formation and metal enrichment history. The model is data-calibrated, i.e. the best values of the free parameters are set to reproduce the observed properties of the MW at $z = 0$ (Section 3). If we assume that all Pop III stars explode as faint SNe, we find that the mass fraction of CEMP stars with [C/Fe] $> +2$ increases at decreasing Galactocentric radii and it is maximum in the Galactic bulge, which is the region of the simulation containing the most ancient stars, which is at odd with observations. However, the *N*-body simulations reveal that the stars dwelling into the present-day bulge form in haloes that experienced the highest mean star formation rate at high redshifts ($z > 8$).

We inferred that the dearth of CEMP stars with [C/Fe] $> +2$ in the Galactic bulge might be linked with the higher star formation rate of its early progenitor haloes, which hosted the first stars. Indeed, star formation rates $> 10^{-2} M_{\odot} yr^{-1}$ in primordial environments might allow the formation of rare very massive Pop III stars (Rossi et al. 2021), which evolve as energetic pair instability SNe (PISNe).

(iii) We thus investigated how the chemical enrichment of the bulge progenitors depends upon the properties of rare Pop III stars that can effectively form in these highly star-forming systems. To this aim, we performed analytical calculations and computed the chemical properties of the ISM in the bulge progenitors after the contribution of Pop III stars. More specifically, we investigated how different Pop III IMFs affect the [C/Fe] and [Fe/H] values of the ISM by assuming various characteristic masses ($m_{ch} = 1, 10, 100 M_{\odot}$), by including the contribution of Pop III stars exploding as PISNe (see Section 4), and by accounting for the incomplete sampling of the Pop III IMF (Rossi et al. 2021). Our results show that very massive Pop III stars can effectively form in the bulge progenitors, and that their contribution to the chemical enrichment as energetic PISNe can partially wash out the distinctive signature of faint SNe, lowering the carbon overabundance down to [C/Fe] $< +2$. In particular, we show that the higher the probability to form very massive Pop III stars, i.e. the larger the $m_{ch}$, the lower is the [C/Fe] value of the ISM after the contribution of Pop III stars. By exploiting the available data we thus tentatively infer $1 M_{\odot} \lesssim m_{ch} \lesssim 10 M_{\odot}$, consistent with the constraints found by Rossi et al. (2021) based on ultrafaint dwarfs. We conclude that the modest [C/Fe] values of CEMP stars identified in the bulge, [C/Fe] $\approx +0.8$, along with the dearth of CEMP stars with [C/Fe] $> +2$ could be an indirect probe of very massive first stars exploding as PISNe, which are extremely rare and hence can only form in highly star-forming progenitors of the MW bulge.

Through careful analysis of observational data and theoretical simulations, we thus suggest that the dearth of CEMP stars in the bulge might be intrinsic, and not only a consequence of systematic observational effects. Furthermore we have shown that the first star-

forming haloes that end up in the bulge have typically significantly higher star formation rates than those in the outskirts of the galaxy. This means that the Pop III IMF is better sampled in these systems, enabling the formation of rare populations such as very massive first stars that explode as PISNe. Ultimately, our analysis showed this to be a very plausible explanation for the dearth of CEMP stars in the bulge, as added PISN contribution will lower the observed [C/Fe]. Furthermore we suggested a new promising method that exploits the lack of CEMP stars to constrain the characteristic mass of the first stars.

# 6 DISCUSSION AND FUTURE OUTLOOK

Our findings provide a plausible explanation for the lack of CEMP stars with high [C/Fe] values in the Galactic bulge, which is not only linked to statistical issues or observational biases, but it is rather a consequence of the different formation path of this ancient but star-rich environment. The idea that Pop III stars exploding as PISNe can efficiently form in the bulge progenitors and lower the [C/Fe] value in their ISM is extremely appealing, since it can also explain the larger fractions of CEMP stars observed in the Galactic halo and in UFDs. The typical star formation rates observationally inferred for UFDs ($\le 10^{-3} M_{\odot} yr^{-1}$, e.g. Salvadori et al. 2014; Gallart et al. 2021) are indeed too low for enabling these galaxies to form very massive Pop III stars (Rossi et al. 2021). The early chemical enrichment stages of these small systems are thus likely dominated by faint SNe, whose explosion energy ($E < 10^{51}$ erg) is small enough to not exceed the halo binding energy (Rossi et al., in preparation) and whose chemical products are characterized by a large amount of carbon and a tiny amount of iron. The progenitors of UFDs have been suggested to build-up the low-Fe tail of the Galactic halo (e.g. Salvadori et al. 2015; Bonifacio et al. 2021), where CEMP-no stars are found, making the overall scenario perfectly consistent.

Other sources that can partially wash-out the chemical signature of faint SN are primordial hypernovae, i.e. Pop III stars with masses $m_{PopIII} \sim (10-100) M_{\odot}$ that experience energetic SN explosions, $E \sim 10^{52}$ erg (Heger & Woosley 2010). In the last years, the descendants of primordial hypernovae have been firstly identified thanks to their sub-solar [C/Fe] values in nearby dwarf spheroidal galaxies and halo stars (Placco et al. 2021; Skúladóttir et al. 2021). However, since Pop III stars exploding as faint SNe and as hypernovae have the same range of masses, we expect them to form with the same efficiency in both the progenitors of the Galactic halo, UFDs, and the bulge. Furthermore, the same reasoning can be applied to the contribution of Pop II stars in the chemical enrichment. In conclusion, this alternative solutions to explain the lower [C/Fe] value of bulge stars, is hard to reconcile with the global scenario. In a forthcoming study, we will further explore the contribution of different Pop III stars to the bulge enrichment by performing a self-consistent calculation, i.e. by following within the *N*-body simulation the chemical elements produced by Pop III stars with different masses and explosion energies (Koutsouridou et al., in preparation). In the same study we will relax the Instantaneous Recycling Approximation, which, although quite robust for $z > 6$, may have affected our results. In fact, taking into account the lifetimes of Pop III stars, the more massive PISNe would have exploded first – raising Fe and keeping C low – and then the faint SNe would have contributed – raising C with constant Fe. Consequently among the older Pop II stars we could have stars with higher [Fe/H] and lower [C/Fe].

Ultimately, the dearth of CEMP stars with high [C/Fe] values in the Galactic bulge might provide an indirect probe for the long-searched very massive Pop III stars evolving as PISNe. A way to prove that primordial PISNe exist is to search for an
underabundance of [Zn/Fe] and [Cu/Fe] in their descendants (Salvadori et al. 2019; Aguado et al. 2023a). Interestingly, stars with sub-solar [Zn/Fe] have been predominantly identified in the bulge (e.g. Barbuy et al. 2015; Duffau et al. 2017) and classical dwarf spheroidal galaxies (e.g. Skúladóttir et al. 2017), which are indeed more massive and star-forming than UFDs. Furthermore, in the nearby future the WEAVE (Dalton 2016), 4MOST (Bensby et al. 2019; Chiappini et al. 2019; Christlieb et al. 2019), and J-PLUS/S-PLUS (Cenarro et al. 2019; Mendes de Oliveira et al. 2019) Galactic Archaeology surveys will allow us to derive a large number of carbon measurements in both halo and bulge metal-poor stars, as well as the satellite dwarf galaxies (4DWARFS; Skúladóttir et al. 2023). Thus, the CEMP-no frequency will be established with unprecedented accuracy, giving new insights into PISN events. Unfortunately, Zn will not be measured by 4MOST. However, if the green grating will be used in the high-resolution WEAVE survey, a large sample of Cu and Zn measurements will be extremely valuable to confirm this scenario and to search for PISNe descendants.

This is the perfect time to deepen our studies on Stellar Ar-
chaeology. The future of the field is indeed promising since great
technology developments and new instruments will provide us
unprecedented data sets to confirm our theoretical findings. In the
coming years, there will be a surge of observations in the inner
Galaxy of ancient, metal-poor stars thanks to the 4MOST Milky Way
Disc and BuLgE Low-Resolution (4MIDABLE-LR; Chiappini et al.
2019) and 4MOST Milky Way Disc And BuLgE High-Resolution
surveys (4MIDABLE-HR; Bensby et al. 2019). This will provide a
completely new view of the central regions of the Milky Way, and,
hopefully, greatly advance our understanding of the first stars.

ACKNOWLEDGEMENTS

SS, MR, DA, IK, and ÁS acknowledge funding from the European Research Council (ERC) under the European Unionion's Horizon 2020 research and innovation programme, project NEFERTITI (grant agreement No 804240). SS also acknowledges funding from the PRIN-MIUR17, The quest for the first stars, prot. n. 2017T4ARJ5. The authors acknowledge the SDSS/APOGEE survey.

DATA AVAILABILITY

The authors confirm that the data analysed in this study are available
within the JINABASE (Abohalima & Frebel 2018), the APOGEE
spectroscopic survey (DR16; Jönsson et al. 2020), and the following
papers: Norris et al. (2010, 2010), Frebel et al. (2010, 2014, 2009),
Spite et al. (2018), Simon et al. (2010, 2015), Lai et al. (2011),
Gilmore et al. (2013), Kirby et al. (2015), Tafelmeyer et al. (2010),
Starkenburg et al. (2014), Jablonka et al. (2015), and Skúladóttir et al.
(2015). PIGS data have been provided via private communication
with Anke Arentsen. The theoretical data that support the findings of
this study are available upon request to the corresponding author GP.

REFERENCES

Abate C., Pols O. R., Izzard R. G., Karakas A. I., 2015, A&A, 581, A22

Abohalima A., Frebel A., 2018, ApJS, 238, 36

Aguado D. S., Allende Prieto C., González Hernández J. I., Rebolo R., Caffau
E., 2017, A&A, 604, A9

Aguado D. S., Prieto C. A., Hernández J. I. G., Rebolo R., 2018, ApJ, 854, L34

Aguado D. S. et al., 2022, A&A, 668, A86

Aguado D. S. et al., 2023a, MNRAS, 520, 866

Aguado D. S. et al., 2023b, A&A, 669, L4

Amsari A. M., Nissen P. E., Skúladóttir Á., 2019, A&A, 630, A104

Aoki W., Beers T. C., Christlieb N., Norris J. E., Ryan S. G., Tsangarides S., 2007, ApJ, 655, 492

Arentsen A., Starkenburg E., Shetrone M. D., Venn K. A., Depagne É.,
McConnachie A. W., 2019, A&A, 621, A108

Arentsen A. et al., 2020, MNRAS, 496, 4964

Arentsen A. et al., 2021, MNRAS, 505, 1239

Arentsen A., Placco V. M., Lee Y. S., Aguado D. S., Martin N. F., Starkenburg E., Yoon J., 2022, MNRAS, 515, 4082

Barbuy B. et al., 2015, A&A, 580, A40

Beers T. C., Christlieb N., 2005, ARA&A, 43, 531

Beers T. C., Preston G. W., Shectman S. A., 1992, AJ, 103, 1987

Bensby T. et al., 2019, The Messenger, 175, 35

Bonifacio P. et al., 2015, A&A, 579, A28

Bonifacio P. et al., 2021, A&A, 651, A79

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

Caffau E. et al., 2011, Nature, 477, 67

Caffau E. et al., 2013, A&A, 560, A71

Caffau E. et al., 2016, A&A, 595, L6

Carollo D., Freeman K., Beers T. C., Placco V. M., Tumlinson J., Martell S. L., 2014, ApJ, 788, 180

Cenarro A. J. et al., 2019, A&A, 622, A176

Cescutti G., Matteucci F., 2011, A&A, 525, A126

Chiappini C. et al., 2019, The Messenger, 175, 30

Christlieb N., 2003, Rev. Mod. Astron., 16, 191

Christlieb N. et al., 2019, The Messenger, 175, 26

Da Costa G. S. et al., 2019, MNRAS, 489, 5900

Dalton G., 2016, ASP Conf. Ser. Vol. 507, Multi-Object Spectroscopy in the Next Decade: Big Questions, Large Surveys, and Wide Fields. Astron. Soc. Pac., San Francisco, p. 97

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

Diemand J., Moore B., Stadel J., 2005, Nature, 433, 389

Duffau S. et al., 2017, A&A, 604, A128

François P. et al., 2018, A&A, 619, A10

Frebel A., Simon J. D., Geha M., Willman B., 2009, ApJ, 708, 560

Frebel A., Kirby E. N., Simon J. D., 2010, Nature, 464, 72

Frebel A., Simon J. D., Kirby E. N., 2014, ApJ, 786, 74

Gallart C. et al., 2021, ApJ, 909, 192

Gehrels N., 1986, ApJ, 303, 336

Gelli V., Salvadori S., Pallottini A., Ferrara A., 2020, MNRAS, 498, 4134

Gilmore G., Norris J. E., Monaco L., Yong D., Wyse R. F. G., Geisler D., 2013, ApJ, 763, 61

González Hernández J. I., Aguado D. S., Allende Prieto C., Burgasser A. J., Rebolo R., 2020, ApJ, 889, L13

Graziani L., Salvadori S., Schneider R., Kawata D., de Bennassuti M., Maselli A., 2015, MNRAS, 449, 3137

Greif T. H., Springel V., White S. D. M., Glover S. C. O., Clark P. C., Smith R. J., Klessen R. S., Bromm V., 2011, ApJ, 737, 75

Hansen T., Andersen J., Nordström B., Beers T., Placco V., Yoon J., Buchhave L., 2016a, A&A, 586, A160

Hansen T., Andersen J., Nordström B., Beers T., Placco V., Yoon J., Buchhave L., 2016b, A&A, 588, A3

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

Heger A., Woosley S. E., 2010, ApJ, 724, 341

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

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

Hosokawa T., Omukai K., Yoshida N., Yorke H. W., 2011, Science, 334, 1250

Howes L. M. et al., 2014, MNRAS, 445, 4241

Howes L. M. et al., 2015, Nature, 527, 484

Howes L. M. et al., 2016, MNRAS, 460, 884

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

Iwamoto N., Umeda H., Tominaga N., Nomoto K., Maeda K., 2005, Science, 309, 451
Jablonka P. et al., 2015, A&A, 583, A67

Jönsson H. et al., 2020, AJ, 160, 120

Keller S. C. et al., 2007, Publ. Astron. Soc. Aust., 24, 1

Keller S. C. et al., 2014, Nature, 506, 463

Kirby E. N. et al., 2015, ApJ, 801, 125

Lai D. K., Lee Y. S., Bolte M., Lucatello S., Beers T. C., Johnson J. A., Sivarani T., Rockosi C. M., 2011, ApJ, 738, 51

Larson R. B., 1998, MNRAS, 301, 569

Lee Y. S., Beers T. C., Kim Y. K., Placco V., Yoon J., Carollo D., Masseron T., Jung J., 2017, ApJ, 836, 91

Lee Y. S., Beers T. C., Kim Y. K., 2019, ApJ, 885, 102

Li H. et al., 2022, ApJ, 931, 147

Lucatello S., Tsangarides S., Beers T. C., Carretta E., Gratton R. G., Ryan S. G., 2005, ApJ, 625, 825

Lucertini F., Monaco L., Caffau E., Bonifacio P., Mucciarelli A., 2022, A&A, 657, A29

Magg M., Klessen R. S., Glover S. C. O., Li H., 2019, MNRAS, 487, 486

Marassi S., Chiaki G., Schneider R., Limongi M., Omukai K., Nozawa T., Chieffi A., Yoshida N., 2014, ApJ, 794, 100

McMillan P. J., 2011, MNRAS, 414, 2446

Meynet G., Ekström S., Maeder A., 2006, A&A, 447, 623

Ness M. et al., 2013, MNRAS, 430, 836

Nomoto K., Maeda K., Mazzali P. A., Umeda H., Deng J., Iwamoto K., 2004, in Stellar Collapse. Springer, p. 277

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

Nordlander T. et al., 2019, MNRAS, 488, L109

Norris J. E., Gilmore G., Wyse R. F. G., Yong D., Frebel A., 2010, ApJ, 722, L104

Norris J. E. et al., 2012, ApJ, 762, 28

Mendes de Oliveira C. et al., 2019, MNRAS, 489, 241

Pacucci F., Loeb A., Salvadori S., 2017, MNRAS, 471, L72

Placco V. M., Frebel A., Beers T. C., Stancliffe R. J., 2014, ApJ, 797, 21

Placco V. M. et al., 2021, ApJ, 912, L32

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

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

Salpeter E. E., 1955, ApJ, 121, 161

Salvadori S., Schneider R., Ferrara A., 2007, MNRAS, 381, 647

Salvadori S., Ferrara A., Schneider R., Scannapieco E., Kawata D., 2010, MNRAS, 401, L5

Salvadori S., Tolstoy E., Ferrara A., Zaroubi S., 2014, MNRAS, 437, L26

Salvadori S., Skúladóttir Á., Tolstoy E., 2015, MNRAS, 454, 1320

Salvadori S., Bonifacio P., Caffau E., Korotin S., Andreevsky S., Spite M., Skúladóttir Á., 2019, MNRAS, 487, 4261

Scannapieco E., Kawata D., Brook C. B., Schneider R., Ferrara A., Gibson B. K., 2006, ApJ, 653, 285

Schneider R., Ferrara A., Natarajan P., Omukai K., 2002, ApJ, 571, 30

Simon J. D., Frebel A., McWilliam A., Kirby E. N., Thompson I. B., 2010, ApJ, 716, 446

Simon J. D., Jacobson H. R., Frebel A., Thompson I. B., Adams J. J., Shectman S. A., 2015, ApJ, 802, 93

Skúladóttir Á., Tolstoy E., Salvadori S., Hill V., Pettini M., Shetrone M. D., Starkenburg E., 2015, A&A, 574, A129

Skúladóttir Á., Tolstoy E., Salvadori S., Hill V., Pettini M., 2017, A&A, 606, A71

Skúladóttir Á. et al., 2021, ApJ, 915, L30

Spite M., Spite F., François P., Bonifacio P., Caffau E., Salvadori S., 2018, A&A, 617, A56

Spite M., Caffau E., Bonifacio P., Spite F., Ludwig H.-G., Plez B., Christlieb N., 2013, A&A, 552, A107

Spite M., Spite F., François P., Bonifacio P., Caffau E., Salvadori S., 2018, A&A, 617, A56

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

Starkenburg E., Shetrone M. D., McConnellie A. W., Venn K. A., 2014, MNRAS, 441, 1217

Starkenburg E., Oman K. A., Navarro J. F., Crain R. A., Fattahi A., Frenk C. S., Sawala T., Schaye J., 2017a, MNRAS, 465, 2212

Starkenburg E. et al., 2017b, MNRAS, 471, 2587

Suda T., Aikawa M., Machida M. N., Fujimoto M. Y., Iben I., Jr, 2004, ApJ, 611, 476

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

Tafelmeyer M. et al., 2010, A&A, 524, A58

Tan J. C., McKee C. F., 2004, ApJ, 603, 383

Tumlinson J., 2009, ApJ, 708, 1398

White S. D. M., Springel V., 2000, in Weiss A., Abel T. G., Hill V., eds, The First Stars. p. 327

Woosley S., Weaver T. A., 1995, Technical Report, The Evolution And Explosion Of Massive Stars II: Explosive Hydrodynamics And Nucleosynthesis. Lawrence Livermore National Lab., CA (United States)

Yong D. et al., 2012, ApJ, 762, 27

Yoon J. et al., 2016, ApJ, 833, 20

Yoon J. et al., 2018, ApJ, 861, 146

Yoon J., Whitten D. D., Beers T. C., Lee Y. S., Masseron T., Placco V. M., 2020, ApJ, 894, 7

Zoccali M., Hill V., Lecureur A., Barbuy B., Renzini A., Minniti D., Gómez A., Ortolani S., 2008, A&A, 486, 177

This paper has been typeset from a *TEX/LATex* file prepared by the author.