Constraining the thermally-pulsing asymptotic giant branch phase with resolved stellar populations in the Small Magellanic Cloud

TL;DR

This study calibrates models of the thermally-pulsing asymptotic giant branch phase using detailed observations of the Small Magellanic Cloud, improving understanding of this uncertain stellar evolution stage.

Contribution

It introduces a calibration method combining resolved stellar observations with detailed simulations to refine TP-AGB models, especially regarding dredge-up and mass loss efficiencies.

Findings

Identified two best-fitting TP-AGB models with different dredge-up and mass loss parameters.

Provided a detailed characterization of the TP-AGB population in the SMC.

Published extensive isochrone tables with improved models.

Abstract

The thermally-pulsing asymptotic giant branch (TP-AGB) experienced by low- and intermediate-mass stars is one of the most uncertain phases of stellar evolution and the models need to be calibrated with the aid of observations. To this purpose, we couple high-quality observations of resolved stars in the Small Magellanic Cloud (SMC) with detailed stellar population synthesis simulations computed with the TRILEGAL code. The strength of our approach relies on the detailed spatially-resolved star formation history of the SMC, derived from the deep near-infrared photometry of the VISTA survey of the Magellanic Clouds, as well as on the capability to quickly and accurately explore a wide variety of parameters and effects with the COLIBRI code for the TP-AGB evolution. Adopting a well-characterized set of observations -- star counts and luminosity functions -- we set up a calibration cycle along which we iteratively change a few key parameters of the TP-AGB models until we eventually reach a good fit to the observations. Our work leads to identify two best-fitting models that mainly differ in the efficiencies of the third dredge-up and mass loss in TP-AGB stars with initial masses larger than about 3 M$_{\odot}$. On the basis of these calibrated models we provide a full characterization of the TP-AGB stellar population in the SMC in terms of stellar parameters (initial masses, C/O ratios, carbon excess, mass-loss rates). Extensive tables of isochrones including these improved models are publicly available.

Tables

Table 1. Table 1: Central coordinates of the VMC tiles used in this work and number of AGB stars identified by SR16 for each tile.

Tile1	N. subreg.	R.A.J2000	Dec.J2000	N. AGB
	used	(deg)	(deg.)
SMC 3_2	4	5.8981	$-$74.1159	136
SMC 3_3	11	11.2329	$-$74.2117	578
SMC 3_4	11	16.588	$-$74.1774	272
SMC 3_5	12	21.8784	$-$74.0137	63
SMC 4_2	11	6.3087	$-$73.0299	205
SMC 4_3	12	11.3112	$-$73.1198	1888
SMC 4_4	12	16.3303	$-$73.0876	1146
SMC 4_5	12	21.2959	$-$72.9339	99
SMC 5_3	12	11.2043	$-$72.0267	484
SMC 5_4	12	16.1088	$-$71.9975	655
SMC 5_5	12	20.7706	$-$71.8633	84
SMC 6_3	9	11.4532	$-$70.9356	57
SMC 6_4	10	15.9581	$-$70.8929	96
SMC 6_5	7	20.3437	$-$70.7697	29
	NTot	Total		NTot
	subreg.	area		AGB
	147	$\approx 21\text{deg}\text{<sup><span>2</span></sup>}$		5792
Notes:
(1) Excluded subregions:
SMC 3_2 G1, G2, G3, G4, G6, G7, G8, G12;
SMC 3_3 G1; SMC 3_4 G4; SMC 4_2 G4; SMC 6_3 G8, G11,
G12; SMC 6_4 G6, G7; SMC 6_5 G5, G6, G9, G10, G11.

Table 2. Table 2: Number counts of RSG/AGB populations classified by SR16. The star counts refer to the selected SMC area used in this work.

Population	Number of stars
C-AGB	1854
O-AGB	2623
X-AGB	343
a-AGB	972
RSG	3150

Table 3. Table 3: Age resolution in TRILEGAL  simulations.

Period	$\log \left(\mathrm{age}/\mathrm{yr}\right)$ interval	Resolution
	(dex)	(dex)
Pre-boosting	6.600 - 9.000	0.020
Boosting	9.000 - 9.400	0.001
Post-boosting	9.400 - 10.170	0.020

Table 4. Table 4: Grid of TP-AGB sets

	Mass Loss							Third Dredge-Up
	Pre-Dust		Dust - driven					Activation	Efficiency
			M-stars		C-stars
SET	Id	$\eta$	Id	$\eta$	Id	$\eta$	${\dot{M}}^{\left(1\right)}$	$\log {(T_b^{\mathrm{dred}}/\left[K\right])}^{\left(2\right)}$	$\lambda$
S_00	SC05	–	BE88	–	BE88	–	(a)	6.40	K02
S_01	CS11	2	BL95	0.05	BL95	0.05	(a)	6.40	K02
S_02	CS11	2	BL95	0.02	CDYN	1	(a)	6.40	K02
S_03	CS11	2	BL95	0.03	CDYN	1	(a)	6.40	K02
S_04	CS11	2	BL95	0.05	CDYN	1	(a)	6.40	K02
S_05	CS11	2	BL95	0.06	CDYN	1	(a)	6.40	K02
S_06	CS11	3	BL95	0.06	CDYN	1	(a)	6.40	K02
S_07	CS11	3	BL95	0.06	CDYN	1	(a)	f1($Z_i$)	K02
S_19	CS11	3	VW93	–	CDYN	1	(a)	f1($Z_i$)	K02
S_12	CS11	3	BL95	0.01	CDYN	1	(b)	f1($Z_i$)	K02
S_08	CS11	3	BL95	0.06	CDYN	1	(a)	f1($Z_i$)	${\lambda }_{\max }$= 0.5
S_09	CS11	3	BL95	0.02	BL95	0.02	(a)	f1($Z_i$)	${\lambda }_{\max }$= 0.5
S_10	CS11	3	BL95	0.01	BL95	0.01	(a)	f1($Z_i$)	${\lambda }_{\max }$= 0.5
S_11	CS11	3	BL95	0.01	CDYN	1	(b)	f1($Z_i$)	${\lambda }_{\max }$= 0.5
S_31	CS11	3	BL95	0.01	CDYN	1	(b)	f1($Z_i$)	${\lambda }_{\max }$= 0.4
									${\lambda }_{\max }^{\ast }$	${\tilde{M}}_c\left[M_{\odot }\right]$	$M_{c,\lambda =0}\left[M_{\odot }\right]$
S_13	CS11	3	BL95	0.01	CDYN	1	(b)	f1($Z_i$)	0.5	0.65	0.95
S_14	CS11	3	BL95	0.01	CDYN	1	(b)	f1($Z_i$)	0.5	0.65	0.85
S_15	CS11	3	BL95	0.01	CDYN	1	(b)	f1($Z_i$)	0.6	0.60	0.85
S_16	CS11	3	BL95	0.01	CDYN	1	(b)	f2($Z_i$)	0.6	0.60	0.85
S_17	CS11	3	BL95	0.01	CDYN	1	(c)	f2($Z_i$)	0.6	0.60	0.85
S_18	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.6	0.60	0.85
S_20	CS11	3	BL95	0.02	CDYN	1	(a)	f2($Z_i$)	0.6	0.60	0.85
S_22	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.7	0.70	0.85
S_23	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.8	0.60	0.85
S_24	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.8	0.60	1.30
S_25	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.8	0.60	1.00
S_26	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.7	0.60	0.85
S_27	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.7	0.60	1.00
S_28	CS11	3	BL95	0.02	CDYN	1	(b)	f2($Z_i$)	0.7	0.60	1.00
S_29	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.7	0.625	1.00
S_30	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.5	0.6	1.00
S_32	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.5	0.5	1.00
S_34	CS11	3	BL95	0.02	CDYN	1	(c)	f2($Z_i$)	0.6	0.6	1.00
S_35	CS11	3	BL95	0.03	CDYN	1	(c)	f2($Z_i$)	0.7	0.60	1.00
Notes:
(1) Actual mass-loss rate:
(a) max(${\dot{M}}_{\mathrm{pre}-\mathrm{dust}}$,${\dot{M}}_{\mathrm{dust}}$)
(b) max(${\dot{M}}_{\mathrm{pre}-\mathrm{dust}}$,${\dot{M}}_{\mathrm{dust}}$), for $C/O>1$ if ${\dot{M}}_{\mathrm{dust}}$is not active, we assume
(c) max(${\dot{M}}_{\mathrm{pre}-\mathrm{dust}}$,${\dot{M}}_{\mathrm{dust}}$) for ; for $C/O>1$ max(${\dot{M}}_{\mathrm{dust}}$(), ${\dot{M}}_{\mathrm{dust}}$($C/O>1$) )
(2) 3DU activation temperature as a function of $Z_i$:
f1($Z_i$): $\log T_b^{\mathrm{dred}}=\max [6.3,T_1+\left(T_2-T_1\right)\left(Z-Z_1\right)/\left(Z_2-Z_1\right)]$ with $T_1=6.3$, $T_1=6.60$, $Z_1=0.001$, $Z_2=0.02$.
f2($Z_i$): $\log T_b^{\mathrm{dred}}=\max [6.2,T_1+\left(T_2-T_1\right)\left(Z-Z_1\right)/\left(Z_2-Z_1\right)]$ with $T_1=6.1$, $T_1=6.75$, $Z_1=0.001$, $Z_2=0.02$.

Table 5. Table 5: Comparison of observed star counts with COLIBRI  models.

	AGB			O-AGB			C-AGB			X-AGB			C/M
	$N_{\mathrm{TOT}}$	$\Delta {N\%}^1$	${\chi }_{\mathrm{LF}}^2$	$N_{\mathrm{TOT}}$	$\Delta {N\%}^1$	${\chi }_{\mathrm{LF}}^2$	$N_{\mathrm{TOT}}$	$\Delta {N\%}^1$	${\chi }_{\mathrm{LF}}^2$	$N_{\mathrm{TOT}}$	$\Delta {N\%}^1$	${\chi }_{\mathrm{LF}}^2$
OBS	5778			3109			2340			329			0.86
S_00	6453	11.7	22.3	4484	44.2	20.5	1655	-29.3	21.6	313	-4.9	4.8	0.44
S_01	5385	-6.8	16.0	3881	24.8	9.1	1335	-42.9	71.9	168	-48.9	6.5	0.39
S_02	7596	31.5	34.0	3756	20.8	7.3	3430	46.6	36.8	409	24.3	7.0	1.02
S_03	7395	28.0	27.2	3690	18.7	6.7	3321	41.9	30.4	383	16.4	5.0	1.0
S_04	6894	19.3	10.7	3682	18.4	5.8	2894	23.7	9.4	317	-3.6	4.8	0.87
S_05	6713	16.2	8.2	3666	17.9	6.2	2755	17.7	7.2	290	-11.9	5.2	0.83
S_06	5728	-0.9	2.4	3155	1.5	2.6	2329	-0.5	3.6	243	-26.1	3.5	0.82
S_07	5805	0.5	2.3	3148	1.3	2.3	2380	1.7	4.0	276	-16.1	3.5	0.84
S_08	5549	-4.0	2.7	3193	2.7	2.6	2095	-10.5	6.7	260	-21.0	2.1	0.74
S_09	5521	-4.4	4.9	3334	7.2	4.9	1990	-15.0	6.7	197	-40.1	5.1	0.66
S_10	5980	3.5	7.4	3387	8.9	6.1	2260	-3.4	3.4	332	0.9	4.9	0.76
S_11	6072	5.1	12.6	3368	8.3	6.2	2366	1.1	9.1	337	2.4	4.5	0.8
S_12	6302	9.1	14.8	3317	6.7	4.3	2599	11.1	11.4	385	17.0	6.9	0.9
S_13	5882	1.8	12.7	3636	17.0	9.3	2014	-13.9	14.9	231	-29.8	5.6	0.62
S_14	5820	0.7	12.5	3774	21.4	13.6	1841	-21.3	19.9	205	-37.7	8.4	0.54
S_15	5843	1.1	10.9	3632	16.8	12.1	1952	-16.6	9.6	258	-21.6	3.4	0.61
S_16	5894	2.0	10.8	3491	12.3	10.5	2117	-9.5	5.2	284	-13.7	3.7	0.69
S_17	6602	14.3	18.2	3486	12.1	10.5	2792	19.3	11.3	323	-1.8	6.2	0.89
S_18	6134	6.2	6.3	3407	9.6	6.2	2398	2.5	2.9	328	-0.3	5.3	0.8
S_19	6542	13.2	29.7	3211	3.3	3.4	2948	26.0	32.2	382	16.1	8.0	1.04
S_20	5991	3.7	10.2	3436	10.5	6.5	2241	-4.2	7.8	314	-4.6	4.2	0.74
S_22	5882	1.8	7.5	4132	32.9	19.2	1508	-35.6	59.6	242	-26.4	4.0	0.42
S_23	6351	9.9	7.8	3289	5.8	5.9	2707	15.7	4.1	354	7.6	6.8	0.93
S_24	6517	12.8	10.9	3230	3.9	3.3	2874	22.8	9.1	412	25.2	8.5	1.02
S_25	6415	11.0	8.2	3223	3.7	4.0	2824	20.7	6.0	368	11.9	6.6	0.99
S_26	6268	8.5	6.7	3331	7.1	5.4	2589	10.6	3.0	347	5.5	4.6	0.88
S_27	6259	8.3	7.4	3307	6.4	5.2	2603	11.2	3.2	347	5.5	4.4	0.89
S_28	5742	-0.6	6.0	3305	6.3	4.9	2111	-9.8	3.2	325	-1.2	5.2	0.74
S_29	6292	8.9	7.2	3325	6.9	4.4	2609	11.5	3.6	356	8.2	4.7	0.89
S_30	6104	5.6	5.6	3483	12.0	6.2	2319	-0.9	4.1	301	-8.5	3.0	0.75
S_31	6098	5.5	11.5	3360	8.1	6.3	2414	3.2	7.6	322	-2.1	5.4	0.81
S_32	6063	4.9	6.8	3482	12.0	7.2	2281	-2.5	4.2	299	-9.1	3.4	0.74
S_34	6201	7.3	6.7	3376	8.6	5.5	2480	6.0	3.0	345	4.9	4.8	0.84
S_35	5905	2.2	4.2	3244	4.3	4.4	2315	-1.1	3.0	345	4.9	4.5	0.82
Notes:
(1) Percentage of the difference between model and observations with respect to the observations: $100\times \left(N_{\mathrm{model}}-N_{\mathrm{obs}.}\right)/N_{\mathrm{obs}.}$

Equations

$$\mbox{$\chi^{2}_{\rm LF}$}=\frac{1}{N}\sum_{i=1}^{N}2\cdot\left[m_{i}-n_{i}+n_{i}\cdot\ln\left(\frac{n_{i}}{m_{i}}\right)\right]$$

