SEMI-EMPIRICAL WHITE DWARF INITIAL–FINAL MASS RELATIONSHIPS: A THOROUGH ANALYSIS OF SYSTEMATIC UNCERTAINTIES DUE TO STELLAR EVOLUTION MODELS

Our reference data set comprises 52 WDs in 10 open clusters, plus Sirius B, for a total of 53 DA objects, as reported in Table 1. We considered the same cluster sample as in the recent work by Ferrario et al. (2005) that contains data for WDs in the Pleiades (Ferrario et al. 2005), Hyades, and Praesepe (Claver et al. 2001; Dobbie et al. 2004), NGC 2516 (Koester & Reimers 1996), NGC 3532 (Koester & Reimers 1993), NGC 2099 (Kalirai et al. 2005), NGC 2168 (Williams et al. 2004; Ferrario et al. 2005). We wish to note that our Pleiades sample contains only the object listed by Ferrario et al. (2005). We did not include data for additional two WDs (GD 50 and PG 0136+251) discussed by Dobbie et al. (2006b), because their membership of the Pleiades is not absolutely certain.

To this cluster sample, we added NGC 6819, NGC 7789, and NGC 1039. For the WDs in the first two clusters, T_eff and g values are taken from Kalirai et al. (2008), and from Rubin et al. (2008) for the last one.

The T_eff, g, and corresponding error estimates for the clusters in common with Ferrario et al. (2005) are the same as in their Table 1, with a few exceptions. Values of T_eff and g for two Praesepe WDs (0837+185 and 0837+218), identified as WD Nos. 13 and 14 in Table 1, have been updated following Dobbie et al. (2006a). The same authors provide data for four additional WDs in the same cluster (identification (ID) numbers from 16 to 19 in Table 1) that we have included in our analysis. And finally, we also use their data for the single Pleiades WD.

As an additional note about the cluster WD sample, our referee (K. Williams) has pointed out that the Hyades WD 0437+138 that appears as WD No. 8 in Table 1, has been incorrectly included as a DA object in (Ferrario et al. 2005) study. In fact, it is a DBA WD, with log(n_H/n_He) = −4.5 in the atmosphere. To assess the errors in the derived mass and cooling age for this object when using DA models, we have made the following test. We have compared log(g) values derived from DB and DA models at the observed T_eff of WD No. 8, for the two extreme cases of 0.55 and 1.0 M_☉ (the mass we obtain from the analysis with DA models is ∼0.77 M_☉). Our adopted DB models have been computed as described in Salaris et al. (2001).

For both masses, the difference in log(g) derived from the DA and DB tracks is smaller than the 1σ error bar on the observed value for WD No. 8 (as listed in Table 1). Also, the difference in cooling ages at the observed T_eff is well below the error bars we derive from our analysis, as reported in Table 5. This implies that the analysis of WD No. 8 with DB models would provide mass and cooling age compatible (within the 1σ errors) with the results obtained from the DA models employed in this analysis. The reason is that the mass–radius relation is not strongly affected by the atmospheric composition, and also at the relatively high T_eff of this object, the evolutionary timescales of DA and DB models are similar. Moreover, given that this WD does not have a pure He envelope, differences will be even more negligible.

In addition to a large sample of cluster WDs, we also include Sirius B in our analysis. The data for Sirius B (also included in Ferrario et al. 2005) in Table 1 are from the recent determination by Barstow et al. (2005). From the estimated mass and radius of the primary MS star (Sirius A) one can derive its age by comparing with stellar evolutionary tracks in the mass–radius plane. We assumed [Fe/H] = 0.0 for Sirius A, following Liebert et al. (2005), with which we associated an 0.10 dex 1σ error. The age of the star is obtained from models bracketing the metallicity of the system (if different from the values provided by the employed isochrone grid). Linear interpolation in [Fe/H] between the two bracketing results gives then the age for the appropriate [Fe/H]. The age of Sirus A provides the age of the system t_sys. Sirius B cooling time t_cool and mass are determined from the spectroscopic T_eff and g estimates, as described before. The difference t_sys-t_cool provides the age of the Sirius B progenitor, hence M_i after using a grid of theoretical lifetimes for solar metallicity.

The sources for the cluster [Fe/H], E(B − V), and their CMDs used for age determinations, are listed in Table 2.

Values for [Fe/H] have been taken mainly from the compilation by Gratton (2000), who presents a collection of [Fe/H] estimates for a large sample of open clusters, recalibrated against results from high-resolution spectroscopy. Alternative estimates have been adopted for NGC 2099, which is not in the Gratton (2000) compilation, NGC 6819 and NGC 1039, for which recent high-resolution spectroscopic measurements have appeared in the literature (see Table 2).

Reddenings are taken from the compilation by Loktin et al. (2001). These are weighted averages of estimates obtained with various methods, all based essentially on color–color relationships. For NGC 6819, we employed a very recent E(B − V) determination based on spectroscopy (Bragaglia et al. 2001), and for the Pleiades reddenings to individual stars have been adopted, as discussed in Percival et al. (2003). The case of NGC 7789 is discussed below.

In case of the Hyades and Praesepe, we have assumed—as usual—zero reddening. Throughout this paper, we use A_V = 3.2E(B − V). We adopted a 0.02 mag uncertainty in E(B − V) when no error is quoted by the source of the data.

To estimate the age of the whole cluster sample, we have followed a two-step procedure where we tried to minimize the inputs from theoretical isochrones. First, cluster distance moduli have been determined using a fully empirical MS fitting to the V − (B − V) CMDs that makes use of the sample of field dwarf with accurate Hipparcos parallaxes and [Fe/H] presented by Percival et al. (2003). Details about the method can be found in the same paper.

