Evolutionary models for ultracool dwarfs

TL;DR

This paper extends an evolutionary code to model ultracool dwarfs, incorporating detailed physics and atmospheric models, and validates the models against observed data, providing a new tool for characterizing these objects.

Contribution

The authors developed and validated ultracool dwarf evolutionary models with detailed physics, including C and O effects, and provided reference tables for parameter estimation.

Findings

Including C and O significantly affects models near the H-burning limit.

Models reproduce observed masses within 1-2.5 sigma.

Systematic errors in models are around 0.0005 solar masses.

Abstract

Ultracool dwarfs have emerged as key targets for searches of transiting exoplanets. Precise estimates of the host parameters (including mass, age, and radius) are fundamental to constrain the physical properties of orbiting exoplanets. We have extended our evolutionary code CLES (Code Li\'egeois d'Evolution Stellaire) to the ultracool dwarf regime. We include relevant equations of state for H, He, as well as C and O elements to cover the temperature-density regime of ultracool dwarf interiors. For various metallicities, we couple the interior models to two sets of model atmospheres as surface boundary conditions. We show that including C and O in the EOS has a significant effect close the H-burning limit mass. The typical systematic error associated with uncertainties in input physics in evolutionary models is $\sim 0.0005 M_\odot$. We test model results against observations for objects whose parameters have been determined from independent techniques. We are able to reproduce dynamical mass measurements of LSPM J1314+1320AB within $1\sigma$ with the condition of varying the metallicity (determined from calibrations) up to $2.5\sigma$. For GJ 65AB, a $2\sigma$ agreement is obtained between individual masses from differential astrometry and those from evolutionary models. We provide tables of ultracool dwarf models for various masses and metallicities that can be used as reference when estimating parameters for ultracool objects.

Figures

Tables

Table 1. Table 1: Stellar evolutionary models for LSPM J1314+1320AB

	CLES	BHAC15
$T_{\mathrm{eff}}$ (K)	$2950\pm 6$	$2950\pm 4$
age (Myr)	$81.7\pm 3.6$	$80.8\pm 2.5$
log $g$	$4.840\pm 0.013$	$4.839\pm 0.009$
$R/R_{\odot }$	$0.1868\pm 0.0021$	$0.1871\pm 0.0016$

Full text

stars: low-mass – stars: late-type

††journal: ApJ††software: CLES (Scuflaire et al., 2008)

1 Introduction

Ultracool dwarfs (UCDs) are cool, small and dim objects lying at the faint, red end of the Hertzsprung-Russell (HR) diagram (M7 and cooler; Kirkpatrick et al. 1995), at the limit of core H-burning. UCD classifications encompasses the very low-masses stars (VLMS) and brown dwarfs (BDs). As these objects emit primarly in the near-infrared (NIR) and have low intrinsic luminosities, it took the development of efficient photometric detectors and spectrographs at NIR wavelengths in the late 1990s to populate the UCD region in the HR diagram (see sect. 2.2 in Reid & Hawley 2005 for a historical review). Today, it is known that UCDs represent at least 15% of the population of star-like objects in the solar neighborhood (see, e.g., Bartlett et al., 2017, and the RECONS collaboration111See www.recons.org).

In parallel with the increased observational access to UCDs, theoretical evolutionary models have improved in the past decades. Key advances include D’Antona & Mazzitelli (1996) who generated models for metal-poor stars, Burrows et al. (1993, 1997, the ”Tucson group”) who included grain opacities in atmosphere models, and Baraffe et al. (1995, 1998, the ”Lyon group”) who coupled interior models to full model atmospheres as boundary conditions. The most recent and commonly used model set is Baraffe et al. (2015) (hereafter, BHAC15), who presented updated evolutionary models for low-mass stars and young BDs that use boundary conditions from BT-Settl model atmospheres (Allard et al., 2012a, b; Rajpurohit et al., 2013, 2018b). These atmospheres include updated molecular opacity linelists and cloud formation, as well as atmospheric convection parameters calibrated with 2D/3D radiative hydrodynamic simulations (Freytag et al., 2010, 2012). BHAC15 models also adopt the solar composition of Asplund et al. (2009), supplemented with abundances from Caffau et al. (2011) for C, N, O, Ne, P, S, K, Fe, Eu, Hf,Os, and Th. The BHAC15 models provide a significant improvement over the Baraffe et al. (1998) models when comparing to observations. In particular, the weaknesses in the Baraffe et al. (1998) colour-magnitude diagrams (e.g., for the optical ($V-I$) colors that were too blue for a given magnitude; and for the NIR colors for relatively old objects) were significantly improved by BHAC15.