$$\mbox{$\mathrm{C}-\mathrm{O}$}=\log(n_{\rm C}-n_{\rm O})-\log(n_{\rm H})+12\,,$$

$$\dot{M}=\mathrm{max}(\mbox{$\eta_{\rm pre-dust}$}\cdot\mbox{$\dot{M}_{\rm pre-dust}$},\,\mbox{$\eta_{\rm dust}$}\cdot\mbox{$\dot{M}_{\rm dust}$})\,,$$

$$\lambda ={\lambda }_{\max }^{\ast }\left[{\left(\frac{\xi -\tilde{\xi }}{\tilde{\xi }-{\xi }_{\lambda =0}}\right)}^2-1\right]\left[{\left(\frac{M_{\mathrm{c}}-{\tilde{M}}_{\mathrm{c}}}{{\tilde{M}}_{\mathrm{c}}-M_{\mathrm{c},\lambda =0}}\right)}^2-1\right],$$

$$\mathrm{K0}=-0.48\times(\mbox{$J\!-\!K_{\rm s}$})+13.022+0.056\mbox{\rm[{\rm M}/{\rm H}]}$$

$$K1$$

$$K2$$

$$[8]=A-11.76\times(\mbox{$J\!-\![8]$}),$$

$$d{A}_J^{\mathrm{RGB}}$$

$$d{A}_J^{\mathrm{RSG}}$$

$$d{A}_J^{\mathrm{AGB}}$$

Full text

††pubyear: 2019††pagerange: Constraining the thermally-pulsing asymptotic giant branch phase with resolved stellar populations in the Small Magellanic Cloud–B

1 Introduction

The thermally-pulsing asymptotic giant branch (TP-AGB) is an advanced stellar evolutionary phase during which stars of low- and intermediate-mass (with initial masses in the range 0.8~{}\mbox{\mathrm{M}{\odot}}\la\mbox{M{\rm i}}\la 6-8~{}\mbox{\mathrm{M}_{\odot}}) experience helium and hydrogen double shell-burning, reach their highest luminosity, synthesise new chemical elements, undergo powerful mass loss via stellar winds, and then finally reach the stage of white dwarfs (Herwig, 2005). The impact of this evolutionary phase spans from chemical evolution and spectral energy distribution of galaxies to dust production in galaxies at low and high redshift (Maraston et al., 2006; Conroy, 2013; Zibetti et al., 2013). Despite its importance in our understanding of galaxy evolution, TP-AGB modelling is still affected by large uncertainties due to the presence of several and interconnected processes – i.e. third dredge-up (3DU), hot-bottom burning (HBB), stellar winds, long period pulsations, reprocessing of radiation by circumsttextllar dust – for which a robust theory is still lacking (Marigo, 2015).

Over the years, TP-AGB evolutionary models have been calibrated using primarily two different kinds of observations: the AGB stars in well-populated star clusters, and the luminosity functions of AGB stars in the fields of nearby galaxies.

The first approach regards mainly the population of M- and C-type AGB stars found in the most populous star clusters in the Magellanic Clouds (Frogel et al., 1990), whose turn-off masses are generally well constrained, e.g. by means of isochrone fitting of deep colour-magnitude diagrams (CMDs). The observed numbers and luminosities of TP-AGB stars can be easily compared to the predictions from single-age, single-metallicity models, providing constraints to lifetimes, surface composition, and core masses of TP-AGB tracks with sub-solar metallicity (e.g. Girardi & Marigo, 2007), as a function of their initial masses. This approach, however, is not completely satisfactory, and possibly fails, for a series of reasons. First, even in the most populous clusters the numbers of TP-AGB stars do not exceed two dozen, and is typically much fewer and/or close to unity. These low number statistics lead to large uncertainties in all results that depend on the numbers of TP-AGB stars, including their lifetimes and their contribution to the integrated light of stellar populations. Second, due to the “AGB boosting” effect (Girardi et al., 2013), there is no guarantee that star counts in Magellanic Clouds’ clusters of ages $\sim\!1.6$ Gyr – the approximate age of the most populous ones – are proportional to the lifetime in the TP-AGB phase. Finally, the discovery of broad main sequence turn-offs in some of the most massive Magellanic Clouds’ clusters, i.e. those with the largest early escape velocity (Goudfrooij et al., 2014), and their possible connection with fast rotation (Brandt & Huang, 2015; Goudfrooij et al., 2017), opens the possibility that the clusters’ population may not be representative of the stars found in galaxy fields.

The second alternative is to use AGB stars in galaxy regions for which the distributions of ages and metallicities are well constrained, comparing their numbers and luminosities with those predicted through the stellar population synthesis approach. Given galaxy regions large enough to contain hundreds of AGB stars, this method overcomes the difficulties related to small-number statistics and the AGB-boosting in clusters. This technique was pioneered by Groenewegen & de Jong (1993), who used the carbon star luminosity function in the Magellanic Clouds to calibrate the third dredge-up efficiency in synthetic AGB models. The same approach was later revised with increasingly detailed models of the AGB evolution (e.g. Marigo et al., 1999, 2003; Izzard et al., 2004; Stancliffe et al., 2005; Marigo et al., 2013; Marigo et al., 2017), and extended to many dwarf galaxies with good estimates of their star formation histories (SFHs) through analysis of Hubble Space Telescope (HST) photometry (e.g. Girardi et al., 2010; Rosenfield et al., 2014, 2016). The obvious shortcoming of this approach, compared to the use of star clusters mentioned above, is that we loose information about the initial stellar masses of TP-AGB stars. Indeed, although the main subtypes of TP-AGB stars broadly separate in CMDs like the $K_{\rm s}$ versus $J\!-\!K_{\rm s}$ one, these sequences represent relatively wide mass ranges, in which tracks of different initial masses partially overlap. Despite the somewhat limited resolution in initial mass, the method benefits from a much better statistics.

Over the last decade, the quality of the data necessary to perform this calibration work using SMC field stars has improved significantly. First, catalogues of candidate TP-AGB stars were extended, from the initial near-infrared samples provided by the Deep Near Infrared Survey of the Southern Sky (DENIS, Cioni et al., 2000) and the Two Micron All-Sky Survey (2MASS Nikolaev & Weinberg, 2000; Cutri et al., 2003), to include the mid-infrared data from Spitzer surveys (Gordon et al., 2011; Boyer et al., 2011; Srinivasan et al., 2016); moreover, a non-negligible fraction of the SMC field TP-AGB stars now have a spectroscopic classification (Ruffle et al., 2015; Boyer et al., 2015). Second, we have now a significantly better description of the space-resolved distribution of stellar ages, metallicities, distances, and mean extinctions across the SMC galaxy (Rubele et al., 2015, 2018), as will be illustrated below. To these positive aspects, we can add the significant advancement in the physical prescriptions adopted in TP-AGB evolutionary models (e.g., gas opacities for varying chemical mixtures; Marigo & Aringer, 2009), and the possibility of performing extremely fast calculations of extended model grids with the COLIBRI code (Marigo et al., 2013).

The combined progress of observations and models opens up the possibility of improving the calibration of the most uncertain processes in TP-AGB models. In this work, we use the population synthesis code TRILEGAL (Girardi et al., 2005) to produce synthetic samples of TP-AGB stars in the SMC. The predicted star counts and luminosity functions (LFs) are compared to the observed ones from the AGB candidate list by Srinivasan et al. (2016, hereafter SR16). Our purpose is to put quantitative constraints on TP-AGB lifetimes as well as on the efficiency of the critical processes of mass loss and third dredge-up.

The structure of the paper is arranged as follows. We first describe the calibration strategy and the adopted sets of observations in Sect. 2. Then, we recall the main properties of our TP-AGB models and their free parameters in Sect. 3. In Sect. 4 we introduce the large grid of new TP-AGB models and the procedure used to identify the best-fitting solutions. Section 5 is devoted to a discussion of the results and the implications derived from the calibration. A summary of the work and the main conclusions are outlined in Sects. 6 and 7.

2 Calibration strategy

The scheme in Fig. 1 illustrates how our calibration cycle operates, and the data and codes involved therein. A fundamental ingredient of our work is the spatially-resolved SFH of the SMC previously derived for many subregions of the SMC (Sect. 2.1). Using the SFH as an input, the TRILEGAL code is used to simulate the stellar populations of every subregion (Sect. 2.3) that is also covered by the TP-AGB catalogues derived from 2MASS and Spitzer observations (Sect. 2.2). These simulations are performed for different sets of TP-AGB models calculated with the COLIBRI code (varying a few key parameters, as detailed in Sect. 3), until a good overall description of the observed $K_{\rm s}$-band LFs and $K_{\rm s}$ vs. J\!-\!\mbox{K_{\rm s}} CMD is reached (Sects. 4 and 5).

2.1 Star Formation History

For many years, the most comprehensive study of the star formation history in the SMC was the spatially-resolved analysis performed by Harris & Zaritsky (2004) for an area of 18 deg2, based on the optical photometry from the Magellanic Clouds Photometric Survey (MCPS). More recently, Rubele et al. (2015, 2018, hereafter R18) recovered the SFH of the SMC using data from the VISTA near-infrared survey of the Magellanic Clouds (VMC, Cioni et al., 2011). Both studies rely on CMD reconstruction methods. The main advantages of using VMC data are i) the lower extinction in the near-infrared passbands, and ii) the photometry reaching the oldest main-sequence turn-off points, ensuring a robust estimate of the ages.

In our study we use the SFH from R18. The most relevant aspects are:

The SFH recovery procedure was applied to 14 SMC tiles, each divided in 12 subregions of $0.143$ deg2 each, for a total area of 23.57 deg2, as illustrated in Fig. 2, and listed in Table 1. 2. 2.

The “partial models”, representing the simple stellar populations used in the SFH recovery, were computed with an updated version of the TRILEGAL code (Marigo et al., 2017), which includes the most recent PARSEC stellar evolutionary tracks111http://stev.oapd.inaf.it/cmd (v1.2S, Bressan et al., 2012, 2015) for evolutionary stages prior to the TP-AGB. The same sets of PARSEC tracks are used in the present work, in addition to the TP-AGB tracks to be discussed later. 3. 3.

For each of the 168 subregions analysed, R18 independently derived the star formation rate (SFR), the age–metallicity relation (AMR), the mean distance modulus and the $V$-band extinction $A_{V}$, together with their $1\sigma$ uncertainties. 4. 4.

The CMD regions used in the SFH recovery are limited to magnitudes \mbox{K_{\rm s}}>12.1 mag and colours -0.52<\mbox{J!-!K_{\rm s}}\ \mathrm{[mag]}<0.88 and -0.82<\mbox{Y!-!K_{\rm s}}\ \mathrm{[mag]}<1.56, therefore they include just a minor number of TP-AGB stars – more specifically, the faint and bluest tails of their distributions in the CMD. Even so, those few TP-AGB stars are largely outnumbered by the RGB stars found in the same dataset at comparable magnitudes. This implies that the different TP-AGB tracks being tested in the present work, would not have affected the derivation of the underlying SFHs in any significant way.

2.2 Observations of AGB stars in the SMC

Our calibration is based on the SMC population of evolved stars identified and classified by Srinivasan et al. (2016, hereafter SR16), based on the program “Surveying the Agents of Galaxy Evolution in the Tidally-stripped, Low Metallicity SMC” (SAGE-SMC, Gordon et al., 2011). The SAGE-SMC photometry includes optical to far-infrared passbands and covers the SMC main body (Bar and Wing) and tail regions. Fig. 2 shows the VMC survey coverage superimposed on the density map of the AGB stars classified by SR16. As can be seen, a large fraction of the area studied by R18 (about 87 per cent) is also covered by the SR16 catalogue.

The first catalogue of evolved stars based on the SAGE-SMC survey was published by Boyer et al. (2011, hereafter B11), followed by the work of SR16 who reconstructed the candidate list using an updated version of the SAGE photometry and optical to mid-IR information from previous studies to refine the candidate selection. Details about their classification method are provided in the Appendix A. In short, SMC stars in SR16 are classified into Red Super-Giant (RSG), Carbon-rich AGB (C-AGB or C-rich), Oxygen-rich AGB (O-AGB or O-rich), anomalous-AGB (a-AGB) and extreme AGB (X-AGB). The a-AGB star sample hosts both C- and O-rich stars that are thought to be low-mass evolved AGB stars (Boyer et al., 2015). The X-AGB sample contains the most dust-enshrouded stars and the majority of them are C-rich stars (see e.g. Sloan et al., 2006). The classification of O-rich and C-rich stars is mainly based on their position in the $K_{\rm s}$ vs. J\!-\!\mbox{K_{\rm s}} CMD, where they draw distinct features. In addition, X-AGB and a-AGB stars are identified with the help of mid-infrared photometry.

The SR16 catalogue also includes the 81 sources spectroscopically classified by Ruffle et al. (2015). In addition, we consider the spectroscopic classification of the 273 sources analysed by Boyer et al. (2015). The star counts for the RSG/AGB populations are summarised in Table 2 and they refer to the area covered by the VMC tiles shown in Fig. 2. Only VMC subregions covered by the SR16 catalogue are considered. As in R18, we also exclude the subregion G6 of tile SMC 6_4, which hosts the Milky Way globular cluster NGC 362. The number of AGB stars present in each tile according to the classification by SR16 is listed in Table 1.

2.3 TRILEGAL simulations

Each subregion covered by both VMC and SAGE-SMC surveys is simulated according to its SFR, AMR, $A_{V}$ and distance as derived from the SFH recovery procedure. Most of the assumptions are identical to those adopted by R18:

• (i)

The adopted initial mass function (IMF) for single stars is that of Kroupa (2001).

• (ii)

Non-interacting binaries are simulated using a binary fraction of 30 per cent and a uniform distribution of mass ratios between the secondary and the primary, ranging from 0.7 and 1.

• (iii)

Our simulations are based on the same set of PARSEC evolutionary tracks used to derive the SFH by R18, for all evolutionary phases prior to the TP-AGB. This ensures an optimal level of consistency between the two works (see also Sect. 2.4).

For the purposes of the present study we build an extended grid of COLIBRI TP-AGB tracks, which are included in all simulations. The main steps for their integration in TRILEGAL are described in Marigo et al. (2017), with some improvements to be detailed in Sect. 3. Our starting point is the set of TP-AGB stellar models described in Rosenfield et al. (2014, 2016), in which the mass loss in COLIBRI was adjusted to reproduce the observed TP-AGB luminosity functions in a sample of nearby galaxies imaged in the optical via the Advanced Camera for Surveys Nearby Galaxy Survey Treasury (ANGST, Dalcanton et al., 2009), and (for a subsample of them) in the near-infrared via a HST/SNAP follow-up campaign (Dalcanton et al., 2012). These works pointed out the need of an efficient pre-dust mass loss, particularly for low-mass stars, in agreement with earlier findings (Girardi et al., 2010). Information about the third dredge-up could not be derived from those data, since the available HST photometry could not discriminate between O- and C-rich stars.

