© The Authors 2025

1 Introduction

The dimension of time broadens our horizon in understanding the processes of formation and evolution of our Galaxy. The main tools for going back in time are stellar ages. The usual method adopted to derive stellar ages consists of comparing observed properties, such as magnitudes and colours, or inferred quantities, such as effective temperature (T_eff) and surface gravity (log g), with the output of stellar evolution models, the so-called isochrone fitting (see, e.g. Mints & Hekker 2018). However, ages of field stars are difficult to obtain with precision, particularly in some phases of stellar evolution, such as the low main sequence or the red giant branch, where there is strong degeneration among isochrones used to derive ages.

To overcome this problem, there are at least two possible alternative approaches: the use of asteroseismology (see, e.g. Chaplin & Miglio 2013) and the use of empirical methods that link some properties of stars with their age. The former is still model dependent, but it allows one to derive more precise ages (with typical errors of 20%, Miglio et al. 2021; or as low as 10%, Montalbán et al. 2021), while the latter needs to be calibrated on stars with known ages. The quantities used to infer stellar ages can be, for example, stellar rotation, magnetic activity, photospheric abundances (see, e.g. Soderblom 2010) and rocky exoplanet composition (Weeks et al. 2024). Several abundances and abundance ratios have been used in the literature to estimate stellar ages, such as, lithium abundance – A(Li) – (e.g. Jeffries et al. 2023), [C/N] (e.g. Masseron & Gilmore 2015; Casali et al. 2019; Spoo et al. 2022), and the abundance ratio between a slow neutron-capture element and an α element, [s/α] (e.g. Feltzing et al. 2017; Delgado Mena et al. 2019; Casali et al. 2020, 2023; Tautvaišienė et al. 2021; Moya et al. 2022; Viscasillas Vázquez et al. 2022; Berger et al. 2022; Ratcliffe et al. 2024; Shejeelammal et al. 2024). Some of these chemical age indicators depend on stellar evolution, such as A(Li) or [C/N], while others are the result of Galactic chemical evolution, such as [α/Fe] or [s/α]. The calibration process is usually based on the measurement of stellar properties in member stars of open star clusters (OCs), whose ages are known from the isochrone fitting of the cluster sequence. In this way, the chemical clocks can be calibrated and the calibrated relations can be applied to field stars. In particular, chemical clocks based on [s/α] have been widely used to derive ‘statistical’ ages of the Galactic stellar populations (e.g. Casali et al. 2020; Manea et al. 2023; Ratcliffe et al. 2024; Boulet 2024).

The production of s-process material comes from rotating massive stars (e.g., Pignatari et al. 2010; Frischknecht et al. 2012, 2016; Limongi & Chieffi 2018) and asymptotic giant branch (AGB) stars (e.g., Gallino et al. 1998; Lugaro et al. 2003; Cristallo et al. 2015; Karakas & Lugaro 2016; Busso et al. 2021). In particular, rotating massive stars are strongly responsible for the production of the first peak s-process elements (Y, Sr, Zr) via the weak s-process; while the main and the strong s-processes are due to low- and intermediate-mass stars (LIMSs) during their AGB phase. The main and the strong s-processes are mainly responsible for the production of elements belonging to the second (Ba, La, Ce) and the third (Pb, Au, Bi) s-process peak, respectively. Because of the longer lifetimes of low-mass stars one should expect an increasing trend of the [s/H] with decreasing stellar ages, particularly for the elements belonging to the second and third s-process peak. Such a trend is observed for the [s/Fe] ratios as well. The increasing trend towards younger ages was first observed by D’Orazi et al. (2009) for the abundance of Ba in a large sample of OCs and then confirmed by a large number of works (e.g. Maiorca et al. 2011, 2012; Jacobson & Friel 2013; Mishenina et al. 2015; Casali et al. 2020; Baratella et al. 2021; Viscasillas Vázquez et al. 2022; Sheminova et al. 2024). Not only Ba but also other neutron-capture elements belonging to both the first and second s-process peaks show a similar increase (Magrini et al. 2018; Frasca et al. 2019, but see also Yong et al. 2012; D’Orazi et al. 2017). The same reasoning is applied to the chemical clocks [s/α]. Since α-elements are mainly produced by massive stars on short timescales while s- process elements are released to the interstellar medium (ISM) at later times, their ratio should be higher at young ages. If such a trend with age is tight, then these abundance ratios can be used to infer stellar ages by means of empirical relations.

The discussion on the universality of the relationship between age and [s/α] derived from the stars in the solar neighbourhood and its dependence on the Galactic stellar populations – thin or thick disc –, on the position across the Galactic disc, and on metallicity remained open for a long time (see, e.g. Feltzing et al. 2017; Delgado Mena et al. 2019). The works based on the large sample of open clusters observed in the Gaia - ESO survey – spanning a wide range in age and Galactocentric distances – confirmed a strong radial variation along the disc and a parameterisable dependence on metallicity (see Casali et al. 2020; Viscasillas Vázquez et al. 2022). A first attempt to understand the variations in the age-[s/α] relationships along the Galactic disc was made in Magrini et al. (2021) by taking into account the inside-out formation of the disc and the metallicity dependence of the yields, particularly for the s-process (Cristallo et al. 2011; Vescovi 2021, 2023). To investigate the different hypotheses, it is important to rely on a self-consistent chemical evolution model.

From a modelling perspective, chemical evolution has encountered difficulties in reproducing the observed rise in the [s/Fe] and [s/α] ratios at very young stellar ages. Models frequently show a flat trend in younger stars, which could result from a balanced delay in enrichment timing between second- peak s-process elements from AGB stars and Fe from Type Ia supernovae (SNe; affecting [s/Fe] ratios) and/or a reduced production of s-process elements during the most recent billion years of chemical evolution (impacting both [s/Fe] and [s/α] ratios). Different solutions have been proposed in the literature to reconcile the model predictions with the observed trend, such as increasing Ba production at higher metallicities (Ratcliffe et al. 2024) or in lower-mass AGB stars (D’Orazi et al. 2009; Maiorca et al. 2012). These solutions represent promising progress in understanding the production and distribution of s-process elements. However, they often depend on comparisons with a single element, Ba, and focus solely on the solar neighbourhood. In this study, we test some of these existing solutions along with new approaches, examining both the production in AGB stars and in rotating massive stars. We compare our models with both Ba abundances (as representatives of the second peak s-process elements) and Y abundances (representing the first peak) across three different Galactocentric regions. We adopt a well-tested chemical evolution model that includes state-of-the-art nucleosynthesis calculations for neutron-capture elements, including contributions from AGB and rotating massive stars for the s- process as well as neutron star mergers and magneto-rotational SNe for the r-process. In fact, it is important to note that most elements typically classified as s- or r-process elements actually have dual (or even more) production pathways. Specifically, the fractions of Y and Ba produced by the s-process at solar metallic- ities are 78% and 89%, respectively (Prantzos et al. 2020), with the remaining fractions primarily accountable to the r-process (and only negligible or no contribution from the p-process). Here, our model results are compared to the open clusters set of the last data release of the Gaia-ESO survey (Randich et al. 2022). Our baseline model has been shown to provide a satisfactory fit to the standard abundance patterns of [El/Fe] versus [Fe/H] for a wide range of neutron-capture elements in different Galactocentric regions observed by the Gaia-ESO survey (see Molero et al. 2023).

The paper is structured as follows: in Section 2, we describe the Gaia-ESO OC sample. In Section 3, we present the chemical evolution framework, including the infall and the nucleosynthesis prescriptions. In Section 4, we present results of models with different prescriptions for infall, AGB, and rotating massive stars and their impact on the observed chemical clocks. Finally, in Sections 5 and 6, we address the limitations of the model, suggest potential improvements, and provide our conclusions.

