STELLAR AGES AND CONVECTIVE CORES IN FIELD MAIN-SEQUENCE STARS: FIRST ASTEROSEISMIC APPLICATION TO TWO KEPLER TARGETS

A detailed comprehension of the physical processes taking place in deep stellar interiors is of paramount importance for an accurate description of stellar populations (e.g., Chiosi et al. 1992). Assumptions and parameterizations used to represent processes like convection, rotation, overshooting and microscopic diffusion directly impact the quantities normally used to characterize different types of stars, such as effective temperatures, colors, surface gravities, and composition (see, e.g., Roxburgh 1978; Vandenberg & Smith 1988; Michaud & Proffitt 1993; Maeder & Meynet 2000, and references therein). As a consequence it is extremely difficult to determine reliable stellar ages of field main-sequence stars, where uncertainties are well above the 20% level (Soderblom 2010).

One aspect of critical relevance for dating main sequence stars is the existence and size of a convective core. It has long been known (e.g., Maeder 1974) that reproducing the observed hook-like feature in the color–magnitude diagram (CMD) of intermediate-age clusters requires that stars at the turn-off phase possess a convective core during their main-sequence lifetime. Since the metallicity of cluster stars is assumed to be fixed and homogeneous, the age of the cluster can be estimated by determining the best fit of isochrones to the observed CMD (e.g., Vandenberg & Bell 1985; Salaris et al. 1997). The great caveat of this technique is that the critical mass at which a convective core appears in evolutionary models, normally around 1.1 M_☉, strongly depends on the input physics (see, e.g., Christensen-Dalsgaard & Houdek 2010; Magic et al. 2010). Moreover, the total size of the mixed central region in stellar models can be controlled by a parameter of mixing efficiency beyond the formal convective boundaries of the core; this is calibrated to match observations of clusters and is again dependent on the rest of the input physics used to construct isochrones (e.g., Maeder & Meynet 1991; Pietrinferni et al. 2004). Therefore, we currently lack a reliable estimation of the size of convective cores in main-sequence stars.

An exciting approach to overcome these impediments comes from the possibility of piercing the surface of stars by studying their pulsations, and the exquisite data currently being obtained by the Kepler mission provide an excellent opportunity to achieve this goal (Christensen-Dalsgaard et al. 2007; Borucki et al. 2009). Several theoretical studies of how asteroseismology can yield information about stellar cores have been carried out (e.g., Audard & Provost 1994; Mazumdar & Antia 2001; Mazumdar et al. 2006; Cunha & Metcalfe 2007; Cunha & Brandão 2011; Silva Aguirre et al. 2011a), but the paucity of data for main-sequence targets prior to the launch of Kepler allowed very few investigations of this type in stars other than the Sun. In fact, all previous data-based asteroseismic studies of stellar cores properties (such as Miglio & Montalbán 2005; Deheuvels et al. 2010; de Meulenaer et al. 2010) relied on ground-based observations of much shorter time coverage and lower quality than the multi-month Kepler observations.

Of the hundreds of main-sequence targets for which the Kepler mission has provided asteroseismic data (Chaplin et al. 2011), individual frequencies have now been determined for dozens of them from time-series spanning several months (e.g., Appourchaux et al. 2012). We carry out here the first investigation using Kepler data on two of these targets with the aim of detecting convective cores, constraining their size, determining accurate stellar parameters, and narrowing down the age uncertainties compared to those obtained by fitting stellar tracks.

To accomplish this, we proceed as follows. In Section 2 we describe the available seismic and spectroscopic observations, which are used in Section 3 to obtain initial stellar parameters from global asteroseismic fitting techniques. The large uncertainties obtained from this analysis, particularly in age, lead us to consider diagnostic tools sensitive to the deep interior of stars, presented in Section 4. Several teams performed detailed modeling of the targets using different evolution and pulsation codes, described in Section 5. The results of this detailed analysis are given in Section 6, together with the derived stellar parameters. We discuss our findings, as well as other possible sources of uncertainties, and give closing remarks in Section 7.

Based on the full ensemble of stars with nine months of Kepler observations available, we initially selected targets with seismic characteristics belonging to the main-sequence phase (no signature of mixed-modes; see Bedding 2011 for a thorough explanation). Since we are interested in detecting convective cores, we further discarded the targets with global seismic properties suggesting a mass value lower than solar, and finally chose two stars from the remaining targets with a large number of reliable individual frequencies determined (see Appourchaux et al. 2012). These are dubbed Perky (KIC 6106415) and Dushera (KIC 12009504).

The power spectrum of solar-like oscillators is modulated by a Gaussian-like envelope, and presents a near-regular pattern of overtone frequencies (e.g., Chaplin et al. 2011). Two parameters can be readily extracted from it, namely the frequency of maximum power of oscillation ν_max and 〈Δν〉, the average separation between the consecutive p-mode overtones

$$\begin{equation} \Delta \nu _{\ell }(n) = \nu _{n,\ell }-\nu _{n-1,\ell }, \end{equation} \tag{ 1 }$$

where ν_{n, ℓ} is the mode frequency of angular degree ℓ and radial order n. These are usually referred to as the global oscillation parameters.

Values of 〈Δν〉 and ν_max were obtained for both targets as described by Verner et al. (2011) using time series prepared from short-cadence data (Gilliland et al. 2010; García et al. 2011). Individual frequencies were extracted in the manner explained by Appourchaux et al. (2012), resulting in 33 oscillation modes for Perky and 34 for Dushera. These are listed in Table 1 and in Table A.53 of Appourchaux et al. (2012) for Perky and Dushera, respectively.

Effective temperature T_eff, surface gravity log  g and iron abundance [Fe/H] were available for both targets from high-resolution spectroscopic observations analyzed with the VWA software (Bruntt et al. 2012 and references therein). The case of Dushera presented a challenge: a set of low-resolution spectra was analyzed with a standard method based on the equivalent widths of iron lines using ARES and MOOG (see, e.g., Sousa et al. 2007). This method, described in detail in Sousa et al. (2008), yielded atmospheric parameter values that were not consistent with those obtained from VWA. Although one might be tempted to disregard the results from the lower quality data, estimations of T_eff for Dushera using the InfraRed Flux Method (IRFM; Casagrande et al. 2010; Silva Aguirre et al. 2012) suggest a value that is consistent with that obtained by ARES+MOOG. Since a thorough comparison of these spectroscopic methods is beyond the scope of this paper, we adopted average values and large uncertainties in the observed parameters that encompassed both determinations.

For the case of Perky only the results from Bruntt et al. (2012) were available and the IRFM T_eff determination is consistent with this value. However, based on our experience with Dushera, we decided to use more conservative uncertainties than those quoted by Bruntt et al. (2012). Since parallax measurements are also available for this star, we estimated its luminosity using the Hipparcos data (Van Leeuwen 2007) and bolometric corrections (Flower 1996; Torres 2010). The complete set of observational constraints considered in this study is given in Table 2.

An initial determination of ages can be made using these spectroscopic parameters and stellar isochrones following a Bayesian analysis (Serenelli et al. 2013). For the values quoted in Table 2, this technique yields age uncertainties of 50% and 60% for Perky and Dushera, respectively. Even when reducing the errors in T_eff to 100 K and considering errors in log  g of 0.1 dex for both targets, the obtained uncertainties only decrease to 36% and 50%. Ages determined with a precision of ∼20% can only be achieved with this technique for targets with accurately determined observational constraints, which are not available for these relatively faint field main-sequence stars (see Section 5.2 in Serenelli et al. 2013). In fact, age results for main-sequence intermediate-mass field stars have often uncertainties of 20%–50% when determined from a single stellar evolution database, that is, before systematics are taken into account (see Section 4.1 in Soderblom 2010).

Using the parameters given in Table 2, it is possible to apply global asteroseismic techniques and determine the stellar properties of the targets. 〈Δν〉 approximately scales as the square root of the mean density (e.g., Ulrich 1986), while ν_max seems to be related to the acoustic cutoff frequency of the atmosphere and thus to the surface gravity and effective temperature of the star (e.g., Brown et al. 1991; Kjeldsen & Bedding 1995; Belkacem et al. 2011). Based on these dependencies, two scaling relations anchored at the solar parameters can be written

$$\begin{equation} \frac{M}{M_\odot } \simeq \left(\frac{\nu _{\mathrm{max}}}{\nu _{\mathrm{max},\odot }}\right)^{3} \left(\frac{\langle \Delta \nu \rangle }{\Delta \nu _\odot }\right)^{-4}\left(\frac{T_\mathrm{eff}}{T_{\mathrm{eff},\odot }}\right)^{3/2}, \end{equation} \tag{ 2 }$$

$$\begin{equation} \frac{R}{R_\odot } \simeq \left(\frac{\nu _{\mathrm{max}}}{\nu _{\mathrm{max},\odot }}\right) \left(\frac{\langle \Delta \nu \rangle }{\Delta \nu _\odot }\right)^{-2}\left(\frac{T_\mathrm{eff}}{T_{\mathrm{eff},\odot }}\right)^{1/2}. \end{equation} \tag{ 3 }$$