The BHAC15 grid is publicly available and is widely used when trying to retrieve the parameters of UCDs (e.g. age, mass) from observations. This grid is given for solar composition with stellar mass from $1.40M_{\odot}$ down to $0.01M_{\odot}$, with the lowest $T_{\rm eff}=1200$ K for 0.01 $M_{\odot}$, which corresponds to an age of 41 Myr. At low effective temperatures, the limit of validity of evolutionary models is set by the limit of validity of model atmospheres (Saumon & Marley, 2008). These models are likely robust for $T_{\rm eff}\gtrsim 2000$ K, as below this temperature there are significant uncertainties associated with cloud formation and condensation of chemical compounds which become important opacity sources in the atmosphere (Tsuji et al., 1996; Marley & Leggett, 2009; Morley et al., 2012).

In this paper, we present adaptations made to our in-house evolutionary code CLES (Code Liégeois d’Evolution Stellaire; Scuflaire et al. 2008) to produce UCD models: relevant equations-of-state (EOS) suitable for the low temperature and high density regime, and relevant model atmospheres used as boundary conditions for the interior. The main motivation behind this work is the SPECULOOS project (Search for Planets EClipsing ULtra-cOOl Stars; Gillon 2017; Delrez et al. 2018), a survey searching for transiting planets in the habitable zone of the nearest and brightest UCDs. The prototype of the SPECULOOS project is the TRAPPIST telescope (TRAnsiting Planet and PlanetesImals Small Telescope; Gillon et al. 2011) which led to the discovery of seven Earth-sized planets transiting the ultracool dwarf TRAPPIST-1 (Gillon et al., 2016, 2017). More generally, transiting surveys are now focused on very cool stars, because of their small radii with close-in habitable zones. Atmospheric characterization of these exoplanets will be within reach of the next-generation telescopes (e.g. James Webb Space Telescope, European Extreme Large Telescope). In this context, precise estimates of the mass, radius, luminosity, effective temperature and age for a host star are important to thoroughly characterize its exoplanets (see, for TRAPPIST-1, Van Grootel et al., 2018).

This paper is organized as follows. Section 2 provides the details of our CLES evolutionary models for UCDs, general properties and systematic error estimate, and comparison to BHAC15 UCD models. Section 3 gives a series of test-cases for model comparison to observations. Section 4 presents our conclusions.

2 CLES models for ultracool dwarfs

2.1 Input physics

The CLES evolutionary code was developed in the early 2000s by the asteroseismology group of the University of Liège, and has been continuously updated. We present here only the main input physics used for UCDs, and refer to Scuflaire et al. (2008) for the main numerical features (see also Buldgen et al., 2016, 2017a, 2017b, for CLES in the context of the Sun and solar-like stars). CLES includes different choices for EOS, opacity and atmosphere tables but here we mention only those relevant for the study of UCDs.

We considered H, He, C, and O for the EOS of UCD objects. We directly adapted tables built for white dwarfs and subdwarf B stars (Brassard & Fontaine, 1994). These tables cover a large domain of the temperature-density plane ($2.10\leq\log T\leq 8.98$ and $-12.0\leq\log\rho\leq 9.0$), which includes UCDs at all evolutionary stages. Four EOS tables are available, one for each element considered. For each table, three regimes are invoked. First, for the low-density region, a network of Saha equations is solved for a mixture of radiation and an almost ideal (including Coulomb corrections), non-degenerate, partially or fully ionized gas composed of a mixture of H, He, C, and O in various proportions. Second, in the partial ionization region where non-ideal and electron degeneracy effects are important (intermediate densities), we used the EOS of Saumon et al. (1995) for H and He, an improved version of the EOS of Fontaine et al. (1977) for C and for O. Third, the high-density domain corresponds to the fully ionized plasma in liquid and, ultimately, solid phases according to the physics described in Lamb (1974) and improved by Kitsikis et al. (2005). The low-density boundary of the second regime matches very smoothly with the high-density boundary of the first regime, thus ensuring that there are no significant jumps in the thermodynamic variables of interest. The connection between the high-density boundary of the partial degeneracy regime with the low-density boundary of the totally ionized domain is made at the location where the electron degeneracy parameter $\eta$ is equal to 20. Care has been taken to insure that thermodynamic consistency is respected as explained in Fontaine et al. (1977). Interpolation in composition is handled following the additive volume prescription of Fontaine et al. (1977). This additive volume prescription is the best viable option for handling mixture in the EOS of individual elements (Vorberger et al., 2007; Wang et al., 2013; Danel & Kazandjian, 2015). The FreeEOS equation of state222See http://freeeos.sourceforge.net/ by A. Irwin. (Irwin, 2012) version 2.2.1. (EOS1 configuration) is also available in CLES. It calculates the EOS without radiation pressure for specified mixtures using an efficient free-energy minimization technique. However, this did not converge for the lowest temperatures/highest densities encountered in H-burning stars below $\sim 0.14$ $M_{\odot}$ at solar composition.