2 Observational data

The Gaia-ESO survey (Randich et al. 2022; Gilmore et al. 2022) is a large public spectroscopic survey. It observed for 340 nights at the Very Large Telescope (VLT) from the end of 2011 to 2018, gathering ~190 000 spectra, for nearly 115 000 targets belonging to all the main Galactic populations. Gaia- ESO observations were carried on at two different resolving powers, R: the medium-resolution sample was observed with GIRAFFE at R~20 000 and the high-resolution sample with UVES at R~47 000. In particular, high resolution spectra, covering a wide spectral range from 480.0 nm to 680.0 nm (U580) or from 420.0 nm to 620.0 nm (U520), provide abundances for about 30 different ions and including numerous neutron-capture elements, such as Y, Zr, La, Ce, Ba, Eu, Nd, Pr, and Sm, discussed in Molero et al. (2023), for example. Gaia-ESO dedicated about 40% of its observation time to the population of open star clusters, used both to determine their specific properties and as calibrators (Pancino et al. 2017; Bragaglia et al. 2022). The sample clusters covered a wide range in age, Galactocentric distance, mass and metallicity (see Randich et al. 2022). Stellar parameters and chemical abundances from Gaia-ESO have been combined with homogeneous ages and distances from Gaia (e.g. Cantat-Gaudin et al. 2020).

In the present work, we adopt the same sample of open clusters used in Magrini et al. (2023b) and in Palla et al. (2024). In this sample of 62 OCs, only clusters older than 100 Myr were considered, to avoid observational biases in the chemical study of young stars. For the distribution in age and distances of this sample clusters, we refer to Viscasillas Vázquez et al. (2022, see their Fig. 1), in which the membership selection is also described. In addition, following the discussion in Magrini et al. (2023b) and in Palla et al. (2024), in each clusters we made a further selection of member stars based on their stellar parameters. From an evaluation of internal trends between abundances and stellar parameters, present in different spectroscopic surveys (see Magrini et al. 2023b, for a discussion), we consider only stars with log 𝑔 > 2.5 and microturbulent velocity ξ < 1.8 km⁻¹ to compute the mean cluster abundances.

We complemented our sample of OCs, with a sample of field stars from the same database, selected as in Viscasillas Vázquez et al. (2022). Due to the wide metallicity range covered by the field stars, we use them to study the evolution in the [El/Fe] versus [Fe/H] planes. However, owing to the very large uncertainties in determining their ages, we do not use them in the age versus abundance ratios diagrams.

The most studied chemical clocks in the literature are [Y/Mg] and [Y/Al] (Tucci Maia et al. 2016; Nissen et al. 2017; Spina et al. 2018; Tautvaišienė et al. 2021; Magrini et al. 2021; Berger et al. 2022; Shejeelammal et al. 2024„ but see also Delgado Mena et al. 2019; Jofré et al. 2020; Ratcliffe et al. 2024). For this work, we selected two s-process elements: Y for the first peak elements and Ba for the second peak. To compute the abundance ratios [s/α] we used Si as representative of the class of the α-elements, due to the larger sampling of OCs relative to Al and the more precise nucleosynthesis prescriptions from massive stars relative to Mg (see discussion in, e.g., Romano et al. 2010; Palla et al. 2022; Jost et al. 2025).

3 The chemical evolution models

In this section, we present the chemical evolution models adopted to study the variations in the [s/α] versus age relations along the Galactic disc. The models are as follows:

• The revised two-infall model proposed by Palla et al. (2020) (see also Spitoni et al. 2019 and Spitoni et al. 2021) and previously adopted in Molero et al. (2023) to study the distribution of the neutron-capture elements along the Milky Way (MW) disc.

• The three-infall model from Spitoni et al. (2023). In particular, here we use the version proposed by Palla et al. (2024), extended to the whole disc.

The revised two-infall model is a variation of the classical two- infall model of Chiappini et al. (1997) (see also Chiappini et al. 2001) designed to explain the dichotomy in α-element abundances observed both in the solar vicinity (e.g. Hayden et al. 2014; Recio-Blanco et al. 2014; Mikolaitis et al. 2017) and at various radii (e.g. Hayden et al. 2015). The model suggests that the initial primordial gas infall event formed the chemically thick disc (corresponding to the high-α sequence), while the second infall event, occurring approximately 3–4 Gyr years later, formed the chemically thin disc (the low-α sequence). It is important to note that the two-infall model used here does not aim to distinguish the thick and thin disc populations based on geometric or kinematic criteria (see Kawata & Chiappini 2016 for further discussion). The three-infall model, on the other hand, is an extension of the two-infall one, designed to replicate the low-α sequence through two distinct gas infall episodes (the most recent of which began approximately 2.7 Gyr ago). This approach aims to account for the recent chemical impoverishment in metallicity with low [α/Fe] values identified by Gaia Collaboration (2023) (reported for the first time by Spina et al. 2017), as well as the recent increase in star formation activity described by Ruiz-Lara et al. (2020) (see also Isern 2019; Mor et al. 2019).

3.1 Basic equations

The basic equations that describe the evolution of the fraction of gas mass in the form of a generic chemical element i, as a function of time, t, and of the Galactocentric distance, R, are $${\dot{G}}_i(R,t)=-\psi (R,t)X_i(R,t)+{\dot{G}}_{i,inf}(R,t)+R_i(R,t),$$(1)

where X_i(R, t) represents the abundance by mass of a given element i. The term ψ(R, t)X_i(R, t) is the rate at which chemical elements are subtracted by the ISM to be included in stars. The star formation rate (SFR) is parameterised following the Schmidt-Kennicutt law (Kennicutt 1998) as $$\psi (R,t)\propto \nu {\sigma }_{gas}{(R,t)}^k,$$(2)

with a law-index k = 1.5 and σ_gas and ν being the surface gas density and the star formation efficiency, respectively.

The term R_i(R, t) is the rate of restitution of matter from the stars with different masses into the ISM in the form of the element i. It takes into account the nucleosynthesis from a variety of stars and phenomena, including stellar winds, SN explosions of all kinds, novae and neutron-star mergers. For the complete expression we refer to Matteucci (2012). Here, we discuss in more detail the nucleosynthesis from the different sources in Section 3.2.

For the two-infall model, the accretion term on the right-hand side of Equation (1) is computed as: $${\dot{G}}_{i,inf}(R,t)=A(R)X_{i,inf}{e}^{-\frac{t}{{\tau }_1}}+\theta \left(t-t_{max}\right)B(R)X_{i,inf}{e}^{\frac{t-t_{max}}{{\tau }_2}},$$(3)

where X_i,inf is the composition of the infalling gas, here assumed to be primordial for both the infall episodes. τ₁ ≃ 1 Gyr is the infall timescale for the first infall event, assumed to be constant for all radii. During the first infall event, the star formation efficiency is set to ν = 2 Gyr⁻¹, fixed at all Galac- tocentric distances as well. On the contrary, for the second gas infall event, the timescale τ₂ is assumed to vary with the radius according to the inside-out scenario (see e.g., Matteucci & Francois 1989; Romano et al. 2000; Chiappini et al. 2001), with longer timescales towards larger distances. The star formation efficiency during the second infall episode is assumed to be a function of the Galactocentric distances as well, with ν ≃ 1 Gyr⁻¹ at R_GC = 8 kpc. In order to correctly reproduce the slope of the abundance gradients of Fe and α-elements, the star formation efficiency assumes larger values towards the inner part of the disc (see e.g., Grisoni et al. 2018). t_max ≃ 3.25 Gyr is the time for the maximum infall on the disc and it corresponds to the start of the second infall episode. θ is the Heavyside step function. The parameters A(R) and B(R) are fixed to reproduce the surface mass density of the MW disc at the present time in the solar neighbourhood as provided by McKee et al. (2015), equal to 47.1 ± 3.4 M_⊙pc⁻². At different Galactocentric radii, as discussed in Palla et al. (2020), the surface mass densities of the chemically thick and thin discs are assumed to follow an exponential profile as $${\Sigma }_{thick;thin}(R)={\Sigma }_{0,thick;thin}{e}^{-R/R_{d,thick;khin}},$$(4)