Here, T_{eff, ☉} = 5777 K, Δν_☉ = 135.1 ± 0.1 μHz, and ν_{max, ☉} = 3090 ± 30 μHz are the observed values in the Sun (Huber et al. 2011). Provided a measurement of T_eff is available, these relations give a determination of mass and radius independent of stellar models (the direct method, see, e.g., Miglio et al. 2009; Chaplin et al. 2011; Silva Aguirre et al. 2011b). It is also possible to search for a best fit to these parameters within pre-calculated sets of evolutionary tracks, where chemical composition can be taken into account (the grid-based method, e.g., Stello et al. 2009; Quirion et al. 2010; Gai et al. 2011; Basu et al. 2012). The latter approach also provides a determination of the stellar age.

We performed a grid-based search using several available sets of evolutionary models to make an initial estimation of the stellar parameters of both targets. Details on how the sets of evolutionary tracks were constructed and the different search methods can be found in the following references: Basu et al. (2010) and Gai et al. (2011) for the YY, YREC, Dotter, and Marigo results; Stello et al. (2009) for RADIUS; Silva Aguirre et al. (2012) for GARSTEC; Creevey et al. (2012) for RadEx10; and Quirion et al. (2010) for SEEK. One grid-based pipeline built using the BaSTI²² isochrones (Pietrinferni et al. 2004) is presented here for the first time.

Briefly, the BaSTI pipeline uses a denser set of isochrones in metallicity space than the publicly available ones, computed specially for the latest revision of the Geneva–Copenhagen Survey (Casagrande et al. 2011) and now extended for asteroseismic analysis. They are built assuming the Grevesse & Noels (1993) solar composition and adopting a value of ΔY/ΔZ = 1.45, in agreement with big bang nucleosynthesis (SBBN) values of primordial helium abundance at low metallicities (Y_P = 0.2482 ± 0.0007; Steigman 2007, 2010). The 22 available chemical mixtures cover a metallicity range between −3.27 < [Fe/H] <+0.5. Models used here do not include overshooting and gravitational settling, and are calculated assuming mass loss according to the Reimers (1975) formulation with a fixed efficiency parameter of η = 0.4. The algorithm to find the best match to the data is that of the GARSTEC grid, based on Monte Carlo realizations to form the probability distribution of each parameter.

Table 3 shows the results obtained with each grid pipeline for the two targets, including their 1σ uncertainties. As mentioned above, some methods return a probability distribution for the stellar parameters characterized by asymmetric error bars. The last column shows, for comparison, the values obtained with the direct method (Equations (2) and (3)). There is an overall good level of agreement in the results for log  g and R/R_☉, explained by the fact that the global asteroseismic parameters are mostly sensitive to ratios between stellar mass and radius (mean density and surface gravity). Due to these strong dependencies, stellar radii and log  g as determined from seismology are very precise (Bedding 2011; Morel & Miglio 2012; Silva Aguirre et al. 2012; Huber et al. 2012).

Table 3. Stellar Parameters as Determined from Each Pipeline for the Two Targets

Star	Parameter	YY	Dotter	Marigo	YREC	RADIUS	GARSTEC	SEEK	RadEx10	BaSTI	Direct
Perky	M(M☉)	$1.06^{+0.09}_{-0.07}$	$1.06^{+0.08}_{-0.07}$	$1.12^{+0.05}_{-0.06}$	$1.09^{+0.07}_{-0.06}$	1.10 ± 0.08	$1.11^{+0.10}_{-0.10}$	$1.03^{+0.10}_{-0.10}$	1.08 ± 0.09	$1.07^{+0.06}_{-0.10}$	1.10 ± 0.07
	R(R☉)	$1.21^{+0.03}_{-0.03}$	$1.21^{+0.03}_{-0.03}$	$1.23^{+0.02}_{-0.02}$	$1.22^{+0.03}_{-0.02}$	1.22 ± 0.03	$1.24^{+0.01}_{-0.04}$	$1.20^{+0.04}_{-0.03}$	1.22 ± 0.03	$1.22^{+0.02}_{-0.04}$	1.23 ± 0.02
	log  g	$4.30^{+0.01}_{-0.01}$	$4.30^{+0.01}_{-0.01}$	$4.31^{+0.01}_{-0.02}$	$4.31^{+0.01}_{-0.01}$	4.30 ± 0.01	$4.30^{+0.01}_{-0.01}$	$4.29^{+0.02}_{-0.02}$	4.30 ± 0.01	$4.30^{+0.01}_{-0.01}$	4.30 ± 0.01
	Age (Gyr)	$5.73^{+2.38}_{-2.44}$	$6.00^{+2.44}_{-2.18}$	$3.71^{+1.80}_{-1.22}$	$6.21^{+2.24}_{-2.15}$	6.77 ± 2.20	$4.21^{+3.36}_{-2.21}$	$5.52^{+4.54}_{-3.66}$	6.37 ± 3.18	$5.80^{+3.70}_{-2.20}$	...
Dushera	M(M☉)	$1.25^{+0.08}_{-0.10}$	$1.21^{+0.09}_{-0.08}$	$1.23^{+0.07}_{-0.05}$	$1.27^{+0.10}_{-0.10}$	1.21 ± 0.10	$1.26^{+0.09}_{-0.09}$	$1.16^{+0.14}_{-0.13}$	1.26 ± 0.09	$1.19^{+0.06}_{-0.08}$	1.22 ± 0.08
	R(R☉)	$1.43^{+0.04}_{-0.04}$	$1.42^{+0.03}_{-0.03}$	$1.42^{+0.03}_{-0.02}$	$1.44^{+0.04}_{-0.04}$	1.41 ± 0.05	$1.44^{+0.03}_{-0.03}$	$1.41^{+0.05}_{-0.05}$	1.44 ± 0.04	$1.41^{+0.02}_{-0.03}$	1.42 ± 0.03
	log  g	$4.22^{+0.01}_{-0.01}$	$4.22^{+0.01}_{-0.01}$	$4.22^{+0.01}_{-0.01}$	$4.22^{+0.01}_{-0.01}$	4.22 ± 0.01	$4.22^{+0.01}_{-0.01}$	$4.21^{+0.02}_{-0.02}$	4.22 ± 0.01	$4.22^{+0.01}_{-0.01}$	4.22 ± 0.01
	Age (Gyr)	$3.21^{+1.83}_{-1.20}$	$3.82^{+1.68}_{-1.36}$	$2.94^{+1.06}_{-1.05}$	$3.78^{+2.13}_{-1.54}$	6.43 ± 2.06	$2.89^{+1.61}_{-1.07}$	$3.75^{+3.83}_{-2.38}$	3.74 ± 1.71	$4.05^{+1.81}_{-1.05}$	...

Notes. Also shown are the values obtained with the direct method. See the text for details.

Download table as:  ASCIITypeset image

On the other hand, asteroseismic mass determinations still require verification from independent methods. The mean values from the grids are 1.08 M_☉ for Perky and 1.23 M_☉ for Dushera. All masses obtained by the nine pipelines are contained within a 10% range from these averages. In spite of that, when taking into account the 1σ uncertainties given by each grid, mass values between ∼0.95–1.20 M_☉ and ∼1.05–1.35 M_☉ are possible for Perky and Dushera, respectively. The large age differences in the results (∼50% or more) are a direct consequence of this mass difference (∼0.3 M_☉). Moreover, the fact that the results are in the vicinity of the ∼1.1 M_☉ limit where the onset of core convection is thought to occur suggests that the existence and extent of a convective core will play a significant role in the internal structure of the star.

These results are to be taken with some caution since the scaling relations given in Equations (2) and (3) are not exact. Moreover, they are of course dependent on the uncertainties in asteroseismic (〈Δν〉 and ν_max) as well as atmospheric (T_eff and metallicity) inputs fed into the grid searches, which we decided to keep at a conservative level as explained in Section 2. In fact, each grid of models has been constructed using different sets of input physics, solar mixture, opacities, equation of state, etc., which accounts for some of the differences seen in the results (see W. J. Chaplin et al. 2013, in preparation, for a thorough discussion of the grid-to-grid scatter). Bearing this in mind it is clear that a detailed asteroseismic study is mandatory to determine masses with higher accuracy, to explore the energy transport mechanisms in the inner structure, and to better constrain the age of the targets. The lack of sensitivity of the global seismic analysis to the finer details of the stellar structure needs to be compensated with tools capable of probing the deep interior of stars.

Different combinations of oscillation frequencies have been proposed in the literature to study the internal structure of stars. The two most commonly used are the large frequency separation (Equation (1)) and the small frequency separation defined as

$$\begin{equation} d_{\ell,\ \ell +2}(n) = \nu _{n,\ell }-\nu _{n-1,\ell +2}\,. \end{equation} \tag{ 4 }$$

The latter quantity is sensitive to the conditions in the stellar core and thus to the central hydrogen content (age) of main-sequence stars (e.g., Ulrich 1986; Gough 1987; Gabriel 1989; Roxburgh & Vorontsov 1994b). In contrast, the large frequency separation is known to be affected by the outer layers of stars (see, e.g., Christensen-Dalsgaard et al. 1988; Christensen-Dalsgaard & Pérez Hernández 1992; Roxburgh & Vorontsov 2003; Roxburgh 2010), where stellar models use simple versions of convection theory instead of properly dealing with turbulent pressure and non-adiabatic effects (e.g., Balmforth 1992; Rosenthal et al. 1999; Ballot et al. 2004; Houdek 2010). For the solar case, differences in Δν around the ν_{max, ☉} region between Birmingham Solar Oscillations Network (Chaplin et al. 1999) data and Model S from Christensen-Dalsgaard et al. (1996) are of the order of 1 μHz, much larger than the typical solar observational uncertainties of σ_Δν ∼ 0.06 μHz.

