The evolution of CNO elements in galaxies

Abstract

After hydrogen and helium, oxygen, carbon, and nitrogen—hereinafter, the CNO elements—are the most abundant species in the universe. They are observed in all kinds of astrophysical environments, from the smallest to the largest scales, and are at the basis of all known forms of life, hence, the constituents of any biomarker. As such, their study proves crucial in several areas of contemporary astrophysics, extending to astrobiology. In this review, I will summarize current knowledge about CNO element evolution in galaxies, starting from our home, the Milky Way. After a brief recap of CNO synthesis in stars, I will present the comparison between chemical evolution model predictions and observations of CNO isotopic abundances and abundance ratios in stars and in the gaseous matter. Such a comparison permits to constrain the modes and time scales of the assembly of galaxies and their stellar populations, as well as stellar evolution and nucleosynthesis theories. I will stress that chemical evolution models must be carefully calibrated against the wealth of abundance data available for the Milky Way before they can be applied to the interpretation of observational datasets for other systems. In this vein, I will also discuss the usefulness of some key CNO isotopic ratios as probes of the prevailing, galaxy-wide stellar initial mass function in galaxies where more direct estimates from the starlight are unfeasible.

1 Introduction and historical perspective

In the late fifties, Burbidge et al. (1957) wrote their monumental paper titled “Synthesis of the Elements in Stars”, hereby laying the foundations of modern stellar nucleosynthesis theory. In that paper, the authors discuss the sequence of processes that lead to the synthesis of all isotopes known in nature and highlight the importance of parallel advances in nuclear physics and astronomy. These, in turn, rest on the constant interplay among observations, experiments and theory. Progresses in nucleosynthesis theory were spurred by Burbidge et al.’s thoughts and recommendations and are recorded in a series of review papers and textbooks (e.g., Truran 1984; Trimble 1991; Pagel 1997; Wallerstein et al. 1997; Clayton 2003; Nomoto et al. 2013). It is beyond the scope of the present work to give a full account of all possible nucleosynthesis paths activating in stars. The aim being to review studies of carbon, nitrogen, and oxygen (CNO) evolution in galaxies, it will suffice to recall that: (i) virtually all $^{14}\hbox {N}$ observed in the universe originates in hydrostatic hydrogen burning in stars; (ii) the synthesis of $^{13}\hbox {C}$, $^{15}\hbox {N}$, and $^{17}\hbox {O}$ is a consequence of hydrogen burning in stars via both the cold and the hot CNO cycles (see Wiescher et al. 2010, for a review); (iii) $^{12}\hbox {C}$, $^{16}\hbox {O}$, and $^{18}\hbox {O}$ have their origin in helium-burning phases (Salpeter 1952; Clayton 2003). Different isotopes are produced in different amounts and on different time scales by stars of various initial masses and chemical compositions (Tinsley 1979). Stellar duplicity also plays a non-negligible role (see De Donder and Vanbeveren 2004; see also Farmer et al. 2021, for a discussion of the effects of binary interactions on carbon production in massive stars). Therefore, in order to follow properly the evolution of the abundances and abundance ratios of different CNO elements in the interstellar medium (ISM), one has to rely on detailed chemical evolution models, relaxing the instantaneous recycling approximation (IRA) and considering the occurrence of binary systems (Matteucci 2021, and references therein).

Oxygen, carbon, and nitrogen follow hydrogen and helium as the most abundant elements in the universe and may aggregate to form the basic bricks of life. Their spectral features are observed in planetary and stellar atmospheres as well as in gaseous matter in a variety of environments (see, e.g., reviews by Wilson and Rood 1994; Tinetti et al. 2013; Maiolino and Mannucci 2019; Randich and Magrini 2021). Presolar dust grains isolated in meteorites feature many CNO compounds (Zinner 2014; Nittler and Ciesla 2016).

Because of their ubiquity and heterogeneous stellar production sites, the CNO elements have been the subject of extensive study. In their pioneering work, Audouze et al. (1975) and Vigroux et al. (1976) used simple Galactic evolutionary models to link the lower than solar $^{12}\hbox {C}/^{13}\hbox {C}$ ratio inferred from molecular data of the local ISM to delayed production of $^{13}\hbox {C}$ from low-mass stars, including novae. They further noted that the hot CNO cycle occurring in the fast explosive phases of novae may be responsible for the production of $^{15}\hbox {N}$ (see also Dearborn et al. 1978) and that the minor CNO isotopes may be not only produced, but also destroyed in stellar interiors (see also Salpeter 1955a), meaning that negative yields¹ must be implemented in chemical evolution codes.

In the eighties, the first numerical models relaxing IRA as well as higher quality molecular data became available. Taking advantage of this, Tosi (1982) analysed the behavior of several CNO isotopic ratios with the galactocentric distance. She concluded that the $^{12}\hbox {C}/^{13}\hbox {C}$ and $^{16}\hbox {O}/^{17}\hbox {O}$ ratios behave as typical primary-to-secondary² element ratios, $^{14}\hbox {N}$ has both a primary and a secondary component and $^{15}\hbox {N}$ is a regular secondary element. She also pointed out that some unknown mechanism is inhibiting the ejection of $^{18}\hbox {O}$ in the ISM over the last few Gyr of Galactic evolution.

The problem of the dual nature of the main N isotope has been a longstanding one. Truran and Cameron (1971) first noticed that some primary $^{14}\hbox {N}$ production from intermediate-mass stars with initial mass in the range 2–4 $M_{\odot }$, resulting from convective mixing of material from the helium-burning region to the hydrogen-burning layer is needed, being the sole secondary production insufficient to account for the solar system abundance of this element. The [N/Fe]³ $\simeq$ 0.0 ratios measured in metal-poor Galactic halo dwarfs by Sneden (1974) then suggested that this mechanism is at work also in massive stars (see also Barbuy 1983). Observations of nitrogen abundances in H II regions of our own as well as nearby spiral and irregular galaxies, comprising low-metallicity dwarf galaxies (e.g., Smith 1975; Peimbert et al. 1978; Lequeux et al. 1979; Garnett 1990), reinforced the idea that $^{14}$N must have a substantial primary component, though the exact stellar mass range appointed to its primary production remained debated (e.g., Edmunds and Pagel 1978; Matteucci and Tosi 1985; Matteucci 1986): hot bottom burning was readily identified as the process leading to primary $^{14}\hbox {N}$ (and $^{13}\hbox {C}$) production in massive asymptotic giant branch (AGB) stars (Renzini and Voli 1981), but we had to wait until the early 2000s to have stellar rotation pinned down as a mechanism able to trigger substantial primary $^{14}\hbox {N}$ production in low-metallicity massive stars (Maeder and Meynet 2000; Meynet and Maeder 2002b).

Coming back to early Galactic chemical evolution (GCE) studies, it is worth mentioning the work of D’Antona and Matteucci (1991) and Matteucci and D’Antona (1991). These authors implemented a self-consistent calculation of the nova rate in their GCE model and concluded that novae can represent a non negligible source of $^{13}\hbox {C}$ and $^{15}\hbox {N}$ in the Galaxy.⁴ In particular, novae would not overproduce $^{13}\hbox {C}$ and $^{15}\hbox {N}$ as long as these elements are synthesized proportionally to $^7\hbox {Li}$ and no major contribution to their primary component comes from intermediate-mass stars. A couple of years later, Henkel and Mauersberger (1993) modeled the evolution of CNO isotopes in the nuclear regions of active galaxies and first suggested that a stellar initial mass function (IMF) skewed in favour of high-mass stars could be responsible for the large $^{18}\hbox {O}/^{17}\hbox {O}$ and small $^{16}\hbox {O}/^{18}\hbox {O}$ ratios observed in starbursts. Indeed, following earlier detections of $^{13}\hbox {CO}$ in five galaxies by Encrenaz et al. (1979), observations of carbon monoxide isotopes at millimeter wavelengths were adding up for extragalactic sources, pointing to physical processes in mergers and starbursts sensibly different than those acting in secularly evolving systems (Aalto et al. 1991; Casoli et al. 1992). The idea of an IMF varying with the environment has been given further support in recent years, as we will see in Sect. 3.2.

The problem of the evolution of C and O isotopes in the Galaxy was then reassessed by Prantzos et al. (1996), who confirmed the puzzling behaviour of $^{18}\hbox {O}$ pointed out in previous studies and the $^{17}\hbox {O}$ overproduction on Galactic scales, possibly due to the low $^{17}\hbox {O}$(p, $\gamma$)$^{18}\hbox {F}$ and $^{17}\hbox {O}$(p, $\alpha$)$^{14}\hbox {N}$ reaction rates adopted in stellar nucleosynthesis calculations. Later on, Romano and Matteucci (2003), Kobayashi et al. (2011), and Nittler and Gaidos (2012) still found the need for a reduction of the stellar yields of $^{17}\hbox {O}$ to explain its solar system abundance. An improved direct measurement of the 64.5 keV resonance strength in the $^{17}\hbox {O}$(p, $\alpha$)$^{14}\hbox {N}$ reaction performed at the Laboratory for Underground Nuclear Astrophysics (LUNA), a unique facility located deep underground in the Gran Sasso National Laboratory, Italy, has recently increased by a factor of two the rate of this reaction at temperatures relevant to shell hydrogen burning in red giant branch (RGB) and AGB stars (Bruno et al. 2016). By adopting the augmented rate, Straniero et al. (2017) have obtained $^{17}\hbox {O}$ yields from 15 to 40 percent smaller for stars with initial mass in the range 2–20 $M_{\odot }$, thus alleviating the discrepancies between GCE model predictions and $^{17}\hbox {O}$ measurements.

When interpreting the CNO abundance data by means of classic GCE models, one relies in essence on a generalization of the time-delay model (Tinsley 1979, see also a discussion in Vincenzo et al. 2016). The time-delay model in its narrower sense explains the [O/Fe] versus [Fe/H] pattern traced by stars in the solar vicinity as due to prompt injection of oxygen from short-lived core-collapse supernovae (CCSNe) and delayed iron pollution from type Ia supernovae (SNeIa), the fundamental point being that elements produced via primary nucleosynthetic processes behave as secondary elements from the point of view of chemical evolution when manufactured by low-mass stars with long lifetimes (Tinsley 1979, see also Matteucci and François 1989). On the other hand, primary elements forged in the outer layers of massive stars may exhibit a pseudo-secondary character because of metallicity dependence driven by stellar mass loss (Maeder 1992). It is to this that no unanimous verdict has been reached yet on whether most of $^{12}\hbox {C}$ we observe today in the solar neighbourhood comes from low-mass (Tinsley 1979; Chiappini et al. 2003; Gavilán et al. 2005; Mattsson 2010) or massive stars (Prantzos et al. 1994; Henry et al. 2000; Akerman et al. 2004; Romano et al. 2020): in fact, depending on the adopted stellar yields, different authors reach opposite conclusions. Possible variations in the stellar IMF complicate further the picture (Mattsson 2010; Romano et al. 2020). Alternative interpretations foresee inhomogeneous chemical enrichment in the framework of hydrodynamical cosmological simulations coupled to alterations of CNO synthesis due to failed SNe at high metallicities (Vincenzo and Kobayashi 2018a, b).

Clearly, our capability of constructing well-sound chemical evolution models that explain self-consistently CNO abundance measurements in various galaxies rests mostly on the robustness of the assumed stellar yields (see Romano et al. 2010; Côté et al. 2017). These must be computed for as dense a grid of initial stellar masses and metallicities as possible. However, until recently, grids of stellar yields dense enough were woefully lacking in the literature. In particular, yields for very metal-poor and/or super-solar metallicity stars were often missing, essentially due to the paucity, or even absence, of corresponding data to calibrate the uncertain parameters of stellar evolution at the lowest/highest metallicities. I will come back to this issue in Sect. 2. Because of this dearth of stellar yield data, uncertain interpolations of the yields in the very metal-poor domain have often been used in chemical evolution models for dwarf spheroidal (dSphs), dwarf irregular (dIrrs), blue compact dwarf (BCDs), and ultrafaint dwarf galaxies (UFDs). Similarly, when dealing with systems forming high fractions of metal-rich stars, such as elliptical galaxies and bulges of spirals, chemical evolution modelers have rather arbitrarily extrapolated solar-metallicity yields to the super-solar metallicity regime.

As regards nearby low-mass galaxies, dSphs and UFDs are currently devoid of gas; their CNO abundances are derived from giant stars, with the obvious drawback that the measured abundances of C and N need corrections for stellar evolutionary effects to obtain information on the ISM composition at the time the stars were born (e.g., Ishigaki et al. 2014; Kirby et al. 2015; Lardo et al. 2016, and references therein). On the other hand, dIrrs and BCDs contain large fractions of gas; thus, their present-day CNO abundances can be obtained quite easily from ultraviolet (UV) and optical nebular emission lines (see Annibali and Tosi 2022, for a review). Models aimed at explaining abundance measurements in star-forming dwarf galaxies considered self-enrichment of H II regions following gasping or bursting star formation (Tosi et al. 1991; Marconi et al. 1995; Tolstoy et al. 2009, and references therein), often taking into account the development of galactic winds ( Matteucci and Chiosi 1983; Pilyugin 1993; Carigi et al. 1995; Pilyugin 1999, among others). In particular, differential outflows powered by SN explosions and preferentially enriched in SN ejecta were invoked to explain the observed low metallicities and relatively high N/O ratios (Marconi et al. 1994; Recchi et al. 2001; Romano et al. 2006). Galactic winds flowing from intensely star-forming regions can significantly impact the structure and subsequent evolution of the host galaxies (Veilleux et al. 2005). Current observational evidence seems to point to a weak feedback in low-mass galaxies (McQuinn et al. 2019; Mancera Piña et al. 2020; Xu et al. 2022) with preferential loss of metals (McQuinn et al. 2015), in agreement with the results of hydrodynamical simulations that resolve the cooling radius of the expanding interstellar bubbles (Romano et al. 2019a, and references therein).

As for metal-rich systems, a distinction must be drawn between those that can be resolved in stars, such as the inner regions of the Milky Way, and those with ensemble properties from light-averaged spectra, like massive early-type galaxies at both low and high redshifts. Recently, Bensby et al. (2021) have measured C and O abundances of 70 microlensed dwarf, turnoff and subgiant stars residing in the Galactic bulge, thus providing the first statistically significant sample of bulge stars with C abundances not altered by stellar evolutionary processes. The almost total lack of C abundances for unevolved bulge stars until 2021 (see Romano et al. 2020, and references therein) was related to the inherent difficulties in obtaining high-resolution spectra of stars that have apparent magnitudes as faint as $V =$ 18–21 at the distance of the bulge. However, during a microlensing event the brightness of the target increases, making it possible to get spectra suitable for detailed abundance analyses with reasonable exposure times on large telescopes (Bensby et al. 2021). Nitrogen abundances for a large number of bulge stars are provided by the APOGEE-2 survey (see Kisku et al. 2021, and references therein) but only for giants, which requires uncertain evolutionary corrections (Grisoni et al. 2021, and references therein). Light-averaged element abundance ratios for 3802 early-type galaxies from the Sloan Digital Sky Survey (SDSS, York et al. 2000) including [C/Fe], [N/Fe], and [O/Fe], are presented by Johansson et al. (2012). In the most massive systems, the [C/O] ratio is slightly higher than the solar value, namely, [C/O]$\sim$ 0.05, while [N/O] is sub-solar. More recently, near-infrared (NIR) and millimeter telescopes have allowed to unravel the properties of distant galaxies, encompassing the determination of the N/O versus O/H and N/O versus stellar mass relations as well as measurements of CNO isotope ratios from observations of the gas phase (e.g., Zhang et al. 2018; Brown and Wilson 2019; Hayden-Pawson et al. 2022, see also Maiolino and Mannucci 2019, and references therein).

Overall, providing a coherent picture of the evolution of all seven stable CNO isotopes in galaxies of different morphological type remains one of the open quests for modern astrophysics. This review aims at giving an overview of recent advances in complementary research fields that will hopefully result in an improved understanding of the evolution of CNO elements in galaxies in the next few years. Section 2 discusses progresses in stellar evolution and nucleosynthesis theory, without neglecting to mention the many uncertainties that still plague stellar yield calculations. Section 3 deals with advances in chemical evolution modeling driven by improvements in the underlying theory as well as instrumental developments that provide more precise and accurate data. In closing the paper, a discussion (Sect. 4) is followed by a list of future challenges and novel perspectives that are opening up (Sect. 5).

2 Synthesis of CNO elements in stars

We have now a firm grasp of the thermonuclear reactions that trigger the synthesis of the CNO elements in stars. We know that $^{12}\hbox {C}$ is synthesised during He burning through the so-called triple-$\alpha$ capture process, whose rate is dominated by capture into the 7.65 MeV Hoyle state in $^{12}\hbox {C}$ (Hoyle 1954; Freer and Fynbo 2014), and that $^{16}\hbox {O}$ is simultaneously produced at the expenses of $^{12}\hbox {C}$ via fusion with another $\alpha$ particle (He nucleus). The cold CNO cycle (von Weizsäcker 1938; Bethe 1939) takes place in the H-burning zones of main sequence⁵ and giant branch stars, leading to the production of $^{14}$N, $^{13}\hbox {C}$, and $^{17}\hbox {O}$ through its different branches; stellar convection can then bring the products of CNO cycle burning to the surface (Fowler et al. 1955; Burbidge et al. 1957; Dearborn and Schramm 1974). The activation of the hot CNO cycle at temperatures exceeding $10^8$ K explains the formation of $^{15}$N and $^{17}\hbox {O}$ in H-rich zones (Audouze et al. 1973; Wiescher et al. 2010, and references therein). Finally, most of $^{18}\hbox {O}$ produced early in He burning through the $^{14}\hbox {N}$($\alpha$, $\gamma$)$^{18}\hbox {F}$($\beta ^+$)$^{18}\hbox {O}$ reaction chain is destroyed by $^{18}\hbox {O}$($\alpha$, $\gamma$)$^{22}\hbox {Ne}$ reactions occurring at higher temperatures in stellar cores. However, some $^{18}\hbox {O}$ may survive in the He shell of massive stars, especially if fresh $^{14}\hbox {N}$ is brought to the shell by rotational mixing (Hirschi 2007; Limongi and Chieffi 2018). Even if these nucleosynthetic paths may seem crystal clear, we are, sadly, far from having a satisfactory account of CNO production in stars.