Model atmospheres based on $T(\tau)$ relations (including grey atmospheres) are inaccurate as boundary conditions (BCs) for interior models of UCDs (Chabrier & Baraffe, 1997). We have implemented in CLES two sets of surface boundary conditions (density, temperature, geometrical depth) from detailed model atmospheres. The first derives from the publicly available BT-Settl model atmospheres (Allard et al., 2012a, b). The second makes use of the model atmospheres of B. Aringer, originally developed for Asymptotic and Red Giant Branch stars (Aringer et al., 2016, hereafter, AR16). In both cases, the transition between interior and atmosphere is performed at the Rosseland mean optical depth $\tau=100$ (similar to the BHAC15 models), as a safe limit to avoid discrepancies at the boundary while treating convection and adiabatic processes (Chabrier & Baraffe, 1997). For BT-Settl model atmospheres and interior structure, we used the models computed with the solar abundances of Asplund et al. (2009), as for the interior structure, and followed the same heavy-element abundance as for the Sun (i.e., $[\alpha/H]=0$). This results in four distinct solar-scaled compositions: $[M/H]=-0.5,0.0,+0.3,+0.5$. The AR16 models are computed for Asplund et al. (2009) solar (and meteoritic when available) abundances as well, supplemented by Caffau et al. (2011) for C, N, O, Ne and Ar. Seven compositions are possible, from $[M/H]=-2$ to 1, by steps of 0.5. AR16 is computed with the COMARCS program which is based on the version of the MARCS code of Gustafsson et al. (2008), and uses the COMA (Copenhagen Opacities for Model Atmospheres) opacity generation code by Aringer (2000). Unlike the BT-Settl models, AR16 models are dust-free and thus limited to $T_{\rm eff}\geq$ 2600 K.

CLES includes opacities from the OPAL (Iglesias & Rogers, 1996) or the OP (Badnell et al., 2005) project. In both cases, these opacities are combined for low temperatures to opacities from Ferguson et al. (2005). The effects of thermal conductivity have been taken into account following Potekhin et al. (1999) and Cassisi et al. (2007).

Nuclear reaction rates for the pp chain come from the review of Adelberger et al. (2011), except for the 7Li$(p,\alpha)^{4}$He reaction which comes from the NACRE II compilation (Xu et al., 2013).

Convection is treated using the mixing length theory (MLT). For UCDs, we generally set $\alpha_{MLT}$ (the ratio between the mixing length and the pressure scale height) to $\sim 2.0$, the value adopted in BT-Settl model atmospheres for such stars (BHAC15). AR16 model atmospheres also follow the MLT formalism, with $\alpha_{MLT}=1.5$, which we adopt for the interior when using these models for surface boundary condition. In Sect. 2.3 we examine the effect of $\alpha_{MLT}$ on the evolution of UCDs.

CLES solar calibration (evolutionary track giving the Sun luminosity and effective temperature at its present age), without diffusion, with OPAL opacities, with our standard H+He+C+O EOS, and BCs from BT-Settl model atmospheres, gives $\alpha_{MLT}=1.8$, $X_{0}=0.729$, and $Z_{0}=0.013$. This is our standard CLES configuration for UCDs (with $\alpha_{MLT}=2.0$). We note that this solar calibration does not depend highly on the chosen input physics, that is model atmospheres, EOS, solar abundances Asplund et al. 2009/Caffau et al. 2011, and opacities.

