THE RADIUS DISCREPANCY IN LOW-MASS STARS: SINGLE VERSUS BINARIES

Low-mass stars (i.e., stars of solar mass or lower) are fascinating objects. During their main sequence (MS) lifetime, these stars have a convectively unstable region located beneath the surface (and extending all the way to the center for M ⩽ 0.3 M_☉), where large scale magnetic fields are generated through the dynamo mechanism (see, e.g., Charbonneau 2013). This interior structure is reminiscent of that of the Sun, with which they share many properties: strong surface magnetic fields (0.01–1 kG; see Reiners 2012); the host of phenomena usually referred to as magnetic activity (e.g., starspots, faculae, flares, ...), a manifestation of a magnetically heated upper atmosphere (Schrijver & Zwaan 2008); and a non-trivial rotational evolution, driven by the braking torque applied to the stellar surface by magnetized winds (Schatzman 1962; Skumanich 1972; Kawaler 1988). Furthermore, M type stars (M ⩽ 0.6 M_☉) are the most abundant objects in the Galaxy, making up almost 75% of its stellar content, and are very promising candidates for the search of exoplanets located in the so-called habitable zone (Kasting et al. 1993; Kopparapu et al. 2013). A solid theoretical understanding of the fundamental parameters of low-mass stars, leading to precise characterization of these objects from the available observables (e.g., mass–magnitude, mass–radius (M–R) relations) is therefore very valuable.

A long-standing issue in the theory of low-mass stars is the disagreement between theoretically derived and observed global parameters. Discrepancies in radius and effective temperature have been known to exist for a very long time (e.g., Hoxie 1973; Lacy 1977; see also Torres et al. 2010 for a review).

The most precise observational constraints on the M–R relation of stars have been traditionally provided by the light curve analysis of detached eclipsing binaries (DEBs). In these systems, masses and radii of the components can be very accurately determined from the light curve analysis (within few percent) and free of additional assumptions. Until recently, only a handful of such very valuable systems hosting stars in the low mass range were known (Torres et al. 2010). However, the situation has improved significantly in the last few years (Feiden & Chaboyer 2012a); moreover, the number of DEBs of known parameters can be expected to further increase in the near future, as a by-product of planet-searching missions, such as Kepler.

In DEBs, inflated radii and cooler effective temperatures, in comparison with the theoretical models, are commonly reported (by 5%–15% and 3%–5%, respectively, see, e.g., Torres 2013 and references therein). The radius and effective temperature discrepancies roughly compensate each other to give the same luminosity of the model, suggesting a surface origin of the phenomenon.

Due to unavailability of theoretical models and/or reliable observational constraints on the chemical composition of the stars, fitting of DEBs has often been attempted in the past using only solar metallicity models. Taking into account the metallicity information, Feiden & Chaboyer (2012a) were able to bring the radius discrepancy down to more moderate levels (around or below 4% for their entire sample).

The radius inflation is usually explained as the consequence of an enhanced level of magnetic activity (see, e.g., López-Morales 2007; Chabrier et al. 2007). Stars in DEB systems can be kept in a regime of fast rotation through the spin–orbit synchronization induced by tidal interactions. The high rotation rate generates strong magnetic fields via dynamo action (Charbonneau 2013), which in turn leads to high star spot coverage of the surface and suppression of convection within the star (Gough & Tayler 1966). Both of these effects limit the flux of energy throughout the star, which expands in response. However, since magnetic activity is enhanced by fast rotation, close-in orbit DEBs, provided they are old enough to have synchronized their rotation periods with the orbital period, should display a larger discrepancy than long-period systems. This is not always the case, as both short-period, non-discrepant systems (e.g., KOI-126, P_orb = 1.77 days; Carter et al. 2011) and long-period, discrepant ones (e.g., Kepler 16, P_orb = 41.1 days; Doyle et al. 2011; Winn et al. 2011) are known.

The constraints on the fundamental parameters of stars coming from DEBs have recently been complemented by the direct interferometric determinations of stellar radii. This technique can now reach a precision on the radius comparable to that of light curve analysis in DEBs (Berger et al. 2006; Boyajian et al. 2012a, 2012b); the mass, however, cannot be obtained independently of other measurements (e.g., empirical mass–magnitude relations). Nevertheless, these measurements can provide very valuable information on whether radius and effective temperature deviations exist in single stars. Since no magnetic activity enhancement mechanism via tidal interaction with a companion can be expected in this case, field stars of solar age and approximate solar mass should show very moderate or no radius discrepancy. A 10% radius deviation with respect to the models of Baraffe et al. (1998) has been found by Berger et al. (2006); however, this result has been challenged by Demory et al. (2009).

We present a grid of models in the mass range 0.1–1.25 M_☉ calculated for various quantities of the metallicity and the mixing length parameter. One of our goals is to extend and update the low mass end of the Yale-Yonsei (Y²) isochrones (Yi et al. 2001; Demarque et al. 2004). Our present calculations rely on up-to-date input physics, capable of accurately reproducing the conditions of stars at the low mass end of the range (e.g., high density/low temperature equation of state (EOS), surface boundary conditions based on realistic model atmospheres).

This work is also intended as the starting point of a series of papers on related topics. One planned application is the study of the rotational evolution of solar-like and very low-mass stars and the physical mechanisms providing rotational coupling in stellar interiors, in particular (see, e.g., Spada et al. 2010, 2011). We have thus calculated rotation-related quantities, such as the moments of inertia of the radiative and convective zones of the stars, as well as the convective turnover timescale. The latter is a global measure of the efficiency of convection and has been often used as a means to connect stellar rotation and activity (Noyes et al. 1984; see also Barnes & Kim 2010; Barnes 2010).

