A Strict Test of Stellar Evolution Models: The Absolute Dimensions of Massive Benchmark Eclipsing Binary V578 Mon
We determine the absolute dimensions of the eclipsing binary V578 Mon, a detached system of two early B-type stars (B0V + B1V, P2.40848 d) in the star-forming region NGC 2244 of the Rosette Nebula. From the light curve analysis of 40 yr of photometry and the analysis of hermes spectra, we find radii of Rsun and Rsun, and temperatures of K and K respectively. We find that our disentangled component spectra for V578 Mon agree well previous spectral disentangling from the literature. We also reconfirm the previous spectroscopic orbit of V578 Mon finding that masses of Msun and Msun are fully compatible with the new analysis. We compare the absolute dimensions to the rotating models of the Geneva and Utrecht groups and the models of Granada group. We find all three sets of models marginally reproduce the absolute dimensions of both stars with a common age within uncertainty for gravity-effective temperature isochrones. However - there are some apparent age discrepancies for the corresponding mass-radius isochrones. Models with larger convective overshoot worked best. Combined with our previously determined apsidal motion of deg cycle, we compute the internal structure constants (tidal Love number) for the newtonian and general relativistic contribution to the apsidal motion, and respectively. We find the relativistic contribution to the apsidal motion of be small . We find that the prediction of of the Granada models fully agrees with our observed .
Detached eclipsing binary stars (dEBs) provide accurate observed stellar masses, radii, effective temperatures, and rotational velocities. See a recent review by Torres et al. (2010) for a discussion of 94 dEBs with accurate masses and radii used to test stellar evolution models. There are only nine total massive dEBs, or equivalently 18 stars whose physical parameters have been determined with an accuracy of better than 3%, making V578 Mon one of only nine EBs with and with sufficient accuracy to be included in the Torres et al. (2010) compilation of benchmark-grade EBs. Figure 1 demonstrates the upper main sequence of all dEBs with and masses and radii determined to 3% (adapted from Torres et al. (2010)). V578 Mon is therefore a benchmark system for testing stellar evolution models of newly formed massive stars. The accurate absolute dimensions of eclipsing binary stars provide a unique opportunity to test stellar evolution models in two ways: the “isochrone test” and the“apsidal motion test”.
The “isochrone test” of stellar evolution models requires that the ages of both components of a dEB predicted from separate stellar evolution tracks to be the same within uncertainty of the absolute dimensions . For the “isochrone test”, we assume that both components of the dEB formed together in the same initial gas cloud. Therefore, both components of a dEB are assumed to arrive at the zero age main sequence (ZAMS) at nearly the same time. Furthermore, their initial chemical compositions must be the same. Finally, we assume that each component of the binary evolves in isolation, where the effects of the companion star on the evolution is small or negligible.
The “isochrone test” is strongest for eclipsing binaries with low mass ratios . For dEBs where component masses , both stars will evolve on the same evolutionary track. This does not allow for strict tests of stellar evolution models unless the chemical composition or effective temperature of the stars is known. Stellar evolution models will predict two stars of the same mass and composition to have the same age. Conversely, the larger the difference in initial mass between the components of the binary star, the larger the difference in main sequence (MS) lifetimes of the two stars. Therefore, the stellar models must have the accurate input physics to correctly predict how quickly stars of different mass evolve relative to each other. The correct input physics in turn yields correct predictions of the observed absolute dimensions of the detached eclipsing binary.
Detached eclipsing binary stars with apsidal motion (precession of the argument of periastron) also allow for the “apsidal motion test” of the stellar internal structure (Claret & Giménez 2010). Physically, the observed apsidal motion rate in an eclipsing binary is a result of the tidal forces of each star on each other. In turn, this tidal force is linked to the internal structure of each star, the star’s separation, their mass ratio and their radii and . The internal structure is quantified by the constant which is the logarithm of twice the tidal Love number (Kramm et al. 2011). The “apsidal motion” test compares the theoretical internal structure constant to the observed internal structure constant . The observed internal structure constant is a function the observed absolute dimensions and apsidal motion of the eclipsing binary. The observed internal structure constant is very sensitive to the radii () - therefore, this test can only be performed with accurate stellar radii. However, including this study of massive dEB V578 Mon, there are only five massive, eccentric eclipsing binaries available for these tests of internal structure (Claret & Giménez 2010).
Here we combine the previous determination of and from Garcia et al. (2011) with a reanalysis of 40 years worth of photometry to re-determine the fundamental properties of V578 Mon. We include the photometry used by the previous light curve analysis (Hensberge et al. 2000). We compare the masses, temperatures and radii of V578 Mon with rotating high mass stellar evolution models by Granada (Claret 2004, 2006), Geneva (Georgy et al. 2013; Ekström et al. 2012), and Utrecht (Brott et al. 2011) groups. We also compare the observed internal structure constant with theoretical using the methods of Claret & Giménez (2010).
2 The Eclipsing Binary V578 Mon in NGC 2244
The photometric variability of the bright (V=8.5), 2.408 day period, eccentric, massive detached eclipsing binary (dEB) V578 Mon (HDE 259135, BD1299), comprising a B1V type primary star and a B2V type secondary star was first identified in the study by Heiser (1977) of NGC 2244 within the Rosette Nebula (NGC 2237, NGC 2246). The identifications, locations and photometric parameters for V578 Mon are listed in Table 1. The absolute dimensions of V578 Mon have been determined from three seasons of Strömgren photometry and one season of radial-velocity data by Hensberge et al. (2000). An analysis of the metallicity and evolutionary status of V578 Mon was undertaken by Pavlovski & Hensberge (2005) and Hensberge et al. (2000). The masses and radii of V578 Mon determined from these data are and , and and for the primary and secondary respectively (Hensberge et al. 2000). V578 Mon was included in the list of 94 detached eclipsing binaries with masses and radii accurate to by Torres et al. (2010). The radii for V578 Mon listed in Torres et al. (2010) were found to be incorrect by Garcia et al. (2013) given system’s eccentric orbit and asynchronous rotation. The apsidal motion and a new eccentricity were determined in Garcia et al. (2011). V578 Mon was observed by MOST (Pribulla et al. 2010).
Given the inclination of V578 Mon, its eclipses are partial, meaning that neither star is fully out of view of Earth. Partial eclipses can translate into a degeneracy between the radii, preventing the component radii and from being individually measured. However, V578 Mon also has an eccentric orbit, meaning that the eclipse durations are not equal, which helps breaks this degeneracy and allows the radii to be determined separately. V578 Mon is observed to not have tidally locked yet. The system has a low mass ratio = as compared to similar systems with well-determined absolute parameters such as V1034 Sco, V478 Cyg, AH Cep, V453 Cyg, and CW Cep (Bouzid et al. 2005; Popper & Etzel 1981; Popper & Hill 1991; Bell et al. 1986; Holmgren et al. 1990; Southworth et al. 2004; Popper 1974; Stickland et al. 1992). Of all of these systems, V578 Mon is also the youngest, making this system a benchmark case for testing stellar evolution models at the youngest ages.
3.1 Johnson UBV and Strömgren Photometry
The available time-series photometry of V578 Mon covers nearly 40 yr and more than one full apsidal motion period. A summary of the various light curve epochs, including filters and observing facilities used, is presented in Table 2. Photometry from Heiser (2010) includes multi-band light curves spanning 1967–2006 from the 16-in telescope at Kitt Peak National Observatory (KPNO) and from the Tennessee State University(TSU) -Vanderbilt 16-in Automatic Photoelectric Telescope (APT) at Fairborn Observatory. The KPNO Johnson light curves comprise 725 data points spanning 1967–1984 with average uncertainties per data point of 0.004 mag computed by Heiser (2010). The APT Johnson light curves span 1994–2006 and consist of 1783 data points with average uncertainties per data point of 0.001 mag for B and 0.002 mag for V (Heiser 2010). Light curves from Hensberge et al. (2000) span 1991–1994 from the 0.5-m Strömgren Automatic Telescope (SAT) at La Silla, with 248 data points in each of the filters and average uncertainty per data point of 0.003 mag (Hensberge et al. 2000). We begin our light curve analysis with the observational errors originally estimated by Heiser (2010) and Hensberge et al. (2000). Table 2 lists these average uncertainties, , as reported by the original authors. However, from our light curve fits (see below) we found that these uncertainties were in most cases underestimated. Thus we also report as in Table 2 the uncertainties that we ultimately adopted for each light curve.
3.2 hermes Spectroscopy
A new series of high-resolution echelle spectra were secured in December 2011 (36 exposures) and February 2012 (8 exposures) with hermes, the fiber-fed high-resolution spectrograph on the Mercator telescope located at the Observatorio del Roque de los Muchachos, La Palma, Canary Islands. hermes samples the entire optical wavelength range (3800-9000 Å) with a resolution of (Raskin et al. 2011). The observations listed in Table 3 cover the orbital cycle uniformly. Groups of two concatenated exposures allow us to obtain a robust estimation of random noise as a function of wavelength, and a check on cosmic ray events surviving the detection algorithm in the data reduction. In total, 44 exposures were obtained at 19 epochs, 16 of which are out of eclipse. One series of six exposures starts near the primary mid-eclipse; One series of two concatenated exposures taken around secondary mid-eclipse has a significantly lower exposure level, but another one consisting of four concatenated exposures starting around secondary mid-eclipse is available.
Exposure times close to 2100 s were used for most spectra, but in case of one out-of-eclipse epoch the exposure time was significantly shorter, 1200 s. The signal-to-noise ratio of the spectra is 50 to 100 at 4000 Å, then rapidly increasing to 120 up to 200 at 5000 Å and remaining close to this level at longer wavelengths. The numbers apply to the sum of two concatenated exposures. The reduction of the spectra has been performed using the heres pipe-line software package. The spectra resampled directly in constant-size velocity bins (ln ), very nearly in size to the detector pixels, were used. Normalization to the continuum is done separately.
The hermes spectra outnumber the caspec spectra used by Hensberge et al.(̃2000), but fall short with regard to signal-to-noise ratio. However, they cover a much larger wavelength region, include epochs in both eclipses and cover the orbit more homogeneously. In the wavelength region covered by both sets, the reconstruction has better signal-to-noise ratio in the caspec set, but the risk of bias due to phase gaps might be higher with the caspec data. Both data sets were obtained in different parts of the apsidal motion cycle.
4.1 Spectral Disentangling & Light Ratio
In V578 Mon binary system the eclipses are partial which causes degeneracy in the light curve solution for the radii of the components. It was checked whether a spectroscopic light ratio has sufficient precision to reduce the degeneracy. This light ratio might be constrained either by the changing line dilution during the eclipse, or/and by constrained fitting of the reconstructed component spectra by theoretical spectra, simultaneously deriving the light ratio as well as the photospheric parameters (Tamajo et al. 2011). In the latter implementation, the light ratio is assumed identical in all observed spectra, hence eclipse spectra are not used.
With partial eclipses of roughly 0.1 mag depth, and less for the secondary eclipse at the epoch of the spectroscopy, line depth in the composite spectrum is affected at the level of 0.5% of the continuum only when the two components have in their intrinsic spectra a line differing by 7% of the continuum depth. The similarity of the components and the rotational broadening in the spectra imply that no metal line approaches this level. Hence, using the changing line dilution to measure the light ratio precisely is challenging. Exceedingly large signal-to-noise ratios would be required to be able to use single or few lines. Including many lines, i.e. large stretches of spectrum offers the opportunity to reduce the requirements on the signal-to-noise ratio. However, bias in tracing the continuum is expected to put an upper limit to the precision with which the light ratio can be measured in a system with components with similar spectra and substantial rotation.
Therefore, we explored the alternative option of constrained fitting,
although it is model-sensitive. Spectral disentangling (Hadrava 1995), further
referred as spd is performed in a spectral
range of about 100 - 150 Å (of the order of 4000 bins) in the
wavelength range 3900 - 5000 Å, centered on prominent lines of
He i, He ii and stronger metal lines. The apsidal motion
study (Garcia et al. 2011) permitted us to fix the eccentricity ,
the longitude of the periastron, , for the epoch of the spectra,
and the time of periastron passage. The spd code used is
spd applied to selected spectral regions of the hermes spectra, well distributed over the full range of Doppler shifts in the orbit (see orbital phases in Table 3), leads to radial velocity amplitudes K1 and K2 compatible with Hensberge et al. (2000) within better than 1 km/s. Thus the spectra are reconstructed using the mean orbital elements (Table 4), now also including regions around and H ( has a broad interstellar band centered on its red wing). For the constrained fitting, optimization was done for hydrogen and helium lines only, and for combinations of them. The reconstructed spectra for both out-of-eclipse and in-eclipse phases are shown in Figures 2 and 3.
The component spectra for different dilution factors can be obtained from a single disentangling computation, followed by an adequate renormalisation. As starting point for the photospheric parameters, K, and K is used, based on the extensive study of Hensberge et al. (2000) and Pavlovski & Hensberge (2005). The surface gravities of the components are fixed to = and = as derived in this paper. This suppresses the degeneracy of line profiles of hot stars in the (temperature, gravity) plane. Calculations for a small grid in has shown that the effect of fixing might produce deviations of about a few tenths of the percentage in determining the light dilution factors.
Optimization of relative light factors includes a search through a grid of theoretical spectra, using a genetic algorithm. A grid of synthetic spectra was calculated assuming non-LTE line formation. The calculations are based on the so-called hybrid approach of Nieva & Przybilla (2007) in which model atmospheres are calculated in LTE approximation and non-LTE spectral synthesis with detailed statistical balance. Model atmospheres are constructed with atlas9 for solar metallicity, [M/H] = 0, and helium abundance by number density, (Castelli et al. 1997). Non-LTE level populations and model spectra were computed with recent versions of detail and surface (Butler & Giddings 1985). Further details on the method, grid and calculations can be found in Tamajo et al. (2011) and Pavlovski et al. (2009).
Depending on the line(s) included, the primary is found to contribute 68 to 72 percent of the total light, with hydrogen lines supporting the larger fractions. Hydrogen suggests a few percent lower temperature for the primary, compared to the starting values. This is compatible with the tendency seen in Figure 7 of Hensberge et al. (2000), that H and He lines for the primary only marginally agree on effective temperature (taking minimum at the relevant gravity, a 1000 K difference in temperature estimation occurs).
The inconsistency between different indicators underlines the importance to develop a more consistent atmosphere model for these stars. One way, following Nieva & Przybilla (2012), is to include more ionization equilibria by analyzing the full wavelength range covered by the new spectra. This work-intensive analysis is out of the scope of the present paper, but probably indispensable to constrain better the degeneracy in the determination of the radii. Its success might be limited by the rotational broadening in the spectra. Another point of attention is the need to take into account temperature and gravity variations over the surface, due to the slightly non-spherical shape of the stars. Our work shows that the purely photometrically estimated light factors (Table 5) lie within the broader range of light factors (primary to total light) derived here from the hermes spectra. However, there are some warnings that improvement is needed - the spectroscopic estimates may be biased as different indicators are not yet fully compatible.
4.2 Light Curve Analysis
We use EB modeling software phoebe (Prša & Zwitter 2005) based on Wilson-Devinney code (Wilson & Devinney 1971; Wilson, 1979) for our light curve analysis. We fit light curves spanning yrs, covering one full apsidal motion cycle and in Johnson and Strömgren photometry.
Figures 4, 5, 6, 7 are the residuals (data-model) for our global best fit model to the light curves for every light curve epoch and filter in Table 2. Overall, the residuals are small - typically mag. The residuals are significantly larger for light curve epochs 1970-1984, since error bars on the photometry data points measured using photometric plates is larger. We explore ranges for our light curve parameters as listed in Table 6. Our global best fit matches observations well - the final light curve parameters , , , and are listed in Table 7.
For our global best fit light curve model, we adopt a square root limb darkening law (Claret 2000), a B1V spectral type for the primary star implying K (Hensberge et al. 2000), no light reflection, and no third light.
We have four light curve parameters of interest - primary potential , secondary potential , inclination and temperature ratio . A “parameter of interest” is defined as a parameter that is varied to compute our confidence intervals. We determine these parameters and their uncertainties by mapping space. Potential is a modified Kopal potential for asynchronous, eccentric orbits (Wilson, 1979). The potential () takes into account contributions from the star itself, its companion, the star’s rotation about its axis, and the star’s rotation in its orbit.
Our fixed parameters are the argument of periastron , eccentricity , apsidal motion , semi-major axis , mass ratio , period , ephemeris HJD, systemic velocity , gravity brightening coefficients & , primary and secondary synchronicity parameters , and albedos , . We fix the argument of periastron , eccentricity and apsidal motion to values determined by a multi-epoch light curve analysis from Garcia et al. (2011). We fix mass ratio , semi-major axis , orbital period , time of minima HJD and systemic velocity to values from Hensberge et al. (2000) analysis of the spectroscopic orbit. As mentioned previously, our hermes spectra analysis derives radial velocity amplitudes and in agreement with the Hensberge et al. (2000) spectroscopic orbit (see Table 4). We adopt gravity brightening coefficients () and surface albedos () to be as appropriate for stars with radiative envelopes. The gravity brightening coefficient for stars with radiative envelopes was first found by von Zeipel (1924). We fixed rotational synchronicity parameters and to values from Hensberge et al. (2000). Our limb darkening coefficients follow the square-root law for hot stars (Claret 2000) and are listed in Table 8.
Our fitting method is adapted from Gómez Maqueo Chew et al. (2014, in prep). We determine our best fit global light curve solution by finding a unique set of light curve parameters , , and that correspond to the minimum chi square in a well mapped grid of parameter space. The chi square is a function of the light curve parameters, . We map parameter space by computing for a grid of unique sets of these light curve parameters. We use our map of parameter space to compute the uncertainties on our light curve parameters using confidence intervals. Plots of vs stellar radii , , temperature ratio , and inclination with confidence intervals are shown in Figure 8.
The step-by-step procedure is as follows:
1. We sample a coarse grid of points defined by a range of potential , potential , inclination and temperature ratio . The parameter ranges and spacings are given in Table 6.
For each grid point, we fit only for the “light levels” in phoebe which is equivalent to the total light contribution from each star in the photometric bandpass. We avoid using the WD2003 differential corrections (DC) fitting algorithm within phoebe to fit our light curve parameters. The DC algorithm can fall into local minima when fitting for many parameters. We compute the total chi square for each light curve fit as the sum of the chi square at each passband and epoch:
Where index corresponds to a unique point in parameter space (, , , ). is the total chi square over all light curves at a unique point . Index corresponds to a unique light curve passband epoch as specified in Table 2. The chi square at specific passband is computed as:
Where is the number of photometry data points minus the number of parameters of interest over all light curve epochs. Each data point has an error bar . Each light curve at a specific epoch and filter has a multiplicative factor which takes into account systematic error. Multiplicative factor is used to normalize the such that or reduced . is the total flux of the binary at an HJD, and flux is the corresponding model. From our coarse grid, we find the minimum total chi square in parameter space.
2. We adjust the error bars of the individual photometry data points for all light curves to take into account any systematic error. For the minimum solution, the passband is computed for each separate light curve epoch and filter using the equation:
Where as in step 1, and is the minimum total of the coarse grid. We choose compute multiplicative factor to weight each light curve such that the minimum reduced chi squared for our global best fit solution. We then rescale the of all other light curve fits using the passband :
Where un-scaled and is the scaled chi square at a unique point in parameter space .
3. We perform steps 1 and 2 for a fine grid of points in parameter space around the location of the minimum . In this way we carefully map out parameter space at the location of the . We use multiple fine grids to precisely find our global best fit minimum. The average grid spacings are , , and respectively for , , and .
We find that the location of the minimum moves slightly, and we recompute the multiplicative factor for each light curve to account for this, again making . Finally, we have a global best fit solution within a finely sampled parameter space. Our global best fit solution listed in Table 7 corresponds to the point in parameter space where scale chi square by such that .
A Comparison of Light Curve Models
In order to ensure our light curve solution is robust and thus our light curve parameters are accurate, we compare our best fit light curve model described above with several other models. As shown in Table 9 we find little effect on our best fit light curve parameters from using different light curve models. All other models are not as favorable due to larger or temperatures that do not agree with the analysis of the component spectra of V578 Mon from spectral disentangling of Hensberge et al. (2000).
For all the tests described below, we start at our best fit solution, then fit all light curves in phoebe for primary potential , secondary potential , temperature ratio , and inclination . Our global best fit uses a fixed primary temperature T K, no light reflection and no third light. Furthermore, our global best fit uses fixed square root law limb darkening coefficients, which are found to work best for hot (K) stars (Diaz-Cordoves & Gimenez 1992; van Hamme 1993). We discuss the different light curve models in the order in which they appear in our summary Table 9:
Fitting for Limb Darkening Coefficients. We test the effect of fitting for square root law limb darkening coefficients, finding a lower chi square due to a larger number of free parameters. We find little effect on , or . However - we do find a much lower . We reject this light curve model since is significantly outside of the acceptable range for from the spectral disentangling of Hensberge et al. (2000). We therefore perform another test: we keep fixed to our best fit value, and fit for the limb darkening parameters, , and . We again find little effect on , or .
Using a different Limb Darkening Law. We test the linear cosine and logarithmic limb darkening laws, finding little effect on our light curve parameters. The linear cosine law has a lower than our best fit model . The light curve model with logarithmic limb darkening has a larger - we therefore reject this model. See Table 8 for a list of the theoretical limb darkening coefficients for each light curve model that we test.
Changing the assumed Primary Star Temperature. We test the effect of changing our adopted primary star effective temperature . Our adopted primary temperature for our best fit solution is K. Once again, we find little effect on , , or .
We start with our best fit global solution, but set K and K, above and below our adopted primary star effective temperature. Fits with lower primary temperature do result in a better - however, K does not agree with the spectral disentangling analysis from Hensberge et al. (2000). This may be due to the fact that the phoebe light curve analysis constrains the temperature ratio and not the individual temperatures themselves. Further light curve tests at lower preferred temperatures and confirmed that changing effective temperatures have little effect on the geometric parameters, , and .
Light Reflection. We fit our light curve model with one light reflection. We find an inclination larger by . However, the is higher than our best fit . We reject this model on this basis.
Third Light. We test the possibility of third light and its effect on our best fit parameters. We fit for a third light parameter starting from our best fit light curve solution. The third light model has a lower due to a larger number of free parameters. We find and to be larger by and respectively from our best fit model.
However, the third light parameter varies on the order of an apsidal period of the system. As shown in Table A Strict Test of Stellar Evolution Models: The Absolute Dimensions of Massive Benchmark Eclipsing Binary V578 Mon, we find at max a small contribution of third light for Johnson B filter of light curve epochs 1967-1984 and 2005-2006. This is likely due to phoebe using the parameter to minimize the small systematic error of mag in the residuals of the 1967-1984 and 2005-2006 light curve epochs. Furthermore, the systemic velocity measured with the hermes spectra and the caspec spectra in Hensberge et al. (2000) does not give any evidence for a large third body in the system that would contribute significantly to the light. This is consistent with the third-light tests performed here.
Uncertainties on Light Curve Parameters
We compute uncertainties on each parameter of interest using confidence intervals as shown in Figure 8. From Press (1988), for four parameters of interest, we find that , , and uncertainties correspond to solutions with confidence intervals of , and respectively. Here, is the minimum of our global best fit solution.
From Figure 8 we see small degeneracies between the geometric parameters - radii , , and inclination . However - as expected we do not see degeneracies between the geometric light curve parameters and the temperature ratio .
Since is not strongly degenerate with these other parameters, we could potentially decrease the number of parameters of interest and in turn decrease the formal parameter uncertainties. Therefore, the uncertainties presented here are possibly conservative, given that we assume all degrees of freedom are parameters of interest (Avni 1976).
The small degeneracies in our parameters leads to uncertainties on potentials and of less than % error - this error already takes into account any systematic error in fitting the light curves, as detailed in §4.2.2. Similarly, the uncertainty on the temperature ratio and inclination are also %.
A source of systematic uncertainty unaccounted for from the confidence intervals and fitting procedure in §4.2.2 is from the comparison of light curve models detailed in §4.2.3 and Table 9. As shown in Table 9, all other light curve models assessed in §4.2.3 with the exception of using linear cosine LD parameters are not as favorable as our best fit model. The linear cosine model has a lower . Nevertheless, The inclination , temperature ratio and secondary potential are all within of our best fit model. However, the primary potential for the linear cosine model as compared to our best fit . Therefore our uncertainty on from our best fit model could be slightly underestimated from these model comparisons.
Consistency of Light Fractions
As mentioned by Torres et al. (2010) an important consistency check of our light curve solution is that the light fraction determined from spectroscopy and photometry agree. Given the small degeneracy between and as seen in Figure 8, we compare our photometrically determined light fraction with the light fraction from the HERMES spectral disentangling and a previous combined light curve and spectral disentangling analysis from Hensberge et al. (2000). We find that all three light fractions agree with each other to within 1.2. A comparison of light fractions is shown in Table 5.
For each of the light curve fits to our 40 yrs of photometry data, we compute the light fraction at each of the passbands Johnson and Strömgren photometry, where and are the contribution of the primary and secondary star to the total light at a specific passband out of eclipse. The distribution of light fractions for light curve models with confidence intervals of and are shown in Figures 9 and 10.
4.3 Comparison with Hensberge et al. (2000)
Hensberge et al. (2000) uses an iterative, combined light curve and spectral disentangling analysis using the Wilson-Devinney light curve modeling program to compute their light cure parameters. We find that from Hensberge et al. (2000) is discrepant from our best fit . We find that our inclination deg is 1.6 discrepant from deg from Hensberge et al. (2000). These discrepancies are likely due to the addition of apsidal motion and an updated eccentricity determined in Garcia et al. (2011). Apsidal motion and eccentricity can affect the potentials and and hence the determination of the radii at a low level. The potential for a non-circular orbit is a function of eccentricity (see Wilson, (1979)). The addition of more light curve epochs may also play a role. Hensberge et al. (2000) only uses the 1991-1994 light curve epoch with Strömgren photometry. As a check, we also recover the Hensberge et al. (2000) light curve solution when we fit only the 1991-1994 light curve epoch. Finally, simply the addition of more photometry data points may play a role. We use 3489 photometry data points in our light curve solution, whereas Hensberge et al. (2000) use 992. Our best fit secondary radius is in agreement with from Hensberge et al. (2000). Our best fit temperature ratio is in agreement with the temperature ratio of from an analysis of the disentangled component spectra (Hensberge et al. 2000).
5 Results: Absolute Dimensions and Apsidal Motion of V578 Mon
The absolute dimensions and other fundamental properties of V578 Mon are compiled in Table 11. Here we detail how each fundamental parameter for V578 Mon is compute in order of which they appear in Table 11:
Orbital Period. We adopt the orbital period d from Hensberge et al. (2000).
Masses. The component masses M and M are determined from the spectroscopic orbit analysis from Hensberge et al. (2000). We do not use RVs from our hermes spectroscopy because the caspec spectra have higher S/N - however, our analysis of the hermes spectroscopy reconfirm the spectroscopic orbit.
Radii. We find precise uncertainties of % for the primary radius and secondary radius from our confidence intervals in Figure 8.
Temperatures. We find a % error on our temperature ratio from our confidence intervals. Combined with the adopted temperature of the primary star T K (Hensberge et al. 2000), our temperature ratio of yields a secondary temperature of T K via propagation of errors.
Rotational Velocities. We compute surface rotational velocities of v km s and v km s using the observed projected surface velocities v km s and v km s from Hensberge et al. (2000) and our inclination of i. The uncertainty on rotational velocities are computed from propagating the error on the inclination and the observed v.
Surface Gravities. We compute the surface gravity from our masses and radii, finding cm s and cm s. We compute the uncertainty on via error propagation:
Where is the uncertainty on the mass and is the uncertainty on the radius.
Luminosities. From our radii and temperatures, we compute compute luminosities for the primary and secondary star of and . We compute the uncertainty on the luminosity using a similar error propagation as above, using errors from the temperature and radii, and .
Synchronicity Parameters. We find the components of V578 Mon to be close but not exactly tidally locked, with and . The synchronicity parameter , where is the rotational velocity at the surface and is the synchronous velocity. We compute the uncertainty via propagation of error from , error on inclination , and error on projected rotational velocities .
Internal Structure Constant. One of us (Dr. Claret) computes the newtonian and general relativistic contributions to the observed internal structure constant, and .
6 The Stellar Evolution Models and Tests
We compare the absolute dimensions of V578 Mon to
the stellar evolution models of three separate groups:
(1) Geneva models of Georgy et al. (2013) and Ekström et al. (2012) hereafter Geneva13;
(2) Utrecht models of Brott et al. (2011) here after Utrecht11
A comparison of the basic input physics of the models is given in Table 12. The models use the same opacity tables of Iglesias & Rogers (1996). The mixing length for all three sets of models differ by only at maximum. The stellar evolution models use similar mass loss treatment from the prescription by Vink et al. (2001). Given the probable young age of V578 Mon due to its location in the open cluster NGC 2244 of the Rosette Nebula, the components of V578 Mon are not expected to have undergone significant mass loss (Vink et al. 2001).
However, all three sets of models differ on the choice of the convective core overshoot parameter . For the H and He burning phases of the convective core, the convective core size of the star is enlarged by , where in units of pressure scale height. The overshoot parameter is designed to accomodate for the non zero velocity of the material moving from the convective core to radiative zone of the star. Observationally, a larger overshoot parameter means longer MS lifetimes for a given star, and thus older ages. The Geneva13 models use a small convective core overshoot of calibrated on width of the main sequence for stars with masses which is characterized by the red most point on the B-V, HR diagram (see figure 8 of Ekström et al. (2012)). The width of the main sequence is defined theoretically by the end of the hydrogen burning phase. The Utrecht11 models use a high convective core overshoot of which is calibrated using the observed width of the main sequence from the VLT-FLAMES survey of B stars (Evans et al. 2005; Hunter et al. 2007). The convective core overshoot parameter is chosen such that a 16 star ends its MS lifetime when . This coincides with the drop in B star rotation rates in a - diagram, which is interpreted as an estimate of the width of the main sequence for B stars. See Brott et al. (2011) for an in depth discussion. The Granada04 models utilize a moderate convective core overshoot , though we performed several tests varying .
Rotationally driven mixing can bring more H and He from the envelope to the core, thus extending the MS lifetime of the star - likewise a larger overshoot parameter extends the size of the core, leading to a longer MS lifetime. The Granada04 models do not incorporate rotational mixing, while the Geneva13 and Utrecht11 models do. However - All three sets of models include rotation. All three sets of models use similar metallicity compositions of near solar. The initial bulk composition for V578 Mon is expected to be close to solar given that Mg surface abundance is within error of the solar surface abundance, despite the fact that several atmospheric abundances such as C, N and O are somewhat metal poor compared to the Sun (Pavlovski & Hensberge 2005). This is because Mg abundance is not expected to be altered from the initial abundance in a star, where as C, N and O atmospheric abundances could vary in V578 Mon due to rotational mixing (Lyubimkov et al. 2005). However, given that the C, N and O atmospheric abundances of V578 Mon may be lower than solar, the metallicity of V578 Mon still remains as a source of systematic error in comparing the evolution models to the observations.
The Granada04 models also compute the internal structure constants and allowing for a test of the internal structure of V578 Mon via apsidal motion. Here we consider only the constant, given that and are very small. For V578 Mon, the tidal Love numbers quantify the deformation for each star’s gravity field due to the companion.
6.1 Isochrone Test for V578 Mon
In Figure 11, we place the primary and secondary star on mass-radius and isochrones for each set of models. For the stellar evolution models to pass the “isochrone test” the models should predict a common age for both components of V578 Mon within uncertainty. Given how different the masses of the primary and secondary star for V578 Mon are the “isochrone test” provides a stringent test of stellar evolution models. We also match all evolution models to the rotational velocities of the primary and secondary star.
We find several Geneva13, Utrecht11 and Granada04 models predict masses, radii and temperatures for the components of V578 Mon that fall within uncertainty of the measured absolute dimensions. Therefore we estimate an age range for each star as shown in Table 12. The age difference for Geneva13, Utrecht11 and Granada04 models is given as the smallest possible difference between the ages of the two stars given age range of each star.
For the Geneva models we use isochrones with
initial rotational velocities and
which allows us
to match the observed rotational velocities for each star. We interpolate the
model evolution tracks for the primary and secondary star using the online
interactive tool provided by the Geneva group
For the Utrecht11 models we use isochrones that match the observed surface velocities of the components of V578 Mon, v km s and km s. The Utrecht11 models are computed at very small steps of mass and initial rotational velocity, such that interpolating between model tracks is unnecessary. From Table 13 we a marginally common age (age difference Myr) for mass-radius isochrones, and a common age of Myr for isochrone. The models were computed at solar metallicity by Dr. Brott (private comm.).
We compute the Granada04 models at the masses of primary and secondary star and chose rotational velocities to match the observed rotational velocities of V578 Mon. We attempt to match the absolute dimensions of V578 Mon to or alternatively mass-radius isochrones for V578 Mon. We find an age gap of 1.5 Myr for mass-radius isochrones, and a marginally common age for isochrones, when both stars have an overshoot of . Again - finding a match on the isochrones is easier given the greater uncertainty in the effective temperatures.
In an attempt to match the ages of the two stars on a mass-radius isochrone, we also compute Granada04 models for and . Figure 12 demonstrates the time evolution of the radii for V578 Mon for these different models. We find a near match on a single mass-radius isochrone with an age difference of only 0.2 Myr - if we assume that the primary star has a convective overshoot and the secondary star has a convective overshoot of . We also find a common age of Myr for the isochrone. This does not mean that an for the primary star is correct for V578 Mon - merely that a higher convective overshoot allows for compatible ages between the two stars. High convective overshoot has been found to work in matching other EBs on a single isochrone (Claret 2007).
In general, we find younger ages by Myr for the Utrecht11 models of V578 Mon and similar ages for the Geneva13 and Granada04 models. This can be attributed to the larger convective overshoot of included in Utrecht11 models than in Geneva13 models (). While the primary star for the Granada04 models does have an even higher convective overshoot of , the models do not include rotational mixing, which also extends the main sequence lifetime of the stars.
6.2 Apsidal Motion Test for V578 Mon
Measurement of apsidal motion in eccentric binary systems allow for a stringent test of the internal structure constant predicted from stellar evolution models (e.g. Claret & Giménez 2010). It is not possible to separate out each individual star’s contribution to the apsidal period U from newtonian apsidal motion.
The apsidal motion for V578 Mon was measured by Garcia et al. (2011). The observed apsidal motion of the V578 Mon, deg cycle, has contributions from both newtonian and general relativity components (Claret & Giménez 2010):
where is given by
We find that which is only of the newtonian apsidal motion .
Both the newtonian and general relativistic observed apsidal motions and have associated observed internal structure constants . The internal structure constant is twice the tidal love number (Kramm et al. 2011), and is related to the density profiles, degree of sphericity, orbital parameters, masses, and rotation rate of both components of a binary star. Specifically, the internal structure constant is related to the solution of the Radau differential equation as in equation 3 of Claret & Giménez (2010). Importantly - constant is one the few ways to directly constrain the internal structure of stars.
From the precise observed apsidal motion, we compute the observed internal structure constant, = , where U is the apsidal period, given by the equations (adopted from Claret & Giménez (2010)):
We compute the internal structure constant due to the newtonian apsidal motion, , and due to general relativity, . The newtonian apsidal motion is much larger than the general relativistic component, and therefore the internal structure constant is also much larger.
We compute the theoretical internal structure constant, using the methods of Claret & Giménez (2010). The theoretical constant was corrected for by rotation (Claret 1999) and dynamical tides (Willems & Claret 2002). The theoretical internal structure constant is a combination of the internal structure constants for both star, such that
which can then be compared to observations.
We find the predicted newtonian apsidal motion to be and consequently the predicted newtonian internal structure constant to be . This is in very good agreement with the observed . From equation 9, the parameter is about 67% larger than . Therefore, the weighted contribution of the primary dominates the theoretical apsidal motion. V578 Mon is a relatively young system - therefore, is almost constant during the early phases of stellar evolution. The “apsidal motion test” is therefore complementary to the “isochrone test”. Claret & Giménez (2010) compile a list of eclipsing binaries with apsidal motion, demonstrating good agreement between observed and predicted apsidal motions. V578 Mon continues this trend of agreement between theoretical and observational internal structure constants. For this relatively young system, matching the radii, temperatures and masses isochrones is key, given that we have so few young massive EBs with non-equal mass ratio.
We have determined the absolute dimensions of the massive, detached eclipsing binary V578 Mon, which is a member of young star forming region NGC 2244 in the Rosette Nebula. We confirm that the the previously published spectroscopic orbit of Hensberge et al. (2000) agree with our current spectroscopic orbit of V578 Mon. From our hermes spectra, we find that our photometric light ratio from the light curve analysis is fully compatible with the disentangled component spectra of V578 Mon.
From 40 yr of Johnson UBV and Strömgren photometry we determine updated radii, measure the temperature ratio and light ratio for the components of V578 Mon. We determine the radii to better than accuracy, and carefully map out parameter space in order to reveal any possible degeneracies. We also compare our global best fit light curve model with models that include different limb darkening parameters, a different assumed temperature for the primary star, light reflection or third light finding little effect on our global model. We do not unambiguously rule out light reflection or a third body, but we confirm that these additional complications to the light curve model will not affect our final solution.
We have compared our observed masses, radii, temperatures and rotational velocities to stellar evolution models of the Geneva, Utrecht, and Granada groups. We find no common match in predicted ages for mass-radius isochrones of the Geneva13 models. We find an age difference of only Myr in predicted ages for the Geneva13 models for isochrones. For the Utrecht11 models, we find a marginally common predicted age with an age difference of only 0.4 Myr for the mass-radius isochrones. For the isochrones we find common ages of Myr for the Utrecht11 models. For the Granada04 models, we find a small age gap of only 0.2 for the mass-radius isochrone, when the primary star has a quite large convective overshoot of . We do not find common ages for the mass-radius isochrone for the Granada04 models when the convective overshoot for both stars is a more moderate .
This work suggests that models with larger convective overshoot predict a closer common age for the components of V578 Mon than models with a more conventional overshoot of pressure scale heights. Evolutionary models with larger convective overshoot extends the size of the convective core for massive stars, thus extending the main sequence lifetime and allowing for isochrones to predict a common age for V578 Mon. However - rotational mixing also can prolong the main sequence lifetime, making the two effects some what degenerate. The radii may in a small way be dependent upon effective temperatures, which are based on imperfect atmosphere models. Furthermore, there are small systematic residuals of mag in the light curve fits which may in a small way affect the radii. Finally, effects of binarity, while likely small, are not taken into account: the side of each star facing the other may be heated and the addition to the potential from the companion is not taken into account into the models. The binarity of V578 Mon may cause single star models explored here to not be applicable.
Given the short apsidal period of V578 Mon of years, our photometry cover one full apsidal motion period. Combined with our precise measurement of the radii of V578 Mon we compute the internal structure constant finding that our observed in agreement with the theoretical internal structure constant .
V578 Mon is a particularly important system for testing stellar evolution models given young age and the difference of in the masses of the primary and secondary component star. B-type detached eclipsing binaries such as V1388 Ori and V1034 Sco have similar differences in mass of 40% and 50% respectively, meaning these systems are also of particular importance to providing constraints on stellar evolution models. However, V578 Mon is unique among such systems by virtue of its young age, thus providing the strongest constraints on the models at the earliest stages of massive stellar evolution.
Future work may include comparing the carefully vetted sample of high mass EBs in the Torres et al. (2010) sample to evolutionary models, include more recent massive EBs such as V 380 Cyg (Tkachenko et al. 2014) and LMC 172231 and ST2-28 (Massey et al. 2012), to see if larger convective overshoot parameters allow for common predictions of age.
|Henry Draper number||HD 259135||Cannon & Pickering (1923)|
|Bonner Durchmusterung||BD +04°1299||Argelander (1903)|
|Hoag number||NGC 2244 200||Hog et al. (1998)|
|06 32 00.6098||Hog et al. (1998)|
|+ 04 52 40.902||Hog et al. (1998)|
|Spectral type||B0V + B1V||Hensberge et al. (2000)|
|8.542||Ogura & Ishida (1981)|
|0.262||Wang et al. (2008)|
|+ 0.165||Ogura & Ishida (1981)|
|0.727||Wolff et al. (2007)|
|+ 0.452||Wang et al. (2008)|
Note. – 16-inch telescope at Kitt Peak (KPNO)
0.5 m telescope at La Silla (SAT)
TSU-Vanderbilt 16-inch telescope at Fairborn University (APT)
|Phase||BJD-2450000.000||Exp Time [s]|
Note. – Time-series Hermes spectroscopy of V578 Mon. Each exposure is less than 0.01 of the orbital period for V578 Mon of days. The time series spectra were obtained to cover the out-of-eclipse, primary eclipse and secondary eclipse phases.
|[km s]||[km s]|
|Hensberge et al. (2000) (LC+Spectroscopy)||0.7078||259.8||183.9||0.0867|
|Hensberge et al. (2000) RV only||0.7050.004||259.80.8||184.4||0.0836|
|HERMES Spectra, fixed||0.710||259.8||184.5||0.07755|
|HERMES Spectra, and fixed||0.709||259.4||184.0||0.07755|
|Light Curve Analysis (this work)||Johnson U, 365|
|Johnson B, 445|
|Johnson V, 551|
|Stromgren u, 365|
|Stromgren v, 411|
|Stromgren b, 467|
|Stromgren y, 547|
|Hensberge et al. (2000)||Stromgren v, 411|
|Stromgren b, 467|
|Stromgren y, 547|
|Parameter||Max||Min||Coarse Grid Spacing||Fine Grid Spacing|
|Primary Surface Potential,||5.36||4.80|
|Secondary Surface Potential,||5.26||4.40|
|Light Curve Parameters||This Work & Garcia et al. (2011)||H2000|
|Primary Surface Potential,|
|Secondary Surface Potential,|
|Angle of Periastron, [deg]|
|Ephemeris, HJD [d]|
|Total Apsidal Motion, [deg cycle]|
|Light Curve Filters||Strömgren , Johnson||Strömgren|
|Total Light Curve Points|
Note. – The uncertainties on light curve parameters , , and are determined from confidence intervals in Figure 8. Light curve parameters , and are taken from Garcia et al. (2011). This work utilizes photometry that span one full apsidal motion period (Uyr). In contrast to the Hensberge et al. (2000) analysis, this work incorporates apsidal motion in the light curve model. Finally, the temperature ratio from Hensberge et al. (2000) is measured from spectral disentangling.
|Square Root Law (adopted)|
Note. – Our best fit model uses the square root limb darkening law. Fits with the linear cosine or logarithmic limb darkening law had little effect on our final light curve solution.
|Fitting for LD coefficients|
Note. – The best fit model uses the square root limb darkening law, a fixed 30000 K, no light reflection, and no third light.
Note. – Our best fit light curve model includes no third light. The small amount of third light varies as a function of epoch.
|Orbital Period, P [d]|
|Effective Temperature, [K]|
|Surface Gravity, [cm s]|
|Surface Velocity, [km s]|
|Apsidal Period, U [yr]|
|Observed Newtonian Internal Structure Constant,|
|Overshoot,||0.10||0.355||0.6 pri, 0.2 sec|
|Opacities||Iglesias & Rogers (1996)||Iglesias & Rogers (1996)||Iglesias & Rogers (1996)|
|Mass loss||Vink et al. (2001)||Vink et al. (2001)||Vink et al. (2001)|
|Model||Primary Age||Secondary Age||Age Diff (lower limit)|
|Granada04||5.05.3||5.56.3||0.2||0.6 pri, 0.2 sec|
|Granada04||4.75.5||4.96.8||Common age||0.6 pri, 0.2 sec|
Note. – The ages for the primary and secondary star are computed from evolutionary tracks at the masses of either star and solar metallicity. The Granada04 models were computed for a high convective overshoot of pressure scale heights for the primary star, which allowed the models to match the observations. It easier to find a common age for the isochrone given the larger uncertainty on the effective temperatures of the stars.
- affiliation: Department of Physics & Astronomy, Vanderbilt University, 6301 Stevenson Center, VU Station B 1807, Nashville, TN 37235; email@example.com
- affiliation: Lowell Observatory Pre-doctoral Fellow, Lowell Observatory, Flagstaff, AZ 86001
- affiliation: Department of Physics & Astronomy, Vanderbilt University, 6301 Stevenson Center, VU Station B 1807, Nashville, TN 37235; firstname.lastname@example.org
- affiliation: Department of Physics, Fisk University, 1000 17th Ave. N., Nashville, TN, USA 37208
- affiliation: Department of Physics, University of Zagreb, Bijeniÿcka cesta 32, 10000 Zagreb, Croatia
- affiliation: Royal Observatory of Belgium, Ringlaan 3, 1180 Brussels, Belgium
- affiliation: Physics Department, University of Warwick, Gibbet Hill Road, Coventry CV4 7AL, United Kingdom
- affiliation: Instituto de Astrofsica de Andaluca, CSIC, Apartado 3004, 18080, Granada, Spain
- The Utrecht Stellar Evolution group is now located in Bonn, Germany
- Argelander, F. W. A. 1903, Eds Marcus and Weber’s Verlag, Bonn (1903), 0
- Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481
- Avni, Y. 1976, ApJ, 210, 642
- Bell, S. A., Hilditch, R. W., & Adamson, A. J. 1986, MNRAS, 223, 513
- Bouzid, M. Y., Sterken, C., & Pribulla, T. 2005, A&A, 437, 769
- Butler, K., & Giddings, J. R. 1985, Newsletter of Analysis of Astronomical Spectra, No. 9
- Brott, I., de Mink, S. E., Cantiello, M., et al. 2011, A&A, 530, A115
- Cannon, A. J., & Pickering, E. C. 1923, Annals of Harvard College Observatory, 98, 1
- Castelli, F., Gratton, R. G., & Kurucz, R. L. 1997, A&A, 318, 841
- Claret, A. 1999, A&A, 350, 56
- Claret, A. 2000, A&A, 363, 1081
- Claret, A. 2004, A&A, 424, 919
- Claret, A. 2006, A&A, 453, 769
- Claret, A. 2007, A&A, 475, 1019
- Claret, A., & Giménez, A. 2010, A&A, 519, A57
- Diaz-Cordoves, J., & Gimenez, A. 1992, A&A, 259, 227
- Evans, C. J., Smartt, S. J., Lee, J.-K., et al. 2005, A&A, 437, 467
- Ekström, S., Georgy, C., Eggenberger, P., et al. 2012, A&A, 537, A146
- Garcia, E. V., Stassun, K. G., Hebb, L., Gómez Maqueo Chew, Y., & Heiser, A. 2011, AJ, 142, 27
- Garcia, E. V., Stassun, K. G., & Torres, G. 2013, ApJ, 769, 114
- Georgy, C., Ekström, S., Granada, A., et al. 2013, A&A, 553, A24
- Gómez Maqueo Chew, Y., Stassun, K. G., Prša, A., et al. 2012, ApJ, 745, 58
- Hadrava, P., 1995, A&AS, 114, 393
- Heiser, A. M. 1977, AJ, 82, 973
- Heiser, A. M. 2010, Journal of the American Association of Variable Star Observers, 38, 93
- Hensberge, H., Pavlovski, K., & Verschueren, W. 2000, A&A, 358, 553
- Hunter, I., Dufton, P. L., Smartt, S. J., et al. 2007, A&A, 466, 277
- Hog, E., Kuzmin, A., Bastian, U., et al. 1998, A&A, 335, L65
- Holmgren, D. E., Hill, G., & Fisher, W. 1990, A&A, 236, 409
- Iglesias, C. A., & Rogers, F. J. 1996, ApJ, 464, 943
- Ilijic, S., Hensberge, H., Pavlovski, K., & Freyhammer, L. M. 2004, Spectroscopically and Spatially Resolving the Components of the Close Binary Stars, 318, 111
- Kramm, U., Nettelmann, N., Redmer, R., & Stevenson, D. J. 2011, A&A, 528, A18
- Lyubimkov, L. S., Rostopchin, S. I., Rachkovskaya, T. M., Poklad, D. B., & Lambert, D. L. 2005, MNRAS, 358, 193
- Massey, P., Morrell, N. I., Neugent, K. F., et al. 2012, ApJ, 748, 96
- Nieva, M.-F., & Przybilla, N. 2012, A&A, 539, A143
- Nieva, M. F., & Przybilla, N. 2007, A&A, 467, 295
- Ogura, K., & Ishida, K. 1981, PASJ, 33, 149
- Pavlovski, K., & Hensberge, H. 2005, A&A, 439, 309
- Pavlovski, K., Tamajo, E., Koubský, P., et al. 2009, MNRAS, 400, 791
- Pribulla, T., Rucinski, S. M., Latham, D. W., et al. 2010, Astronomische Nachrichten, 331, 397
- Popper, D. M. 1974, ApJ, 188, 559
- Popper, D. M., & Hill, G. 1991, AJ, 101, 600
- Popper, D. M., & Etzel, P. B. 1981, AJ, 86, 102
- Prša, A., & Zwitter, T. 2005, ApJ, 628, 426
- Press, W. H., Teukolsky, S. A., Vetterling, W. T. and Flannery, B. P., 1988, Numerical Recipes in C, Second Edition, Cambridge University Press, NY, pg 697.
- Raskin, G., van Winckel, H., Hensberge, H., et al. 2011, A&A, 526, A69
- Song, H. F., Maeder, A., Meynet, G., et al. 2013, A&A, 556, A100
- Stickland, D. J., Koch, R. H., & Pfeiffer, R. J. 1992, The Observatory, 112, 277
- Southworth, J., Maxted, P. F. L., & Smalley, B. 2004, MNRAS, 351, 1277
- Sterne, T. E. 1939, MNRAS, 99, 451
- Tamajo, E., Pavlovski, K., & Southworth, J. 2011, A&A, 526, A76
- Tkachenko, A., Degroote, P., Aerts, C., et al. 2014, MNRAS, 438, 3093
- Torres, G., Andersen, J., & Giménez, A. 2010, A&A Rev., 18, 67
- Wang, J., Townsley, L. K., Feigelson, E. D., et al. 2008, ApJ, 675, 464
- Wilson, R. E. & Devinney, E. J. 1971, ApJ, 166, 605
- Wilson, R. E. 1979, ApJ, 234, 1054
- Willems, B., & Claret, A. 2002, A&A, 382, 1009
- Wolff, S. C., Strom, S. E., Dror, D., & Venn, K. 2007, AJ, 133, 1092
- van Hamme, W. 1993, AJ, 106, 2096
- Vink, J. S., de Koter, A., & Lamers, H. J. G. L. M. 2001, A&A, 369, 574
- von Zeipel, H. 1924, MNRAS, 84, 665