Here, we just recall that for each cluster an MS fiducial line has been derived by plotting color histograms in V-magnitude bins, and considering the mode of the resulting distribution. In this way, one minimizes the impact of unresolved binaries. The resultant points are fitted to a cubic function that provides the cluster fiducial line. After applying (empirical) color corrections, which account for the actual cluster metallicity and for the reddening E(B − V), the field dwarf template sequence is shifted in magnitude to match the MS fiducial line, and the apparent distance modulus is obtained by χ² minimization. For objects like NGC 2168 or the Pleiades, with an extremely well-defined and thin (in color) MS, the separation of the binary sequence is very evident, and a cubic fit directly to the single-star MS is performed in order to derive the fiducial line.

In case of NGC 3532, we could not find available data with good photometry in the MS magnitude range of our template field dwarf sample. The distance modulus of this cluster has been obtained by comparison with Praesepe that has the same [Fe/H] within errors (and also a comparable age). The larger (as compared to our other estimates) error bar on its distance (hence on its age) reflects this less direct determination of the distance modulus.

As a final remark, we wish to note that, as discussed, e.g., in Percival et al. (2003), a simultaneous MS fitting to both VI and BV data can provide both the cluster distance modulus and reddening. This procedure has been applied (using the same field dwarf sample that has homogeneous BVI photometry) to NGC 7789, employing additionally the VI data by Gim et al. (1998). Our adopted E(B − V) for NGC 7789, used also in Percival & Salaris (2003), has been determined in this way. This method to determine the reddening could not, however, be applied to the whole cluster sample because of the lack of suitable BV and VI photometry for all clusters.

For a few clusters, we have used distances already published in the literature, obtained with the same technique and the same Hipparcos field dwarf sample. More in detail, for the Hyades, we have used the result by Percival et al. (2003) that is in perfect agreement with the Hipparcos distance by Perryman et al. (1998). For the Pleiades, we have used the distance modulus determined by Percival et al. (2005), while for NGC 7789, we have employed the result by Percival & Salaris (2003). The cluster NGC 6791, which is also part of the investigation by Kalirai et al. (2008) was not considered in our analysis because its metallicity [Fe/H] ∼ +0.45 (Gratton et al. 2006) is higher than the upper [Fe/H] range of the field dwarf sample ([Fe/H] = +0.30), and therefore we could not apply this empirical method to determine its distance. Since we did not want to introduce a second, inconsistent method, we disregarded this cluster.

Second, isochrone fits to the TO region (in the V − (B − V) plane), employing the distance moduli derived in the previous step, determine the cluster age.

The final ages for various choices of the isochrone grid, as well as the adopted cluster [Fe/H] and E(B − V) values, are reported in Table 3. Age uncertainties are due to both distance (which is in turn affected by the assumed [Fe/H] and E(B − V) values and their errors) and [Fe/H] uncertainties. Figure 1 compares the ages obtained with the BaSTI overshooting isochrones (that throughout the paper will be considered as our reference set of models to determine the IFMR; see Section 2.3) with the ages adopted by Ferrario et al. (2005) for the clusters in common (plus Sirius A). In spite of the dishomogeneity in terms of methods, data and isochrones employed, the two sets of ages compare well, within the respective error bars. The only exception occurs for NGC 2099 (M37), for which our derived distance modulus is 0.4 magnitudes larger than implicitly assumed in Ferrario et al. (2005)—who adopt the age estimate by Kalirai et al. (2001b)—leading to a smaller cluster age (320 Myr as compared to 520 Myr). However, if we adopt a lower metallicity and reddening as suggested by Kalirai et al. (2005) and discussed in Section 2.4 (see Table 3), we derive an age of 550 Myr for NGC 2099, in close agreement with Ferrario et al. (2005). Obviously, at this time, the uncertainties in the cluster parameters dominate the determination over those from using different theoretical isochrones, which constitutes our first result.

Zoom In Zoom Out Reset image size

Our reference IFMR is defined by the use of the following set of models. We refer the reader to the original papers quoted below for details of the models and calculations.

• 1.  
  Pietrinferni et al. (2004) models from the MS to the AGB with the inclusion of core overshooting during the central H-burning phase (hereafter BaSTI models with overshooting), to estimate in a consistent way both t_clus and t_prog.
• 2.  
  Salaris et al. (2000) WD cooling models (hereafter S00 WD models) to estimate t_cool. The code employed to compute these WD models is the same as for the progenitors. Some major sources of input physics in the WD models (see S00) are different from what has been employed in the BaSTI models, most notably equation of state and low-temperature opacities. This is somewhat unavoidable, given the different range covered by many physical parameters throughout the WD structures, as compared to their progenitors. It is very important to remark here about the internal CO stratification of the WD models. The CO profile in the inner core of a WD is built up during the central He-burning phase, whereas the more external part of the WD core is processed during the early AGB phase and the thermal pulses. Current uncertainties in the treatment of mixing in stellar interiors and in some key reaction rate (i.e., the C¹² + α reaction) introduce additional uncertainties in the predicted CO profiles. S00 models employ a CO stratification obtained from evolutionary models of solar metallicity progenitors, from the MS until the completion of the first thermal pulse (Salaris et al. 1997). The core mass and the associated CO profile at the end of the pulse were taken as being representative for the final WD object. When estimating the total error budget in the IFMR, we have accounted for the uncertainties in the CO profiles by using WD models with different C/O mass fraction ratios in their cores. This will be described in detail in Section 3. We anticipate, however, that these uncertainties are very small, since all the WD stars in our sample have not entered yet the crystallization phase, where the detailed core composition could have a larger effect.

We have studied systematic effects on the IFMR relation due to the following.

• 1.  
  Different model grids. We redetermine the IFMR, employing Girardi et al. (2000) models from the MS to the AGB (hereafter Padua models; they also include core overshooting but with a different formalism compared to BaSTI), keeping the WD models unchanged (S00).
• 2.  
  Effect of overshooting. We redetermine the IFMR, employing Pietrinferni et al. (2004) models from the MS to the AGB without core overshooting (hereafter BaSTI non-OV models), again keeping the WD models unchanged (S00).