2.2 Properties of CLES models

2.2.1 H-burning limit mass

The formal H-burning limit mass ($M_{HBL}$), i.e., the mass where the fractional contribution to the total luminosity due to hydrogen fusion $L_{nuc}/L_{total}$ never exceeds 50% (Reid & Hawley, 2005), is slightly above 0.078 $M_{\odot}$ in standard CLES configuration (Fig. 1). $M_{HBL}$ is about 0.07 $M_{\odot}$ in the BHAC15 models according to their public grids, and 0.073 $M_{\odot}$ for Burrows et al. (1993, 1997) models. Figure 1 shows that stars with mass slightly higher than 0.08 $M_{\odot}$ achieve stable luminosities and temperatures for many Hubble times, while $0.080M_{\odot}$ and $0.079M_{\odot}$ objects are transition objects, able to sustain fusion for several hundreds of million years, but eventually cooling degenerately as brown dwarfs. The region where $L_{nuc}/L_{total}>0.5$ at early ages for all masses in Fig. 1 corresponds to the short-lived phase of D-burning. The $M_{HBL}$ depends on the chemical composition of the star such that the mass limit increases as metallicity decreases (see also section 3.4 in Reid & Hawley, 2005). For example, in the CLES models the $M_{HBL}$ is 0.074, 0.078 and 0.080 $M_{\odot}$ at $[M/H]=+0.5,+0.0$ and $-0.3$, respectively.

2.2.2 Abundance of light elements

The fusion of light elements, namely Li, Be, B, can occur in the interior of UCDs depending on their mass and age. Presence/absence of light elements in the spectra provide a powerful diagnostic to distinguish UCDs of various ages and masses, including BDs (Chabrier & Baraffe 1997). For example, observing the lithium doublet at $\lambda LiI=6708\AA$ was proposed by Rebolo et al. (1992) as the famous lithium test that confirmed the first BD candidates (Rebolo et al., 1995; Basri et al., 1996).

We computed the minimum burning mass for 7Li, 9Be and 10B in CLES models to be 0.053 $M_{\odot}$, 0.065 $M_{\odot}$ and 0.079 $M_{\odot}$, respectively, in agreement with literature values, e.g., 0.055 $M_{\odot}$, 0.065 $M_{\odot}$ and 0.08 $M_{\odot}$, respectively, for Chabrier & Baraffe (1997), and of about 0.055 $M_{\odot}$, 0.07 $M_{\odot}$ and 0.09 $M_{\odot}$, respectively, for Burrows et al. (1997). Fig. 2 shows the boundary for lithium depletion for when $90\%$ of lithium has been consumed for different masses.

2.3 Systematic error in evolutionary models

2.3.1 The choice of EOS

In our standard EOS, C and O are proxies for all metals, with the same C/O proportion as in Asplund et al. (2009). BHAC15 (and older generation) models gather all metals, including C and O, in an increased fraction of He.333The mixture of H EOS and He EOS is also handled by the additive volume prescription in BHAC15 and older generation models. This is only valid when ionic pressure is negligible compared to electronic pressure, which is true through most of a UCD interior (Chabrier & Baraffe, 1997). Figure 3 quantifies this effect for the first time. We compare stellar luminosity and effective temperature as a function of age, at solar composition, for 0.075 $M_{\odot}$ (below $M_{HBL}$), 0.08 $M_{\odot}$ (close to $M_{HBL}$), and 0.09 $M_{\odot}$ (above $M_{HBL}$) stars with the CLES standard configuration EOS (H+He+C+O; solid lines) and the EOS used in the BHAC15 models (H+He; dashed lines). Figure 3 shows that including the EOS for C and O has a strong effect at and close to $M_{HBL}$, with diverging evolutionary tracks at 0.08$M_{\odot}$. Close to the $M_{HBL}$ (0.079 and 0.081 $M_{\odot}$), the effect on the luminosity and effective temperature is $+4\%$ and $+1\%$, respectively. Whereas farther from $M_{HBL}$ (0.09$M_{\odot}$ and higher masses; 0.075$M_{\odot}$ and lower), the impact is smaller, about $+1\%$ in luminosity and $+0.3\%$ in effective temperature. The conclusion of this experiment is that assimilating all of the metals into He is generally a valid hypothesis, but care must be taken close to the $M_{HBL}$. Some modern stellar evolution codes (e.g. DARTMOUTH, Feiden 2016, and references therein; and PARSEC, Chen et al. 2014), used to model very low-mass stars, provide public grids down to $\sim$ 0.10 $M_{\odot}$, based on the FreeEOS EOS. We found it impossible to make FreeEOS converge for densities/temperatures corresponding to stellar masses below $\sim$ 0.14 $M_{\odot}$. It is plausible that the differences in minimum mass found by different groups come from choosing different FreeEOS configurations when computing its tables, in particular regarding the numerous ionization states of various elements. In Fig. 4, we show the evolution in luminosity of a 0.14 $M_{\odot}$ star with FreeEOS and our standard EOS. Differences observed at canonical age are about $2\%$ in luminosity and $0.2\%$ in effective temperature.