Model quantities are converted into the relevant photometry by means of extensive tables of bolometric corrections, derived from the stellar spectral libraries by Aringer et al. (2009) for C-rich stars, and Castelli & Kurucz (2003) plus Aringer et al. (2016) for O-rich stars.

The effect of the circumstellar dust in mass-losing stars is taken into account following the work of Marigo et al. (2008), to which the reader should refer for all details. The approach couples radiative transfer calculations across dusty envelopes (Groenewegen, 2006; Bressan et al., 1998), with the scaling formalism introduced by Elitzur & Ivezić (2001) and a few key relations from the dust-growth model originally developed by Ferrarotti & Gail (2006). Tables of spectra are computed for a few dust mixtures and then interpolated to produce tables of bolometric corrections. The spectra used for the present simulations are computed for a fixed dust mixture consisting of amorphous carbon (85 per cent) and SiC (15 per cent) for C-rich stars, and silicates for O-rich stars (Groenewegen, 2006) In this work we revise some input prescriptions of the dust treatment: we adjust the abundances of some elements (Si, S, Fe) to follow the scaled-solar pattern of Caffau et al. (2011) in line with the evolutionary tracks, and we compute the condensation degree of carbon dust based on the results from dynamical atmosphere models by Eriksson et al. (2014), in place of a fitting relation based on the results of Ferrarotti & Gail (2006). This modification improves the consistency of our simulations since the same data presented in Eriksson et al. (2014) are also used to predict the mass-loss rates of C-stars during the dust-driven regime (see Sect. 3.2).

To handle the effect of the AGB boosting period at $\sim\!1.6$ Gyr (Girardi et al., 2013), our simulations are split into three age intervals for each subregion (see Table 3). In the younger and older age intervals, the simulations are not affected by the boosting effect, therefore a resolution of 0.02 dex in log(age) suffices to represent the SFR for each subregion. During the boosting period, however, the log(age) resolution is improved to 0.001 dex, in order to sample the fast variations that occur in the production rate of TP-AGB stars. This boosting interval is intentionally kept very wide (0.4 dex), to account for the variation of its mean age with metallicity.

As a more general note, we recall that TRILEGAL is a quite complex population synthesis code able to generate catalogues covering stars over very wide ranges of stellar parameters, and including many physical processes specific to TP-AGB stars, like the temperature and luminosity variations driven by thermal pulses, changes in surface chemical surface and spectral type, and reprocessing of radiation by circumstellar dust shells (see Marigo et al., 2017). However, there are two important effects still not considered in the present work. The first is the long-period variability that characterizes the photometry of most TP-AGB stars. Its effect is quite limited in the infrared filters considered here. Indeed, even Mira variables with large-amplitude variations in the optical have their $K$-band variations limited to just a few tenths of magnitude (see figure 2.44 in Lattanzio & Wood, 2004). The second potentially important effect is the evolution in interacting binaries, which depends on additional parameters like the distributions of orbital separations and mass ratios. While interacting binaries might be important to explain some classes of photometric variables and chemical anomalies in red giants (see De Marco & Izzard, 2017, for a review), in the present work we assume that the fraction of stars that interact before the end of the TP-AGB evolution is negligible. This approach can be justified by the current stage of the TP-AGB model calibration, in which modellers still struggle to have uncertainties in model lifetimes reduced to within the $\lesssim\!20$ per cent level that characterizes other main evolutionary phases. Once this first step of the calibration is successfully completed, interacting binaries will have to be considered.

2.4 SFH uncertainties

First of all, the CMD reconstruction algorithm applied by R18 produces simulations that quite nicely fit the original VMC observations for every subregion of the SMC, as illustrated in their figure 4. Differences between star counts in the data and in their best-fitting models are in general $\la\!10$ per cent across the VMC colour-magnitude diagrams, and concentrated in the immediate vicinity of the red clump. Although such differences could be reflecting small uncertainties in stellar evolutionary models for pre-AGB phases, they are negligible compared to the sizeable heterogeneity in the TP-AGB lifetimes among various models in the literature. Therefore, we can safely conclude that our simulations, using the same pre-AGB models and SFH as in R18, represent a robust starting point for the subsequent calibration of the TP-AGB models.

We verify how the uncertainties in the best-fitting SFR, AMR, distance and extinction impact the simulated LFs and number counts. To this aim, we select the subregions 5, 8, and 9 of the SMC tiles 4_3 and 5_4, which can be considered as representative of most of the SMC regions. For each subregion, we compute two distinct sets of 10 simulations. The first set is computed using the best-fitting values of the SFH. For the second one, the value of the SFH at each age bin is extracted from a Gaussian distribution, centred at the best-fitting value, and with standard deviation equal to the uncertainty associated with the best-fitting solutions222 We use a single Gaussian distribution to describe the symmetric uncertainties in the values of distance and reddening. To describe the asymmetric uncertainties in the SFR and AMR, we use two Gaussian distributions, having the same mean but different standard deviations, $\sigma_{\rm low}$ and $\sigma_{\rm up}$. .

Figure 3 compares the mean $K_{\rm s}$-band LFs for the two sets of simulations. We see that: (i) the standard deviation and the Poisson error in each luminosity bin are equivalent, both in the case of a fixed SFH (as expected) and in the case with randomly extracted SFH, and (ii) the difference in number counts between the two cases is within the Poisson error in each luminosity bin. Additionally, we find that the uncertainties in the SFH produce smaller than 10 per cent variations in the total number of AGB stars. These results justify the use of a simplified approach, in which, for each SMC subregion, the LF and the associated uncertainty in each luminosity bin are obtained as the mean value and the standard deviation from 10 TRILEGAL simulations obtained with the best-fitting values of the SFH. Given the large number of TP-AGB sets we tested, this approach allows us to significantly speed-up the calibration procedure.

2.5 Model selection criteria

To compare our results to the observed catalogue, we need to properly identify the different classes of stars in the synthetic catalogue. One possibility is to use the same photometric criteria as in B11, also applied in SR16. In this case, however, small offsets in the model colours and magnitudes with respect to the observations can significantly mix different kinds of stars, affecting the comparison with the observed star counts in the different AGB populations. On the other hand, the catalogues of synthetic stars contain the surface chemical composition that allows the straightforward separation into different chemical types. For these reasons, we used a ‘hybrid’ selection approach: whenever possible we use theoretical parameters – i.e. C/O ratio, evolutionary stage – to distinguish different populations, complementing the selection with photometric criteria when necessary (i.e. for the X-AGB stars). Our selection criteria are detailed in Appendix A.2.

In addition, we have identified a systematic colour shift between models and observations, which appears to widen at increasing luminosity of the O-rich stars (including RSG and O-AGB stars). In Appendix B we discuss this feature and describe the way models are corrected to reduce its impact on the comparison between data and simulations.

2.6 Metrics for identifying the best-fitting model

Before proceeding with the model calibration, we should better establish what exactly defines a good TP-AGB model in our case: essentially, it means producing the correct lifetimes and luminosity distributions for the main sub-types of TP-AGB stars. In this perspective, we use as primary calibrators the star counts in the $K_{\rm s}$ vs. J\!-\!\mbox{K_{\rm s}} CMD and the $K_{\rm s}$-band LFs. The basic requirement is that the models simultaneously reproduce these quantities for the O-, C-, and X-AGB star samples identified by SR16. Among the many possible criteria to quantify the performance of each simulation, we choose the $\chi^{2}$ of the LF distributions, $\chi^{2}_{\rm LF}$. It is computed following the $\chi^{2}$-like statistics defined by Dolphin (2002), i.e. assuming that the stars are distributed into Hess diagram cells following a Poisson statistics:

[TABLE]

where $n$ and $m$ are the number of simulated and observed stars, respectively, in the $i$-th cell of the Hess diagram or the $i$-th bin of the LF, and $N$ is the number of cells or bins in which both $n$ and $m$ are different from zero.

For each of the three classes of AGB stars and for the entire sample, we compute the value of the $\chi^{2}_{\rm LF}$. The performance of each set of models is evaluated using a total $\overline{\chi^{\,2}_{\rm\,LF}}$, that is the average of the four $\chi^{2}_{\rm LF}$ values weighted by the respective number of observed stars. To identify the best-fitting model, we first select the sets that simultaneously reproduce the total number of observed AGB, O-, C- and X-AGB within 3-$\sigma$. The one with the lowest $\overline{\chi^{\,2}_{\rm\,LF}}$ value is considered as the best-fitting solution. This criterion is applied to the two series of models presented in Sects. 4.3 and 4.4.

3 TP-AGB models

The evolution prior to the TP-AGB phase is computed with the PARSEC code (Bressan et al., 2012), while the TP-AGB phase is followed by the COLIBRI code (Marigo et al., 2013). We refer the reader to this latter paper for a detailed explanation of all TP-AGB calculations. Below we mention just a few details, aimed to provide a quick picture of what is being considered.

3.1 The PARSEC and COLIBRI codes: basic assumptions

PARSEC v1.2S models take into account overshooting from stellar cores, following the Bressan et al. (1981) formalism. The overshooting parameter, $\Lambda_{\mathrm{ov}}$, is assumed to increase from null to its maximum value, 0.5 pressure scale heights ($H_{P}$; measured across the convective border), over a 0.3-$\mathrm{M}_{\odot}$ wide interval that starts at the maximum mass of stars with radiative cores on their main sequence (which is 0.95 to 1.1 $\mathrm{M}_{\odot}$, depending on metallicity). This maximum value of overshooting implies an overshooting region of $\sim 0.25\,H_{P}$ above the convective border (and, equivalently, a $f_{\mathrm{ov}}\sim 0.022$ parameter in diffusive methods), which is among the typical values in the current literature (see Claret & Torres, 2016, 2017, and references therein). This overshooting efficiency determines the maximum masses of the TP-AGB models, referred to as $M_{\mathrm{IM}}$ in Bressan et al. (2012). In practice, PARSEC models of smaller masses develop at least one thermal pulse in the helium shell, after forming an elecron-degenerate carbon-oxygen core. The COLIBRI evolution starts at this point. The $M_{\mathrm{IM}}$ limit depends on metallicity and is generally located between 5 and 6 $\mathrm{M}_{\odot}$ (see figure 3 in Bressan et al., 2012). Stars with slightly higher masses (up to $\sim 9$ $\mathrm{M}_{\odot}$, in models with overshooting), build a degenerate oxygen-neon core before developing any significant thermal pulse. These stars are expected to evolve through a short-lived super-AGB phase (e.g. Siess, 2010), and will not be considered in this paper.

PARSEC models are computed assuming the mixing length theory with an $\alpha_{\mathrm{MLT}}=1.74$ parameter, as determined from the calibration of the solar model (Bressan et al., 2012). COLIBRI adopts the same description for the external convection. The physical process most affected by this choice is the HBB, whose energetics and nucleosynthesis are computed by means of a detailed nuclear network (Marigo et al., 2013). In Sect. 5.4 below will discuss the impact of this model assumption in our calibration.

COLIBRI includes the ÆSOPUS code (Marigo & Aringer, 2009) as an internal routine, to compute the equation of state on-the-fly (for $\simeq 800$ species) and Rosseland mean opacities in the outer regions of the star. Moreover, it accounts for the hot-bottom burning nucleosynthesis and energetics through a complete nuclear network coupled to a diffusive description of convection. Given the uncertainties in the mass loss and 3DU, these processes are treated with parametrized descriptions. The related free parameters need to be calibrated to reproduce the observations (see Sects. 4.3 and 4.4). For any combination of input prescriptions, the code allows for a very quick computation of large grids of TP-AGB tracks. Such a feature makes COLIBRI a suitable tool for the calibration cycle depicted in Fig. 1.

COLIBRI usually starts from a stellar configuration (defined by luminosity, total mass, core mass, and envelope composition) extracted from a PARSEC model just before its first significant thermal pulse (1TP). The subsequent TP-AGB evolutionary track is then joined with the PARSEC one, inside TRILEGAL, to form a continuous track spanning from the pre main-sequence to the end of the TP-AGB. As illustrated in Marigo et al. (2013), there is in general a quite good continuity in evolutionary properties as a given track moves from PARSEC to COLIBRI.

However, PARSEC tracks are computed at constant mass, whereas in this work we allow for the possibility that mass loss occurs before the first thermal pulse. Therefore, we adopt an additional step before joining PARSEC and COLIBRI tracks: along the final sections of the PARSEC tracks (corresponding to Early-AGB stages), we compute the mass-loss rate according to the adopted mass-loss formalism. At each time step, the effective temperature along the PARSEC track is re-evaluated, using the envelope integration routines in the COLIBRI code, together with the corrected stellar mass. The entire procedure is justified by the fact that, as the TP-AGB is approached, the luminosity and the rate of brightening are mainly controlled by the core mass, and are much less sensitive to the envelope mass. The effect of accounting for mass loss during the Early-AGB is discussed in Sect. 4.3.

In the following we outline the main ingredients that are relevant to the physical calibration.

3.2 Mass loss

The treatment of mass loss in AGB models is critical as it controls the duration of the phase, which normally terminates when the whole envelope is ejected into the interstellar medium. The abundant literature on the topic proves that over the years substantial efforts have been devoted to understand the physical mechanisms responsible for mass loss in AGB stars, and more generally in red giants (see Höfner & Olofsson, 2018; Willson, 2000, for extensive reviews). While the key roles of stellar pulsation and dust growth are today universally acknowledged, a comprehensive theoretical framework is still missing. As has been assumed in earlier studies (Girardi et al., 2010; Rosenfield et al., 2014, 2016), we adopt a scheme that considers two regimes of mass loss, which we shortly refer to as pre-dust mass loss (with a rate $\dot{M}_{\rm pre-dust}$) and dust-driven mass loss (with a rate $\dot{M}_{\rm dust}$), respectively.

In our scheme, the pre-dust mass loss occurs during red-giant stages, characterized by relatively low luminosities and high effective temperatures, which set unfavourable conditions to drive a stellar wind by radiation pressure on dust grains. During these phases mass loss should instead be produced by another, still not clearly identified, mechanism. In this work we implement a routine based on the model developed by Cranmer & Saar (2011, hereafter CS11), in which the wind is assumed to be driven by Alfvén waves and turbulence that originate in cool extended chromospheres (see also Schröder & Cuntz, 2005). To this aim we adopt all the equations described in the paper by Cranmer & Saar (2011) related to the cold wave-driven mass loss case (their section 3.2), except for the photospheric parameters (e.g. opacity, mass density, gas pressure) that are obtained directly from the COLIBRI code.