• 3.  
  Different WD model grids. We redetermine the IFMR, employing the WD models computed with the LPCODE (hereafter LPCODE models; details of the LPCODE are given in Althaus et al. 2003), keeping the progenitor models unchanged (BaSTI with overshooting). This tests not only different codes, but also—in an unsystematic way—variations of the physics employed and of the WD chemical stratification.
• 4.  
  Systematic variation in WD input physics and chemical profiles. We redetermine the IFMR, employing the LPCODE WD models, by varying the CO-core stratification, H-envelope thickness, electron conduction opacities, and neutrino energy loss rates. In this case, the progenitor models are unchanged throughout (BaSTI with overshooting).

We present in this section our reference IFMR, while we defer to subsequent sections the discussion of the uncertainties involved in the IFMR determination. We use the S00 WD models to derive the WD (final) masses and cooling ages. Using our reference cluster ages (fifth column in Table 3) and the BaSTI models with overshooting, we construct the IFMR by deriving the progenitor (initial) masses. Figure 2 displays this reference semiempirical IFMR. Error bars include all uncertainty sources we consider. We defer a detailed discussion to Sections 3 and 4, but note here that the most relevant sources of uncertainty are due to the following effects: (1) errors in the WD T_eff and g values; (2) uncertainties in the cluster [Fe/H]; these affect the age through their influence on the derived cluster distance (mainly) and on the isochrones to be used for the age determination; (3) errors in the reddening E(B − V), which influence the age determination through the effect on the cluster distance. In case of one star (WD No. 38 in NGC 2099), the progenitor age results to be negative. We did not consider this object in the analysis of the IFMR.

Zoom In Zoom Out Reset image size

We have performed at this stage a simple check of our reference IFMR, comparing the distribution of M_i and M_f values with Ferrario et al. (2005)—disregrading the attached error bars—by means of a bi-dimensional Kolmogorov–Smirnov (KS) test (Press et al. 1992). We considered only clusters (plus Sirius) in common between the two studies. As customary, we accept the existence of a significant difference between the two IFMRs if the probability P that the two samples of (M_i, M_f) pairs were drawn from different distributions was higher than 95%. We obtain P values smaller than 50%, highlighting the statistical agreement between our reference IFMR and Ferrario et al. (2005) results for the clusters (plus Sirius B) in common.

We have also performed analytical fits to the data. In doing this, we have left out the data from NGC 2099, because our results for this cluster show some anomalies, e.g., a progenitor with negative age (WD No. 38 mentioned above), a 0.45 M_☉ WD with a progenitor mass above 8 M_☉ (WD No. 37), and some other issues to be discussed later in the paper.

In Figure 2, we show two analytic fits to the data. The first one, linear, gives

$$\begin{equation} M_{\rm f}= 0.084\, M_{\rm i} + 0.466 \end{equation} \tag{ 1 }$$

and a reduced chi-squared χ²_r = 3.7. Approximate values for the 1σ rms dispersion about this relation are: 0.075 M_☉ below 2 M_☉, 0.12 M_☉ between 2.7 and 4 M_☉, 0.09 M_☉ between 4 and 6 M_☉, and then it increases almost linearly up to 0.20 M_☉ at around 8.5 M_☉.

Guided by theoretical models, which show a break in the slope of the IFMR at around 4 M_☉ (see Figure 3), we have also performed a piecewise linear fit with a pivot point at M_i = 4 M_☉ that yields

$$\begin{equation} M_{\rm f} = \left\lbrace \begin{array}{@{}llll} 0.134 M_{\rm i} + 0.331 & & & 1.7 \;{M}_\odot \le M_{\rm i} \le 4\;{M}_\odot \\ 0.047 M_{\rm i} + 0.679 & & & 4\;{M}_\odot \le M_{\rm i}. \end{array} \right. \end{equation} \tag{ 2 }$$

This fit gives χ²_r = 2.7 at the expense of introducing only one more free parameter (we force the fit to be continuous at the pivot point). In this case, the dispersion is reduced to 0.05 M_☉ below 2 M_☉ and to 0.09 M_☉ above 6 M_☉, while it remains basically unchanged anywhere else. The reader should keep in mind that these cluster data do not constrain the IFMR at masses lower than about 1.7 M_☉. As a consequence, the present analytic relations, particularly the piecewise linear fit, should not be extrapolated toward lower initial masses.

Zoom In Zoom Out Reset image size

We could attempt more complicated fits to the data, but in our opinion the uncertainties in mass determinations, and possibly the existence of an intrinsic spread in the IFMR to be discussed below, render such efforts probably unnecessary. In fact, we have found that polynomial fits beyond second order do not produce any relevant improvement in the goodness of fit.

As an additional check, we performed a direct comparison of our fits with the polynomial relationship from Ferrario et al. (2005) that was derived by combining an assumed star formation rate and initial mass function for the Galactic Disk, with the observed field WD mass distribution. For masses between 2.5 and 6.5 M_☉—the range of initial masses covered in their study—at a given M_i, differences in the M_f values are less than 0.04 M_☉, well below the dispersions given above.

Another point, which we wish to address, is the presence of a possible intrinsic spread in the reference IFMR, as noted by Reid (1996) in an analysis of Praesepe WD population. With intrinsic spread, we mean a signature of differential mass loss among stars with the same mass and same initial chemical composition. Catalán et al. (2008a) have investigated very recently the spread around their mean semiempirical IFMR (determined from cluster WDs and WDs in common proper motion pairs) and found no clear cut indication that is either due to a spread in mass loss or to some other individual stellar properties (e.g., the thickness of the envelope layers) that vary between objects.