The relative rates of the triple-$\alpha$ and $^{12}\hbox {C}$($\alpha$, $\gamma$)$^{16}\hbox {O}$ reactions as functions of temperature and density determine the carbon-to-oxygen abundance ratio at the end of He burning, thus affecting the abundances of later nucleosynthesis products (in massive stars) and the composition of the stellar remnants. Unfortunately, the $^{12}\hbox {C}$($\alpha$, $\gamma$)$^{16}\hbox {O}$ cross section in the Gamow window corresponding to temperatures for hydrostatic He-burning in stars is very small, of the order of a few $10^{-17}$ barn, well below the current direct measurement sensitivity (Kunz et al. 2002). Laboratory measurements are thus performed at higher energies and the results extrapolated to the astrophysically relevant energy range, which makes the stellar rate of the $^{12}\hbox {C}$($\alpha$, $\gamma$)$^{16}\hbox {O}$ reaction one of the most important unsettled rates in nuclear astrophysics (see, e.g., Buchmann and Barnes 2006; Aliotta et al. 2022). As for the triple-$\alpha$ reaction, a new measurement of the radiative width of the Hoyle state (after more than 40 years!) by Kibédi et al. (2020) has resulted in a 34% increase in the reaction rate; its consequences on full stellar nucleosynthesis calculations still need to be explored (see, however, Fields et al. 2018). Regarding H burning through the CNO cycles, it is worth mentioning the results obtained over 30 years by the LUNA Collaboration as well as ongoing experiments for many key reactions such as, for instance, $^{14}\hbox {N}$(p, $\gamma$)$^{15}\hbox {O}$, $^{15}\hbox {N}$(p, $\gamma$)$^{16}\hbox {O}$, $^{17}\hbox {O}$(p, $\alpha$)$^{14}\hbox {N}$, $^{17}\hbox {O}$(p, $\gamma$)$^{18}\hbox {F}$, $^{18}\hbox {O}$(p, $\alpha$)$^{15}\hbox {N}$, $^{18}\hbox {O}$(p, $\gamma$)$^{19}\hbox {F}$, $^{12}\hbox {C}$(p, $\gamma$)$^{13}\hbox {N}$, and $^{13}\hbox {C}$(p, $\gamma$)$^{14}\hbox {N}$ (e.g., Ananna et al. 2022, and references therein).

Updated rates for these reactions impact the production of CNO elements in stars (Herwig et al. 2006; Fields et al. 2018) and, thus, our comprehension of the origin of stardust grains (Palmerini et al. 2011; Lugaro et al. 2017), mixing processes in giant stars (Straniero et al. 2017, and references therein), stellar age dating (Imbriani et al. 2004; Salaris et al. 2015; Casali et al. 2019), and chemical evolution of galaxies (Romano and Matteucci 2003; Kobayashi et al. 2011; Romano et al. 2017, 2019b).

2.1 Single stars

Even if the rates of all the reactions involved in stellar nucleosynthesis calculations were known with the necessary accuracy and precision, we would still be struggling with several major problems. First, the integration of nuclear reaction networks is a difficult task that requires special care in the choice of the integration method to avoid numerical instability outbreaks (see, e.g., discussions in Longland et al. 2014; Yagüe López et al. 2022).

Second, different treatments of convection⁶ in one-dimensional (1D) stellar evolution codes—including overshoot of convective eddies⁷—rely on parametric recipes to approximate physical phenomena that would require a hydrodynamical non-local description, which leads to several shortcomings (Renzini 1987; Kupka and Muthsam 2017; Arnett et al. 2019). Other mixing processes have been identified and incorporated in the computations (see the review by Salaris and Cassisi 2017), among which semiconvection (Schwarzschild and Härm 1958; Langer et al. 1985; Merryfield 1995), thermohaline (or fingering) convection (Charbonnel and Zahn 2007b; Eggleton et al. 2008; Cantiello and Langer 2010; Denissenkov 2010), atomic diffusion (Michaud 1970; Michaud et al. 2015, and references therein) magnetic buoyancy (Busso et al. 2007; Denissenkov et al. 2009), meridional circulation and shear-induced turbulence⁸ (Maeder 2009, and references therein). Besides taking into account the contribution of these and additional processes generally not included in stellar models—such as, for instance, gravity waves or internal magnetic fields—it would be desirable to have their couplings and interactions properly modelled. For example, Charbonnel and Zahn (2007a) have shown that stellar magnetic fields may inhibit thermohaline mixing, while the rotation might contrast the effects of atomic diffusion (see, e.g., Deal et al. 2020) in low-mass stars.⁹ Stellar mass loss (see next paragraph) can suppress microscopic diffusion in the outer stellar layers as well (e.g., Vauclair and Charbonnel 1995) and Maeder et al. (2013) have clearly shown that some well-known instabilities must be treated jointly in rotating stars, as they interact with each other.

Mass loss prescriptions deeply impact the outcome of stellar evolution models, including the predicted yields (e.g., Chiosi et al. 1978; Maeder 1992; Stancliffe and Jeffery 2007; Beasor et al. 2021). Stars lose mass at widely varying rates during their lifespan, with a strong dependence of the wind strength on stellar mass, metallicity and molecular opacities (see, e.g., Marigo 2002; Vink and de Koter 2005; Ventura and Marigo 2010). At present, mass loss rates can not be deduced from first principles in 1D stellar models that, hence, rely on empirical relations or theoretical predictions. Only time-averaged rates are implemented, although it is clear that both episodic and eruptive events can occur in real stars. Low- and intermediate-mass stars (objects with initial masses between about 1 and 9 times the mass of our Sun) lose a large fraction of their mass before leaving the AGB. Mass loss rates ${\dot{M}}_{\mathrm {AGB}} \sim 10^{-8}$–$10^{-5}\,M_{\odot }\, \mathrm {yr}^{-1}$ are quite typical of AGB stars; rates as high as $10^{-4}\,M_{\odot }\, \mathrm {yr}^{-1}$ are observed at the end of the AGB phase¹⁰. Mass loss rates measured in massive stars span about 8 orders of magnitude. Stars with initial mass in the range 10–25 $M_{\odot }$ become red supergiants (RSGs), while stars above 40 $M_{\odot }$ become Wolf-Rayet (WR) stars, with some overlapping in-between. At variance with low- and intermediate-mass stars, for which wind mass loss is relatively unimportant until the most advanced phases of the evolution, massive stars above $\sim$20 $M_{\odot }$ expel significant fractions of their envelopes already on the main sequence; mass loss in the post-main sequence determines the final fate of the star, i.e., the kind of SN explosion (see Smith 2014). I refer the reader to the excellent reviews by Höfner and Olofsson (2018), Decin (2021), and Vink (2022)—focusing, respectively, on low- and intermediate-mass stars, AGB and RSG stars, and stars more massive than 25 $M_{\odot }$—for exhaustive reports on mass-loss rate measurements in different wavelength bands, recent advances in the understanding of the mechanisms underpinning and sustaining the outflows, the role of wind clumping and porosity, and any remaining issues in assessing the dependence of the mass-loss rates on the fundamental stellar parameters.

In the following, I review studies of stellar evolution and nucleosynthesis, focusing mostly on those that have produced grids of yields adequate for use in GCE studies. Moreover, I concentrate on CNO element production.

2.1.1 Low- and intermediate-mass stars

The evolution and nucleosynthesis of single low- and intermediate-mass stars, including super-AGB stars, are nicely reviewed in Karakas and Lattanzio (2014) and Doherty et al. (2017). In this section, I summarize briefly the evolutionary paths of these stars in relation to CNO element production. In addition, I report on the latest theoretical developments.

Stars with zero-age main-sequence masses in the range $\sim$1–9 $M_{\odot }$ have lifetimes varying from about 12 Gyr to $<30$ Myr, and enrich the ISM with significant quantities of $^{12}\hbox {C}$, $^{13}\hbox {C}$, $^{14}\hbox {N}$, and $^{17}\hbox {O}$ (only to mention the elements that are the focus of this review). They are, furthermore, efficient sinks of $^{15}\hbox {N}$ and $^{18}\hbox {O}$. The outer stellar layers are expelled by pulsation-enhanced dust-driven winds when the stars are climbing the giant branches.

During the ascent of the RGB, when the star is burning H in a shell surrounding the contracting He core, material previously exposed to partial H burning in the stellar interior is conveyed to the surface (first dredge-up episode). Low-mass stars ($m \le 3\,M_{\odot }$, with the limiting mass partly dependent on the initial chemical composition) may experience non-standard (also referred to as extra) mixing processing in the upper part of the RGB that further modifies the chemical composition of the envelope (see Stancliffe et al. 2009; Lagarde et al. 2019, for the impact of such processes on C and N atmospheric abundances). At this stage, the stellar envelope is enriched in $^{13}\hbox {C}$, $^{14}\hbox {N}$, and $^{17}\hbox {O}$, while the abundances of $^{12}\hbox {C}$, $^{16}\hbox {O}$, and $^{18}\hbox {O}$ decrease. The RGB phase ends when He ignites in the core.¹¹ Following core He exhaustion, the star develops a complex structure, with a degenerate CO core surrounded by a thermally unstable He-burning shell, separated from an outer H-burning shell by a He-rich region that contains the ashes of H burning. On top of the H-burning shell, a thin radiative layer¹² is overlaid by a deep convective envelope. When He firstly ignites in the shell following the contraction of the core, a substantial amount of energy is released, which leads to an inwards expansion of the envelope. In stars more massive than about 4 $M_{\odot }$, the H-burning layer is extinguished and the convective envelope can penetrate the inner regions: this is the second dredge-up episode that brings to the surface the products of complete H burning, leading to a further increase in $^{14}$N and decrease in $^{12}\hbox {C}$, $^{13}\hbox {C}$, and $^{15}$N. In the most massive intermediate-mass stars C ignites off center in the core in conditions of partial degeneracy, burning in an inward-propagating flame that eventually produces a CO-Ne or ONe core (Garcia-Berro and Iben 1994; Farmer et al. 2015). Objects climbing the AGB with ONe cores are called super-AGB stars.¹³

Now, the thermally-pulsing AGB evolution begins, in which long periods of quiescent H-burning in shell, known as interpulse phases, are interspersed with shell He flashes/pulses powering effective convective motions that homogenize the intershell region; as the pulse fades, convection in the intershell retreats. The energy produced by the pulse makes the envelope to expand (while the star is overall contracting), which switches off the H-burning shell. The base of the outer convective envelope can thus reach regions previously mixed by intershell convection and carry to the surface matter enriched in He-burning products (third dredge up): the star becomes enriched in $^{12}\hbox {C}$. In stars with mass above 3–5 $M_{\odot }$ (the exact threshold value depending on the metallicity and on the details of convection modelling) the profiles of the thermodynamic variables in the regions above the degenerate core are so steep that the convective envelope is partially overlapped with the H-burning shell. The innermost regions of the external mantle of these stars are therefore exposed to H-burning nucleosynthesis, a phenomenon usually referred to as a hot bottom burning (HBB): most $^{12}\hbox {C}$ is converted to $^{14}\hbox {N}$ and $^{13}\hbox {C}$, moreover, large amounts of $^{17}\hbox {O}$ are produced and quickly mixed into the envelope (Renzini and Voli 1981; Wood et al. 1983). The comparison between model predictions and abundances observed in AGB stars seems to suggest that some extra mixing is required in low-mass, low-metallicity AGBs (Boothroyd et al. 1995), but the nature of the process driving the mixing is unclear at present. Regarding solar-metallicity AGB stars, Karakas et al. (2010) have demonstrated that the inclusion of extra mixing on the RGB impacts the subsequent evolution, leading to theoretical C and N abundances consistent with those measured in Galactic carbon-rich stars.

Yields from low- and intermediate-mass stars Low- and intermediate-mass stars pollute the ISM through wind mass ejection. The net yield of element j of a star of initial mass m, that is, the amount of newly-produced element j ejected in the ISM during the star lifespan, is

$$\begin{aligned} m_j^{\mathrm{{new}}} = m p_j = \int _0^{\tau (m)}[X_j(t) - X_j^{\mathrm{{init}}}] \, {\dot{M}}(t) \, \mathrm{d} t, \end{aligned}$$

(1)

where $\tau (m)$ is the stellar lifetime, $X_j(t)$ is the abundance of element j in the ejecta at any time, $X_j^{\mathrm{{init}}}$ is the abundance of element j in the gas out of which the star is born, and ${\dot{M}}(t)$ is the mass loss rate. A negative yield simply means that element j is destroyed in the stellar interior rather than produced. The total amount of mass ejected in the form of element j is obtained by adding to $m_j^{\mathrm{{new}}}$ the mass originally in the form of element j that is ejected without any nuclear processing in the star:

$$\begin{aligned} m_j^{\mathrm{{tot}}} = m_j^{\mathrm{{old}}} + m_j^{\mathrm{{new}}} = (m-m_{\mathrm{{remn}}}) \, X_j^{\mathrm{{init}}} + m p_j, \end{aligned}$$

(2)

where $m_{\mathrm{{remn}}}$ is the mass of the remnant; $m_j^{\mathrm{{tot}}}$ is, obviously, always positive (it can be zero in practice for deuterium and $^6\hbox {Li}$). Allowing for episodic/eruptive mass loss events rather than for a smooth mass loss rate during a star’s evolution could lead to large variations in the predicted yields, especially if the eruptive events coincide with peaks in the surface abundances. Therefore, while many uncertainties (in nuclear reaction rates, molecular opacities, treatment of mixing processes, boundaries of mixed regions, etcetera) affect the calculation of the yields, perhaps the largest uncertainties arise from the assumed mass loss recipes.

Renzini and Voli (1981) computed parameterized synthetic evolutionary models for stars with initial mass in the range 1–8 $M_{\odot }$ and initial metallities $Z =$ 0.004 and 0.02, from the main sequence to the planetary nebula ejection (or carbon ignition in the core), and produced the first grid of metallicity-dependent yields for low- and intermediate-mass stars. The effects of dredge ups, HBB and mass loss were included in the computations and, interestingly, the primary and secondary $^{13}\hbox {C}$ and $^{14}$N components provided separately. Since then, improved treatments of the AGB phase and exploration of wider ranges of parameters have made available finer grids of stellar yields suitable for use in GCE studies (e.g., van den Hoek and Groenewegen 1997; Marigo 2001; Izzard et al. 2004; Karakas 2010). The first yields for super-AGB stars spanning a useful range of initial metallicities, $Z =$ 0.0001$-$0.04, have appeared more recently (Siess 2010).

However, chemical evolution modellers have the annoying habit of asking always for more, with more refined and extended grids of yields usually being on top of the list (closely followed by ‘more data’). Indeed, to avoid pernicious interpolations and extrapolations (see Romano et al. 2010, their Sect. 2.3), yields for use in GCE codes should cover, as homogeneously and densely as possible, the mass and metallicity ranges of stars that contribute to the chemical enrichment of a variety of stellar systems, both ‘here and now’ and ‘there and then’. In recent years, there has been a considerable improvement in this direction. In the last decade, Ventura and collaborators (Ventura et al. 2013, 2014, 2018, 2020, 2021) have put forward a homogeneous yield set for low- and intermediate-mass stars, including super-AGB stars, covering a large grid of initial masses (from about 1 to 8 $M_{\odot }$) and metallicities ($Z = 3 \times 10^{-7}$, $3 \times 10^{-5}$, $3 \times 10^{-4}$, 0.001, 0.004, 0.008, 0.014, 0.03, and 0.04). All the models are computed with the stellar evolutionary code aton (Ventura et al. 1998, and references therein); the yields are currently provided for light elements up to Fe and are being complemented with the ones for s-process elements using a newly developed post-processing code coupled to the aton models (Yagüe López et al. 2022).

Fig. 1

Stellar yields of $^{12}\hbox {C}$ (dashed black lines), $^{13}\hbox {C}$ (dash-dotted green lines), $^{14}\hbox {N}$ (dotted red lines), $^{16}\hbox {O}$ (solid blue lines), and $^{17}\hbox {O}$ (dash-dot-dot yellow lines) for a homogeneous set of low- and intermediate-mass stellar models (Ventura et al. 2013, 2014, 2018, 2020, 2021) computed by means of the aton numerical code (Ventura et al. 1998, and references therein). The yields of $^{13}\hbox {C}$ and $^{17}\hbox {O}$ are multiplied by magnification factors ($\mathrm {M} = 100$, $\mathrm {N} = 1000$ for $Z < Z_\odot$ and $\mathrm {M} = 10$, $\mathrm {N} = 100$ for $Z \ge Z_\odot$) to make them clearly visible. The yields do not include the amount of mass ejected without nuclear processing; hence, negative yields indicate that a given element is destroyed inside the star rather than produced

Figure 1 shows the full grid of yields for the stable CNO isotopes from the work of Ventura and collaborators. The yields of $^{15}\hbox {N}$ and $^{18}\hbox {O}$—mostly negative, especially at high metallicities, owing to the more efficient CNO cycling—are not shown, in order to avoid overcrowding. It is immediately seen that low-mass stars are the main $^{12}\hbox {C}$ producers, which is readily explained by the fact that they experience only the third dredge up, while their higher-mass counterparts are exposed to HBB. In the low-mass domain,¹⁴ the $^{12}\hbox {C}$ yield increases with increasing the stellar mass, owing to the larger number of third dredge-up episodes experienced at higher masses. When HBB turns on, the $^{12}\hbox {C}$ yields turn negative. The $^{16}\hbox {O}$ yield follows a similar qualitative behaviour. As regards $^{14}\hbox {N}$, the (scarce) production in the low-mass domain is due to the occurrence of the first dredge up, while significant production takes place when HBB is activated. It is worth emphasizing that the effects of HBB are more pronounced at low metallicities and that the minimum mass for HBB increases with Z. Overall, the largest $^{14}\hbox {N}$ yields are found at high masses and high metallicities ($Z \ge 0.008$), where secondary production dominates (see a similar figure in Vincenzo et al. 2016, where the primary and secondary contributions to $^{14}\hbox {N}$ synthesis are shown separately). Likewise $^{14}\hbox {N}$, $^{13}\hbox {C}$ is brought to the surface by the first dredge up and its surface abundance thereby increases in stars of all masses; however, unlike $^{14}\hbox {N}$, its surface abundance may decrease as a consequence of the second dredge up. A further enhancement can result from HBB if creation via proton capture on $^{12}\hbox {C}$ prevails over destruction through the reaction $^{13}\hbox {C}$(p, $\gamma$)$^{14}\hbox {N}$. The spikes in the predicted $^{13}\hbox {C}$ yield may depend on some fine tuning of the HBB process, as also noticed by Marigo (2001). The only oxygen isotope produced in noticeable amounts by the models is $^{17}\hbox {O}$. For it to be produced, it is necessary that the reaction chain $^{16}\hbox {O}$(p, $\gamma$)$^{17}\hbox {F}$($\beta ^+$ $\nu$)$^{17}\hbox {O}$ overwhelms the sink reactions $^{17}\hbox {O}$(p, $\gamma$)$^{18}\hbox {F}$ and $^{17}\hbox {O}$(p, $\alpha$)$^{14}\hbox {N}$. Interestingly, the astrophysical rate of the latter reaction has been revised upwards lately (Bruno et al. 2016). We refer the reader to the original papers by Ventura et al. (2013, 2014, 2018, 2020, 2021) for a more exhaustive discussion of the dependencies of the yields on the model parameters.

