© ESO, 2014

1. Introduction

The primary goal of NASA’s forthcoming exoplanet mission, Transiting Exoplanet Survey Satellite (TESS), is to detect small transiting planets around bright, nearby host stars. Due to the increased signal-to-noise ratio (S/N) in the spectroscopic observations obtained for a bright host star, it is much easier to derive both the mass of an orbiting planet using the radial velocity technique and to measure the spectroscopic signatures from the planet’s atmosphere. To identify bright planet hosting stars over the whole sky, TESS will reside in a High-Earth Orbit and continuously monitor each field of target stars for 27 consecutive days. While stars that reside in regions where the fields overlap will have longer duration light curves, the majority of newly discovered TESS planets will have relatively short orbital periods (<27 days). The consequence of this observing strategy is that potentially habitable worlds where liquid water could exist on the surface will only be found around cool, very late-type M dwarf stars. Luckily, intrinsically faint M dwarf stars also have relatively high contrast ratios between the star and the orbiting planet; thus, these systems will be optimal targets for detecting atmospheric signatures from the planet itself.

As always, to fully understand the planets, it is paramount to accurately characterise their host stars. The mass of a planet hosting star, which directly determines the derived planet mass, is typically obtained by comparing measurable star properties (e.g., colours, T_eff, luminosity) to theoretical stellar evolution models (e.g., Dotter et al. 2008) and/or empirical relationships (e.g., Torres et al. 2010). Although M dwarfs comprise the majority of stars in the Galaxy, our understanding of the relationships between their masses, radii, temperatures, and metallicities is still incomplete – particularly at the very bottom of the main sequence. For example, recent analysis of a newly discovered eclipsing M dwarf (KIC 1571511) with a mass of ~0.2 M_⊙ suggests its temperature to be hotter than theoretical models predicted by ~900 K (Ofir et al. 2012). If this is typical, despite the accurate radii benefiting from Gaia parallaxes, the mass of the planets identified by TESS that are orbiting stars in the habitable zone will not be accurately characterised. A temperature that is hotter than expected for planet hosting M dwarfs would mistakenly imply a higher stellar mass and, consequently, a higher mass for the planet. This erroneous characterisation is not unique to the transiting method but would also affect other kinds of systems, such as those discovered astrometrically by Gaia (Perryman et al. 2001).

For low mass stars, the majority of existing knowledge necessary to calibrate relationships between fundamental properties has come from two types of systems: (1) nearby, single stars with interferometric radii measurements (Ségransan et al. 2003; Berger et al. 2006; Demory et al. 2009; Boyajian et al. 2012, and references therein) and (2) M+M dwarf eclipsing binaries (EBs; Lacy 1977; Leung & Schneider 1978; Metcalfe et al. 1996; Torres & Ribas 2002; Delfosse et al. 1999; Ribas 2003; López-Morales & Ribas 2005; Morales et al. 2009; Nefs et al. 2013). However, the number of very late-type stars (M_⋆≤ 0.25 M_⊙) for which these analyses can be done is extremely small. In the literature, there are only 18 measurements of stellar mass and radius of very low mass stars; of those, only 7 have temperatures (Nefs et al. 2013; Zhou et al. 2014, and references therein). Furthermore, deficiencies are inherent in these techniques that contribute to our insufficient knowledge of the properties of very late-type M dwarfs. Specifically, single interferometric systems allow for the measurement of radius, effective temperature (T_eff), and sometimes metallicity (with large uncertainties), but not mass. Furthermore, M+M EBs provide direct measurements of the masses and radii of both components and their relative temperatures, but individual temperatures and metallicities are difficult to derive from the complex spectra of two unresolved M dwarfs.

To address these issues, we have an ongoing program designed to measure significant numbers of masses and radii of very low mass stars with accurate metallicity and temperature determinations by analyzing M dwarfs in eclipsing systems with higher mass F, G, or K stars (Triaud et al. 2013). Hereafter, we refer to these systems as EBLMs. In an EBLM, the primary star dominates the light allowing an accurate temperature and metallicity to be determined from the relatively simple and well-understood spectrum of the F/G/K star. The mass of the primary is derived from its stellar parameters and then used in combination with the radial velocity curve and light curve to find the mass of the M dwarf, the radius of both components, and the M dwarf temperature. Finally, the metallicity derived for the primary star is adopted for the M dwarf secondary by assuming the close binary pair formed from the same parent molecular cloud.

Our program is based around a large sample of EBLM systems (~150) discovered in the SuperWASP survey. SuperWASP (Pollacco et al. 2006) is a dedicated, ultra-wide field sky survey, continuously monitoring stars of V ~ 8–15 mag over a quarter of the sky every clear night. The SuperWASP archive is a rich source of new eclipsing binaries and, in particular, systems with low mass secondaries are found in the course of exoplanet candidate selections. These objects are bright and their light curves have thousands of data points obtained over multiple years show clear eclipse signatures with well-defined periods. Although EBLMs are sources of false alarm detections when searching for exoplanets, they are ideal objects to use for determining fundamental parameters of very low mass stars.

Here, we present the study and characterisation of a newly discovered EBLM from the SuperWASP survey, J0113+31. The binary is composed of a G0–G2 V primary and an M₂ = 0.19 M_⊙ secondary. The system is eccentric with an orbital period of 14.277 days. Like KIC 1571511, we find the temperature of the M dwarf to be significantly hotter than stellar evolution models predict for a star of this mass and metallicity. The following paper describes our analysis of this interesting system. It is structured as follows: in Sect. 2 we describe the photometric and spectroscopic observations utilised in the analysis of the eclipsing binary; in Sect. 3, we describe our data analysis and our approach to the modelling of the system, and in Sect. 4, we discuss our results and set them in context of other M-dwarf measurements. Finally, we draw our conclusions in Sect. 5.

2. Observations and data reduction

In this section, we describe the data utilised to derive the physical properties of the orbit and the stellar components of the eclipsing binary. Although J0113+31 was discovered from its SuperWASP time-series photometry, the survey-quality light curve (Sect. 2.1; Fig. 1) does not have sufficient precision to model the eclipses. Thus, we obtained higher precision light curves of the primary eclipse at optical wavelengths with the 0.4-m Near Infra-red Transiting Exoplanet Survey (NITES) telescope at La Palma (Sect. 2.2), the 0.6-m telescope at the Observatori Esteve Duran (OED; Sect. 2.3), and the 0.91-m telescope at the BYU West Mountain Observatory (Sect. 2.4). We also acquired a light curve of the secondary eclipse in the near-infrared 2.1-m telescope using the instrument FLAMINGOS at the Kitt Peak National Observatory (KPNO; Sect. 2.5). Additionally, we characterised the reflex motion of the primary star around the system’s centre of mass via radial velocity measurements (Sects. 2.6–2.7).

Fig. 1 SuperWASP light curve. We show the time-series photometry that led to the discovery of J0113+31, as part of the SuperWASP transiting planet survey, overplotted on the final model light curve from the EB modelling (Sect. 3.5) in the V-band, which is used as a proxy for the non-standard SuperWASP filter. The primary eclipse of the M dwarf blocking the light of the more massive star has a transit-like shape (shown in inset). The error bar to the bottom left represents the median uncertainty of the SuperWASP photometric points, which are excluded for clarity.

2.1. Photometry-optical: SuperWASP

The SuperWASP North telescope is located in La Palma (ORM – Canary Islands). The telescope consists of eight Canon 200 mm f/1.8 focal lenses coupled to e2v 2048 × 2048 pixel CCDs, which yield a field of view of 7.8 × 7.8 square degrees, and a pixel scale of 13.7′′ (Pollacco et al. 2006).