To investigate this issue with our data, we made the following test. We considered objects in a single cluster, namely Praesepe, with M_i values in agreement within their respective 1σ error bars. There are eight stars with M_i between 3.045 and 3.543 M_☉, satisfying this constraint, and we assume that—within the present error bars—they share the same M_i, equal to the average value 〈M_i〉 = 3.207 M_☉. For these objects, the average WD mass is 〈M_f〉 = 0.769 M_☉ with 1σ rms equal to 0.045 M_☉. To test whether this observed dispersion in M_f is due exclusively to the errors in the M_f estimates, we built a synthetic sample by drawing randomly 10,000 M_f values for each of the eight objects, according to a Gaussian distribution centered on 〈M_f〉 = 0.769 M_☉, with a 1σ spread given by their individual error bars. The synthetic distribution has obviously 〈M_f〉 = 0.769 M_☉, but a 1σ rms equal only to 0.028 M_☉. We compared statistically the σ values of the two distributions by means of an F-test that returned a probability of P = 96.3% that σ of the synthetic sample is different from (smaller than) the observed one. By fixing, as customary, a threshold of P = 95% to accept the existence of a real difference between the two σ values, our simple test seems to point to the existence of a spread in the IFMR, at least in this cluster. Adding five more objects in the Hyades and in NGC 3532 (the [Fe/H] estimates for these two clusters plus Praesepe span a range of only ∼0.10 dex, and overlap within their associated ∼2σ error bars, see Table 3), with initial masses within the same overlapping M_i range, and repeating the previous analysis, does not change this conclusion. In fact, for this enlarged sample, 〈M_f〉 = 0.729 M_☉, with a 1σ rms equal to 0.10 M_☉. The probability that the corresponding synthetic sample—constructed in the same way as for the eight Hyades WDs—has a different σ than observed is P = 96.9%.

However, this evidence for an intrinsic spread (i.e., due to a significantly varying mass loss at a fixed mass and chemical composition) in the IFMR is weakened if we take into account theoretical expectations for near-solar metallicities, displayed in Figure 3. Although the M_i values estimated for these 13 WDs (or just the eight Praesepe objects discussed above) overlap within the 1σ error bars, the same error bars suggest that they may also cover a range of values. Just considering the range covered by the central M_i estimates for these 13 WDs, theoretical IFMRs (see also another completely independent theoretical IFMR for solar metallicity in Figure 3 of Prada Moroni & Straniero 2007) generally display a steep rise of M_f with M_i. Over the M_i range between 3.05 and 3.5 M_☉, the theoretical IFMRs of Figure 3 predict an increase in M_f by about 0.1 M_☉, which corresponds to ∼2σ of the observed M_f spread in the eight Praesepe objects, and to ∼1σ of the observed spread in the 13-object enlarged sample. Without accounting for errors in the estimate of M_f, the predicted change within the M_i range allowed by the error bars of the estimated initial masses can potentially explain the observed M_f spread. This would imply that the observed spread is not the result of star-to-star mass-loss variations, but due to, as predicted from theory, the steep increase of M_f with M_i in the relevant mass range.

Moving this analysis to different M_i ranges can hardly provide additional constraints on the existence of an intrinsic spread in the IFMR, because of either smaller numbers of objects in comparable M_i intervals (when moving to lower M_i values) or a combination of much larger error bars on M_i and small numbers (when moving to higher M_i values). It may, however, be interesting to note the situation at M_i ∼ 1.8 M_☉, where we have five objects (belonging to NGC 7789 and NGC 6819) with formally the same M_i, small errors in both M_i and M_f and theoretical predictions for an IFMR essentially flat (see Figure 3). By applying the same statistical analysis described before, we do not find any statistically significant spread in the IFMR for this sample.

Before concluding this section, we discuss briefly another test of our reference IFMR. In the case of NGC 2099, we have considered the alternative value [Fe/H] = −0.20 and reduced reddening (E(B − V) = 0.23) used by Kalirai et al. (2005) and given by Deliyannis et al. (2002). Note, however, that Deliyannis (private communication) considers this as a possible but as yet unconfirmed value for [Fe/H]. Using this alternative metallicity and reddening, we derive an age for NGC 2099 of 550 Myr which, compared to 320 Myr derived with [Fe/H] = 0.09, has a dramatic effect in the inferred initial masses (Figure 7). As expected from the relatively large change in the inferred cluster age, changes in M_i due to the different cluster parameters generally exceed the uncertainties in M_i determination.

This test, based on estimates that are admittedly not confirmed, is just to exemplify in a quantitative way how crucial the assessment of the cluster [Fe/H] and E(B − V) values is.

We will return to this issue in Section 4.

We have tried to cover all major sources of uncertainties in the WD mass and cooling age determinations that are going to be discussed in this section. While we have used S00 models as our standard set of WD models, all the WD evolutionary sequences to test the effect of WD input physics and chemical stratification have been computed with the LPCODE described below. This should not invalidate our procedure, since we are using LPCODE models only to determine differential changes in the WD properties (e.g., core composition), not absolute values. In what follows, X represents the logarithmic values of either the WD mass or the cooling age.

• 1.  
  Observational uncertainties. The central values and 1σ uncertainties for log g and T_eff are listed in Table 1. We leave them unmodified as they are given in the original papers quoted in Section 2.1.
• 2.  
  Systematic uncertainties from WD cooling tracks. In addition to our standard choice of WD cooling tracks (S00), we have specifically computed a completely independent set of WD cooling models with the LPCODE. The aim is to estimate the systematic uncertainties in the inferred WD properties that arises from using different WD models. We define the 1σ uncertainty by

  $$\begin{equation} \sigma _{\rm syst}(X)=|X_{\rm S00}-X_{\rm LP}|, \end{equation} \tag{ 3 }$$

  where X_S00 is derived from S00 models and X_LP from models computed with the LPCODE. WD models computed with the LPCODE have been used in Serenelli & Fukugita (2007), where the evolution of a set of stellar models from 1 to 8 M_☉ and solar metallicity has been consistently followed from the main-sequence phase up to the WD phase, including the thermal pulse (TP) AGB phase. It should be noted that, in addition to differences in numerics, σ_syst(X) also accounts for differences in input physics between the two sets of models. While S00 and LPCODE models have used the same neutrino energy loss rates (Itoh et al. 1996), conductive opacities in the CO cores (Itoh et al. 1993), and very similar hydrogen envelope thickness and core composition, they differ in other relevant aspects such as the equation of state and low-temperature opacities (see Althaus et al. 2003 and Salaris et al. 2000 for details about the LPCODE and S00 input physics respectively) and the inclusion of gravitational settling in the model envelopes (S00 models mimic the effects of gravitational settling by adopting a pure hydrogen envelope and fixed chemical H/He and He/C/O interfaces, while in LPCODE the chemical profiles are evolved according to time-dependent equations of element diffusion).