In order to study the deep interior of the star we need combinations of frequencies that isolate the signal arising from the stellar center as much as possible from surface contamination. Roxburgh & Vorontsov (2003) proposed to use ratios of small to large separations, which are mainly determined by the inner structure of the star (Roxburgh & Vorontsov 2003; Otí Floranes et al. 2005; Roxburgh 2009). They are constructed as

$$\begin{equation} r_{02}(n)=\frac{d_{0,2}(n)}{\Delta \nu _{1}(n)} \end{equation} \tag{ 5 }$$

$$\begin{equation} r_{01}(n)=\frac{d_{01}(n)}{\Delta \nu _{1}(n)},\quad r_{10} = \frac{d_{10}(n)}{\Delta \nu _{0}(n+1)}, \end{equation} \tag{ 6 }$$

where d₀₁(n) and d₁₀(n) are the smooth five-points small frequency separations:

$$\begin{equation} d_{01}(n)=\frac{1}{8}(\nu _{n-1,0}-4\nu _{n-1,1}+6\nu _{n,0}-4\nu _{n,1}+\nu _{n+1,0}) \end{equation} \tag{ 7 }$$

$$\begin{equation} d_{10}(n)=-\frac{1}{8}(\nu _{n-1,1}-4\nu _{n,0}+6\nu _{n,1}-4\nu _{n+1,0}+\nu _{n+1,1})\,. \end{equation} \tag{ 8 }$$

The ratios of degree 0 and 1 (Equation (6)) have already been used to identify the location of the convective envelope in the Sun (Roxburgh 2009), and also in other main-sequence targets (e.g., Lebreton & Goupil 2012; Mazumdar et al. 2012). Inspection of Equations (7) and (8) shows that r₀₁(n) and r₁₀(n) are alternatively centered in ν_{n, 0} and ν_{n, 1}, respectively. Therefore we can consider the ratios r₀₁(n) and r₁₀(n) as one unique set of observables and as a function of frequency which we call r₀₁₀:

$$\begin{eqnarray} r_{010}&=&\lbrace r_{01}(n),r_{10}(n),r_{01}(n+1),r_{10}(n+1),r_{01}(n+2),\nonumber\\ &&\ r_{10}(n+2),\ldots \rbrace. \end{eqnarray} \tag{ 9 }$$

Existence of a growing convective core in stellar models produces an abrupt discontinuity in the chemical composition at the boundary between the core and the radiative envelope. The formal location of this boundary is determined by a stability criterion, most commonly the Schwarzschild criterion. When additional mixing such as overshooting is included in stellar models, the transition between the homogeneously mixed center and the radiative zone does not occur at the formal edge of the convective core but at a certain distance from it determined by the adopted mixing efficiency. In either case, the discontinuity in chemical composition (and thus adiabatic sound speed) produces a shift in the frequencies relative to those in a corresponding star with a smooth sound-speed profile that is an oscillatory function of the frequency itself, and whose period is related to the location of the discontinuity (e.g., Roxburgh & Vorontsov 1994a, 2007; Cunha & Metcalfe 2007).

It is customary to represent the radial coordinate in stellar interiors by the wave travel time from the surface to the center, called the acoustic depth τ, or its alternative representation named the acoustic radius:

$$\begin{eqnarray} \tilde{t} & = & \int _{0}^{r}\frac{dr}{c_s}, \end{eqnarray} \tag{ 10 }$$

where c_s is the adiabatic sound speed. This quantity measures the travel time from the center toward the surface. Considering that the total acoustic radius is given by $\tilde{t}_\mathrm{tot} = \tilde{t}(R_{\rm tot})$, then $\tau = \tilde{t}_\mathrm{tot} - \tilde{t}$. If the sound-speed discontinuity is located at, say, r₁ in radial coordinates, and that same position is represented in acoustic radius and acoustic depth by $\tilde{t}_1$ and τ₁ respectively, the periods of the oscillation in frequency space induced by this sharp variation are 1/(2$\tilde{t}_1$) and 1/(2τ₁), according to the considered seismic coordinate (e.g., Montgomery et al. 2003).

We stress the fact that the frequency combinations will be sensitive to the position of the transition between the homogeneously mixed and radiative region, regardless of whether or not this coincides with the formal boundary of the convective core. If the ratios were only affected by the location of the outer edge of the mixed region, information about its location could be directly extracted from the asteroseismic data. Nonetheless, the small acoustic radius characterizing the size of the convective core introduces an oscillatory signal with a period much larger than the observed frequency range. Moreover it has been shown that, although the ratios are mostly sensitive to the presence and size of a convective core, they are also affected by the central hydrogen content (Brandão et al. 2010; Silva Aguirre et al. 2011a) and still contain the signal of the base of the convective envelope (Roxburgh & Vorontsov 2003; Roxburgh 2010).

There are several other ingredients of stellar evolution that play more or less important roles in the determination of the sound speed. Use of different opacities, equation of state, and convective theory applied affects the internal temperature stratification in stellar models. When using frequency ratios one must bear in mind that their behavior reflects the interplay between evolutionary stage (through the central hydrogen content) and a central convective region (when it is present), while the rest of the microphysics affects them to a much lesser (but not negligible and still not fully quantified) extent. A detailed study of these issues goes beyond the scope of this paper and will be presented elsewhere (V. Silva Aguirre et al. 2013, in preparation).

Figure 1 shows the observed large frequency separations and ratios r₀₁₀ for the two targets considered in this study. Inspection of the lower panel shows that the stars have different mean values and slopes in their ratios. Moreover, the oscillatory component introduced by the position of the convective envelope is clearly visible in the ratios of Perky.

Zoom In Zoom Out Reset image size

The behavior of the frequency ratios in both targets shows a sudden increase at high frequencies that occurs regardless of the length of the data used and is not predicted by models. We have compared individual frequencies extracted from data sets of different lengths and confirmed that, due to the low signal-to-noise ratio and the larger line width at the high-frequency end, this behavior is the result of realization noise. For this reason, we will not consider the frequency ratios on this regime when comparing to stellar models. In the lower panel of Figure 1 we mark the frequency range that will be taken into account to study r₀₁₀ (the linear range, see Silva Aguirre et al. 2011a). The range considered to fit the ratios r₀₂ has also been restricted following the same criterion.

Together with 〈Δν〉 the average small frequency separation is used to construct an asteroseismic H-R diagram known as the C-D diagram (Christensen-Dalsgaard 1984, 1988). A slightly modified version of this diagram is constructed by using the ratios r₀₂ instead of the small frequency separations because it is a more effective indicator of stellar age (White et al. 2011). In Figure 2 we compare the position of Perky and Dushera in this diagram with evolutionary tracks at solar metallicity, showing that for a given composition and input physics the age of the targets is well constrained by r₀₂. Also plotted in the figure are tracks with a metallicity ∼0.2 dex lower, showing that the offset due to composition has a strong impact in the derived ages from this diagram. This is particularly important for Dushera, as its position is closer to where the tracks collapse in this diagram when entering the subgiant phase of evolution.

Zoom In Zoom Out Reset image size

Asteroseismic modeling of individual targets has experienced a leap forward thanks to the accurate frequencies extracted from CoRoT and Kepler data (e.g., Barban et al. 2009; Ballot et al. 2011; Campante et al. 2011; Appourchaux et al. 2012). The usual strategy applied in these types of studies consists of calculating different sets of one-dimensional hydrostatic stellar models. Using these as the equilibrium structures for a code capable of computing frequencies of oscillations, it is possible to produce a theoretical set of the spectroscopic and asteroseismic observables. Finally, a criterion must be applied that defines which set of calculations best reproduces the available data (e.g., Metcalfe et al. 2010; Brandão et al. 2011; Mathur et al. 2012; Doğan et al. 2013).

There are many ingredients in the computations of stellar models that are still not fully constrained, both in the microphysics (e.g., equation of state, opacities, nuclear reactions) and in the macrophysics (e.g., convective energy transport) of stellar evolution (see, e.g., Cassisi 2009; Kippenhahn et al. 2013). Even when the same assumptions are used in different evolutionary codes, the resulting models can differ in their global properties, internal structure (see, e.g., Lebreton et al. 2008; Montalbán et al. 2008), and also in the theoretical frequencies obtained (Moya et al. 2008).

For these reasons, several teams performed detailed modeling to match the stellar parameters presented in Table 2 and individual frequencies, using the stellar evolution and pulsation code of their choice (see Section 5.1 below). This allows us to obtain an at least qualitative measure of the systematic uncertainties inherent in the full modeling process, arising from the different input physics and implementations of search methodologies. A quantitative investigation of the individual influence of the various differences between the models presented here is beyond the scope of this paper.