Our object J0113+31 was observed from 2004 June 11 through to 2011 January 16 with a total of seven cameras and 29 324 photometric points covering 357 primary eclipses with a total of 5607 points during eclipse. The SuperWASP data were first processed with the custom-built reduction pipeline described in Pollacco et al. (2006). The resulting light curves were analysed using our implementation of the box least-squares fitting and SysRem de-trending algorithms (see Collier Cameron et al. 2006; Kovács et al. 2002; Tamuz et al. 2005) to search for signatures of planetary transits. Once the candidate planet is flagged, a series of multi-season, multi-camera analyses are performed to strengthen the candidate detection. In addition, different de-trending algorithms (e.g., TFA, Kovács et al. 2005) were used on one season and multi-season light curves to confirm the transit-like signal and the physical parameters of the putative planet candidate. These additional tests allow a more thorough analysis of the system parameter derived solely from the SuperWASP data thus helping to identify the best candidates, and to reject possible spurious detections. At this time, J0113+31 was flagged as a candidate exoplanet. Subsequently, during eye-balling of the SuperWASP light curves, their morphology and the characteristics of the photometric signal lead us to revise J0113+31 as an EBLM. The SuperWASP light curve folded on the ephemerides derived in Sect. 3.2 is shown in Fig. 1.

2.2. Photometry-optical: NITES

The NITES telescope is an f/10 0.4 m Meade LX200GPS with advanced coma-free optics. The CCD camera is a Finger Lakes Instrumentation (FLI) Proline 4710 with a back illuminated 1024 × 1024, 13 μm/pixel, deep-depleted CCD made by e2v (for more details see McCormac et al. 2014).

Two partial light curves were observed with the NITES telescope on La Palma on the nights of 2011 September 18 and 2012 December 18. The telescope was defocused to 5′′ and 2385 images of 5 s exposure time were obtained during the September observations. On the December observations, the telescope was defocused to 4.6′′ and 1105 images with 5 s exposure time were obtained; data on both nights were taken without a filter. The data were bias subtracted, corrected for dark current, and flat-fielded using the standard routines in IRAF; and aperture photometry was performed using DAOPHOT (Stetson 1987). The stars GSC 2295-0229 and GSC N322212139 (01:13:35.3 +31:49:01 J2000) were used as comparison stars. Other stars in the field were rejected to minimise the scatter in the light curve. Shown in the bottom panel of Fig. 2, this data was only used in Sect. 3.2 to derive the best ephemeris but was excluded from the EB analysis from which the properties of the M dwarf are derived (Sect. 3.5).

Fig. 2 Follow-up photometry of the primary eclipse. On the top, we show the time series photometric points acquired at the OED (Sect. 2.3), phased to the ephemeris from Sect. 3.2, and overplotted with the best fit EB model described in Sect. 3.5 (red line). We also show directly below the residuals to the fit and the uncertainty of the individual photometric data points. Similarly in the middle panel, we show the data from BYU and the residuals to the fit. In the bottom panel, we show the partial primary eclipse data acquired with NITES. Overplotted in grey is a model light curve in the R-band, shown only for guidance and comparison to the upper panels. The NITES data were taken without a filter, and as such, we only utilise the NITES phtometry to derive the ephemeris, excluding it from the EB analysis (Sect. 3.5).

2.3. Photometry-optical: OED

Because of the long period of the system, it is challenging to acquire a full eclipse on a single night from a given location. Thus, photometry from Observatori Esteve Duran (OED) was performed on three nights with a ST-9XE CCD camera attached to the 0.6 m Cassegrain telescope and a 512×512 pixels CCD. Each pixel is 20 μm × 20 μm, giving an image resolution of 1.37″ × 1.37″ on the CCD plane. An Optec I-band Bessel filter was used for the observations and the exposure time was between 15 and 30 s. Photometry was extracted by performing synthetic aperture photometry. The star GSC 2295-0229 was used as comparison. The flux of our target was divided by the flux of the comparison star and differential magnitudes was then obtained. The final light curve (shown in Fig. 2) was obtained after binning of the original photometric points with a typical magnitude scatter of ~2–3 mmag. The features in the OED photometry (at the end of the ingress and at mid-transit) are most likely due to the normalisation of the three individual partial eclipses that comprise this light curve. Based on the residuals between the data and the model (Fig. 2), these features in the light curve are not systematically affecting the fit to the light curve.

2.4. Photometry-optical: BYU

The Brigham Young University (BYU) 0.91-m telescope is a f/5.5 system located at the West Mountain Observatory in Utah, USA. It is fitted with a Finger Lakes PL-09000 CCD camera with a 3056 × 3056, 12 μm/pixel array that gives a field of view (FOV) of 25.2′×25.2′ (Barth et al. 2011).

A near-complete I-band light curve of the primary eclipse was obtained on 2012 October 8. Pre-eclipse, ingress, flat bottom, and egress was obtained. The resulting I-band light curve consists of 2200 observations with a 17 s cadence and is presented in Fig. 2. The raw data were processed in the standard way. After the instrumental signatures were removed, source detection and aperture photometry were performed on all science frames using the Cambridge Astronomical Survey Unit (CASU) catalogue extraction software (Irwin & Lewis 2001). The optimal aperture radius was chosen through empirical testing of several different sizes. The median seeing of our observations was 3−4 pixels. We adopted a circular aperture with a 14 pixel radius, which is several times the FWHM of our relatively bright targets, to obtain the final instrumental magnitude measurements. Two nearby stars of comparable brightness (GSC 2295-0229 and GSC N322212139) in the FOV of the detector were chosen as comparison objects for deriving differential photometry. Other stars were considered as comparisons but were excluded to minimise the scatter in the photometry.

The total flux enclosed in the photometry aperture for the reference objects was divided by the instrumental flux of the target for each data point and then converted to magnitudes. All the measurements were then normalised using the out-of-eclipse portions of the light curve.

2.5. Photometry-NIR: FLAMINGOS

We observed J0113+31 with the Florida Multi-Object Imaging Near-Infrared Grism Observational Spectrometer (FLAMINGOS) in its imaging mode mounted on the 2.1 m telescope at Kitt Peak National Observatory on the nights of 2012 October 28 and 29. Continuous near-infrared (NIR) time-series photometry of the target with the Barr J-band filter was obtained on both nights: the first night was dedicated to obtain out-of-eclipse photometry, and during the second night, the secondary eclipse of the target was observed. The FLAMINGOS field-of-view is 20′× 20′, which allowed the placement of the brightest reference star in the same quadrant as the target. Our observing strategy consisted in keeping the stars in the same position on the detector with no dithering. We defocused the telescope at the beginning of the night, and then actively changed the focus to keep the detector from saturating as the temperature of the telescope and airmass of the target changed its FWHM. The 5 s science exposure times ensured that no shutter correction was needed during the reduction process of the images. Calibration frames were obtained on both nights. Dome flats were obtained at the beginning of each night with an exposure time of 8 s, and a series of night sky flat field images (150 s). Dark frames of 5, 8, and 150 s were also obtained.

We reduced the 1700 target frames in the standard manner: we first subtracted the dark current frames from both the flat-field frames and the target exposures. We then created a master flat-field frame and normalised its average counts to 1. Then, all the target frames were divided by that master flat field frame, and in a further step, bad pixels and cosmic ray hits were deleted in all the target exposures.

We extracted the stellar counts of the target star and two comparison stars (GSC 2295-0229 and GSC N322212139) by adopting aperture photometry. We chose the two comparison stars because of their brightness and because no intrinsic flux variations in our photometric data were found in the two stars. We used a self-written code for properly centring the stars in the 30 pixel wide apertures. To determine the average sky background value, we measured the count rates in a ring around the star. We paid particular attention at masking out faint stars present in those sky background rings. In the next step, we subtracted the mean sky background values from the stellar count rates and then calculated the final light curve by dividing the target count rates by the sum of the fluxes of both comparison stars. In the final step, we de-trended the light curve, thereby normalising the out-of-eclipse flux to one, and calculated the Heliocentric Julian date. The middle of the secondary eclipse was affected by the non-linearity of the detector reached while the target was a its lowest airmass of 1.0. When we exclude the affected section of the light curve the measured depth of the secondary eclipse is consistent with the depth measured including all photometric points.

Fig. 3 Follow-up J-band photometry of the secondary eclipse. We show the time-series photometry acquired at the Kitt Peak 2.1 m telescope with FLAMINGOS, which is superimposed is the best fit model (Sect. 3.5; red line) from which we derive a temperature of 3922±42 K for the M dwarf. Below, we show the residuals to the fit and the individual photometric uncertainties. On the top panel, we also show a model light curve (green, dashed line) considering a cooler temperature for the M dwarf (3350 K), as expected from theoretical models (see discussion in Sect. 4). This cooler temperature model light curve is significantly shallower than our best fit model (shown in red) and does not fit the depth of the secondary eclipse.