• 3.  
  Neutrino energy loss rates. For the WD evolutionary sequences, we are considering in this work, plasma neutrinos always contribute at least 90% of the total neutrino emission. Uncertainties in the theoretical calculation of plasma-neutrino-emission rates, within the framework of the standard model of particle physics, are at most a few percent (M. Fukugita 2007, private communication), and thus we do not expect them to represent an important source of uncertainty in determining WD cooling ages. Several additional cooling mechanisms that could operate in WD interiors have been proposed, e.g., enhanced neutrino cooling due to a non-null neutrino magnetic dipole moment, axion production or other particles like weakly interacting, massive particles (WIMP; see Raffelt 1996 for a complete discussion of these mechanisms). It is not our aim here to consider the effects on the WD properties—particularly on the cooling rate—of all such mechanisms, especially because they remain purely hypothetical⁷. We do not want, however, ignore the possibility that some additional cooling might be present. Consequently, we have adopted a more conservative point of view and have computed additional WD cooling models where standard neutrino energy loss rates (Haft et al. 1994; Itoh et al. 1996) have been multiplied by factors of 2 and 0.5 respectively to take into account the possibility of such additional cooling mechanisms. These large factors were chosen to represent our conservative approach. We define the WD mass (and cooling age) 1σ uncertainties due to neutrino energy losses as

  $$\begin{equation} \sigma _\nu (X)=\frac{|X_{\rm 2}-X_{\rm 0.5}|}{2}, \end{equation} \tag{ 4 }$$

  where the subindices refer to quantities corresponding to the WD cooling models with the modified neutrino-emission rates as described above.
• 4.  
  Conductive opacity. For the physical conditions and composition of a representative set of WD models, we have computed conductive opacities according to both Itoh and collaborators (Itoh et al. (1993) and references therein—our references set of conductive opacities) and Potekhin and collaborators (see Cassisi et al. 2007 for the most up to date calculations). We find differences between both sets of calculations to be usually below 20% for the range of temperatures and densities relevant to the WD models used in this paper. In analogy with the neutrino-emission rates, we have estimated the effect of electron conduction uncertainties by computing additional WD models where the conductive opacities from Itoh and collaborators have been multiplied by factors of 1.25 and 0.8 and have used these models to define the WD mass and cooling age 1σ uncertainties due to errors in conductive opacities calculations as

  $$\begin{equation} \sigma _\kappa (X)=\frac{|X_{\rm 1.25}-X_{\rm 0.8}|}{2}. \end{equation} \tag{ 5 }$$
• 5.  
  WD core composition. The mass fraction of the central oxygen X_O in S00 WD models ranges from 0.65 up to 0.80, depending on the mass of the WD. The internal CO stratification is, however, subject to uncertainties due to the uncertainty on the C¹² + α reaction rate (see, e.g., S00) and, very important, the treatment of mixing (semiconvection, breathing pulses and the mechanism to suppress them) during the central He-burning phase (Straniero et al. 2003). These effects are larger than the expected change in CO profiles due to the small metallicity range covered by the cluster sample. In our analysis, we have been conservative and have defined the 1σ uncertainty due uncertainties in the core composition as

  $$\begin{equation} \sigma _{\rm core}(X)=\frac{|X_{\rm 0.9}-X_{\rm 0.3}|}{2}, \end{equation} \tag{ 6 }$$

  where X_0.9 and X_0.3 are derived from the WD cooling tracks, where the original chemical profiles have been rescaled to have X_O = 0.9 and 0.3 at the center, respectively.
• 6.  
  Hydrogen envelope thickness (M_H). As shown by, e.g., Provencal et al. (1998) the H-envelope thickness of DA WDs displays a range of values. To include this effect in our error analysis, in addition to the standard thick envelopes considered in our LPCODE calculations (usually 10⁻⁴ M_☉, but larger values were used for WD masses below 0.6 M_☉) we have computed also cooling sequences with thinner H-layers of 10⁻⁶M_☉. Similar to the above procedures, we have defined the corresponding 1σ errors

  $$\begin{equation} \sigma _{\rm env}(X)=\frac{|X_{\rm thick}-X_{\rm thin}|}{2}, \end{equation} \tag{ 7 }$$

  where X_thick and X_thin are derived from the WD models with thick and thin hydrogen envelopes respectively. We have performed an additional test to explore very thin hydrogen envelopes (as in Catalán et al. 2008a) by computing a 0.6 and a 1 M_☉ WD models with hydrogen envelopes of M_H ∼ 10⁻¹⁰ M_☉. When compared to models with thick envelopes, we have found that differences in the cooling ages never exceeded 11% in both cases, and are usually at the level of 7% for the 0.60 M_☉ model and 5% for the more massive one in the range of effective temperatures of interest to this work. Asteroseismology studies of ZZ Ceti stars disclose a wide range of hydrogen envelope thickness, ranging between ∼10⁻⁴M_☉ and ∼10⁻⁹ M_☉ (see, e.g., Castanheira & Kepler 2008). In this regard, the changes in the cooling ages due to uncertainties in the H-envelope thickness quoted above can probably be taken as robust upper limits.