where R_d;thick = 2.3 kpc and R_d;thin = 3.5 kpc are the disc scale lengths for the high-α and for the low-α disc, respectively.

For the three-infall case, the accretion term in Equation (1) has the following form: $$\begin{array}{cc}{\dot{G}}_{i,\text{\hspace{0.17em}inf\hspace{0.17em}}}(R,t)=\hfill & A(R)X_{i,inf}{e}^{-\frac{t}{{\tau }_1}}+\theta \left(t-t_{max,1}\right)B(R)X_{i,inf}{e}^{\frac{t-t_{max,1}}{{\tau }_2}}\hfill \\ \hfill & +\theta \left(t-t_{max,2)}\right)C(R)X_{i,inf}{e}^{\frac{t-t_{max,2}}{{\tau }_3}},\hfill \end{array}$$(5)

where τ₃ = 1 Gyr is the timescale of the third infall event and t_max,2 = 11 Gyr is the Galactic time of the start of the third accretion event. The coefficient C(R) is set in order to fit the present-day total surface density of the third accretion phase, Σ₃. Since the three-infall model splits the low-α sequence in two gas accretion sequences, the sum between Σ₂ and Σ₃ is equal to the density profile of the low-α sequence of the two-infall model. In particular, the ratio Σ₂/Σ₃ as a function of the Galactocen- tric distance as well as all the variables of the third-infall event are fixed following the prescriptions of the best model of Palla et al. (2024). All the other parameters for the first and the second infall phases are as in the two-infall model case. In Spitoni et al. (2023), the three-infall model has been first introduced to predict the recent chemical impoverishment characteristic of the young, massive stellar populations observed in Gaia DR3. Similarly, in Palla et al. (2024), the model was employed to account for the properties of young OCs of the Gaia-ESO survey.

Comparisons of the evolution of some important quantities predicted by the three-infall model to present-day observations are reported in Figure 1. The predicted SFR, surface densities of stars and gas are computed in the solar neighbourhood and compared with present-day estimates as suggested by Prantzos et al. (2018). Rates of Type Ia, Type II SNe, MR-SNe and MNS are averaged over the whole disc and compared with the observational values estimated of Cappellaro et al. (1999) and Abbott et al. (2021).

Fig. 1 Time evolution of the rate of Type Ia SNe, Type II SNe, MR-SNe, and MNS (upper-left panel); SFR (upper-right panel); surface densities of gas (lower-left panel); and stars (lower-right panel) as predicted by the three-infall model. Predictions of the model are compared to present day values from Cappellaro et al. (1999) (for SNe rates), from Abbott et al. (2021) (for MNS rate) and from Prantzos et al. (2018) (for SFR and surface gas and stars densities).

3.2 Nucleosynthesis prescriptions

In this work, we adopt yield prescriptions similar to those used by Molero et al. (2023), where a detailed description can be found. Here, we provide a brief recap.

Yields for low- and intermediate-mass stars (LIMSs; with initial masses 1.0 ≤ M/M_⊙ ≤ 8.0) are taken from the FRUITY database (Cristallo et al. 2009, 2011, 2015). We adopt the non- rotational set which provides yield grids for progenitors with 8 initial stellar masses, from 1.3 M_⊙ to 6.0 M_⊙, and 12 values of metallicities from Z = 4.8 × 10⁻⁵ to Z = 2.0 × 10⁻². Yields for the 8 M_⊙ have been obtained extrapolating the FRUITY yields for all metallicity values to reduce the gap between LIMSs and massive stars yields. Both the main and the strong s-process component are produced in AGB stars. However, as underlined in Magrini et al. (2021) (see Vescovi et al. 2021), isotopic ratio measurements in presolar SiC grains showed that the neutron density in FRUITY models is likely overestimated (Liu et al. 2018). Consequently, s-process yields from FRUITY AGB models are often reduced by a factor of two or more in chemical evolution studies (e.g., Rizzuti et al. 2019, 2021; Molero et al. 2024, see Cescutti & Matteucci 2022 for a recent review). In this work, we relax this reduction across all metallicities unless stated otherwise.

Massive stars (with initial masses 8 ≤ M/M_⊙ ≤ 120) are assumed to explode either as normal core-collapse SNe (CC- SNe) or as magneto-rotational SNe (MR-SNe). For normal CC-SNe, yields are taken from Limongi & Chieffi (2018) set R. We run models with either a constant initial rotational velocity of 150 km s⁻¹ or a velocity distribution. About 20% of massive stars in the 10 – 25 M_⊙ range are assumed to explode as MR- SNe, representing a source of r-process material with yields from Nishimura et al. (2017) (model L0.75).

The second r-process site is represented by merging neutron stars (MNS). They are computed as systems of two neutron stars with masses of 1.4 M_⊙, originating from progenitors in the 9 – 50 M_⊙ range. Mergers follow the delay-time distribution of Simonetti et al. (2019) (see also Greggio et al. 2021 for a detailed discussion), and their rate is constrained to reproduce the latest estimation from Abbott et al. (2021). Yields of r-process material from MNS are obtained by scaling to solar the yield of Sr measured in the re-analysis of the spectra of the kilonova AT 2017gfo by Watson et al. (2019) (after having considered uncertainties in its derivation; see Molero et al. 2021 for details). The Eu yield obtained in this manner is consistent with the theoretical calculation of Korobkin et al. (2012) and with estimates from Matteucci et al. (2014).

It is, however, important to stress that some of the main uncertainties in the production of r-process elements in chemical evolution models are the adopted prescriptions for MR-SNe. While for MNS we can rely on observations for the yields (Watson et al. 2019), for the delay-time distributions (from observations of short gamma-ray bursts; e.g. Ghirlanda et al. 2016) and for the rate (Kalogera et al. 2004; Abbott et al. 2021), for MR-SNe we have to arbitrarily choose (i) the fraction of normal CC-SNe which can die as a MR-SNe, (ii) the mass range and (iii) the yields (since different authors often obtain different results). In Molero et al. (2023), we show some possibilities in which these free parameters can be fixed in a self-consistent manner. Nevertheless, ambiguities persist, and the possibility of a substantial level of degeneracy remains a significant consideration. In particular, while for normal CC-SNe yield grids for a wide number of different progenitors exist, this is not the case for MR-SNe for which the prescription for one single stellar mass (e.g., 35 M_⊙) is interpolated to all the range of masses considered. Therefore, in our chemical evolution model, MR-SNe are included solely as sources of r-process material, excluding their contribution to Fe and other elements.

Finally, for Type Ia SNe, we assume the single-degenerate scenario for the progenitors, with stellar yields from Iwamoto et al. (1999) (model W7).

As a sanity check for the nucleosynthesis prescriptions of r-process, Fe and α-elements, we show the evolution of the Si, Y, Ba and Eu abundances in the typical [X/Fe] versus [Fe/H] diagnostic diagrams (see Figure 2) predicted by the three-infall model for the three different Galactocentric distances considered. From a nucleosynthesis point of view, here prescriptions are identical to the one adopted in Molero et al. (2023), with reduced s-process AGB production. Note that, following Palla et al. (2024) (see also Spitoni et al. 2024), we accounted for observational uncertainties in the model predicted abundances by adding a random error at each time step, t: $${[X/H]}_{\mathit{\text{new}}}(t)={[X/H]}_{\mathit{\text{old}}}(t)\pm N\left(\left[0,{\sigma }_{[X/H]}\right]\right),$$(6)

