© The Authors 2025

1. Introduction

The standard big bang nucleosynthesis (BBN) framework connects cosmology with the standard model of particle physics, describing the formation of light elements (D, ³He, ⁴He, ⁷Li, and ⁷Be) during the first minutes after the big bang. The only free parameter in standard BBN is the baryon-to-photon ratio (η), which is now best determined through measurements of the cosmic microwave background (CMB) radiation by instruments aboard the Wilkinson Microwave Anisotropy Probe (WMAP) and Planck satellites (Spergel et al. 2003; Planck Collaboration VI 2020; Planck Collaboration VI 2021). In particular, CMB measurements of baryon density serve as inputs for primordial nucleosynthesis calculations, enabling predictions of elemental abundances that can be directly compared with observational data (see, e.g., Fields 2011 and the references therein).

The combination of BBN and CMB predictions shows excellent agreement with high-redshift observational data for D and ⁴He. However, the same is not true for ⁷Li: BBN+CMB predictions suggest a primordial abundance of (⁷Li/H)_{BBN + CMB} = (4.94 ± 0.72)×10⁻¹⁰ (Yeh et al. 2021; Fields & Olive 2022), whereas metal-poor near-main-sequence stars in the Spite plateau indicate an observed value¹ of (⁷Li/H)_obs = (1.6 ± 0.3)×10⁻¹⁰ (Fields & Olive 2022 and references therein). This discrepancy, with the theoretical value exceeding the observed one by a factor of ∼3 − 4, is known as the cosmological lithiumproblem.

Many studies have attempted to address the lithium problem, proposing various hypotheses to explain the discrepancy between theoretical predictions and observations. One straightforward explanation is that primordial lithium was depleted through stellar processes post-formation. Stellar models that incorporate diffusion and lithium transport mechanisms have been developed to quantify lithium abundances in old stars, partially explaining the observed levels (see, e.g., Korn et al. 2007; Gruyters et al. 2016; Deal et al. 2021). Additional processes such as rotational mixing, turbulence, and magnetic field effects can further reduce surface lithium over time, particularly in low-mass stars, helping explain Spite plateau lithium levels despite variations in metallicity and stellar temperature among Population (Pop) II stars (Borisov et al. 2024). In particular, a compilation of key studies on lithium abundances in metal-poor near-main-sequence stars over the past few decades can be found in Norris et al. (2023).

An alternative explanation posits that primordial nucleosynthesis predictions require modification due to physics beyond the standard model. The existence of exotic particles or modifications in baryon–photon interactions in the early Universe could have altered lithium production during primordial nucleosynthesis (Fields 2011, see also Perivolaropoulos & Skara 2022 for a recent review).

Similarly, supersymmetry theories introduce particles that could have influenced lithium evolution by acting as catalysts for nuclear reactions or interfering with beryllium decay, thereby affecting lithium production (Cyburt et al. 2009). Dark matter (DM) interactions have also been explored, with hypotheses suggesting that DM particles decayed during or shortly after primordial nucleosynthesis, releasing energy and modifying nuclear reaction rates relevant to lithium abundances (Pospelov & Pradler 2010). Another possibility, motivated by models involving dynamic scalar fields, such as Guts or dilaton theories, suggests that the fine-structure constant was larger than the current value by 2 − 20 ppm during the BBN period (see Clara & Martins 2020; Martins 2021; Deal & Martins 2021). While these alternative scenarios are theoretically well motivated, they introduce additional complexities when confronted with observational data.

The Spite plateau was first identified in low-metallicity Galactic stars (Spite & Spite 1982), which are considered good indicators of the primordial lithium abundance due to their minimal contamination from later nuclear processes. The Sagittarius (Sgr) dwarf galaxy was among the first extragalactic systems studied for Spite plateau lithium (Mucciarelli et al. 2014a,b). Low-metallicity stars in Sgr exhibit lithium abundances consistent with the Spite plateau observed in the Milky Way, though slightly lower, potentially reflecting differences in the star formation history and chemical evolution of these two galaxies.

The Sculptor dwarf galaxy is another extragalactic system in which the lithium abundance in metal-poor stars has been studied, with observations confirming the presence of a Spite plateau in its oldest stars (Hill et al. 2019; Matteucci et al. 2021). More recently, Molaro et al. (2020) analyzed stars belonging to the Gaia-Enceladus galaxy and found a Spite plateau similar to that of the Milky Way. The detection of a lithium plateau across different environments, particularly in nearby dwarf galaxies, suggests that primordial lithium production was largely uniform across galaxies. Small variations observed in some systems may offer valuable insights into evolutionary processes influencing lithium in distinct galactic environments.

In this study we incorporated lithium into the semi-analytical model proposed by Corazza et al. (2022), originally developed to explore the cosmological evolution of metallicity in the Universe. This model accounts for contributions from Pop III and Pop II stars to the synthesis of 11 elements: Fe, Si, Zn, Ni, P, Mg, Al, S, C, N, and O.

In Sect. 2 we review the model from Corazza et al. (2022), hereinafter referred to as the CMW model, which couples the cosmic star formation rate (CSFR) with chemical evolution equations, and discuss how lithium is integrated into this framework. In Sect. 3 we examine the role of Pop III and Pop II stars in lithium evolution, emphasizing the significance of nova systems and Galactic cosmic rays (GCRs) for ⁷Li production when metallicity² exceeds [Fe/H]>  − 1.0. Conclusions are presented in Sect. 4.

2. Methodology

The cosmological model used by CMW assumes that DM halos form at redshift z ∼ 20 with masses on the order of 10⁶ M_⊙. More massive structures form later as the redshift decreases. Using the Press-Schechter formalism (Press & Schechter 1974), the number of DM halos per comoving volume at a given redshift z and with masses in the interval [M, M + dm] can be determined (see in particular Sect. 2.1 of CMW).

Once formed, DM halos generate potential wells that draw primordial gas into their interiors. In this way, the chemical elements synthesized by primordial nucleosynthesis (such as hydrogen, helium, and lithium) are gradually incorporated into the halos through an infall process. The physical parameter associated with this process is the baryon accretion rate, given by