2.3.2 The choice of model atmospheres as BCs

Figure 5 compares luminosity and effective temperature as a function of age, at solar composition, for 0.09 $M_{\odot}$ and 0.13 $M_{\odot}$ with BCs extracted from BT-Settl model atmospheres and AR16 models. Systematic error is small in this case, with maximum differences (in H-burning phase) of $1-2\%$ in luminosity, $0.4-0.7\%$ in effective temperature, and $0.2-0.4\%$ in radius. Models built with BCs based on BT-Settl atmospheres, which include grain formation below $T_{\rm eff}<$ 2600 K, are used in a more extended range of stellar mass (including young BDs) than those based on AR16 models, but the latter is more extended in terms of metallicities.

2.3.3 Other sources of systematic error

We computed evolutionary tracks with OPAL and OP opacity tables for various masses and various ages at solar composition. OP and OPAL tables give very similar results, i.e., no systematic error is observed by choosing OP rather than OPAL.

We also computed evolutionary tracks with Caughlan & Fowler (1988) nuclear reaction rates, which were implemented in some older versions of CLES. We found a typical $+1\%$ increase in luminosity for 0.09 $M_{\odot}$ at 10 Gyr using Caughlan & Fowler (1988) instead of Adelberger et al. (2011).

It is difficult to give estimates linked to the ’imperfect’ description of convection by the mixing length theory. We note that for UCD objects, the exact value of $\alpha_{MLT}$ (1.5, 2.0 or 2.2) has very little influence on model global parameters such as luminosity, radius, and effective temperature, implying that systematic error linked to MLT theory is negligible. Chabrier & Baraffe (2000) also reports that varying $\alpha_{MLT}$ value between 1 and 2 in the interior has no impact for stars below 0.60 $M_{\odot}$ (see also Montalbán et al. 2000). We observed the same in our models of UCDs. Only very small values of $\alpha_{MLT}$ significantly changes UCD parameters (see section 3.2).

2.3.4 Systematic error in models: conclusions

In order to provide estimates of the typical systematic error associated to evolutionary models, we carried out a series of experiments reproducing the traditional use of evolutionary models by observers, which consists in inferring a mass from luminosity measurements (determined from parallax and spectral energy distributions typically; see Filippazzo et al. 2015). In order to ’translate’ systematic error in luminosity by changing input physics (EOS, model atmospheres, nuclear reaction rates) into a systematic error in mass, we identified which shift in mass (for models of a given age) would give a shift in luminosity that is typical for the uncertainty in the measurement of luminosity. We found that all systematic shifts observed in luminosity correspond to a modest change in mass, typically between 0.0001 and 0.0005 $M_{\odot}$. We propose that the typical systematic error in mass from evolutionary models, in the range of mass of UCD objects, is 0.0005 $M_{\odot}$.

2.4 Comparison to existing UCD models

The BHAC15 models are commonly used when characterizing ultracool dwarfs. To compare CLES with BHAC15 models, we have adapted, as far as possible, identical input physics at solar metallicity: Grevesse & Noels (1993) abundances for the interior models and Asplund et al. (2009) supplemented by Caffau et al. (2011) for some elements for BCs from BT-Settl model atmospheres; EOS for H and He only, with an increased fraction of He for assimilating metals; identical initial composition (I. Baraffe, priv. comm.); OPAL opacities, no diffusion, $\alpha_{MLT}=1.6$. A comparison between CLES models and those of BHAC15 with identical input physics has already been presented in Van Grootel et al. (2018) with evolutionary tracks for 0.08, 0.09 and 0.10 $M_{\odot}$ stars. Both models have similar luminosity and differences in effective temperature and stellar radius are of about 1% and 3%, respectively. We complement this study by also including tracks for 0.130 and 0.06 $M_{\odot}$ (the latter below $M_{HBL}$), seen here in figure 6.