where N is a random generator with a normal distribution, and σ_[X/H] represents the average observational error of the open cluster sample. The model well agrees with the overall abundance patterns, in particular in the last 3 Gyr of evolution, as shown by the synthetic model results. Only in the case of Eu in the outer regions the model underproduces the OCs data in the age range of interest, reproducing their pattern too early in the Galactic evolution. In fact, the synthetic model, plotted only for Age ≤ 3 Gyr, falls below the observed OCs data. On the other hand, the track of the chemical evolution model, corresponding to Age ≤ 10.45 Gyr, aligns with the OCs trend. As a consequence, although the model may initially appear to fit the OCs trend well, it actually underestimates the data within the corresponding age range. This discrepancy occurs however exclusively for [Eu/Fe] at high metallicities and it is confined to the outer regions. The lower metallicity data described by the field stars is correctly reproduced for all elements and in all the region of interest, except for Y in the inner regions where data are overestimated. In the case of Eu, in particular, the agreement with the trends observed for field stars in the low metallicity regime ensures that the prescriptions for the MR-SNe nucleosynthesis are correctly fine-tuned.

Fig. 2 [X/Fe] versus [Fe/H] abundance patterns for Si, Y, Ba, and Eu as predicted by the three-infall model at RGC = 6, 8, 12 kpc compared to the sample of OCs (colored dots) and field (grey dots) stars in the inner, solar, and outer regions, respectively. Results of the synthetic model, as well as OCs, are plotted only for Age ≤3 Gyr. In the lower right panel, for the prediction of the [Eu/Fe] versus [Fe/H] in the outer regions, the filled magenta circles indicate the [Eu/Fe] at a given age.

4 Results

In this section, we use stellar ages as indicators of the passage of time, instead of look-back time. Thus, for relationships between ages and abundance ratios, positive slopes indicate relationships that decrease with the passing of time, while negative slopes indicate relationships that increase. Here we show results from our chemical evolution model to be compared to the [s/H], [s/Fe] and [s/α] versus age trends observed in the OCs dataset at different Galactocentric distances. The focus is mainly on the chemical clock trends and on how these are influenced by different infall events and nucleosynthesis prescriptions from AGB and massive stars.

4.1 The effect of primordial gas inflows

As outlined in Sect. 1, different slopes for the [X/Y] versus age relations are indications of a different contribution from different nucleosynthesis sources of the two elements X and Y. It is well-known that, while the observed relations of [α/Fe] versus age have a positive slope, the one of the [s/Fe] has negative ones. The positive slopes of the α-elements represent their production on a shorter timescale with respect to Fe. Type II SNe are indeed producing α-elements on short timescales (~10⁻² Gyr) as well as Fe. However, since the bulk of Fe comes from Type Ia SNe over longer timescales, the observed [α/Fe] ratio will decrease towards younger ages, producing a positive slope. On the other hand, the observed negative slope of the [s/Fe] versus age relations is attributed to the delayed contribution to the s-process elements production from AGB stars. Nevertheless, although AGB stars are the main producers of s- process elements, rotating massive stars play an important role in interpreting the [s/Fe] versus age relations. In Figure 3, we show the observed [s/Fe] versus age relations for Y, as representative of the first s-process peak elements, and for Ba, as representative of the second peak. The relation is plotted for the three different Galactocentric regions of interest together with a logarithmic fit of the OC sample and the results from the two- and three-infall models described in the previous sections (Model 1 and 2, respectively; see Table 1 for reference). For both Y and Ba, it is possible to see the increase of the [s/Fe] ratio towards younger ages in the OC set discussed above. However, it is interesting to note that, in the case of Y, the rise shown by the logarithmic fit in the inner and outer regions is not significantly pronounced. Notably, the abundances observed in the youngest OCs remain within the uncertainty range of those of the older ones, resulting in an overall flat trend. This long plateau is most probably due to the production of Y by rotating massive stars which happens on the same timescales as that of Fe. The flat trend, in fact, is less evident in the case of Ba, for which the contribution from massive stars is expected to be weaker. The logarithmic fit in the solar zone is computed excluding the young (Age ≃ 0.2 Gyr) OC NGC 6709, which exhibits a high abundance of s-process elements. NGC 6709 is enriched in Y, Ba, and Nd but shows no significant overabundance in Fe ([Fe/H] = −0.025 dex) or α-elements such as Mg (e.g. [Mg/H] = −0.038 dex). However, due to the presence of only two confirmed members and the significant uncertainties in the measured s-process element abundances, NGC 6709 is excluded from the computation of the logarithmic fit, which, as a consequence, displays an overall flat trend. Interestingly, for Ba in particular, the trend appears to increase towards younger ages up to Age ≃ 1 Gyr.

The two- and three-infall models shown in Figure 3 both fail to reproduce the observed growth in the [s/Fe] versus age relations. Instead of producing an increasing trend towards younger ages, the models produce either a decreasing or a flat trend. The decrease is more pronounced in the case of Ba than Y, indicative of faster Y enrichment due to the contribution of massive stars. The flat trend predicted by the chemical evolution model for Y is more consistent with observations in the outer regions compared to the solar and inner zones. Similarly, for Ba, the model’s trend in the outer region aligns better with observed data. However, while the observations reveal a clear increasing trend, none of the model predictions capture this feature in any of the regions of interest. The flat trend predicted for Y is a consequence of the similar production mechanism for Y and Fe, both in terms of timescales and overall amounts. Both elements are produced on comparable timescales, and their ratio remains relatively stable over time. Notably, although metallicity increases from the outer to the inner regions, the [Y/Fe] ratio predicted by the models generally exhibits a flat or only slightly decreasing trend towards younger ages in the three different regions. In contrast, the predicted [Ba/Fe] ratio consistently decreases, with the steepest decline observed towards the inner regions, possibly reflecting a reduced contribution to the Ba production at higher metallicities. The trends predicted by the two- and three-infall models are similar, though the three-infall model slightly mitigates the decline of the [Ba/Fe] at very young ages. The third infall event brings in fresh, metal-poor gas that dilutes the existing metal-rich ISM. As a consequence, it also introduces a new episode of star formation that may alter the chemical evolution dynamics. However, although this newly formed stars have the potential to slow the decrease in [Ba/Fe], it is insufficient to reverse the overall trend. The additional gas inflow impacts both s-process elements and Fe similarly, thereby failing to create the observed increasing trend in the [s/Fe] ratio. Thus, while the third infall event introduces new complexities to the Galaxy chemical evolution, it does not provide a solution to the discrepancy between model predictions and observed [s/Fe] trends. An inversion of the predicted trend might therefore be obtained only from a nucleosynthetic point of view, by modifying either the timescales or the amount (or both) of production of s-process material.

The three-infall model has nevertheless been demonstrated by Palla et al. (2024) to successfully reproduce the recent dilution at intermediate age (1 < Age/Gyr < 3) in Fe and α-elements observed in the same set of OCs used in this work. Specifically, Palla et al. (2024) attribute the observed dilution in Fe to a third, recent gas accretion event of primordial (or slightly enriched) material. This approach allows the authors to accurately replicate the age-metallicity relation observed in the current OC dataset across the entire disc, as well as the present-day metallicity gradient. A similar dilution is observed in the s-process elements Y and Ba, as illustrated in Figure 4, in which we show both the observed and predicted trends for the [s/H] versus age relations.

