Free Access
Issue
A&A
Volume 624, April 2019
Article Number L3
Number of page(s) 5
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/201935094
Published online 01 April 2019

© ESO 2019

1. Introduction

Classical Wolf-Rayet (WR) stars are evolved massive objects that are in advanced nuclear burning stages, have shed their outer hydrogen-rich layers, and possess strong stellar winds with enhanced abundances of nuclear processed elements (Maeder & Meynet 2000; Crowther 2007; Langer 2012). Two basic channels are believed to give rise to these WRs. The first is the binary interaction channel, whereby the originally more massive star fills its Roche lobe as it leaves the main sequence and transfers mass to its companion, thereby shedding its outer layers and exposing the regions close to the nuclear burning core, as described by van den Heuvel (1976). The second channel is the single-star channel, whereby the outer layers are shed through a combination of stellar wind and violent mass outflow events occurring near the end of the core hydrogen-burning phase (Conti 1975). The second channel is generally favored for stars lacking eclipses or a periodic radial velocity variation curve, even though several of these are extremely variable and have long been suspected to be in binary systems, the most vexing example of which is EZ CMa.

EZ CMa (HD 50896 = WR6, van der Hucht 2001) is the brightest WR star that can be observed from the northern hemisphere. The variability in its emission lines and light intensity were first reported by Wilson (1948) and Ross (1961), respectively. Subsequent studies disclosed possible variability periods of 1 d (Kuhi 1967; Lindgren et al. 1975) and 13 d (Schmidt 1974), but clear periodicity eluded observers until Firmani et al. (1978) discovered the now well-established 3.76 d period, confirming it shortly thereafter (Firmani et al. 1980). The authors concluded that EZ CMa is a binary system containing a ∼1.3 M companion, which, they proposed, could be a collapsed companion such as a neutron star because such WR + cc systems are predicted within the scenario proposed by van den Heuvel (1976).

The collapsed-companion scenario was first challenged by Stevens & Willis (1988) and Pollock (1989), who argued that the X-ray luminosity is too low for the expected accretion rate onto a neutron star, unless some mechanism is present that inhibits accretion. However, Skinner et al. (2002) analyzed a scenario in which the companion is a non-degenerate low-mass star and showed that the hard component of the observed X-rays could be explained by the shock that is formed when the WR wind collides with the surface of such a companion.

A significantly more serious challenge for the binary scenario is the loss of coherence in the 3.76 d periodic variations when observations covering timescales longer than about two weeks are analyzed. This has prompted scenarios invoking a rotating inhomogeneous wind or disk (Underhill & Yang 1991; Schulte-Ladbeck et al. 1991) or one with two oppositely directed outflows of enhanced density (Matthews et al. 1991), which eventually morphed into models involving corrotating interaction regions, CIRs (Morel et al. 1998; St-Louis et al. 2018; Moffat et al. 2018). According to these scenarios, the 3.76 d period corresponds to the rotation of the WR core, with a yet to be explained instability producing the streams or inhomogeneities in the wind that cause the observed variability. One common feature of these purported instabilities is that their location drifts over time with respect to the underlying rotating stellar surface, hence leading to the incoherence when observed over long timescales.

In this paper we solve the puzzle of the non-coherence in the variability of EZ CMa by showing that it is caused by an exceptionally fast apsidal motion of the eccentric binary orbit. This result is made possible through the recent observations obtained by the BRITE satellite constellation (Moffat et al. 2018), which provides a nearly continuous data train over 136 days.

2. Interpretation of the photometric variability

Moffat et al. (2018) observed EZ CMa with the BRITE satellite constellation from 2015 November 15 to 2016 March 29. We reproduce their published data in our Fig. 1. Adopting the model in which EZ CMa is a binary system with a WR primary and a low-mass companion as proposed by Skinner et al. (2002), we interpret this light curve as an overall ∼3.6 d variation that is twice eclipsed, yielding basically two minima per orbit. The overall variation is due to brightening of the shocked WR wind during the time of periastron passage. The eclipses are produced when either the WR or the companion occult a portion of the shocked wind emission. Below we present clear orbital signatures by analyzing the timings of brightness minima.

