Mixing Uncertainties in Low-Metallicity AGB Stars: The Impact on Stellar Structure and Nucleosynthesis

: The slow neutron-capture process ( s -process) efﬁciency in low-mass AGB stars (1.5 < M/M (cid:12) < 3) critically depends on how mixing processes in stellar interiors are handled, which is still affected by considerable uncertainties. In this work, we compute the evolution and nucleosynthesis of low-mass AGB stars at low metallicities using the MESA stellar evolution code. The combined data set includes models with initial masses M ini /M (cid:12) = 2 and 3 for initial metallicities Z = 0.001 and 0.002. The nucleosynthesis was calculated for all relevant isotopes by post-processing with the NuGrid mppnp code. Using these models, we show the impact of the uncertainties affecting the main mixing processes on heavy element nucleosynthesis, such as convection and mixing at convective boundaries. We ﬁnally compare our theoretical predictions with observed surface abundances on low-metallicity stars. We ﬁnd that mixing at the interface between the He-intershell and the CO-core has a critical impact on the s -process at low metallicities, and its importance is comparable to convective boundary mixing processes under the convective envelope, which determine the formation and size of the 13 C-pocket. Additionally, our results indicate that models with very low to no mixing below the He-intershell during thermal pulses, and with a 13 C-pocket size of at least ∼ 3 × 10 − 4 M (cid:12) , are strongly favored in reproducing observations. Online access to complete yield data tables is also provided.


Introduction
Asymptotic giant branch (AGB) stars are important contributors to the chemical enrichment of our galaxy. In particular, low-mass AGB stars (1.5 < M/M < 3) are the main site of both the main and strong component of s-process ('slow' neutron-capture process [1][2][3][4]), i.e., the nucleosynthesis processes mainly responsible for around half of the neutron-capture element abundances between Zr and Bi in the solar system. The main component produces s-isotopes in the range of 90 < A ≤ 204, and takes place in low-mass AGB stars with around solar metallicity [5,6]. The strong component explains about half of the solar 208 Pb and is hosted by low-metallicity ([Fe/H] ≤ −1), low-mass AGB stars (e.g., [7][8][9]). Throughout the present study, we will adopt the square-bracket notation to express elemental abundance ratios, which are defined as follows: [X/Y] = log((X * /Y * )/(X /Y )). (1) The main difference between models, adopting these different convective boundary mixing (CBM) formalisms, is the remarkable difference in the 13 C-pocket size and the relative abundances of 13 C and 14 N within the pocket, which directly impacts the total amount of s-process elements produced. Indeed, 14 N is a major neutron poison that competes with neutron captures on iron-group and heavier nuclei, reducing the s-process efficiency [29][30][31][32]. For example, the average 13 C-pocket size in [26] models was ∼10 −4 M , while [22] obtained a pocket size between 2 and 3 × 10 −3 M . Moreover, Ref. [33] computed multi-D simulation of the He-flash in low-mass AGB star (M = 2 M and Z = 0.01), showing how CBM is actually active under the intershell during TPs. Including this CBM under the He-intershell is indeed essential to reproduce the observed surface abundances of post-AGB H-deficient stars of the PG1159 class [34] and to match the highest observed [hs/ls] (with hs and ls being the average abundance of s-process elements at the second and first peak, respectively) on the surface of carbon stars. For these reasons it was included in all the models presented in [26]. On the other hand, this additional mixing was not included in [22].
The different methods to include mixing processes in 1D stellar models, reflected in the large range of s-process results presented in the literature, is indicative of the large stellar modeling uncertainties that are still affecting our understanding of AGB nucleosynthesis. Additionally, some of the CBM formalisms mentioned before were validated only for around solar metallicity models. The main reason for this is the large number of observations available for AGB stars with close to solar metallicities, which includes presolar grain measurements (e.g., [35]). However, even if only tested for near-solar metallicity models, these formalisms were then applied to the whole metallicity range, introducing more systematic uncertainties when investigating the s-process at low metallicities.
In this work, we compute the evolution and nucleosynthesis of low-mass AGB stars at low metallicities using the MESA stellar evolution code and the NuGrid mppnp code for the full nucleosynthesis [36]. The combined data set includes models with an initial mass of M ZAMS /M = 3 for Z = 0.001. The main goal is to show the impact of the uncertainties affecting the main mixing processes, such as convection and mixing at convective boundaries, on heavy-element nucleosynthesis.
This work is organized as follows. In Section 2, we describe the stellar code and postprocessing nucleosynthesis tools. In Section 3, the stellar models are presented, together with our nucleosynthesis results, before finally presenting our conclusions in Section 4.