For neutrino-emission rates, conductive opacities, core composition and hydrogen-envelope thickness we have assumed that uncertainties distribute normally. In the case of systematic uncertainties from different WD cooling tracks we have assumed uniform distributions of uncertainties (f_X) which are defined, in each case, by $\int _{-\sigma _{\rm syst}}^{+\sigma _{\rm syst}}{f_X dX}=0.683$. We tabulate σ_ν(X), σ_κ(X), σ_core(X), σ_env(X), and σ_syst(X) in (log g,T_eff) grids for using them in the MC simulations. We summarize the results for the WD uncertainties, both for mass and cooling age, in Table 4, where for a small subset of (log g, T_eff) values, we give the 1σ fractional uncertainties in WD mass and cooling age as defined by Equations (3)–(7), but transformed to linear scale to facilitate interpretation of the results.

Progenitor masses M_i are obtained from the estimated progenitor lifetime and stellar models. In addition to the uncertainty sources discussed in the previous section for t_cool, the determination of M_i is affected by

• 1.  
  Cluster age. Progenitor lifetimes are obtained simply as t_prog = t_clus − t_cool, so uncertainties in t_clus directly affect the determination of M_i. Central values and adopted 1σ uncertainties for the cluster ages determined for this paper are listed in Table 3. As our standard choice, we have adopted uniform distributions f_t(t_clus) for cluster age uncertainties, where f_t is a constant defined by the simple relation $\int _{-\sigma _{t {\rm clus}}}^{+\sigma _{t {\rm clus}}} f_t \,dt_{\rm clus}=0.683$.
• 2.  
  Stellar metallicity (Δ[Fe/H]). Δ[Fe/H] affects the determination of M_i in two different ways. The first and most important relates to the determination of the cluster age (mainly through the change in the distance modulus obtained via MS fitting) and has already been taken into account in the uncertainty of t_clus. Additionally, Δ[Fe/H] introduces uncertainty in the determination of the progenitor mass because, for different metallicities, a fixed t_prog corresponds to different progenitor masses. We treat this effect by adopting the central values and uncertainties in the metallicities listed in Table 3 to generate stellar metallicity distributions (normal distributions are our standard choice). For a fixed t_prog, we then compute M_i for the corresponding metallicity distribution.
• 3.  
  Different stellar models and isochrones. In order to estimate the systematic uncertainties arising from using different sets of stellar models and isochrones, we have proceeded in a similar fashion as in the case of the WD cooling tracks and define as the 1σ uncertainty

  $$\begin{equation} \sigma _{\rm syst}(M_{\rm i})=|M_{\rm i, BaSTI}-M_{\rm i, Padua}|, \end{equation} \tag{ 8 }$$

  where M_i,BaSTI represents the progenitor masses derived with the Basti stellar models and isochrones and M_i,Padua those derived with the Padua models and isochrones. We have adopted a uniform distribution for systematic uncertainties. It should be noted here again that consistency in the calculations requires the same set of models and isochrones be used to derive M_i and t_prog, i.e., the cluster age. It is in principle not correct to adopt a cluster age derived with, let us say the BaSTI isochrones, and then the Padua stellar models to derive M_i from the calculated t_prog. As already stated in the introduction, such inconsistencies are ubiquitous in the current literature and can introduce noticeable systematic changes in the derived M_i, particularly for M_i ≳ 5 M_☉.

We close this section by stressing that both codes, BaSTI and Padua, are employing up-to-date physics input, and therefore are, at least nominally, quite similar in their treatment of stellar evolution. Systematic uncertainties arising from code-to-code differences are, therefore, possibly underestimated in the above error analysis, if they arise from different physical assumptions, as, for example, the inclusion or neglect of convective overshooting, which we did not include in the discussion of this section (we provide separate results for models without core overshooting).

When considering models without core overshooting, the systematic effect of using different models and isochrones was taken into account with a slightly different approach because Padua calculations are only available for one chemical composition, i.e., Z = 0.019 and Y = 0.273. This pair of (Y, Z) values corresponds to [Fe/H] = 0.04 when considering the Grevesse & Noels (1993) solar mixture (that provides Z/X = 0.0245 ± 0.005 for the present Sun) employed in those stellar model calculations. Also, the BaSTI models employ the Grevesse & Noels (1993) solar mixture, and the (Y, Z) pair closest to the Padua one around the solar metallicity is Z = 0.0198, Y = 0.2734 that corresponds to [Fe/H] = 0.06. In this case, and just for estimating this systematic uncertainty, we have computed all progenitor masses for both BaSTI and Padua non-overshooting models assuming a fixed [Fe/H] = 0.04 for Padua and [Fe/H] = 0.06 for BaSTI. As above, the difference between the masses has been adopted as a measure of the systematic uncertainty, and we have assumed that this difference does not depend on the actual [Fe/H] value. This is, in fact, a reasonable approximation and, in addition, systematic uncertainties from isochrones and models turn out to be a small contribution to the overall uncertainty budget, so our simplified treatment should not affect our results and conclusions.

For each star in our sample, the MC simulations consist of the following steps.

• 1.  
  Generate {log g_j, T_eff,j}^N_j=1 distributions according to observational uncertainties and use S00 WD models to get the initial WD mass {M⁰_f,j}^N_j=1 and cooling age {t⁰_cool,j}^N_j=1 distributions. Here N(=10⁵) is the total number of trials in the MC simulations and quantities with the j suffix denote results of each one of the trials.
• 2.  
  For each trial j, include additional WD mass and cooling age uncertainties according to the details given in Section 3.1. Each mass in the the final WD mass distribution {M_f,j}^N_j=1 (analogously for the WD cooling age) is given by

  $$\begin{equation} M_{f,j}= M^0_{f,j} + \sum _{k=1,5} \delta M_{j,k}, \end{equation} \tag{ 9 }$$

  where each δM_j,k, the individual contributions to the WD mass global uncertainty, is randomly drawn from distributions constructed as described in Section 3.1. The final WD cooling age distribution {t_cool,j}^N_j=1 is computed in the same way.
• 3.  
  Generate the cluster age distribution {t_clus,j}^N_j=1 and then the progenitor age distribution {t_prog,j}^N_j=1, where for each trial t_prog,j = t_clus,j − t_cool,j.