Techniques to find the best-fit model to a given set of observables usually include minimization of a reduced χ² (e.g., Deheuvels & Michel 2011; Brandão et al. 2011; Gruberbauer et al. 2012) defined as

$$\begin{equation} \chi ^2 = \frac{1}{N} \sum _k \left(\frac{x_{\mathrm{obs}}(k) - x_{\rm model}(k)}{\sigma (x_{\rm obs}(k))} \right)^2,\ \end{equation} \tag{ 11 }$$

where k runs over the observable parameters, x_obs(k) represent each observable with its uncertainty σ(x_obs(k)), x_model(k) are the theoretically produced values, and N is the total number of observations. In the analysis of the models that best fit the data, we will consider χ² values for individual frequencies and ratios. It is clear from Equations (7) and (8) that the differences entering the ratios r₀₁₀ are strongly correlated. We take this into account by including their covariance term in the χ²:

$$\begin{equation} \chi ^2 = \frac{1}{N} \left(\boldsymbol {x}_{\mathrm{obs}} - \boldsymbol {x}_{\rm model}\right)^{T} \, {\bf C}^{-1} \, \left(\boldsymbol {x}_{\mathrm{obs}} - \boldsymbol {x}_{\rm model}\right),\ \end{equation} \tag{ 12 }$$

where C is the covariance matrix calculated after perturbing the observed frequencies using Monte Carlo simulations. To calculate the goodness of fit for the ratios r₀₁₀ we use Equation (12), which reduces to Equation (11) for the individual frequencies and ratios r₀₂ as they are independent.

As mentioned in Section 4, theoretical frequencies of oscillations are not able properly to reproduce the observed ones due to our poor understanding of the physical processes taking place in the outer layers of stars (see, e.g., Christensen-Dalsgaard et al. 1996, for the solar case). There are different formulations to correct for this offset (see, e.g., the supplementary materials for Carter et al. 2012 and references therein), the most commonly used is that of Kjeldsen et al. (2008). It consists of an empirical power law derived using helioseismic data and standard solar models, which can be written as

$$\begin{equation} \nu _\mathrm{obs}(n)-r\nu _\mathrm{ref}(n)=a\left(\frac{\nu _\mathrm{obs}(n)}{\nu _0}\right)^b\,. \end{equation} \tag{ 13 }$$

Here ν_obs(n) correspond to the observed frequencies and ν_ref(n) are the theoretical ones. To apply this surface correction one initially has to specify the values of the exponent b and the reference frequency ν₀. Then r is determined after calculating an observed and theoretical average of the large frequency separation, obtained using a linear fit to the individual ℓ = 0 modes as a function of radial order n. The value of a is determined from the average theoretical and observed ℓ = 0 frequencies (see Equations (6) and (8)–(10) in Kjeldsen et al. 2008). Although most modeling teams use this formulation in their analysis, they determine the parameters involved in the correction in slightly different manners (see Section 5.2 below).

Another point where discrepancies can arise from is the estimation of the surface iron abundance in stellar models, to be compared with the [Fe/H] value determined from spectroscopic observations. Although the issue of the solar surface abundances is still a matter of debate (see, e.g., Grevesse et al. 2010; Caffau et al. 2010), the observed [Fe/H] in stars barely depends on the solar reference considered because the largest discrepancies between the different solar mixtures come from the C, N, O, and Ne elements (e.g., Serenelli et al. 2009). On the other hand, in stellar models it is assumed that

$$\begin{equation} \log (Z/X) - \log (Z/X)_\odot \simeq \mathrm{[Fe/H]}, \end{equation} \tag{ 14 }$$

where X and Z are respectively the mass fraction of hydrogen and all elements heavier than helium, and log (Z/X)_☉ is the surface solar value. Equation (14) clearly depends on the chosen solar mixture, and assumes that the fraction of each element comprising Z is distributed according to the selected set of solar abundance ratios. The commonly adopted values of (Z/X)_☉ vary for different mixtures between 0.0245 (Grevesse & Noels 1993) and 0.0181 (Asplund et al. 2009). This means that for models with the same (Z/X) at the surface but computed with a different set of solar abundances, the theoretically determined [Fe/H] can vary by ∼0.13 dex, comparable to the uncertainties quoted in Table 2 from spectroscopic determinations.

Due to the aforementioned issues, and the large error bars in the spectroscopically determined [Fe/H] and T_eff (see Table 2), we will not include these constraints in the calculations of goodness of fit via χ². These observables were used by all teams as guidelines when exploring the parameter space (except for models G; see Section 5.1.7 below), and will be considered as reference for comparison that only in extreme cases of discrepancy could be used to discard models.

5.1. Evolutionary Calculations and Search Algorithms

In the following sections we describe the evolutionary code, input physics, and fitting techniques adopted by each modeling team to find the best match for both targets.

5.1.1. Models A

Models A were computed using the Garching Stellar Evolution Code (GARSTEC; Weiss & Schlattl 2008). We adopted the 2005 version of the OPAL equation of state (EOS; Rogers et al. 1996; Rogers & Nayfonov 2002) complemented by the MHD equation of state for low temperatures (Hummer & Mihalas 1988), low-temperature opacities from Ferguson et al. (2005) and OPAL opacities for high temperatures (Iglesias & Rogers 1996), the Grevesse & Sauval (1998) solar mixture, and the NACRE (Angulo et al. 1999) compilation for thermonuclear reaction rates with the updated cross section for ¹⁴N(p, γ)¹⁵O from Formicola et al. (2004). Due to the large uncertainty in the surface composition of Dushera (see Table 2), diffusion of helium and heavy elements was included only for the models of Perky using the Thoul et al. (1994) prescription. The atmospheric stratification is given by Eddington's gray plane-parallel formalism.

Convective zones are treated with the mixing-length theory (MLT) as described by Kippenhahn et al. (2013), using the solar calibrated value of the convective efficiency α_mlt = 1.791. Overshooting is implemented in GARSTEC as a diffusive process consisting of an exponential decline of the convective velocities within the radiative zone (Freytag et al. 1996). The diffusion constant is given by

$$\begin{equation} D_{\mathrm{ov}}\left(z\right) = D_0 \ \exp \ \left(\frac{-2z}{\xi H_p}\right), \end{equation} \tag{ 15 }$$

where ξ corresponds to an efficiency parameter, H_p is the pressure scale height, z the distance from the convective border, and the diffusion constant D₀ is derived from MLT-convective velocities. Using open clusters, the calibrated value for the overshooting efficiency is defined as ξ = 0.016. To prevent the survival of the pre-main sequence convective core due to overshooting in the solar model, a geometric factor for efficiency restriction is used when convective cores are small (Schlattl & Weiss 1999; Magic et al. 2010). When a convective core was present, we explored overshooting efficiencies between ξ = 0.004–0.035.

Oscillation frequencies of the specific models analyzed were computed using the Aarhus Adiabatic Oscillation Package (ADIPLS; Christensen-Dalsgaard 2008a). Models were selected by matching the observed large frequency separation and minimizing the χ² fit to both frequency ratios r₀₁₀ and r₀₂. No surface correction was applied to define the best fit to the data.

5.1.2. Models B

Models B were computed using the Aarhus Stellar Evolution Code (ASTEC; Christensen-Dalsgaard 2008b). The models are calculated using the NACRE compilation of nuclear reaction rates and the Grevesse & Noels (1993) solar mixture. The EOS is the OPAL 2005 version, while the opacities are those of OPAL for high temperatures and Ferguson et al. (2005) for low temperatures. Convection is treated using the MLT as described by Böhm-Vitense (1958), and no microscopic diffusion is considered.

Core overshooting is implemented as an extension of the central mixed region using the radiative temperature gradient, and controlled by an efficiency parameter α_ov

$$\begin{equation} d_\mathrm{ov}=\alpha _\mathrm{ov} \mathrm{min}(r_\mathrm{cc},H_p), \end{equation} \tag{ 16 }$$

where r_cc is the radial size of the convective core. This prescription restricts overshooting if convective cores are smaller than H_p. Values of α_ov = 0.0, 0.1, and 0.2 were explored in the calculations of models B, similar to the calibrated efficiencies of other codes with the same prescription (see models F below).

The atmospheric stratification is given by the T–τ_op relation of model C in Vernazza et al. (1981). Oscillation frequencies are computed using ADIPLS to find the best fit to the data by minimizing the χ² values of individual frequencies and frequency ratios after applying the Kjeldsen et al. (2008) surface correction. There is one subtle difference in the way the surface correction is applied in models B, which will turn out to be of particular importance in the case of Perky (cf. Equation (13)). To obtain the average large frequency separation needed to compute r, a fit is made to all available frequencies instead of only the radial modes. Similarly, radial and non-radial modes are used in the average frequencies used to obtain the a parameter.

5.1.3. Models C

Models C were also computed using ASTEC with the same input physics as described for models B, but no overshooting was included. A grid of models was calculated covering the following range in the mass-initial composition space: M = 1.0–1.6 M_☉, Y = 0.24–0.32, and Z/X = 0.01–0.07, while the mixing length parameter α_MLT was fixed to 1.8. Frequencies were computed with ADIPLS for the models that have [Fe/H] and log  g within 3σ, and T_eff within 1σ from the values given in Table 2, since the uncertainties on T_eff were relatively conservative. The best-fitting models were selected based on minimization of the χ² value for only the individual frequencies after applying the Kjeldsen et al. (2008) empirical surface correction.