3 Test cases

We test our CLES models by comparing theoretical results with observations for objects for which stellar parameters (luminosity, dynamical mass, age and/or inferred radius from interferometry) have been estimated from independent techniques: the M7 binary LSPM J1314+1320AB (Dupuy et al., 2016, hereafter D16), and the spectroscopic twins GJ 65AB and Proxima Centauri (Kervella et al., 2016, hereafter K16). We also compare these results with those from BHAC15 models.

3.1 LSPM J1314+1320AB

LSPM J1314+1320AB is a pre-main sequence, nearby binary which has precise measurements of the total dynamical mass and integrated luminosity ($M_{tot}=0.1761\pm 0.0015M_{\odot}$ and $\log(L_{tot,bol}/L_{\odot})=-2.322\pm 0.009$, respectively), as well as the individual masses and luminosities within $\lesssim 1\%$ precision, obtained through spatially-resolved absolute and relative astrometric monitoring and optical and NIR photometry and spectroscopy (D16). The metallicity of LSPM J1314+1320AB has been estimated by D16, using the calibration of Mann et al. (2014), to [Fe/H]=0.04$\pm$0.08, thus consistent with solar metallicity.

In a first test, we used the individual masses and luminosities to model LSPM J1314+1320AB, assuming solar composition and our standard CLES configuration (see Sect. 2.1). Given the A and B components are nearly equal in mass and luminosity, we took the averaged values $\langle M\rangle=0.0881\pm 0.0008M_{\odot}$ and $\langle\log(L/L_{\odot})\rangle=-2.623\pm 0.010$. Results are presented in Table 1. The quoted errors are from propagating errors on measured mass and luminosity. More precisely, we computed various evolutionary tracks by varying observational constraints ($L$, $M$) within their given 1-$\sigma$ range and computed the respective 1-$\sigma$ confidence interval for output quantities $T_{\rm eff}$, age, $R$, and log $g$. Table 1 also shows results obtained with BHAC15 models, as determined by D16, which are similar. Note that, the derived $T_{\rm eff}$ is about 180 K hotter than that given by the spectral type-$T_{\rm eff}$ relation of Herczeg & Hillenbrand (2014) used by D16, which gives $T_{\rm eff}=$ 2770 $\pm$ 100 K for LSPM J1314 A and B.

In a second test, we reproduced the test proposed by D16 that mimicked the typical application of models by observers, i.e., that we have access to a measurement of the luminosity and to an effective temperature based on a spectral type$-T_{\rm eff}$ relation, but not to stellar masses. To carry out this test, we used our in-house Levenberg-Marquadt optimization algorithm (Press et al. 1992, see Van Grootel et al. 2018 for details). We first assumed solar composition and obtained a stellar mass $M=0.049\pm 0.017M_{\odot}$ (51$\pm 18M_{\rm Jup}$) with age $=25\pm 15$ Myr (Fig. 7, dashed black curve). Here again, quoted errors simply come from error propagation on $T_{\rm eff}$ and $L$. No systematic error (see Sect. 2.3) was included. This result is consistent with the values from D16 interpolated from BHAC15 models: $M=50^{+20}_{-13}M_{\rm Jup}$ and age $=25^{+10}_{-17}$ Myr, still showing the discrepancy between mass of stellar models and direct measurements. Adjusting the metallicity of LSPM J1314 to [Fe/H]=0.12 (+1$\sigma$) yields a higher mass $M=0.057\pm 0.019M_{\odot}$, closer to the measured value. To fully reconcile models with observations, we found that we have to increase the metallicity to +2.5$\sigma$, i.e., [Fe/H]=0.24 (Fig. 7, dot-dashed green curve). Increasing the metallicity to recover stellar parameters was also found in the study of TRAPPIST-1 (Van Grootel et al., 2018). In the same line, Lindgren & Heiter (2017) and Rajpurohit et al. (2018a) recently observed a typical average deviation of 0.2 to 0.4 dex in [Fe/H] between direct spectral fitting and calibration-based techniques. One possibility is that metallicity is underestimated, another is residual error in models that can be ’adjusted’ by artificially increasing metallicity. A third possibility is that the spectral type$-T_{\rm eff}$ relation depends on the scale in use: for a same spectral type, Herczeg & Hillenbrand (2014), Rajpurohit et al. (2013) and Passegger et al. (2018) show differences of 10% in $T_{\rm eff}$ in the M-dwarf regime.