When applied to these relations, the extended three-infall model from Palla et al. (2024) (Model 2) achieves a reasonable agreement with the data, maintaining a generally increasing trend. Unlike the two-infall scenario, this model successfully reproduces the dilution observed at Age ≃ 2 Gyr in the solar and inner regions. However, in some cases, the predicted relations are higher than the observations. In particular, Model 2 nicely fits the [Y/H] in the observed uncertainty ranges, but fails to reproduce the [Ba/H] in the solar regions. In the inner region, while the model overestimates the general trend indicated by the logarithmic fit, it still falls within the data uncertainty ranges.

The overproduction of s-process elements is a well-known issue in Galactic chemical evolution modelling and is often attributed to an excessive production of s-process material from the AGB stars yields set of the FRUITY database. Typically, this issue is addressed by applying a reduction factor to the FRUITY yields across all elements, set independently of the progenitor star (as first done by Rizzuti et al. 2019). This approach implicitly assumes that the overproduction is constant across mass, metallicity, and elements. However, from Figure 4 it seems that this is not the case. Elements belonging to the first s-process peak, such as Y, are generally less overestimated than those belonging to the second s-process peak, such as Ba. Furthermore, the different discrepancies between our model and the observed trend across different disc regions, show that the necessary reduction factor should not be uniform across the disc, indicating that it should vary with metallicity. Age-related differences also emerge; for instance, in the solar region, the model matches the [Ba/H] of a few OCs at Age ≃ 1 Gyr, while it overestimates this ratio for both older and younger ages.

These discrepancies propagate into the [s/Fe] versus age relations, as previously discussed, and ultimately affect the interpretation of chemical clocks, complicating the understanding of the observed trends. The observed pattern of the [s/Si] versus age relations shown in Figure 5, shares many similarities with that of [s/Fe]. The logarithmic fit shows a general increase towards younger ages, with a steeper growth for Ba compared to Y. Similar to the [s/Fe], the solar region exhibits an overall flat trend (with an increase only up to Age ≃ 1 Gyr) rather than the increasing trend observed in other regions. We remind that the OC NGC 6709 is excluded from the computation of the fit. The models, both in the case of the two- and three-infall scenarios, struggle to reproduce the observed patterns. Similar to the issues seen with the [s/Fe] ratios, these models exhibit an inverse trend compared to the observed one for [Ba/Si] in the inner and local regions. The trend is better reproduced in the outer region, but the models still show an overproduction of Ba, further highlighting the challenges in accurately modelling the chemical evolution of s-process elements. The models agree only with the data of the [Y/Si] in the outer regions. These discrepancies suggest that the current models may not fully capture the complexities of the s-process nucleosynthesis, particularly in terms of the contributions from different stellar populations.

Fig. 3 [s/Fe] versus Age observed trends for Y (upper panel) and Ba (lower panel) divided into the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. (See Table 1 for reference.)

Fig. 4 [s/H] versus age observed trends for Y (upper panel) and Ba (lower panel) divided in the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. (See Table 1 for reference.)

4.2 The influence of asymptotic giant branch stars

With the current set of nucleosynthesis prescriptions, the issues in reproducing the observed [s/Fe] versus age and chemical clock trends are twofold: i) the models tend to overestimate the observed patterns in some regions of interest, and ii) they almost always predict an inverted or flat trend compared to the observed increasing ones.

Concerning the first point, predictions from our model appear to be more compatible with the observed patterns in the outer region than in the solar and most inner regions. This observation aligns with findings by Casali et al. (2020), who noted that the content of neutron-capture elements belonging to the first s-process peak, is lower than expected from chemical evolution models. As a result, they found that the [Y/α] ratios in clusters located in the inner region are lower than in clusters of the same age in the solar vicinity. Their direct conclusion was that stellar dating relations between abundance ratios and ages, based on samples of stars in the solar neighbourhood, cannot be universally applied across different regions of the disc (see also Viscasillas Vázquez et al. 2022; Casali et al. 2023; Ratcliffe et al. 2024). According to the authors, one should expect less s- process elements to be produced at high metallicity and tested a set of empirical yields from AGB stars, in which the supersolar metallicity yields were depressed by a factor of ten. This adjustment allowed them to successfully reproduce the chemical abundances of the young OCs located at R_GC ≃ 6 kpc. However, when comparing this approach to the more extensive dataset used in this work, it becomes clear that a set of reduced AGB yields at super-solar metallicity does not provide a satisfactory fit across the board. This is illustrated by the dashed line in Figure 6 (Model 3; see Table 2), where following Casali et al. (2020), we reduced the yields of AGB stars by a factor of ten for Z ≥ 1.4 × 10⁻² (namely, the solar metallicity as computed by Asplund et al. 2009). Although this adjustment improves the model’s match with OCs in the inner regions for ages between approximately 1.0 and 2.0 Gyr, it accentuates the inverted trend compared to the observed data. Moreover, the reduction in AGB yields does not impact the model’s result in the solar vicinity and in the outer regions, as expected. This suggests that while reducing the AGB yields at high metallicities may help to some extent in specific cases, it does not provide a comprehensive solution and may introduce new issues in reproducing the observed trends across different Galactic regions.

The second issue, in our opinion, is the most crucial. The inability of the models to capture the correct observed directional trend suggests that fundamental aspects of the chemical evolution, particularly the timing and contribution of s-process elements relative to Fe and α-elements, may be inaccurately modelled. This highlights the need for a reassessment of the nucleosynthesis yields. This aspect has emerged also from the recent work of Ratcliffe et al. (2024), where the non-universality of the chemical clock [Ba/Mg] versus age relation is studied. Once radial migration has been included, the chemical evolution model adopted by Ratcliffe et al. (2024) is able to capture the overall trend of the [Ba/Mg] radial gradient as a function of time. However, the same model fails in reproducing the [Ba/Mg] versus age relation since, similar as our case, it predicts a significantly decrease of the [Ba/Mg] abundance instead of the observed steep increase with time. According to the authors, a possible explanation for this decrease is that the amount of Ba produced in the adopted FRUITY models of high metallicity AGB stars is not high enough. A promising solution proposed by Ratcliffe et al. (2024) was to replace the AGB yields for Z > 0.01 with those from models at Z = 0.01, which, in the case of Ba, are approximately two times larger than those at higher metallicities. This small adjustment to the high metallicity tail of the AGB yields, allowed the authors to mostly resolve the discrepancies observed in the first place. However, here we note that this solution is not universally effective. While it corrects the Ba yields, it fails to address issues with Y. The yields for Y at Z = 0.01 are actually lower – by about 1.3 times, depending on the progenitor mass – than those at higher metallicities. This reduction would create a similar discrepancy for Y as that seen in the Model 3 in the [Y/Si] versus age relation.

It is important to highlight the complexity of finding a universal solution in terms of AGB yields across different metallicities. The production of s-process elements in AGB stars is influenced by several factors, including the quantity of seed nuclei available and the number of thermal pulses the star undergoes during the AGB phase which, in turn, depend on the metallicity. Moreover, these dependencies do not behave uniformly across all s-process elements. As expected, elements belonging to different s-process peaks respond differently to changes in metallicity which, together with an inside-out scenario of disc formation and a radial variable SF efficiency, makes it challenging to create a single model that accurately predicts the behaviour of the chemical clocks across the disc.