• 4.  
  From {t_prog,j}^N_j=1 get the initial progenitor mass distribution {M⁰_i,j}^N_j=1 using the BaSTI stellar models.
• 5.  
  The final progenitor mass distribution $\lbrace M_{i,j} \rbrace _{\rm {j=1}}^{\rm N}$ is obtained by adding the star metallicity and stellar models systematic uncertainties, i.e., M_i,j = M⁰_i,j + δM_i,[Fe/H] + δM_i,syst as discussed in Section 3.2.

Table 5 displays the results about cooling age and mass for all the WDs in our sample. The central values have been obtained using S00 WD models and the uncertainties are the result of the MC simulations described in the previous section and thus include all sources of uncertainty considered in this work. We also provide, as additional information, the bolometric luminosity of the WDs, as obtained from the interpolation among the WD model grid in log g and T_eff. It is important to note that none of the objects have yet entered the crystallization phase. This minimizes the effect of uncertainties in the interior CO profiles that play a very important role during the phase separation associated with the crystallization process (S00).

Our reference IFMR relation (BaSTI with core overshooting plus S00 models) is displayed in Table 6 (labelled "ov"). For the WD masses, the quoted uncertainties are the average of σ₋(M_f) and σ₊(M_f) given in Table 5. For the progenitor masses, we give the results of our MC simulations that show strongly asymmetric error bars arising from the highly nonlinear relation between evolutionary lifetime and mass for MS stars. Results for the BaSTI models without core overshooting are also shown (labelled "non-ov").

As already discussed, our estimates of the global uncertainties in the IFMR include systematic effects due to the use of different sets of stellar models and isochrones, and different WD cooling models. It is instructive, however, to try to single out the influence of the stellar models on the final error. To this aim, Figure 4 shows a comparison between the progenitor masses obtained with two different sets of stellar models and isochrones (BaSTI and Padua). In both cases, we used S00 WD models. When compared with the error bars on M_i shown in Figure 3, Figure 4 leads us to conclude that systematic uncertainties from using different sets of isochrones and stellar models do not contribute appreciably to the global error budget on M_i, as long as the assumptions concerning the stellar physics are similar. We emphasize again that to achieve internal consistency the same set of isochrones and stellar models must be used to determine the cluster ages and to obtain progenitor masses from the progenitor lifetimes.

Zoom In Zoom Out Reset image size

We next study the changes in the IFMR from using different WD models. Figure 5 displays our reference IFMR (S00 + BaSTI with overshooting) and that resulting from using the LPCODE WD models. Although WD masses derived from both sets of WD models are very similar, the effect on the progenitor mass is much larger, especially for progenitor masses above ∼5 M_☉. The changes are due to differences in the WD cooling ages obtained with the different WD models. Both S00 and LPCODE models give WD ages consistent within the observational errors in log g and T_eff, however, the short lifetime of massive stars makes the determination of M_i very sensitive to small differences in t_cool for higher progenitor masses.

Zoom In Zoom Out Reset image size

Figure 6 gives an overview of the relative contributions of observational and model (cooling tracks, including systematic effects due to variations of input physics and chemical stratification) uncertainties to the total error budget for the final WD masses and ages. We display also the relative contribution of cluster age (due to isochrones, cluster metallicity, reddening, distance) and WD age (due to cooling tracks) uncertainties to the total error budget for the progenitor masses and ages. These fractional contributions are expressed in terms of (σ_i/σ_tot)² for each WD, denoted by its number (as in Table 1), whilst vertical dashed lines separate the various clusters. The uncertainty of the WD mass appears to be completely dominated by the observational errors in log g and T_eff, while in case of WD ages the systematic uncertainties of the cooling tracks dominate for some clusters (NGC 3532 and NGC 2099). For the progenitor ages and masses (lower panels of Figure 6) cluster ages are the largest contributors to the total error, except mainly for NGC 2099.

Zoom In Zoom Out Reset image size

In closing this section, we comment on the possibility that systematic uncertainties (to be added to the quoted errors) may be present in the empirical determination of WDs T_eff and log g(see, e.g., Napiwotzki et al. 1999). We have not included this source of uncertainty in our global error assessment simply because it is not quantified in the literature. We have performed, however, a simple test to understand how such hypothetical uncertainties may influence the determinations of WD masses and cooling ages. We have considered that T_eff determinations have systematic uncertainties amounting up to ±1000 K that we represent by means of a uniform distribution. Analogously, in the case of log g we have assumed an additional systematic uncertainty of ±0.15 dex. Note that these are probably very generous systematic errors, because we assume all WDs are similarly affected. Our results show that uncertainties in log g have a more noticeable effect than those on T_eff. This is because WD masses are more robustly predicted by WD models than cooling ages, and consequently the error budget is dominated by observational uncertainties (as already shown in Figure 6). For about 40% of our sample, uncertainties in WD masses are about a factor of 2 larger when systematic uncertainties as described above are included; these are the stars that show the smallest observational uncertainties in Table 1. For the rest of the stars in our sample, uncertainties are increased by not more than a factor of 1.5. Cooling ages, on the other hand, are much less affected by systematic uncertainties, because the total error budget here is heavily influenced by input physics and systematics in WD models. As a consequence, initial mass determinations are only very mildly affected. Indeed, the largest increase in positive uncertainties for M_i is only 50%, with only four objects having changes above 40%. There are 40 objects for which positive uncertainties increase by less than 25%. Negative uncertainties, which are always smaller, are less affected, and for 49 objects, increases are smaller than 25%. We, thus, conclude that adding systematic uncertainties to WD observational parameters would only have a very modest influence on the inferred IFMR. Based on these results and also considering that any attempt to quantify systematic uncertainties remains, at this point, arbitrary, we have chosen not to include these effects neither in our reference IFMR nor in the discussion of our results.