2.6. Radial velocities: NOT+FIES

We obtained follow-up spectroscopic observations to determine the EBLM’s orbital and stellar parameters. Our object J0113+31 was initially observed using the FIbre-fed Echelle Spectrograph (FIES) mounted on the 2.5 m Nordic Optical Telescope (NOT). In total, 15 usable spectra were obtained between 2011 August 23 and 2012 August 14. The FIES was used in medium resolution mode (R = 46 000) with interlaced ThAr calibrations, and observations were conducted using exposure times of 1200 s (S/N ~ 50) that covers the wavelength range between 3630 and 7260 Å. We used the bespoke data reduction package FIEStool¹ to extract the spectra. An IDL cross-correlation routine was used to obtain radial velocities (RVs) by fitting Gaussians to the cross-correlation functions (CCFs) of the spectral orders and taking the mean. A template spectrum was constructed by shifting and co-adding the spectra, against which the individual spectra were cross-correlated to obtain the final velocities. The template was cross-correlated with a high S/N spectrum of the Sun to obtain the absolute velocity to which the relative RVs were shifted. The RV uncertainty is given by rms(ν)/$\sqrt{\mathrm{(}}\mathit{N}\mathrm{)}$, where ν is the RV of the individual orders and N is the number of orders.

Figure 4 shows a phase-folded orbital fit to the NOT (red diamonds) and HET observations (blue asterisks). The modelling to obtain the orbital parameters is explained in Sect. 3.2, and the best-fit values are listed in Table 4.

Fig. 4 Primary radial velocity curve. We show the radial velocities describing the motion of the primary star around the system’s centre of mass. The best fit model (Sects. 3.3 to 3.5) is shown by the continuous (black) line with the systemic velocity given by the dotted line; the RV measurements acquired at the NOT are denoted by red diamonds and the HET data are shown with blue asterisks. The bottom panel shows the residuals to the fit and the uncertainties of the individual RVs.

2.7. Radial velocities: HET+HRS

Additional spectroscopic observation were obtained between October and November 2011 using the High Resolution Spectrograph (HRS; Tull et al. 1998) mounted on the 9.2 m Hobby-Eberly Telescope (Ramsey et al. 1998). The 316g5936 cross-disperser setting with a 2′′ optical fiber and a slit that provides resolving power of R ~ 30 000 were used. Each observation was a 120 s integration, which yielded S/N ~ 150 per resolution element. Science observations were bracketed before and after with a ThAr hollow-cathode lamp exposure for wavelength calibration. The HRS data was extracted, reduced, and wavelength-calibrated using a custom optimal extraction pipeline written in IDL (Bender et al. 2012).

We computed the CCFs for each epoch by the cross-correlation of fully reduced and calibrated spectra with a weighted G2 stellar template mask (Pepe et al. 2002; Baranne et al. 1996) which was created by using an NSO FTS solar atlas (Kurucz et al. 1984). Radial velocities are determined by fitting the CCF with a Gaussian (Wright et al. 2013). We find that the HRS RV measurements of σ Dra, a well-known stable star, which was observed with the HET–HRS at the same setting over a period of three months give an rms of ~56 m s^-1. We therefore quote this as our formal measurement error. We note that this error, though higher than the error bars from photon noise, is more appropriate to use, since the velocity precision is limited by instrumental effects, not by photon noise.

3. Data analysis

In this section, we describe the procedure by which we derive the physical properties of the eclipsing system. Combined analysis of the observed radial velocity curve and the photometric light curve of an eclipsing binary provides the individual masses and radii of the component stars, as well as their temperature ratio and the absolute dimensions of the system. Given that J0113+31 is a single-lined eclipsing binary (see Fig. 4), one additional parameter, the mass of the primary star in this case, is required to complete the solution (e.g., Kallrath & Milone 2009).

The analysis of a single-lined eclipsing binary is similar in many ways to a transiting planet system. However, in the case of the single lined EB, the mass of the secondary component is not negligible. Therefore, the complete solution requires an interative approach to the analysis. Derivation of the primary star radius from the eclipse light curve depends directly on the primary star mass through the semi-major axis. However, the primary star mass is derived bye comparing stellar evolutionary models using its radius, its temperature, and metallicity. Below, we describe our iterative analysis procedure. which has been applied to J0113+31 until convergence of all the parameters was achieved.

3.1. Stellar characterisation: primary star temperature and metallicity

To derive the temperature, metallicity, and gravity of J0113+31, we performed stellar characterisation using the Spectroscopy Made Easy (SME, Valenti & Piskunov 1996) spectral synthesis code. At the core of SME is a radiative transfer engine that generates synthetic spectra from a given set of stellar parameters. Wrapped around the core engine is a Levenberg-Marquardt solver that finds the set of parameters (and corresponding synthetic spectrum) that best matches observed input data in specific regions of the spectrum. The basic parameters used to define a synthetic spectrum are temperature (T_eff), gravity (log  g₁), metallicity ([M/H]), and iron abundance ([Fe/H]). To match an observed spectrum, we solve for these four parameters, and the rotational broadening (v sin i) of the star. Our implementation of SME (described in Cargile et al., in prep.) varies from the SME described in Valenti & Fischer (2005) in several ways. First, we used the ACCRE High-Performance Computing Center at Vanderbilt University to run SME 300 times on the same input data but with a range of initial conditions. By allowing SME to find a best-fit synthetic spectrum from a large distribution of initial guesses, we explore the full χ-squared space and find the optimal solution at the global minimum. In addition, we apply a line list based on Stempels et al. (2007) and Hebb et al. (2009) that includes more lines for synthesis, including the gravity sensitive Mg b triplet region. Furthermore, we use the MARCS model atmosphere grid in the radiative transfer engine, and we obtained the microturbulence (v_t) from the polynomial relation defined in Gómez Maqueo Chew et al. (2013).

We applied our SME pipeline to the two independent spectra obtained with the different telescope and instrument configurations (Sects. 2.6 and 2.7). We shifted and stacked four NOT observations obtained in the 2011 observing season to generate a single moderate S/N spectrum for stellar characterisation. We also shifted and stacked the 15 HET observations into a single, high S/N spectrum for use in an independent stellar characterisation analysis. We applied our pipeline to each combined spectra while allowing all five parameters to vary freely. The results of the SME analysis are presented in Table 1. The two datasets have different resolutions and S/N properties, but we derive the same stellar parameters from both spectra. However, we were only able to derive a robust vsini from the NOT spectrum given its higher resolution and signal-to-noise ratio, which we found to be 5.87 km s^-1. For the final stellar parameters of J0113+31, we adopt a weighted mean of the results from the two independent analyses.

The uncertainties on each set of parameters include both statistical and systematic uncertainties added in quadrature. The formal 1σ errors are based on the Δχ² statistic for five free parameters derived from our distribution of 300 final solutions. To derive the systematic uncertainties, we compared the results of four independent stellar characterisation analyses and report the mean absolute deviation of our results for the case of the planet-hosting star WASP-13 (Gómez Maqueo Chew et al. 2013). The systematic errors we derived are σ = ±48 K in T_eff, σ = ±0.07 in log  g₁, and σ = ±0.03 in [Fe/H].

3.2. Ephemeris, orbital properties, and K₁