Second, we will investigate the effect of variable convective efficiency in stellar one-dimensional models. In conventional stellar models, convection is described by means of the classical mixing length theory (MLT; Böhm-Vitense 1958). The so-called mixing length parameter α, the ratio between the mixing length itself and the local pressure scale height, is a measure of the efficiency of convection, and, being a free parameter in the MLT, is typically calibrated on the Sun. There is increasing evidence, however, that the use of a solar-calibrated mixing length is not warranted for other stars, and that α should be tied to stellar properties instead (Bonaca et al. 2012). Theoretical investigations of these effects, with the goal of formulating prescriptions that could be incorporated in one-dimensional models, are well under way (e.g., Tanner et al. 2013).

A reduction of convective efficiency has also been proposed as the main global effect on stellar structure of the presence of magnetic fields. This effect has been taken phenomenologically into account by, e.g., Chabrier et al. (2007), using an effective value of α, reduced in comparison with its solar-calibrated value. A more quantitative approach, where the reduction of α is connected with the vertical magnetic field intensity according to the formulation of Gough & Tayler (1966), has been adopted by MacDonald & Mullan (2012, 2013). Finally, the formalism developed by Lydon & Sofia (1995) in an early attempt at a self-consistent one-dimensional modeling of the Sun with magnetic fields was recently applied to the DEB system EF Aquarii by Feiden & Chaboyer (2012b). Further theoretical efforts along these lines are a natural continuation of our present investigation on the radius discrepancy problem.

Our grid will also be useful in those studies where knowledge of a precise mass–radius theoretical relationship is crucial. For example, in the characterization of exoplanet-hosting systems, the accuracy in the determination of the mass and the radius of the planets is contingent on that of the hosting star (e.g., Henry 2004).

In the present paper, we investigate the (theoretical versus observational) relationships among the fundamental parameters of low-mass stars. Our aim is to shed light on the nature of the physical mechanisms producing the radius and effective temperature discrepancies in low-mass stars, by disentangling the intrinsic disagreement from that caused by other effects (e.g., age and/or metallicity uncertainties), and by ascertaining whether the observational data on binaries and single stars (from DEBs and interferometry, respectively) can be reconciled within a consistent paradigm. To this end, we compare our models with the two largest samples currently available that contain low-mass stars and satisfy a high accuracy criterion (i.e., the compilation of DEBs by Feiden & Chaboyer 2012a and the interferometric measurements by Boyajian et al. 2012b, respectively).

2.1. Description of the Code

All of the models were calculated using the Yale Rotational stellar Evolution Code (YREC) in its non-rotational configuration. Compared with the standard version of YREC (see, e.g., Demarque et al. 2008), the code used in this work contains some major improvements. The most important ones concern the treatment of the EOS and the use of outer boundary conditions based on updated non-gray atmospheric models (see Section 2.1.2 below). These improvements are particularly significant when modeling low mass objects, i.e., for M ⩽ 0.6 M_☉. Our present calculations, therefore, supersede the Y² models in this range of mass.

2.1.1. Standard Input Physics in YREC

The details of the input physics are as follows. We used Ferguson et al. (2005) opacities at low temperatures, and the OPAL Rosseland opacities at high temperatures (Iglesias & Rogers 1996); the EOS is the OPAL 2005 EOS (Rogers & Nayfonov 2002). The energy generation rates are calculated according to the prescription of Bahcall & Pinsonneault (1992); diffusion of helium and heavy elements is taken into account, with the diffusion coefficients calculated according to the prescription of Thoul et al. (1994). Convection is described with the MLT (Böhm-Vitense 1958); the value used for the MLT parameter α (the mixing length scaled over the pressure scale height) is discussed below (see Section 2.2). In all calculations, we adopt the Grevesse & Sauval (1998) value of the solar metallicity, (Z/X)_☉ = 0.0230 and an α-enhancement [α/Fe] = 0. Convective core overshooting is not necessary in modeling stars in the mass range considered here, and therefore, it is not taken into account.

2.1.2. Equation of State and Surface Boundary Conditions for Low-Mass Models

The standard treatment of the atmospheric boundary conditions in YREC relies on the specification of a temperature–optical depth (T–τ) relationship (typically, the one of the classical Eddington gray atmosphere model or the empirical one from Krishna Swamy 1966). However, although this approximation is adequate for stars of solar mass or higher, it can lead to significant errors in the determination of the stellar radius and effective temperature for stars of mass ≲ 0.3 M_☉ (Chabrier & Baraffe 1997; see also Spada & Demarque 2012). Since the mass range we wish to investigate extends below this limit, we incorporated in the standard version of YREC, described so far, the improvements discussed in Sills et al. (2000): the use of the SCVH EOS (Saumon et al. 1995) in the low temperature–high density regime and the surface boundary conditions derived from PHOENIX model atmospheres (Hauschildt et al. 1999).

In the current YREC implementation, the OPAL and SCVH EOS are ramped at T = 5000 K. For the surface boundary conditions, the value of the gas pressure at the stellar photosphere is interpolated from a table compiled from the PHOENIX atmosphere models. To achieve overall consistency within the grid, we used the surface boundary conditions derived from PHOENIX throughout the entire mass range covered. The original input tables used by Sills et al. (2000), generated from NextGen models (Hauschildt et al. 1999), were updated to include the new BT-Settl models (Allard et al. 2011), available from F. Allard's Web site,⁵ to account for non-solar metallicities.

As the BT-Settl models are only available for either the Grevesse & Noels (1993) or the Asplund et al. (2009) chemical compositions, complete consistency cannot be achieved with the composition we adopted for the interior calculations (i.e., that of Grevesse & Sauval 1998). We have used the BT-Settl atmosphere models with Grevesse & Noels (1993) mixture, which is the closest to our choice for the interior. The interior and atmosphere models are matched according to the value of [Fe/H] (i.e., they will not have exactly the same value of Z/X or of Z). The relative difference between the Grevesse & Noels (1993) and the Grevesse & Sauval (1998) value of (Z/X)_☉ is however, quite small (about 6%). Moreover, the difference in Z will only have second order effects on the stellar models, as a change in the atmospheric Z will result in a change in the value of the photospheric pressure used in the outer boundary conditions, and this, in turn, will affect the stellar model.