At later stages along the AGB, we assume that the chromospheric mass loss is quenched and the outflow is initially triggered through stellar pulsations (McDonald et al., 2018) and then fully accelerated by radiation pressure on dust grains thanks to the dust–gas dynamical coupling (Höfner & Olofsson, 2018). Several relations for $\dot{M}_{\rm dust}$ as a function of stellar parameters (current mass $M$, luminosity $L$, effective temperature $T_{\rm eff}$, initial metallicity $Z_{\rm i}$) are available in the literature. Here we investigate a few among the most popular ones, i.e. Vassiliadis & Wood (1993, hereafter VW93), which relates the efficiency of mass loss to the pulsation period333The pulsation periods are computed as in the original paper of VW93., and Bloecker (1995, hereafter BL95), which is characterized by a significant dependence on the luminosity. In addition, we use the recent results of dynamical atmosphere models for carbon stars (Mattsson et al., 2010; Eriksson et al., 2014) to predict $\dot{M}_{\rm dust}$ as a function of stellar parameters: mass, luminosity, effective temperature, and carbon excess, $\mathrm{C}-\mathrm{O}$. This latter is defined as:

[TABLE]

where $n_{\rm C}$, $n_{\rm O}$, and $n_{\rm H}$ denote the number densities of carbon, oxygen and hydrogen, respectively. During the TP-AGB evolutionary calculations, if not otherwise specified, the current mass-loss rate is taken as:

[TABLE]

where $\eta_{\rm pre-dust}$ and $\eta_{\rm dust}$ are adjustable efficiency parameters.

Mass loss on the RGB of low-mass stars, instead, is taken into account by simply removing the amount of mass inferred from the Reimers (1975) formula with a multiplicative factor $\eta_{\mathrm{Rei}}=0.2$ (see Miglio et al., 2012), as we pass from the RGB-tip to the initial stage of core-helium burning. At SMC metallicities, only stars with \mbox{M_{\rm i}}\lesssim 1.2 $\mathrm{M}_{\odot}$ turn out to have their RGB masses reduced by more than 5 per cent. Therefore, RGB mass loss does not impact on the TP-AGB populations here simulated in any significant way.

3.3 The third dredge-up

Given the uncertainties and the heterogeneous results that characterize convective mixing in full TP-AGB models, we need to treat the third dredge-up following a parametric scheme as detailed in Marigo et al. (2013). We need to specify three main characteristics of the third dredge-up: i) the onset, ii) the efficiency, and iii) the chemical composition of the intershell material. We now discuss each of these prescriptions in turn.

3.3.1 Onset of the 3DU

Following earlier studies (Marigo et al., 1999; Marigo & Girardi, 2007), in COLIBRI the occurrence of a mixing event is determined by a temperature criterion based on the parameter $T_{\rm b}^{\rm dred}$, which is the minimum temperature that should must be reached at the base of the convective envelope for a mixing event to take place. Typical values lie in the range \mathrm{6.3\la\log(\mbox{T_{\rm b}^{\rm dred}}/\mathrm{K})\la 6.7}, as first shown by Wood (1981). At each thermal pulse, the COLIBRI code performs complete envelope integrations to check whether the temperature criterion is satisfied or not at the stage of the post-flash luminosity maximum. An alternative choice, equivalent to $T_{\rm b}^{\rm dred}$, is the classical parameter $M_{\rm c}^{\rm min}$ that defines the minimum core mass for the onset of the sequence of third dredge-up events. Typical values are usually within the range 0.54\la\mbox{M_{\rm c}^{\rm min}}\,(\mbox{\mathrm{M}_{\odot}})\la 0.60.

In this work we rely on the $T_{\rm b}^{\rm dred}$ parameter. The simplest approach is to adopt a constant value (e.g. \log(\mbox{T_{\rm b}^{\rm dred}}/\mathrm{K})=6.4) in all models. This assumption may be relaxed to include a dependence on metallicity, so that the minimum temperature increases with $Z_{\rm i}$, as indicated by Marigo & Girardi (2007) following the results of their calibration. The effect is to increase the minimum initial mass for the formation of carbon stars at higher metallicity. In our models we test both the case of constant $T_{\rm b}^{\rm dred}$ as well as metallicity-dependent relations.

3.3.2 Efficiency of the 3DU

The efficiency of each mixing event is usually described by the parameter $\lambda=\Delta M_{\rm dred}/\Delta M_{\rm c}$, the fraction of the increment of the core mass during an inter-pulse period that is dredged up during the next thermal pulse. Typical values are in the range $\lambda=0$ (no dredge-up) and $\lambda=1$ (no net increase of the core mass), even though cases with $\lambda>1$ are predicted in models that assume efficient convective overshoot (e.g. Herwig, 2004).

We aim to devise a parametric description of $\lambda$ that a) is qualitatively consistent with full TP-AGB models calculations, and b) keeps the necessary flexibility, through the adoption of free parameters, to perform a physically-sound calibration based on observations. Fig. 4 provides an example of the 3DU efficiency predicted by complete TP-AGB models from a few different authors. During the TP-AGB evolution of a given star, $\lambda$ first increases following the strengthening of thermal pulses, reaches a maximum value $\lambda_{\rm max}$, and then decreases as the envelope mass is reduced by stellar winds, until the 3DU is eventually quenched. The peak value $\lambda_{\rm max}$ is not fixed, but varies from star to star depending on the initial stellar mass. We note, however, that this behaviour is not homogeneously described by existing models. For instance, the models of Karakas et al. (2002, hereinafter also K02), and, in general, all TP-AGB models computed with the Mount Stromlo Evolutionary code (e.g. Karakas et al., 2002; Karakas, 2010; Fishlock et al., 2014, and references therein), reach a high 3DU efficiency ($\lambda\simeq 1$) for $M_{\rm i}\gtrsim 3\,{\rm M}_{\odot}$, regardless of metallicity. In contrast, TP-AGB models by Cristallo et al. (2011, the FRUITY database); Cristallo et al. (2015, the FRUITY database) and by Ventura & D’Antona (2009) are characterized by lower values of $\lambda$, and by a decreasing trend of $\lambda_{\rm max}$ with the initial stellar mass.

To account for this behaviour of the 3DU efficiency, we express $\lambda$ as a function of the ratio between the current total mass and the current core mass, $\xi=M/M_{\rm c}$, a dimensionless parameter that is expected to decrease during the TP-AGB evolution as a consequence of mass loss and core mass growth. A scheme of the new parametric formalism is shown in Fig. 5. We assume a parabolic dependence for $\xi$, that allows to capture the bell-shaped evolution of $\lambda$. The controlling parameters for this relation are the centre $\widetilde{\xi}$ of the parabola, the value $\xi_{\lambda=0}$ at which $\lambda$ drops to zero, and the peak value $\lambda_{\rm max}$. The latter is not treated as a free parameter, but rather expressed as a parabolic function $\lambda_{\rm max}(M_{\rm c})$ of the core mass, to reflect the increase in 3DU efficiency with initial stellar mass, while also allowing for a decrease towards the high-mass range. The controlling parameters, equivalent to the previous case, are $\widetilde{M}_{\rm c}$, $M_{{\rm c},\lambda=0}$, and $\lambda_{\rm max}^{\ast}$.

The explicit form of our parametrization is:

[TABLE]

that is used to assign the 3DU efficiency at each thermal pulse as a function of $M$ and $M_{\rm c}$, provided the $T_{\rm b}^{\rm dred}$ criterion for the occurrence of the 3DU is fulfilled. The physical interpretation of the parameters is the following.

For a given star, $\widetilde{\xi}$ is the value of $M/M_{\rm c}$ at which the maximum $\lambda$ is reached, while $\xi_{\lambda=0}$ is the value at which the 3DU is quenched due the decrease of the envelope mass caused by mass loss. The adopted ranges for these parameters ($\widetilde{\xi}\simeq 3$ – $4$, $\xi_{\lambda=0}\simeq 1.5$ – $2$) are based on full evolutionary calculations available in the FRUITY database (Cristallo et al., 2015). 2. 2.

$\lambda_{\rm max}^{\ast}$ is the maximum efficiency of the 3DU among all TP-AGB stars, that is attained at $M_{\rm c}=\widetilde{M}_{\rm c}$, while $M_{{\rm c},\lambda=0}$ is the value of core mass beyond which the 3DU does not occur. These are the key parameters that we aim at constraining with the aid of the observed AGB stars in the SMC.

3.3.3 Chemical composition of the interhsell

In the COLIBRI code the pulse-driven nucleosynthesis is computed coupling a synthetic description of the intershell with a complete nuclear network which includes the most important $\alpha$-capture reactions. In addition, we have the possibility to explore the impact of various efficiencies of the convective overshoot (refer to Sect. 7.5.5 in Marigo et al., 2013). In the present work we limit to consider the standard case for the intershell composition (i.e. no convective overshoot, oxygen-poor case), with typical abundances (in mass fraction) of helium, carbon and oxygen of $\mathrm{{}^{4}He/^{12}C/^{16}O\approx 0.70-0.75/0.25-0.20/0.005-0.01}$. This pattern may be significantly altered by the inclusion of convective overshoot at the bottom of the pulse-driven convective zone (Herwig, 2004), the main effect being the enhancement of the oxygen abundance. We plan to extend our investigation about the impact of an oxygen-rich intershell to a future study.

4 Calibration along the sequence of TP-AGB sets

In this section we first describe the results of our starting model. We then move to illustrate our calibration presenting the large grid of computed models and their performance in comparison with observations. In a first series of models we explore the effect of different mass-loss prescriptions for both O- and C-rich stars. Starting from the resulting best fitting-model of this series, we calculate a second series of models in which we focus on the calibration of the 3DU.

4.1 The starting set of TP-AGB models

As anticipated in Sect. 2, we start by simulating the SMC population using the TP-AGB tracks described in Rosenfield et al. (2014, 2016) and presented in Marigo et al. (2017). We will refer to this initial set as S_00. The mass-loss prescription for the pre-dust driven wind is a modified version of the Schröder & Cuntz (2005, hereafter SC05) prescription (Rosenfield et al., 2016). The dust-driven phase is described with a mass-loss formalism similar to Bedijn (1988), for both O-rich and C-rich stars. The 3DU onset is controlled by a constant temperature parameter \log(\mbox{T_{\rm b}^{\rm dred}}/\mathrm{K})=6.40 and the efficiency parameter $\lambda$ follows the relations provided by K02.

Fig. 6 compares the observed $K_{\rm s}$-band LFs for each class of AGB stars with the simulated ones based on the set S_00. As shown in the top panel, the total number of AGB stars is reproduced to within 10 per cent, a result that supports the earlier calibration performed by Rosenfield et al. (2016). However, from a closer inspection we note that the model underpredicts the number of stars at magnitudes brighter than \mbox{K_{\rm s}}\approx 11.6 mag, and overpredicts the number of stars at fainter magnitudes. The deficit of predicted stars at brighter magnitudes corresponds mainly to an underestimation ($\approx 30$ per cent) of C-rich stars (third panel from top). The excess of predicted stars at fainter magnitudes corresponds to an overestimation ($\approx 40$ per cent) of O-rich AGB stars with initial stellar mass in the range 1–1.5 $\mathrm{M}_{\odot}$ and low metallicity, \mbox{Z_{\rm i}}\approx 0.001 (second panel from top). A slight excess is also present at the bright end of the LF of O-rich stars. Our analysis indicates that such stars are experiencing HBB, with initial masses larger than $\approx$ 3 $\mathrm{M}_{\odot}$ and initial metallicity \mbox{Z_{\rm i}}\simeq 0.003–0.004.

The total number of X-AGB stars is in agreement with the observations within $\approx 15$ per cent (bottom panel), but the shape of the LFs shows discrepancies around \mbox{K_{\rm s}}\approx 12 mag and \mbox{K_{\rm s}}\approx 11 mag. In particular the peak of the simulated LF is located at fainter $K_{\rm s}$ magnitudes compared to the observed one. The difference is likely a consequence of the adoption of the Vassiliadis & Wood (1993) prescription for the mass loss during the super-wind phase. With that prescription, mass-loss rates reach a maximum ($10^{-5}$–10^{-4}~{}\mbox{\mathrm{M}{\odot}}\rm{yr}^{-1}) and then remain almost constant during the final thermal pulses, causing significant dust extinction in the $K_{\rm s}$-band and an accumulation of stars at \mbox{J!-!K{\rm s}}\ga 4 mag and \mbox{K_{\rm s}}\approx 12 mag.

Before leaving this section we consider it worth commenting on the sample of red super-giants present in the SMC data. The comparison between observed and simulated LFs is shown in Fig. 7. The overall shape of the observed LF of the RSG stars is well reproduced by the starting model. However, there is an evident lack of stars at faint magnitudes. Such discrepancy, already pointed out by Melbourne et al. (2012), can be due to a number of reasons: (i) a residual foreground contamination in the data up to $\approx$ 35 per cent (Boyer et al., 2011, their Sect. 2.2); (ii) an underestimation of the SFR for these young populations; (iii) an underestimation of predicted lifetimes for core-helium burning stars with intermediate and high initial masses. Since the present calibration does not rely on this class of stars, we will only focus our analysis on AGB stars, postponing the investigation of the RSGs to a future work.

4.2 The new grid of TP-AGB tracks

Starting from the results obtained with the initial S_00 set, we put in operation a calibration cycle along which the mass-loss prescriptions and the third-dredge up parameters are varied until a satisfactory reproduction of the observations is achieved. The calibration procedure relies on a large grid of TP-AGB tracks computed with COLIBRI, adopting different descriptions for the third dredge-up and mass loss. Each combination of the adopted parameters corresponds to a set of tracks. Each set covers a range of initial masses from about 0.5 $\mathrm{M}_{\odot}$ to 5–6 $\mathrm{M}_{\odot}$ (typically 70 values) and spans a wide metallicity interval, from \mbox{Z_{\rm i}}=0.0005 to \mbox{Z_{\rm i}}=0.03 (typically 10 values). Therefore, each set includes roughly 700 TP-AGB tracks, computed up to the complete ejection of the envelope. For the purposes of this work the grid of evolutionary TP-AGB models consists of 33 sets ($\sim$ 23000 TP-AGB tracks). Such a demanding computational effort can only be achieved with a flexible and fast code such as COLIBRI. The full grid of TP-AGB tracks is summarised in Table 4.

The observed and simulated number counts for each class of AGB stars, the corresponding $\chi^{2}_{\rm LF}$ and the C/M ratio, i.e. the number of C-type (C-rich and X-AGB) to M-type (O-rich) AGB stars, are given in Table 5. All the LFs calculated for each set are presented in Appendix C available online. In Fig. 8 we show a summary of the performance of the calculated models and we illustrate the fundamental steps of the calibration in the next section.