Nevertheless, despite variations in efficiency at different Galactocentric distances, the general increasing trend in the [s/α] towards younger ages is observed consistently across different regions of the Galactic disc. These similarities in the timing of this trend imply that modifications of the s-process yields from AGB stars may depend more on the initial stellar mass than on the metallicity of the progenitor star. In the context of the FRUITY model database, the lowest progenitor mass included is a star with an initial mass of 1.3 M_⊙. The lifetime of such a star (assuming solar chemical composition) is approximately 5.19 Gyr, meaning it would have enriched the ISM with s-process elements around 8 Gyr ago. Therefore, increasing the yields from a 1.3 M_⊙ star would cause the model to predict higher [s/Si] ratios than observed at each age of interest, implying that stars with even lower masses may be responsible for the observed increases in [s/α] ratios (see also D’Orazi et al. 2009). In Figure 6, we show results of models in which we artificially include yields from a 1.1 M_⊙ star. Yields are assumed to be 10–20 times greater than those of a 1.3 M_⊙ (Models 4 and 5). Adding yields from a lower-mass progenitor has different effect on Y and Ba. Specifically, since the Y production is less sensitive to the contribution of AGB stars, the inclusion of the 1.1 M_⊙ progenitor does not significantly impact the model, especially in the outer regions of the Galaxy, where different models produce nearly identical results. In the case of the [Ba/Si], on the other hand, the inclusion of the lower-mass progenitor does have a notable effect, particularly in the solar region, where the models are able to reverse the declining trend in Ba abundance at younger ages. However, while this adjustment leads to an increase in Ba at younger ages, it still falls short of fully invert the modelled trend in accordance to the observational data. This, however, raises questions about whether lower-mass stars (below 1.3 M_⊙) contribute more significantly to the production of s-process elements than currently accounted for in models. According to D’Orazi et al. (2009) and Maiorca et al. (2012), this might be possible in a scenario where the efficiency of the extra-mixing processes producing the neutron source 13C is anti-correlated with the initial mass of the star. It is also possible that the intermediateprocess (i-process) nucleosynthesis takes place during the early AGB phase of low-mass (1.0 M_⊙) stars. However, as shown in Choplin et al. (2024) (see also Choplin et al. 2021, 2022), models are in favor of an i-process operating in AGB stars up only to metallicity [Fe/H] ≃ −1 dex.

Fig. 5 [s/Si] versus age observed trends for Y (upper panel) and Ba (lower panel) divided in the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. (See Table 1 for reference.)

Fig. 6 [s/Si] versus age observed trends for Y (upper panel) and Ba (lower panel) divided into the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of different assumptions concerning the production of s-process material from AGB stars are compared to the OC sample. (See Table 2 for reference.)

4.3 The influence of rotating massive stars

In this section, we shift focus from AGB stars to rotating massive stars as a key nucleosynthetic site for the production of s-process elements in the model. While the previous section demonstrated the challenges in finding a universal model for the chemical clocks based on the s-process yields from AGB stars (due to variations in metallicity and progenitor mass), here we explore the contribution of rotating massive stars to the synthesis of s-process elements. Rotating massive stars play indeed a crucial role in the Galactic chemical evolution, particularly for s-process elements at lower metallicities (see e.g., Frischknecht et al. 2012, 2016; Limongi & Chieffi 2018; Prantzos et al. 2018; Rizzuti et al. 2019, 2021), and their initial rotation speeds influence the production of these elements.

In our reference model (Model 2), we adopted a fixed initial rotational velocity for massive stars of 150 km s⁻¹. Even though such a nucleosynthesis prescription allowed us in Molero et al. (2023) to correctly reproduce the majority of the abundance patterns of the neutron capture-elements observed by Gaia-ESO, the adoption of a distribution of rotational velocities is usually favoured. A more realistic approach involves adopting a distribution of rotational velocities rather than a single fixed value. Prantzos et al. (2018) firstly implemented the idea of such a distribution in the context of Galactic chemical evolution models, based on the massive star yields from Limongi & Chieffi (2018). This distribution accounts for different initial rotational velocities and adjusts them as a function of [Fe/H], to match the observed behaviour of primary nitrogen (14N) at low metal-licities and to avoid an overproduction of s-process elements at higher metallicities. The obtained distribution favors higher rotational velocities at lower [Fe/H] and lower or negligible velocities at solar and supersolar metallicity. This is in agreement with also other distributions in literature which where parametrised later (e.g. Romano et al. 2019; Rizzuti et al. 2019, 2021; Molero et al. 2024).

Here, we test two different distributions: DIS 3 from Molero et al. (2024) and the one from Prantzos et al. (2018) (Model 6 and 7, respectively; see Table 3 for reference). DIS 3, similar to the distribution adopted by Romano et al. (2019) for studying CNO isotopes, assumes that massive stars rotate at 150 km s⁻¹ for Z < 3.236 × 10^—3 ([Fe/H] ≲ −1.0) and rotation becomes negligible beyond the metallicity threshold of Z = 3.236 × 10⁻³. The results of such a distribution is reported by the dashed lines in Figure 7.

In the case of the [Y/Si], the absence of high rotation in massive stars reduces the production of Y. Si production is however also reduced, and since both Y and Si are reduced by similar amounts, the net effect is negligible in the [Y/Si] ratio across the Galactic disc. In fact, despite this reduction, the overall pattern of the chemical clock remains similar to models where massive stars consistently have higher rotational velocities, showing that even with reduced rotation, the general trend of the chemical clock [Y/Si] can still resemble higher-velocity models due to balanced effects on multiple elements. In the case of [Ba/Si], the impact of adopting a distribution of rotational velocities for massive stars differs from the results seen with [Y/Si]. In all the regions both Ba and Si production are reduced due to the absence of rotation. However, Si is reduced more than Ba, leading to a higher [Ba/Si] ratio compared to the model with constant faster rotation.

The trends predicted by models with a distribution of rotational velocities for massive stare might be, in general, highly dependent on how the model transitions between populations of rotating and non-rotating massive stars at a fixed metallicity. In Model 6, this transition occurs abruptly at a specific metallicity threshold. A more realistic approach would account for the coexistence of different populations of rotating massive stars across different metallicities, including a small fraction of highly rotating stars (e.g. with initial rotational velocities of 300 km s⁻¹) towards lower Z. This approach is shown by Model 7 (dotted lines in Figure 7, representing models in which the Prantzos et al. 2018 distribution is adopted). With this distribution, Model 7 yields similar [s/Si] results to Model 6. For the [Ba/Si], results are slightly under the one of Model 6, as expected, because of the lower percentage of stars with null rotational velocities. Differences are minor, but noticeable. Overall, none of the tested distributions successfully capture the increasing trend observed in chemical clocks towards younger ages. Additionally, in the case of [Ba/Si], the models tend to worsen the overestimation of the observed pattern. Concerning this last point, it appears clear that the general idea according to which reducing the rotational velocity of massive stars would help in reproducing the lower level of [s/α] observed in the inner region (due to reduced s-process material production), is challenged. The different contribution to different elements from different progenitor stars are more complex, and their effect are not linear and, as a consequence, it appears difficult to lead back the issue to a specific choice in the distribution or a specific progenitor star without performing a more extensive parametric study.

Fig. 7 [s/Si] versus age observed trends for Y (upper panel) and Ba (lower panel) divided in the three Galactocentric regions of interested. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region) and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. Predictions of the chemical evolution model in the case of different assumptions concerning the production of s-process material from massive stars are compared to the OC sample. See Table 3 for reference.

5 Discussion on the missing amount of Ba

In this study, we have investigated the influence of star formation history/gas infall events, AGB stars, and massive stars on the production and distribution of s-process elements across the MW disc, in the context of chemical clocks. We focused on the [Y/Si] and [Ba/Si] versus age relationships across different Galactic regions, which, as observed in the OC datasets from the Gaia-ESO survey used here, exhibit an increasing trend towards younger ages. This trend is more pronounced for second s-process peak elements compared to first peak elements. The discrepancies between our model predictions and observed data – particularly the overestimation of certain patterns and the prediction of inverted trends – highlight the challenges of developing a universal Galactic chemical evolution model for elements belonging to different s-process peaks.