We applied a Markov-chain Monte Carlo (MCMC) analysis to simultaneously model all our available data: the SuperWASP photometry, the higher precision, optical photometry from NITES, OED, and BYU, and the KPNO NIR secondary eclipse light curve, with the NOT and HET radial velocity measurements. A detailed description of the method is given in Collier Cameron et al. (2007) and Pollacco et al. (2008), as is typically applied for transiting planetary systems. The star–planet system represents an example of a single-lined eclipsing binary with an extreme mass ratio (q ≪ 1), and our model follows that described by Mandel & Agol (2002), which assumes that the mass of the secondary component (e.g., the planet or lower mass star) is significantly lower than the mass of the primary star (M₂ ≪ M₁). Because this assumption on the mass ratio is no longer valid in the case of J0113+31, which has an secondary mass ~0.2 M_⊙, we do not utilise this MCMC analysis to measure the absolute dimensions and masses of the binary stellar components. Instead, we use this well-tested code to measure the properties that can be directly measured from the complete dataset with robust uncertainties, namely the orbital period P_orb, the time of mid-primary eclipse T₀, the eccentricity e, the argument of periastron ω, the radial velocity semi-amplitude K₁, the centre-of-mass velocity of the system γ, and depth of the secondary eclipse (ΔF_sec). In the case of the orbital geometry, e and ω are derived from the Lagrangian elements $\sqrt{\mathit{e}}\cos \mathit{\omega }$ and $\sqrt{\mathit{e}}\sin \mathit{\omega }$ that are tightly constrained by the well-sampled RV curve and the secondary eclipse light curve. It must be noted that the derived value of e = 0.3098 is very close to the eccentricity (0.308) at which the rotational angular motion of the star at periastron is highest as compared to the mean orbital motion, as defined by Hut (1981, after Eq. (48)). Furthermore, γ was allowed to vary independently for the different RV datasets to allow for systematic offsets between the instruments/telescopes/observing runs. All the values derived from the MCMC analysis are given in Table 2.

3.3. Primary star mass determination

The spectroscopically-determined stellar parameters derived above (Sect. 3.1), and the radius of the primary star are compared with subsolar metallicity Yonsei-Yale (Y²) stellar evolutionary models (see Fig. 5; Demarque et al. 2004) to constrain the mass of the primary star. We generate a highly dense grid of mass tracks by interpolating the Y² models in temperature, radius, and metallicity. The spectroscopic properties of the primary star (e.g., T_eff and [Fe/H]) are independent from the primary mass, and as such, are kept fixed. Since our target is an eclipsing binary, we are able to measure the radius of the star from the combined analysis of the radial velocity curve and light curve (Sect. 3.5). Therefore, we compare the radius of the star instead of the luminosity or spectroscopic gravity, as is typically done for single stars, because the stellar radius is a more precise and direct observational quantity for our system.

Fig. 5 Mass tracks from Y2 models of Demarque et al. (2004) interpolated to [Fe/H] = − 0.4 (solid lines) for five different masses. The grey bands around each track indicate the range of properties at those masses, which are consistent with the uncertainty on the metallicity. The values for the radius and temperature of the primary star are plotted as the solid circle with uncertainties. Through this comparison with the stellar evolution models, we find the mass of the primary star to be 0.945 ± 0.045M⊙.

We then find the range of stellar masses that are consistent with the derived stellar parameters including their uncertainties. Figure 5 shows a modified Hertzprung-Russel diagram that compares the derived properties of the primary star to theoretical mass tracks interpolated to [Fe/H] = − 0.4 dex. The figure shows a series of mass tracks at exactly the most probable metallicity of the star, while the grey bands around each track indicate the range of properties at those masses, which are consistent with the uncertainty on the metallicity.

While the spectroscopically derived T_eff and [Fe/H] are independent of the mass and are kept constant, it is important to note that the radius derived from the eclipsing binary model is not. Thus, we derive the primary star mass from the stellar models and then perform the eclipsing binary model several times in an iterative fashion until all the quantities are consistent and agree within their uncertainties. We iterate until the primary mass derived from the theoretical models using a particular primary radius is consistent with the primary mass that is input into the EB model from which that same radius is derived. We first estimate the radius of the star based on a light curve solution using the transiting planet approximations (see Sect. 3.2). Although the planet approximations are not valid in the case of EBLMs, we use this preliminary radius (1.273±0.028 R_⊙) only as a starting point for our analysis deriving an initial primary star mass. This resulting primary mass is then fed back into a full eclipsing binary model (Sect. 3.5) to derive a new primary star radius. After two iterations of this procedure, we achieved convergence. In addition, our independently derived spectroscopic gravity was found to be consistent with the gravity derived from the mass and radius of the primary star after convergence. The final primary star mass is given in Table 4. The uncertainties in the primary mass are obtained from the range of Y² mass tracks that are consistent with the derived properties. We derive an approximate age for the system from the models of ~9.5 Gyr. This is reasonable given the low gravity and low metallicity.

However, as a sanity check on the relatively old age, we calculate the rotation period of the primary star and compare it to the expected rotation period if the system has had enough time to sychronise. Since the binary is eccentric, this is called the pseudo-psynchronous period. We find the period of the primary star to be, P_rot = 11.88 d, by using the spectroscopically-determined vsini and the primary radius and by assuming that the spin axis of the primary star is aligned with the orbital axis of the secondary star. While not always a reasonable assumption for transiting planet systems, this is reasonable for binary stars. We then calculated an estimate of the expected pseudo-synchronous period as defined by Hut (1981, Eq. (45)) for eccentric stellar binaries by assuming a constant configuration for the binary and its stellar components to be P_pseudo = 11.42 d. Considering the uncertainties in the measurements, the rotation period of the primary is consistent with being synchronized at periastron. In addition, a close binary such as J0113+31 is expected to have synchronized its rotation in ~4 Gyr, following Zahn (1977). Thus these additional constraints indicate that the age of J0113+31 is older than ~4 Gyr, which agrees with the age derived from the evolutionary models.

3.4. Mass ratio and semi-major axis

Using the orbital elements, K₁ derived in Sect. 3.2, and M₁ derived in Sect. 3.3, we solve for K₂ from M₁sin³i = 1.036149 × 10^-7(1 − e²)^3/2(K₁ + K₂)²K₂P_orb, which includes updated values for the Solar radius and for the heliocentric gravitational constant (Torres et al. 2010, and references therein). The units of the primary mass are solar mass; the K₁ is in km s^-1, and the orbital period P_orb is in days. The 1σ uncertainties are given by the extremes allowed by the errors in the parameters from which K₂ is derived. Once K₂ has been computed, the mass ratio is directly obtained from q = K₁/K₂; the semi-major axis of the orbit (in solar radii) as a function of the inclination is derived from asini = 1.976682 × 10^-2(1 − e²)^1/2(K₁ + K₂)P_orb. Table 2 contains these derived values and their 1σ uncertainties.

3.5. Eclipsing binary modelling: inclination, stellar radii and secondary temperature

Because the MCMC model assumes (a) a secondary component with a negligible mass that has (b) an opaque dark surface, and requires (c) the secondary radius be less than 10% of the stellar radius (Mandel & Agol 2002), it is necessary to use standard eclipsing binary modelling tools to primarily derive the inclination of the orbit, the stellar radii, and the temperature of the secondary component. We applied two techniques both based on the widely used Wilson-Devinney code (hereafter WD; Wilson & Devinney 1971). The WD code allows us to simultaneously fit the light and the radial velocity curves of eclipsing binaries deriving consistent parameters for all the data. For this analysis, we used the OED, BYU, and KPNO light curves, and the radial velocities of the primary. The WASP and NITES photometry were excluded from these analyses because WD calculates each light curve based on model atmospheres and observed pass-bands; this photometry was obtained with non-standard filters. We adopted the parameters derived in the previous sections and kept them fixed in the EB modelling, namely: q, e, ω, P_orb, T₀, T_eff,1, and [Fe/H]. It must be noted that the EB model depends directly on the primary mass (Sect. 3.3), as q is a needed input. Thus, the derived primary radius resulting from the EB modelling is fedback iteratively into the primary mass determination until the derived primary mass is consistent with the primary mass that is used in the EB model. Spin-orbit synchronisation was assumed for the two components. Furthermore, rotational effects for such a long period binary are not expected to be significant. The reflection albedos were fixed at the value 0.5, appropriate for stars with convective envelopes (Ruciński 1969), and the bolometric gravity darkening exponents were set to 0.4 for the primary and 0.2 for the secondary, following Claret (2000).

As shown in Table 3 from each independent EB analysis, we derived the inclination angle of the orbit (i), the fractional radii (r_j = R_j/a), and the temperature ratio (T_eff,1/T_eff,2). We then combined these two sets of results via a weighted mean to obtain their adopted values. These parameters depend solely on the light curves, and, as such, are independent from the mass ratio and asini. As mentioned above, we iterated until convergence because the primary star mass depends on the stellar radius: updating the stellar radius in the mass determination, and then deriving updated mass ratio and orbital separation, and consequently, new stellar radii.