Most recently, Cinquegrana and Karakas (2022) have applied the post-processing nucleosynthesis code monsoon developed at Monash University (Cannon 1993; Lugaro et al. 2012) to stellar evolutionary models computed by Karakas et al. (2022) and published a grid of yields for stars with initial mass in the range 1–8 $M_{\odot }$ and initial metallicities $0.04 \le Z \le 0.10$. It is worth noticing that there are no other yields in the literature for metallicities exceeding twice the solar value. This study shows that the efficiency of the third dredge up is strongly reduced as the metallicity increases. Furthermore, metal-rich stars have higher opacities, are less dense and, thus, reach lower temperatures at the base of their convective envelopes, which allows for only partial hydrogen burning, if any. Nevertheless, irrespective of the temperature at the base of the envelope, the metal-rich intermediate-mass stars produce large amounts of $^{13}\hbox {C}$, $^{14}\hbox {N}$, and $^{17}\hbox {O}$, while efficiently destroying $^{15}\hbox {N}$, $^{16}\hbox {O}$, and $^{18}\hbox {O}$, which suggests that the nucleosynthetic output is determined by a combination of the effects of the first and second dredge ups.

Given the importance of the adopted mass loss prescriptions in determining how and when the AGB phase ends and, hence, in setting the stellar yields, it should be noted that the mass loss rates implemented in the super-solar metallicity models of Cinquegrana and Karakas (2022) and Ventura et al. (2020)—which are taken, respectively, from Vassiliadis and Wood (1993) and Bloecker (1995)—have been calibrated against metallicity in a limited range. In fact, Vassiliadis and Wood (1993) have evolved models with $0.89 \le m/\hbox {M}_\odot \le 5$ and $0.001 \le Z \le 0.016$, while Bloecker (1995) considered stars with initial masses between 1 and 7 $M_{\odot }$ and initial composition $(X, Y, Z) = (0.739, 0.24, 0.021)$. It is, thus, unclear if the formulae employed still hold in the case of very metal-rich stars.

The yield set by Cinquegrana and Karakas (2022) complements previous grids computed with the same post-processing tool for stars of different initial chemical composition presented in Lugaro et al. (2012, $0.9 \le m/ \hbox {M}_\odot \le 6$ and $Z = 1 \times 10^{-4}$), Fishlock et al. (2014, $1 \le m/M_{\odot } \le 7$ and $Z = 1 \times 10^{-3}$), Karakas et al. (2018, $1 \le m/$ $\hbox {M}_\odot \le 7$ and $Z = 0.0028$), and Karakas and Lugaro (2016, $1 \le m/M_{\odot } \le 8$ and $Z =$ 0.007, 0.014, and 0.03). Coupling with the results of computations by Doherty et al. (2014a, 2014b), who investigated five metallicities from $Z = 1 \times 10^{-4}$ to 0.02, allows to extend the mass grid up to $\sim 9\,M_{\odot }$ (the exact value of the limiting mass depending on metallicity) and offers a homogeneous grid of yields suitable for chemical evolution modelling applications.¹⁵

Other yield sets well-suited for chemical evolution studies can be retrieved by querying the Frascati Raphson Newton Evolutionary Code (FRANEC) Repository of Updated Isotopic Tables and Yields (FRUITY) database.¹⁶ The FRUITY web interface allows to download yield tables of all elements from H to Bi for stars with initial masses in the range 1.3–6 $M_{\odot }$, for a fine grid of initial metallicities ($Z = 2 \times 10^{-5}$, $5 \times 10^{-5}$, $1 \times 10^{-4}$, $3 \times 10^{-4}$, 0.001, 0.002, 0.003, 0.006, 0.008, 0.01, 0.014, and 0.02). Moreover, it is possible to explore the effects of different stellar rotational velocities ($v_{\mathrm {rot}} = 0$, 10, 30, or 60 km $\hbox {s}^{-1}$) and standard versus extended $^{13}\hbox {C}$ pockets on the final yields. The relevant stellar models are thoroughly discussed in the accompanying series of papers (see Cristallo et al. 2009, 2011, 2015, and references therein).

Last but not the least, the NuGrid stellar data set ought to be mentioned. Stellar models are computed with the mesa (Modules for Experiments in Stellar Astrophysics) stellar evolution code (Paxton et al. 2011) and nucleosynthesis is calculated in post-processing using the NuGrid mppnp code (see Pignatari et al. 2016). At present, the public database¹⁷ includes complete yields from H to Bi for $m =$ 1, 1.65, 2, 3, 4, 5, 6, 7, 12, 15, 20, and 25 $M_{\odot }$ and $Z = 1 \times 10^{-4}, 0.001, 0.006, 0.01,$ and 0.02 (for a detailed description of the models and yields, see Pignatari et al. 2016; Ritter et al. 2018; Battino et al. 2019, 2021). The use of the same stellar evolution and post-processing codes as well as the adoption of a common nuclear reaction network ensure a high degree of internal consistency. Relevant to the subject of this review, in low-mass and massive models at low metallicities H-ingestion in the He shell leads to substantial $^{15}\hbox {N}$ production (see Pignatari et al. 2015). Although this database has clear advantage of covering homogeneously different stellar mass realms, it has a major shortcoming in that it does not provide stellar yields for stars above 25 $M_{\odot }$.

2.1.2 Notes on the fate of the most massive AGB stars

Figure 2 illustrates the endpoints of the evolution of intermediate-mass stars with zero-age main-sequence masses between 5 and 10 $M_{\odot }$ and initial metallicities in the interval $Z =$ 0.0001$-$0.02 as reported by Doherty et al. (2015, top panel) and in the range $Z = 3 \times 10^{-7}$$-$0.04 after Ventura et al. (2013, 2014, 2018, 2020, 2021, and P. Ventura, private communication, bottom panel). The two sets of models are in reasonable qualitative agreement; both show that the minimum masses for the onset of carbon burning, formation of neutron stars, and explosion as SN increase with metallicity. However, in the aton models the ONe white dwarf (WD) progenitors are confined to a narrower range of initial masses.

Broadly speaking, stars in a narrow initial mass range¹⁸ climb the AGB with ONe cores after igniting carbon off-center under partially degenerate conditions (see previous section). Simulations can now follow the evolution of this class of objects all the way through the entire thermally pulsing super-AGB phase. This is a rather demanding computational task, mostly because of the extremely short time step, of the order of $10^{-3}$–$10^{-2}$ yr, dictated by H burning during the interpulse phases (Siess 2010), and the large number of thermal pulses expected (Doherty et al. 2017). The evolution of super-AGB stars resembles that of their lower-mass counterparts, however, owing to the massive ONe cores, they develop weaker thermal instabilities and may reach very high temperatures at the base of the convective envelope (up to $140 \times 10^7$ K, depending on the efficiency of convective energy transport, Siess 2010). This leads to very efficient HBB. Therefore, super-AGB stars may become powerful $^{13}\hbox {C}$, $^{14}\hbox {N}$, and $^{17}\hbox {O}$ forges, with yields heavily dependent on the temperature at the base of the convective envelope and on the rate of consumption of the mantle (hence, on the treatment of convection, see Siess 2010; Ventura and D’Antona 2011; Doherty et al. 2017).

Fig. 2

Top panel: final fates of intermediate-mass stars as a function of metallicity after Doherty et al. (2015). Each circle represents a model star and the different colours tell if the remnant is a CO WD, a CO-Ne WD, an ONe WD, or an ECSN. Bottom panel: same as the top panel, but for the stellar models computed by Ventura et al. (2013, 2014, 2018, 2020, 2021, and P. Ventura, private communication)

The fate of the most massive intermediate-mass stars is highly uncertain, as the velocity of erosion of the convective mantle and the rate of core growth deeply influence the outcome of the models (Poelarends et al. 2008; Ventura and D’Antona 2011; Woosley and Heger 2015; Leung et al. 2020). Stars with ONe core masses in excess of $\sim$1.35 $M_{\odot }$ (Schwab et al. 2016) ignite neon and end their lives as electron-capture SNe (ECSNe). Objects with smaller ONe cores die quietly, leaving cooling ONe WDs as remnants (e.g., Leung et al. 2020). Primordial (metal-free) stars in this mass interval could grow their degenerate cores up to the Chandrasekhar limit owing to the inefficient mass loss and explode as type 1.5 SNe (Gil-Pons et al. 2007; Poelarends et al. 2008). In summary, the competition between core growth (influenced by the third dredge-up efficiency) and envelope removal by mass loss (influenced by dredge-ups that change the envelope opacity and, hence, the strength of the stellar winds) dictates whether a star dies as an ONe WD or an ECSN. Poelarends et al. (2008) estimate that from about 3 to 20% of all SNe might be ECSNe. Doherty et al. (2015) and Jones et al. (2019) both point to low percentages, 1–5%. Since ECSNe are characterized by explosion energies and nickel production factors an order of magnitude lower than those of typical Fe CCSNe (Kitaura et al. 2006; Wanajo et al. 2009; Hiramatsu et al. 2021; Zha et al. 2022), a reliable estimate of their expected rate would be highly welcome as it would improve the modelling of the energetic and chemical feedback in galaxies.

2.1.3 High-mass stars

In the next paragraphs, I briefly recap the main facts regarding the evolution of massive stars that evolve in isolation. The section closes with some considerations about CNO element production from single stars that enter the main sequence with masses in excess of 10 $M_{\odot }$.

It is common wisdom by now that OB stars with initial masses between 8–12 and 20–25 $M_{\odot }$ evolve through the RSG phase, which is characterized by the presence of an extended convective envelope, and explode as type II-P SNe, leaving neutron stars as remnants and dominating the rate of explosions in the nearby universe (Smartt 2009; Langer 2012; Chieffi and Limongi 2013, and references therein). The evolutionary routes of stars above 20–25 $M_{\odot }$, instead, are still shrouded in mistery and there is no unanimous consensus on their final destiny (e.g., Langer 2012). Stars with initial masses in the range 40–120 $M_{\odot }$ possibly evolve through the WR phase sheding most of their H-rich envelopes and explode as type Ibc SNe, leaving behind either neutron stars or black holes (see Heger et al. 2003). Of particular interest are the evolution and fate of the most massive objects ($m > 80\,M_{\odot }$ and up to 350 $M_{\odot }$) that may die as pair-instability SNe (PISNe, also called pair creation SNe, Langer et al. 2007; Yusof et al. 2013; Leung et al. 2019). For a star to die as a PISN, it is necessary that it is left with a massive H-rich envelope at the pre-SN stage. Since the higher the metallicity, the more efficient the erosion of the outer layers by mass loss, there is a metallicity threshold above which essentially no PISNe are expected (but see Georgy et al. 2017, and references therein). Langer et al. (2007) estimate that at redshift $z = 5$ one should see one PISN per 100 SNe, while this fraction drops to 1/1000 in the local universe. As we will see in Sect. 3, including PISN nucleosynthesis in chemical evolution models might have important consequences for our understanding of the earliest phases of galactic evolution, especially with regard to early CNO element enrichment: the ejecta of low-metallicity very massive stars, in fact, are expected to be very rich in C and O (Goswami et al. 2021) and underabundant in N (Salvadori et al. 2019).

Single massive stars with initial masses between about 10 and 120 $M_{\odot }$ undergo all stages of nuclear burning to the point of Fe collapse. In the case of a successful conversion of the collapse into an explosion, the mantle experiences shock heating and explosive burning. In spherically symmetric 1D simulations, the transition from collapse to explosion is usually forced, e.g., by locating a piston¹⁹ somewhere in between the Fe core edge and the base of the O shell and setting the kinetic energy of the ejecta at infinity to a value inferred from observations: the piston is first moved inward, then rebounced generating a shock wave (see, e.g., Woosley and Weaver 1995; Woosley and Heger 2007). In the last decade, three-dimensional (3D) radiation hydrodynamic simulations of the collapse have succeeded in producing powerful explosions without any artifice, starting from 3D progenitor models. Spherical symmetry breaking, neutrino heating behind a stalled shock, and neutrino-driven turbulence²⁰ are the key drivers of the explosion of most CCSNe (Burrows and Vartanyan 2021), while the thermonuclear (e.g., Jones et al. 2019) and magnetically-driven (e.g., Burrows et al. 2007; Eisenberg et al. 2022) channels underlie a minor fraction of the global events, respectively, at the low- and high-mass ends.

Yields from massive stars Following Maeder (1992), the net yield of element j of a massive star of initial mass m consists of two contributions:

$$\begin{aligned} m_j^{\mathrm{{new}}} = m p_j = m p_j^{\mathrm{{wind}}} + m p_j^{\mathrm{{SN}}}. \end{aligned}$$

(3)

The first term on the right-hand side of Eq. (3) refers to the contribution from stellar winds,

$$\begin{aligned} m p_j^{\mathrm{{wind}}} = \int _0^{\tau (m)}[X_j(t) - X_j^{\mathrm{{init}}}] \, {\dot{M}}(t) \, \mathrm{d} t, \end{aligned}$$

(4)

while the second one generally accounts for the effects of explosive nucleosynthesis,