$$\begin{array}{c}\hfill a_{\mathrm{b}}(t)={\mathrm{\Omega }}_{0,\mathrm{b}}{\rho }_{0,\mathrm{c}}{\left(\frac{\mathrm{d}t}{\mathrm{d}z}\right)}^{-1}\left|\frac{\mathrm{d}{f}_{\mathrm{b}}}{\mathrm{d}z}\right|,\end{array}$$(1)

where Ω_0, b is the baryonic density parameter at z = 0, ${\rho }_{0,\mathrm{c}}=3H_0^2/8\pi G$ is the critical density of the Universe (H₀ = 100 h km⁻¹ Mpc⁻¹ is the Hubble parameter at the present time), and f_b(z) is the fraction of primordial baryons incorporated into the halos as determined from the Press-Schechter formalism. The cosmological framework determining the function a_b(t) is similar to that used by other authors (e.g., Daigne et al. 2006; Pereira & Miranda 2010; Tan et al. 2016; Gribel et al. 2017; Vangioni et al. 2018).

As gas gradually accumulates inside virialized halos via a_b(t), conditions are created for stars to begin forming within these structures. Using Schmidt’s law (Schmidt 1959), which characterizes the conversion of gas into stars, and Salpeter’s law (Salpeter 1959), which characterizes the initial mass function (IMF), the CSFR represented by ${\dot{\rho }}_{\star }(z)$ can be determined (see Sect. 2.2 of Corazza et al. 2022). This characterizes the cycle associated with star formation, which is primarily related to the removal of gas to form stars, the return of chemically processed gas by stars to the environment, and the formation of stellar remnants at the end of the stars’ life.

The mechanisms regulating ${\dot{\rho }}_{\star }(z)$ are described by the star formation efficiency ⟨ε_⋆(z)⟩ and the star formation timescale τ_s(z), which represents the gas depletion timescale. As discussed in CMW, at high redshifts, the τ_s(z) function is dominated by the atomic gas depletion timescale, while for z < 1.0, τ_s(z) is dominated by the molecular depletion time. On the other hand, the ⟨ε_⋆(z)⟩ function exhibits nearly constant behavior for z > 3.5, with values between ∼0.15 and 0.40, depending on the IMF exponent. For z < 3.5, the star formation efficiency gradually decreases with redshift, reaching values close to 0.01 − 0.02 at z = 0.

Additionally, the cosmological scenario in this study incorporates the results of Gribel et al. (2017), in which the authors demonstrated the connection between the CSFR and local (Galactic) star formation (LSF). In particular, Gribel et al. (2017) identified that turbulence has a dual role in star formation: it induces the star formation efficiency to remain high and almost constant until a critical Mach number (M_crit) is reached, after which it induces a rapid decline in the star formation efficiency.

The unified model of Gribel et al. (2017) consistently derives Larson’s first law (Larson 1981) associated with the velocity dispersion (⟨V_rms⟩) in local star-forming regions. Finally, the most important characteristic of Gribel et al. (2017), relevant to this work, is that part of the star formation taking place inside a halo with mass ∼10⁶ M_⊙ occurs in substructures (or subhalos) with masses similar to those of globular clusters. It is within this cosmological framework that the analysis of the lithium abundance is performed.

2.1. Lithium yields of Pop III and Pop II stars

Once star formation begins at z ∼ 20, conditions are created for the chemical enrichment of the Universe. The first stars to form have zero metallicity (Pop III) and gradually enrich the environment. As discussed in CMW, upon reaching the transition metallicity Z = 10⁻⁶, conditions are created for the formation of the first generation of Pop II stars. As a result, Pop III stars are no longer formed, but due to the mass that governs the lifetime τ_m (Z = 0), some zero-metallicity stars can coexist with Pop II stars with a metallicity of 10⁻⁶, as long as τ_m (Z = 0) > t₁ (where t₁ is the moment in time when the metallicity of the Universe reaches 10⁻⁶).

When the metallicity of the Universe reaches Z = 10⁻⁴, the second generation of Pop II stars is formed. This generation can coexist with some of the stars from previous generations (Pop III and Pop II of the Z = 10⁻⁶ branch), as long as the conditions τ_m (Z = 0) > t₂ and τ_m (Z = 10⁻⁶) > t₂ are satisfied (where t₂ is the moment in time when the metallicity of the Universe reaches 10⁻⁴). This process repeats itself at each metallicity transition in the Universe determined by the values 10⁻³, 4 × 10⁻³, 8 × 10⁻³, and 2 × 10⁻². Thus, the stars of Pop II and Pop III are situated in different classes of evolutionary trajectories dictated by the metallicity of the Universe and the lifetime of the stars (τ_m), which is itself a function of the metallicity (see CMW and references therein). The reference t₀ is determined by the age of the Universe at redshift z = 20.

Tables 1 and 2 show the ranges of stellar masses and metallicities from which the yields were chosen. The Pop II branches shown in Table 2 were chosen to provide the best combination of masses and metallicity ranges, enabling the construction of a consistent chemical model encompassing transitions across stellar populations from metallicity 10⁻⁶ to 2 × 10⁻². Additionally, the Pop II stars in Table 2 share the same parametric basis – MONSTAR code for stellar evolution, OPAL opacities, and compatible mass loss models (see Corazza et al. 2022 and references therein) – ensuring that the chemical scenario is built on a common physical foundation.

2.2. Contribution of primordial lithium to the structures in formation

The function a_b(t) given in Eq. (1) supplies primordial gas to the structures formed, and it can be written as

$$\begin{array}{c}\hfill a_{\mathrm{b}}(t)=a_{\mathrm{b}-\mathrm{H}}(t)+a_{\mathrm{b}-\mathrm{He}}(t)+a_{\mathrm{b}-\mathrm{Li}}(t),\end{array}$$(2)

where a_b − H(t) is the accretion rate of primordial hydrogen, a_b − He(t) accounts for the accretion of primordial helium, and a_b − Li(t) is the corresponding accretion rate of primordial lithium. It is important to highlight that primordial gas is continuously injected into the star-forming environment through the function a_b(t). This term, in turn, depends on the halo mass function, involving halos of different masses that decouple from the Universe’s expansion at different redshifts. Once injected into the star-forming environment, part of the primordial gas is used to form stars, which is quantified, at each time t, by the parameter ⟨ε_⋆(z)⟩ – the star formation efficiency (Corazza et al. 2022). The return of (enriched) gas to the environment occurs more slowly and is regulated by the stellar lifetimes τ_m, which in turn depend on both stellar mass and metallicity.