5.1.4. Models D

Models D were calculated using the CESAM2K evolutionary code (Morel 1997; Morel & Lebreton 2008), with the OPAL 2005 EOS and opacities, the Grevesse & Noels (1993) solar mixture, and microscopic diffusion following Michaud & Proffitt (1993). Convection is modeled using the Canuto et al. (1996) formulation, which for this case gives a solar calibrated value for the convective efficiency of α_conv = 0.64. Nuclear reactions are those from the NACRE compilation, while the atmosphere is described by Eddington's gray law and connected to the stellar interior at an optical depth of τ_op = 10. Core overshooting is implemented using an adiabatic temperature stratification and efficiency as given by Equation (16).

Oscillation frequencies are computed using the Liège Oscillation code (LOSC; Scuflaire et al. 2008). The best-fit model is found via optimizations following the Levenberg–Marquardt algorithm, using spectroscopic constraints and matching individual frequencies after applying the Kjeldsen et al. (2008) surface correction. The free parameters adjusted during the procedure are mass, age, initial helium abundance and metallicity, and the convective and overshooting efficiencies. The same procedure was used by Miglio & Montalbán (2005) using the frequency ratios r₀₂ instead of individual frequencies.

5.1.5. Models E

Models E were computed using the Asteroseismic Modeling Portal (AMP; e.g., Metcalfe et al. 2012; Mathur et al. 2012). It is based on a genetic algorithm for minimization (Metcalfe & Charbonneau 2003) that matches the observed individual frequencies and spectroscopic constraints to those produced by a stellar evolution and pulsation codes. Its current implementation uses ASTEC for the evolutionary calculations, with similar input physics as models B but low-temperature opacities from Alexander & Ferguson (1994) and nuclear reactions from Bahcall et al. (1995). No overshooting is included.

Theoretical frequencies are computed with ADIPLS and corrected using the Kjeldsen et al. (2008) formulation. The best-fit model is chosen by minimizing the sum of two separate χ² values: one determined for the individual frequencies and another one for the spectroscopic constraints.

5.1.6. Models F

Models F were computed using the Yale Rotation and Evolution Code (Demarque et al. 2007) in its non-rotating configuration. The code uses low-temperature opacities from Ferguson et al. (2005) and OPAL opacities for high temperatures, OPAL EOS, the Grevesse & Sauval (1998) solar mixture, and Adelberger et al. (1998) nuclear reactions using the updated ¹⁴N(p, γ)¹⁵O cross section from Formicola et al. (2004). Models include helium and heavy-element diffusion (Thoul et al. 1994), and overshooting as given by Equation (16) using a fixed value of α_ov = 0.2 and an adiabatic temperature gradient when a convective core exists. This value of overshooting efficiency is motivated by the fit to the open cluster NGC 2420, and described in Demarque et al. (1994).

The atmospheric stratification is given by an Eddington T–τ relation, and frequencies calculated with the pulsation code used by Antia & Basu (1994). In order to find the best fit to the data, sequences of models were calculated with masses between 0.94–1.12 M_☉ for Perky and 1.08–1.24 M_☉ for Dushera. For each mass, we considered three different values of the mixing-length parameter, α_mlt = 1.55, 1.7, and 1.826 (the solar calibrated value). The best fit was found using the reduced χ² determined for the individual frequencies and frequency ratios assuming the solar surface term described in the supplementary material of Carter et al. (2012).

5.1.7. Models G

Models G were constructed using ASTEC with the same input physics as models B, exploring different values of overshooting efficiency α_ov = 0.0, 0.1, 0.2 and metallicity −0.5 < [Fe/H] < 0.5. No diffusion of elements was included in these calculations. A grid of models spanning a wide region of the Hertzsprung–Russell diagram was considered, with masses varying from 1.0 to 1.6 M_☉. ADIPLS was used to compute the oscillation frequencies and the best models were selected solely based on their asteroseismic characteristics, in particular, their r₀₂ and r₀₁₀ ratios, the frequency derivatives (slopes) of the ratios and the large frequency separations. In short, the range of possible models was first reduced based on the following criteria: that the slopes of r₀₂ and r₀₁₀ were within the 1σ observed values and that the large frequency separation laid between the observed value ±3 μHz. The identification of the best model within this subset was based on the minimum of the χ² computed by comparing the ratios r₀₂ and r₀₁₀ with their observed values.

5.1.8. Models H

Models H were extracted from a grid of models calculated using GARSTEC, described in detail in Silva Aguirre et al. (2012). Oscillation frequencies were computed using ADIPLS, and the best fits to the data were obtained looking for the minimum χ² of the individual frequencies after applying the Kjeldsen et al. (2008) surface correction. No further optimization beyond the pre-computed grid was made. The input physics is the same as for models A with the following exceptions. The EOS is the FreeEOS (A. W. Irwin,²³ see also Cassisi et al. 2003), hydrogen burning nuclear reaction rates are from Adelberger et al. (2011), and α_mlt = 1.811 is fixed from the calibration of a standard solar model. Overshooting is implemented as in Equation (15) with a fixed efficiency of ξ = 0.02 (also constrained geometrically in the case of small convective cores and thin shells) and microscopic diffusion of helium and metals has been included in all these models.

5.2. Applying the Surface Correction

In the previous sections we showed that search strategies and criteria for defining which model best reproduces the data are not unique and vary among teams. Some consider only individual frequencies while others also take into account frequency combinations.

Figure 3 depicts the effect of applying the surface correction in one of the models of Perky (cf. Equation (13)). The upper panel shows the difference between corrected and uncorrected frequencies, presenting the same functional form as the solar case to which it is calibrated. The differences in the region of observed modes can be of up to ∼20 μHz, a factor of 100 or more larger than the typical uncertainties in frequencies.

Zoom In Zoom Out Reset image size

A different picture arises when looking at the difference between the ratios r₀₁₀ of corrected and uncorrected frequencies (lower panel in Figure 3). Discrepancies are at the 10⁻⁶ level, three orders of magnitude below the average uncertainties, showing that this combination effectively cancels out the contribution from the outer layers (Roxburgh & Vorontsov 2003). The value of the reduced χ², when considering only individual frequencies, is thus heavily affected by the use of the surface correction, while the same is not true for the ratios.

Those teams applying the empirical surface correction of Kjeldsen et al. (2008) consider somewhat different implementations of it. Before discussing how the parameters involved in the correction are actually obtained by each modeling team, a few words on the assumptions behind this formulation are necessary. It requires the value of r to be close to unity, meaning that the mean stellar density of the model actually reproduces that of the real star. It considers only radial modes when finding the parameters for the corrections, assuming that non-radial frequencies are equally affected by the surface layers. Finally, usage of a single power law is justified by considering the offsets in the mode frequencies to be caused solely by the properties of the outer stellar layers, and thus dependent on surface gravity, effective temperature, and surface composition. As a consequence of these assumptions models found using this empirical treatment are expected to reproduce the global stellar properties (mass and radius), but not necessarily the core structure nor provide reliable age determinations (see Section 5 in Kjeldsen et al. 2008). Considering also non-radial modes in the fit is therefore necessary to ensure a proper representation of the deep stellar interior.

Obtaining the parameters involved in the surface correction deserves some comments. In principle, for each combination of evolutionary and pulsation code the power-law exponent b should be determined from a solar calibration. This would also define the input physics to be used when modeling the targets, including the mixing-length parameter. In practice, variations of b do not affect the obtained mean density of the models (Kjeldsen et al. 2008), and one unique value can be used to model all stars. The b parameters applied by the modeling teams in this study have been obtained either via a solar calibration of their own codes or taken from others. Therefore the adopted values vary among the results.

Another point where differences arise is in the calculation of the scaling parameter r. The best model for each star should have an r value very close to unity. In the original formulation, this is determined from theoretical and observed frequencies of radial modes only. Although it is the most commonly adopted manner, some modeling teams do not consider this r-scaling, or they estimate r using also the non-radial modes (models B in Section 5.1.2 above). Finally, the reference frequency ν₀ is arbitrarily chosen.

Since the surface correction is an empirical law based on the observed departure between the frequencies of the Sun and a solar model, there is no unique way of applying it. The net result is that χ² values for individual frequencies determined by each team are not directly comparable due to these different implementations. Thus, we have uniformly applied the surface correction to all results to homogeneously compare them. As suggested by Kjeldsen et al. (2008), we fixed b = 4.9, calculated r only with ℓ = 0 modes, and considered the observed frequency of maximum power ν_max from each star as the reference frequency ν₀.

5.3. Convective Core Location

In Section 4 we described how the presence of a convective core produces an abrupt change in the adiabatic sound speed. The sharpness and shape of that transition will depend on the mixing processes taking place and the way overshooting is treated, if considered at all. In Figure 4 some examples of these transition regions are shown for models with convective cores, where distinctive features can be observed in different cases.

Zoom In Zoom Out Reset image size