Figure 3 compares the inferred IFMR with that predicted by the same BaSTI models (we chose to display the BaSTI results for [Fe/H] = 0.06) used in the determination of t_prog and M_i. We display results obtained from progenitor models and isochrones with and without MS core overshooting.

Cluster ages determined from models without overshooting are smaller (see Table 3), hence t_prog is decreased, resulting in larger progenitor masses.

As a consequence, although the associated error bars are large, models without overshooting predict for some stars progenitors masses well above the maximum possible values allowed by theory. Examples are the WD in the Pleiades and two of the WDs in NGC 2516 (see Table 6), for which the estimated progenitor masses are between ∼12 and ∼16 M_☉. Theoretical stellar evolution models predict an upper limit to the progenitors of carbon/oxygen WDs in the range between ∼6 and ∼8 M_☉, the exact values depending on the chemical composition of the models and on the treatment of mixing in the cores; lower upper limits are obtained when core overshooting during the central H-burning phase is included. In several other cases in Table 6 even negative progenitor ages were obtained.

The relationship between M_i and the core mass at the onset of the thermal pulses (M_c1TP) is also displayed in Figure 3. For models without core overshooting, the theoretical IFMR lies largely above the semiempirical determination, especially for M_i < 5 M_☉, while the theoretical M_i − M_c1TP relationship is in better agreement with the data. This would imply that the progenitor stars did not experience any thermal pulses at all.

On the other hand, for core overshooting models the theoretical IFMR follows much better the semiempirical results, and no very high initial masses are typically derived in this case. The theoretical M_i − M_c1TP relation constitutes a lower envelope to the data, consistent with general expectations about AGB evolution. This would lead us to conclude that core overshooting during the MS results in a better agreement with the derived IFMR, because of the internal consistency about the end-product of stars with inferred M_i.

Before moving forward in our comparison with theory, we wish to comment briefly on NGC 2099 that represents an interesting case and a warning about the role played by the cluster [Fe/H] and E(B − V) estimates. As already mentioned, we have considered two possible values for its metallicity and reddening. Results for both choices are given in Table 6 and are also plotted in Figure 7. It becomes apparent from the figure that data points fall systematically below even the theoretical M_i − M_c1TP (models and IFMR with overshooting) relation for [Fe/H] = 0.09. The same result is obtained if models without core overshooting are used instead. This inconsistency is strongly alleviated if the alternative [Fe/H] and E(B − V) used by Kalirai et al. (2005; see Table 3) are adopted. While the aim of this paper is not to use the IFMR to improve cluster properties, this example demonstrates how informative it can be, when systematic effects are well controlled.

Zoom In Zoom Out Reset image size

Focusing again on the core overshooting IFMR (Figure 3, top panel), we can finally try to obtain some constraints on TP-AGB evolution by comparing the inferred IFMR with that resulting from sequences with different TP-AGB modeling. To do this, we first compare in Figure 3 the M_i − M_c1TP relations from BaSTI and LPCODE stellar models (dotted and dashed lines respectively). Both relations are very similar, regardless of differences in the treatment of convective and additional mixing episodes. We recall that in BaSTI instantaneous mixing in the overshooting region beyond the Schwarzschild boundary of the H-burning convective core is assumed, as well as semiconvective mixing during He-core burning, which is briefly described in Pietrinferni et al. (2004), where also the appropriate references are given. On the other hand, no extra mixing is taken into account at the bottom of convective envelopes. LPCODE models include diffusive mixing due to overshooting at all convective boundaries. We also point out that evolutionary timescales until the onset of the TPs are very similar in the BaSTI models with overshooting and the LPCODE models. Since TP-AGB lifetimes are very short and their contribution to the the total progenitor age can be neglected, we can safely use semiempirical data obtained from the BaSTI models for comparisons with the LPCODE theoretical IFMR.

Having assessed the similarity of the M_i − M_c1TP relations, differences in the IFMRs predicted by BaSTI and LPCODE models can be traced to differences during the TP-AGB phase. While BaSTI models do not include extra mixing during the thermal pulses, LPCODE includes, as mentioned above, diffusive overshooting at all convective boundaries. The theoretical IFMRs are also shown in Figure 3 with thick (BaSTI) and thin (LPCODE) solid lines and do indeed differ. This is unrelated to the choice of mass loss during the AGB phase, and it is mainly due to the presence of very strong third dredge-up episodes in LPCODE sequences, which inhibit the growth of the H-free core for models with M_i ≳ 3 M_☉. These strong third dredge-up episodes can be traced back mainly to the presence of overshooting at the base of the thermal pulse driven convective zone. This was demonstrated by Herwig (2000) and confirmed by us running test models where convective overshooting was included only at separately selected convective boundaries during the AGB evolution. Due to strong third dredge-up events, the final mass of LPCODE sequences with M_i ≳ 3 M_☉ is very similar to M_c1TP, the core mass at the first thermal pulse (Figure 3, top panel). On the other hand in BaSTI sequences the final mass is mainly set by the length of the TP-AGB phase and, thus, by the adopted mass-loss rate. The apparent failure of the theoretical LPCODE IFMR to match the final masses in the range 2.75 M_☉ ≲ M_i ≲ 3.5 M_☉ is an indication that such strong dredge-up episodes do not take place in real stars, and therefore the adopted exponentially decaying overshooting at the base of the thermal pulse driven convective zone in the LPCODE models is too strong. This allows us to conclude that the overshooting parameter f (see Herwig et al. 1997, for a definition of f ) must be significantly smaller than f = 0.016 employed at the convective core of upper MS stars. This is in line with recent hydrodynamical simulations of AGB thermal pulses which point to a value of f ≲ 0.01 at the bottom of the pulse driven convective zone (Herwig et al. 2007). On the other hand, given the large error bars of the data points, the existence of strong third dredge-up episodes for stars with M_i ≳ 5 M_☉ cannot be excluded by the same arguments. The IFMR may therefore teach us that overshooting varies a lot between different convective layers and may also depend on stellar mass.