As an ansatz, we assumed that the following relations hold and remain valid throughout the entire process of star formation within the halos:

$$\begin{array}{cc}\hfill a_{\mathrm{b}-\mathrm{H}}(t)& =0.7512a_{\mathrm{b}}(t),\hfill \end{array}$$(3)

$$\begin{array}{cc}\hfill a_{\mathrm{b}-\mathrm{He}}(t)& =0.2487a_{\mathrm{b}}(t),\hfill \end{array}$$(4)

$$\begin{array}{cc}\hfill a_{\mathrm{b}-\mathrm{Li}}(t)& =4.94\times {10}^{-10}{a}_{\mathrm{b}}(t).\hfill \end{array}$$(5)

Once the parameters and functions are characterized, it is possible to write the differential equation for the mass density of the element lithium as

$$\begin{array}{c}\hfill \frac{\mathrm{d}{\rho }_{\mathrm{Li}}}{\mathrm{d}t}={\displaystyle {\int }_{m(t)}^{m_{\mathrm{s}}}}[(m-m_{\mathrm{r}})Z_{\mathrm{Li}}(t-{\tau }_{\mathrm{m}})+P_{Z_{\mathrm{Li},\mathrm{m}}}]{\dot{\rho }}_{\star }(t-{\tau }_{\mathrm{m}})\\ \hfill \phi (m)dm-Z_{\mathrm{Li}}{\dot{\rho }}_{\star }(t)+a_{\mathrm{b}-\mathrm{Li}}(t),\end{array}$$(6)

where the term (m − m_r) Z_Li (t − τ_m) accounts for the amount of lithium incorporated when the star was born, which later returns to the interstellar medium (m_r Z_Li (t − τ_m) is the part of the lithium retained in the stellar remnant), and the term Z_Li = ρ_Li/ρ_g (where ρ_g is the total gas density). The P_ZLi, m parameter is the lithium mass produced by a star of mass m, the term $Z_{\mathrm{Li}}{\dot{\rho }}_{\star }(t)$ takes into account the removal of part of the lithium element to form a new generation of stars, and a_b − Li(t) is the primordial lithium incorporated into structures.

The total gas density contained in the term Z_Li = ρ_Li/ρ_g is obtained by solving the equation

$$\begin{array}{c}\hfill {\dot{\rho }}_{\mathrm{g}}=-\frac{d^2{M}_{\star }}{\mathit{dVdt}}+\frac{d^2{M}_{\mathrm{ej}}}{\mathit{dVdt}}+[a_{\mathrm{b}-\mathrm{H}}(t)+a_{\mathrm{b}-\mathrm{He}}(t)],\end{array}$$(7)

where the first term on the right-hand side describes the removal of gas to form stars, the second term describes the return of gas to the halos through material ejected at the end of the stars’ life cycle, and the last term represents the infall of primordial gas.

There are some differences between Eq. (6) of this study and that used by CMW. Those authors wrote Eq. (6) for each of the i metals they studied, satisfying the initial condition (at z = 20) ρ_g, i = 0. Thus, the term a_b(t) was included through the total gas density at Z_i = ρ_g, i/ρ_g. In the present study, primordial lithium is continuously inserted into the structures, an effect described by the last term of Eq. (6). In turn, the contribution of primordial gas to the total density, according to Eq. (7), is mainly due to hydrogen and helium, as represented by the term [a_b − H(t)+a_b − He(t)] in Eq. (7). With these conditions established, the numerical integration of Eq. (6) gives the mass density ρ_Li contained in the structures (DM halos). This allows us to determine the quantity ⁷Li/H and compare it with observational data.

In Sect. 3 we analyze the role of the a_b − Li function in the ⁷Li/H ratio, as well as the contributions of Pop III and Pop II stars to the results. For comparison with the model, the following observational datasets were selected: 7 extremely metal-poor stars (EMPSs) from Bonifacio et al. (2018), 41 plateau stars from the sample of Meléndez & Ramírez (2004), 49 Gaia-Enceladus sources from the sample of Molaro et al. (2020), the lithium abundance in the Small Magellanic Cloud (SMC) as determined by Howk et al. (2012) and Izzo et al. (2022), and the ⁷Li/H abundance in Solar System meteorites presented in Asplund et al. (2009). We also included two extremely iron-poor stars, namely J0023+0307, an extremely dwarf main-sequence (MS) star (Aguado et al. 2019), and SMSS J0313–6708, considered one of the oldest stars in the Galaxy (Keller et al. 2014). Finally, we took the primordial limits for the ratio (⁷Li/H)_{BBN + CMB} from the results of Yeh et al. (2021) and Fields & Olive (2022).

3. Results and discussion

The parameters used for modeling and characterizing the cosmological framework are the same as those used by CMW. These parameters are summarized in Table 3.

3.1. The formation of the Spite plateau

The first result of this investigation is the confirmation that the model described and characterized in the previous section accurately reproduces the lithium abundance observed in the stars distributed across the Spite plateau. As shown in Fig. 1, the gray shaded area represents the typical Spite plateau, extending in the range −3.7 ≤ [Fe/H] ≤ −1.7. The gray dashed line marks the extension of the Spite plateau up to [Fe/H] ∼ −8.0, allowing us to evaluate the behavior of the model for the ⁷Li/H versus [Fe/H] ratio. It is important to note that, without adding any new parameters to the CMW model, and by simply applying the ansatz defined by Eqs. (3), (4), and (5), we can confirm that the ⁷Li abundance predicted by this model fits very well with the usual Spite plateau. Moreover, the model suggests that the Spite plateau extends up to a metallicity of [Fe/H] ≲ −8.0.

Fig. 1. Lithium abundance as a function of [Fe/H]. The gray shaded area represents the usual Spite plateau. The primordial infall term in association with the formation of Pop III and Pop II stars causes the Spite plateau to be extended to [Fe/H] ≲ −8.0.

The yellow shaded band indicates the limits for the primordial ⁷Li abundance determined from the combination of BBN and CMB data. The difference between (⁷Li/H)_{BBN + CMB} ∼ 5.0 × 10⁻¹⁰ and the model’s predicted value of (⁷Li/H)∼1.8 × 10⁻¹⁰ is primarily due to the capture of primordial lithium during the birth of Pop III stars within structures formed at z ∼ 20. The high star formation efficiency ⟨ε_⋆⟩∼0.34 is responsible for removing primordial gas from the system, keeping the ⁷Li abundance within the Spite plateau.