thumbnail Fig. 1.

Copy of Fig. 1 published by Moffat et al. (2018), displaying the photometric light curve of EZ CMa as observed by the Toronto BRITE satellite from 2015 November 14 until 2016 March 29. The measured times of well-defined minima are indicated by green dots. The calculated times of upper and lower panels conjunctions of the orbit, with the parameters given in Table 2 that are described in the text, are marked by plus signs. The uneclipsed photometric light curve from a shocked zone, which is assumed to vary proportionally to the square of the orbital separation, is given by the orange curve. The light from the shock zone is eclipsed twice each orbit, and the light curve resulting from a simple Gauss-shaped eclipse model is shown by the blue line.

Open with DEXTER

Light from a shocked wind zone has been identified by Lamberts et al. (2017), who resolved the WR binary γ2 Vel and estimated that the brightness contribution from the shock zone amounts to about 10% of the system’s total brightness. Richardson et al. (2017) showed that the periodic photometric brightening is produced by spectral line emission and that it varies proportional to the distance between the stars. Similarly for EZ CMa, we modeled the basic variability as the brightness contribution of the shocked WR wind, which varies with the distance between the objects. Below we derive the eccentricity of the system to be about 0.1, which implies that the emission from the shocked region is variable by ±20%, if it varies with the square of the separation. The observed variability in the BRITE satellite filter is on the order of ±0.03 mag, which means that for EZ CMa, the shock zone contributes approximately 15% to the total light. The computed uneclipsed light curve is shown by the orange curve in Fig. 1.

The observed light curve then means that the light from the brightest part, that is, the most variable part, is eclipsed twice each orbit by either the WR wind or by the as yet undetected companion. Thus, there is a minimum twice in each orbit, similarly as in close binary systems, at both the upper and lower conjunction. We find that basically only the excess emission is occulted because the eclipse depth typically reaches the minimum value observed out of eclipse. As a consequence, the eclipses are most clearly seen when the conjunctions coincide with periastron passages, and less so in the other orbital configurations.

3. Determining the apsidal motion

In the BRITE light curve we measured the timings of clearly visible minima, which are marked by green dots in Fig. 1. We note that the intervals between the minima are not constant but always alternate systematically between short and long duration, as shown in Fig. 2. This interval pattern is a distinctive signature of an apsidal motion in an eccentric binary orbit. Random drifting of the minima timings on timescales of fractions of the period have been observed in low-mass contact binaries due to reasons other than an apsidal motion (Tran et al. 2013). However, in the case of EZ CMa, we observe that the argument of periastron systematically drifts over more than 360°, which challenges such explanations. Additional evidence for an apsidal motion comes from radial velocities that do not follow the anomalistic phase but agree with the sidereal phase (Koenigsberger, in prep.).

thumbnail Fig. 2.

Measured intervals between consecutive minima. The red dots mark the differences in even number minus odd number minima, and green dots show the differences in odd minus even minima. The green and red lines are sine curve fits to the measurements. The blue crosses connected by the blue dotted line show the calculated intervals from the final orbital solution (see text).

Open with DEXTER

A double sine curve fit yields the orbital parameters given in Table 1. The position of the periapsis for the time T is given by the equation

Table 1.

Orbital parameters derived from a double sine curve fit to the observed minima intervals (Fig. 2).

(1)

where is the rate of periapsis advance, is the epoch when the periapsis angle is 90°, and Pa is the anomalistic period, the time between two passages at periastron. The apsidal motion completes one 360° cycle within the apsidal period U. The mean sidereal period, Ps, is the mean time that elapses between consecutive conjunctions, upper or lower, respectively. As there are two eclipses for each orbit, there is half the mean sidereal period equal to the mean interval value. Half the sidereal period and the apsidal period result directly from the sine curve fit, whereas the anomalistic period needs to be calculated from the relation between the three periods,