To facilitate the presentation of the calibration it is convenient to split the TP-AGB sets into two groups, depending on the assumptions adopted for the two key processes under examination, namely mass loss (see Sect. 4.3) and 3DU (see Sect. 4.4).

4.3 Changing mass-loss

The first series of models (S_01 – S_07) shares the same prescriptions for the 3DU. In particular, the efficiency $\lambda$ follows the formalism introduced by K02, which assigns quite high values ($\lambda\simeq 1$) to stars with initial masses larger than $\simeq 3$ $\mathrm{M}_{\odot}$.

The initial step is to reduce the excess of O-rich stars present in our starting set S_00. This goal requires increasing the efficiency of mass loss both during the pre-dust and the dust-driven regimes. This is partly achieved with the adoption of the CS11 formalism to describe the cool chromospheric mass loss during the early stages when dust is expected not to be the main driver of stellar winds. In fact, the pre-dust mass loss is already active prior to the 1TP, i.e. during the brightest part of the early-AGB. For instance, in the case of a 1 $\mathrm{M}_{\odot}$ model with \mbox{Z_{\rm i}}=0.001 the current mass at the beginning of the TP-AGB phase is reduced by $\approx 15$ per cent. The net effect is the shortening of the duration of the TP-AGB phase, mostly for low-mass O-rich stars, thus reducing their numbers. All S_01 – S_07 sets exhibit a reduction in the number of O-rich stars compared to the set S_00, improving the comparison with observations.

To explore the impact of mass loss in the dust-driven regime, the set S_01 assumes a very efficient mass loss, adopting BL95 with \mbox{\eta_{\rm dust}}=0.05 for all models. However, the reduction of O-rich stars is accompanied by a further deficit in the predicted number of C-rich stars, which worsens the comparison with observations compared to S_00.

In all sets from S_02 – S_07 we vary the mass-loss laws in the dust-driven regime, depending on the surface C/O ratio: for C/O $<1$ we use the relation proposed by BL95, while for C/O $>1$ we use a routine based on the recent results of dynamical atmosphere models (Eriksson et al., 2014; Mattsson et al., 2010, hereinafter CDYN). All these models adopt the same CDYN prescription for $\dot{M}_{\rm dust}$ when C/O $>1$, while the efficiency of the BL95 relation, active for C/O $<1$, is increased from \mbox{\eta_{\rm dust}}=0.02 (for S_02) to \mbox{\eta_{\rm dust}}=0.06 (for S_07). Interestingly, this variation results in a relatively modest effect on the total number of O-rich stars, but has a significant impact on the predicted C-rich LF. This effect can be appreciated in Fig. 8 which shows that the total number of O-rich stars do not vary significantly along the sequence from S_02 to S_05, whereas the total number of C-rich (and X-AGB) stars progressively decreases (see also Figs. C.3, C.4, C.5 and C.6 in the online Appendix C).

The reason of the small effect on the number of O-rich stars is that the bulk of them have low-mass progenitors, i.e. \mbox{M_{\rm i}}\lesssim 2 $\mathrm{M}_{\odot}$, for which the duration of the AGB phase is mainly controlled by the cool chromospheric mass loss active before the TP-AGB phase. Due to its high luminosity dependence, the BL95 relation mainly affects the evolution of more massive (and brighter) AGB stars that experience HBB. These represent a small contribution to the total number of O-rich stars (see Sect. 5.4 below).

On the other hand, the sizeable impact on the number of C-rich stars can be explained as follows. The activation of the CDYN $\dot{M}_{\rm dust}$ takes place once specific physical thresholds are met, i.e. sufficiently low effective temperatures, enough carbon excess to be condensed into dust grains, and suitably large L/M ratios. In general, increasing the efficiency of mass loss during the O-rich phase tends to favour larger mass-loss rates also during the stages with C/O $>1$, since TP-AGB stars will enter the C-rich phase with lower current masses and lower effective temperatures.

However, we should also note that extremely efficient mass loss during the O-rich stages may even inhibit the formation of C-stars by either anticipating the termination of the TP-AGB phase, or preventing the star from attaining the minimum temperature, $T_{\rm b}^{\rm dred}$  for the onset of the 3DU if the envelope mass has been significantly reduced. Indeed, this circumstance is met in the evolutionary models of the set S_07, so that the maximum initial mass for the formation of C-stars at \mbox{Z_{\rm i}}=0.004 is found to be around 2.6 $\mathrm{M}_{\odot}$ (see Sect. 5.2).

All LFs computed with an efficiency \mbox{\eta_{\rm pre-dust}}=2 (as in the sets S_02 – S_05) show an excess ($\approx 20$ per cent) of low-mass O-rich stars that populate the faintest bins (\mbox{K_{\rm s}}\gtrsim 11.8 mag). This leads us to increase the $\eta_{\rm pre-dust}$ in the CS11 prescription, passing from \mbox{\eta_{\rm pre-dust}}=2 to \mbox{\eta_{\rm pre-dust}}=3, as in sets S_06 and S_07. As can be appreciated in Table 5 and in Fig. 8, this modification is successful and the residual excess of faint O-rich stars is removed. Models S_06 and S_07 only differ in the choice of the $T_{\rm b}^{\rm dred}$ parameter, which increases linearly with the initial $Z_{\rm i}$ in the case of S_07 (see notes of Table 4). This assumption is equivalent to reducing the formation of C-stars at higher metallicity, as indicated by observations of the C/M ratio in galaxies (e.g. Boyer et al., 2013; Cioni & Habing, 2003; Groenewegen, 2002), and by earlier population synthesis studies (Marigo & Girardi, 2007; Marigo et al., 1999). Anyhow, we note that at SMC metallicities the differences in the $T_{\rm b}^{\rm dred}$ parameter are negligible, becoming significant only at lower and higher $Z_{\rm i}$.

In addition to the BL95 relation for mass loss, we run a set of models (S_19) in which we adopt the formalism of VW93 for the efficiency of dust-driven winds during the O-rich stages, while keeping the CDYN prescription for C-stars. The LF of O-rich stars is recovered very well, with a relative difference in the number counts within $\sim 3$ per cent. However, the same model leads to an excess of C-stars (in the brighter bins, \mbox{K_{\rm s}}\la 10.5 mag) of about 25 per cent.

Among all the sets of this first series we identify the best-fitting model as the set S_07, for which we obtain the lowest $\chi^{2}_{\rm LF}$ values. The resulting LFs are shown in Fig. 9. This set recovers quite well both the LFs of all classes of AGB stars and the number counts. The optimal combination of parameters includes: CS11 pre-dust mass loss with efficiency \mbox{\eta_{\rm pre-dust}}=3, dusty regime with BL95 and \mbox{\eta_{\rm dust}}=0.06 for \mbox{{\rm C/O}}\leq 1 and the CDYN for \mbox{{\rm C/O}}>1; 3DU efficiency described following K02 (see Table 4).

4.4 Changing the third dredge-up

The first part of our work shows that, adopting our standard prescriptions for the 3DU, based on Karakas et al. (2002) models, the AGB star population in the SMC can be well reproduced with suitable choices for the mass-loss regimes. However, at this stage we cannot consider the calibration work complete. Given the lack of a robust treatment of stellar convection and mixing in stellar evolution codes it is not sufficient to limit our investigation to a single 3DU prescription.

As already mentioned in Sect. 3.3 and shown in Fig. 4, the predicted properties of the 3DU are quite heterogeneous in the literature. The largest differences show up at larger stellar masses, typically for \mbox{M_{\rm i}}\ga 3 $\mathrm{M}_{\odot}$. In this mass range the K02 formalism predicts a very efficient 3DU, $\lambda\simeq 1$, for all models with $M>3$ $\mathrm{M}_{\odot}$, irrespective of metallicity. TP-AGB evolutionary calculations by other authors predict lower efficiencies of the 3DU, with typical values $\lambda\la 0.5$ for $M_{\rm c}>0.8$ $\mathrm{M}_{\odot}$ (Cristallo et al., 2015; Ventura & D’Antona, 2009).

We should note that hints for a reduced efficiency of the 3DU in intermediate-mass AGB stars ($\lambda<1$ for \mbox{M_{\rm i}}>3-4\,M_{\odot}) come also from studies on the initial-to-final mass relation (IFMR) derived from white dwarfs (WD) data in Galactic open clusters (see Fig. 12). As already pointed out by Kalirai et al. (2014), current TP-AGB models tend to predict WD masses that are too low for initial masses $\ga 3$ $\mathrm{M}_{\odot}$. A significant improvement is obtained by reducing the efficiency of the third dredge-up during the TP-AGB phase of progenitor stars. In fact the 3DU has a direct impact on the evolution of the core mass reducing its net growth (e.g. Kalirai et al., 2014; Marigo et al., 2013; Herwig, 2004). Therefore, lowering the 3DU in intermediate-mass stars may help to obtain more massive WDs.

In the following we will focus our discussion on those sets that share the same description for the dust-driven mass-loss regime. This relies on two alternatives depending on the C/O ratio, namely: the BL95 relation for O-rich stars and the CDYN option for C-rich stars. This latter represents our preferred option as it is based on well-tested and state-of-the-art dynamical atmosphere models for C-stars (Eriksson et al., 2014; Mattsson et al., 2010).

Given these premises, we now analyse the impact of varying the 3DU description in our simulations. As an initial step we explore the effect of keeping the K02 relations for $\lambda$ while lowering $\lambda_{\mathrm{max}}$  down to a fixed value of 0.4-0.5 (S_08, S_11, S_31), i.e. we impose that the 3DU efficiency does not exceed the selected $\lambda_{\mathrm{max}}$in all models. By analysing the C-rich and X-AGB LFs we find that, unless we assume a very efficient mass loss (\mbox{\eta_{\rm dust}}=0.06 as in S_08) during the O-rich dusty regime, the models predict an excess of bright C-rich stars (\mbox{K_{\rm s}}\lesssim 10 mag), most evident in S_11 and S_31.

To carry out a more systematic exploration on the effects of varying the 3DU efficiency, we relax our standard choice, based on K02, and adopt the new formalism, introduced in Sect. 3.3 (sets from S_13 to S_35). Its parametric formulation enables us to investigate various choices of the 3DU law, including also the possibility that the efficiency of the mixing episodes become lower at larger stellar masses. Since we have already explored the cases of extremely efficient mass loss for O-rich stars (first series of models in Sect. 4.3), we now opt to reduce $\eta_{\rm dust}$ down to 0.01-0.03, values that are also adopted in other widely used AGB models (e.g. Karakas et al., 2018a; Ventura & D’Antona, 2009).

In Fig. 10 we show the resulting carbon star luminosity functions (CSLFs) for a few sets, together with the evolution of $\lambda$ as a function of the core mass $M_{\mathrm{c}}$ for some selected TP-AGB evolutionary tracks at \mbox{Z_{\rm i}}=0.004, representative of the initial metallicity for the bulk of C-rich stars. The emerging picture is quite complex due to the interplay of the 3DU with mass loss, still we can extract a few key general indications that help the interpretation of the results.

• (i)

The faint tail contains information about the lowest-mass stars that become C-stars. This limit is controlled by the temperature parameter $T_{\rm b}^{\rm dred}$, as well as by $\widetilde{M}_{\mathrm{c}}$. As illustrated in panel (a) of Fig. 10, lowering $T_{\rm b}^{\rm dred}$, i.e. moving from set S_15 to set S_16, favours an earlier formation (i.e. at lower $M_{\mathrm{c}}$ and $L$) of C-stars with low masses and therefore increases the predicted number of C-stars in the faintest bins of the CSLF. Decreasing $\widetilde{M}_{\mathrm{c}}$ has a similar effect on the faint tail, as shown in panel (b) of Fig. 10.

• (ii)

The position of the peak is mainly affected by the choice of $\widetilde{M}_{\mathrm{c}}$ as shown in panel (b) of Fig. 10. At lower $\widetilde{M}_{\mathrm{c}}$, i.e. moving from set S_22 to set S_26, the CSLF peak shifts towards fainter magnitudes. At the same time, lowering the value of $\widetilde{M}_{\mathrm{c}}$ results in an increase of the total number of C-stars, in particular those that populate the faint tail of the CSLF.

• (iii)

The amplitude of the peak is mainly controlled by the choice of $\lambda_{\rm max}^{\ast}$. As shown in panel (c) of Fig. 10 increasing the maximum 3DU efficiency, i.e. moving from set S_18 to set S_23, produces a higher peak in the CSLF.

• (iv)

The bright end contains information about the highest-mass stars that become C-stars, and is controlled by $M_{\mathrm{c},\lambda=0}$ together with the efficiency of mass loss. As illustrated in panel (c) of Fig. 10, increasing the value of $M_{\mathrm{c},\lambda=0}$ , i.e. moving from set S_23 to set S_02, causes the appearance of C-stars at brighter magnitudes. It essentially reflects the larger $\lambda$ attained by stars with \mbox{M_{\rm i}}\ga 2.5 $\mathrm{M}_{\odot}$in the present models.

We may conclude that with a very efficient 3DU in more massive AGB stars ($\lambda\approx 1$, as in K02) models tend to largely over-predict the number of bright C-stars (rightmost plot of Fig.10, panel (c), for the set S_02). As already discussed in Sect. 4.3, such excess can be removed only by invoking a powerful mass loss during the preceding O-rich stages (e.g., as for S_07 in the first series of sets).

Among all sets of the second series, the best-fit is provided by S_35 which is found to reproduce well all LFs (see Fig. 11) and the corresponding number counts within 5 per cent. The optimal combination of parameters includes: CS11 pre-dust mass loss with efficiency \mbox{\eta_{\rm pre-dust}}=3, dusty regime with BL95 and \mbox{\eta_{\rm dust}}=0.03 for \mbox{{\rm C/O}}\leq 1 and the CDYN for \mbox{{\rm C/O}}>1; $\lambda_{\rm{max}}^{\ast}=0.7$, $\widetilde{M}_{\mathrm{c}}$ = 0.6, (see Table 4 and panel (d) of Fig. 10).

In the next section we discuss the results of the two best-fitting models S_07 and S_35 and the main differences with respect to the starting set S_00.

5 The best-fitting models

Our final results are summarised in Table 5 and Fig. 8. Models S_07 and S_35 stand out clearly by their smaller $\chi^{2}_{\rm LF}$ and indeed provide an excellent description of the observed CMDs (see Fig. 17) and LFs (see Figs. 9 and 11). In the following we discuss the main implications we may derive from the analysis of the results. The main differences between the two sets deal with the chemical composition of the ejecta and the final core masses. While a proper analysis on the predicted chemical yields is postponed to a follow-up work (Marigo et al., in prep.), we briefly discuss the IFMR in the next Sect. 5.1.