A comparison of the evolutionary tracks in the Hertzsprung–Russell diagram (HRD), obtained using the surface boundary conditions from NextGen and from BT-Settl models, is shown in Figure 1. Although the effects on the isochrone shown are rather modest, the differences in the pre-main sequence (PMS) tracks are non-negligible, especially for the lowest masses.

Zoom In Zoom Out Reset image size

2.1.3. Standard Solar Model Calibration

Using the input physics described so far, we calculated a standard solar model, i.e., we calibrated the initial helium content Y₀ and the MLT parameter α in order for a 1 M_☉ model to match the solar radius and the solar luminosity at the age of 4.57 Gyr, having set the photospheric value of (Z/X)_☉ = 0.0230 (Grevesse & Sauval 1998). Our best solar model has the initial composition Y_{0, ☉} = 0.27783, Z_{0, ☉} = 0.018811, and α_☉ = 1.8750; other details are listed in Table 1.

Table 1. Standard Solar Model

Parameter	Adopteda	Model
Age (Gyr)	4.57	⋅⋅⋅
Mass (g)	1.9891 × 1033	⋅⋅⋅
R (cm)	6.9598 × 1010	log (R/R☉) = 5 × 10−8
L (erg s−1)	3.8418 × 1033	log (L/L☉) = 4 × 10−8
RBCZ/R☉	0.713	0.715
(Z/X)surf	0.0230	0.0230

Note. ^aSee, e.g., Basu & Antia (2008).

Download table as:  ASCIITypeset image

2.2. Parameters of the Grid

Our grid of models covers the mass range 0.10 M_☉ ⩽ M ⩽ 1.25 M_☉ with increments of 0.05 M_☉. This suffices to encompass the entire range of solar-like stars (i.e., stars with a subsurface convection zone and an inner radiative region), as well as stellar objects below the transition to the fully convective regime. Modeling the brown dwarf regime is beyond the scope of the present work.

The initial chemical composition of the models is assigned as follows. We have calculated models for five different values of the metallicity (scaled to the solar value): [Fe/H] = +0.3, 0.0, −0.5, −1.0, −1.5, where (Z/X)_* = (Z/X)_☉ 10^[Fe/H]. The helium mass fraction Y is scaled with the heavy elements mass fraction Z according to the relation Y = Y_p + (ΔY/ΔZ) Z, where Y_p = 0.25 is the primordial helium abundance (Cyburt et al. 2008).⁶ In all of the calculations, the enrichment parameter (ΔY/ΔZ) is equal to 1.48, based on our standard solar model calibration (see Section 2.1.3).

The initial composition of the models, given by

$$\begin{eqnarray*} &&X = \frac{1 - Y_p}{1 + \left(\frac{Z}{X} \right)_*\left[1 + \frac{\Delta Y}{\Delta Z} \right]} \; \ \ \ \ Z = \left(\frac{Z}{X}\right)_* \, X, \end{eqnarray*}$$

is listed in Table 2.

Table 2. Initial Compositions of the Models in the Grid

[Fe/H]	Y = 0.25 + 1.48 Z	X	Z
+0.3	0.29573	0.67336	0.03090
0.0	0.27415	0.70952	0.01631
−0.5	0.25793	0.73671	0.00535
−1.0	0.25253	0.74574	0.00171
−1.5	0.25080	0.74865	0.00054

Download table as:  ASCIITypeset image

The solar-calibrated value of the MLT parameter α is equal to α_☉ = 1.875 (see Section 2.1.3); however, in order to investigate its impact on the structure and on the evolutionary tracks, we have also calculated models with α = 0.50, 1.00, and 3.00. In general, reducing the value of α results in models with larger radii and cooler effective temperature; the sensitivity to α decreases approximately exponentially with decreasing mass. An effective value of the MLT parameter α has been used in the past by many authors to mimic the reduced convection efficiency induced by magnetic fields (e.g., Chabrier et al. 2007). A detailed analysis of this effect in the context of the radius discrepancy problem, discussed in Section 3, is deferred to future work.

The initial models have homogeneous composition and a polytropic structure. They are based on the birth line models constructed by Dr. Sydney Barnes at an early stage of the Y² project (see Yi et al. 2001). Each track is evolved through the PMS and MS, until either the bottom of the red giant branch or the age of 13 Gyr is reached.

The evolutionary tracks in the HRD for a representative subset of our grid ([Fe/H] = 0.0, α = α_☉) are shown in Figure 2.

Zoom In Zoom Out Reset image size

2.3. Description of the Tracks and Isochrones

2.4. Comparison with Other Models

In this section, we compare our models with those by the Lyon, Dartmouth, Y², and Padova groups (Baraffe et al. 1998; Dotter et al. 2008; Demarque et al. 2004; Bressan et al. 2012, respectively). The comparison between tracks (in our mass range, 0.1 ≲ M/M_☉ ≲ 1.25) and isochrones, all calculated for solar metallicity and solar-scaled MLT parameter, is shown in Figure 3.