Model B in Figure 4 does not include overshooting and shows simply a sharp vertical transition at the edge of the core. Model G on the other hand applies a chemically homogeneous extension of the convective region when applying overshoot, and thus presents a similar structure beyond the formal limits of the convective core. In contrast model F shows a smoother profile and a more extended transition region, produced by the interplay between overshooting and microscopic diffusion at the edge of the core.

Since different assumptions and evolutionary codes produce different sound-speed profiles, it is necessary to determine one criterion defining the position of the convective core for all models regardless of the considered input physics. As the p-modes propagate according to the sound speed, they are most sensitive to the regions where rapid variations in it take place. Therefore, we use the maximum in the adiabatic sound-speed derivative as a common definition for the boundary of the homogeneously mixed central region. This basically defines the extent of the mixed zone, which is equivalent to the formal convective core size unless overshooting is included.

In the following sections, we compare the best-fitting models found by each modeling team as described in Section 5.1. It is worth mentioning that, when comparing radius determinations, we decided to use the radius measured at the photosphere in order to avoid systematics introduced by the different prescription for atmospheric stratification used in each evolutionary code.

In all cases, the values of $\chi ^2_\nu$ for frequencies and $\chi ^2_{r02}$ and $\chi ^2_{r010}$ for frequency ratios (Equations (5) and (6)) given in the results below are calculated using frequencies corrected as explained in Section 5.2. For comparison purposes we also give the results of the fit to the individual frequencies using the original frequencies as corrected by each modeling team, which we denote $\chi ^2_{\nu \, \mathrm{orig}}$. Since the frequency correction does not affect the computed ratios, we will sort models according to the average between $\chi ^2_{r010}$ and $\chi ^2_{r02}$. We remind the reader that the correlations in the ratios r₀₁₀ are taken into account when computing $\chi ^2_{r010}$, as mentioned in Section 5.

6.1. Perky

The main characteristics of the best-fit model found by each team for the target Perky are given in Table 4. The first feature to notice is that the agreement in surface gravity between all models is better than 0.02 dex. This level of precision is already better than what was found using the direct and grid-based methods (see Table 3), and the obtained mass range from detailed modeling is now restricted to values between 1.05–1.17 M_☉. This decrease is partly responsible for the reduction in the age difference compared to the grid results, now ranging from ∼3.8 to ∼5.5 Gyr. The age spread is constrained within a level of ∼18%, which although not extremely good is much better than that obtained from the global asteroseismic fit (see Section 3). Effective temperatures are also contained within the values given in Table 2. The overall agreement in radius is of the order of 0.05 R_☉, only slightly better than that found with the global analysis.

Table 4. Main Characteristics of the Models Calculated for Perky from Each Method as Described in Section 5.1

Model	M	Age	log  g	Teff	RPhot/R☉	[Fe/H]a	Zi	Yi	log  (L/L☉)
	(M☉)	(Gyr)		(K)
Model A	1.09	4.62	4.299	6074	1.2257	−0.09	0.0181	0.273	0.263
Model B	1.11	4.36	4.305	6122	1.2284	+0.00	0.0175	0.267	0.279
Model C	1.05	5.53	4.294	5970	1.2093	−0.09	0.0145	0.260	0.222
Model D	1.17	4.97	4.306	6016	1.2570	−0.01	0.0210	0.240	0.268
Model E	1.17	4.94	4.307	5951	1.2580	+0.01	0.0197	0.226	0.250
Model F	1.12	4.88	4.301	6100	1.2370	−0.09	0.0179	0.264	0.279
Model G	1.06	4.47	4.299	6071	1.2088	−0.10	0.0139	0.271	0.250
Model H	1.13	3.80	4.305	6134	1.2382	−0.07	0.0188	0.269	0.289

Note. ^aUsing the Grevesse & Noels (1993) solar mixture; see Equation (14).

Download table as:  ASCIITypeset image

Models D and E have initial helium abundances below the currently accepted SBBN value. This problem occurs regularly when models are selected solely on the base of best reproducing the observed individual frequencies after applying the surface correction (e.g., Metcalfe et al. 2010; Mathur et al. 2012). Interestingly enough, models D and E have basically the same characteristics. Both cases present the highest mass values among the results, hinting that the correlation between mass and helium abundance through the mean molecular weight is reflected in the asteroseismic fits (see also Figure 5 in Metcalfe et al. 2009). The initial metallicities of these models are also the highest ones, probably to compensate for the low helium abundance when achieving a similar luminosity and effective temperature. This might point toward deficits in the use of opacity or EOS data. We discuss this further in Section 7.

Zoom In Zoom Out Reset image size

As we aim to make a better characterization of the target properties we analyze the spectrum of oscillations in detail. In Table 5 the internal structure and asteroseismic properties are given for the best-fit models of Perky. The values of $\chi ^2_{r010}$ have been calculated using the frequency range marked in Figure 1. Interestingly, the models with the lowest $\chi ^2_{\nu }$ value for the individual frequencies are not the ones that best reproduce the observed ratios, and are in fact those with the lower than SBBN initial helium abundances (D and E). Model F presents a large difference between the original surface correction $\chi ^2_{\nu \, \mathrm{orig}}$ and the one we applied ($\chi ^2_{\nu }$). This model was the only one selected using a prescription different from the Kjeldsen et al. (2008) one (see Section 5.2), reflecting the large impact of the assumed recipe for surface correction and warning us about selecting models solely based on fits to individual frequencies.

Table 5. Structure and Seismic Characteristics of the Models for Perky

Model	Xc	$\tilde{t}_\mathrm{cc}$	$\tilde{t}_\mathrm{Phot}$	$\tilde{t}_\mathrm{Tot}$	$\tilde{t}_\mathrm{cc}/\tilde{t}_\mathrm{Phot}$	$\chi ^2_{\nu \, \mathrm{orig}}$	$\chi ^2_{\nu }$	$\chi ^2_{r010}$	$\chi ^2_{r02}$	Symbol
		(s)	(s)	(s)	(%)
Model A	0.124	...	4549.5	4577.4	...	...	21.28	0.79	0.35	Filled circle
Model B	0.249	95.0	4541.9	4621.7	2.09	9.15	8.96	1.07	0.43	Open square
Model C	0.078	...	4557.5	4593.7	...	11.03	15.40	1.52	1.12	Filled cross
Model D	0.085	...	4565.3	4691.6	...	2.37	3.25	1.67	1.21	Filled square
Model E	0.120	...	4573.0	4633.3	...	5.04	6.05	2.23	1.39	Downward triangle
Model F	0.085	...	4553.4	4725.7	...	3.50	10.27	1.04	3.62	Open diamond
Model G	0.122	...	4546.3	4626.1	...	...	21.78	2.22	3.56	Upward triangle
Model H	0.195	...	4547.6	4613.6	...	19.19	20.13	1.97	5.46	Filled diamond

Notes. Values of $\chi ^2_{\nu \, \mathrm{orig}}$ were obtained using the frequencies as corrected by each modeling team, while $\chi ^2_{\nu }$ was calculated using the standard version of the surface correction. Last column gives the symbol used to plot each model in all figures.

Download table as:  ASCIITypeset image

A more detailed picture can be made when looking at Figures 5 and 6, where we show the comparison to the data of models A, E, F, G, and H. From Figure 5 it can be seen that the ratios r₀₁₀ of model E (G) are systematically below (above) the observed ones, allowing us to discriminate against them. Ratios r₀₂ in Figure 6 have larger error bars and less data, but those of model G are also in disagreement with the observations. Following the same line of reasoning models F and H can also be discarded due to their large discrepancy with the observed values of r₀₂.

More interesting than discarding models is trying to understand the structural differences that produce this behavior. As discussed in Section 4, the ratios r₀₂ are a proxy for age and thus central hydrogen content of the star. Inspection of Table 5 does not reveal a correlation between X_c and the goodness of fit to r₀₂. In fact, models with the same central hydrogen content (D and F) present very different values of $\chi ^2_{r02}$, while model B gives an excellent fit to these ratios with more than twice the amount of hydrogen remaining in the core than any other model. Because of the systematic differences in the input physics and evolutionary codes, we no longer find a correlation between r₀₂ and X_c.

The case of model B presents an interesting challenge to our current understanding of asteroseismic diagnosis. With the exception of this one, all other results point toward a star with no convective core. Nonetheless, this model was selected using a modified version of the Kjeldsen et al. (2008) surface correction that includes also the non-radial modes (see Section 5.2), naturally resulting in a $\chi ^2_{\nu }$ value that is slightly different than what would be obtained using the radial modes only. Interestingly enough, the algorithm applied to select model B returns a different best-fit case when the "normal" surface correction is applied (using only ℓ = 0 modes). This new solution, which we call model B', gives values of $\chi ^2_{\nu \, \mathrm{orig}}=9.27$ and $\chi ^2_{\nu }=8.18$. When compared to the results of model B given in Table 5 of $\chi ^2_{\nu \, \mathrm{orig}}=9.15$ and $\chi ^2_{\nu }=8.96$ we see that model B' was not selected because it used the modified surface correction, while it would have been chosen if the standard one was applied instead.