Note that the Pop III branch is divided into three classes, as shown in Table 1. Once formed, the stars follow a distribution ∝m^{−(1 + x)}, where x is the IMF exponent. Stars from the CL08 (low-mass) and HW10 (intermediate-mass) classes of the Pop III branch return some of the primordial lithium they initially incorporated at the end of their lives. During their evolution, stars in these classes also produce lithium, which is partially returned to the system when they reach the end of their life cycle. On the other hand, stars from the HW02 class end their lives as pair-instability supernovae (PISNe), which removes ⁷Li but does not return it. Additionally, the HW02 group leaves no remnants after the disruptive process is complete.

When the metallicity of the Universe exceeds the threshold Z = 10⁻⁶, the model assumes that no new Pop III stars can form. However, stars formed within the metallicity range [0 − 10⁻⁶] will continue their evolutionary trajectories as dictated by their lifetimes τ_m (see Sect. 2.1). In other words, the threshold Z = 10⁻⁶ marks the beginning of Pop II stars’ contribution to the chemical enrichment of the Universe, starting with the CL04 group at Z = 10⁻⁶. The evolutionary paths of the other Pop II branches are described in Sect. 2.1.

We analyzed three different values for the IMF exponent x: 0.85, 1.00, and 1.35 (the Salpeter exponent). The results are only weakly dependent on the value of x. Figure 1 illustrates the importance of the term a_b − Li in establishing the Spite plateau when the model includes contributions from both Pop III and Pop II stars.

Excluding the a_b − Li term from Eq. (6), as shown by the red curve for x = 1.35, results in the ⁷Li/H ratio reaching the Spite plateau only when [Fe/H] > −2.0. For metallicities above [Fe/H] > −2.0, the curves with (black) and without (red) the a_b − Li term show similar behavior. Thus, at first glance, the Spite plateau could be seen as a signature of primordial infall. Although the zero-metallicity stars from the CL08 and HW10 branches incorporate the primordial abundances of H, He, and Li, the fraction of ⁷Li ejected by these stars at the end of their lives is generally smaller than their initial compositions. In contrast, the HW02 branch removes ⁷Li during the star formation process and does not return any ⁷Li at the end of these stars’ lives. Therefore, the term a_b − Li is crucial for establishing the Spite plateau at high redshifts (or low metallicities as represented by [Fe/H] < −2.0).

Once the Universe begins to receive contributions from Pop II stars’ mass ejections, particularly from the Z = 10⁻³ branch, the ⁷Li abundance no longer strongly depends on the a_b − Li term. This occurs when [Fe/H] > −2.0. In Sect. 3.4 we assess this analysis by removing the HW02 branch from the model. Finally, the model accurately reproduces the ⁷Li abundance in the SMC.

Thus, by incorporating the ansatz (Eqs. (3), (4), and (5)) into the CMW model, we can obtain the observed ⁷Li/H ratio of the Spite plateau and align it with the primordial abundance. Furthermore, as discussed earlier, the Spite plateau is naturally generated and can extend to very low metallicities.

It is important to highlight that Fig. 1 shows a large dispersion of the lithium abundance for [Fe/H] < 4.0, which corresponds to the EMPS sample. These seven stars were analyzed by Bonifacio et al. (2018), and as discussed by those authors, stars with T_eff > 6000 K and [Fe/H] = −4.0 lie on the Spite plateau, whereas EMPSs with T_eff < 6000 K and [Fe/H] < −3.5 always lie below it. According to the authors, this suggests that the T_eff threshold at which lithium destruction becomes important increases as metallicity decreases. This is likely one of the factors involved in the observed spread relative to the Spite plateau. Other factors likely contribute to this spread as well, such as turbulence and atomic diffusion.

However, the model does not reproduce the lithium abundance observed in Solar System meteorites. This discrepancy may indicate the need to introduce a new source of lithium not accounted for in the CMW model.

3.2. Lithium production by nova systems

As discussed in Cescutti & Molaro (2019), the detection of lithium in classical nova outbursts suggests that these systems could play a significant role in the production and enrichment of ⁷Li in galactic systems. In particular, only binary systems formed by stars with 0.8 < M/ M_⊙ < 8 can develop nova systems capable of producing ⁷Li. As a result, the total mass (M_B) of the ⁷Li−producing binary system lies in the range 1.6 < M_B/ M_⊙ < 16. The distribution function for the primary and secondary stars in the binary system is given by (see Greggio & Renzini 1983; Cescutti & Molaro 2019 for details)

$$\begin{array}{c}\hfill f(\mu )=2^{1+\gamma }(1+\gamma ){\mu }^{\gamma }(0<\mu \le 1/2),\end{array}$$(8)

where μ = M_sec/M_B, with M_sec being the mass of the secondary star, and γ = 2.

Additionally, for ⁷Li production to occur in a nova system, the primary star must evolve into a white dwarf and accumulate enough material from the secondary star to reach the ignition conditions for the nova outburst. This implies that binaries with stars of the same mass cannot develop nova systems. The lifetimes of the primary and secondary stars must satisfy the condition τ_sec − τ_prim ≥ τ_nova, where τ_nova is the delay time.

The ⁷Li production in nova systems can be determined by

$$\begin{array}{cc}\hfill {\dot{\rho }}_{\mathrm{N}-\mathrm{Li}}& ={\mathcal{N}}_{\mathrm{Li}}{\displaystyle {\int }_{M_{\mathrm{Bm}}}^{M_{\mathrm{BM}}}}\varphi (M_{\mathrm{B}})\hfill \\ \hfill & \left\{{\displaystyle {\int }_{{\mu }_{{\mathrm{m}}_{\min }}}^{1/2}}f(\mu ){\dot{\rho }}_{\star }(t-{\tau }_{\sec })d\mu \right\}d{M}_{\mathrm{B}}{\mathcal{Y}}_{\mathrm{N}-\mathrm{Li}},\hfill \end{array}$$(9)