3.2 Proxima Centauri and GJ65 AB

Proxima Centauri (GJ 551, M5.5Ve) is the nearest star after the Sun. It most likely part of a gravitationally bound triple system with $\alpha$ Centauri AB, with an orbital period of about 550 000 yr (Kervella et al., 2017). GJ65 AB (BL+UV Ceti, M5.5Ve+M6Ve) is a pair of red dwarfs in a sufficiently large binary to neglect gravitational and magnetic interactions, but sufficiently close to resolve the orbital motion within a human lifetime, with $P_{\rm orb}\sim 26.3$ yr (K16). These three nearby stars are among the most studied red dwarfs, and as such constitute cornerstones for models of very low mass stars. Although they are not strictly speaking UCDs, we can still confront our models to these unique benchmarks.

Proxima and GJ65 AB all have three direct radii measurements from interferometry with, respectively, $0.141\pm 0.007R_{\odot}$ (Demory et al., 2009), $0.165\pm 0.006R_{\odot}$ and $0.159\pm 0.006R_{\odot}$ (K16). By orbital modeling, the total mass of GJ65 is accessible, as well as individual masses based on differential astrometry, indicating $0.1225\pm 0.0043M_{\odot}$ and $0.1195\pm 0.0043M_{\odot}$ for A and B components, respectively (K16). No direct mass measurement of Proxima exists, but it can be estimated through mass-absolute magnitude relations, $0.1239\pm 0.0032M_{\odot}$ using Mann et al. (2019) at solar-metallicity (distance and $K_{S}$ magnitude in Cutri et al. (2003)), $0.123\pm 0.006M_{\odot}$ using Delfosse et al. (2000) relations, and $0.118\pm 0.011M_{\odot}$ using Henry et al. (1999). The mass-luminosity relations obtained with CLES models are shown in fig. 8 for solar-metallicity, giving us an estimated mass of $0.120\pm 0.001M_{\odot}$ for Proxima at 5 Gyr assuming a luminosity $L=0.00155\pm 0.00002L_{\odot}$ (Boyajian et al., 2012). Directly determining metallicities from high-resolution spectroscopy of red dwarfs is a delicate task. For Proxima, iron abundance measurements range from [Fe/H]$=-0.07\pm 0.14$ (Passegger et al., 2016) to [Fe/H]$=+0.16\pm 0.20$ (Neves et al., 2014). We could also consider that its membership in the $\alpha$ Cen system indicates a common origin from the same formation cloud, and thus the same age and initial composition. The metallicities of $\alpha$ Cen A and B have been accurately determined from atmospheric abundances of many elements (e.g. Porto de Mello et al., 2008), and stellar modeling and asteroseismology allows us to obtain the initial composition and age, with $X_{0}\sim 0.70$, $Z_{0}\sim 0.025$ and an age about 6 Gyr (Bazot et al., 2016). It is also possible that Proxima has been captured by $\alpha$ Cen, and thus does not share a common composition. Whatever the actual metallicity of Proxima, K16 adopted a differential approach to Proxima to determine the metallicity of GJ65 AB. Adopting [Fe/H]$=+0.05\pm 0.20$ for Proxima, gave [Fe/H]$=-0.03\pm 0.20$ and [Fe/H]$=-0.12\pm 0.20$ for GJ65 A and B.

K16 used the BHAC15 models (with solar composition) to model Proxima and GJ65 AB. They found that, while Proxima fits with the expected mass-radius relation, GJ65 AB both appear inflated, exceeding model expectations by 14$\pm$4% and 12$\pm$4%, respectively. This radius inflation could naturally be explained by a young age (about 250 Myr), but the GJ65 velocity vector likely indicates a star of at least 1 Gyr, and possibly much older. K16 carefully examined possible sources of discrepancy, and concluded that the most likely explanation is a reduced convection efficiency for GJ 65AB by a strong internal magnetic field linked to the relatively fast rotation of both stars ($v\sin i\sim$ 30 km s*-1*). Proxima, on the contrary, is a slow rotator ($v\sin i\sim$ 2 km s*-1*).