The final EB solution indicates that the two stars are spherical, as expected for long period binaries. The primary eclipse is annular, while the secondary is total, and the totality phase lasting about 3.1 h. This provides the opportunity to constrain the effective temperature and metallicity of the primary component from high S/N spectra that are not contaminated by the secondary component. In conjunction with the results from Sects. 3.1 to 3.4, we derived the final physical properties of the orbit and stellar components, as given in Table 4.

Fig. 6 Comparison to stellar evolution models and other direct measurements. On the left, our measurements of stellar radius (top), effective temperature (middle), and luminosity (bottom) for both components of J0113+31 are shown as black diamonds against the stellar evolution models of Baraffe et al. (1998, in red) and Dotter et al. (2008, in blue), as interpolated at a metallicity of −0.4 dex for different ages (solid line for 1 Gyr; dashed line for 5 Gyr; dash-dotted line for 9.5 Gyr). The insets highlight the area around the measurements for the M dwarf for clarity. On the right, we show the same measurements of J0113+31 (cyan-filled diamonds) and compare them to other measurements of M dwarfs from double-lined eclipsing binaries (black-filled circles), from Kepler eclipsing binaries with circumbinary planets (green-filled triangles), and with interferometric radii (brown, open squares). Apart from J0113+31, the only other M dwarf in an single-lined EB with a measured temperature (KIC 1571511) is marked by the yellow diamond. See the text for the references of the measurements included in this figure. The model evolutionary tracks are the same as on the left panel, except for the black lines, which show the models of Baraffe et al. (1998, solid) and Dotter et al. (2008, dashed) for solar metallicity at 9.5 Gyr for comparison. A colour version can be found in the online version of the paper.

3.5.1. Wilson-Devinney (2010) modelling

The first EB modelling of the light and radial velocity curves was performed using the 2010 version of the WD code. The main parameters adjusted were the orbital inclination (i), the pseudo-potentials (Ω₁ and Ω₂), the temperature ratio (T_eff,1/T_eff,2), the semi-major axis (a), the systemic velocity (γ), the primary luminosity (L₁) for each bandpass (i.e., the light ratio), and the time of periastron passage (to account for variations of the times of eclipses). Initially, we tried to also fit e and ω; however, the solution diverged given that the parameters are highly correlated. Based on the method of multiple subsets described by Wilson & Biermann (1976), we fit e and ω independently of the other parameters and obtained values for e and ω that were consistent with those derived in Sect. 3.2. Emergent intensities used in the program were taken from model atmospheres described by van Hamme & Wilson (2003), and limb-darkening coefficients were computed from van Hamme (1993) as implemented in the WD code. The coefficients were dynamically adjusted, according to the current effective temperatures and surface gravities of the stars at each iteration.

In the case of this first method, the convergence in the final fit was considered to have been achieved when the corrections to the elements were smaller than the internal errors in three consecutive iterations. This procedure was repeated five times, and the solution with the smallest χ² value was chosen as our best solution. Observational weights for each light curve were adjusted, according to their residuals using a preliminary fit and scaled to get a χ² ~ 1 for the light curves. The observational weights for the radial velocities were also scaled accordingly. To provide more realistic uncertainties for the geometric and radiative parameters than the internal errors estimated by WD, we continued the iterations in our adopted fit beyond convergence for another 200 steps, and we examined the scatter of those 200 solutions. We adopted the larger of the estimates for the parameters. The third light was also tested for but was always found to either converge toward negative (unphysical) values, or be roughly consistent with zero. The depth of the secondary eclipse measured from this WD analysis is 0.00741 ± 0.00001, which is consistent with both the MCMC and phoebe measurements.

3.5.2. PHOEBE modelling

The second EB analysis was performed using the WD-based code phoebe (Prša & Zwitter 2005). Firstly, we made use of the scripter capability of phoebe and cluster computing to sample a large parameter space. The parameters sampled were the inclination of the orbit (i) and the stellar radii via the potentials (Ω₁, Ω₂), as defined for eccentric orbits by Wilson (1979).

Utilising the Vanderbilt University Advanced Computing Center for Research and Education (ACCRE), we randomly selected values between 15.5 ≤ Ω₁< 23.5, 22.0 ≤ Ω₂< 32.0, and 88.2 ≤ i< 90.0. For each combination of i, Ω₁, and Ω₂, we calculated a from asini (Sect. 3.4) and the chosen i. To make this computationally expensive process more efficient while mapping the space around the minimum sufficiently well after the first 10 000 points, we narrowed the sampled parameter range around the solution with the lowest χ². The total number of different parameter combinations sampled was 36972. The calculation of χ² includes the radial velocity curve and the follow-up light curves. It must be noted that the sampled parameters (i, Ω₁, and Ω₂) depend solely on the light curves. However, the stellar radii depend on the absolute dimension of the binary orbit determined by a (semi-major axis) which is derived from asini (see Sect. 3.4) that depends on the RV curve and the light curves via i. This allowed us to create a multi-dimensional χ² map and ensure that our best solution (i.e., with the lowest χ² value) corresponds to the global minimum. Furthermore, we are able to define the 1σ uncertainties around the best solution from the χ² map and to assess the correlations between the parameters (e.g., Gómez Maqueo Chew et al. 2012).

To derive the temperature ratio (T_eff,1/T_eff,2), we included the secondary eclipse light curve (Sect. 2.5) in the dataset modelled with phoebe. We fitted primary and secondary eclipse light curves simultaneously to derive the best solution. We obtained fractional radii and an inclination that agree within 1σ with the result from the best ACCRE solution and a temperature ratio of 1.509±0.012, which combined with the T_eff,1 derived in Sect. 3.1 renders T_eff,2 = 3950 K. The depth of the secondary eclipse derived from our phoebe analysis is 0.0074 in flux, which agrees with the measurements from the WD and MCMC analyses.

4. Discussion

Table 4 presents our final physical properties for the system and its stellar components which in Fig. 6 are compared to two sets of stellar evolution models: BCAH98 (Baraffe et al. 1998) and DSEP (Dotter et al. 2008). Both sets of evolutionary models have been interpolated to the metallicity of the EBLM (−0.4 dex). It must be noted that the BCAH98 models for the estimated age of our system (9.5 Gyr) only go up to masses of 0.9 M_⊙, because stars of higher mass at this metallicity beyond this age have evolved off the main sequence. The differences between both set of models could be due to the different mixing-length parameters used in the readily available theoretical calculations (α_BCAH98 = 1.0; α_DSEP = 1.938). In general, a more efficient convection (i.e., higher mixing-length α) causes the overall radius to decrease and the temperature to increase to maintain the same luminosity.

As shown on the left side of Fig. 6, the 9.5 Gyr DSEP model reproduces our measurements for the primary radius, temperature, and thus luminosity for the derived mass of the primary very well. Given that the primary mass and the system’s age are derived directly from the Y² models (Sect. 3.3), the intersection of the primary radius, temperature, and luminosity with the 9.5 Gyr DSEP isochrone denotes good agreement between these models in this mass regime. More specifically, our measurements of the primary radius and luminosity fall within their uncertainties on the 9.5 Gyr DSEP isochrones in the mass-radius and mass-luminosity planes, respectively, and are not consistent with the ones at 1 and 5 Gyr. In the temperature-mass plane, our measurement of T_eff,1 is consistent within uncertainties with all three DSEP isochrones and with the 5 Gyr BCAH98 isochrone.

While both sets of models are able to reproduce the radius of the M dwarf given its measured mass, its measured T_eff,2 is significantly hotter than that predicted by the theoretical isochrones. Consequently, the observed luminosity is also well above the predicted one. Due to the slow evolution of low-mass stars, our measurement of the secondary radius is consistent with both sets of models from 1 to 9.5 Gyrs. At these low masses, the BCAH98 and the DSEP models are indistinguishable in the mass-radius plane. Both sets of models predict a star with the mass of the secondary to have a temperature of ~3350 K, which is ~600 K cooler than our measured T_eff,2.