(2)

where the apsidal motion in degrees per anomalistic period is given by .

When is negative, which is the case for EZ CMa, the subtraction in Eq. (2) changes to an addition in the relation

(3)

if the apsidal period U is given as a positive number. The fact that the periastron angle is retreating yields the sidereal period Ps = 3.77 d, which is longer than the anomalistic period, Pa = 3.64 d. This can be verified by counting in Fig. 1 that there is one more light maximum than eclipse-pairs over the 100 days of the apsidal period U.

4. Interpretation of the O-C diagram

For a clear detection of an apsidal motion we have used the intervals between successive minima (Fig. 2). For the purpose of determining an apsidal motion, a O-C diagram may be used as well, which illustrates the difference between an observed minimum epoch and one computed with a linear ephemeris. In Fig. 3 we display the O-C diagram for the measured minima, which are indicated in Fig. 1. The apsidal motion is also clearly visible in the O-C diagram, but in addition, an overall parabolic trend is indicated by the quadratic fit. This dependence indicates a changing period. From the measured minima timings, it cannot be determined which period is changing. The O-C diagram can be fitted equally well with a change in sidereal period as with a change of the anomalistic period.

thumbnail Fig. 3.

Observed epochs of minima, marked with red dots for even number minima and with green dots for odd minima, minus computed minima epochs calculated with the relation T(E) = E0 + PE, with P = 1.88 d, which is half the sidereal period. E is the number count of the minima starting with zero. The dashed curve is a quadratic fit to all measured minima. The blue crosses mark the calculated minima from the final orbital solution relative to the linear ephemeris.

Open with DEXTER

If the changing period were the anomalistic period, then the observed trend seen in Fig. 3 would require a = −5.7 × 10−4 d/d. This rate is too fast given the reported periods for EZ CMa in the literature, which are all of about 3.7 d (see Georgiev et al. 1999). An estimate of a possible period decrease, by including historic data, yields a much lower value, a = −1 × 10−5 d/d, which is a considerably fast period change but not sufficient to reproduce the observed shape in the O-C diagram.

We therefore fit the O-C measurements of Fig. 3 allowing for a change in apsidal motion. Equation (1) is now modified to read

(4)

where accounts for the change in rate of the periapsis retreat.

We solved the binary motion for the full time interval of Fig. 1 in time steps of 0.1 days with the changing argument of periastron as given by Eq. (4). From this resulting table of projected positions in the sky, (x(t),y(t),z(t)), we determined the epochs of conjunction by searching for the sign change in the projected x-coordinate and interpolated between the two epochs of the adjacent grid points. The fit comprises six parameters, which were used to minimize the difference squared between the calculated minima epochs and the measurements as well as the difference between the calculated and observed intervals shown in Fig. 2. Table 2 lists the best-fit orbital parameters.

Table 2.

Orbital parameters derived from fitting the measurements of minima epoches shown in Fig. 3 and simultaneously the intervals shown in Fig. 2.

5. Discussion

We have derived the orbital parameters exclusively from the timing of the epochs of mid-eclipse, and did not use the shape of the observed light curve. We now construct a synthetic light curve based on our eclipse model by subtracting a Gaussian function at each eclipse from the shocked wind light curve. The free parameters are the depth and the width of the Gaussian, and the inclination of the system. In Fig. 1 we have chosen the inclination i = 74° and adjusted the depth. Visually, however, the same result is obtained with i = 60° and a larger depth. We find that basically the two eclipses are identical, that is, when the shock zone is eclipsed by either the WR or the unseen object. The resulting light curve is plotted in Fig. 1 with the blue line. The agreement with the observations is surprisingly good given the simplicity of our eclipse and wind collision models.