• Lyon models (α = 1.900, Z = 0.0172). The models by Baraffe et al. (1998) cover the PMS phase, but do not go all the way through the subgiant branch. The largest difference with our models is found for the 0.1 M_☉ track; however, this mass is close to the lower limit of validity of our current configuration of YREC. These authors used the SCVH EOS and the NextGen version of the PHOENIX atmosphere models in their surface boundary conditions. The SCVH EOS is a pure H and He EOS, specifically designed to take into account non-ideal effects due to high density in low-mass stars; it has been shown by Chabrier & Baraffe (1997) to be in very good agreement with the MHD EOS (Hummer & Mihalas 1988; Mihalas et al. 1988; Daeppen et al. 1988) in the common range of validity (for a comparison of MHD and OPAL EOS, see Trampedach et al. 2006). Note that our version of YREC uses the SCVH EOS only in the low temperature/high density regime. The effects of the different versions of the PHOENIX atmospheric models (i.e., NextGen versus BT-Settl), on the other hand, have been discussed in Section 2.1.2.
• Dartmouth models (α = 1.938, Z = 0.0188). The Dartmouth Stellar Evolution Program (DSEP) is a distant relative of YREC. The two most important differences between the Dartmouth models and our calculations are in the choice of the EOS and the atmospheric boundary conditions. The latest version of DSEP (see Dotter et al. 2007, 2008 and references therein) uses NextGen PHOENIX atmospheres (see Section 2.1.2), and the FreeEOS.⁷ This set of models is probably the one in closest agreement with our calculations.
• Y² models (α = 1.743, Z = 0.0200). One of the aims of this work was to revise the Y² isochrones at low masses. The Y² models (Yi et al. 2001; Kim et al. 2002; Demarque et al. 2004) cover the mass range from 0.4 to 5.0 M_☉; convective core overshooting is taken into account. The models were constructed with a pre-Grevesse & Sauval (1998) solar mixture (i.e., with (Z/X)_☉ = 0.0245), and the OPAL 2001 EOS. Moreover, only helium diffusion is present in Y² models. However, the most significant difference with our current calculations is in the treatment of the atmosphere, which was based on a T–τ relation, unsuitable for M ⩽ 0.6 M_☉ (their lowest mass track has M ⩽ 0.4 M_☉). The large differences in Figure 3 for the tracks with lowest masses are therefore not unexpected, while a much better agreement can be observed moving toward the high mass end of the range.
• Padova models (α = 1.740, Z = 0.0190). In this comparison we use the latest version of the PARSEC models (Bressan et al. 2012). Their calculations make use of the FreeEOS; surface boundary conditions are based on the Eddington gray T − τ relation; convective core overshooting is taken into account. Both the tracks and the isochrone shown here are in good agreement with ours for masses larger than 0.6 M_☉. At lower masses, the use of gray atmospheres results in effective temperatures systematically larger than ours, as expected (see discussion in Section 2.1.2).

Zoom In Zoom Out Reset image size

In conclusion, the comparisons shown here demonstrate the sensitivity of the low mass models to the surface boundary conditions and, to a lesser extent, to the EOS. This sensitivity lessens moving toward higher masses. Our isochrones and tracks are thus in good agreement with those of Baraffe et al. (1998) and Dotter et al. (2008), who adopted similar choices for these pieces of input physics.

2.5. AI Phoenicis: a Stellar Evolution Benchmark

Binary systems with accurate determination of the metallicity and of the physical parameters of both components (masses, radii, effective temperatures) provide very valuable benchmarks for stellar evolution models. One of the best examples is the evolved system AI Phe (Andersen et al. 1988; see also Torres et al. 2010).

The metallicity of the system has been determined as slightly below solar ([Fe/H] = −0.14 ± 0.1). The masses of the two components, known with a relative accuracy of 0.4% (see Table 6), are sufficiently separate to give rise to a significant difference in the evolutionary timescales. As a consequence, a non-trivial test for stellar models is whether the best-fitting models consistently give the same age for both components.

Table 6. Measured Physical Parameters of the System AI Phe, from Andersen et al. (1988)

	A (hotter)	B (cooler)
M/M☉	1.1954 ± 0.0041	1.2357 ± 0.0045
R/R☉	1.816 ± 0.024	2.930 ± 0.048
Teff (K)	6310 ± 150	5010 ± 120
log L/L☉	0.67 ± 0.04	0.69 ± 0.04
[Fe/H]	−0.14 ± 0.1
Distance (pc)	173 ± 11

Download table as:  ASCIITypeset image

Figure 4 shows how our models perform in this test. The tracks shown in the figure were calculated using the observed metallicity and masses (listed in Table 6), our solar-calibrated α parameter, and the initial helium content fixed by our enrichment parameter ΔY/ΔZ = 1.48. As in all of the models discussed in the rest of the paper, core overshooting is not taken into account.

Zoom In Zoom Out Reset image size

We estimate the age of each component as the age of the model whose radius coincides with its observed value. We obtain ages of 4.44 Gyr and 4.54 Gyr, respectively, for the A and B component, compatible with each other at the 2% level. This gives an age of the system in between those obtained by Torres et al. (2010) using the Y² isochrones (Demarque et al. 2004) and the Victoria models (VandenBerg et al. 2006).

2.6. Empirical Mass–Magnitude Relations

Empirical mass–luminosity relations (MLR) in infrared and optical wavelengths have been established by Henry & McCarthy (1993) using data from speckle interferometry of binaries and from the analysis of eclipsing binary systems (see also Henry 2004; Henry et al. 2006). The MLRs were obtained through a fitting procedure which does not take into account metallicity and age differences. The data and the empirical MLRs are compared with the synthetic mass–magnitude relations derived from our models in Figure 5. The figure shows 5 Gyr isochrones with [Fe/H] = 0.0, [Fe/H] = −0.5, [Fe/H] = +0.3, and α = 1.875; the M_V and M_K absolute magnitude were obtained using the Lejeune et al. (1998) conversions. Notably, the MLR in the K band is much less sensitive to metallicity compared to the one in the optical band (see also Baraffe et al. 1998).

Zoom In Zoom Out Reset image size

The K-band MLR is particularly important in the following, since it has been used by Boyajian et al. (2012b) to estimate the mass for the stars in their sample. They assigned a fixed 10% error to the masses determined in this way; this is the largest source of uncertainty in the interferometric sample.