We modelled Proxima and GJ65 AB using the same optimization procedure as in the previous section. First, we placed GJ65 AB and Proxima in our CLES mass-radius (M-R) plot with 0.1-5.0 Gyr isochrones, for solar composition (Fig. 8). The masses reference for Proxima plotted in Fig. 8 are the values from Delfosse et al. (2000) and Henry et al. (1999) (red and orange, respectively) in which we also include our result (blue). We agree with the K16 conclusions: while Proxima reasonably fits the model expectations, GJ65 AB is either a young star (about 300 Myr), or somewhat inflated compared to its mass. In more detail, we found that models slightly overestimate the stellar radius of Proxima, opposite to the usual trend (see, e.g., section 3.1 and Van Grootel et al. 2018). The situation would be worse considering Proxima to have the same metallicity as $\alpha$ Cen. We tested the suggestion of K16 to decrease the convection efficiency for GJ65 and found that we had to turn $\alpha_{\rm MLT}$ down to 0.03 (and even slightly lower for GJ65 A) in order to reconcile the measured radii with model expectations. If we play with metallicity, we have to increase it by a great amount (more than [Fe/H]=+0.5) to reconcile interferometric and model radii. This is also unlikely given the differential spectroscopy carried out by K16, which indicates that Proxima and GJ65 have similar metallicities.

We reversed the problem by identifying which stellar mass would fit the radii measured from interferometry, assuming several Gyrs stars in both cases. For Proxima we found a mass of 0.113$\pm$0.007 $M_{\odot}$. This is within 2$\sigma$ of the value from Delfosse et al. (2000) relations, and within 1$\sigma$ of the Henry et al. (1999) ones. Proxima is at the limit of validity for applicability of such relations (11, 10, and 9.5 for the $JHK$ absolute magnitudes, respectively, for Delfosse et al. 2000). Indeed, applying such relations for GJ65 (also at the validity limit), K16 found that predicted masses are 12% and 17% lower than the measured masses, respectively. Repeating the same exercise with GJ65, we found that we are within 2$\sigma$ of the dynamical value: a mass of 0.137$\pm$0.006 $M_{\odot}$ is needed to account for the radius of 0.165$\pm$ 0.006 $R_{\odot}$ (GJ 65A), and a mass of 0.132$\pm$0.006 $M_{\odot}$ is needed to account for the radius of 0.159$\pm$ 0.006 $R_{\odot}$ (GJ 65B).

In conclusion, given uncertainties in stellar models in metallicity and $\alpha_{MLT}$ for stars with inhibit convection (e.g., fast rotating or strong magnetic field, see review in Brun & Browning 2017), and the uncertainties that are likely to occur in mass-absolute magnitude relations, we posit that an agreement within 2$\sigma$ is acceptable.

4 Conclusion

We have presented new evolutionary models for $T_{\rm eff}\geq$ 2000 K UCD objects, which encompass the very low-mass stars and young, still contracting brown dwarfs. These models are based on our in-house evolutionary code CLES adapted to produce UCD models. In particular, we included a relevant EOS that includes H, He, C and O elements, and appropriate boundary conditions from two sets of model atmospheres. We presented some properties of our models and investigated their systematic error associated with uncertainties in input physics. We showed a comparison with the reference BHAC15 models in the UCD regime. Finally, a series of test-cases was carried out to visualize the strengths and limits of our models.

Tables for our CLES models for various UCD masses and metallicities can be found at http://www.astro.ulg.ac.be/ASTA/cles-models-UCDs/

The authors thank the anonymous referee for an helpful review. We warmly thank Didier Saumon for various advices when developing our UCD models. We thank the whole SPECULOOS team (http://www.speculoos.uliege.be/) for motivation, advice and encouragement when developing these models. We thank Uffe G. Jørgensen for his comments after careful reading of the manuscript. C.S.F. is funded by an Action de Recherche Concertée (ARC) grant financed by the Wallonia-Brussels Federation. V.V.G. is F.R.S.-FNRS Research Associate. The research leading to these results has received funding from the ARC grant for Concerted Research Actions, financed by the Wallonia-Brussels Federation. B.A. was supported by the ERC Consolidator Grant funding scheme (project STARKEY, G.A. n. 615604). G.F. acknowledges the contribution of the Canada Research Chair Program.