where μ_mmin is the minimum mass fraction contributing to the nova rate at time t, ϕ(M_B) is the IMF for the total mass of the binary system ($\varphi (M_{\mathrm{B}})\propto M_{\mathrm{B}}^{-(1+x)}$), ${\dot{\rho }}_{\star }(t-{\tau }_{\sec })$ is the CSFR at the time when the secondary star is formed, 𝒴_N − Li is the amount of ⁷Li produced in each nova event, and 𝒩_Li is a fixed parameter to reproduce the rate of nova outbursts in the Galaxy at the present time.

The nova rate is determined following the formalism presented in Matteucci et al. (2003). As discussed by these authors (and also by Cescutti & Molaro 2019), the nova outburst rate is related to the white dwarf formation rate in binary systems. Each nova system can produce an average of 10⁴ nova outbursts. These points are considered in determining the parameter 𝒩_Li, which is adjusted to reproduce a nova rate of 50 yr⁻¹, in agreement with recent estimates of the Galactic nova rate, yielding a value of ${50}_{-23}^{+31}{\mathrm{yr}}^{-1}$ (Shafter 2017). We determined the yield parameter 𝒴_N − Li to match the observed ⁷Li abundance in meteorites. Finally, we adopted τ_nova = 1 Gyr for the delay time (Matteucci et al. 2003; Cescutti & Molaro 2019).

In Fig. 2 we present the contribution of nova outbursts to the lithium yield by incorporating Eq. (9) into the right-hand side of Eq. (6). Nova systems become effective at producing ⁷Li only for [Fe/H] > −1.0. This delay is due to the time required for the primary star to evolve into a white dwarf, as well as the time necessary for the white dwarf to cool enough for a nova outburst to occur.

Fig. 2. Lithium abundance as a function of [Fe/H]. The figure shows the contribution of nova binary systems to the lithium abundance when [Fe/H] > −1.0.

The parameters 𝒩_Li and 𝒴_N − Li required to produce a Galactic nova rate of 50 yr⁻¹ and return the lithium abundance observed in meteorites are shown in Table 4. We used τ_nova = 1 Gyr in this model, consistent with the fiducial model of Cescutti & Molaro (2019). It is important to highlight some differences between this work and that of Cescutti & Molaro (2019). While they explored the values 0, 1, 2, and 5 Gyr for τ_nova, they considered 1 Gyr as the best value for the delay time. Another key difference concerns the parameter 𝒩_Li. They found that a value of 0.03 provided a nova burst rate compatible with 20 − 30 yr⁻¹. However, we adopted the most recent value (see Shafter 2017) of ${50}_{-23}^{+31}{\mathrm{yr}}^{-1}$, which is approximately 2 − 3 times larger than the rate used by Cescutti & Molaro (2019).

In our model, for a nova outburst with x = 1.35, we obtain 𝒩_Li = 0.066, a result consistent with the analogous factor of Cescutti & Molaro (2019). Finally, Cescutti & Molaro (2019) found that the lithium abundance in meteorites could be matched with 𝒴_N − Li = 4.14 × 10⁻⁵ M_⊙, while we obtain agreement with 𝒴_N − Li = 5.4 × 10⁻⁵ M_⊙ (for x = 1.35). Assuming each nova system undergoes ∼10⁴ outbursts during its lifetime, the required lithium production per nova outburst is approximately 5.4 × 10⁻⁹ M_⊙ in our model.

3.3. Lithium from the spallation of CNO by Galactic cosmic rays

Type II supernovae are natural sites for the generation and acceleration of high-energy cosmic rays, which can, in turn, produce elements such as lithium, beryllium, and boron through spallation of CNO nuclei in the Galactic environment (Ryan et al. 2000; Prantzos 2012; Matteucci et al. 2021). As discussed by Prantzos (2012), GCRs may play an important role in ⁷Li production for [Fe/H] > −3.0, possibly complementing the contribution of novae. To include GCRs in our model, we obtained a fit consistent with the curve shown in Fig. 1 of Prantzos (2012) for [Fe/H] > −3.0. Combining GCRs and novae, it is possible to reach the meteoritic ⁷Li abundance with 9.3% from GCRs and 90.7% from novae. This allows for a 16.7% reduction in the total lithium yield in the model with IMF slope x = 1.35 (fiducial model from Table 4), resulting in 𝒴_N − Li = 4.5 × 10⁻⁵ M_⊙. This value is closer to the estimate of Cescutti & Molaro (2019). The yield reduction factor of 16.7% also applies to the other models in Table 4, particularly with 𝒴_N − Li = 1.67 × 10⁻⁴ M_⊙ for x = 1.00 and 2.50 × 10⁻⁴ M_⊙ for x = 0.85.

3.4. Lithium production without the HW02 branch of Population III stars

Between zero and 10⁻⁶ metallicity, Pop III stars from branches HW02 and HW10 contribute predominantly to the chemical enrichment of the environment, with HW02 not producing lithium but instead removing it from the system. Thus, the lithium content in the star-forming environment is slightly higher – specifically, 1.002 (or ∼0.2% higher) – when only the HW10 branch is present. The presence or absence of the HW02 branch has a more significant impact on the amount of iron. Without HW02, iron injection into the star-forming environment occurs more slowly. The amount of iron injected jointly by HW02 and HW10 is approximately 3.2 times greater than that injected by HW10 alone during the metallicity phase from zero to 10⁻⁶. Therefore, removing HW02 Pop III stars slightly increases lithium but substantially reduces iron in the metallicity range between zero and 10⁻⁶. Since the reference metallicity of 10⁻⁶ marks the beginning of Pop II star formation, the subsequent chemical evolution in the star-forming environment becomes linked to the contribution of Pop III stars.

In particular, a Pop III star with an initial mass of ∼200 M_⊙ (HW02 branch) returns ∼100 M_⊙ in metals and ∼10 M_⊙ in iron (Heger & Woosley 2002). PISNe end their lives without leaving any remnants due to a complete disruptive process. Given the characteristics of these stars, the results for [Fe/H] depend on the presence of the PISN branch. Moreover, these stars remove ⁷Li from the environment and do not return it at the end of their lives.