In the following, we discuss in detail the theoretical relations among global stellar properties derived from our models and their agreement with the observational data.

3.1. The Data

3.2. The M–R Relation

A comparison between our theoretical M–R relationship and both data samples is shown in Figure 6. The isochrones plotted have [Fe/H] = −0.5, 0.0, 0.3, and an age of 5 Gyr. Note that the theoretical M–R relation is insensitive to metallicity for M ≲ 0.7 M_☉. A larger deviation from the theoretical isochrones is immediately apparent for the short-period DEBs (plotted as open squares in the figure). In the interferometric sample, the largest disagreement is for M ≈ 0.4 M_☉.

Zoom In Zoom Out Reset image size

A quantitative measure of the discrepancy between observed radii and the model predictions can be defined as follows:

$$\begin{equation*} \frac{\delta R}{R_{\rm obs}} = \frac{R_{\rm obs} - R_{\rm mod}}{R_{\rm obs}}, \end{equation*}$$

where R_obs is the observed radius and R_mod is its theoretical counterpart, derived from the models.

3.2.1. The M–R Relation, DEB Sample

The tightest observational constraint on the M–R relationship is provided by the DEBs, due to the high accuracy in the determination of both the mass and the radius. The discrepancy δR/R_obs for the DEB sample is plotted in Figure 7.

Zoom In Zoom Out Reset image size

In all three panels of the figure, the value of R_mod for each star in the DEBs sample has been determined by interpolating our theoretical M–R isochrones with the observed mass of the star.

In the upper panel, we used theoretical isochrones of the same age and metallicity reported to be the best-fitting by Feiden & Chaboyer (2012a). We obtain an average δR/R of 2.4%, consistent with the findings of these authors. No correlation seems to be present between δR and the metallicity of the star (represented by the color of the symbols, see the caption of the figure). Interestingly, systems with orbital periods shorter than 1.5 days (open symbols) are significantly more discrepant than the others, as can be seen in the middle panel of Figure 7, where δR/R_obs is plotted as a function of P_orb (see also the average δR/R_obs for these two subsamples, shown in the legend). This was already found by Kraus et al. (2011), although these authors caution that this effect could be a spurious result of the light curve analysis.

The availability of model-independent constraints on either the age or the metallicity of DEBs systems is very scarce, leaving these two parameters essentially freely adjustable (with a few exceptions, e.g., KOI-126; Carter et al. 2011). To investigate the impact of the assumed age and metallicity on the radius discrepancy, we have recalculated the δR using values of the theoretical radii R_mod interpolated from a fixed solar metallicity, fixed 5 Gyr M–R isochrone (lower panel of Figure 7). The resulting average δR/R_obs is about 3%, still comparable with the value found by Feiden & Chaboyer (2012a). Varying the age within the range 1–9 Gyr, on the other hand, leads to "age error bars" (shown in blue in the figure) of the same size of or smaller than the observational uncertainty for all but three stars in the sample (which have M ≳ 0.9 M_☉), thus not significantly affecting the overall conclusions. Note that the correlation between discrepancy and short orbital periods is robust to these different choices of age and metallicity (see the legend in the lower panel of Figure 7).

3.2.2. The M–R Relation, Inteferometric Sample

A similar analysis was performed for the interferometric sample. Values of R_mod for each star have been calculated by interpolating in mass from isochrones of 5 ± 4 Gyr (the same for the whole sample) and metallicities equal to those quoted by Boyajian et al. (2012b). The results are shown in Figure 8; blue bars represent the uncertainty arising from the age range, while the colors of the symbols encode the metallicity information, as before. It should be emphasized that the comparison with the results on the M–R relationship for the DEBs sample, discussed previously, is hindered by the much larger errors on R_mod, which are in turn due to the larger errors in the masses derived from a mass–magnitude relation.

Zoom In Zoom Out Reset image size

For the interferometric sample, the average δR/R_mod ≃ 3% is compatible with that of the DEBs sample. Unfortunately, for the majority of the stars, this value is comparable with the uncertainty in R_mod. A notable exception are the eight stars of M ≲ 0.42 M_☉ in the sample, which appear to have inflated radii to the level of several sigma (see also Boyajian et al. 2012b). Note that five of these eight stars are reported to have sub-solar metallicity. We checked whether an erroneous metallicity determination could explain the anomalous behavior of these stars; using a solar value of Z for the theoretical isochrone, their discrepancy is lessened, but it remains significant.

Our analysis so far suggests that the empirical M–R relations of both the DEBs and the interferometric samples are, on average, in agreement with their theoretical counterpart within 2%–3%. Higher values of the discrepancy exist, but they are more the exception than the rule. This conclusion appears to be rather insensitive to the choice of the (poorly constrained) age and metallicity. However, the eight lowest mass stars in the interferometric sample challenge this claim, with an average discrepancy of the order of 8.5%. This result is robust to changes in the age or the metallicity of these stars. Notably, no such trend seems to be present at the low mass end of the DEBs sample.

3.3. The M–T_eff Relation

Radius deviations are usually observed to be associated with a concurrent effective temperature deviation, in such proportions to keep the luminosity constant (e.g., Torres et al. 2010). It is therefore interesting to compare the M–R relation with the mass–T_eff relation. Effective temperature measurements for DEBs typically have large errors, of the order of 100 K. Nevertheless, we compare the theoretical and empirical M–R and M–T_eff relations, for both the DEBs and the interferometric samples, in Figure 9.

Zoom In Zoom Out Reset image size

As we have already seen, the M–R relation of the two samples are, in general, compatible with each other. The largest deviations from our theoretical isochrones are observed in the mass ranges 0.6–0.8 M_☉ (both samples) and below about 0.4 M_☉ (most prominently in the interferometric sample). Figure 9 shows, on the other hand, that the DEBs have significantly cooler temperatures than the models for masses ≲ 0.55 M_☉, while such a trend is much less evident among the single stars.