Model B', has exactly the same initial composition as model B in Table 4, but with a mass of 1.10 M_☉ and no convective core present. Its age is 4.72 Gyr, T_eff = 6105 K, log  g = 4.300, and luminosity log  (L/L_☉) = 0.275, in good agreement with the results of model A and similar to model B. However its central hydrogen abundance of X_c = 0.03 is much lower than found in the other solutions, and gives fits to the frequency ratios of $\chi ^2_{r010}=3.55$ and $\chi ^2_{r02}=0.43$.

Figure 7 shows the comparison to the data of models A, B, and B'. Despite the very different central hydrogen contents (0.115, 0.249, and 0.03 for models A, B, and B', respectively), the ratios r₀₂ show no significant sensitivity to this parameter. Inspection of the ratios r₀₁₀ allows us to discriminate against model B', but it is clear that for this particular case our current diagnostic tools are not capable of discriminating against the presence of a convective core (i.e., between models A and B). We remind the reader that using the standard formulation of the surface correction from Equation (13) would have resulted in model B' as the best fit found by this team, a model clearly not compatible with the observed ratios r₀₁₀.

To explore this issue in more detail, we reproduced the results of model B with the evolutionary code and input physics of model A. For the same initial metallicity and mass it was possible to sustain a convective core beyond the pre-main sequence phase only by removing the geometrical restriction for small convective cores (see Section 5.1.1). The resulting model has an age of 4.71 Gyr at T_eff = 6148 K, log  g = 4.303, and luminosity log  (L/L_☉) = 0.288. Convective core size is slightly smaller than the one found in model B, placed at an acoustic radius of $\tilde{t}=81.2$ s that corresponds to 1.77% of the total photospheric acoustic radius. The fit to the frequency data is also comparably good as in model B, with $\chi ^2_{r010}=1.45$ and $\chi ^2_{r02}=0.32$.

The overall picture we are left with is the following. Based on the analysis described in the previous paragraphs, we can discard models E, F, G, and H and obtain from the rest of the results a mean and standard deviation pointing toward a star of R = 1.23 ± 0.02 R_☉, M = 1.11 ± 0.05 M_☉, and 4.87 ± 0.50 Gyr. This corresponds to uncertainties of 1.6%, 4.5%, and 10.4% in radius, mass, and age, respectively; the level of precision obtained in the stellar parameters by this method of asteroseismic analysis is extremely good, particularly in age. However, the experience of model B shows that slightly different models with equally good or even better fits to the data can still be lurking among our results, and have not been picked out due to our assumptions on the fitting algorithms. Although the mass, radius and age are not particularly affected by this issue, there could be other sources of systematic errors arising from unexplored input physics. In any case, we expect the true uncertainties to be only slightly larger given the heterogeneous set of codes and input physics used in this analysis.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We are left with the unsolved issue of whether or not a convective core exists and, in relation to that, the amount of hydrogen still remaining for central burning. This has a direct impact on the inferred age at the turn-off: while model A will remain in the main sequence for ∼1 Gyr, model B will spend another ∼2 Gyr in this evolutionary phase.

6.2. Dushera

For our second target of interest, Dushera, we performed a similar analysis to the one made for Perky. In Table 6 we give the main characteristics of the best-fit model found by each team, ranked again by their average fit to the frequency ratios. In this case no model gives a result with an initial helium abundance below the SBBN value. Agreement in log  g values is again better than 0.02 dex, while the radii of the results are contained within a range of 0.06 R_☉.

Table 6. Main Model Characteristics of the Best Fits to the Data of Dushera

Model	M	Age	log  g	Teff	RPhot/R☉	[Fe/H]a	Zi	Yi	OV
	(M☉)	(Gyr)		(K)
Model F	1.12	3.53	4.208	6318	1.3803	−0.15	0.0162	0.310	αov = 0.2
Model A	1.17	4.06	4.215	6305	1.3990	−0.12	0.0135	0.261	ξ = 0.03
Model D	1.26	3.09	4.223	6305	1.4378	−0.03	0.0202	0.266	αov = 0.015
Model E	1.19	3.44	4.217	6262	1.4083	−0.11	0.0147	0.263	αov = 0
Model C	1.20	3.34	4.218	6094	1.4119	+0.09	0.0210	0.280	αov = 0
Model G	1.14	4.28	4.208	6040	1.3926	+0.00	0.0173	0.278	αov = 0.1
Model B	1.20	3.67	4.216	6205	1.4157	+0.00	0.0175	0.267	αov = 0
Model H	1.19	4.00	4.214	6150	1.4119	−0.08	0.0188	0.269	ξ = 0.02

Notes. Models are sorted according to the value of their $\chi ^2_{r010}$ as given in Table 7. The last column gives the value of the overshooting efficiency parameter used. ^aUsing the Grevesse & Noels (1993) solar mixture, see Equation (14).

Download table as:  ASCIITypeset image

The mass range of the models is similar to the one obtained via the grid-based method, but with an age spread of the order of ∼15%. All analyses for this target returned models where a convective core exists, but there are differences in the size of the centrally mixed region that we aim to characterize using the frequency combinations. Table 7 presents the internal structure characteristics of each best-fit model, as well as the values of the fits to individual frequencies and ratios.

Table 7. Internal Structure and Seismic Characteristics for the Models of Dushera

Model	Xc	$\tilde{t}_\mathrm{cc}$	Mcc	Rcc	$\tilde{t}_\mathrm{Phot}$	$\tilde{t}_\mathrm{Tot}$	$\tilde{t}_\mathrm{cc}/\tilde{t}_\mathrm{Phot}$	$\chi ^2_{\nu \, \mathrm{orig}}$	$\chi ^2_{\nu }$	$\chi ^2_{r010}$	$\chi ^2_{r02}$
		(s)	(%)	(%)	(s)	(s)	(%)
Model F	0.236	129.5	8.32	7.10	5367.1	5577.5	2.41	7.09	11.81	0.92	0.93
Model A	0.262	131.4	7.94	7.08	5368.1	5401.2	2.45	...	20.3	1.44	0.47
Model D	0.256	125.1	5.89	6.57	5387.4	5546.1	2.32	2.04	1.90	1.60	0.53
Model E	0.213	110.1	4.77	5.79	5276.0	5347.2	2.09	4.78	4.67	1.74	0.72
Model C	0.236	116.0	5.01	6.09	5390.8	5433.9	2.15	9.93	10.31	2.06	0.64
Model G	0.214	126.2	6.94	6.63	5396.9	5498.1	2.34	...	9.16	0.81	2.91
Model B	0.173	117.3	5.28	5.99	5319.5	5417.1	2.21	4.72	3.87	1.58	2.30
Model H	0.243	122.0	6.66	6.67	5373.7	5430.9	2.27	20.98	16.93	1.92	4.10

Note. Values of $\chi ^2_{\nu \, \mathrm{orig}}$ were obtained using the frequencies as corrected by each modeling team, while $\chi ^2_{\nu }$ was calculated using the standard version of the surface correction.

Download table as:  ASCIITypeset image

Once again, inspection of Table 7 reveals that models best fitting the individual frequencies after surface correction (D and B) do not necessarily give the best match to the frequency combinations. Just as in the case of Perky, the data warn us to be cautious when relying on results solely based on fits to individual frequencies. Also, comparison of the values given in Table 7 shows no significant correlation between the central hydrogen content X_c and r₀₂ among these models built with different input physics.

We can further constrain the properties of the target by analyzing the fits to the ratios, plotted in Figure 8 for the three models presenting the largest average of $\chi ^2_{r010}$ and $\chi ^2_{r02}$. Model H is systematically below the observed ratios r₀₁₀, while model B shows a slope that is not consistent with the data. Moreover, when looking at the ratios r₀₂, all models are systematically lower than the observed values.

Zoom In Zoom Out Reset image size

In Figures 9 and 10 we plot the ratios r₀₁₀ and r₀₂ of the first five models in Table 6. The slopes of models C and E are clearly different from the observed one, with model C (E) showing lower (higher) values than the data in the low (high) frequency range of the considered ratios. Model D presents a similar behavior as model C.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Models A and F are the only two that consistently fit all the available data of Dushera, reproducing also the slope of the observed ratios r₀₁₀ and r₀₂. Considering only these results, our analysis points to a star with R = 1.39 ± 0.01 R_☉, M = 1.15 ± 0.04 M_☉, and 3.80 ± 0.37 Gyr corresponding to uncertainties of 1.0%, 3.1%, and 9.9% in radius, mass and age, respectively. As in the case of Perky we caution the reader that these values are to be taken as a reference, since systematics can arise from unexplored input physics and from models with slightly different properties that can also fit the data correctly but are overlooked by the fitting algorithms. We have shown that models fitting all asteroseismic combinations are restricted to a very small stellar parameter space and thus we expect the total uncertainties to be only slightly larger than those given above.

The analysis of Dushera not only reveals the presence of a convective core but it also puts constraints on its size. Our diagnostics are sensitive to the total extent of the mixed region, which in this case was reproduced by the inclusion of overshooting. It is clear from Table 7 and Figure 9 that only models with a mixed core comprising more than 2.35% of the total acoustic radius are able to reproduce the ratios. The overall results leave us with a homogeneously mixed region of ∼2.4% $\tilde{t}/\tilde{t}_\mathrm{Phot}$, and in any case extending beyond the formal convective boundary due to the inclusion of overshooting. Thus, we measure here for the first time mixing beyond the formal convective boundaries in the core of a Kepler field main-sequence star.