On the other hand, a Pop III star from the HW10 branch with ∼100 M_⊙ returns ∼9 M_⊙ in metals, ∼10⁻¹² M_⊙ in iron, and ∼10⁻⁶ M_⊙ in lithium, assuming the standard model with an energy explosion of 1.2 B (see Heger & Woosley 2010 for details). As a result, the [Fe/H] versus redshift relationship depends on the presence (or absence) of HW02 branch stars. Figure 3 illustrates the behavior of the [Fe/H] versus ⁷Li/H relationship without the PISN branch. A comparison of the results shown in Figs. 2 and 3, with a_b − Li, reveals virtually no difference. There is only a slight increase in the amplitude of ⁷Li/H for the model sans PISNe (Fig. 3) compared to the model with PISNe (Fig. 2). This happens because stars from the HW10 and CL08 branches synthesize ⁷Li, while stars from the PISN branch (HW02) only remove this element from the environment.

The main difference between the models in Figs. 2 and 3 is observed when we compare the red curves for the models without the primordial infall term. In the model with PISNe, the transition from Pop III to Pop II (with metallicity Z = 10⁻⁶) occurs when [Fe/H] ∼ −2.3. This transition happens very quickly due to the large amount of iron injected by the HW02 stars. Only when [Fe/H] > −2.3 do the stars of the HW02 branch stop forming and removing ⁷Li from the system. During this phase, lithium enrichment in the structures primarily comes from synthesis by stars in the HW10 branch, while stars in the HW02 branch are responsible for determining [Fe/H].

Fig. 3. Lithium abundance versus [Fe/H]. In this case, only the CL08 (low-mass) and HW10 (intermediate-mass) branches of Pop III stars are included.

For the sans PISNe model (Fig. 3), during the phase characterized by the chemical enrichment of the Universe from zero metallicity up to Z = 10⁻⁶, stars from the HW10 branch are responsible for enriching both ⁷Li and [Fe/H]. The transition to Pop II is slower without the HW02 branch, occurring when [Fe/H] ∼ −3.1 (corresponding to Z = 10⁻⁶). Although Pop III stars no longer form when Z > 10⁻⁶, the lower-mass stars from the HW10 branch, along with those from the CL08 branch, continue to contribute to ⁷Li injection into the system up to [Fe/H] ∼ −0.84, at which point they end their lives, and lithium synthesis now depends solely on Pop II stars from the various branches listed in Table 2.

Despite the fact that the sans PISNe model achieves a ⁷Li/H ∼ 1.20 × 10⁻¹⁰ without the contribution of the infall term, this value represents the lower limit of the usual Spite plateau. Only by including the primordial infall term can ⁷Li/H ∼ 1.81 × 10⁻¹⁰ be maintained, which is closer to the central value of the Spite plateau. From Fig. 3, we can see that without the term a_b − Li, in the sans PISNe model, the extended Spite plateau breaks down for [Fe/H] < −6.7, whereas with the primordial infall term (black lines in Fig. 3), the Spite plateau is maintained up to [Fe/H] ≲ −8.0.

3.5. Lithium without the contribution of the Pop III stars

Figure 4 shows the behavior of ⁷Li/H considering only Pop II stars. In this case, similar to the study by CMW, it is assumed that zero-metallicity stars have yields identical to Pop II stars with Z = 10⁻⁶ (CL04 branch, as indicated in Table 2). These stars are responsible for synthesizing lithium, iron, and other elements, chemically enriching the Universe up to Z = 10⁻⁴, when stars from the subsequent branches such as K10, D14b, and CL04 begin to form, as described in Sect. 2.1.

Fig. 4. Lithium abundance versus [Fe/H] considering the contribution of Pop II stars and nova systems working together with GCRs. The Spite plateau is established by the primordial infall term. Even with the infall term, only Pop II stars, novae, and GCRs do not return the lithium abundance inferred for Solar System meteorites (see the main text).

Stars from the CL04 branch (with Z = 10⁻⁶) return 3.5 M_⊙ in metals, 10⁻¹ M_⊙ in iron, and ∼3.0 × 10⁻¹¹ M_⊙ in lithium at the end of their lives. The cycle characterized by the formation of these stars continues until the metallicity of Z = 10⁻⁴ is reached, at which point the formation of stars with this metallicity begins. Although the iron production of stars from the CL04 (with Z = 10⁻⁶) and HW10 branches is similar, Pop II stars with Z = 10⁻⁶ produce one-third of the total amount of metals and much less lithium compared to HW10 stars. As a result, the values of [Fe/H] and ⁷Li/H increase much more slowly over time.

It can be seen in Fig. 4 that the Spite plateau is established through the primordial infall term (a_b − Li). The lithium production by Pop II stars with Z = 10⁻⁶ is ∼10⁵ times lower than the average yield of Pop III stars from the HW10 branch. The next class of Pop II stars with Z = 10⁻⁴ contributes very weakly to the production of ⁷Li. Most stars in the 1 − 7 M_⊙ range have negative yields, contributing to the removal of ⁷Li from the system. Stars in the 7.5 − 35 M_⊙ range, on average, synthesize five times less lithium than stars in the Z = 10⁻⁶ class.

The stars of the next metallicity class (Z = 10⁻³) return an average of ∼10⁻⁸ M_⊙ of ⁷Li, which is greater than the contribution from stars in the less enriched Pop II classes. The contribution from novae occurs simultaneously with that of stars from the Z = 10⁻³ class, increasing the value of ⁷Li/H when −0.13 > [Fe/H] > −0.70. Pop II stars from subsequent classes (i.e., Z = 4 × 10⁻³ and Z = 8 × 10⁻³) mostly have negative yields, leading to the incorporation and destruction of ⁷Li in the environment. During this phase, the enrichment of lithium in the structures comes from the infall of primordial gas and nova systems together GCRs. Additionally, stars from these two metallicity classes produce almost no iron. This explains thebehavior of the ⁷Li/H versus [Fe/H] ratio, where the abundance of lithium increases while the availability of iron in the environment decreases.

Although the inclusion of the term a_b − Li produces the Spite plateau, it is important to emphasize that by considering only the yields from Pop II stars and novae working together with GCRs, the relation between ⁷Li/H and [Fe/H] observed in Solar System meteorites cannot be recovered. This is due to the characteristics of stars with metallicities of 4 × 10⁻³ and 8 × 10⁻³, as discussed earlier. On the other hand, the removal of the term a_b − Li (red curve in Fig. 4) does not provide a physically adequate model for the ⁷Li/H versus [Fe/H] relation.