It is worth noting that Boyajian et al. (2012b) reported a systematic difference between the effective temperature of stars in their sample and those from DEBs. The same authors have also stressed that the DEBs measurements could be affected by systematic errors and careful validation of these data is required before any conclusion can be drawn from them.

3.4. The L–R–T_eff Relations

The high precision in the determination of radius, luminosity, and effective temperature for the interferometric sample allows us to test the theoretical L–R–T_eff relations as well.

Figure 10 shows the comparison between our 5 Gyr isochrones (calculated with [Fe/H] = 0.0, − 0.5, and +0.3 and solar-calibrated α) and the interferometric data in the L–R and L–T_eff planes. Note that the models display a prominent metallicity dependence which was absent in the M–R relation (Figure 6). Remarkably, even if the range of the measured [Fe/H] is comparable to that of the theoretical isochrones plotted, the data are less scattered, in both planes.

Zoom In Zoom Out Reset image size

As was done in the previous section, we define the radius and effective temperature deviations as:

$$\begin{equation*} \frac{\delta R}{R_{\rm obs}} = \frac{R_{\rm obs} - R_{\rm mod}}{R_{\rm obs}}\ \; \ \ \ \ \frac{\delta T}{T_{\rm obs}} = \frac{T_{\rm obs} - T_{\rm mod}}{T_{\rm obs}}, \end{equation*}$$

where R_mod and T_mod are calculated interpolating in luminosity (with the observed value L_obs) the L–R and L–T_eff theoretical relations, respectively.

For each star in the sample, 5 Gyr isochrones with the observed metallicity were used in the interpolation; the results are shown in Figure 11. A significant trend emerges of the observed radii being inflated (by about 7%, on average) and the effective temperatures being cooler (by about 4%) than predicted by theoretical models. The radius and effective temperature deviations are significant for R ≲ 0.7 R_☉ and T_eff ≲ 5000 K. The R and T_eff discrepancies are not correlated with the stellar activity (as measured by the X luminosity proxy L_X/L_bol). These results are consistent with those of Boyajian et al. (2012b).

Zoom In Zoom Out Reset image size

Clearly, for the same interferometric sample, we obtain quite a different picture for the radius discrepancy depending on whether the theoretical radius R_mod is calculated by interpolating in mass or in luminosity. In particular, the average discrepancy <δR/R > is more than a factor of two larger in the latter case. This is, however, not surprising, since the theoretical M–R relation and MLR have very different properties from a stellar evolution standpoint. First of all, the L–R and L–T_eff relations have a much stronger dependence on metallicity than the M–R relation (see Figure 10 and the lower panel of Figure 6, respectively). Moreover, the uncertainty in age plays a very different role for large enough masses, as illustrated in Figure 12. In the mass–radius plane, a star moves along a vertical track as it ages; this effect was represented by the blue "age bars" in Figures 7 and 8. In the luminosity–radius plane, on the contrary, stars move along the isochrone as they age. This can lead to compensation effects between the age and the mechanism responsible for the radius discrepancy (whatever it may be), which cannot be disentangled without independent information on the stellar age. The age shift along the isochrone starts to be significant above 0.7 M_☉, and it is very likely to be the cause of the threshold for non-significant δR and δT_eff, shown in Figure 11.

Zoom In Zoom Out Reset image size

As has been suggested by Boyajian et al. (2012b), the R and T_eff deviations themselves could be explained if the effect of metallicity on the stellar parameters is exaggerated by theoretical models. This would lead to a larger discrepancy when interpolating in luminosity than in mass, while the residual, more modest, discrepancy in the M–R relation could be due to an entirely different mechanism.

3.5. Effect of the MLT Parameter Variation

The MLT parameter, α, is a measure of the efficiency of convection. In models of solar-like stars, the global radius is quite sensitive to the value of α, its reduction resulting in larger radii (see Christensen-Dalsgaard 1997). In the case of fully convective stars, the impact of α on the radius is much more modest (Baraffe et al. 1998). Although the standard practice in stellar modeling is to calibrate α on the Sun, lower effective values of this parameter have been used in the past as a means to phenomenologically account for the additional physics (e.g., magnetic fields) responsible for the radius discrepancy (e.g., Chabrier et al. 2007).

We compare the mass–radius and mass–T_eff isochrones for various values of α with the data from the DEB and the interferometric samples in Figure 13. As α is freely adjusted in this comparison, we restrict ourselves to qualitative remarks, without calculating the δR and δT_eff. Figure 13 is consistent with the standard results discussed above (larger radii for reduced-α models of solar-like stars, very weak sensitivity for fully convective models). It also shows that the data are compatible with moderately reduced values of the MLT parameter, i.e., 1 ≲ α ≲ 1.875, while α = 0.5 seems to be too extreme.

Zoom In Zoom Out Reset image size

We have calculated models of low-mass stars in the range 0.1–1.25 M_☉ for various choices of the composition and of the MLT parameter α. In our calculations, we have used an improved version of the YREC stellar code, optimized to provide a realistic treatment of the EOS and atmospheric boundary conditions for low-mass stars (M ≲ 0.6 M_☉). The present models complement and update the low mass end of the Y² isochrones (Yi et al. 2001; Demarque et al. 2004).

We have compared the theoretical relations among fundamental stellar parameters predicted by our models with the best currently available data from DEBs and single stars interferometric measurements. We find that radius deviations are present, on the average, at the ≈3% level in both samples when the theoretical isochrones are interpolated in mass, consistent with the results of Feiden & Chaboyer (2012a). Among the binaries, the most discrepant stars are those members of short orbital period systems. Among the single stars, very large deviations (10%–20%) are found for stars of mass ≲ 0.4 M_☉. This very discrepant subgroup has no counterpart in the DEB sample; notably, it is composed of predominantly metal-poor stars. No correlation is observed with the activity indicator L_X/L_bol. The values adopted for the metallicity and the age of the stars do not significantly affect these results, since the theoretical M–R relation is almost insensitive to metallicity.