We showed that neither the two-infall nor the three-infall model can accurately reproduce the increasing trend of the chemical clocks. However, the three-infall model more successfully reproduced the observed depression in the [s/H] versus age trends at ≃2 Gyr.

While AGB stars play a significant role in s-process nucleosynthesis, current models with state-of-the-art nucleosynthesis prescriptions struggled to accurately reproduce the observed chemical clock trends across the Galactic disc. Our model predictions aligned better with the observed patterns in the outer regions than in the solar and inner regions and, in general, with the [Y/Si] than the [Ba/Si]. Attempts to adjust AGB yields for super-solar metallicities, as originally proposed by Casali et al. (2020), provided some improvement in reproducing the observed data for certain OCs, particularly those in the inner disc aged 1 to 2 Gyr. However, this approach did not yield a satisfactory fit across all regions, and instead introduced new discrepancies and failed to provide a universal solution for the entire disc. Suggestions from Ratcliffe et al. (2024) are also discussed, highlighting the difficulties that the models face in simultaneously reproducing elements belonging to different s-process peaks. A potential solution to reproduce the increasing trend in chemical clocks is the inclusion of lower-mass stars (below 1.3 M_⊙) as contributors to s-process element production, as first suggested by D’Orazi et al. 2009. However, though the s-process yield from such low-mass stars might be substantial, it is likely unrealistic. The potential contribution from the i-process was briefly considered and has been tentatively ruled out given its current inability to produce neutron-capture elements above [Fe/H] ≳ −1 dex (Choplin et al. 2024).

We recognized the essential role of massive stars in the production of s-process elements, especially at lower metallicities. The adoption of a fixed initial rotational velocity for massive stars in our starting model (Model 2) proved insufficient for capturing the complexities of the chemical clocks. Testing the distributions of rotational velocities, as first introduced by Prantzos et al. (2018), revealed unexpected impacts on the chemical clock patterns, particularly in the case of the [Ba/Si] versus age. A distribution of initial rotational velocities for massive stars can lead to an overproduction of [Ba/Si], even when null rotational velocities are included, due to the reduced contribution to Si production. The net effect of adopting a rotational velocity distribution strongly depends on the parameter choices for the distribution itself. Overall, adopting a distribution of rotational velocities either worsens the agreement with observations or leaves it unchanged.

The different models tested in this work show a better agreement with the [Y/Si] rather than with the [Ba/Si]. It is interesting to quantify the overproduction of the overall observed trend shown by our model as well as to compute the amount of Ba, in terms of its surface gas density, that our standard model (Model 2) fails to produce in the few most recent Gyr of chemical evolution, as it causes the trend to be inverted with respect to the observations. We provide here the calculation for the inner, solar, and outer regions. In Figure 8, we show the logarithmic fit of the observed chemical clock [Ba/Si] versus age together with the predictions from Model 2. As already pointed out in Section 4.1, the reduction in the nucleosynthesis prescriptions of s-process elements seem not to be constant among the disc. Indeed, reducing the production of s-process material, independently by the source, by a factor of 2.8, 1.8 and 0.4, leads to the results shown by the dashed black lines in Figure 8 for the inner, solar, and outer regions, respectively. Implying that, while in the solar region the nucleosynthesis is overestimated approximately by 44% (consistent with prior estimates by Rizzuti et al. 2019), in the outer region this percentage is reduced to 28% and in the inner region shows the most dramatic value, around 64%. The obtained reduced results are in agreement with the observations for all the age ranges of interest in the solar and outer regions. On the other hand, in the inner region, the agreement holds for Age ≳ 3.0 Gyr. For younger stellar ages, the model exhibits the insufficient growth trend. To address this, we computed the amount of Ba that our model should produce in the last few Gyr in order to correct its inverted/flat trend with respect to the observations. We computed the residuals between the observed logarithmic fit and our reduced model results, which are shown in the Figure by the yellow dashed lines. In the solar and outer regions, as expected, the residual curve is consistent with the zero, indicating the nice fit between the reduced model and observations. The residuals become positive at very young ages in these regions (Age ≃ 2 Gyr in the local and Age ≃ 1 Gyr in the outer disc), reflecting the flattening of the curve towards younger ages with the consequent underestimation of the observations. The discrepancies are however not dramatic. In the inner region, on the other hand, the residual curve becomes positive for Age ≲ 3 Gyr, indicating that the model increasingly underestimates the observed [Ba/Si] over time when moving towards younger ages. To quantify this discrepancy, we can estimate the missing amount of Ba, both in terms of absolute abundance and of surface mass density (expressed in M_⊙ pc⁻²), provided with the assumption that predictions for Si are correct (see Figure A.1 for a brief discussion). We write, by definition, $$X{(Ba)}_{\text{O}}-X{(Ba)}_{\text{M}}=X{(Si)}_{\text{M}}\times \left({10}^{{[Ba/Si]}_{\text{O}}+\text{Sun}}-{10}^{{[Ba/Si]}_{\text{M}}+\text{Sun}}\right)$$(7)

where $\text{Sun}=\log {\left(\frac{\text{X(Ba)}}{\text{X(Si)}}\right)}_{\odot }$ (taken from Asplund et al. 2009 solar abundances); X(Si)_M is the absolute abundance of Si predicted by the model during the relevant time range; and [Ba/Si]_O and [Ba/Si]_M are the observed and predicted [Ba/Si] ratios, respectively. The obtained missing abundance of Ba is shown in the lower panels of Figure 8 (dashed red lines) for the reduced model in the three regions of interest. It is possible to see that in the inner region, nearly 50% more of the actual predicted amount of Ba, on average, should be produced by the model in order to reproduce the high rise observed towards young ages. In the outer and solar regions, on the other hand, the missing amount of Ba is not significantly highs and this is reflected by the current calculation. In particular, after multiplying by the surface gas density predicted by the model, we can compute the average local surface mass density of Ba that our model fails to produce in the last 3 Gyr of chemical evolution, which is equal to ∑_Ba = 5.71 × 10⁻⁸ M_⊙pc⁻².

6 Conclusions

In summary, we conclude that:

• Two- and three-infall chemical evolution models with state-of-the-art nucleosynthesis prescriptions for neutron-capture elements show difficulties in reproducing the high rise observed in the chemical clocks [s/α] versus age towards younger ages. The deviations between models and data are more evident towards the inner part of the disc and for elements belonging to the second s-process peak (here represented by Ba) than for those of the first peak (here represented by Y). The three-infall model captures the observed trend better, particular the decrease around ≃ 2 Gyr;

• Modifying the production of s-process material from AGB stars, as suggested by previous studies, does not improve the agreement with the data. Modifications that improve the fit for elements of the second s-process peak often fail to match observations for elements of the first peak. An increased s-process contribution from low-mass stars (~1.1 M_⊙) could potentially help reconcile model predictions with observed abundances, yet such an increase is not supported by current nucleosynthesis models, even with the inclusion of the i-process. Future developments in the nucleosynthesis of AGB that take into account, for example, magnetic buoyancy, might lead to improvements in the agreement between models and observations because they affect the production of the s-elements of the first and second peaks differently in terms of mass and metallicity (Magrini et al. 2021; Vescovi 2021);

• Assuming a distribution of rotational velocities for massive stars, which favors high velocities at low metallicities and no velocity at higher ones, does not provide any improvement for elements of the first s-process peak due to the balanced effect on multiple elements. In the case of second peak s-process elements, it accentuates the overproduction across the disc due to the lower production of α-elements;

• Once the residual between the observed trend shown by the OC dataset and our model results is accounted for, it is possible to estimate the missing amount of Ba predicted by our model in the last 3 Gyr in terms of its surface mass density. After scaling the observed trend to our model result, the average local missing Ba is Σ_Ba = 5.71 × 10⁻⁸ M_⊙pc⁻². Namely, the production of Ba in the inner part of the disc should be approximately 50% more than the current level during the last few Gyrs of chemical evolution.