It is crucial to emphasize, as highlighted in Sect. 2.1, that the physical basis for the Pop II yields used in this study includes the MONSTAR code for stellar evolution, OPAL opacities, and compatible mass loss models, providing a common physical foundation for chemical modeling across stellar populations from metallicity 10⁻⁶ to 2 × 10⁻² (Corazza et al. 2022). In particular, MONSTAR treats the transport of chemical elements by convection only, not accounting for diffusion or turbulence. On the other hand, codes such as MESA+NuGrid (see, e.g., Ritter et al. 2018) incorporate microscopic diffusion, rotation, and magnetic fields.

Comparison between the ⁷Li yields from this work and those from MESA+NuGrid is not direct for the same stellar mass and metallicity ranges presented in Table 2. Within the possible comparisons for Pop II stars at specific metallicities Z ∼ 10⁻⁶ and 10⁻⁴, the present work shows higher ⁷Li yields for the mass range 5–25 M_⊙, while above 25–60 M_⊙ MESA+NuGrid with diffusion gives higher lithium yields. Integrating the total ⁷Li yield from 5 to 60 M_⊙, this work produces approximately 1.5 times more lithium than the diffusion and rotation inclusive models. This result agrees with the analysis by Richard et al. (2005). Thus, within the possible comparisons to be made between MONSTAR and MESA+NuGrid, the removal of the a_b − Li term would not produce a scenario capable of generating the Spite plateau with only Pop II stars.

3.6. J0023+0307 and SMSS J0313–6708

While this work does not aim to evaluate the ⁷Li/H versus [Fe/H] ratio for individual objects, it is still worthwhile to make some observations about J0023+0307 and SMSS J0313–6708 in light of the model presented here. The inferred abundances [Fe/H] ≲ −6.1 ([Fe/H] ≲ −7.1 for SMSS J0313–6708) classify these objects as extremely iron-poor stars.

J0023+0307 is a dwarf star with A(Li) = 2.02³ (Aguado et al. 2019), placing it at the lower limit of the Spite plateau. Since no iron lines were detected in the spectrum of this star, the upper limit of [Fe/H] = −6.1 was inferred from the resonant Ca II lines, assuming [Ca/Fe] ≥ 0.40 (Aguado et al. 2019). The chemical abundance analysis of J0023+0307 by Frebel et al. (2019), using supernova yield models from HW10, suggested that J0023+0307 formed in a re-collapsed minihalo. Specifically, using nonlocal thermodynamic equilibrium abundances, the authors identified four best-fit models from the entire HW10 database. Using the StarFit algorithm⁴, we obtained the same set of models as Frebel et al. (2019). The model with the lowest χ² corresponds to a star with a mass of 50 M_⊙, an explosion energy of 3.0 B, a mixing parameter of 3.98 × 10⁻³, and a remnant black hole mass of ∼3.6 M_⊙.

The condition [Fe/H] = −6.1 is met in our model when the total metallicity of the Universe is Z = 2.3 × 10⁻¹⁰ in the scenario that includes PISNe (Z = 9.8 × 10⁻¹⁰ without PISNe) and assuming an IMF exponent of x = 1.35. J0023+0307 would have formed approximately 4.4 × 10⁵ yr after the explosion of the first PISN or ∼1.3 × 10⁶ yr after the explosion of the first star in the HW10 branch (model without PISNe).

The statistical analysis by Frebel et al. (2019) suggests that the inferred dilutions of hydrogen masses (M_dil) from their fits range from 10^4.5 − 10^5.5 M_⊙, which is consistent with the idea that J0023+0307 formed in a re-collapsed minihalo. Since the present study, along with the CMW model, incorporates the structure of the unified model developed by Gribel et al. (2017), it is also plausible that J0023+0307 could have formed in a subhalo (a globular cluster-like substructure) with a mass of ∼4.6 × 10⁴ M_⊙, contained within one of the first halos with a mass of 10⁶ M_⊙. While these scenarios are distinct, there is a similarity between our findings and the proposal of a re-collapsed minihalo by Frebel et al. (2019) to explain the chemical abundances of this extremely iron-poor star.

SMSS J0313–6708 is a red giant with [Fe/H] < −7.1 and A(Li) = 0.7 (Keller et al. 2014; Frebel & Norris 2015). In particular, Keller et al. (2014) compared the abundances of SMSS J0313–6708 with the yields of Pop III supernovae and obtained a good fit to a progenitor star with a mass of 60 M_⊙ and primordial composition. Using the HW10 database and the StarFit algorithm in this study, we reached a similar conclusion: a progenitor star with 60 M_⊙, an explosion energy of 1.8 B, a mixing parameter of 1.58 × 10⁻³, and a remnant black hole mass of 25 M_⊙.

The condition [Fe/H] = −7.1 is satisfied in our model when the total metallicity of the Universe is Z = 2.3 × 10⁻¹¹ in the PISN scenario (Z = 2.3 × 10⁻¹⁰ without PISNe). The [Fe/H] = −7.1 abundance is reached approximately 2.2 × 10⁵ yr after the first PISN event or ∼4.4 × 10⁵ yr after the explosion of the first star in the HW10 branch (model without PISNe).

As SMSS J0313–6708 is in the red giant phase, its lithium abundance must have decreased during its evolution from the MS to the upper red giant branch (RGB), due to the destruction of lithium in the convective envelope during the RGB phase. As discussed by Frebel & Norris (2015), it is possible that the lithium abundance was depleted by approximately 1.3 dex, which would imply that the star initially had A(Li)_MS ∼ 2.0, placing it near the lower boundary of the Spite plateau.

4. Conclusions

In this work we address the cosmological lithium problem by employing a chemical model that incorporates yields from stars with metallicities ranging from zero to solar. The present scenario follows the unified CSFR-LSF model proposed by Gribel et al. (2017), which consistently derives Larson’s law, among other results. The cosmological structure formation framework used here mirrors that of Pereira & Miranda (2010), Gribel et al. (2017), and Corazza et al. (2022) and is based on the formalism originally developed by Press & Schechter (1974) to describe the formation and evolution of large-scale structures in the Universe.