These results are consistent with the commonly accepted explanation of the R and T_eff discrepancies as manifestations of enhanced magnetic activity (see, e.g., Chabrier et al. 2007). If a binary system is old enough (e.g., age ≳ 1 Gyr for P_orb ⩽ 5 days, see Figure 9 of Meibom & Mathieu 2005), its components will have synchronized their rotation periods with the orbital period via tidal interaction. Close-in DEB systems are thus expected to be more magnetically active due to the faster rotation rate powering their dynamo (Charbonneau 2013). The interferometric stars, on the other hand, are single and likely old enough (≈5 Gyr) that their level of activity has already declined (e.g., Wright et al. 2011). Theoretical efforts are under way to incorporate the effects of magnetic fields in one-dimensional stellar models with a self-consistent approach (based on the prescription of Lydon & Sofia 1995; see also Feiden & Chaboyer 2012b). We plan to tackle this intriguing problem in a forthcoming paper.

The very precise determinations of luminosity and effective temperature, along with radius, available for the interferometric sample, allow us to perform a complementary analysis on these stars, to test the L–R–T_eff relations. In this analysis, the luminosity is held constant, as is (approximately) expected to be the case even for discrepant stars, where the radius and effective temperature deviations should compensate each other (Torres et al. 2010). We obtain much larger discrepancies in this case (of the order of 7 and 4%, on the average, for radii and effective temperatures, respectively). A trend of vanishing deviations above 0.7 M_☉ is also found. However, rather than attribute the deviation as a feature of the underlying physical mechanism causing the radius discrepancy, we interpret the lack of observed deviation for higher mass stars as the consequence of compensation effects due to the age scatter, since evolutionary tracks are almost parallel to the theoretical isochrones in the L–R plane. For stars of mass below the threshold of 0.7 M_☉, the larger discrepancies suggest that the effect of metallicity on the global stellar parameters is overestimated by the theoretical models in comparison with what is observed in the data (see also Boyajian et al. 2012b).

In DEBs, the existence of the radius discrepancy phenomenon seems to be established beyond a doubt and its extent and physical origin to be understood at least at the most basic level. A major source of concern for this picture, however, is the possible presence of unaccounted systematic errors. Morales et al. (2010) have suggested that systematic errors in the radii determined from the light curve analysis of DEBs could arise due to the variable star spots distribution and size. Systematic errors of the order of 3% would be sufficient to completely mask the signal (if any) of the genuine radius discrepancy (Feiden & Chaboyer 2012a).

Stellar radii measured by interferometry, with precision comparable to DEBs, are very valuable to ascertain whether the radius discrepancy affects single and binary stars in the same way. Currently, the main differences between these two groups are in the M–T_eff relations and in the mild metallicity dependence at the low mass end of the M–R relation of single stars, which is absent in binaries. However, both effective temperatures of DEBs and metallicity determinations (in both samples, especially for M dwarfs) are notoriously difficult measurements and are plagued by large errors (≈100 K and ≈0.2 dex, respectively). More observational evidence needs to be gathered before attempting to interpret these differences.

F.S. acknowledges support from the Leibniz Institute for Astrophysics Potsdam (AIP) through the Karl Schwarzschild Postdoctoral Fellowship. F.S. is grateful to Dr. Sydney Barnes for reading the manuscript and providing useful comments, and to Professor Alessandro Bressan for sending the latest version of the evolutionary tracks of the Padova group. P.D. thanks Dr. Tabetha Boyajian for interesting discussion.

One of our aims in this work is to provide a grid of models suitable for studies of rotational evolution of low-mass stars, their magnetic activity, and the properties of their convective envelopes, which are believed to be the seat of the stellar dynamo (see, e.g., Charbonneau 2013). We thus include in our evolutionary tracks rotation-related parameters, such as moments of inertia and global and local convective turnover timescales.

For illustrative purposes, we show in Figure 14 the moments of inertia of the convective envelope and of the radiative zone (if present) for our solar metallicity models.

Zoom In Zoom Out Reset image size

The convective turnover timescale, τ_c, represents a characteristic time for the rise of a convective element through the stellar convection zone. Following Kim & Demarque (1996), we have calculated a "global" and a "local" estimate of τ_conv:

$$\begin{eqnarray*} \tau _{\rm c,global} &=& \int _{R_{\rm BCZ}}^{R_*} \frac{dr}{v_{\rm conv}}; \\ \tau _{\rm c,local} &=& \frac{\ell ^*}{v_{\rm conv}^*}; \end{eqnarray*}$$

the local definition refers to a distance of one half of the mixing length above the bottom of the convective envelope, i.e., ℓ* = αH_P/2 and $v_{\rm conv}^* = v_{\rm conv}(R_{\rm BCZ}+\ell ^*)$.

Even though the MLT is an admittedly crude approximation of real convection in the atmosphere (e.g., Kim et al. 1996; Ludwig et al. 1999; Tanner et al. 2013), it is a remarkably reliable approximation of stellar convection in the deep layers (Chan & Sofia 1989). The parameter τ_conv has the advantage of providing a simple parameterization of how the properties of convection scale with, e.g., stellar mass. An example of that is the Rossby number, a dimensionless parameter calculated from τ_conv and the surface rotation period:

$$\begin{equation*} {\rm Ro} = \frac{P_{\rm rot}}{\tau _{\rm c}}. \end{equation*}$$