Computational Methods
We used the stellar code MESA (revision 3709, [37]) to compute all the stellar models presented in this work. For the initial composition, we used the alpha-enhanced solar distribution based on [38], which implies a solar metallicity Z = 0.018. Additionally, in Table 1 we provide the initial carbon, nitrogen and oxygen (CNO) abundances adopted. The modeling assumptions are the same as in [39], except for CBM modeling, which is included the same way as in [25,26]. For the simulations, the MESA nuclear network agb.net is used, including 18 isotopes from protons to 22 Ne. Table 1. Initial CNO abundances in mass fractions adopted for the stellar models presented in Tables 2 and 3.

Specie
Initial Abundance Z = 0.001 12  Full nucleosynthesis simulations are obtained by using a post-processing code and the precalculated stellar structure. The post-processing code mppnp is described in detail in [36]. The stellar structure evolution data are computed and saved with MESA for all zones at all time steps, then used as inputs and processed with mppnp. This means that the stellar structure and the full nucleosynthesis are computed separately, hence requiring less computing time and resources. In order to maintain consistency between stellar and nucleosynthesis calculations, MESA and mppnp adopt the same nuclear reaction rates relevant for energy generation and, therefore, for the evolution of the star. The nuclear reaction rates adopted are the same as in [26]. Exceptions relevant for this work are the 22 Ne(α,n) 25 Mg and the 22 Ne(α,γ) 26 Mg reaction rates, for which we use [27].

Description of the Stellar Models
Our stellar models are listed in Tables 2 and 3, giving details of our M ini = 3 M and M ini = 2 M models, respectively. Each model corresponds to a different setting of mixing processes inside the star. The same initial metallicity (Z = 0.001) is adopted for all models in Table 2. All models' names start with a 'm3z1m3'. The first 'm3' means that this is a 3 M model, and 'z1m3' is to be read as Z = 1 × 10 −3 , where the second 'm3' means 'minus three', referring to the exponent. The suffix abbreviations of the models' names indicate the internal mixing setting adopted. "−mlt" means that "mixing length theory" parameters different from other models have been tested, "−bigpoc" indicates that the model has a "bigger pocket of 13 C", while "−no f TP" means that the CBM parameter f adopted under the intershell during a TP is zero. Tables 2 and 3 list key global features like core masses and lifetimes for all the models, which were all computed with the same stellar code and input physics of [39] (hereafter RI18), but with the CBM model of [25] during TDUs. The f 1 and f 2 parameters shown in both tables enter the CBM equation presented in [25] in order to mimic the mixing induced by internal gravity waves (IGWs). Table 2. Main properties of our 3 M asymptotic giant branch (AGB) models: H-free core mass at the beginning and the end of the AGB phase (in solar masses), total lifetime, convective boundary mixing parameters during thermal pulses (TP), and third dredge-up events (TDU) and the adopted mixinglength parameter α are given. Core masses, total lifetimes, and mixing parameters in RI18 for models with the same mass/metallicity combinations are also presented. MLT-mixing-length theory.  Table 3. Same as in Table 2, but for our 2 M AGB models. m2z1m3-bigpoc and m2z2m3-bigpoc have an initial metallicity of Z = 0.001 and Z = 0.002, respectively. For a distance from the Schwarzschild boundary dr < dr 2 , the diffusion coefficient profile is given by the exponentially-decaying diffusion mixing of [40]:

Name
where dr is the geometric distance to the convective boundary. The term f 1 × H p 0 identifies the scale height of the CBM regime. The values D 0 and H p 0 are, respectively, the diffusion coefficient D and the pressure scale height at the convective boundary. A second, slowerdecreasing mixing coefficient is included in the same way as in [25], in order to also consider the IGW contribution discussed by [10]. IGWs are the dominating CBM mechanism at a distance dr 2 from the Schwarzschild boundary, when the mixing coefficient defined in Equation (2) falls below a value D 2 defined as with length scale f 2 × H p0 , which is adopted for distances dr > dr 2 . Therefore, for dr > dr 2 : In our models, we adopt the same D 2 value as calibrated and adopted in [26], i.e., D 2 = 4.3 × 10 11 cm 2 s −1 .
We also included the M ini = 3 M , Z = 0.001 model from RI18 in Table 2 as a reference. Notice that the full nucleosynthesis of the RI18 model has been recomputed, adopting the newest 22 Ne(α,n) 25 Mg and 22 Ne(α,γ) 26 Mg nuclear reaction rates from [27]. The grid of models presented explores the impact of different combinations of mixing parameters. In particular, one of the biggest uncertainties in 1D AGB models is the treatment of convection, which would require 3D modeling [41][42][43]. The most adopted treatment in 1D stellar evolution relies on mixing-length theory (MLT [44,45]). The free mixing-length parameter α MLT specifies the mean free path (i.e., the mixing-length) of a convective blob in units of the pressure scale height, and is usually calibrated on the Sun. Ref. [46] claimed that a mixing-length dependent on metallicity is required to match the surface temperatures of red giants in the APOKASC catalog [47], reporting a decreasing mixing-length when decreasing metallicity of about 0.16 per metallicity dex. Based on [46] results, m3z1m3mlt adopts a lower α MLT parameter compared to the one adopted for solar metallicity models in [26], while all other models in Table 2 keep the solar-calibrated value. Indeed, a remarkable impact on the initial H-free core mass (which is then reflected on the final core mass) is visible when considering model m3z1m3-mlt. This is due to the fact that reducing the α MLT directly impacts convection efficiency in the star, in particular, the second-dredge up, which is what stops the H-free core growth after the end of core He-burning. On the other hand, Ref. [48] derived a mixing-length parameter from their 3D simulations of the external part of the AGB envelope convection. They found α MLT = 2.6, larger than our solar-calibrated value, in contrast with the value for m3z1m3-mlt. Ref. [48] results point to a strong evolutionary phase dependence of the mixing-length, which is not included in the simpler metallicity-dependence relation derived by [46].
Another relevant impact, visible from Table 2, comes from the treatment of CBM under the He-intershell during TPs. Both m3z1m3-nofTP and m3z1m3-bigpoc do not include this CBM. As a consequence, they both show a higher final H-free core mass compared to m3z1m3, where CBM during TPs is included. The reason is shown in Figure 1. The upper panel shows the Kippenhahn diagram of m3z1m3-nofTP zoomed in the Heintershell. The lower panel presents the same kind of diagram, but it shows the structure evolution of m3z1m3. When a He-flash occurs in the intershell, it develops a pulse-driven convective zone (PDCZ) mixing the whole intershell, which is then followed by a TDU event, where the hydrogen convective envelope penetrates into the H-free core. The major difference between the two stellar structures is the lower reduction of the H-free core size immediately after a TDU event, which is the main cause of the difference in the final core size, in particular, after the first five TDUs. This is due to the lower He-flash luminosity of m3z1m3-nofTP compared to m3z1m3, as shown in Figure 2. A lower luminosity of the He-flashes leads to a smaller expansion of the convective envelope, and consequently, it will cool down less. This will cap the opacity increase in the envelope, limiting the increase of convection efficiency in the envelope. This will directly impact the TDU penetration into the He-intershell, which is visibly reduced in m3z1m3-nofTP compared to m3z1m3. The origin of this difference in the luminosity of He-flashes lies in the different He abundance in the intershell, as shown in Figure 3, which is up to 40% less abundant in m3z1m3. This is a direct consequence of the different treatment of CBM under the PDCZ, since through CBM, C-rich material can be dredged-up into the intershell from the core, decreasing the He mass-fraction. In turn, a lower abundance of He causes the higher temperature necessary to trigger the He-flash.   The slightly lower core mass of m3z1m3-bigpoc, compared to m3z1m3-nofTP, is due to the higher CBM efficiency under the convective envelope during TDUs, which also has a strong direct impact on heavy element nucleosynthesis (see next chapter). CBM parameters in m3z1m3-nofTP during TDUs are the same ones adopted in [26], guided by the internal gravity waves mixing profile in solar metallicity AGB stars obtained by [10]. However, as already discussed, this may not apply to the different structure of low-Z AGB stars. Hence, we explore the effects of a stronger CBM under the convective envelope. Indeed, this causes a deeper diffusion of hydrogen from the envelope, forming a bigger 13 C-pocket, as shown in Figure 4. In particular, the average 13 C-pocket size in m3z1m3-bigpoc is between 3 and 4 × 10 −4 M , about a factor of three larger than the pocket size formed in m3z1m3-nofTP.   13 C-pocket at the end of the oxygen-rich phase from m3z1m3. Lower panel: 13 C-pocket at the end of the oxygen-rich phase from m3z1m3-bigpoc. The comparison shows a much larger pocket formed by m3z1m3-bigpoc compared to m3z1m3 (about a factor of three) as a consequence of the slower decay of the mixing coefficient adopted for this model (see Table 2 for the mixing parameters adopted).