Figure 3 compares our best fit light curve model for a 3922 K M dwarf and that of a 3350 K, as predicted by stellar evolution models. Our data is not consistent with the cooler temperature model light curve. Below, we discuss different phenomena/effects that could impact our measurement of the secondary effective temperature by either affecting our measurement or by providing a source of heating, that is additional to fusion in the M-dwarf interior. Any additional heating mechanism must contribute ~10³¹ erg/s, which is about 70% of the energy produced by fusion for the 3350 K version of this star.

• Treatment of model atmospheres in EB modelling: using the most up-to-date model atmospheres for cool stars (PHOENIX; Husser et al. 2013), we investigate and verify that our measurement of the T_eff,2 from the EB modelling is consistent with the latest stellar atmospheres, since the underlying model stellar atmospheres used by both PHOEBE and WD are those of Kurucz (van Hamme & Wilson 2003, and references therein). These PHOENIX atmospheres extend to very low effective temperatures (2300 K) and include the wavelength range of our secondary eclipse (~1.1−1.4 μm). We were not able to use the synthetic spectra of the AMBRE project which based on the MARCS model atmospheres, because their cut-off wavelength is 1.2 μm (de Laverny et al. 2012). Because of the discrete grid points calculated by Husser et al. (2013) for this test we chose the closest values to our derived physical properties. Thus, we utilised [Fe/H] = − 0.5 dex for both components (with no alpha element enhancement; α = 0.0 dex), T_eff,1 = 6000±100 K and log  g₁ = 4.0 dex for the primary, and log  g₂ = 5.0 dex for the secondary. To compare our measurement of the eclipse depth to that expected by model atmospheres, we integrated the stellar atmospheres over the transmission function of the Barr J-band filter of our NIR observations and scaled the model atmospheres by the measured stellar radii. As shown in Fig. 7, we tested temperatures for T_eff,2 between 2300 and 5000 K in 100 K intervals and confirmed that the measured secondary eclipse depth (0.00737±0.00024 in flux) is consistent with the PHOENIX model atmospheres for a secondary star of ~3900 K.

• Determination of primary temperature and metallicity: given that T_eff,1/T_eff,2 is measured from the relative depth of the eclipses and the spectroscopic T_eff,1 and that we adopt the metallicity of the primary as that of the M dwarf, we discuss the possibility that the stellar characterisation of the primary star is not accurate. If T_eff,1 were cooler, then T_eff,2 would so consequently. However, with the measured temperature ratio, the primary temperature would need to be ~900 K lower to reconcile T_eff,2 with evolutionary models. Similarly, the metallicity of the primary star (and thus of the M dwarf) would need to be between −1.5 and −2.0 dex for our measured T_eff,2 to be consistent. Both of the effects (i.e., cooler T_eff,1 and/or a more metal poor system) would imply the primary star to be lower in mass than what we have estimated in Sect. 3.3, and thus the derived T_eff,2 in turn would still likely be higher than expected. The observed spectra are not consistent wtih such cool temperatures nor such a low metallicity for the primary. Furthermore, our spectroscopic analysis of the two available spectra were done independently and render results within 1σ of each other (see Table 1 and Sect. 3.1). Although some studies have found that SME measurements show strong correlations between effective temperature, metallicity, and log  g₁ (e.g., Torres et al. 2012), our primary star is cooler than 6000 K, above which these systematics start to appear. Because our spectroscopically-determined log  g₁ is consistent with the log  g₁ that is derived from the EB modelling, this indicates that the spectroscopically-derived temperature and metallicity are not subject to the systematics problems outlined in Torres et al. (2012). Thus, we conclude that the T_eff,2 discrepancy is not likely due to an inaccurate stellar characterisation of the primary star and that the system is not sufficiently metal-poor to reconcile T_eff,2 with the evolutionary models.

• Alpha element enhancement: our stellar characterisation of the primary star and the model atmospheres used in the EB modelling assume solar abundances of the alpha elements. Thus, as done similarly above, using the latest PHOENIX model atmospheres (Husser et al. 2013), we explored the possibility of alpha element enhancement as the cause of the deep secondary eclipse (i.e., instead of due to a high effective temperature). As shown in Fig. 7, we tested the upper and lower limit of α-element enchancement considered by the stellar atmospheres, −0.2 and + 1.2 dex, respectively. We find that even an enhancement of + 1.2 dex is not sufficient to reconcile the observed depth of the secondary eclipse with that of a secondary M dwarf with a temperature as expected from theoretical evolutionary models (~3350 K; Baraffe et al. 1998; Dotter et al. 2008).

  Fig. 7 Expected secondary eclipse depth based on the latest PHOENIX stellar atmospheres from Husser et al. (2013). The magenta-lined area marks the measured depth of the secondary eclipse from our NIR light curve, 0.00737±0.00024 in flux. The red crosses denote the expected secondary eclipse depth when assuming [Fe/H] = −0.5 dex, Teff,1 = 6000 K, log  g1 = 4.00 dex, log  g2 = 5.00 dex, and no alpha element enhancement, and the grey area denotes an uncertainty in the primary temperature of 100 K. The triangles are similar to the crosses, except that they represent different alpha enhancement levels for both components: the blue triangles are for α = + 1.2 dex, and the green triangles are for α = −0.2 dex. The red dotted line is the expected eclipse when using black bodies instead of the model atmospheres. Even the highest alpha element enhancement considered by the PHOENIX model atmospheres is not enough to reconcile the simulated secondary eclipse depth for the Teff,2 that is expected from evolutionary models with our measurement of the secondary eclipse.

• Contamination by tertiary component or background/foreground star: an unresolved component would affect the measured relative depth of the eclipses and, thus, the derived T_eff,1/T_eff,2. If there were an unresolved low-mass star, then the primary eclipse depth would not be severly affected by the third light contamination due to the high luminosity ratio, while the secondary eclipse in the NIR would be shallower. This effect would cause the measured T_eff,2 to be cooler than expected (i.e., not hotter as we observe). A more exotic blend, such as that caused by a white dwarf, would make the primary eclipse shallower in the I-band and would not significantly affect secondary eclipse depth in the NIR, causing the secondary star to appear hotter. The radial velocity measurements would not be affected in the case of contamination of a background/foreground star. In the case of a physically-associated tertiary component, either a low-mass star or a white dwarf in a long enough period (relative to our RV sensitivity and timespan of our observations) would not cause a significant RV motion. However, the low-mass component (~0.14 M_⊙) of the only other EBLM-type object with a measurement of the M-dwarf temperature was discovered from its Kepler light curves in which the secondary eclipse was evident (KIC 1571511; Ofir et al. 2012). This object was also found to be much hotter than expected by the DSEP models for its measured mass (ΔT_eff ~ 900 K). In this case, either blend scenario would dilute both eclipse depths equally because they are observed in the same filter; the temperature ratio for this object would not be affected by either a low-mass star blend or by an exotic blue blend. In the case of J0113+31, contamination by a low-mass star is not able to explain the observed hot T_eff,2, while a blend with a white dwarf could. However, the M-dwarf component of KIC 1571511 is also found to be significantly hotter than expected, and in that case, no blend scenario is able to explain it. Thus, we conclude that a blend scenario is not likely the cause of the hot secondary of J0113+31; more measurements of the temperature of M dwarfs would be able to confirm this high temperature trend.

• Magnetically-induced starspots: low-mass M dwarfs are thought to be very spotted; magnetically-induced spots are expected to be cooler than the photosphere, reducing the average temperature while increasing the stellar radius (e.g., Chabrier et al. 2007). Neither of these two effects of starspots would explain the measured, lower-than-expected temperature ratio (i.e., the high secondary temperature), given that T_eff,1/T_eff,2 is independent from the stellar radii, and a cooler effective temperature of the secondary would lead to an increase of T_eff,1/T_eff,2. Although the presence of hot spots (plages) due to activity is not certain on M-dwarfs, several studies have shown that some eclipsing binary light curves are better reproduced with hot spots (e.g., López-Morales & Ribas 2005; Morales et al. 2009). However, either a cool spot on the primary component or a bright plage on the secondary would need to be too large to explain the temperature ratio of the system or the secondary eclipse depth. On the contrary, the light curves do not show a significant modulation due to spots. Therefore, we conclude that starspots are not likely the cause of the lower-than-expected temperature ratio.