This work gives strong indications that may guide future studies for s-process nucleosynthesis in both AGB and massive stars. Indeed, the theoretical understanding of these elements is important for many aspects of both nucleosynthesis and stellar evolution, as well as the chemical dating of stars, which is becoming increasingly important in view of new large spectroscopic surveys (e.g. Bensby et al. 2019; Jin et al. 2024) and new instrumentation (e.g. Zhang et al. 2018; Magrini et al. 2023a; Mainieri et al. 2024).

Fig. 8 Upper panels: [Ba/Si] versus age trend in the inner, solar, and outer regions. The observational data trend is represented by the logarithmic fit with the associated 1σ error bar curve. The result of Model 2 is shown by the black line. The dashed black and yellow lines represent Model 2 with reduced s-process production and its associated residual curve. Bottom panel: absolute abundance of Ba as predicted by the reduced Model 2 (dashed black line) together with the computed missing abundance of Ba in the age ranges of interest (see text for details).

Acknowledgements

This work was supported by the Deutsche Forschungsge-meinschaft (DFG, German Research Foundation) – Project-ID 279384907 – SFB 1245, the State of Hessen within the Research Cluster ELEMENTS (Project ID 500/10.006). M.M. and F.M. thank INAF. for the 1.05.12.06.05 Theory Grant – Galactic archaeology with radioactive and stable nuclei. This research was partially supported by the Munich Institute for Astro-, Particle and BioPhysics (MIAPbP) which is funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy – EXC-2094 – 390783311. M.M. thanks Almudena Arcones for the useful discussions which improved the paper. L.M. and S.R. thank INAF for the support (Large Grant EPOCH), the Mini-Grants Checs and PILOT (1.05.23.04.02), and the financial support under the National Recovery and Resilience Plan (NRRP), Mission 4, Component 2, Investment 1.1, Call for tender No. 104 published on 2.2.2022 by the Italian Ministry of University and Research (MUR), funded by the European Union –NextGenerationEU– Project ‘Cosmic POT’ Grant Assignment Decree No. 2022X4TM3H by the Italian Ministry of Ministry of University and Research (MUR). M.P. acknowledges financial support from the project “LEGO – Reconstructing the building blocks of the Galaxy by chemical tagging” granted by the Italian MUR through contract PRIN2022LLP8TK_001. CVV acknowledges funding from the Research Council of Lithuania (LMTLT, grant no. P-MIP-23-24).

Appendix A The Si production

The overall conclusions of this work remain robust, even when analyzing different chemical clocks based on various s-process and α-elements. However, it is important to note that some discrepancies may still arise. These uncertainties can stem not only from the detailed nucleosynthesis of s-process elements but also from that of α-elements. For instance, different α-elements exhibit varying production rates from Type Ia SNe, which can influence the steepness of the chemical clock trends as they evolve towards younger ages. While the general conclusions are strong, deviations may occur, particularly in the final computed missing Ba production. This missing Ba production has been estimated under the assumption that our model accurately reproduces the observed Si in the adopted OC dataset. Indeed, our model provides a good fit for the [Si/Fe] versus [Fe/H] trend, as shown and discussed in Figure 2 in Section 3.2. Figure A.1 shows the [Si/H] versus age trend for Model 2, Model 6, and Model 7, the latter two incorporating different massive star nucleosynthesis prescriptions. Model 2, used to compute the missing Ba production, shows a generally good agreement with the observed trend, particularly in the solar and outer Galactic zones. In particular, in the solar zone, the depression around ∼2 Gyr is well-reproduced. The solar Si abundance predicted by our model is 7.44, which compares reasonably with the value of 7.51 ± 0.03 from Asplund et al. (2009). Model 6 and 7, which both assume a distribution of initial velocities of massive stars, predict a lower trend compared with Model 2 as a greater percentage of stars with lower or null rotational velocity. This effect, as discussed in Section 4.3 worsen the agreement with chemical clocks based on Ba and left it unchanged in the case of Y.

Fig. A.1 [Si/H] versus age observed trends divided into the three Galactocentric regions of interest. Predictions of the chemical evolution model in the case of Model 2, Model 6, and Model 7 are compared to the OC sample. (See Table 3 for reference.)

Appendix B The open clusters sample

Abundance ratios of the elements studied in this work, together with Galactocentric distances and ages of the sample of 62 OCs described in section 2 are reported in Table B.1. See Magrini et al. (2023b); Palla et al. (2024) for further discussion.

References

All Tables

All Figures

Fig. 1 Time evolution of the rate of Type Ia SNe, Type II SNe, MR-SNe, and MNS (upper-left panel); SFR (upper-right panel); surface densities of gas (lower-left panel); and stars (lower-right panel) as predicted by the three-infall model. Predictions of the model are compared to present day values from Cappellaro et al. (1999) (for SNe rates), from Abbott et al. (2021) (for MNS rate) and from Prantzos et al. (2018) (for SFR and surface gas and stars densities).

In the text

Fig. 2 [X/Fe] versus [Fe/H] abundance patterns for Si, Y, Ba, and Eu as predicted by the three-infall model at RGC = 6, 8, 12 kpc compared to the sample of OCs (colored dots) and field (grey dots) stars in the inner, solar, and outer regions, respectively. Results of the synthetic model, as well as OCs, are plotted only for Age ≤3 Gyr. In the lower right panel, for the prediction of the [Eu/Fe] versus [Fe/H] in the outer regions, the filled magenta circles indicate the [Eu/Fe] at a given age.

In the text

Fig. 3 [s/Fe] versus Age observed trends for Y (upper panel) and Ba (lower panel) divided into the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. (See Table 1 for reference.)

In the text

Fig. 4 [s/H] versus age observed trends for Y (upper panel) and Ba (lower panel) divided in the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. (See Table 1 for reference.)

In the text

Fig. 5 [s/Si] versus age observed trends for Y (upper panel) and Ba (lower panel) divided in the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. (See Table 1 for reference.)

In the text

Fig. 6 [s/Si] versus age observed trends for Y (upper panel) and Ba (lower panel) divided into the three Galactocentric regions of interest. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region), and magenta (outer region) curves. Predictions of the chemical evolution model in the case of different assumptions concerning the production of s-process material from AGB stars are compared to the OC sample. (See Table 2 for reference.)

In the text

Fig. 7 [s/Si] versus age observed trends for Y (upper panel) and Ba (lower panel) divided in the three Galactocentric regions of interested. The logarithmic fits of the OC samples in the three regions are shown as blue (inner region), green (solar region) and magenta (outer region) curves. Predictions of the chemical evolution model in the case of a two- (dashed lines) and of a three-infall (solid lines) scenario are compared to the OC sample. Predictions of the chemical evolution model in the case of different assumptions concerning the production of s-process material from massive stars are compared to the OC sample. See Table 3 for reference.

In the text

Fig. 8 Upper panels: [Ba/Si] versus age trend in the inner, solar, and outer regions. The observational data trend is represented by the logarithmic fit with the associated 1σ error bar curve. The result of Model 2 is shown by the black line. The dashed black and yellow lines represent Model 2 with reduced s-process production and its associated residual curve. Bottom panel: absolute abundance of Ba as predicted by the reduced Model 2 (dashed black line) together with the computed missing abundance of Ba in the age ranges of interest (see text for details).

In the text

	Fig. A.1 [Si/H] versus age observed trends divided into the three Galactocentric regions of interest. Predictions of the chemical evolution model in the case of Model 2, Model 6, and Model 7 are compared to the OC sample. (See Table 3 for reference.)
In the text