From the results of the hydrodynamic simulations made for the γ2 Vel system, it is clear that the reality of eclipsing the shock zone is much more complex than can be represented by a simple Gaussian curve. It is also clear that the emission from the shock zone, which basically consists of emission lines, is variable (Richardson et al. 2017; Lamberts et al. 2017). Therefore, the calculated light curve may be shifted up or down, depending on whether a part of the zone is in a higher or lower state. Such a misfit usually persists for several orbital periods at the same phase. It is interesting to note that the two segments of data for which the simple model has the most difficulty in reproducing the minima occur ∼100 days apart; or in other words, at the same orbital geometry as defined by ω.

We adopted for the mass of the WR star MWR ≃ 20 M, determined from the helium-star mass-luminosity relation (Schaerer & Maeder 1992), where the luminosity was derived independently of distance from the energy needed to accelerate the WR wind (Schmutz 1997)1. For the companion we may follow the proposal of Skinner et al. (2002) that the companion is a low-mass star (LMS), and we adopt MLMS ≈ 1.5 M. The combined mass and the anomalistc period yield a separation of the two objects of a ≃ 0.13 AU. This value can be compared with the optically thick (electron scattering) radius of the WR star given in Fig. 3 of Schmutz (1997), RWR ≈ 5 R*, where R* = 3.5 R, which yields RWR ≈ 0.08 AU, or 0.6 a. This ratio compares favorably to the value of the width of the Gauss curve, w ≈ 0.6 to 0.65 a, which is needed to describe the eclipse light curve in Fig. 1.

The eclipse by the unseen object has about the same duration as that by the WR star. This is possible if the eclipsed part of the shock zone, that is, the brightest variable part, is very close to the other object such that the object covers about a quarter of the sky seen from the shock. This would be indeed the case if the other object were a low-mass star, as speculated above. The shock zone is a direct result from the collision of the wind with this star, or, somewhat speculatively, with an accretion disk or magnetopause.

In order to determine the orbital phases and orientation of the orbit during observations obtained in the past (e.g., those of XMM-Newton from 2001), a much more precise set of orbital and precession parameters is needed, which requires more observations. Moreover, the very fast change rate that we find in the apsidal motion, , cannot be sustained for very long. This also requires further investigation, as does the physical mechanism that could produce the fast apsidal motion. The primary assumption is that it is caused by a third body. A rough estimate of the minimum period of the outer body, using the equations given by Borkovits et al. (2015), neglecting tidal and relativistic effects, yields Pout ≈ 22 d if m3 ≃ mWR or 31 d if m3/m123 ≈ 1, where m3 and m123 are the mass of the third body and the total system mass, respectively (personal communication by the referee). Such a short period may be revealed by searching for a variable inner systemic velocity in the radial velocity curves of the WR star. Further investigation is also needed to refine the wind collision model that provides the periastron times and to search for the signatures of the associated emission lines in optical spectra.

6. Conclusions

The long-duration light curve of 136 days allowed us to identify the minima epochs and to detect a regular pattern that is characteristic for an apsidal motion. The timing measurements of the minima were used to determine six parameters of the orbit of the system. Based on these parameters, and using a simple Gauss-shaped eclipse model and a simple wind-collision model producing enhanced emission at periastron, we calculated a light curve that reproduces the observed complicated light curve shape over the full observational time-span. No doubt remains that EZ CMa is a binary.

Many eccentric binaries are known to exhibit an apsidal motion (Wolf et al. 2006, 2008), which in extreme cases is dominated by the perturbations of a third body. The high value of and the fact that instead of an advancing periapsis it is receding are very unusual, however. In principle, the presence of a third body could also explain the curvature in the O-C diagram by the light travel time in the outer orbit. If this were the case, then the radial velocity variations of the inner system would have to be greater than 100 km s−1, and we expect that this would have been noted in previous investigations.

Our findings have implications for the determination of the stellar parameters of the WR star because an analysis has to take into account that there is light from a shock zone and possibly also from other objects in the system. The binary scenario also changes the view of Morel et al. (1997) and St-Louis et al. (2018). What they interpreted as additional emission from corotating interaction regions are the brightening episodes of the shock zone. As illustrated in Fig. 1, the large majority of the observed light maxima are explained by the corresponding binary phase.