Tests were made on the impact of the core size by running model A with higher and lower overshooting efficiency values to produce different sizes of mixed zones. These models were not able properly to reproduce the ratios r₀₁₀ and show that, unlike the case of Perky, these ratios are mostly sensitive to the size of a convective core when this is already well developed during the main-sequence phase. We also performed tests using the model F set of evolutionary calculations on the impact of reducing the considered number of ratios r₀₁₀ by two, and found no significant difference in the obtained solutions.

Asteroseismic analysis aims to provide accurate stellar parameters of pulsating stars by theoretically reproducing their frequency spectrum. This has usually been done by matching the observed oscillation frequencies to theoretically calculated ones after applying an empirical correction to account for our poor modeling of stellar outer layers. Our results confirm that best-fit models found by this method do not always reproduce frequency combinations, and that different criteria for defining the best-fit model can lead to solutions with similar stellar parameters but very different interior structures. As originally stated by Kjeldsen et al. (2008), the surface correction only ensures an accurate reproduction of the stellar mean density by assuming that the frequency offsets are caused solely by the properties of the outer layers. Usage of non-radial modes and frequency combinations is necessary to match properly the deep interior and provide a reliable estimate of stellar age.

Of particular importance is the implementation of the surface correction itself, which has a direct impact in the goodness of fit via the reduced $\chi ^2_\nu$ value. The values of the corrections are much larger than the errors in the frequencies, making the definition of the best-fit model sensitive to the way these corrections are applied. Thus, the $\chi ^2_\nu$ criterion cannot be used objectively among best-fits found using different implementations of the surface correction. We have instead proposed to characterize the goodness of fit to asteroseismic data by using the frequency ratios r₀₁₀ and r₀₂. These combinations are thought to be mostly determined by the central conditions of the star and almost unaffected by the outer layers, providing a more robust approach than matching individual frequencies. The masses and radii obtained by both techniques agree within 2σ, but the resulting models can have very different structural properties. Thus, although the overall results are not hugely affected by the fitting techniques and the issue of the surface correction, using frequency ratios instead of just reproducing the individual frequencies of oscillation is a more reliable way to yield accurate stellar parameters, in particular age (as was also suggested by Miglio & Montalbán 2005).

The frequency ratios have been used as proxies for different internal properties of stars. For instance, a tight correlation between r₀₂ and central hydrogen content (age) holds for a given set of input physics, as shown by the modified C-D diagram in Figure 2. However, no correlation is found between these parameters among the set of models that fit the seismic data. The reason for this is that the equivalence between central hydrogen content and r₀₂ is only valid for a fixed set of input physics, which is naturally not the case in this study (Monteiro et al. 2002; White et al. 2011).

In a similar way the ratios r₀₁₀, which were thought to be sensitive to the existence and size of the convective core as well as a proxy of the central hydrogen abundance, were found to be affected by other parameters not yet fully identified. The latter sensitivity is particularly important in the case of Perky, where analysis of r₀₁₀ did not allow us to discriminate between models with or without a convective core. The stellar parameters obtained at this stage are not significantly affected by the presence of a convective core; thus we can still derive a reliable mass, radius and age for this star. Nevertheless the remaining central hydrogen content of the target is not well constrained, and therefore neither is the age it will have at the end of the main sequence.

More importantly, the experience of Perky is a warning flag to our fitting algorithms and the use of the surface correction. It is easy to imagine that if all teams used a different implementation of the surface correction, different solutions from the ones presented here would have been found. In fact, the restriction to the amount of overshooting in the case of small convective cores might be overestimated in cases like Perky. It has been shown for the solar case that restricting the amount of pre-main sequence overshooting is necessary to reproduce helioseismic data (Schlattl & Weiss 1999), but for stars slightly more massive than the Sun this needs to be further explored.

Related to fitting techniques, it remains fully to understand why in some cases low initial helium abundances are obtained when only individual frequencies are matched. Although the results for Perky suggest that mass is the major contributor to this feature, the overall picture is certainly more complex. While there is no reason to believe a universal value of ΔY/ΔZ applies for all stars, and small variations in helium are therefore expected, sub-SBBN values are clearly a challenge to stellar models. Such a behavior already appeared from the analysis of the main-sequence broadening in field K dwarfs (Casagrande et al. 2007), and simple analytical considerations suggested that model opacities could have been a likely culprit (Portinari et al. 2010).

Mathur et al. (2012) performed a uniform asteroseismic analysis of 22 Kepler targets using AMP, the same code and fitting algorithm of models E described in Section 5.1.5. Their final results show that 10 of their 22 best-fit models have a mean value of initial helium below that of SBBN (see their Table 5). The models in Mathur et al. (2012) are not directly comparable to ours since they are based on analysis of only one-month-long observations, but the overall outcome clearly points out to the need for further studies on this topic. These are currently underway (A. C. S. Jørgensen et al. 2013, in preparation). Understanding these discrepancies is important not only for the sake of stellar models, but also for the possibility of using these stars as tracers of the Galactic ΔY/ΔZ and possibly one day for multiple populations in stellar clusters.

Stars like Perky lie at the limiting mass where the onset of core convection is thought to occur. Although such targets have been used in open clusters as a test for the solar abundances (Vandenberg et al. 2007), it has been shown that the existence or absence of a convective core in stellar models heavily depends on our assumptions about the input physics. The upcoming releases of Kepler data will give us more frequencies with higher accuracy and we can certainly hope to solve this issue.

In the case of Dushera, we have made the first direct detection of a convective core in a Kepler main-sequence target. Use of fits to the ratios r₀₁₀ and r₀₂ allow us to discard outliers and determine the characteristics of the star with an excellent accuracy. The size of the central mixed region is estimated to be ∼2.4% of the photospheric acoustic radius, showing that mixing beyond the formal Schwarzschild convective boundary exists. In this study the extra mixing has been modeled using different overshooting prescriptions; in reality it could also be produced, for instance, by rotational mixing or a combination of more physical processes.

For both targets, masses and radii obtained by our asteroseismic analysis have a precision of ∼4% and ∼1.5%, respectively. Regarding stellar ages, we determine them with an encouraging ∼10% level of precision. This is probably optimistic due to our fitting algorithms potentially overlooking models that might better reproduce the data, and also since our models do not explore all possible sets of input physics. Thus, additional sources of uncertainties can exist that have not been explored yet. In fact, systematics arising from, for instance, the unknown input physics can easily account for an extra ∼7% in age at the end of the main-sequence phase (Valle et al. 2013). Although our analysis confirms that we are certainly doing better than the ∼20%–50% in age obtained with the usual isochrone fitting from a single stellar evolution database (e.g., Takeda et al. 2007; Soderblom 2010), asteroseismology has the potential to deliver even better ages once we fully understand our current diagnostic tools and pursue further theoretical studies to develop even better ones.

The authors thank the anonymous referee for useful suggestions and comments that improved the quality of the paper, and also Martin Asplund and Nicolas Grevesse for clarifying discussions about spectroscopic determination of stellar parameters. Funding for this Discovery mission is provided by NASA's Science Mission Directorate. The authors wish to thank the entire Kepler team, without whom these results would not be possible. We also thank all funding councils and agencies that have supported the activities of KASC Working Group 1. We are also grateful for support from the International Space Science Institute (ISSI). The authors acknowledge the KITP staff of UCSB for their warm hospitality during the research program "Asteroseismology in the Space Age." This KITP program was supported in part by the National Science Foundation of the United States under grant No. NSF PHY0551164. Funding for the Stellar Astrophysics Centre is provided by The Danish National Research Foundation (grant agreement No. DNRF106). The research is supported by the ASTERISK project (ASTERoseismic Investigations with SONG and Kepler) funded by the European Research Council (grant agreement No. 267864). V.S.A. received financial support from the Excellence cluster "Origin and Structure of the Universe" (Garching). S.B. acknowledges NSF grant AST-1105930. The National Center for Atmospheric Research (NCAR) is sponsored by the National Science Foundation. I.M.B., M.S.C., and S.G.S. acknowledge the support from the Fundacão para a Ciência e Tecnologia (Portugal) through the grants SFRH/BD/41213/2007, SFRH/BPD/47611/2008 and SFRH/BPD/84810/2012 and the ERC, under FP7/EC, through the project FP7-SPACE-2012-312844. T.S.M. acknowledges support from NASA grant NNX09AE59G. A.M.S is partially supported by the European Union International Reintegration Grant PIRG-GA-2009-247732, the MICINN grant AYA2011-24704, by the ESF EUROCORES Program EuroGENESIS (MICINN grant EUI2009-04170), by SGR grants of the Generalitat de Catalunya and by the EU-FEDER funds. W.J.C. acknowledges support from the UK Science and Technology Facilities Council (STFC). This research was carried out while O.L.C. was a Henri Poincaré Fellow at the Observatoire de la Côte d'Azur, partially funded by the Conseil Général des Alpes-Maritimes. R.A.G. acknowledges the funding from the European CommunityÕs Seventh Framework Programme (FP7/2007-2013) under grant agreement No. 269194 (IRSES/ASK).