Many studies on stellar activity have shown that the Rossby number is a critical parameter, marking the transition between different regimes in the chromospheric and X emission (Noyes et al. 1984; Wright et al. 2011). Although our estimates of τ_c are based on non-rotating models, this should not significantly affect the results when calculating Ro, as the effects of rotation on the structure of low-mass stars are quite moderate if the rotation rate is not extreme (i.e., far from break-up speed; see Sills et al. 2000). The time evolution and mass dependence of τ_c for our solar composition subgrid are shown in Figure 15. Note that the value of τ_c for a star of given mass is approximately constant for the whole MS lifetime.

Zoom In Zoom Out Reset image size

Theoretical calculations of the convective turnover timescale have been reported by Kim & Demarque (1996), Landin et al. (2010), and Barnes & Kim (2010). The time evolution of τ_c calculated by these authors agrees quite well with our results (compare the left panel of Figure 15 with, e.g., Figure 3 of Kim & Demarque 1996 or Figure 3 of Landin et al. 2010). When comparing the absolute values of τ_c at a given age, however, a ≈10% scatter emerges. Other choices of the input physics being equal, a strong dependence of τ_c on the value of the MLT parameter α can of course be expected. To illustrate the influence of other parameters, we compare, in Figure 16, the various sets of τ_c for 0.6, 0.7, ..., 1.2 M_☉ at 0.5 Gyr, scaling them over the value corresponding to the 1 M_☉ model. Once a constant (i.e., mass-independent) scale factor is accounted for, our calculations are in remarkably close agreement with those of Barnes & Kim (2010). For the other two sets of τ_c, the scale factor depends on mass. This is most likely the consequence of different choices in the input physics to which stellar models have different sensitivity according to their mass, e.g., the treatment of the atmospheric boundary conditions. It should also be noted that both the Kim & Demarque (1996) and the Landin et al. (2010) models take stellar rotation into account, which more strongly affects massive models than less massive ones (Sills et al. 2000).

Zoom In Zoom Out Reset image size

To construct synthetic color–magnitude diagrams (CMDs), we need to transform from the theoretical variables (log g, T_eff) to colors. Various such transformations are reported in the literature, either constructed purely from theoretical atmosphere models or established on the basis of semi-empirical calibrations. We have compared the following color transformations.

• 1.  
  Synthetic transformations derived from PHOENIX atmosphere models, available in various photometric systems from F. Allard's Web site.⁸ This is the natural choice for our models, ensuring consistency (whose importance has been stressed by Baraffe et al. 1998) with the surface boundary conditions used (see Section 2.1.2).
• 2.  
  Semi-empirical transformations in the (UBV)_J(RI)_CJHKLL'M by Lejeune et al. (1998), extending the grid produced by the same authors (Lejeune et al. 1997) to include M dwarfs. Starting from a hybrid grid of theoretical stellar spectra from many literature sources, they developed a semi-empirical calibration of the spectral continua; color–T_eff relations were then calculated from the corrected spectra.
• 3.  
  Semi-empirical transformation in the (BV)_J(RI)_C system by VandenBerg & Clem (2003), obtained from a grid of theoretical stellar spectra, applying empirical corrections based on constraints from clusters and field stars.

All of the color–T_eff relations above come in tabular form as a function of (log g, T_eff) for various metallicities; bolometric corrections and color indices are obtained through interpolation in these variables.

As an example, we compare our theoretical isochrones with data for the Pleiades in (B − V) and (V − I_C) CMDs in Figure 17. The data and the empirical single star sequence are from the compilation by Stauffer et al. (2007); the BVI_C data are from various literature sources, among which the authors have favored photoelectric photometry when possible (see references in the paper). We adopt the same fundamental parameters for the cluster used by these authors, i.e., solar metallicity, a distance modulus m − M = 5.62 (corresponding to the ground-based distance estimate, i.e., 133 pc), an age of 100 Myr, and the following values of the reddening: A_B = 0.16, A_V = 0.12, A_I = 0.07 (see Stauffer et al. 2007).

Zoom In Zoom Out Reset image size

The theoretical isochrones shown in the figure were constructed from our solar metallicity models by applying the PHOENIX, Lejeune et al. (1998), and VandenBerg & Clem (2003) color transformations. All of the isochrones agree very well with the data for M_V ≲ 8, which correspond to M ≳ 0.6 M_☉. At lower masses, however, the purely theoretical PHOENIX transformations produce the largest discrepancy in both CMDs. In the (B − V) CMD, both semi-empirical transformations perform equally well in providing a good match of the data over the whole mass range considered here (i.e., M ≳ 0.30 M_☉). In the (V − I) CMD, on the other hand, all the relations fail to reproduce satisfactorily the fainter magnitude end (M_V ≳ 8.5). The VandenBerg & Clem (2003) transformation results in a reasonable fit up to M_V ≲ 10.

It is worth noting that Lejeune et al. (1998) caution that their empirical colors between 4500 and 2000 K have higher uncertainties than the ones at higher temperatures, due to the difficulties of defining a continuum for these stars. These results are also compatible with the discussion in Section 3.4 of VandenBerg & Clem (2003), who argued that the failure of the Baraffe et al. (1998) model to match the observations at the faint end of the (V − I) CMD is not due to a flaw intrinsic to the models, but to the color–T_eff relations used (for the M_V magnitude in particular), inaccurate for low-mass stars. A similar result has been reported by Stauffer et al. (2007) about the models of both Baraffe et al. (1998) and Siess et al. (2000).

The disagreement between models and observations for T_eff ≲ 4000 K is usually ascribed to the existence of sources of opacity unaccounted for in atmosphere and/or interior models. Moreover, for a young cluster as the Pleiades, and up to about 500 Myr, stars in this range of T_eff are still in the PMS/early zero-age main sequence phase. The missing opacities problem, combined with the notorious uncertainties of PMS models, can also lead to a disagreement between cluster ages estimated on the basis of MS and PMS isochrones (Bell et al. 2012). A detailed analysis of PMS evolution is deferred to future work.