• Hot spot due to irradiation from the primary star: it is also not likely that the deeper than expected secondary eclipse is due to a hot spot on the secondary caused by the irradiation from the primary component. Because of the eccentric orbit, the components are not rotating synchronously with the orbital motion, and thus there is not one side of the secondary that is always being irradiated by the primary. From the measurement of vsini and the stellar radius, the rotation period of the primary component is consistent with the pseudo-synchronous rotation (~12 d; Sect. 3.3). However, we have considered the most conservative case in which both components are rotating synchronously and the mutual irradiation effects between the components are taken into account in the EB modelling. There is no evidence of a hot spot in the observed and model light curves. Furthermore, when modeling the light curve with a range in albedo from 0.0 to 1.0 for either star, the depth of the secondary eclipse remains constant to within 0.001% in flux. Thus, it is unlikely that a hot spot due to irradiation from the primary star can be the cause of the larger than expected depth of the secondary eclipse.

• Residual heat from formation: it is unlikely that the hot M dwarf temperature we observe is due to residual heat from its formation (i.e., the M dwarf is younger than we estimate). Firstly, there are no youth signatures in the spectrum of the primary, and close binaries are generally thought to be formed at the same time (e.g., Prato et al. 2003). Moreover, M dwarfs of 0.2 M_⊙ are thought to stay at roughly a constant temperature from the first several Myrs until the end of their main-sequence life (e.g., Baraffe et al. 2002), and thus, it is never expected during its evolution before and through the main sequence to have such a high temperature.

• Mass transfer and/or accretion: the stellar components of J0113+31 are well detached and inside their Roche lobes. They are not transfering mass and, thus, are not interacting. There is also no signature in the light curves or spectra of the presence of a circumbinary or circumstellar disk from which the stars could be accreting new material. Right after formation, the radii of the stars could have been large enough to interact, and they could have surrounding material from which to accrete. However, the effects of episodic accretion and/or mass transfer at these young ages disappears after a few Myr (Baraffe et al. 2009). Thus, given J0113+31’s old age and that mass transfer and accretion are not currently occuring, accretion and/or mass transfer are not likely the cause of the temperature difference.

• Tidal heating: tidal heating cannot explain the apparent temperature discrepancy. Given the orbital period and significant eccentricity, we examined the possibility that tidal heating could contribute extra energy to the M dwarf and raise its temperature to the observed value. We calculated the present rate of tidal heating with the “equilibrium tide model,” which treats the star as a deformed spheroid and determines dissipation by just a single parameter, such as the tidal Q (e.g. Zahn 1975; Ferraz-Mello et al. 2008; Leconte et al. 2010). In particular, we use the “constant-phase-lag” (CPL) model (Greenberg 2009) as described in Barnes et al. (2013), (see also Heller et al. 2011) and refer the reader to the former reference for a complete description of this model. In the CPL framework, the rate of dissipation is inversely proportional to Q, which probably has a value in the range 10⁴–10⁹ in the case of stars (e.g. Mardling & Lin 2002; Jackson et al. 2009; Matsumura et al. 2010; Adams et al. 2011). Tidal heating is also a function of the rotation rate, which is unknown for the M dwarf, and, thus, we treat it as a free parameter. We considered a range of tidal Qs from 10⁴ and 10⁷ and rotation periods from 1 to 100 days and calculated the tidal heating rate. The largest tidal heating occurs at the smallest period and Q and is about 4 × 10²⁸ erg/s, which is about 2.5 orders of magnitude too low. We also examined non-zero obliquities and found that heating from obliquity is negligible. We conclude that tidal heating cannot explain the apparent temperature discrepancy.

• Missing physics in atmosphere and/or evolutionary models: although stellar atmosphere and evolutionary models are able to reproduce some of the direct measurements of M dwarfs (e.g., Fig. 6), as discussed below, the comparison of other mass, radius, metallicity, and temperature measurements for M dwarfs are marred with different assumptions for the distinct kinds of systems (EBLMs, M+M EBs, or single M dwarfs). There is always the possibility that either or both the atmospheric and evolutionary models are missing relevant physics for systems, such as J0113+31 (e.g., in unequal-mass binaries, at low metallicities, and/or with a low-mass M dwarf). A thorough assessment of the models is beyond the scope of this paper.

On the right side of Fig. 6, we compare the known measurements for stars with masses <1.0 M_⊙ with the BCAH98 and DSEP evolutionary models. Typically, these comparisons are done against the double-lined EBs from which masses, radii, and temperatures are measured, and thus, luminosities are derived (Bass et al. 2012; Birkby et al. 2012; Blake et al. 2008; Brogaard et al. 2011; Carter et al. 2011; Coughlin 2012; Creevey et al. 2005; Fleming et al. 2011; Hartman et al. 2011; Hebb et al. 2006; Hełminiak & Konacki 2011; Hełminiak et al. 2011, 2012; Irwin et al. 2009, 2011; Kraus et al. 2011; Liakos et al. 2011; López-Morales 2007; Morales et al. 2009; Rozyczka et al. 2009; Thompson et al. 2010; Torres et al. 2010; Torres & Ribas 2002; Vaccaro et al. 2007; Young et al. 2006). We also compare the Kepler EBs with circumbinary planets (Doyle et al. 2011; Orosz et al. 2012; Welsh et al. 2012; Schwamb et al. 2013), and KIC 1571511 (Ofir et al. 2012), which is the only other single-lined EB with measurements of the mass, radius and temperature for the M-dwarf companion (apart from J0113+31). The direct radius measurements from EBs typically find the stars to be larger than predicted by evolutionary models (e.g., López-Morales 2007; Morales et al. 2008; Torres et al. 2010). Direct measurements of the radii via interferometry for nearby stars are also used for comparison (Ségransan et al. 2003; Berger et al. 2006; Demory et al. 2009; Boyajian et al. 2012) and exhibit the same trend. It is thought this is due to magnetic activity and/or reduction of convection efficiency (e.g., Morales et al. 2010; Mullan & MacDonald 2001). The increase in radius is compensated with a decrease in effective temperature to maintain the stellar luminosity. The spread in the measurements of radius and temperature as a function of stellar mass in Fig. 6 is not only due to the metallicity but age becomes important for more massive stars where the evolution is much faster than for the lower mass stars (≲0.6−0.7 M_⊙).

In the case of M dwarfs, their temperatures as a function of mass and metallicity remain uncertain. Most M dwarf temperatures are derived from the analysis of spectra line indices and/or broadband SED modelling (e.g., Rojas-Ayala et al. 2012) for nearby, single M dwarfs with interferometric radii (see brown open squares in Fig. 6). However, for single stars, it is not possible to obtain a dynamical measurement of their mass, so their mass estimates rely on comparing observed properties to empirical and/or theoretical scales (e.g., Rajpurohit et al. 2013). In the case of the EBs that exhibit primary and secondary eclipses, as shown in this paper, the temperature is typically measured from the temperature ratio (i.e., from the relative depth of the eclipses) in conjunction with a measure of one of the individual temperatures (e.g., from spectral type or stellar characterisation of the primary). Typically, for double-lined EBs composed of two M dwarfs, the temperature ratio is measured from the light curve(s) and the integrated temperature from colour indices or SED models (e.g., Torres & Ribas 2002). Getting the system’s metallicity is particularly challenging because of the two sets of complex spectral features of the two unresolved M dwarfs. Although it is this double-lined nature that allows the direct mass determination for these M+M EBs, it hinders the spectroscopic determination of individual temperatures and metallicity. In the case of EBLMs, the method to derive the temperatures is the same as for EBs. However, it is more precise because the stellar characterisation of our primary star is well-understood in the solar-type regime and because of the high luminosity contrast between the components allows for high-quality (single-lined) spectra of the primary to be acquired. Systems in the EBLM sample, like J0113+31, in which the primary is a solar-type star and the secondary is an M dwarf, will provide a large number of measurements of the mass, radius, temperature, metallicity, and age for M dwarfs.

5. Summary