5.1 The IFMR relation

Two recent semi-empirical IFMR together with model predictions from the two best-fitting sets are shown in Fig. 12.

While the general trend is satisfactory for both sets S_07 and S_35, at larger initial masses, $M>3M_{\odot}$, the S_07 models are only marginally consistent with the semi-empirical IFMRs. The WD masses for intermediate-mass stars (\mbox{M_{\rm i}}>3\,\mbox{\mathrm{M}_{\odot}}) appear to be underestimated by the models, in comparison with both relations presented by Cummings et al. (2018) and El-Badry et al. (2018). The predicted IFMR based on the set S_35 improves the comparison, mostly with respect to the initial masses larger than about 3 $\mathrm{M}_{\odot}$. The larger WD masses predicted by these models, compared to those of the set S_07, result from reducing the efficiency of both the third dredge-up and mass loss during the TP-AGB phase of the progenitors. We also note that the agreement becomes quite good with the semi-empirical relation based on Gaia DR2 data (El-Badry et al., 2018), whereas a substantial discrepancy still affects the comparison with the IFMR from Cummings et al. (2018). This latter runs systematically above the mean relation from Gaia data. The reason for such a large difference between the two semi-empirical relations is not known at present. More work is required to clear up this point, given its importance for the calibration of the AGB stellar evolution models.

5.2 The initial masses of C-stars

Figure 13 shows the predicted ranges of initial masses and bolometric magnitudes of C-stars for the starting set S_00 and the two best fitting sets S_07 and S_35.

First, it is striking that in both S_07 and S_35 models, C-stars are expected to form only in a limited interval of initial masses, between 1.4 and 2.8 $\mathrm{M}_{\odot}$, at SMC-like metallicities (i.e. \mbox{Z_{\rm i}}=0.004). In particular, the upper limit for C-stars formation at around 3 $\mathrm{M}_{\odot}$ is mainly constrained by the observed deficit of C-rich stars for magnitudes brighter than \mbox{K_{\rm s}}\simeq 10 mag. Our calibration indicates that C-stars with initial masses $\gtrsim 3$ $\mathrm{M}_{\odot}$ would form a bright red tail that is not actually observed in the 2MASS CMD.

The two best-fitting sets converge to the same threshold mass following somewhat different paths. In the case of S_07, the formation of massive C-stars is prevented because of a high mass-loss efficiency that drastically shorten the lifetimes, despite the efficient 3DU based on K02 formalism. In the case of S_35, the reduced 3DU efficiency – assumed for TP-AGB models of higher masses – contributes to confining the formation region of C-stars to lower initial masses. In both sets of models the upper mass limit for C-stars formation is controlled by the efficiency of mass loss and 3DU, rather than by the onset of HBB.

The results of our calibration need to be compared with the predictions available in the literature. Past works that try to reproduce the CSLF in the SMC present much wider ranges of masses for the formation of C-stars. For instance, Groenewegen & de Jong (1993) and Marigo & Girardi (2007) find C-stars still forming at initial masses as high as 4 $\mathrm{M}_{\odot}$. Also, recent AGB models, more focused on the chemical yields, predict that the maximum initial mass for C-star formation at SMC-like metallicity lies between between 4.0 and 4.5 $\mathrm{M}_{\odot}$ (Karakas et al., 2018b, \mbox{Z_{\rm i}}=0.028), or up to 6 $\mathrm{M}_{\odot}$ (Cristallo et al., 2015, \mbox{Z_{\rm i}}=0.003). On the other hand, AGB models presented by Ventura & D’Antona (2009) and Dell’Agli et al. (2015) form carbon stars with initial masses up to \simeq 3\,\mbox{\mathrm{M}{\odot}} at \mbox{Z{\rm i}}=0.004, a value that is in close agreement with our calibration.

It is also interesting to compare these predictions with the cluster data in the Magellanic Clouds. According to the classical compilation by Frogel et al. (1990), ‘the youngest clusters in which C-stars are found have an age of about 100 Myr implying a maximum initial mass for these stars of 3–5 $\mathrm{M}_{\odot}$’. This conclusion appears in contradiction with our results; however, it was based on quite crude age estimates for the clusters, and uncertain memberships for the C-stars. According to the revised compilation by Girardi & Marigo (2007), the youngest LMC cluster to contain a C-star in its central region is NGC 1850, which has isochrone ages of $\sim\!60$–$90$ Myr (Correnti et al., 2017) and hence turn-off masses $\gtrsim 5.5$ $\mathrm{M}_{\odot}$. This high turn-off mass appears compatible with AGB stars that experience efficient HBB, a process that usually prevents the formation of C-stars; therefore the C-star in NGC 1850 is usually regarded as a chance alignment of a field LMC star. The second youngest LMC clusters to contain C-stars are NGC 1987 and NGC 2209, with one and two C-stars, respectively. NGC 1987 has estimated ages between $\sim\!1$ (Goudfrooij et al., 2017) and 1.3–1.5 Gyr (Milone et al., 2009), while NGC 2209 has ages of $\sim\!1.15$ Gyr (Correnti et al., 2014). Both clusters therefore have turn-off masses slightly below 2 $\mathrm{M}_{\odot}$. We remark that NGC 2209 was regarded as a much younger cluster (ages between 120 and 370 Myr, hence turn-off masses between 4 and 2.7 $\mathrm{M}_{\odot}$) in the work by Frogel et al. (1990). In the SMC, the revised cluster ages by Glatt et al. (2008) indicate that the youngest such cluster is the $\sim\!1.5$-Gyr old NGC 419, with its impressive population of ten C-stars. The presence of a double red clump in this cluster (Girardi et al., 2009) firmly points to turn-off masses of $\sim\!1.75$ $\mathrm{M}_{\odot}$ (Girardi et al., 2013). Overall, the revised cluster data do not appear to indicate the presence of C-stars with progenitors more massive than about 2 $\mathrm{M}_{\odot}$. On the other hand, Magellanic Clouds star clusters with turn-off masses between 2 and 3 $\mathrm{M}_{\odot}$ are rare and more scarcely populated than older clusters, so that the maximum turn-off mass giving origin to C-stars cannot be identified with confidence from these data alone.

Fig. 13 shows other aspects that are worth note. The faintest transition to the C-star domain happens at \mbox{M_{\rm i}}\simeq 1.8~{}\mbox{\mathrm{M}{\odot}}. This value coincides (at \mbox{Z{\rm i}}=0.004) with the initial mass boundary between stars that develop degenerate He-cores after the main sequences and the more massive ones that skip electron degeneracy. As discussed by Lebzelter et al. (2018), the presence of a well-defined minimum in the M-to-C transition luminosity at increasing stellar mass explains the splitting of the data along two branches (labelled (a) and (c) in fig. 3 of Lebzelter et al., 2018) in the diagram that combines 2MASS and Gaia photometry for long period variables.

We expect that some fraction of TP-AGB stars may appear below the RGB-tip (at \mbox{M_{\rm bol}}\approx-3.5 mag), mostly during the relatively long-lived low-luminosity dip that follows the occurrence of a thermal pulse (see e.g., Boothroyd & Sackmann, 1988). These faint stars may be either O-rich or C-rich. The O-rich stars below the RGB-tip span a range of initial masses 1.0~{}\mbox{\mathrm{M}{\odot}}\la\mbox{M{\rm i}}\la 2.2~{}\mbox{\mathrm{M}{\odot}} and cover a fraction of $\approx$ 26-27 per cent of the entire O-rich TP-AGB population, for both sets S_07 and S_35. The C-rich stars below the RGB-tip should correspond to initial masses in the range 1.5~{}\mbox{\mathrm{M}{\odot}}\la\mbox{M_{\rm i}}\la 2.0~{}\mbox{\mathrm{M}_{\odot}}. Our simulations indicate these faint C-stars represent $\approx$ 7 per cent of the entire C-rich population, for both sets S_07 and S_35.

5.3 Lifetimes

TP-AGB lifetimes are relevant to quantify the contribution of low- and intermediate-mass stars to i) the chemical enrichment of the interstellar medium, and ii) the integrated light of galaxies (Marigo, 2015). In fact, the duration of the TP-AGB phase, $\tau_{\rm TP-AGB}$, determines the level of chemical enrichment of the ejecta, which is set through the number of thermal pulses and 3DU episodes experienced by a star and/or the period in which HBB is operating. At the same time $\tau_{\rm TP-AGB}$ controls the energy emitted during the phase through the relation $E_{\rm TP-AGB}=\int_{\tau_{\mathrm{TP-AGB}}}L(t)dt$, where $L(t)$ is the stellar luminosity at time $t$.

In Fig. 14, we present the comparison between the predicted lifetimes of M- and C-stars as derived from the starting set of TP-AGB models S_00 and the two best-fitting models S_07 and S_35. We consider only models with luminosity higher than $\log(L/L_{\odot})=3.3$, i.e. brighter than the RGB tip.

As in other models in the literature (e.g. Weiss & Ferguson, 2009; Karakas, 2010), the predicted TP-AGB lifetimes peak at \mbox{M_{\rm i}}\simeq 2 $\mathrm{M}_{\odot}$, close to the value below which stars develop a degenerate helium core after the main sequence, the exact value depending on the initial metallicity and other model details. The general trends of the lifetimes for O-rich stars are similar in the three set of models, yet some differences are present for \mbox{M_{\rm i}}\gtrsim 3 $\mathrm{M}_{\odot}$. Compared to set S_00, the best-fitting models S_07 and S_35 predict shorter lifetimes by a factor of $\approx 2$. As to C-stars, we note an increase of the lifetimes in the new models compared to set S_00, reaching a factor of 2 near the peak at \mbox{M_{\rm i}}\sim 2 $\mathrm{M}_{\odot}$.

Within the metallicity range covered by the calibration (0.001\la\mbox{Z_{\rm i}}\la 0.004), we expect shorter O-rich stages at decreasing $Z_{\rm i}$. This is the consequence of the earlier transition to the C-rich regime. The lifetimes for C-rich stars also follow a decreasing trend with lower $Z_{\rm i}$, particularly evident for initial masses close to the peak and larger than $\approx 2$ $\mathrm{M}_{\odot}$. The reason is linked to the fact that, at a given initial mass, a TP-AGB star of lower metallicity is characterized by a larger core mass. As a consequence, the TP-AGB phase proceeds at higher luminosities and with shorter inter-pulse periods, i.e. thermal pulses and 3DU episodes take place more frequently with consequent quicker increase of the carbon excess during the C-star stages. Both factors (related to luminosity and carbon excess) tend to favour an earlier onset of the dust-driven wind in C-stars of lower metallicity, in the range under consideration.

5.4 The overall impact of HBB

As already mentioned, our models assume a constant and solar-calibrated mixing-length parameter, with $\alpha_{\mathrm{MLT}}=1.74$, which turns out to affect the efficiency of HBB in the models. Since there are significant claims that such a convection treatment could be inaccurate and that $\alpha_{\mathrm{MLT}}$ (or the equivalent parameter coming from different approaches) may vary along the sequences of red giants (e.g. Tayar et al., 2017, and references therein), let us verify the impact of such an assumption in our results.

First of all, it turns out that adopting $\alpha_{\mathrm{MLT}}=1.74$ in our best-fitting S_07 and S_35 sets, stars with initial masses larger than 3 $\mathrm{M}_{\odot}$  represent a little contribution to the entire AGB star population of the SMC, $\simeq 8$ and $\simeq 9$ per cent respectively. Therefore we expect that, overall, stars in this mass range should have a small impact in our calibration.

To test this point in detail we run an additional set of AGB models over the relevant range of metallicity, from \mbox{Z_{\rm i}}=0.001 to \mbox{Z_{\rm i}}=0.008, keeping the same input prescriptions as in the set S_35, except for the mixing-length parameter which is changed from $\alpha_{\mathrm{MLT}}=1.74$ to $\alpha_{\mathrm{MLT}}=2.00$ in the models with initial masses >3~{}\mbox{\mathrm{M}{\odot}}. The main effects of increasing $\alpha_{\mathrm{MLT}}$ are those of i) lowering the minimum stellar mass for the activation of HBB and ii) making stars with HBB hotter and more luminous, as can be seen in Fig. 15. Given the significant luminosity dependence of the B95 mass-loss formula, models with $\alpha_{\mathrm{MLT}}=2.00$ have a shorter duration of the AGB phase. At \mbox{Z{\rm i}}=0.004 the reduction of the TP-AGB lifetimes does not exceed $\simeq 30$ per cent in most cases. The impact on star counts is relatively modest, so that AGB stars with \mbox{M_{\rm i}}>3~{}\mbox{\mathrm{M}_{\odot}} make an overall contribution of $\simeq 6$ with $\alpha_{\mathrm{MLT}}=2.00$, in place of $\simeq 9$ predicted with $\alpha_{\mathrm{MLT}}=1.74$. Therefore, this test confirms that uncertainties in the efficiency of HBB, due to reasonable changes of the mixing-length parameter, do not affect the results of the present calibration.

5.5 Integrated luminosities

One of the most important motivations for calibrating TP-AGB evolutionary tracks, resides in their contribution to the integrated light, or spectral energy distribution (SED), of unresolved galaxies. Indeed, integrated spectra are often the only source of information about the total stellar masses, mean ages and metallicities of galaxies. Since Maraston (2005), TP-AGB stars have been recognized as a crucial contributor to the light of galaxies, especially at (rest-frame) infrared wavelengths, and for galaxies containing $\sim 1$-Gyr populations, at which this contribution is expected to peak. A few subsequent works attributed a much smaller importance to TP-AGB stars (e.g. Conroy & Gunn, 2010; Kriek et al., 2010; Zibetti et al., 2013), lowering down their contribution to values similar to those presented in the classical population synthesis models by Bruzual & Charlot (2003). Entering into the details and consequences of this ongoing debate (e.g. Maraston, 2013) is far beyond the scope of this paper; anyhow let us check what amount of TP-AGB contribution is favoured by our calibration.

Fig. 16 displays the expected contribution of TP-AGB stars to the total integrated bolometric luminosity and the $K$-band luminosity as a function of age for a stellar population of $Z_{\rm i}$= 0.004 for the initial set S_00 and the two best-fitting models S_07 and S_35. We also show the contribution of TP-AGB stars as derived by Bruzual & Charlot (2003, private communication) and Maraston (2005, private communication) for the same metallicity444Maraston (2005) does not include the results for metallicity $Z_{\rm i}$= 0.004. We used a linear interpolation between the models with $Z_{\rm i}$= 0.001 and $Z_{\rm i}$= 0.02..