The results of this work show that the initial assumption introduced by the ansatz, expressed in Eqs. (3), (4), and (5), is consistent with explaining the lithium abundance over cosmic time. In particular, the introduction of the ansatz shows that the so-called cosmological lithium problem can be addressed using a standard chemical evolution model coupled to the hierarchical structure formation scenario. Additionally, the model presented in this work does not introduce any additional parameters to the scenario described by Gribel et al. (2017) and Corazza et al. (2022), thereby preserving all the results obtained by those authors. This strengthens the interpretation that the ansatz considered here is valid throughout the star formation process.

Additionally, diffusion processes, which reduce the lithium abundance, are evaluated in studies that complement the scenario proposed here. Using the Spite plateau reference established in this work, these processes could decrease lithium abundances by a factor of roughly 1.6–2.0, as determined by Richard et al. (2005). This reduction might partly account for the dispersion observed in the lithium abundance versus [Fe/H] among the samples presented in Figs. 1 through 4.

Our main conclusions are summarized as follows:

• There is no so-called cosmological lithium problem, i.e., there is no discrepancy between the primordial abundance of ⁷Li determined from combined BBN+CMB predictions and the inferred abundance in metal-poor stars near the MS.

• The use of the ansatz enables the generation of the Spite plateau for the ⁷Li/H versus [Fe/H] relation, in agreement with the BBN+CMB predictions and inferred abundances for several samples: EMPSs (Bonifacio et al. 2018), plateau stars (Meléndez & Ramírez 2004), Gaia-Enceladus (Molaro et al. 2020), and the SMC (Howk et al. 2012; Izzo et al. 2022). The term a_b − Li, together with processes associated with stellar formation and the yields of Pop III and Pop II stars, produces an extended Spite plateau within the range −8.0 ≲ [Fe/H] ≲ −2.0. The inclusion of a_b − Li allows us to obtain the Spite plateau with ⁷Li/H ∼ 1.81 × 10⁻¹⁰, which is consistent with the result found by Fields & Olive (2022) and references therein.

• The enrichment of ⁷Li over cosmic time is weakly dependent on the IMF. However, x = 1.35 is the exponent that best brings the parameters of our scenario in line with the fiducial parameters used by Cescutti & Molaro (2019) in their study of nova binary systems. In our model, the abundance of ⁷Li in meteorites is obtained from the combined contribution of novae and GCRs.

• J0023+0307 could have formed 4.4 × 10⁵ yr after the first PISN explosion or 1.3 × 10⁶ yr after the first star of the HW10 branch (model without PISNe). The results in Sect. 3.6 are similar to the analysis performed by Frebel et al. (2019), particularly their suggestion that the hydrogen mass dilutions (M_dil) inferred from their fits could indicate that J0023+0307 formed in a re-collapsed minihalo. This idea is consistent with the results of Gribel et al. (2017), who propose that part of the star formation at high redshifts occurred in DM subhalos with masses similar to those of globular clusters. In this study, J0023+0307 could have formed in a subhalo-like structure of ∼4.6 × 10⁴ M_⊙, contained in one of the first halos formed with a mass of 10⁶ M_⊙ at redshift 20.

• SMSS J0313–6708 could have formed 2.2 × 10⁵ yr after the first PISN event or 4.4 × 10⁵ yr after the explosion of the first star in the HW10 branch (model without PISNe). As SMSS J0313–6708 is in the red giant phase, its lithiumabundance must have decreased during its evolution from the MS to the RGB. As discussed by Frebel & Norris (2015), it is possible that this star initially had A(Li)_MS ∼ 2.0 and was depleted by 1.3 dex during the RGB phase. Therefore, SMSS J0313–6708 could have initially had a lithium abundance close to the Spite plateau.

Although this work demonstrates that the Spite plateau can be generated by three different scenarios, the formation of large-scale structures in the Universe, as described by the Press-Schechter-like formalism, suggests that the a_b − Li term cannot be omitted. The DM halos that decouple from the expansion of the Universe, collapse, and virialize are immersed in the primordial environment, which consists of elements synthesized during primordial nucleosynthesis. It is in this environment that the potential wells of the halos create conditions for baryonic matter to flow into the structures, agglomerate, and form the first stars. Thus, the physically appropriate scenarios are those that include the a_b − Li term, operating together with Pop III stars in the low-metallicity regime.

The remaining question to resolve is whether the HW02 branch could have contributed to the chemical enrichment of the Universe. Although the analyses performed here for J0023+0307 and SMSS J0313–6708 support the results of other authors, particularly that the chemical abundances of these stars are compatible with progenitors of the HW10 branch and not the HW02 branch, it is possible that, if part of the star formation at high redshifts occurred in subhalos (see Gribel et al. 2017) or in re-collapsed minihalos (Frebel et al. 2019), the contribution of PISNe to the establishment of the extended Spite plateau would not be ruled out.

Finally, the cosmological lithium problem can be resolved without any changes to standard BBN. The hierarchical structure formation scenario, as described by Press & Schechter (1974), with the appropriate characterization of the CSFR (Daigne et al. 2006; Pereira & Miranda 2010; Gribel et al. 2017; Vangioni et al. 2018) and yields from Pop III and II stars (Corazza et al. 2022), along with Schmidt’s law for star formation (Schmidt 1959), the Salpeter IMF (Salpeter 1959), nova systems (Greggio & Renzini 1983; Matteucci et al. 2003; Shafter 2017; Cescutti & Molaro 2019), and GCRs (Ryan et al. 2000; Prantzos 2012; Matteucci et al. 2021), produces a consistent scenario for the observed ⁷Li abundance over cosmictime.

References

All Tables

All Figures

	Fig. 1. Lithium abundance as a function of [Fe/H]. The gray shaded area represents the usual Spite plateau. The primordial infall term in association with the formation of Pop III and Pop II stars causes the Spite plateau to be extended to [Fe/H] ≲ −8.0.
In the text

	Fig. 2. Lithium abundance as a function of [Fe/H]. The figure shows the contribution of nova binary systems to the lithium abundance when [Fe/H] > −1.0.
In the text

	Fig. 3. Lithium abundance versus [Fe/H]. In this case, only the CL08 (low-mass) and HW10 (intermediate-mass) branches of Pop III stars are included.
In the text

	Fig. 4. Lithium abundance versus [Fe/H] considering the contribution of Pop II stars and nova systems working together with GCRs. The Spite plateau is established by the primordial infall term. Even with the infall term, only Pop II stars, novae, and GCRs do not return the lithium abundance inferred for Solar System meteorites (see the main text).
In the text