On a more general level, we note that the interpretation of the observed X-ray spectrum of EZ CMa by Skinner et al. (2002) has lead to the correct conclusion. They argued that the presence of hard X-rays implies the presence of a companion. When we extend this reasoning to other apparently single nitrogen-type WR stars with X-ray detections, then we find that all of them show an admixture of cool (kT1 <  1 keV) and hot (kT2 >  2 keV) X-ray plasma and prominent emission lines from ions such as Si XIII and S XV (Skinner et al. 2010). This indicates that they might have as yet undetected companions. The most extreme candidate after EZ CMa is WR 134. This in turn would mean that the binary scenario for forming WR stars is more important, if not dominant, and that current evolution models of massive single stars assume too large an amount of mass is lost by single massive stars in order to force the simulations to reach a WR stage.


1

The spectroscopic analyses of EZ CMa published by Schmutz (1997) adopted a distance of d ≃ 1.8 kpc as determined by Howarth & Schmutz (1995). ESA’s Gaia satellite measures π = 0.41 ± 0.05 mas, yielding d ≃ 2.4 kpc, and based on this, the system is about 20 to 50% more distant than assumed in his analysis. However, allowing for 15 to 20% of the light to come for the shock zone, and allowing another 10% for other objects in the system, we conclude that the correction needed for a larger distance is cancelled and that Schmutz (1997) assumed the correct order of magnitude for the absolute brightness of the WR component. Thus, the stellar parameters of the WR star given by Schmutz (1997) can be considered valid estimates. This holds in particular for the electron-scattering radius of the WR wind that is used in the discussion section.

Acknowledgments

We thank an anonymous referee for very constructive and helpful comments. Based on data collected by the BRITE Constellation satellite mission, designed, built, launched, operated and supported by the Austrian Research Promotion Agency (FFG), the University of Vienna, the Technical University of Graz, the University of Innsbruck, the Canadian Space Agency (CSA), the University of Toronto Institute for Aerospace Studies (UTIAS), the Foundation for Polish Science & Technology (FNiTP MNiSW), and National Science Centre (NCN). GK acknowledges support from CONACYT grant 252499.

References

All Tables

Table 1.

Orbital parameters derived from a double sine curve fit to the observed minima intervals (Fig. 2).

Table 2.

Orbital parameters derived from fitting the measurements of minima epoches shown in Fig. 3 and simultaneously the intervals shown in Fig. 2.

All Figures

thumbnail Fig. 1.

Copy of Fig. 1 published by Moffat et al. (2018), displaying the photometric light curve of EZ CMa as observed by the Toronto BRITE satellite from 2015 November 14 until 2016 March 29. The measured times of well-defined minima are indicated by green dots. The calculated times of upper and lower panels conjunctions of the orbit, with the parameters given in Table 2 that are described in the text, are marked by plus signs. The uneclipsed photometric light curve from a shocked zone, which is assumed to vary proportionally to the square of the orbital separation, is given by the orange curve. The light from the shock zone is eclipsed twice each orbit, and the light curve resulting from a simple Gauss-shaped eclipse model is shown by the blue line.

Open with DEXTER
In the text
thumbnail Fig. 2.

Measured intervals between consecutive minima. The red dots mark the differences in even number minus odd number minima, and green dots show the differences in odd minus even minima. The green and red lines are sine curve fits to the measurements. The blue crosses connected by the blue dotted line show the calculated intervals from the final orbital solution (see text).

Open with DEXTER
In the text
thumbnail Fig. 3.

Observed epochs of minima, marked with red dots for even number minima and with green dots for odd minima, minus computed minima epochs calculated with the relation T(E) = E0 + PE, with P = 1.88 d, which is half the sidereal period. E is the number count of the minima starting with zero. The dashed curve is a quadratic fit to all measured minima. The blue crosses mark the calculated minima from the final orbital solution relative to the linear ephemeris.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.