We compute the integrated luminosities assuming an age interval of $\log(\mathrm{age/yr})=7.5-10.5$, with a step of 0.001. The spike visible at $\log(\mathrm{age/yr})\sim 9.2$ is caused by the TP-AGB boosting effect, already mentioned in Section 2.3 and fully described by Girardi et al. (2013).

The results for the two best-fitting models are very similar. They differ mainly in the contribution of the youngest TP-AGB stars ($\log(\mathrm{age/yr})\lesssim 8.5$), which is slightly lower in S_07. The lower TP-AGB contribution for $\log(\mathrm{age/yr})\lesssim 8.8$ and $\log(\mathrm{age/yr})\gtrsim 10$ predicted by S_07 and S_35 with respect to S_00 is a consequence of the reduced lifetimes of the more massive stars (\mbox{M_{\rm i}}\gtrsim 3~{}\mbox{\mathrm{M}_{\odot}}) and the low-mass O-rich stars respectively.

In general, our calibration favours a mild TP-AGB contribution that does not exceed $\simeq$ 50 per cent in the $K$-band luminosity ($\simeq$ 20 per cent in bolometric luminosity) for ages between $\sim 0.5$ and $\sim 2$ Gyr. These results are closer to the BC03 models, whereas M05 models show a larger contribution of the TP-AGB, up to 80 per cent in the $K$-band luminosity ($\simeq$ 30 per cent in bolometric luminosity). We note that a sizable reduction of the TP-AGB contribution in M05 models is also suggested by Noël et al. (2013), based on a re-analysis of the integrated colours of Magellanic Clouds’ clusters.

5.6 Characterization of AGB star population

Given our detailed population synthesis simulations we can provide a full characterization of the AGB population in the SMC in terms of stellar parameters. A comparison between the best-fitting models S_07 and S_35 and the observations in the $K_{\rm s}$ vs. J\!-\!\mbox{K_{\rm s}} CMD is shown in Fig. 17. Both sets of models give essentially the same results in terms of the predicted stellar photometry. Most of the stars classified as X-AGB are C-rich (in agreement with the results from the SED fitting by SR16), with a small contamination of O-rich stars undergoing HBB. The stars populating the upper part of the CMD, with $K_{\rm s}$ magnitudes brighter than 10 mag, are exclusively O-rich and most of them are experiencing HBB. The inset of Fig. 17 clearly shows that C-rich and O-rich stars cannot be accurately separated by using classical $J\!-\!K_{\rm s}$ criteria, in agreement with Boyer et al. (2015).

A more quantitative analysis can be performed with the help of the $K_{\rm s}$-band LFs decomposed in bins of selected stellar parameters. In Fig. 18, we show some of the most relevant ones for the O-rich and the C- and X-AGB populations (i.e. ${\rm C/O}$, $\mathrm{C}-\mathrm{O}$, $\dot{M}$ and $M_{\rm i}$) as predicted by the best-fitting set S_35. The bulk of O-rich stars contains low-mass stars with \mbox{M_{\rm i}}\lesssim 2 $\mathrm{M}_{\odot}$ and $K_{\rm s}$ magnitudes fainter than $\approx 11$ mag. Most massive TP-AGB stars mainly populate the bright-end of the LF (see panel (a) of Fig. 18). As shown in panel (b) of Fig. 18, the predicted values of the C/O ratios are between 0.2 and 0.6 for the $\approx 80$ per cent of the stars. The contribution of stars that have a lower surface abundance of carbon as a consequence of the HBB is visible at \mbox{K_{\rm s}}\approx 9.5 mag. As for the mass-loss rates attained by O-rich stars, the models predict a rate around $10^{-7}$ $\mathrm{M}_{\odot}\,\mathrm{yr}^{-1}$ for a half of them, and the higher mass-loss rates, i.e. \mbox{\dot{M}}\gtrsim 10^{-5} $\mathrm{M}_{\odot}\,\mathrm{yr}^{-1}$, are predicted for the most massive and luminous stars (see panel (c) of Fig. 18). A fundamental parameter for the C-rich stars is the $\mathrm{C}-\mathrm{O}$ as it is the main parameter that determines the activation of the CDYN $\dot{M}_{\rm dust}$, the minimum $\mathrm{C}-\mathrm{O}$ value being $\approx 8.2$. The bulk of C-stars has \mbox{\mathrm{C}-\mathrm{O}}\gtrsim 8 and values larger than 9 are attained by the X-AGB stars, as can be appreciated in panel (d) of Fig. 18. The predicted values of the C/O ratio are also shown in panel (e) of Fig. 18 and they essentially reflect the distribution of $\mathrm{C}-\mathrm{O}$ values. Finally, the predicted mass-loss rates for the C-rich stars are similar to those of the O-rich, with the X-AGB stars reaching the higher mass-loss rates (see panel (f) of Fig. 18).

The distributions of the predicted mass-loss rates for the C-, X- and O-AGB synthetic population are shown Fig.19. The predicted mass-loss rates for both O-rich and C-rich stars are in agreement with the mass-loss rates estimates derived from the spectral energy distribution (SED) fitting performed by Nanni et al. (2018) for the same sample of C- and X-AGB stars used in this work. In particular, the separation between C- and X-AGB stars at \mbox{J!-!K_{\rm s}}\ \approx 2 occurs when the mass-loss rate reaches $\approx 10^{-6}$ $\mathrm{M}_{\odot}\,\mathrm{yr}^{-1}$ (see Fig. 17), in agreement with the value found by Nanni et al. (2018). The predicted ranges of mass-loss rates are also in agreement with the results of three additional SED fitting studies of AGB stars in the SMC by SR16, Groenewegen & Sloan (2018) and Goldman et al. (2018).

It is also interesting to compare the distributions of the bolometric magnitudes ($M_{\rm bol}$) as predicted by the best-fitting model S_35 with the results of SR16, B11 and Nanni et al. (2018). SR16 estimated the luminosities of the observed AGB stars in the SMC using the results of their best-fitting SED models, whereas B11 performed a trapezoidal integration from the observed optical $U$-flux through the mid-infrared 24 $\mu$m-flux. Nanni et al. (2018) performed a similar study on the sample of C-rich and X-AGB stars identified by SR16. The predicted $M_{\rm bol}$-LFs as derived from the average LFs of the 10 TRILEGAL simulations with the TP-AGB set S_35 are shown in Fig. 20. The predicted peaks of the $M_{\rm bol}$ distributions, i.e. the median values, of the O-rich, C-rich and X-AGB LFs are $-4.11$, $-4.49$ and $-4.86$ mag, respectively. The combined C- and X-AGB LF peaks at $-4.54$ mag. The values reported by B11 are $-4.59$ mag for the C-AGB and $-4.63$ mag for the C- and X-AGB samples combined, and the values estimated by SR16 are $-4.2$ mag for the O-rich and $-4.5$ mag for the C-rich stars (including the stars classified as X-AGB). The bolometric luminosity distribution of the C-rich and X-AGB stars estimated by Nanni et al. (2018) are very similar and are well reproduced by our best-fitting model. The predicted bolometric luminosities are in excellent agreement with both the results of B11 and SR16, given that small discrepancies can be due to the fact that both authors adopted a fixed distance modulus across the SMC.

5.7 2MASS and Spitzer luminosity functions

Figure 21 shows the comparison between our best-fitting model S_35 and the LFs in the 2MASS and Spitzer bands. We find a general satisfactory agreement between our best-fitting model and the observed LFs, which supports the validity of our treatment of the circumstellar dust based on an improved version of Marigo et al. (2008). A few discrepancies exist, in particular affecting the LFs of X-AGB class in some photometric bands. This is the case of the Spitzer filters [4.5], [5.8] and [8.0]; the predicted X-AGB LFs are shifted to brighter magnitudes with respect to the observed ones. At the same time, the predicted C-AGB LFs in the Spitzer filters [5.8] and [8.0] appear slightly shifted towards magnitudes fainter than observed. Overall, these differences are of small entity, and more importantly, we expect they do not impact on the results of our TP-AGB calibration, since the X-AGB sample includes less than 6 per cent of the total number of AGB stars.

6 Summary

In this study we present and discuss a detailed calibration of the TP-AGB phase, focusing on resolved stellar populations in the Small Magellanic Cloud. To this aim, we couple high-quality observational data (from the optical to the mid-infrared) to detailed models for AGB stars. Several ingredients are put in place, namely: the spatially resolved SFH and AMR of the SMC derived from the deep near-infrared photometry of the VMC survey (Rubele et al., 2018); the complete census of TP-AGB stars in the SMC identified in the 2MASS and Spitzer pass-bands, together with their classification (O-rich, C-rich, extreme-AGB stars) contained in the catalogues by Boyer et al. (2011) and Srinivasan et al. (2016); large grids of stellar evolutionary tracks computed with the PARSEC code (Bressan et al., 2012), corrected by the early-AGB mass loss, and coupled with detailed TP-AGB evolutionary tracks computed with the COLIBRI code (Marigo et al., 2013); state-of-the-art stellar spectral libraries for O-rich (Aringer et al., 2016) and C-rich (Aringer et al., 2009) stars; a treatment to account for the effect of circumstellar dust on the photospheric emission based on radiative transfer calculations (an improved version of Marigo et al., 2008); new pulsation models for long-period variables (Trabucchi et al., 2019); extended grids of stellar isochrones with rich information about the properties of AGB stars (modulation of luminosity and effective temperature during the thermal pulse cycles, surface chemical abundances) following the implementation described in Marigo et al. (2017); a well-tested stellar population synthesis code (Girardi et al., 2005, TRILEGAL) that includes a proper estimation of photometric errors and completeness.

We produced a large grid of TP-AGB models with the aim of identifying suitable combinations of mass loss and third dredge-up laws that best fit a set of observables at the same time, namely: the $K_{\rm s}$-band luminosity functions of the observed O-rich, C-rich and extreme-AGB stars. Among the many computed sets, two are found to recover the observational constraints with comparable performance. Both sets share the same prescriptions for the pre-dust mass loss (Cranmer & Saar, 2011), and the dust-driven mass loss of C-rich stars (Eriksson et al., 2014; Mattsson et al., 2010). The main differences between them refer to the adopted efficiency of the dust-driven winds during the O-rich stages, and the description of the 3DU.

The first best-fitting set describes the efficiency $\lambda$ of the 3DU following the formalism introduced by Karakas et al. (2002) on the basis of their full TP-AGB models. Given the high efficiency ($\lambda\simeq 1$) predicted by this formalism for \mbox{M_{\rm i}}\ga 3 $\mathrm{M}_{\odot}$, the only way to avoid a sizeable excess of C-stars at magnitudes \mbox{K_{\rm s}}<10 mag is to invoke a strong mass loss during the O-rich stages that precede the formation of C-stars (Bloecker, 1995, with \mbox{\eta_{\rm dust}}=0.06). In this way, for stars with 3~{}\mbox{\mathrm{M}{\odot}}\la\mbox{M{\rm i}}\la 4~{}\mbox{\mathrm{M}{\odot}} (at \mbox{Z{\rm i}}=0.004) the complete ejection of the envelope takes place before the 3DU is able to increase the C/O ratio above unity. For \mbox{M_{\rm i}}\ga 4~{}\mbox{\mathrm{M}_{\odot}} the lack of C-stars is also due to the concomitant occurrence of HBB which consumes carbon in favour of nitrogen.

The second best-fitting set adopts a new formalism for the 3DU. Starting from the results of full AGB models in the literature (Ventura & D’Antona, 2009; Cristallo et al., 2015), we build a parametrization in which $\lambda_{\mathrm{max}}$ (the maximum value attained during the TP-AGB phase) first increases with the mass of the star, reaches a maximum value and then decreases in the most massive AGB stars. Our calibration indicates that at \mbox{Z_{\rm i}}=0.004 the efficiency of the 3DU increases up to about \mbox{\lambda_{\mathrm{max}}}\simeq 0.6 in stars with \mbox{M_{\rm i}}\simeq 2.0-2.5 $\mathrm{M}_{\odot}$, and then decreases at larger masses, down to \mbox{\lambda_{\mathrm{max}}}\la 0.2 for \mbox{M_{\rm i}}>5 $\mathrm{M}_{\odot}$. We note that such decline of $\lambda_{\mathrm{max}}$, for \mbox{M_{\rm i}}>4 $\mathrm{M}_{\odot}$in particular, is at odds with the 3DU properties of the TP-AGB models by Karakas et al. (2002), while appears in line with the trends predicted by the TP-AGB models of (Ventura & D’Antona, 2009) and Cristallo et al. (2015). Thanks to the reduction of $\lambda_{\mathrm{max}}$ at larger masses, the second best-fitting set avoids naturally the excess of bright carbon stars, so that the O-rich stages can be suitably described with a lower efficiency parameter for the dust-driven mass loss (Bloecker, 1995, with \mbox{\eta_{\rm dust}}=0.02-0.03).

The two sets produce similar star counts for all classes of AGB stars in the SMC, recover their $K_{\rm s}$-band luminosity functions quite well, and yield comparable contributions to the integrated light of simple stellar populations as a function of the age. The comparison with the Spitzer observations in the mid-infrared bands is also satisfactory, indicating that the treatment of the bolometric corrections for the circumstellar dust around mass-losing AGB stars is reasonable.

The similarly good performance of both sets suggests that the final results are affected by a subtle degeneracy. However, this latter can be resolved by considering complementary predictions for AGB stars, namely the IFMR and the AGB chemical ejecta. As to the IFMR, we find that, compared to the S_07 set, the second best-fitting set S_35, characterized by a reduced efficiency of the 3DU for \mbox{M_{\rm i}}>3 $\mathrm{M}_{\odot}$, at solar metallicity predicts higher white dwarf masses, getting in better agreement with the semi-empirical relations newly derived by El-Badry et al. (2018) and Cummings et al. (2018). Based on this result, we may conclude that the set S_035 is to be preferred.

The second aspect that is relevant to further probe the models is related to the chemical yields of AGB stars. In particular, a careful examination of the surface elemental abundances (e.g., lithium, CNO elements, sodium, magnesium, aluminum, s-process elements) in high-luminosity AGB stars that experience HBB may help to discriminate between the two sets. We postpone a detailed investigation to a follow-up work.

7 Conclusions

The main conclusions we draw from the results of the present work can be outlined as follows:

• •

We confirm the importance of the pre-dust mass loss for the low-mass (with initial masses around 1 $\mathrm{M}_{\odot}$), low-metallicity O-rich stars, as already found by earlier works (Girardi et al., 2010; Rosenfield et al., 2014, 2016).

• •

The inclusion of a mass-loss prescription based on detailed dynamical models for C-stars (Eriksson et al., 2014; Mattsson et al., 2010) is of key importance to correctly describe the population of this class of TP-AGB stars, given its dependence on stellar parameters, in particular the carbon excess.

• •