In this second paper of the EBLM Project, we derive the orbital parameters of J0113+31and the fundamental properties of its stellar components. We present the first full analysis of an EBLM in our sample of ~150 systems discovered from their WASP light curves, thereby defining the project’s methodology. Our object J0113+31 is an old and metal-poor system, as determined by the large radius and the spectrum of the solar-type primary star with an eccentric and long-period orbit. The secondary radius of the low-mass M dwarf is consistent with stellar evolution models for its given mass, but its temperature is measured to be ~600 K hotter than expected. We discard different sources of possible error in our measurement of the M dwarf temperature, including the treatment of model atmospheres in the eclipsing binary model, the stellar characterisation of the primary, and contamination by an unresolved star. We also discuss different physical processes that could have an impact on the M dwarf by affecting its effective temperature (e.g., hot or cold spots, younger age) or by providing an additional source of energy (e.g., tidal heating, mass transfer, accretion). These scenarios are not able (or are not very likely) to account for such a large difference between the temperature expected by the stellar evolutionary models and the one measured from the secondary eclipse. Until the relationship between the mass, radius, and temperature for M dwarfs as a function of their metallicity is well defined, caution must be taken when deriving M dwarf masses from luminosities, temperatures, and/or colours. The EBLM Project will be able to provide these empirical constraints which will be crucial, for example, when deriving physical properties of planets around M dwarfs discovered by TESS and/or Gaia.

Acknowledgments

The authors would like to thank Isabelle Baraffe, Ignasi Ribas and Barry Smalley for helpful discussions. The research leading to these results has received funding from the European Community’s Seventh Framework Programme (FP7/2007-2013) under grant agreement number RG226604 (OPTICON). L.H.H. acknowledges funding support from NSF grant, NSF AST-1009810. A.H.M.J. Triaud received funding from the Swiss National Science Foundation in the form of an Advanced Mobility post-doctoral fellowship (P300P2-147773). E.G.M. was supported by the Spanish MINECO project AYA2012-36666 with FEDER support. R.D. and S.M. acknowledge funding support from the Center for Exoplanets and Habitable Worlds. The Center for Exoplanets and Habitable Worlds is supported by the Pennsylvania State University, the Eberly College of Science, and the Pennsylvania Space Grant Consortium. The Hobby-Eberly Telescope (HET) is a joint project of the University of Texas at Austin, the Pennsylvania State University, Stanford University, Ludwig-Maximilians Universitat München, and Georg-August-Universität Göttingen. The HET is named in honour of its principal benefactors, William P. Hobby and Robert E. Eberly. This research is based on observations made with the Nordic Optical Telescope, operated by the Nordic Optical Telescope Scientific Association at the Observatorio del Roque de los Muchachos, La Palma, Spain, of the Instituto de Astrofisica de Canarias, as well as from Kitt Peak National Observatory, National Optical Astronomy Observatory, which is operated by the Association of Universities for Research in Astronomy (AURA) under cooperative agreement with the National Science Foundation. FLAMINGOS was designed and constructed by the IR instrumentation group (PI: Elston) at the University of Florida, Department of Astronomy, with support from NSF grant AST97-31180 and Kitt Peak National Observatory. The BYU West Mountain Observatory 0.91 m telescope was supported by NSF grant AST0618209 during the time these observations were secured. This work was conducted in part using the resources of the Advanced Computing Center for Research and Education at Vanderbilt University, Nashville, TN.

References

All Tables

All Figures

Fig. 1 SuperWASP light curve. We show the time-series photometry that led to the discovery of J0113+31, as part of the SuperWASP transiting planet survey, overplotted on the final model light curve from the EB modelling (Sect. 3.5) in the V-band, which is used as a proxy for the non-standard SuperWASP filter. The primary eclipse of the M dwarf blocking the light of the more massive star has a transit-like shape (shown in inset). The error bar to the bottom left represents the median uncertainty of the SuperWASP photometric points, which are excluded for clarity.

In the text

Fig. 2 Follow-up photometry of the primary eclipse. On the top, we show the time series photometric points acquired at the OED (Sect. 2.3), phased to the ephemeris from Sect. 3.2, and overplotted with the best fit EB model described in Sect. 3.5 (red line). We also show directly below the residuals to the fit and the uncertainty of the individual photometric data points. Similarly in the middle panel, we show the data from BYU and the residuals to the fit. In the bottom panel, we show the partial primary eclipse data acquired with NITES. Overplotted in grey is a model light curve in the R-band, shown only for guidance and comparison to the upper panels. The NITES data were taken without a filter, and as such, we only utilise the NITES phtometry to derive the ephemeris, excluding it from the EB analysis (Sect. 3.5).

In the text

Fig. 3 Follow-up J-band photometry of the secondary eclipse. We show the time-series photometry acquired at the Kitt Peak 2.1 m telescope with FLAMINGOS, which is superimposed is the best fit model (Sect. 3.5; red line) from which we derive a temperature of 3922±42 K for the M dwarf. Below, we show the residuals to the fit and the individual photometric uncertainties. On the top panel, we also show a model light curve (green, dashed line) considering a cooler temperature for the M dwarf (3350 K), as expected from theoretical models (see discussion in Sect. 4). This cooler temperature model light curve is significantly shallower than our best fit model (shown in red) and does not fit the depth of the secondary eclipse.

In the text

Fig. 4 Primary radial velocity curve. We show the radial velocities describing the motion of the primary star around the system’s centre of mass. The best fit model (Sects. 3.3 to 3.5) is shown by the continuous (black) line with the systemic velocity given by the dotted line; the RV measurements acquired at the NOT are denoted by red diamonds and the HET data are shown with blue asterisks. The bottom panel shows the residuals to the fit and the uncertainties of the individual RVs.

In the text

Fig. 5 Mass tracks from Y2 models of Demarque et al. (2004) interpolated to [Fe/H] = − 0.4 (solid lines) for five different masses. The grey bands around each track indicate the range of properties at those masses, which are consistent with the uncertainty on the metallicity. The values for the radius and temperature of the primary star are plotted as the solid circle with uncertainties. Through this comparison with the stellar evolution models, we find the mass of the primary star to be 0.945 ± 0.045M⊙.

In the text

Fig. 6 Comparison to stellar evolution models and other direct measurements. On the left, our measurements of stellar radius (top), effective temperature (middle), and luminosity (bottom) for both components of J0113+31 are shown as black diamonds against the stellar evolution models of Baraffe et al. (1998, in red) and Dotter et al. (2008, in blue), as interpolated at a metallicity of −0.4 dex for different ages (solid line for 1 Gyr; dashed line for 5 Gyr; dash-dotted line for 9.5 Gyr). The insets highlight the area around the measurements for the M dwarf for clarity. On the right, we show the same measurements of J0113+31 (cyan-filled diamonds) and compare them to other measurements of M dwarfs from double-lined eclipsing binaries (black-filled circles), from Kepler eclipsing binaries with circumbinary planets (green-filled triangles), and with interferometric radii (brown, open squares). Apart from J0113+31, the only other M dwarf in an single-lined EB with a measured temperature (KIC 1571511) is marked by the yellow diamond. See the text for the references of the measurements included in this figure. The model evolutionary tracks are the same as on the left panel, except for the black lines, which show the models of Baraffe et al. (1998, solid) and Dotter et al. (2008, dashed) for solar metallicity at 9.5 Gyr for comparison. A colour version can be found in the online version of the paper.

In the text

Fig. 7 Expected secondary eclipse depth based on the latest PHOENIX stellar atmospheres from Husser et al. (2013). The magenta-lined area marks the measured depth of the secondary eclipse from our NIR light curve, 0.00737±0.00024 in flux. The red crosses denote the expected secondary eclipse depth when assuming [Fe/H] = −0.5 dex, Teff,1 = 6000 K, log  g1 = 4.00 dex, log  g2 = 5.00 dex, and no alpha element enhancement, and the grey area denotes an uncertainty in the primary temperature of 100 K. The triangles are similar to the crosses, except that they represent different alpha enhancement levels for both components: the blue triangles are for α = + 1.2 dex, and the green triangles are for α = −0.2 dex. The red dotted line is the expected eclipse when using black bodies instead of the model atmospheres. Even the highest alpha element enhancement considered by the PHOENIX model atmospheres is not enough to reconcile the simulated secondary eclipse depth for the Teff,2 that is expected from evolutionary models with our measurement of the secondary eclipse.

In the text