$$\begin{aligned} m p_j^{\mathrm{{SN}}} = \int _{m_{\mathrm{{remn}}}}^{{\widetilde{m}}(\tau )}[X_j(m') - X_j^{\mathrm{{init}}}] \,\mathrm{d} m'. \end{aligned}$$

(5)

In the above formulae, $\tau (m)$, ${\dot{M}}(t)$, $X_j^{\mathrm{{init}}}$, $X_j(t)$, $X_j(m')$, $m_{\mathrm{{remn}}}$, and ${\widetilde{m}}(\tau )$ are, respectively, the stellar lifetime, the assumed mass loss rate, the initial abundance of element j, the abundance of element j at the surface at any time, the abundance of element j at the Lagrangian mass coordinate $m'$, the mass of the remnant, and the remaining mass at age $\tau$. Again, to obtain the total amount of mass ejected in the form of an element j one needs to add to $m_j^{\mathrm{{new}}}$ the mass originally in the form of element j that is injected in the ISM without having suffered any nuclear processing inside the star [see Eq. (2)].

Quantitative estimates of the yields from numerical simulations of massive star evolution performed in systematic studies that span a large range of initial masses and chemical compositions are exactly what chemical evolution modellers need. Several teams, adopting different (1D) stellar evolutionary codes and different input microphysics, have provided homogeneous grids of massive star yields for use in GCE models. All these works have important limitations. In the following, we comment on a few that have been widely used in GCE studies.

Back in the nineties, Maeder (1992) investigated the effects of metallicity on the yields for two different initial chemical compositions, $Z =\,0.001$ and 0.02. His models include mass loss, but the evolutionary tracks of massive stars stop at the end of central C-burning. Moreover, only a few representative elements are considered. Woosley and Weaver (1995) computed the evolution and explosion of massive stars with initial metallicities $Z/Z_\odot = 0, 0.0001$, 0.01, 0.1, and 1 and masses in the range 11–40 $M_{\odot }$ and provided detailed yield tables for several isotopes. However, their nuclear reaction network is truncated at Zn, convective core overshooting is not taken into account and, more important for the purposes of the present discussion, mass loss is not considered. Thielemann et al. (1996) used a larger reaction network, but investigated only solar-metallicity stars, had a lower number of mass points in their grid and, again, did not consider the effects of mass loss during the pre-SN evolution. Limongi and Chieffi (2006, 2012) computed yields for solar-metallicity and metal-free stars with initial masses between 11 and 120 $M_{\odot }$ (a slightly reduced mass grid was employed for zero-metallicity stars). They took mass loss into account (see Limongi and Chieffi 2006, for the adopted mass loss prescriptions as well as other details of the code). In the meantime, the Geneva group published the first grids of yields including the effects of both mass loss and rotation (Meynet and Maeder 2002a; Hirschi et al. 2005; Hirschi 2007; Ekström et al. 2008). Their models are limited to the pre-SN evolution, but this should not be a serious issue for galactic CNO evolution studies (see below). More recently, Nomoto et al. (2013) presented a compilation of extant, homogeneous yields for CCSNe and hypernovae²¹ (HNe) from the Japanese team, to which they added unpublished, consistent yields for super-solar metallicity stars (see Nomoto et al. 2013, their Sect. 4.2). Nomoto et al. ’s models include the effects of mass loss, but not those of rotation. The upper mass limit is set to 40 $M_{\odot }$, apart from the grid of zero-metallicity models that extends to 100 $M_{\odot }$. Finally, Limongi and Chieffi (2018) presented new yields from a grid of models with initial metallicities [Fe/H] = $-3$, $-2$, $-1$, and 0 (corresponding to $Z = 3.236 \times 10^{-5}$, $3.236 \times 10^{-4}$, $3.236 \times 10^{-3}$, and $1.345 \times 10^{-2}$), initial masses $m =$ 13, 15, 20, 25, 30, 40, 60, 80, and 120 $M_{\odot }$, and equatorial velocities at the beginning of the main sequence $\upsilon =$ 0, 150, and 300 km $\hbox {s}^{-1}$; the model stars are susceptible to mass loss and some may explode as PISNe.

Fig. 3

Image reproduced with permission from Hirschi et al. (2005), copyright by ESO

IMF-weighted stellar yields (see text) of different elements for non-rotating (left-handed panel) and rotating (right-handed panel) solar-metallicity stellar models as a function of the initial mass, piled up on top of each other. The dotted areas highlight the contribution of the wind. Note that $^{16}\hbox {O}$ has negative SN yields for $m > 50\,M_{\odot }$ ($m > 30\,M_{\odot }$) in the non-rotating (rotating) case.

At this point, it is worth recalling that CNO elements in massive stars are synthesized mostly during hydrostatic burning and seldom suffer modification of their yields during explosive burning (e.g., Woosley and Heger 2007). Therefore, the yields depend crucially on the pre-SN evolution and, hence, on the assumed nuclear reaction rates, convection treatment, mixing processes, and mass loss rates, while having weak or no dependence on other uncertain parameters of the models, such as the mass cut.²² At high metallicities, in the presence of mass loss large amounts of He and C are vented out and escape further processing to oxygen and its descendants; as a consequence, high-metallicity high-mass stars have large C, but low O yields (Maeder 1992). In hypernovae, oxygen burning takes place at higher temperatures in more extended regions, which results in lower C and O abundances in the ejecta (Umeda et al. 2002). On the other hand, rotation may significantly increase the yield of oxygen. Figure 3 shows the stellar yields [computed according to Eq. (3)] from the solar-metallicity, non-rotating (left-handed panel) and rotating (right-handed panel) models of Hirschi et al. (2005), weighted by an extrapolated Salpeter ’s (1955b) IMF and multiplied by an arbitrary constant. Rotation clearly boosts the production of oxygen and, to a lesser extent, that of carbon, a result qualitatively confirmed by Limongi and Chieffi (2018). At low metallicities, mass loss plays a minor role, but rotation-induced instabilities trigger the interplay between the core He-burning and the H-burning shell that results in a significant primary production of $^{13}\hbox {C}$, $^{14}\hbox {N}$, $^{17}\hbox {O}$, and $^{18}\hbox {O}$ (Hirschi 2007, see also Limongi and Chieffi 2018). This, as we will see in Sect. 3, has a dramatic impact on the predictions of chemical evolution models.

In closing this section, it is worth mentioning the possibility that a large fraction of massive stars with primordial, or nearly primordial, metallicity explode with lower explosion energies, i.e., as faint SNe. In this case, higher [(C + N)/Fe] and [(C + N)/Mg] ratios, corresponding to smaller ejected Fe masses and to larger compact remnant masses, have to be expected (e.g., Heger and Woosley 2010; Tominaga et al. 2014) with important measurable consequences for the chemical composition of the stars at extremely low metallicities. More details on this issue will be given in Sect. 3.1.2.

2.2 Binary systems

So far, I have discussed the evolution of stars that spend their lives in isolation. However, most stars form in clusters and associations (Lada and Lada 2003), where interactions are frequent. In this section, I focus on binary systems in which close interaction may boost CNO element production. A thorough review of binary classes and their pathways can be found in De Marco and Izzard (2017).

Most young Galactic stars with $m \ge$ 15 $M_{\odot }$ interact with a companion before exploding as CCSNe (Sana et al. 2012), which modifies their evolutionary paths, nucleosynthesis, SN appearance, and remnant properties (see, e.g., Langer 2012; Modjaz et al. 2019; Woosley 2019; Schneider et al. 2021). Podsiadlowski et al. (2004) investigated how the presence of a close companion affects the final core structure and terminal evolution of a massive star and found, among other effects, that ECSN explosions are favoured in binary systems.

The nucleosynthetic outputs of massive binary interactions haven’t been scrutinized in so much detail as to allow implementation in GCE codes yet, which is understandable. Each binary component, in fact, brings with it all the uncertainties on the treatment of convection, extra mixing, mass loss, and microphysics discussed in the previous sections. Duplicity adds further ambiguity linked, for example, to the details of mass transfer and leakage by Roche-lobe overflow or to the effects of tides on rotation velocity evolution and internal mixing (Song et al. 2016). The separation and the mass ratio between the primary and the secondary component of the system must also be specified. All of this makes the parameter space to expand significantly.

Fig. 4

Image adapted from Farmer et al. (2021), copyright by the authors

Net $^{12}\hbox {C}$ yield including the contributions from stellar winds and SN explosions as a function of the initial stellar mass. Blue lines and areas denote the results from single-star models, while red lines and areas represent binary models. The $^{12}\hbox {C}$ yields are dominated by SN explosions in the range 11–22 $M_{\odot }$, whereas winds dominate above 35 $M_{\odot }$.

Farmer et al. (2021) have calculated the C yields of solar-metallicity massive stars that either evolve alone or are stripped in binary systems during H burning in a shell (notice that this is only one of the many possible fates of massive stars with close companions). In binary-stripped stars the convective He-burning core does not grow (it may rather retreat) and the outermost layers disconnect forming a $^{12}\hbox {C}$-rich pocket that never reaches temperatures high enough for $^{12}\hbox {C}$ burning. In single stars, instead, these layers get mixed into the growing core where $^{12}\hbox {C}$ is readily destroyed by $\alpha$ captures first, and C burning later. As a consequence, binaries are more efficient in enriching the ISM in $^{12}\hbox {C}$ than single stars (see Fig. 4). However, the authors themselves caution that different choices of wind mass-loss rates could have a not negligible impact on the results. The $^{12}\hbox {C}$ yields are, instead, remarkably robust against changes in the explosion parameters (see also Sect. 2.1.3). The main source of uncertainty (in both single and binary models) appears to be the treatment of convective boundary mixing. Yields of other species should become available in the future, which will open new avenues to our knowledge of how chemical enrichment proceeds in the cosmos taking stellar duplicity into account.

Regarding lower mass systems, it has to be emphasized that many solar-type and intermediate-mass stars have companions, but the interaction is less frequent than in the case of more massive stars (see, e.g., Kouwenhoven et al. 2005; Raghavan et al. 2010; Moe and Di Stefano 2017). This notwithstanding, some low-mass close binaries may lead to spectacular events that have been, and still are, the subject of a great deal of theoretical work. In the remainder of this section, I focus on classical novae that prove to be efficient manufacturers of some rare CNO isotopes.

Classical nova explosions occur in close binary systems where a WD primary accretes H-rich matter from a secondary, less evolved star that overflows its Roche lobe. Once the temperature of the deepest layers accreted on top of the WD exceeds about $10^8$ K, a thermonuclear runaway is initiated (see Gurevitch and Lebedinsky 1957, and references therein). The resulting outburst does not disrupt the WD and, after the ejection of the envelope, the process is re-initiated. On average, the typical nova experiences $10^4$ outbursts during its lifetime (Bath and Shaviv 1978; Ford 1978).

Fig. 5

Image adapted from Romano et al. (2017), copyright by the authors

Yields of $^{13}\hbox {C}$, $^{15}\hbox {N}$, and $^{17}\hbox {O}$ from hydrodynamic models of the nova outburst for CO (triangles) and ONe WDs (circles). Results from representative 1D models (José and Hernanz 1998, coloured symbols; Yaron et al. 2005, black symbols; Starrfield et al. 2009, small empty circles) are shown together with yields from simulations that switch from one to three dimensions, and viceversa (123–321 models, José et al. 2020, large empty symbols). In the latter suite of models, the metal enhancement of the envelope is not forced but naturally produced. The big $\times$ signs represent the average yields adopted in chemical evolution simulations by Romano et al. (2017).

In their pioneering hydrodynamical studies, Starrfield et al. (1972, 1978), Prialnik et al. (1978), and Prialnik (1986) developed models of the nova outburst and showed that for a successful outburst to develop, the envelope must be enriched in CNO elements. Moreover, only massive WDs produce outbursts: Kovetz and Prialnik (1985) find that models with $m_{\mathrm{{WD}}} = 0.6\,M_{\odot }$ do not lead to a nova outburst, but instead settle into a thermal equilibrium configuration with a luminosity typical of a red giant. In the models of Starrfield et al. (1972), the CNO cycle is always far out of equilibrium in the outburst and the $\beta ^+$-unstable nuclei $^{13}\hbox {N}$, $^{14}\hbox {O}$, $^{15}\hbox {O}$, and $^{17}\hbox {F}$ become very abundant. Convection carries these unstable nuclei to the surface where they decay to $^{13}\hbox {C}$, $^{14}\hbox {N}$, $^{15}\hbox {N}$, and $^{17}\hbox {O}$ and the temperatures are too low for further p-captures, which affects the composition of the ejecta. José and Hernanz (1998) performed a thorough investigation of the role of classical novae as ISM polluters, adopting the 1D implicit hydrodynamic code SHIVA to follow the course of nova outbursts from the onset of accretion up to the ejection of the processed envelope. Fourteen models were computed, spanning WD masses from 0.8 to 1.35 $M_{\odot }$ and mixing levels between the accreted envelope and the underlying WD from 25% to 75%. Both CO WDs and ONe WDs were considered. To limit the parameter space, the rate of mass accretion and the initial luminosity were kept fixed, at $2 \times 10^{-10}\,M_{\odot }\, \mathrm {yr}^{-1}$ and $10^{-2}$ $\hbox {L}_\odot$, respectively. From these simulations, novae turn out to be significant producers of $^{13}\hbox {C}$, $^{15}$N, and $^{17}\hbox {O}$, possibly contributing most of their Galactic abundances (Romano and Matteucci 2003).

Over the last twenty years, multi-dimensional studies of mixing at the core-envelope interface during outbursts have demonstrated that Kelvin-Helmholtz instabilities develop, naturally resulting in an effective mixing of the outermost WD material into the accreted layers. This makes the envelope metallicity to achieve values in agreement with the observed ones, namely, $Z \sim 0.20$–0.30 (e.g., Glasner et al. 1997; Casanova et al. 2010, 2016; José et al. 2020), which ensures powerful nova outbursts. Simulations implementing—or guided by—these results have started to produce more coherent and self-consistent nova yields (José et al. 2020; Starrfield et al. 2020).

Figure 5 compares the masses ejected in the form of $^{13}\hbox {C}$, $^{15}$N, and $^{17}\hbox {O}$ in a single outburst from different nova models. We display both results from 1D simulations (José and Hernanz 1998; Yaron et al. 2005; Starrfield et al. 2009) and results from models that combine 1D and 3D hydrodynamics (123–321 models, José et al. 2020). Yields computed with different codes and specifics may differ by up to four orders of magnitude. The $\times$ signs denote the values adopted in the GCE models of Romano et al. (2017). A detailed discussion of nova nucleosynthesis implementation in GCE models is deferred to Sect. 3.1.1.

3 Modeling the evolution of CNO elements in galaxies

From all the above it is clear that the stellar yields bring with them considerable uncertainties. Nonetheless, they are fundamental ingredients of GCE models (actually, their most important inputs; see Romano et al. 2010; Côté et al. 2017).

GCE models follow the variation of the chemical composition of the ISM in space and time in galaxies taking into account various processes, such as gas inflows and outflows, star formation, feedback from stars, and radial motions of gas and stars (see Matteucci 2021, and references therein). They rely on empirical formulae and semi-analytical, yet physically-motivated laws to describe complex processes. This makes the computational cost affordable even in the case of a full screening of the parameter space. This kind of models can thus be used, among other things, to efficiently test different grids of stellar yields, providing useful hints to more elaborate simulations—it would be rather unpleasant to run expensive chemodynamical simulations first, then discover that the adopted stellar nucleosynthesis prescriptions are in need of a radical overhaul!

The theoretical abundances and abundance ratios likewise depend on the adopted galaxy-wide stellar IMF (gwIMF), which is used to weigh the yields. A broad consensus has emerged over the years that the first stars – called, somewhat counterintuitively, Pop III stars – formed with a top-heavy IMF; chemical enrichment by the first SNe then suddenly turned the star formation mode from high-mass to low-mass dominated (e.g., Karlsson et al. 2013). While such claims have solid theoretical foundations, they still lack an unambiguous observational confirmation. In fact, it is found that low-mass stars can form at very low metallicity ($Z < 7 \times 10^{-7}$; Caffau et al. 2011a). Also, the average abundance patterns of Milky Way halo stars—as well as those of some weird objects—do not necessitate a peculiar gwIMF to be reproduced (e.g., Limongi et al. 2003). From the point of view of chemical evolution, Pop III stars, if ever existed, produced negligible effects in the earliest phases of Galactic evolution (Ballero et al. 2006). Possible systematic variations of the gwIMF in dependence on a galaxy’s metallicity and star formation rate have been proposed (Jerabkova et al. 2018; Yan et al. 2020, 2021, and references therein) that deserve further investigation.²³

In particular, it has been suggested that the abundance ratios of some specific CNO isotopes may provide an effective litmus test of possible gwIMF variations in starburst galaxies, once reliable measurements are analysed by means of GCE models (Henkel and Mauersberger 1993; Papadopoulos et al. 2014; Romano et al. 2017; Zhang et al. 2018, see Sect. 3.2). This procedure, however, requires an a-priori careful calibration of the models against local CNO abundance data. This is the subject of the following section.

3.1 The local universe as a benchmark

In this section, we deal with the evolution of the CNO elements and their isotopes in the Milky Way and other selected galaxies in and beyond the Local Group, up to about 20 Mpc distance from us. The latter are grouped in galaxies with (Sects. 3.1.3 and 3.1.4) and galaxies without (Sect. 3.1.2) a detectable gas component at present.

3.1.1 The Milky Way galaxy

In a number of galaxies in the Local Group it is possible to resolve single stars. In particular, in our own Galaxy high-resolution spectra of unevolved stars can be analysed, which provides detailed, high-precision CNO abundances reflecting the chemical composition of the gas out of which the stars were born. Low-mass stars formed at different epochs and still alive today can thus be used to trace the evolution of CNO elements in the Galaxy from its earliest formation phases to the present time.

CNO abundances in stars can be determined from permitted atomic lines of C I and O I, forbidden [C I] and [O I] lines, and molecular bands of $\hbox {C}_2$, CH, NH, OH, CN, and CO in the UV, optical and NIR (Barbuy et al. 2018; Randich and Magrini 2021, and references therein). Significant efforts are ongoing to obtain accurate and precise abundances for, ideally, millions of stars in the Galaxy representing the variety of its stellar populations, from the inner bulge to the disc outskirts (see Jofré et al. 2019, and references therein). Coupling these abundance measurements with the exquisite astrometric data provided by the European Space Agency (ESA) Gaia mission (Perryman et al. 2001; Gaia Collaboration et al. 2016) and with precise stellar ages (Gallart et al. 2019; Miglio et al. 2021) allows an unprecedented multidimensional mapping of stars within about 10 kpc from the Sun. These high-precision positions, motions, abundances, and ages can be effectively used to constrain galaxy formation and evolution scenarios as well as stellar evolution and nucleosynthesis theory (Freeman and Bland-Hawthorn 2002; Nissen and Gustafsson 2018).

In the last fifteen years, large ground-based spectroscopic surveys, such as the Gaia-ESO Survey (GES, Randich et al. 2013, 2022; Gilmore et al. 2022), the Apache Point Galactic Evolution Experiment (APOGEE, Majewski et al. 2017), the Large sky Area Multi-Object fiber Spectroscopic Telescope (LAMOST) survey (Zhao et al. 2012), and the GALactic Archaeology with HERMES survey (GALAH, De Silva et al. 2015), have measured detailed properties of several hundred thousand stars in the Milky Way discs, finding that the stars populate two distinct sequences in the [$\alpha$/Fe]–[Fe/H] space (Mikolaitis et al. 2014; Nidever et al. 2014; Recio-Blanco et al. 2014; Hayden et al. 2015; Yu et al. 2021). The observed dichotomy confirms and extends to a larger portion of the disc original findings by Fuhrmann (1998), who found the thin- and thick-disc components to be chemically disjunct from the analysis of a sample of solar neighbourhood stars (see also Gratton et al. 1996).

Pure chemical evolution models for the Milky Way disc have been developed assuming either a sequential or a parallel formation of its components (for a thorough review, see Sect. 5 of Matteucci 2021). The two-infall model of Chiappini et al. (1997, 2001), for instance, assumes two main episodes of gas infall, out of which the thick and thin discs form in a sequence. Spitoni et al. (2019) fine-tuned this model by leveraging recent constraints on stellar ages from asteroseismology to reproduce the [$\alpha /\hbox {Fe}$] versus [Fe/H] relationship at different Galactic epochs, the age-metallicity relation, and the metallicity and age distributions of the stars populating the high-$\alpha$ and low-$\alpha$ sequences. A crucial assumption of their revised two-infall model is the longer time delay imposed for the start of the second gas infall event, i.e., $\sim 4$ Gyr as opposed to 1 Gyr in the original version. The parallel approach developed by Grisoni et al. (2017, 2018) envisages instead the formation of the thin and thick discs out of two distinct accretion episodes on different timescales. Also in this case, the parallel sequences in the [$\alpha /\hbox {Fe}$]–[Fe/H] plane and the stellar metallicity distributions of the thin- and thick-disc populations can be reproduced. Chemically distinct discs are also predicted in the context of chemodynamical simulations of the formation and evolution of Milky Way-sized galaxies. In some models the two components emerge from separate evolutionary pathways (e.g., Grand et al. 2018; Mackereth et al. 2018; Agertz et al. 2021; Khoperskov et al. 2021), whereas others propose a continuous star formation history and stellar migration to generate the double sequence (Sharma et al. 2021, and references therein). To fine tune the free parameters relevant to the chemical enrichment process, the stellar metallicity distributions and [$\alpha /\hbox {Fe}$]–[Fe/H] trends are commonly used. The evolution of the CNO elements and their isotopes are rarely investigated. This is partly due to the inherent difficulties in obtaining trustworthy CNO element abundances that reflect the ISM enrichment processes.

The Sun It is sobering, in fact, that even the solar photospheric abundances of the CNO elements have been revised several times over the years. Lambert (1978) suggested A(O)$_\odot = 8.92 \pm 0.10$²⁴, $A\hbox {(C)}_\odot = 8.67 \pm 0.10$, and $A\hbox {(N)}_\odot = 7.99 \pm 0.10$, using as primary abundance indicators the [C I] $\lambda$ 8727 Å line, the CH A-X and $\hbox {C}_2$ Swan systems, the N I lines and the forbidden atomic oxygen lines at 6300 and 6363 Å. In their prominent review, Anders and Grevesse (1989) recommended a value of A(O)$_\odot = 8.93 \pm 0.04$. The solar O abundance has undergone major downwards revisions since, culminating with the low value of A(O)$_\odot = 8.66 \pm 0.05$ suggested by Asplund et al. (2004). Later on, Ayres et al. (2006) pointed again to a high value, $A\hbox {(O)}_\odot = 8.84 \pm 0.06$, which was revised downwards to $A\hbox {(O)}_\odot = 8.76 \pm 0.07$ by Caffau et al. (2008). Recommended literature values for the solar photospheric C abundance similarly oscillate from high to low figures—for instance, $A\hbox {(C)}_\odot = 8.592 \pm 0.108$ in Holweger (2001), $A\hbox {(C)}_\odot = 8.39 \pm 0.05$ in Asplund et al. (2005), and $A\hbox {(C)}_\odot = 8.50 \pm 0.06$ in Caffau et al. (2010). The determination of the solar abundance of nitrogen also suffered from important downwards revisions: Anders and Grevesse (1989) list $A\hbox {(N)}_\odot = 8.05 \pm 0.04$, Holweger (2001) recommends $A\hbox {(N)}_\odot = 7.93 \pm 0.11$, while Caffau et al. (2009) favour $A\hbox {(N)}_\odot = 7.86 \pm 0.12$, to cite only a few relevant references. Incorrect or incomplete modelling of the different lines considered in different studies offers a likely explanation for the discrepancies (see, e.g., Caffau et al. 2013).

The chemical portrait of the Sun by Asplund et al. (2021), based on 3D spectral line formation calculations also taking into account departures from local thermodynamic equilibrium (non-LTE), ends up with a set of present-day photospheric abundances, $A\hbox {(O)}_\odot = 8.69 \pm 0.04$, $A\hbox {(C)}_\odot = 8.46 \pm 0.04$, and $A\hbox {(N)}_\odot = 7.83 \pm 0.07$, that nicely agree within the mutual errors with previous 3D-based studies ( Caffau et al. 2008, 2009, 2010, see the above-quoted values). In their most recent work, Magg et al. (2022) recommend $A\hbox {(O)}_\odot = 8.77 \pm 0.04$, $A\hbox {(C)}_\odot = 8.56 \pm 0.05$, and $A\hbox {(N)}_\odot = 7.98 \pm 0.08$, again in excellent agreement with the results by Caffau and collaborators. The concordance of these latest studies is very encouraging.

Concerning the minor isotopes, values of $91.4 \pm 1.3$, $2738 \pm 118$, and $511 \pm 10$ (at 1 $\sigma$) are inferred for the $^{12}\hbox {C}$/$^{13}\hbox {C}$, $^{16}\hbox {O}$/$^{17}\hbox {O}$, and $^{16}\hbox {O}$/$^{18}\hbox {O}$ ratios, respectively, from IR (2–6 $\mu \hbox {m}$) rovibrational bands of CO in the solar spectrum²⁵ using 3D convection models (Ayres et al. 2013). The nitrogen isotopic ratio presents outstanding variability on the Solar System scale, which can be explained by invoking fractionation processes (Füri and Marty 2015). For the protosolar nebula, a low ratio, $^{14}\hbox {N}/^{15}\hbox {N}$ $\simeq 441 \pm 6$, is suggested by Marty et al. (2011).

Other stars Modern 3D line formation models including non-LTE effects provide solar photospheric abundances consistent with helioseismic constraints (Magg et al. 2022, see also Pinsonneault and Delahaye 2009; Caffau et al. 2011b). Consideration of 3D non-LTE effects is even more crucial in the low-metallicity regime. Amarsi et al. (2019) derived homogeneous 3D non-LTE carbon, oxygen, and iron abundances for a sample of 187 F and G dwarfs belonging to the Galactic discs and halo. The inferred abundance trends are tighter than—and, at low metallicities, deviate markedly from—those obtained in the classical 1D LTE approximation. Thick-disc stars display [C/Fe] and [O/Fe] ratios higher than those of their thin-disc counterparts, while the opposite is true in the [C/O]–[O/H] plane (see also Nissen et al. 2014, but see Bensby and Feltzing 2006, for different results). Notably, the upturn of [C/O] for metallicities below [O/H] $\simeq -$1 spotted out in 1D LTE analyses (Akerman et al. 2004) is not recovered by 3D non-LTE line formation models (see Amarsi et al. 2019, their Fig. 13). This has profound implications for the characterization of the stellar factories responsible for the synthesis of C and O at early times.

Apart from the inclusion or omission of 3D and/or non-LTE effects in the derivation of the abundances, other systematics may affect the data, impacting the comparison between observations and theoretical predictions and, hence, the conclusions of GCE models. As an example, in Fig. 6 we display different datasets in the [C/Fe]–[Fe/H] and [C/Mg]–[Mg/H] planes (top and bottom panels, respectively). The red and blue symbols refer to abundance measurements from high-resolution, high signal-to-noise ratio (SNR) spectra of single thin- and thick-disc stars, respectively, presented in Bensby and Feltzing (2006), Bensby et al. (2014), Nissen and Schuster (2010), Nissen et al. (2014), Amarsi et al. (2019). The solid red and blue lines represent the average trends derived by Franchini et al. (2020) from, respectively, about 1300 thin-disc and nearly 100 thick-disc dwarf stars from the fifth internal data release of GES (GES iDR5). The dashed red and blue lines have similar meanings, but refer to the analysis of more than 12,000 stars from the second data release of GALAH (GALAH DR2), with the stars chemically divided into high-Ia and low-Ia sequences based on the higher/lower Fe contribution from SNeIa in the thin-/thick-disc phases (Griffith et al. 2019). Only unevolved stars with C abundances not altered by evolutionary processes are considered. To remove at least partly systematic offsets, zero-point corrections are applied to each dataset, so that thin-disc stars with [X/H] $\simeq 0$ also have [C/X] $\simeq 0$ (where X stands for either Fe or Mg).

Fig. 6

Image reproduced with permission from Romano et al. (2020), copyright by ESO

Carbon-to-iron (top panel) and carbon-to-magnesium (bottom panel) abundance ratios as functions of [Fe/H] and [Mg/H], respectively. Red and blue symbols are abundance estimates for thin- and thick-disc unevolved stars (triangles: Nissen and Schuster 2010; Nissen et al. 2014; squares: Bensby and Feltzing 2006; Bensby et al. 2014; circles: Amarsi et al. 2019). Solid red and blue lines are abundance ratio trends of, respectively, $\sim 1,300$ thin- and $\sim 100$ thick-disc dwarf stars selected on the basis of chemical criteria from GES iDR5 (Franchini et al. 2020). Dashed red and blue lines are abundance ratio trends of more than 12,000 high-Ia and low-Ia stars, respectively, from GALAH DR2 (Griffith et al. 2019). First-order zero-point offsets (see text) are applied to the data. Pale red and blue areas show the spread in thin- and thick-disc GES iDR5 abundance data within boundaries corresponding to the 10th and 90th percentiles. Red and blue dashed areas show the analogous spread in the GALAH DR2 thin- and thick-disc datasets (bottom panel only).

From Fig. 6, top panel, it appears that, in general, the agreement between the offset [C/Fe] versus [Fe/H] trends traced by GES and GALAH data is very good. Moreover, there is an overall agreement with the area covered by individual data points from the other high-resolution spectroscopic studies, especially if the spread in the survey data (here only for the GES dataset—pale red and blue areas for thin- and thick-disc stars, respectively) is taken into account. The rather large extension of the two scatter areas, $\pm 0.2$ dex, is likely due to a combination of uncertainties in the derived abundance ratios—some of the GES spectra have relatively low SNR—and cosmic scatter. The fact that the individual points from the selected high-resolution studies, which are affected by much smaller uncertainties, span similar wide areas indicates that the cosmic scatter is significant for both thin- and thick-disc stars. Figure 6, bottom panel, with Mg as the metallicity tracer, points to a less favourable figure: the GES [C/Mg] trends are flatter than the GALAH ones for both discs and agree better with the loci of the individual observations, especially for the thick-disc component. The spreads in GES and GALAH data are comparable in size (when also taking the different number of stars in the two samples into account) and resemble that of the individual points, thus supporting the presence of a cosmic scatter also in [C/Mg]. In this respect, it is worth noting that examination of [O/Fe] abundance ratios of Galactic disc stars with H-band spectra from APOGEE reveals a star-to-star cosmic variance at a given metallicity of about 0.03$-$0.04 dex in both the thin- and thick-disc components that is tentatively associated with the wide range of galactocentric distances spanned by APOGEE targets (Bertran de Lis et al. 2016).

It is, therefore, interesting to analyse the behaviour of stars that span a narrower distance range. Delgado Mena et al. (2021) have derived homogeneous C abundances for a volume-limited sample (most targets are within 60 pc distance from us) of 757 FGK dwarfs with $T_{\mathrm {eff}} > 5200$ K and $-1.2<$ [Fe/H] $< 0.6$. The stars have been observed at very high resolution ($R \sim$ 115,000) with the High Accuracy Radial velocity Planet Searcher (HARPS) spectrograph (Mayor et al. 2003) mounted on the ESO La Silla 3.6 m telescope. The sample does not contain any fast rotators and active stars. Two well-known atomic C lines (C I $\lambda$ 5052 $\text{\AA}$ and C I $\lambda$ 5380 $\text{\AA}$) have been subjected to a standard LTE analysis. The two lines provide consistent C abundances in the considered $T_{\mathrm {eff}}$ range (see Fig. 1 of Delgado Mena et al. 2021). Oxygen abundances have also been derived using either the O I $\lambda$ 6158 $\text{\AA}$ line or the [O I] $\lambda$ 6300 $\text{\AA}$ line. The usage of the first O indicator results in a tighter [C/O] versus [Fe/H] relation, suggesting that the O abundances derived from the permitted line at 6158 $\text{\AA}$ are more reliable. For the hot dwarf stars analysed by Delgado Mena et al. (2021), the forbidden line at 6300 $\text{\AA}$ provides systematically higher carbon-to-oxygen ratios.

Fig. 7

Circles: carbon-to-iron (top panel) and carbon-to-oxygen ratios (bottom panel) as functions of metallicity for neighboring stars ($d \le$ 60 pc) from Delgado Mena et al. (2021). Carbon abundances are from two well-known atomic lines; oxygen abundances have been derived from the permitted O I $\lambda$ 6158 Å line (see text). Different colours indicate different Galactic populations (legend in the top panel). The observational sample comprises only stars for which the C and O abundances are measured reliably, i.e., dwarf stars hotter than 5200 K. Solid lines are theoretical predictions for the solar vicinity from Romano et al. (2020), their model MWG-11 rev. The empty squares along the theoretical tracks highlight the abundance ratios at the specified ages. Stars older than 9.4 Gyr belong to the thick-disc component (see Spitoni et al. 2019). Data and models are normalized to solar abundances from Delgado Mena et al. (2021)

In Fig. 7, the dataset of Delgado Mena et al. (2021) is compared to the predictions of one of the GCE models for the solar vicinity discussed in Romano et al. (2020), i.e., their model MWG-11 rev. The observed stars are assigned to different groups—thin-disc, thick-disc, halo, and high-$\alpha$ metal-rich stars—based on chemical criteria (see Delgado Mena et al. 2021, and references therein). The model displayed in the figure is a two-infall model, assuming that the two disc components form in a sequence. In particular, it assumes a delay of 4.3 Gyr for the beginning of the second accretion event that originates the thin disc, based on recent constraints on the stellar ages from asteroseismology (Spitoni et al. 2019). The adopted stellar yields²⁶ guarantee a satisfactory fit to the data. From Fig. 7, top panel, it is seen that the model reproduces very well the separation between thin- and thick-disc stars in the [C/Fe]–[Fe/H] plane, the decreasing [C/Fe] ratio at super-solar metallicities, and the subsolar [C/Fe] ratios for solar-metallicity stars (also pointed out by other authors; see, e.g., Botelho et al. 2020; Franchini et al. 2020). Furthermore, the model naturally explains the [C/Fe] versus [Fe/H] pattern of the high-$\alpha$ metal-rich stars, in the hypothesis that they represent the metal-rich tail of the thick disc as proposed by Bensby et al. (2014). A metal-rich tail of thick-disc stars (according to a kinematic classification) has been also found in other high-resolution studies (e.g., Amarsi et al. 2019). Some shortcomings are present though: model MWG-11 rev does not reproduce the observed [C/Fe] bump for thick-disc stars in the $-0.8<$ [Fe/H] $< -0.4$ range, it predicts [C/O] ratios higher than observed (by 0.10 dex at high metallicities and up to 0.4 dex in the metal-poor regime), and it does not reproduce the flattening²⁷ of the carbon-to-oxygen ratio at super-solar metallicities (Fig. 7, bottom panel). In this respect, we note that the adopted massive star yields (Limongi and Chieffi 2018) are computed for four initial metallicity values, [Fe/H] $= -$3, −2, −1, and 0. Therefore, when the ISM metallicity computed by the GCE code exceeds the solar value, one has to extrapolate arbitrarily the solar-metallicity grid. The inclusion of yields tailored to super-solar metallicity stars could improve the model performances. In particular, the implementation of suitable mass loss rates is crucial in determining the final outcome of the models, including the yields (see Sect. 2.1). Our results may imply that a milder mass loss is needed at super-solar metallicities, resulting in lower (higher) yields for carbon (oxygen). In the metal-poor domain, the inclusion of HNe (Limongi and Chieffi 2018 consider only normal CCSNe) could improve the agreement with the data. Another intriguing possibility would be to take into account the existence of very massive stars exploding as PISNe in the early Galaxy. Such a possibility is explored by Goswami et al. (2021) for a number of elements comprising oxygen, but not for carbon.

In the parallel Galaxy formation scenario (Grisoni et al. 2017, 2018), thin- and thick-disc stars are predicted to populate distinct sequences in the [C/Fe]–[Fe/H] and [C/O]–[O/H] spaces (see Romano et al. 2020, their Fig. 4) but, at variance with the sequential framework discussed above, they may have a common range of ages. While precise stellar ages for larger samples of stars and an exact quantification of the effects of radial migration are required to allow discriminating among different evolutionary pathways, pure GCE models can still be used to provide some indications of the robustness of the adopted stellar yields. In fact, when considering the run of abundance ratios with metallicity the effects of changes in several free parameters of the GCE models tend to cancel out, while those of different nucleosynthesis prescriptions are maximized (Tosi 1988).

To that end, it is important to notice that different C diagnostics produce different trends (e.g., Bensby and Feltzing 2006; Suárez-Andrés et al. 2017). It has become clear by now that it is not a good idea to mix different kinds of data and try to fit the grossly averaged trend; rather, efforts must be made to understand any systematics affecting the data, to define meaningful trends and better constrain the models (see, e.g., discussions in Jofré et al. 2019; Romano et al. 2020; Delgado Mena et al. 2021).

The above considerations apply, of course, also to other elements. Nitrogen abundances in both the [N/Fe] versus [Fe/H] and [N/O] versus [O/H] planes are characterized by a large scatter at low metallicities, even if only unmixed stars²⁸ are considered (Spite et al. 2005). The dispersion could be due to inhomogeneous enrichment during the early Galaxy assembly, either from massive stars with different initial rotational velocities (Matteucci 1986; Chiappini et al. 2006; Prantzos et al. 2018; Romano et al. 2019b) or from AGB stars (Kobayashi 2022, and references therein). Spite et al. (2005) further notice that N determinations from the NH band and from the CN band at 388.8 nm in giants show systematic differences of 0.4 dex that they attribute to physical parameters of the NH band. The scatter in N measurements sensibly reduces at disc metallicities, where reliable N abundances are derived homogeneously for statistically significant samples of dwarf stars from the NH molecular band at 336 nm (e.g., Suárez-Andrés et al. 2016). We note that, at present, high-resolution spectra of dwarf stars including the NH feature can be obtained, basically, only with one instrument, the Ultraviolet and Visual Echelle Spectrograph (UVES, Dekker et al. 2000) mounted on the ESO Very Large Telescope (VLT).

Atomic lines suitable to O abundance derivation are scarce in late-type stars. Moreover, each of them presents specific challenges. The strong O I IR triplet lines at $\lambda$ 777 nm are easily accessible and free of blends, but highly sensitive to 3D non-LTE effects (e.g., Cavallo et al. 1997; Caffau et al. 2008; Amarsi et al. 2019). Other features (the O I $\lambda$ 6158 Å and [O I] $\lambda$ 6363 Å lines) are safe from non-LTE effects, but are weak and difficult to measure. High-resolution, high SNR spectra permit the derivation of trustworthy abundances, though (Bertran de Lis et al. 2015; Delgado Mena et al. 2021). The forbidden [O I] $\lambda$ 6300 Å line is very weak at low metallicities; it can also be contaminated by telluric features and severely blended by two isotopic Ni I lines (Lambert 1978). The two forbidden [O I] lines also contain weak CN blends (e.g., Skúladóttir et al. 2015). The picture is no better as regards the numerous OH bands in the near UV and NIR. Although they are strong enough to be measured in main-sequence stars down to metallicities [Fe/H] $\simeq -$3 (Israelian et al. 1998, 2001; Boesgaard et al. 1999), the extraction of the abundances is not that straightforward, because of poorly-understood physical processes affecting the lines.

Fig. 8

Oxygen-to-iron ratio as a function of metallicity in the solar vicinity predicted by model MWG-11 rev (solid black line) and inferred from observations of atomic lines (green circles) and molecular bands (red circles) in high-resolution stellar spectra. The empty squares along the theoretical track highlight the abundance ratios at the labeled ages (13, 12, 9.4, 7.5, 4,5, and 1 Gyr). Data and models are normalized to solar abundances from Amarsi et al. (2019)

This notwithstanding, it has been claimed that GCE models can reproduce fairly well the main features of the [O/Fe] versus [Fe/H] diagram of solar neighbourhood stars (e.g., Romano et al. 2010; Prantzos et al. 2018; Kobayashi et al. 2020). In Fig. 8, we show the run of [O/Fe] as a function of [Fe/H] predicted by model MWG-11 rev of Romano et al. (2020). The theoretical trend is compared to O abundance determinations for two samples of unevolved stars, resting on a set of atomic oxygen lines the one (Amarsi et al. 2019) and on the UV molecular OH band the other (Israelian et al. 1998). The atomic lines have been subjected to a full 3D, non-LTE analysis. While at relatively high metallicities, [Fe/H] $> -1.0$, the two datasets almost overlap and the theoretical predictions are in reasonable agreement with the observations, in the low-metallicity regime the linear increase of the [O/Fe] ratio with decreasing [Fe/H] drawn from the OH data sharply contrasts with the flattening of the ratio traced by the 3D, non-LTE analysis of the atomic lines. The model predictions during the earliest phases of Galactic evolution (it takes only $\sim 350$ Myr to reach [Fe/H] $= -1$ in the local disc) are marginally consistent with both the flat trend and the steep increase of [O/Fe] with decreasing [Fe/H] suggested by Amarsi et al. (2019) and Israelian et al. (1998), respectively. Actually, if we were to consider the data all together, we would be led to happily conclude that the gross [O/Fe] versus [Fe/H] trend is well reproduced; on closer inspection, however, the inconsistency of the two datasets and the shortcomings of the model for [Fe/H] $< -1.5$ dex become apparent. Accounting for the existence of HNe and/or PISNe (both not included in the MWG-11 rev model shown in Fig. 8) would decrease the predicted [O/Fe] ratio at low metallicities (see Fig. 8 of Romano et al. 2010 and Fig. 11 of Goswami et al. 2021, respectively, for the effects of HNe and PISNe on the model predictions). Moreover, some low-metallicity stars might be not rotating fast (our reference model in Fig. 8 assumes that all massive stars below [Fe/H] $= -1$ are spinstars), which would again bring down the predicted [O/Fe], because of the lower O production expected from non-rotating massive star models (Limongi and Chieffi 2018; Prantzos et al. 2018).

Things get even more complicated when it comes to the minor isotopes. The local enrichment history of the secondary CNO isotopes can be unlocked only by measuring their abundances in unevolved stars. In fact, since the abundances of the minor CNO isotopes are severely altered in the advanced phases of stellar evolution (see Sect. 2), measurements in giant stars do not reflect the original composition of the ISM out of which the stars were born. The problem is that measuring CNO isotopic ratios in dwarf stars is extremely challenging, if not impossible. Spite et al. (2006) determined the $^{12}\hbox {C}$/$^{13}\hbox {C}$ ratio of a sample of 35 extremely metal-poor stars on the lower RGB and found it to be very close to 30 in a subsample of stars that likely suffered minimal mixing processes. Chiappini et al. (2008) adopted a GCE model especially deviced to the Galactic halo including gas infall and outflow and demonstrated that, if fast rotators are common in the early universe, $^{12}\hbox {C}$/$^{13}\hbox {C}$ ratios between 30 and 300 are obtained for [Fe/H] $< -3$ dex. This result implies that some stellar mixing is needed to explain Spite et al. ’s (2006) data, but by far less extreme than what would be necessary if low-metallicity massive stars were not rotating at all. Stellar rotation, in fact, triggers the production of primary $^{13}\hbox {C}$ in the first generations of massive stars, significantly lowering the C isotopic ratio in the ISM of the early Galaxy (see Chiappini et al. 2008, their Fig. 1). Measurements of $^{12}\hbox {C}$/$^{13}\hbox {C}$ ratios are available also for carbon-enhanced metal-poor (CEMP) stars at the turnoff (e.g., Sivarani et al. 2006), but in CEMPs these ratios may mirror the chemical composition of material transferred from AGB companions (e.g., Aoki et al. 2006), or localized enrichment from faint SNe or spinstars (e.g., Hansen et al. 2016). The first and still unique measurement of the carbon isotopic ratio in an ancient, metal-poor, C-normal star that just passed the turnoff rests on the analysis of exquisite spectra obtained with ESPaDOnS at the Canada-France-Hawaii Telescope (CFHT), ESPRESSO at the VLT, and HARPS at the ESO 3.6 m telescope and returns $^{12}\hbox {C}$/$^{13}\hbox {C}$ = $33^{+12}_{-6}$ (Spite et al. 2021), very close to previous estimates for unmixed giants (Spite et al. 2006).

Fig. 9

Expected $^{12}\hbox {C}$/$^{13}\hbox {C}$ ratio as a function of age (top panel) and metallicity (middle panel) according to the GCE models for the solar vicinity of Romano et al. (2021). The models implement nova nucleosynthesis assuming that the nova WD progenitors are either in the range 1–8 $M_{\odot }$ (solid black lines) or in the range 3–8 $M_{\odot }$ (solid blue lines). The evolution of the $^{16}\hbox {O}$/$^{18}\hbox {O}$ ratio is also shown (bottom panel). The entry model is a parallel model. Here we show only the evolution of the thin-disc component. Data are from Spite et al. (2021, HD 140283), Crossfield et al. (2019, GJ 745 AB) Botelho et al. (2020, 63 solar twins, colour-coded from light-yellow to dark-orange going from metal-poorer/older to metal-richer/younger objects), Ayres et al. (2013, the Sun), and Zhang et al. (2021, 2M0355, light blue cross)

Useful constraints to the evolution of the C isotopic ratio in the local disc are set by Botelho et al. (2020), who estimate $^{12}\hbox {C}$/$^{13}\hbox {C}$ ratios from CH A-X and CN B-X features in HARPS spectra of 63 solar twins spanning a wide range of isochrone ages (0.5–9 Gyr). Crossfield et al. (2019) add further information from their detailed analysis of a system of two low-mass M dwarfs orbiting in a wide binary observed at high resolution ($R =$ 70,000) with the iSHELL spectrograph (Rayner et al. 2016) on the NASA Infrared Telescope Facility (IRTF). The stars are fully convective, hence chemically homogeneous, hence reflective of the chemical composition of the ISM at the time of their formation. Spectroscopy of carbon monoxide (CO) at the 4.7 $\mu$m fundamental and 2.3 $\mu$m first-overtone rovibrational bandheads reveals multiple rare isotopologues. For the first component, GJ 745 A, $^{12}\hbox {C}$/$^{13}\hbox {C}$ = $296 \pm 45$ and $^{16}\hbox {O}$/$^{18}\hbox {O}$ = $1220 \pm 260$, while for the second component, GJ 745 B, $^{12}\hbox {C}$/$^{13}\hbox {C}$ = $224 \pm 26$ and $^{16}\hbox {O}$/$^{18}\hbox {O}$ = $1550 \pm 360$. While observations of the CO fundamental rovibrational band at high resolution present several advantages over previous medium-resolution observations (see discussion in Crossfield et al. 2019), the huge amount of observing time per target and the scarceness of spectrographs capable of performing high-resolution, M-band observations make the collection of a statistical sample an unattainable task for now. Last but not the least, we mention the determination of the C isotopic ratio in a young (age $\sim$ 125 Myr) brown dwarf: Zhang et al. (2021) have obtained $^{12}\hbox {CO}/^{13}\hbox {CO}$ = $97^{+25}_{-18}$ for the nearby ($d = 9.1 \pm 0.1$ pc) L dwarf 2MASS J03552337+1133437 (hereinafter 2M0355) from an analysis of archival K-band (2.03$-$2.38 $\mu \hbox {m}$) spectra taken with the spectrograph NIRSPEC at the Keck II 10 m telescope.

In Fig. 9, we compare the stellar measurements discussed in previous paragraphs to the predictions of recent GCE models including nova nucleosynthesis (Romano et al. 2021, and Romano et al., in prep.). Nova systems are expected to inject non-negligible amounts of $^{13}\hbox {C}$, $^{15}\hbox {N}$, and $^{17}\hbox {O}$ in the ISM, with yields varying widely in dependence on the assumptions underlying any specific model of the outburst (see Sect. 2.2 and, especially, Fig. 5). The current nova rate in the Galaxy is also largely uncertain. Recent studies report rates of $43.7^{+19.5}_{-8.7}$ $\hbox {yr}^{-1}$ (De et al. 2021), $26 \pm 5$ $\hbox {yr}^{-1}$ (Kawash et al. 2022), and $28^{+5}_{-4}$ $\hbox {yr}^{-1}$ (Rector et al. 2022). The past Galactic nova rate is basically unknown, though some population synthesis studies suggest the existence of an anticorrelation between metallicity and the number of novae produced: for instance, according to Kemp et al. (2022), the number of novae at $Z = 0.03$ is roughly half that at $Z = 10^{-4}$. Clearly, there is a high degeneracy between the nova rate and the nova yields assumed by any GCE model. Moreover, the models discussed here only accommodate average yields (see, again, Sect. 2.2). From all the above, it is clear that we are treading on thin ice: extreme caution is required when discussing the impact of novae on the Galactic enrichment.

This notwithstanding, it is common wisdom that a large fraction of the meteoritic lithium ($^7\hbox {Li}$) abundance comes from low-mass stars, with novae being the favored candidates (see Romano et al. 1999, 2021, and discussions therein). Resting on the analysis of a sample of the open cluster (OC) and field stars with $^7\hbox {Li}$ abundance measurements from GES iDR6, Romano et al. (2021) have suggested that a better match between theoretical predictions about $^7$Li evolution and observations is possible if the lower mass limit of primary stars entering the formation of nova systems is increased from $\sim$1 $M_\odot$ (a standard assumption of GCE models including novae; see D’Antona and Matteucci 1991; Romano and Matteucci 2003) to 2–3 $M_\odot$. This is to say that systems hosting small WDs are not effective in developing the outbursts (see also Kovetz and Prialnik 1985). Interestingly, while the standard assumption leads to a continuous, mild increase in the predicted nova rate in time, increasing the mass threshold for WD progenitors that potentially lead to nova outbursts results in a nova rate that peaks at early times and declines afterwards, in agreement with more general findings that the fraction of close binary systems decreases with increasing [Fe/H] (e.g., Gao et al. 2014). In Romano et al.’s (2021) model, this happens because the more massive the allowed WD progenitors, the more closely the theoretical nova rate follows the star formation rate (which is overall decreasing in time in the Milky Way), due to the shorter time lag introduced by the lifetimes of more massive stars (although the peak in the nova rate is always delayed relative to that of the star formation rate, because of the time required for the WDs to cool at a level that ensures strong enough nova outbursts; see Romano et al. 2021, and references therein).

A model assuming a mass range for WD progenitors in nova systems intermediate between the two extreme cases discussed here would agree better with the data in the $^{12}\hbox {C}$/$^{13}\hbox {C}$ versus age diagram (see Fig. 9, top panel). However, the increase of the C isotopic ratio with increasing metallicity spotted out in Botelho et al.’s (2020) sample (Fig. 9, middle panel) remains difficult to decipher. A viable explanation could be recent dilution of the local ISM by poorly-processed gas, which would lower the metallicity whilst the $^{12}\hbox {C}$/$^{13}\hbox {C}$ ratio would remain almost constant or slightly decrease in time. The possibility that some of the solar twins are intruders, i.e., stars that were born at inner radii and migrated to their current positions, should be investigated as well. The chemical properties of any intruders, in fact, should not be compared to the predictions of models for the solar neighbourhood (e.g., Guiglion et al. 2019; Romano et al. 2021). Also shown in Fig. 9 is the run of the $^{16}\hbox {O}/^{18}\hbox {O}$ ratio with [Fe/H] (bottom panel). Both $^{16}\hbox {O}$ and $^{18}\hbox {O}$ are produced mostly by massive stars, without any significant contribution from novae, so the two models described above overlap in the O isotope plot. The $^{16}\hbox {O}/^{18}\hbox {O}$ ratio first decreases because of primary $^{18}\hbox {O}$ production triggered by stellar rotation (Limongi and Chieffi 2018). Then, the trend is reversed when the yields from non-rotating stars are adopted. Finally, at higher metallicities secondary $^{18}\hbox {O}$ production becomes dominant and the isotopic ratio decreases again. The model predictions are marginally consistent with Crossfield et al.’s (2019) data as well as with the solar value.

Planet hosts Modeling planet formation and evolution requires detailed information about the chemical composition of the host stars (e.g., Santos et al. 2017, and references therein). In particular, it is important to ascertain the role of CNO elements in the planet formation process, with special regard to C and O, which enter the composition of complex organic molecules that are the precursors of biogenic compounds (Caselli and Ceccarelli 2012). A thorough recap of this and related issues would require a dedicated review and is beyond the scope of the current discussion. However, for the sake of completeness we report below a few facts about the determination of CNO abundances in planet-harboring solar-type stars.

In their pioneering work, Meléndez et al. (2009) used high-resolution abundances for a dozen solar twins and the Sun itself and highlighted that our star has a low ratio between refractory and volatile elements²⁹ relative to the majority of the solar twins. They propose the sequestration of refractory elements in terrestrial planets as a possible explanation, further suggesting that the planets may imprint the signature of their existence in the [X/Fe] versus elemental condensation temperature relation. However, the refractory-to-volatile element ratios may be related to the stellar ages and birth radii: Adibekyan et al. (2014) point out that old stars and stars with an inner disc origin display lower refractory-to-volatile element ratios.

The tendency of stars with giant planets to have a metallicity excess with respect to stars without planets was first spotted out by Gonzalez (1997) and has been confimed by many subsequent studies (see, e.g., the review by Gonzalez 2006). However, this correlation seems to disappear if the planets have masses lower (higher) than 0.05 (4) times the mass of Jupiter (e.g., Buchhave et al. 2012; Santos et al. 2017). As regards the CNO element abundances, we stress again that detailed abundances of CNO elements are difficult to determine in stars. Delgado Mena et al. (2021) find tentative evidence from high-resolution, high-quality spectra that stars hosting low-mass planets have [C/Fe] ratios higher than their counterparts without planets for [Fe/H] $< -0.2$. For higher metallicities, there is no C-enhancement associated with the presence of planets. Oxygen seems to follow a similar trend, but given the larger errors affecting O abundance determinations, no strong conclusion can be drawn (see also Nissen et al. 2014). Homogeneous nitrogen determinations for 74 solar-type stars by Suárez-Andrés et al. (2016) show that planet hosts (42 stars, 57% of the sample) are N-rich with respect to single stars. A likely explanation resides in the fact that their planet hosts are metal-rich and, according to GCE models, higher [N/Fe] ratios are expected at higher [Fe/H].

We conclude that a good understanding of the evolution of the CNO elements in different Galactic components is mandatory to improve our knowledge of exoplanet formation and evolution. High-accuracy, high-precision abundances of statistically significant control samples of stars and increasingly refined GCE models are necessary tools to reach this goal.

Fig. 10

Carbon-to-iron (top panel) and oxygen-to-iron (bottom panel) ratios as functions of [Fe/H] for microlensed unevolved stars in the bulge (Bensby et al. 2021, circles). Red, orange, yellow, green, and blue colours indicate stars from older to younger. Predictions from models for the solar neighbourhood (solid black curves) and for the Galactic bulge (solid and dotted gold curves) from Romano et al. (2020) are also shown. The solar neighbourhood model assumes a gwIMF slope $x = 1.7$ for massive stars, while the bulge models assume the Salpeter slope, $x = 1.35$. The solar neighbourhood model assumes the yields by Limongi and Chieffi (2018) for fast-rotating massive stars below [Fe/H] $= -1$, while above this threshold the ones for non-rotating stars are used. The models for the bulge are computed assuming either that all massive stars are fast rotators (solid curve), or that all massive stars do not rotate (dotted curve), independently of metallicity

Fig. 11

Same as Fig. 10, but for [C/O] versus [O/H]

The Galactic bulge Bensby et al. (2021) have determined the abundances of C and O of 70 microlensed dwarf, turnoff, and subgiant stars found in the Galactic bulge today, thereby providing the first statistically significant sample of bulge stars with C abundance measurements free from the effects of stellar evolution. To the best of our knowledge, before their study C abundances of the unevolved bulge stars were available only for three objects (Johnson et al. 2008; Cohen et al. 2009), making it hard to sensibly constrain carbon enrichment scenarios for the inner Galaxy (Romano et al. 2020, but see Cescutti et al. 2009). To measure carbon, Bensby et al. (2021) rely on the spectral synthesis of six C I lines around 910 nm, while the O I triplet at 777 nm is used to derive oxygen abundances. The abundances are corrected for non-LTE effects.

The trends of [C/Fe] and [O/Fe] as functions of [Fe/H] obtained by Bensby et al. (2021) are similar to the corresponding trends observed in the local thin and thick discs, with bulge stars older than 8 Gyr behaving as thick-disc members and stars younger than 8 Gyr following the thin-disc trend (Bensby et al. 2021). In Fig. 10, Bensby et al.’s data are compared to the predictions of model MWG-11 rev for the solar vicinity³⁰ and to the outputs of two models tailored to the bulge from Romano et al. (2020). In the high-mass domain, a gwIMF flatter than the one used for the solar neighbourhood is assumed for the bulge, to fit better the stellar metallicity distribution function and the run of [$\alpha /\hbox {Fe}$] with [Fe/H] observed for bulge stars (Matteucci et al. 2019, and references therein). The models for the bulge assume two extreme prescriptions about stellar rotation—massive stars are all rotating fast in one case, or do not rotate at all in the other, independently of metallicity. The new data seem consistent with a scenario in which the bulk of the bulge stars owes its origin in the secular evolution of the disc component and with the existence of a metallicity threshold for stellar rotation. The latter hypothesis allows to explain also the N abundances of a sample of giant bulge stars selected from Data Release 16 (DR16) of the APOGEE-2 survey (Kisku et al. 2021), after some corrections for stellar evolutionary effects are applied (see Grisoni et al. 2021).

The [C/O] ratio in the bulge is found to steeply increase with increasing [O/H] (Fig. 11). This contrasts with the almost flat behaviour of [C/O] at super-solar metallicities found by Delgado Mena et al. (2021) for nearby stars. The different abundance diagnostics used in the two studies may be at the origin of the discrepancy.

As a final word of caution, it is worth reminding that the inner Milky Way region is very complex. It hosts different stellar populations that coexist close to the Galactic midplane (Queiroz et al. 2021): an inner thin- and thick-disc component, a bar, a population of counter-rotating stars, and a spheroidal component. The latter experienced an early rapid, vigorous star formation, in agreement with the predictions of classical GCE models (Matteucci et al. 2019). When modeling the inner Galactic regions, the coexistence of these distinct components has to be kept in mind (see, e.g., Haywood et al. 2018).

Fig. 12

Image reproduced with permission from Arellano-Córdova et al. (2021), copyright by the authors

Gradients of oxygen (top panel), nitrogen (middle panel) and N/O (bottom panel) at the present time from H II regions. Orange squares are results from optical emission lines (Arellano-Córdova et al. 2020, 2021). Burgundy circles are results from FIR observations (Rudolph et al. 2006). Solid and dotted lines are the fits to the two datasets. The Sun symbols show the solar photospheric values after Lodders (2020).

The Galactic gradient Nowadays, samples of stars (both in the field and in OCs) from large surveys spanning wide ranges of distance and metallicity can be used to undisclose the full chemical enrichment histories of different Galactic components, revealing the details of the Milky Way assembly process over an unprecedented Galactic volume (e.g., Hayden et al. 2015; Xiang et al. 2019; Gaia Collaboration et al. 2022; Randich et al. 2022). Additional information on the dynamical processes that regulate the distribution of chemical species along the disc comes from the abundance gradients.

Recent observations of H II regions, classical Cepheids, B-type stars, and OCs provide sound estimates of present-day gradients of CNO elements along the disc. Arellano-Córdova et al. (2020, 2021) put forth a homogeneous analysis of deep spectra of 42 H II regions with galactocentric distances $R_{\mathrm {GC}} \simeq$ 4–17 kpc. All the nebulae have direct determinations of the electron temperatures, $T_{\mathrm {e}}$, and revised distances based on Gaia parallaxes of associated stars. Moreover, the most appropriate ionization correction factor (ICF) scheme is chosen on an element-by-element basis. Resting on their careful analysis, Arellano-Córdova et al. (2020, 2021) conclude that the present-day ISM is fairly well mixed, at least in the disc quadrant covered by their observations. Figure 12 illustrates the improvement obtained with respect to previous investigations. The new determinations are consistent with linear radial abundance gradients for oxygen and nitrogen with slopes of $-0.042 \pm 0.009$ dex $\hbox {kpc}^{-1}$ and $-0.057 \pm 0.011$ dex $\hbox {kpc}^{-1}$, respectively (Arellano-Córdova et al. 2021). A flattening of the O gradient in the inner disc (e.g., Esteban and García-Rojas 2018) is ruled out. The N/O ratio gradient is flatter, with a small negative slope of $-0.015 \pm 0.007$ dex $\hbox {kpc}^{-1}$ (Arellano-Córdova et al. 2021). For C/H and C/O, Arellano-Córdova et al. (2020) report slopes of $-0.072 \pm 0.018$ dex $\hbox {kpc}^{-1}$ and $-0.037 \pm 0.016$ dex $\hbox {kpc}^{-1}$, respectively.³¹

Bragança et al. (2019) find a steeper gradient of $-0.07$ dex $\hbox {kpc}^{-1}$ for oxygen in the outer disc ($8.4 \le R_{\mathrm {GC}} \le 15.6$) from a sample of 28 young OB stars with high-resolution spectra obtained with the spectrograph MIKE on the Magellan Clay 6.5-m telescope in Las Campanas. Stanghellini and Haywood (2018) use planetary nebulae (PNe) as tracers of the radial O gradient. By separating their sample in PNe with progenitors younger (older) than 1 Gyr (7.5 Gyr), they find that young objects define an O gradient, $-0.015$ dex $\hbox {kpc}^{-1}$, shallower than that derived from older objects, $-0.027$ dex $\hbox {kpc}^{-1}$. Magrini et al. (2018) present N and O abundances in a sample of 17 GES OCs spanning wide ranges of galactocentric distances (4.5–15 kpc) and ages (0.1–7 Gyr), while Donor et al. (2020) discuss the O abundance gradient based on a high-quality sample of 71 OCs with IR data from APOGEE DR16; both find general agreement with data published in the literature. Finally, the gradient from Cepheids that is derived by Luck (2018) based upon the analysis of 1127 spectra for 435 stars and adopting Bayesian distance estimates based on Gaia data from Bailer-Jones et al. (2018) has slopes of $-0.0665 \pm 0.0038$ dex $\hbox {kpc}^{-1}$, $-0.0470 \pm 0.0040$ dex $\hbox {kpc}^{-1}$, and  $-0.0429 \pm 0.0023$ dex $\hbox {kpc}^{-1}$ for C, N, and O, respectively. It is particularly reassuring that updated gradients from gaseous and stellar diagnostics resting on homogeneous analyses of large samples of objects agree very well with one another within the errors (see Table 1).

Table 1 Slopes (in dex $\hbox {kpc}^{-1}$) of current X/H abundance gradients (X = C, N, O) from recent observations of H II regions and Cepheids. Theoretical predictions from selected GCE models are listed in the last four columns for comparison

Fig. 13

Image adapted from Palla et al. (2020), copyright by the authors

Oxygen (top panels) and nitrogen (bottom panels) gradients predicted by GCE models. The models on the right adopt a constant star formation efficiency and different prescriptions about the intensity of radial gas flows along the disc. The models on the left adopt a star formation efficiency that decreases along the disc without (blue lines) or with (magenta and cyan lines) radial gas flows. The black curves represent the reference model with constant star formation efficiency and no gas flows. Different symbols denote different datasets for H II regions, Cepheids, young OCs, and PNe (see Palla et al. 2020, for references).

The effects of different model assumptions on the predicted current radial abundance gradients are nicely summarized in Palla et al. (2020). We report their results for O and N in Fig. 13 and Table 1. Besides an inside-out formation of the thin disc (Larson 1976; Matteucci and François 1989), a variable star formation efficiency (see Boissier and Prantzos 1999, their Sect. 2.1.3) and radial gas flows (Lacey and Fall 1985; Portinari and Chiosi 2000; Spitoni and Matteucci 2011) must be considered to reproduce the observed gradient. In the early work of Lacey and Fall (1985), radial gas flows arise as a consequence of the viscosity of the gas, the angular momentum difference between the gas accreted by infall and that already present in the disc, and the presence of non-axisymmetric density patterns that modify the angular momentum of the gas. Nowadays, high-resolution, magnetohydrodynamical cosmological simulations of the formation of Milky Way-like galaxies can be used to study the radial transport of material within the disc planes (e.g., Okalidis et al. 2021), which provides scaling laws that can be usefully employed to parametrize radial gas flows between concentric rings in GCE models. The tight relations between abundance and position along the disc obtained through improved data analyses (e.g., Luck 2018; Arellano-Córdova et al. 2020, 2021, see Table 1) have the potential to significantly constrain the models and, hence, the relative importance of the physical mechanisms underlying the emergence of Galactic gradients. Finally, PN and OC samples can be separated in subgroups of different ages, thus enabling the study of the temporal evolution of the gradient. This is especially important since it allows discriminating between models that predict a steepening (e.g., Chiappini et al. 2001) or a flattening of the gradients in time (e.g., Prantzos and Boissier 2000; Vincenzo and Kobayashi 2018a).

Fig. 14

Image adapted from Romano et al. (2019b), copyright by the authors

Clockwise from top left: present-time gradients of C, N, and O isotope ratios across the Milky Way disc. The solid lines are the predictions of model MWG-11 of Romano et al. (2019b), including primary production of $^{13}\hbox {C}$, $^{14}\hbox {N}$, $^{17}\hbox {O}$, and $^{18}\hbox {O}$ from massive fast rotators at low metallicity and $^{13}\hbox {C}$, $^{15}\hbox {N}$, and $^{17}\hbox {O}$ production from novae. Data are from Wilson and Rood (1994, filled triangles), Milam et al. (2005, empty triangles) Polehampton et al. (2005, large empty circles), Wouterloot et al. (2008, small circles), Li et al. (2016, diamonds), Colzi et al. (2018, small pentagons), Boogert et al. (2000, empty squares), and a preliminary analysis of proprietary data (Romano et al. 2019b, large pentagons).

Derivation of CNO isotopic ratio gradients along the disc involves measurements of $^{12}\hbox {C}/^{13}\hbox {C}$, $^{14}\hbox {N}/^{15}\hbox {N}$, $^{16}\hbox {O}/^{18}\hbox {O}$, and $^{18}\hbox {O}/^{17}\hbox {O}$ from various CNO-bearing molecules, such as CO, CN, $\hbox {CH}^+$, $\hbox {H}_2\hbox {CO}$, HCN, HNC, OH, and their isotopologues. A summary of observations can be found in Romano et al. (2017, 2019b). Nowadays, the increased sensitivity and high spectral resolution of millimeter telescopes allow to reach far outer Galactic disc regions and probe their low-metallicity ISM. This provides unique information on the synthesis of the minor CNO isotopes in metal-poor environments (see Romano et al. 2017, 2019b; Zhang et al. 2020; Fontani et al. 2022). Spectroscopy of unevolved stars, in fact, provides only fragments of the local enrichment histories of $^{12}\hbox {C}/^{13}\hbox {C}$ and $^{16}\hbox {O}/^{18}\hbox {O}$ (see Fig. 9), with basically no hope to get any information about $^{14}\hbox {N}/^{15}\hbox {N}$. Observations of the $^{17}\hbox {O}$/$^{18}\hbox {O}$ ratio in stars are limited to a few giants, where the ratio is known to be altered by stellar evolutionary processes (e.g., De Nutte et al. 2017). Therefore, molecular data are a perfect complement to the scarce stellar data. It is necessary to bear in mind, though, that switching from a given isotopologue abundance ratio to its corresponding isotope abundance ratio is potentially complicated by local chemical fractionation effects (e.g., Roueff et al. 2015), like isotope-selective photodissociation (e.g., Casoli et al. 1992). Neglecting local fractionation processes might result, for instance, in highly scattered diagrams for $^{14}\hbox {N}/^{15}\hbox {N}$ (Colzi et al. 2019; Bergner et al. 2020; Evans et al. 2022).

Fig. 14 shows the radial gradients of several CNO isotopic ratios predicted by a specific GCE model (model MWG-11 of Romano et al. 2019b) compared to the observations. The current gradients of $^{12}\hbox {C}/^{13}\hbox {C}$ and $^{16}\hbox {O}/^{18}\hbox {O}$ (Fig. 14, left-handed panels) are in satisfactory agreement with the observed trends, although the latter are characterized by a significant dispersion. As for nitrogen, the model predicts an increase of the isotopic ratio when moving from the inner Galaxy to the solar position and a decrease afterwards, when moving from solar to outer radii (Fig. 14, top right-handed panel). Such a prediction still waits to be proven, or disproven, by the observations. Finally, the theoretical $^{18}\hbox {O}/^{17}\hbox {O}$ gradient presents a moderate variation from a value of approximately 3 at $R_{\mathrm {GC}} =$ 4 kpc to $^{18}\hbox {O}/^{17}\hbox {O}$ $< 5$ at $R_{\mathrm {GC}} =$ 16 kpc (Fig. 14, bottom right-handed panel), which fits only in part the data.

We close this section with a word of caution regarding the model predictions discussed here for the inner Galaxy: the adopted yields for massive stars do not extend to super-solar metallicities, meaning that arbitrary extrapolations are needed when the computed metallicity exceeds solar. It is hard to foresee how this may affect the GCE model predictions, but secondary elements such as $^{18}\hbox {O}$ are likely to be the most affected.

3.1.2 Dwarf spheroidal and ultrafaint dwarf galaxies

The local universe contains a large number of dwarf galaxies, the faintest of which are still being found lurking in and around the haloes of the Milky Way and Andromeda, after the first discovery of an UFD by Willman et al. (2005). The properties of stellar populations in local dwarf galaxies are summarized in excellent reviews by Tolstoy et al. (2009); Simon (2019); Annibali and Tosi (2022), to which we address the interested reader.

Being concerned here with the evolution of the CNO elements, it comes natural to divide the dwarf galaxies into two loose groups, the one comprising galaxies that still have large amounts of gas and active star formation, the other including galaxies with no detectable (or only traces of) gas at present. The first is the subject of the next section, while the second is dealt with in the following paragraphs.

Many of the closest dwarf galaxies host old stellar populations and have lost their gas, hence their ability to form stars. Being placed relatively nearby, their stellar populations can be studied on a star-by-star basis with current observational facilities. Their low mean metallicity gives us the possibility to investigate how star formation and chemical enrichment proceeded in the early universe (e.g., Frebel and Norris 2015). High-resolution abundance measurements, however, are possible only for giant stars. A fair comparison between the predictions of chemical evolution models and the observations thus implies uncertain corrections to the derived abundances of C and N, which hampers the achievement of firm conclusions, at least as far as these elements are concerned (e.g., Lardo et al. 2016, and references therein).

Fig. 15

Image reproduced with permission from Jeon et al. (2021), copyright by the authors

A(C) and [C/Fe] versus [Fe/H] diagrams (left-hand and right-hand panels, respectively) of individual stars observed in UFDs (coloured crosses, top panels) and in dSphs (small violet $\times$ symbols, bottom panels), compared to simulated UFD stars (Jeon et al. 2021) assuming contributions from normal (cyan circles) and faint SNe (grey circles). Simulated dSph stars from the work of Jeon et al. (2017) are shown as pink squares.

It is, however, the existence of a non-negligible fraction of CEMP stars in these galaxies that has attracted major attention, both in an attempt to understand the origin of their anomalously high C abundances and in connection to the formation of the Galactic halo: if the dwarf galaxies orbiting the Milky Way today are the remnants of the hierarchical merging process that built up the Galactic halo many Gyr ago, the halo and the surviving satellites should contain similar fractions of CEMPs.

Metal-poor stars having [C/Fe] $> +1$ (CEMPs)³² are classified in subcategories taking the relative abundances of heavy neutron-capture elements into account: CEMP-s and CEMP-r are CEMPs characterized by s- and r-process element enhancement, respectively, CEMP-r/s present both kinds of enhancements, and CEMP-no no enhancement at all (Beers and Christlieb 2005). For [Fe/H] $< -1$, many stars in the Galactic halo and in dwarf satellites are CEMPs, with relative fractions with respect to C-normal stars increasing with decreasing the metallicity (e.g., Aoki et al. 2007; Yong et al. 2013; Salvadori et al. 2015), reaching 70% below [Fe/H] $= -4$ (Yoon et al. 2018). For [Fe/H] $> -3$, CEMP-s stars are more numerous; their peculiar chemical composition is likely due to mass transfer from a close AGB companion (e.g., Aoki et al. 2006). On the other hand, CEMP-no stars with [Fe/H] $< -2$ might either have a nucleosynthetic origin and reflect the properties of the elusive first (Pop III) stars, or be again the result of mass transfer from an AGB star in a binary system (Bonifacio et al. 2015). With their low-metallicity stellar populations, [Fe/H] $< -2$ (Simon 2019), UFDs are the ideal laboratory to study CEMPs and the effects of pristine chemical enrichment from Pop III stars (Salvadori et al. 2019).

In their recent study based on cosmological hydrodynamic zoom-in simulations of isolated UFDs, Jeon et al. (2021) find that most CEMP-no stars in the UFDs and the Milky Way halo can be explained by normal Pop III and Pop II SNe (with explosion energy of the order of $10^{51}$ ergs). Faint SNe (with explosion energy of the order of $0.6 \times 10^{51}$ ergs) provide a viable path to CEMP-no stars with [C/Fe] exceeding a value of 2, though the highest C overabundances remain unexplained. Figure 15 summarizes the results obtained by Jeon et al. (2021) for both UFDs and dSphs. The loci occupied by the observations are, generally, very well recovered by the models, but there are simulated stars that fall in unpopulated regions, and viceversa, densely-populated regions that are missed by the models. The same problem is faced by Komiya et al. (2020), using their stochastic chemical evolution model in cosmological context that implements the yields from two sets of mixing and fallback models for faint SNe. Overall, this calls for a refinement of the nucleosynthesis prescriptions underlying the simulations. Indeed, detailed analyses of high-resolution CEMP spectra providing the abundances of many elements suggest that different kinds of Pop III polluters must be in action (Skúladóttir et al. 2021, see also Bonifacio et al. 2015).

3.1.3 Dwarf irregular and blue compact dwarf galaxies

Local dIrrs and BCDs contain large fractions of gas and are subject to ongoing star formation. This, joint to the low metallicities, made people think of them as convenient, readily accessible counterparts of the basic galactic building blocks predicted numerous in the distant universe (White and Rees 1978). However, the physical conditions within galaxies and their surroundings vary with cosmic time (Kewley et al. 2013; Shapley et al. 2019), which challenges the primeval galaxy interpretation (even leaving aside the effects of internal chemical enrichment). The high spatial resolution of the Hubble Space Telescope (HST) photometry has allowed us to resolve local (up to $\sim 20$ Mpc distance from us) dIrrs and BCDs in stars and to build deep colour-magnitude diagrams (CMDs) from which detailed star formation histories (SFHs) have been extracted (see Tosi et al. 1991; Gallart et al. 2005; Tolstoy et al. 2009; Cignoni and Tosi 2010; Annibali and Tosi 2022, and references therein): local dIrrs and BCDs prove to be definitively old objects and are found to have formed stars as far back in time as the observations can reach. The disparate SFHs revealed via CMD analyses for objects that share similar mass and metallicity suggest that many different physical processes—gas accretion, mergers, interactions with halo hosts and close companions, stellar chemical and energetic feedback—combine in a complex way to lead to the specific configuration observed for each system.

With the notable exceptions of the Large and Small Magellanic Clouds (hereinafter, LMC and SMC, respectively), for which statistically significant samples of (giant) stars with high-resolution spectra are available (Pompéia et al. 2008; Lapenna et al. 2012; Nidever et al. 2020), chemical abundances in dIrrs and BCDs are mainly measured from optical and NIR spectroscopy of bright H II regions (see Annibali and Tosi 2022, and references therein) and UV nebular emission lines (e.g., Berg et al. 2019). Thus, for these galaxies we do not have a complete fossil record of the past enrichment history, but only the current snapshots (see, e.g., Chiappini et al. 2003).

Fig. 16

Image reproduced with permission from Berg et al. (2019), copyright by AAS

Theoretical tracks in a log(C/O) versus log(O/H) + 12 plane (coloured curves) predicted by chemical evolution models assuming a series of star formation bursts of different duration as well as various fractions of oxygen expelled by galactic winds, as given in the legend (Berg et al. 2019). Different symbols are data for Damped Lyman $\alpha$ systems, Milky Way halo and disc stars, and H II regions in dwarf galaxies (references in Berg et al. 2019). The panel to the right illustrates the physical processes responsible for the sawtooth behaviour of the theoretical tracks (see text).

In Fig. 16, taken from Berg et al. (2019), determinations of the C/O ratio in different systems are compared to the predictions of generic chemical evolution models. The sawtooth behaviour (see Pilyugin 1992) seen in the theoretical curves arises from the assumption of a bursty star formation: each burst makes the metallicity increase and the C/O ratio decrease due to the short time delay with which oxygen, produced by short-lived massive stars, is injected into the ISM relative to carbon, which owes a large fraction of its production to low- and intermediate-mass stars acting on longer timescales (Tinsley 1979). Differential galactic outflows, venting out of the system mainly the products of CCSN nucleosynthesis, contrast both the increase of O/H and the decrease of C/O. Nitrogen follows a qualitatively similar path. Interestingly, from their analysis Berg et al. (2019) conclude that the C/N ratio of the low-metallicity sample is constant with metallicity, but has a significant scatter: they suggest that this may hint at the dominant production of C from stars in a different mass range relative to that responsible for most of N. The broad variations in the observations can be reproduced by heuristic chemical evolution models that vary solely the efficiency of star formation, its duration, and the fraction of ejected oxygen (see Berg et al. 2019, their Fig. 12). However, as pointed out by Romano et al. (2006), if a name is given to the generic point in the C/O, N/O, C/N diagrams, and a present-day gas fraction, a stellar mass, and a detailed SFH can be associated to it, it may become difficult to reproduce the observed abundances. This is because a better-constrained GCE model often reveals all the limitations of its simple approach.

While a galactic outflow that expels preferentially metal-rich matter processed by massive stars is a convenient way of reducing the effective yield of a stellar population, the very low star formation rates observed in some metal-poor and extremely metal-poor dwarf galaxies would actually point to a mild outflow effect. In fact, star formation rates lower than $1\,M_{\odot }\, \mathrm {yr}^{-1}$ imply a deficit of massive clusters able to generate massive stars (Yan et al. 2017), hence a less effective SN feedback (both in terms of metal and energy production). Yan et al. (2020) apply the integrated galactic IMF theory that allows to compute the gwIMF as a function of a galaxy star formation rate and metallicity, and find a mildly bottom- and top-light gwIMF for the UFD Boötes I, which leads to a good agreement of their GCE model predictions with the data. We note that a gwIMF in which oxygen production from massive stars is depressed, whilst nitrogen production from intermediate-mass stars is enhanced, could provide a compelling explanation for the high N/O ratios observed in some extremely metal-poor, gas-rich dwarf galaxies, such as Leo P, DDO 68, and Leoncino (AGC 198691) (Skillman et al. 2013; Annibali et al. 2019b; Aver et al. 2022), which lie significantly above the very narrow low-metallicity plateau at $\log \mathrm{(N/O)} = -1.60 \pm 0.02$ suggested by Izotov and Thuan (1999).

Another ingredient with which GCE models can play a bit to reduce the ISM metallicity of a galaxy is the rate of infall of primordial gas. The effect of infall is maximum on O/H but tends to cancel out when considering the abundance ratios (except for a second order effect on the yields, due to their metallicity dependence). In a recent study by Pascale et al. (2022), N-body hydrodynamical simulations are used to investigate whether the extremely low metallicity of DDO 68 (Annibali et al. 2019b) can be explained by infall of gas of primordial chemical composition. The models, tailored to reproduce the main observed kinematic and structural properties of this galaxy (Annibali et al. 2016, 2019a), show that the infalling mass necessary to dilute the oxygen that is produced according to the CMD-inferred SFH and a canonical gwIMF is inconsistent with any possible dynamical solution, unless the measured metallicity refers to the accreted material and not that of the main galaxy.

The roles of all the different physical mechanisms mentioned above in setting the current chemical properties of dwarf galaxies have been recognized and discussed long ago (Matteucci and Chiosi 1983). Nowadays, also thanks to the parallelization of algorithms that boosts the computational power of supercomputers, hydrodynamical simulations can be coupled with chemical evolution models to better explain the observations of gas-rich small galaxies and increase the predictive power of both (e.g., Recchi et al. 2001; Recchi and Hensler 2007; Emerick et al. 2020, among others). Because of the dramatic increase of the necessary computational resources when introducing a new chemical species, however, usually a choice has to be made and only a few elements can be selected. We stress that it would be important to include at least the main CNO elements in the simulations.

Turning to the minor CNO isotopes, abundance estimates of CNO isotopes in dwarf galaxies are available only for the Magellanic Clouds. From a spectral line survey of the prominent star-forming region N 113 in the LMC covering parts of the frequency range from 85 GHz to 357 GHz, ratios of $^{12}\hbox {C}$/$^{13}\hbox {C}$ $\simeq 49 \pm 5$, $^{16}\hbox {O}$/$^{18}\hbox {O}$ $\simeq 2000 \pm 250$, and $^{18}\hbox {O}$/$^{17}\hbox {O}$ $\simeq 1.7 \pm 0.2$ are estimated by Wang et al. (2009). Atacama Large Millimeter/submillimeter Array (ALMA) observations in the $^{12}$CO(J = $2-1$), $^{13}$CO(J = $2-1$), $\hbox {C}^{18}\hbox {O}$(J = $2-1$), $^{12}$CO(J = $3-2$), and $^{13}$CO(J = $3-2$) lines toward the active star-forming region N83C in the SMC are presented by Muraoka et al. (2017), from which the $^{12}\hbox {C}$/$^{13}\hbox {C}$ and $^{16}\hbox {O}$/$^{18}\hbox {O}$ isotopic ratios can be derived. We will comment on this in Sect. 3.2.

3.1.4 The Andromeda galaxy

We conclude this roundup of local galaxies with Andromeda (M 31), a giant spiral galaxy that is also the most massive member of the Local Group. At variance with the disc of the Milky Way, that suffered a last major merger about 10 Gyr ago, leading to the formation of the thick disc (Helmi et al. 2018), and evolved peacefully since, M 31 is representative of a lively history of mergers, with clear signs of recent encounters (McConnachie et al. 2009). Bhattacharya et al. (2022) provide oxygen abundance measurements from direct detection of the [O III] $\lambda$ 4363 Å line for a magnitude-limited sample of 205 PNe in the disc of M 31. The sample is divided into high- and low-extinction PNe that are associated, respectively, with ages $\le 2.5$ Gyr and $\ge 4.5$ Gyr, namely, with thin- and thick-disc populations. The sample as a whole defines a flat gradient, $0.001 \pm 0.003$ dex $\hbox {kpc}^{-1}$, consistent with zero as also found by previous studies based on PNe (e.g., Peña and Flores-Durán 2019). However, it is also seen that older PNe are metal-poorer and trace a mild positive gradient of $0.006 \pm 0.003$ dex $\hbox {kpc}^{-1}$, while younger PNe are metal-richer and define a steeper gradient of $-0.013 \pm 0.006$ dex $\hbox {kpc}^{-1}$. The latter is marginally consistent with the O abundance gradient measured from other young indicators, namely, H II regions ($-0.023 \pm 0.002$ dex $\hbox {kpc}^{-1}$, Zurita and Bresolin 2012).

N-body simulations of galaxy interactions show that mergers tend to flatten the radial metallicity gradient (e.g., Rupke et al. 2010), while classic GCE models emphasize the role of the gas threshold and efficiency of the star formation process in steepening or flattening the gradient (Marcon-Uchida et al. 2010). Ad hoc simulations specifically tailored to reproduce the observed characteristics of Andromeda, like the ones presented by Yin et al. (2009) and Robles-Valdez et al. (2014), are necessary to quantify the relative importance of all the processes at play.

3.1.5 Other systems

Fig. 17

Image reproduced with permission from Vincenzo et al. (2016), copyright by the authors

Heuristic chemical evolution model from (Vincenzo et al. 2016, black solid line) reproducing the main features of local galaxies drawn from the SDSS DR7 in the log(N/O)–log(O/H) + 12 plane. The observed distribution is shown as a 2D histogram, with a bin size along both axes of 0.025 dex. Different colours represent the number of galaxies in each bin. The changes in slope of the predicted relation are linked to different physical properties of the model, as indicated.

For local elliptical galaxies, indications about the abundances of CNO elements come essentially from absorption features in integrated starlight, while emission lines in UV and optical spectra are used for gas-rich star-forming galaxies. A comprehensive review of these topics can be found in Maiolino and Mannucci (2019). In the following, we report a couple of examples focusing on the theoretical implications.

Table 2 $^{12}\hbox {CO}/^{13}\hbox {CO}$, $^{12}\hbox {CO}/\hbox {C}^{18}\hbox {O}$, and $^{13}\hbox {CO}/\hbox {C}^{18}\hbox {O}$ abundance ratios for external galaxies from the literature (adapted from Table 1 of Romano et al. 2017, copyright by the authors)

Conroy et al. (2014) model stacked spectra of early-type galaxies from the SDSS (see also Johansson et al. 2012) as a function of velocity dispersion ($\sigma$). The abundances of C, N, and O measured from the integrated light reveal that [O/Fe] increases from about 0 for $\sigma =$ 90 km $\hbox {s}^{-1}$ to 0.25 for $\sigma =$ 300 km $\hbox {s}^{-1}$. In the same velocity dispersion range, C and N track O and their ratios relative to Fe exceed 0.2 at the high $\sigma$ end. The observed ratios are better explained if the gwIMF of elliptical galaxy progenitors is skewed in favour of massive stars (see Romano et al. 2020, see also Yan et al. 2021). The idea of a flatter IMF slope in the high stellar mass range for elliptical galaxies relative to that assumed for the solar neighbourhood is not new: it has been proposed long ago to fit the [Mg/Fe] abundance ratios of ellipticals derived from Lick indices (Matteucci 1994).

The log(N/O) versus log(O/H) + 12 diagram of star-forming galaxies in the local universe has been studied by Vincenzo et al. (2016). These authors use [O II] $\lambda \lambda$ 3726, 28, [O III] $\lambda$ 5007, H $\beta$, H $\alpha$, [N II] $\lambda$ 6584, and [S II] $\lambda \lambda$ 6717, 31 line fluxes in emission to infer gas-phase O and N abundances of galaxies selected from SDSS Data Release 7 spanning a wide metallicity range. The effects of dust depletion are taken into account. The main characteristics of the diagram are well reproduced by the heuristic GCE model sketched in Fig. 17 that includes: (i) primary N production from fast-rotating massive stars at low metallicities, (ii) a major contribution to N production (mainly of secondary origin) from intermediate-mass stars at higher metallicities, and (iii) the effects of galactic outflows. The outflows are assumed to eject preferentially the products of massive star nucleosynthesis (oxygen), which explains the steepening of the relation relative to the flatter slope obtained at previous times, due solely to the late injection of N from low- and intermediate-mass stars.

Local gas-rich galaxies can also be targeted for isotopic CO transition measurements. While the limits of past observational capabilities severely hampered the investigations, which were restricted to metal-rich galaxies (see Encrenaz et al. 1979), ALMA now allows to reach much fainter objects. In Table 2, that is a slightly modified version of Table 1 of Romano et al. (2017), determinations of $^{12}\hbox {CO}/^{13}\hbox {CO}$, $^{12}\hbox {CO}/\hbox {C}^{18}\hbox {O}$, and $^{13}\hbox {CO}/\hbox {C}^{18}\hbox {O}$ in a number of local and high-redshift objects (see Sect. 3.2.1) are reported. The abundance ratios of isotopes that originate from stars in different mass intervals are seen to vary widely. In particular, in starburst galaxies the ratio of two generic isotopes L and M, where M is produced mainly by massive stars on short timescales and L is produced mainly by low- and intermediate-mass stars on longer timescales (or has a mixed origin) is lower than in secularly evolving objects (Henkel and Mauersberger 1993; Romano et al. 2017; Zhang et al. 2018; Brown and Wilson 2019). We will come back to this in Sect. 3.2.1.

3.2 Moving beyond the local universe

3.2.1 The high-redshift universe

Fig. 18

Images reproduced with permission from Hayden-Pawson et al. (2022), copyright by the authors

Top panel: nitrogen-to-oxygen abundance ratio versus gas-phase metallicity (as traced by oxygen) for the SDSS (grey contours) and KLEVER (magenta points) galaxy samples discussed by Hayden-Pawson et al. (2022). The white points are the median binned SDSS data with corresponding standard deviations. The solid red curve is the best fit to local data. The yellow star refers to the average position of the KLEVER galaxies. The dot-dashed vertical line indicates the metallicity above which secondary N production from low- and intermediate-mass stars becomes dominant. The histogram in the top-left corner shows the deviations of the SDSS sample and the KLEVER sample from the best-fitting line to the local data. Bottom panel: nitrogen abundance as a function of the stellar mass for the SDSS and KLEVER samples. Contours, line and symbols have the same meanings as in the top panel.

Fig. 19

Image reproduced with permission from Zhang et al. (2018), copyright by the authors

Lines: run with time of the theoretical $^{13}\hbox {C}$/$^{18}\hbox {O}$ isotopic ratio in the ISM (bottom panel) for a reference galaxy that forms stars in an early burst lasting 1 Gyr (top panel). The different colours refer to different choices of the gwIMF. The IMF slopes for $m > 1\,M_\odot$ are $x = 1.7$ for the green and yellow curves (changing the prescriptions for the low-mass stars), $x = 1.1$ for the red curve, and $x = 0.95$ for the black one, where $x = 1.35$ is the Salpeter slope. Red symbols indicate the SMGs in Zhang et al.’s (2018) sample. Blue symbols refer to the $^{13}\hbox {CO}/\hbox {C}^{18}\hbox {O}$ ratios measured in a few representative galaxies.

The properties of local star-forming galaxies can be usefully compared to those of their higher-redshift counterparts, providing important clues on the evolution of the CNO elements, in particular the N/O ratio. Hayden-Pawson et al. (2022) compare the trends of log(N/O) as functions of metallicity and stellar mass emerging from a sample of 37 $z \sim$ 2 galaxies drawn from the KMOS Lensed Emission Lines and VElocity Review (KLEVER) Survey to the corresponding patterns of a sample of local star-forming galaxies selected from the SDSS (see Fig. 18). Despite the large scatter in the data (see also Strom et al. 2018), a modest enhancement of $+0.1$ dex is apparent when comparing the KLEVER sample to local galaxies in the log(N/O)–log(O/H) + 12 plot, whereas at a fixed stellar mass a depletion is observed, at the level of $-0.35$ dex. For local galaxies, a strong anticorrelation is found between N/O and star formation rate in the N/O–stellar mass plane, resembling that between O/H and star formation rate in the mass-metallicity relation. Such anticorrelation is used to construct a fundamental N relation analogous to the fundamental metallicity relation. The KLEVER galaxies are consistent with both. Therefore, the depletion of N/O in the KLEVER sample at a fixed stellar mass can be interpreted as an effect of the redshift evolution of the mass-metallicity relation coupled to a redshift-invariant N/O–O/H relation (Hayden-Pawson et al. 2022), a scenario consistent with the models discussed in Vincenzo et al. (2016).

Moving to even higher redshifts, z up to 8–9, the [C II] and [O III] gas emission can be detected in normal galaxies thanks to ALMA superb sensitivity (Capak et al. 2015; Combes 2018; Hodge and da Cunha 2020). Moreover, the most powerful starbursts, thought to be the progenitors of the current population of quenched, passively-evolving, massive early-type galaxies (Toft et al. 2014), can be observed through their enshrouding dust curtains. Zhang et al. (2018) have observed the rotational transitions of the $^{13}\hbox {CO}$ and $\hbox {C}^{18}\hbox {O}$ isotopologues in the cold molecular gas of four submillimeter galaxies (SMGs) at $z \sim$ 2–3 with ALMA and determined their $^{13}\hbox {C}$/$^{18}\hbox {O}$ abundance ratio. Importantly, the determination of this ratio is immune to the effects of dust (Zhang et al. 2018).

Figure 19 shows the predictions of some reference chemical evolution models in comparison to Zhang et al. ’s dataset as well as other data from the literature. The effects of different assumptions about the gwIMF on the predicted $^{13}\hbox {C}$/$^{18}\hbox {O}$ ratio are explored. It is seen that, assuming star formation rates that oscillate between 100 and 1000 $M_\odot$ $\hbox {yr}^{-1}$ for the observed SMGs, the high-redshift measurements are best explained by a gwIMF skewed towards massive stars relative to that assumed to fit the Milky Way properties. This is due to the fact that $^{13}\hbox {C}$ and $^{18}\hbox {O}$ come mainly from intermediate- and high-mass stars, respectively, so a $^{13}\hbox {C}$/$^{18}\hbox {O}$ ratio as low as $\sim$1, as observed, requires a higher fraction of massive stars. Although the requirement of a top-heavy IMF in Zhang et al. (2018) is based on a set of yields that do not take the effects of stellar rotation into account, it has been demonstrated that by adopting the new yields for fast rotators from Limongi and Chieffi (2018) the conclusions are qualitatively the same (Romano et al. 2019b). This is better understood since rotation produces its largest effects at low-metallicity regimes, while SMGs reach soon solar metallicities because of the vigorous star formation regimes they experience.

4 Building the bridge from “here and now” to “there and then”

In the Local Group, reliable CNO element abundances can be derived from high-resolution stellar spectra in the Milky Way and its closest satellites, with recent large surveys performed with multifiber spectrographs having led to major advances in the field (see Randich and Magrini 2021, for a review). High-resolution spectroscopy provides also the necessary anchor to the development of new data-driven techniques that allow the determination of precise chemical abundances from low-resolution spectra (see, e.g., Ting et al. 2018). In this way, low-resolution ($R =$ 1800) LAMOST spectra of millions of stars are fully exploited to recover the abundances of many elements, including the CNO ones (Wheeler et al. 2020, see also Xiang et al. 2019).

Gaseous diagnostics are also routinely analysed to explore disc gradients and to unravel the properties of the ISM in metal-poor, gas-rich dwarf galaxies leading, generally, to a good agreement with stellar indicators (e.g., Bhattacharya et al. 2022; Esteban et al. 2022). Together with estimates of the star and gas masses and distributions, the rates of star formation, the rates of CCSN and type Ia SN explosions, the nova outburst rates, and full SFHs from CMDs, these abundances provide significant constraints to GCE models of our own and other nearby galaxies (Matteucci 2021).

Exploring the universe beyond our backyard add invaluable information about the temporal evolution of fundamental relations, such as the mass-metallicity relation (e.g., Maiolino et al. 2008; Wang et al. 2022), or invariance thereof (e.g., Hayden-Pawson et al. 2022), and populate the galactic zoo with more exotic objects (e.g., Cardamone et al. 2009), allowing us to explore different physical conditions for star formation and chemical enrichment. As we have seen in previous sections, CNO elements are fundamental to this kind of investigation. GCE models for external galaxies can be constructed, but their predictive power is weakened by the shortage of some of the fundamental observational constraints listed above. For high-redshift galaxies a minimal list would include well-determined star formation rates, stellar masses, and ISM metallicities, as well as any information on the outflows. In this respect, we note that molecular outflows, which were largely out of reach till only a few years ago, can now be studied with the most sensitive mm/sub-mm interferometers, such as ALMA and NOEMA, though large integrations are needed to unravel their signatures through the CO line emission (Cicone et al. 2018).

Before building models for external galaxies, moreover, it would be highly desirable to have the adopted stellar nucleosynthesis prescriptions tested against observations for the Milky Way, the best observed among all the local systems. In general, all the uncertainties underlying the model ingredients (which lead to degeneracies among the possibile solutions) must be kept well in mind and, as we have seen in Sect. 2, the stellar yields come with large associated uncertainties.

5 Open issues and future directions

The above discussion leads to a list of desiderata that I translate and split into two lists of theoretical and observational challenges for the years to come in the remaider of this section.

5.1 Theoretical challenges

There are several improvements to the theory that would be highly welcome. In the following, a minimal personal selection of the most crucial items is given.

• The secular evolution of stars is still out of reach of 3D hydrodynamic simulations. At present, grids of yields covering a range of initial masses and metallicities adequate for chemical evolution studies are obtainable only via 1D models that rely on empirical prescriptions and parametric recipes to implement complex physical phenomena, such as convective motions and mass loss. However, small portions of the stars as well as small fractions of their lifetimes can be studied through 3D simulations. HPC resources should be devoted to these studies since they provide fundamental hints to 1D modeling.

• Important physical processes leading to modification of the yields, such as stellar rotation, magnetic fields, and interactions with a companion in close binary systems, are still largely unexplored and deserve further efforts.

• GCE models must include the effects of stellar duplicity on nucleosynthesis. Besides type Ia SNe and novae, which originate from low- and intermediate-mass stars, massive binaries must be implemented. Stellar duplicity is particularly relevant to CNO element evolution, as we discuss in this review.

• Full coupling of GCE models with 3D hydrodynamical simulations is not possible yet. However, it is crucial to implement the evolution of all CNO elements (and, possibly, of their minor isotopes) in the simulations.

• Galactic outflows are critical ingredients of GCE models. High-resolution 3D hydrodynamical simulations are needed to provide information on their effectiveness not only in removing metals but also in entraining the neutral ISM in the flow. This allows to reduce sensibly the free parameter space of GCE models, hence reducing the degeneracy of the proposed solutions.

5.2 Observational challenges

Observations provide the necessary tests to model predictions and a guidance towards the development of more refined/new theories. In the following, a list of observations that may potentially have a huge impact on the implementation of GCE models is given.

• Related to the second item in the previous list, more time should be spent to improve the statistics of the binary fraction, period, and mass ratio distributions. Particular attention should be devoted to high-mass stars. In fact, massive stars often interact with companions that influence their lives at every stage, but more so during the giant phases in which the stars lose copious mass.

• Although considerable progress has been made in recent years, further observational efforts are needed to understand which CNO diagnostics are the best to use in stars in dependence on their physical parameters. It is also important to keep in mind the usage that has to be made of the abundance determinations. For instance, Delgado Mena et al. (2021) point out that, since the O I $\lambda$ 6158 Å line has an high excitation potential, as the C I lines, in relatively warm dwarf stars in the disc, it is better to use O I $\lambda$ 6158 Å lines to evaluate the C/O ratio in those stars. Since using different indicators may give different trends, the information and physical justification provided by Delgado Mena et al. (2021) is very useful to GCE modelers.

• The N I atomic line at 7468 Å is very weak also in solar-type stars, which forces to study the very crowded molecular band at 3360 Å, where continuum placement is challenging. High resolution, of the order of $R = 80,000$, and SNR = 150 are necessary for useful NH analysis in stars. This should be considered in building the next generation of instruments that will probe the UV band (e.g., CUBES, HRMOS).

• Better characterization of the external regions of the Galactic disc is needed to improve our understanding of how CNO element enrichment proceeds at low metallicities. The project CHEMOUT (Fontani et al. 2022), now in its infancy, will address this point. Other surveys dedicated to the outer Galaxy, exploring different quadrants, would be highly welcome.

• PISNe and faint SNe seem necessary to improve the predictions of GCE model on CNO evolution at low metallicities: the James Webb Space Telescope (JWST) promises to provide important clues on the nature of these objects at high redshifts (e.g., Hartwig et al. 2018).