The minimum and maximum initial masses for the formation of C-stars are around 1.4 $\mathrm{M}_{\odot}$ and 3 $\mathrm{M}_{\odot}$  respectively, at SMC-like metallicity (i.e. \mbox{Z_{\rm i}}=0.004) that is representative for the bulk of C-stars. The lack of C-stars for \mbox{K_{\rm s}}<10 mag clearly indicates that the third dredge-up in TP-AGB stars with \mbox{M_{\rm i}}\ga 3~{}\mbox{\mathrm{M}_{\odot}} should not take place over many thermal pulses, nor with extreme efficiency. The minimum mass threshold is found to increase with the initial metallicity, in agreement with previous studies (e.g. Marigo & Girardi, 2007). The maximum mass threshold also shows a positive correlation with metallicity, though less pronounced. In general, these trends would lead the C/M ratio to decrease in host systems of increasing metallicity (depending also on their SFH), a well-known fact in the literature (e.g. Boyer et al., 2013).

• •

The number counts of the brightest O-rich stars experiencing HBB can be reproduced by invoking a mass loss with a significant luminosity dependence (Bloecker, 1995, with \mbox{\eta_{\rm dust}}=0.02-0.03 for a shallow 3DU, or with \mbox{\eta_{\rm dust}}=0.06-0.07 for an extremely deep 3DU). We find that a good reproduction of the O-rich LF is also attained with the empirical relation between mass-loss rate and pulsation period proposed by Vassiliadis & Wood (1993).

• •

A direct outcome of our calibration, based on resolved stellar populations, is the full characterization of the AGB population in the SMC in terms of stellar parameters (i.e. initial mass, mass-loss rate, C/O, carbon excess, luminosity). The distribution of mass-loss rates for C-rich and extreme-AGB stars agree well with the one derived by Nanni et al. (2018) through the technique of SED fitting applied to the same observed catalogue adopted in this work. We also find a good agreement with the mass-loss rate estimates derived by Groenewegen & Sloan (2018) and Goldman et al. (2018). The main features (the position of the peak and the width) of the predicted distributions in bolometric magnitude for the O-rich, C-rich and extreme-AGB star classes are in line with studies based on the SED-fitting methodology (Srinivasan et al., 2016), as well as with the work of Boyer et al. (2011).

• •

Considering the characteristics of the method and the data employed in the calibration, it is fair to state that some conclusions might be more robust than others. The most reliable constraints apply to the lifetimes and the bolometric luminosities of the different TP-AGB subtypes (especially Figs. 13 and 14), at SMC metallicities. These properties are directly related to the observational LFs being fitted, with some uncertainty coming mainly from the limited resolution in initial masses inherent to our method. Lifetimes and bolometric luminosities are also linked to the contribution of TP-AGB stars to the integrated SED of galaxies. Our results, obtained at SMC metallicities, appear to support a TP-AGB contribution intermediate between those favoured by the popular models from Bruzual & Charlot (2003) and Maraston (2005) – although pending towards the “light contribution” side given by Bruzual & Charlot (2003) models. In this respect additional indications may be derived in the near future, through the analyses of other, more metal-rich galaxies (namely the LMC and M31) using the same method.

While the present calibration covers the metallicity range defined by the AMR of the SMC ($0.0005\lesssim Z\lesssim 0.008$), we extend the calculation of the TP-AGB tracks to cover the interval $0.0005<Z<0.03$ ($-1.5\lesssim\mbox{\rm[{\rm Fe}/{\rm H}]}\lesssim+0.3$). Extended PARSEC+COLIBRI isochrones derived from these tracks (of both sets S_07 and S_35) can be retrieved from our web interfaces http://starkey.astro.unipd.it/cmd and http://stev.oapd.inaf.it/cmd, for over 50 different photometric systems including major photometric surveys and instruments (e.g. 2MASS, SDSS, Spitzer, AKARI, HST/ACS, HST/WFC3, WISE, VISTA, Gaia, TESS), and future datasets (e.g. JWST, LSST, Euclid). The isochrones also include pulsation properties of long-period variables according to the prescriptions of Trabucchi et al. (2019). A general description of the isochrone characteristics and data format is provided in Marigo et al. (2017).

Acknowledgements

This work is supported by the ERC Consolidator Grant funding scheme (project STARKEY, G.A. n. 615604). We thank the entire VMC team for producing the space-resolved SFH used in this work; J. Cummings and J. Kalirai for kindly providing us with their IFMR data. MRC acknowledges funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No. 682115). Many thanks go to C. Maraston, S. Charlot, G. Bruzual for providing us with their stellar population synthesis models. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation. This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration et al., 2013, 2018) and matplotlib, a Python library for publication quality graphics (Hunter, 2007).

Appendix A Criteria for classifying AGB stars

A.1 Observational criteria

SMC stars in B11 and SR16 are classified into red Supergiants (RSG), Carbon-rich AGB (C-AGB), Oxygen-rich AGB (O-AGB), anomalous Oxygen-rich AGB (a-AGB) and extreme AGB (X-AGB). C-AGB and O-AGB stars are selected in the $J\!-\!K_{\rm s}$ CMD, following Cioni et al. (2006). All C-AGB and O-AGB stars are brighter than K0 line, defined as:

[TABLE]

where the value of the metallicity is $\mbox{\rm[{\rm M}/{\rm H}]}=\log(Z_{\mathrm{SMC}}/Z_{\odot})$. The adopted values of metallicity and distance of the SMC are $Z_{\mathrm{SMC}}=0.2$ $Z_{\odot}$ and $d_{\mathrm{SMC}}=61$ kpc. Two additional lines separate between C- and O-AGB stars:

[TABLE]

C-AGB stars have $J\!-\!K_{\rm s}$ colors redder than the K2 boundary while O-AGB stars lie between K1 and K2. To minimize the contamination from RGB stars, sources fainter than the TRGB level in both $K_{\rm s}$ and 3.6$\mu$m passbands – estimated to be at 12.7 mag and 12.6 mag, respectively – are excluded from the sample.

The class of the reddest AGB stars is denoted as extreme AGB. The majority of them are probably in the “superwind” phase with high mass-loss rates. Since they are obscured at optical wavelengths by thick dusty envelopes, the selection also includes mid-IR photometry. The sources classified as X-AGB stars are brighter than the 3.6$\mu$m TRGB and redder than \mbox{J!-![3.6]}>3.1 mag, or redder than \mbox{[3.6]!-![8]}>0.8 mag if the $J$-band detection is not available. Two additional criteria allow to minimize the contamination from YSOs and unresolved background galaxies (see 3.1.2 in Boyer et al., 2011).

The classification criteria adopted for selecting the RSGs are as follows: (1) they are bluer than the K1 line; (2) they have $K_{\rm s}$ brighter than the TRGB and, (3) they are redder than the KRSG line which is parallel to the K1 line but shifted by \Delta(\mbox{J!-!K_{\rm s}})=-0.2. These criteria minimize the contamination from O-rich stars, from RGB stars and foreground/background sources.

Boyer et al. (2011) identified a new feature in the $[8]$ vs. $J\!-\![8]$ CMD suggesting the presence of a class of stars distinct from C- and O-AGB. They referred to these stars as anomalous Oxygen stars (aO-AGB), as they were originally classified as O-rich by Cioni et al. (2006). These stars were selected from the original O-AGB sample if they are redder than the line ($aO$):

[TABLE]

with $A=27.95$ mag, and with an 8$\mu$m absolute magnitude fainter than $M_{8}=-8.3$ mag. Stars in the original C-AGB sample fall in this population if they are bluer than the line ($aC$) defined by Eq. 8 with $A=31.47$ mag. In a more recent work, Boyer et al. (2015) carried out a detailed analysis of chemistry, pulsation properties, stellar parameters and dust production of these stars and concluded that they are low-mass dusty AGB stars at the very end of their evolution, with median current stellar masses of about 0.94 $\mathrm{M}_{\odot}$ and initial masses $M<1.25$ $\mathrm{M}_{\odot}$. The spectral classification has been performed for a sample of 273 aO-AGB stars and resulted in 122 C-rich, 100 O-rich, 23 S-type and 28 unknown spectra. Since a high fraction ($\sim 45$ per cent) of the aO-AGB stars turned out to be C-rich, Boyer et al. (2015) refer to aO-stars simply as a-AGB stars. Their results suggest that a-AGB stars should be photometrically selected using the $J\!-\![8]$ colour, which is more reliable with respect to $J\!-\!K_{\rm s}$ colour for the most evolved stars.

Figure 22 shows the $K_{\rm s}$ vs. J\!-\!\mbox{K_{\rm s}} CMD of the observed sample of evolved stars in the SMC area used in this work. The grey lines are the photometric criteria adopted to separate RSG from AGB stars and O-rich from C-rich stars. The dashed line is the tip of the RGB in the $K_{\rm s}$ band. Stars are colour-coded according to the evolutionary stage as shown in the legend of the figure.

A.2 Model selection criteria

The X-AGB stars can only be selected using photometric criteria, as they do not correspond to a well defined set of intrinsic stellar parameters. Therefore, a star in the synthetic catalogue is identified as X-AGB following the same criteria as B11: 1) \mbox{J!-![3.6]}>3.1 mag or \mbox{[3.6]!-![8]}>0.8 mag if $J>19.69$ mag; 2) $[3.6]<[3.6]_{\mathrm{TRGB}}$. To identify the RSG sequence we use again the photometric criteria: 1) Ks<\mbox{K_{\rm s}}_{\mathrm{TRGB}}; 2) $[3.6]<[3.6]_{\mathrm{TRGB}}$ and 3) K_{\mathrm{RSG}}<=\mbox{K_{\rm s}}<=K1.

A simple comparison between models and observations suggests that the sequence of O-AGB stars is populated by both TP-AGB stars and early-AGB stars, which can be easily distinguished in the synthetic catalogues (looking at their evolutionary stages) but not in the data. To properly compare the synthetic and observed O-AGB sequences we first select them as in B11: 1) $J\!-\!K_{\rm s}$ redder than $K1$; 2) \mbox{K_{\rm s}}<\mbox{K_{\rm s}}_{\mathrm{TRGB}}; 3) $[3.6]<[3.6]_{\mathrm{TRGB}}$ . We then considered the value of the ${\rm C/O}$ ratio and selected only stars with \mbox{{\rm C/O}}<1.

The C-stars are selected according to the ${\rm C/O}$ ratio, i.e. if \mbox{{\rm C/O}}>1. We exclude from the C-AGB sequence the stars fainter than both \mbox{K_{\rm s}}_{\mathrm{TRGB}} and $[3.6]_{\mathrm{TRGB}}$ and the already selected sample of X-AGB stars.

The observed class of a-AGB stars could be identified in the synthetic catalogues by using the photometric criteria, but this class hosts both C- and O-rich AGB stars that cannot be clearly separated using photometry. The only information we have about their chemical type is the spectral classification performed by Boyer et al. (2015). They concluded that in the SMC the ratio between C-rich and O-rich a-AGB stars is $\approx 50$ per cent. When comparing the simulated and observed $K_{\rm s}$-band LFs, we take into account the contribution of the a-AGB stars by weighting the observed C- and O-AGB LFs with the a-AGB LF, according to the result of Boyer et al. (2015). Similarly, the total number counts of the C- and O-AGB stars are corrected to include the contribution of the aAGB stars. Therefore, the total number of observed C-AGB and O-AGB changes from 1854 to 2340 and 2623 to 3109 respectively.

Appendix B The colour shift of giants and supergiants

By comparing the observed and simulated $K_{\rm s}$ vs. J\!-\!\mbox{K_{\rm s}} CMD we found an evident discrepancy in the colour and slope of the RSG and O-rich TP-AGB sequences. The left panels of Fig. 23 show the CMD of the SMC from the 2MASS data and the TRILEGAL simulation performed with the TP-AGB set S_00. By visual inspection it is immediately evident that the simulated RSG and O-rich AGB sequences are redder than the observed ones.

To quantify such a shift in colour, we calculate the fiducial lines of both the observed and simulated RSG, O-rich and RGB sequences. The fiducial lines are obtained by fitting a Gaussian distribution in the case of the RGB and a double Gaussian for the RSG and AGB to the distributions of the $J\!-\!K_{\rm s}$ colours in different magnitude bins. We use a linear fit to infer the value of the shift as a function of the $K_{\rm s}$ magnitude, for the RGB stars (\mbox{K_{\rm s}}>12.7), for the AGB stars (11.8<\mbox{K_{\rm s}}<12.7) and for the RSG stars (\mbox{K_{\rm s}}<12.7). The linear fits are shown as solid lines. Figure 24 shows the $J\!-\!K_{\rm s}$ difference between the observed and simulated sequences as a function of $K_{\rm s}$. The colour of the RGB sequence is reproduced within 0.025 mag, while the RSG and the AGB are shifted toward redder colours, with differences spanning from 0.05 to 0.1 mag.

For the purpose of this work, it is important to have the right slope and colours for the different evolutionary sequences because of the photometric selections that have been performed on the data. This is particularly important for the selection of RSG and O-AGB stars for which we rely on the photometric criteria of SR16 and B11.

To reconcile the differences in colours, we correct the synthetic photometry a posteriori by assuming that the $J\!-\!K_{\rm s}$ colour difference is due to an extinction-term correction $dA_{\lambda}$. We first derive the $dA_{J}$ correction from the $J\!-\!K_{\rm s}$ colour difference as a function of the $K_{\rm s}$ magnitude for the RGB, RSG and AGB sequences, by assuming a zero correction for the $K_{\rm s}$ band. We also calculate the value of $dA_{\lambda}$ as a function of $K_{\rm s}$ magnitude for all the photometric bands bluer than $K_{\rm s}$, considering the corresponding extinction coefficients. The corrections for the $J$-band, $dA_{J}$, are the following:

[TABLE]

We find that similar shifts occur also at shorter wavelengths, and that they can be likewise mitigated by applying such corrections.

In Figure 23 we show the simulated CMD after correcting the sequences and the resulting $J\!-\!K_{\rm s}$ colour differences. The differences for the RGB and RSG sequences are reduced to below 0.01 mag, while a slightly larger difference (up to 0.05 mag) is still present in the AGB sequence.

Such residual difference is due to the fact that in the synthetic catalogues, we apply the correction derived for the RSG sequence to the stars brighter than \mbox{K_{\rm s}}=12.7 mag and the correction derived for the RGB sequence to the stars fainter than \mbox{K_{\rm s}}=12.7 mag. This choice is motivated by the fact that RSG and RGB models are computed with the PARSEC code and they are not the subject of the present calibration.

Given their wavelength-dependent behaviour, it is unlikely that these shifts are due to some abundance or opacity effect. They could well be caused by small temperature offsets. Independently of their possible cause, we emphasize that these shifts do not significantly affect the results of our calibration.