Postprocessing Nucleosynthesis Calculations
The s-process nucleosynthesis in low-mass AGB stars heavily depends on internal mixing processes. The resulting heavy element production of all the models presented in this work is shown in Figure 5. In the upper panel, RI18 exhibits the lowest s-process production, globally between a factor of two and three lower compared to m3z1m3, whose only difference from RI18 is a more efficient CBM under the convective envelope, which produces a 13 C-pocket about three times larger. The same difference is visible when comparing m3z1m3 with m31m3-bigpoc in the lower panel, since in this case, the typical 13 C-pocket also differs by a factor of three between the two models. Surface abundances of the Ba star HD 123396 [49] and CEMP-s star HD 26 [50] for two first-peak elements (Y and Zr), two second-peak elements (Ce and Nd), and Pb are also plotted for comparison. The two Pb abundances are only from HD 26, with the higher one being determined including non-LTE correction. m3z1m3-bigpoc is the only model consistent with first, second, and third s-process peak elements at the same time. Unfortunately, Pb abundances are not available for most of the sample stars we selected, and the determination of Pb abundances has some difficulties, often introducing big uncertainties [51][52][53]. Considering elements lighter than first-peak ones, m3z1m3 shows high enrichment, as shown in more detail in Figure 6, zoomed in the 30 < A < 41 region. The [Rb/Sr] ratio is visibly higher in m3z1m3 compared to both m3z1m3-nofTP and m3z1m3-bigpoc. In particular, we recall that the only difference between m3z1m3 and m3z1m3-nofTP is the inclusion of CBM processes under the He-intershell in m3z1m3, while they are not considered in m3z1m3-nofTP. As we already discussed, this causes more luminous and hotter TPs, which eventually results in opening both the 85 Kr and 86 Rb branching points, efficiently producing both the neutron magic 86 Kr and 87 Rb isotopes. This high enrichment in Kr and Rb during the TP is then converted into Sr, Y, and Zr by neutron captures during the radiative burning of the 13 Cpocket. It is worth noticing that it is harder for the s-process nucleosynthesis flux to go past the first peak passing through 86 Kr and 87 Rb, compared to the lower neutron-density case where the two branching points are less open and 86 Sr is produced. Indeed, in high-neutron density conditions, caused by an efficient 22 Ne(α,n) 25 Mg in hot TPs, the neutron-capture chain will need to take place on four to five neutron magic nuclei, compared to only three in the case of colder TPs. This will cause substantial over-production of first-peak elements, and therefore, a lower [Ce/Y] in m3z1m3 compared to both m3z1m3-nofTP and m3z1m3-bigpoc. This last aspect is confirmed by Figure 7. The [Ce/Y] ratio was used by Cseh et al. [49] as a representative of the [hs/ls] ratio. These elements were chosen based on the reliability and accuracy of the determined abundances of the sample Ba stars in their study, since Ba abundances were not available and [La/Fe] showed unexpectedly high values for some of the sample stars. [hs/Y] and [hs/Zr] show very similar ratios; here, we use [Ce/Y] to match their choice. The colored area between the two vertical dashed lines in Figure 7 indicates the range within which 80% of our low-Z AGB abundance observations are found. Despite the fact m3z1m3-nofTP is not able to reproduce the full [Ce/Y] range, it is still the best performing model among those presented in the figure when compared to observations. This suggests that models with no CBM under the He-intershell are largely favored, as this is apparently an essential condition to reproduce the observed high [Ce/Y] on low-Z AGB surfaces.  [49] and CEMP-s star HD 26 [50] for two first-peak (Y and Zr) and second-peak elements (Ce and Nd); Pb are also plotted for comparison. The two Pb abundances are only from HD 26, with the higher one being determined including non-LTE correction. Lower panel: Same as in the upper panel, but for m3z1m3, m3z1m3-nofTP, and m3z1m3-bigpoc.  On the other hand, Ref. [33] performed multi-D simulations of the TP in a Z = 0.01 AGB star. They found that mixing across convective boundaries is significant for He-shell flash convection. Guided by these findings, Refs. [25,26] implemented a consistent CBM under the He-intershell during TPs in their around solar-metallicity models. However, the Z = 0.001 models presented here are remarkably different from their Z = 0.01 counterparts with the same mass. For example, m3z1m3 from this work has an H-free core mass about 0.15 M larger than m3z1m2 (M = 3 M , Z = 0.01) in [26]. This means that the conditions found in m3z1m3 are probably significantly different from those explored by [33], questioning the validity of their results in low-Z AGBs. In particular, a significantly massive core is also more compact, hence, possibly inhibiting CBM processes at the interface with the He-intershell.
In Figure 8, we show again our predicted [Ce/Fe] vs. [Ce/Y], but including m3z1m3bigpoc, m2z1m3-bigpoc, m2z2m3-bigpoc, and the observed abundances from different objects enriched with s-process elements from AGB stars: 10 post-AGBs [54,55], CEMP-s [50] and Ba stars [49,[56][57][58]. Names, references and [Fe/H] for each star are given in Table 4. CEMP-s and Ba stars are members of a binary system, in which the accreted material from the former AGB companion polluted the now-observed star. Thus, these stars retain the abundance pattern of their AGB companion, and are appropriate for comparison with AGB nucleosynthesis models. The straight line across the observational data is a linear regression showing the general trend of observed abundances. Comparing m3z1m3-bigpoc and m3z1m3-nofTP, we conclude that a very low to no CBM under the He-intershell is a necessary yet insufficient condition to reproduce the bulk of the observed s-process values. Indeed, a 13 C-pocket as large as in m3z1m3-bigpoc, i.e., typically 3 × 10 −4 M , is an additional essential condition to be in agreement with the observed abundances. Hence, any CBM process under the convective envelope during TDUs and included in stellar models should be efficient enough to produce a 13 C-pocket as large as at least 3 × 10 −4 M . Additionally, m3z1m3-bigpoc is consistent with both models with the same mass and metallicity from the FRUITY database [59,60], which did not include CBM under the He-intershell and form a 13 C-pocket comparable in size to m3z1m3-bigpoc ones. Both FRUITY and [60] models are also shown in Figure 8.  Here, we also include results from the FRUITY database and the Monash group models as a comparison. Observed abundances from 10 post-AGB [55,61], CEMP-s [50], and Ba stars [49,56,58] are also included. The straight line across the observational data is a linear regression showing the general trend of observed abundances.
Finally, we notice how adding our 2 M models helps with covering the observed range of both [Ce/Fe] and [Ce/Y]. As expected, the higher metallicity of m2z2m3-bigpoc causes less compact and a hotter interior, leading to less-efficient Rb production during TPs. As already discussed (e.g., Figure 6), this translates into lower first-peak element abundances, increasing the [Ce/Y] ratio. In particular, m2z2m3-bigpoc successfully matches the two most-enriched s-process stars in the sample, with a surface [Ce/Y] > 1. This means that at least part of the observed spread in s-process efficiencies can be explained by a variation of mass and metallicity. Consistent with our models, the FRUITY 2 M , Z = 0.002 model shows both a [Ce/Fe] and [Ce/Y] higher than its 3 M , Z = 0.001 counterpart. On the other hand, several nuclear and stellar evolution uncertainties, such as those related to stellar binarity, can still have an impact. For instance, one could speculate that a mass transfer event from a massive AGB companion, which typically produces an envelope composition with [Ce/Y] < 0 (e.g., [39,62]), might be able to explain those stars with both high [Ce/Fe] and [Ce/Y] close to zero. We will explore this hypothesis in detail in a follow-up study (Battino et al. 2021b in preparation).

Conclusions
In this work, we computed the evolution and nucleosynthesis of low-mass AGB stars at low metallicities (Z = 0.001 and 0.002), testing the impact of the uncertainties affecting the main mixing processes, such as convection and mixing at convective boundaries, on heavy element nucleosynthesis.
We found that in order to be consistent with observed heavy-element abundances, AGB models should include no CBM under the He-intershell, though a CBM process under the convective envelope, efficient enough to form a large 13 C-pocket of at least 3 × 10 −4 M , is required. In this way, our models show a final surface composition with 1.16 < [Ce/Fe] < 2.08 and 0.61 < [Ce/Y] < 1.20, consistent with the bulk of observations. We will verify in a forthcoming study if internal gravity waves are actually capable of producing such a 13 C-pocket in low-Z AGB stars.
Uncertainties related to the treatment of convection through parametric MLT impact the final [Ce/Y] by about +/−0.1 dex. On one hand, this is minor compared to the role of CBM at convective boundaries, yet, enough to make future improvement of the handling of convection in 1D stellar models highly desirable, especially considering the visible impact on core masses.
Additionally, Pb abundances (third peak of the s-process) are important to understand the nucleosynthetic processes leading to the overabundance of this element in s-process enriched stars. Unfortunately, Pb abundances are not available for most of our sample stars and the determination of Pb abundances has some difficulties: the Pb I line at 4057.8 Å, detectable in the optical wavelength, is often blended with CH and needs high-resolution, high signal-to-noise blue spectra for the derivation of the abundance. Additionally, the inclusion of non-LTE corrections is necessary to the most accurate comparison of the models with Pb and other elements [51][52][53]. Pb abundances from CEMP-s star HD 26 are indeed reproduced by our m3z1m3-bigpoc model, but larger statistical modeling is still prevented by the observational difficulties just mentioned.
Further uncertainties affecting the production of heavy elements, yet only briefly mentioned in this work, originate from nuclear physics and stellar binarity, which could in principle explain the spread in observed [Ce/Y] values, in particular, the lowest ones. These additional sources of uncertainties will be investigated in detail in a forthcoming paper (Battino et al. 2021b in preparation).

Data Availability Statement:
The data presented in this study are openly available in NuGRID repository, as specified at https://nugrid.github.io/content/data.