Long Secondary Periods in variable red giants Free

Abstract

We present a study of a sample of Large Magellanic Cloud red giants exhibiting Long Secondary Periods (LSPs). We use radial velocities obtained from VLT spectral observations and MACHO and OGLE light curves to examine properties of the stars and to evaluate models for the cause of LSPs. This sample is much larger than the combined previous studies of Hinkle et al. and Wood, Olivier & Kawaler.

Binary and pulsation models have enjoyed much support in recent years. Assuming stellar pulsation, we calculate from the velocity curves that the typical fractional radius change over an LSP cycle is greater than 30 per cent. This should lead to large changes in T_eff that are not observed. Also, the small light amplitude of these stars seems inconsistent with the radius amplitude. We conclude that pulsation is not a likely explanation for the LSPs. The main alternative, physical movement of the star – binary motion – also has severe problems. If the velocity variations are due to binary motion, the distribution of the angle of periastron in our large sample of stars has a probability of 1.4 × 10⁻³ that it comes from randomly aligned binary orbits. In addition, we calculate a typical companion mass of 0.09 M_⊙. Less than 1 per cent of low-mass main-sequence stars have companions near this mass (0.06–0.12 M_⊙) whereas ∼25–50 per cent of low-mass red giants end up with LSPs. We are unable to find a suitable model for the LSPs and conclude by listing their known properties.

1 INTRODUCTION

A subset of Long-Period Variable stars (LPVs) was found several decades ago to show a secondary period of variation, in addition to their primary pulsation (e.g. Payne-Gaposchkin 1954; Houk 1963). These Long Secondary Periods (LSPs) exceed the primary period in length by approximately one order of magnitude, a fact more recently confirmed by Wood et al. (1999). Wood et al. showed that LPVs fall on four distinct period-luminosity sequences A–D, and they found an additional sequence, E, of red giant binaries. These multiple sequences have been confirmed in subsequent studies, and a splitting of sequence B into two sequences has since been discovered (Ita et al. 2004; Soszyński et al. 2004a; Fraser et al. 2005; Soszynski et al. 2007). LSPs occupy sequence D, the sequence corresponding to variations with the longest period. The primary pulsation of these stars is usually found on sequence B (Wood et al. 1999).

Approximately 25–50 per cent of LPVs show an LSP (Wood et al. 1999; Percy et al. 2004; Soszyński et al. 2004b; Soszynski et al. 2007; Fraser, Hawley & Cook 2008). At present, there is no accepted explanation for LSPs. Wood et al. (1999) initially proposed a model wherein the sequence D stars are semi-detached red giant binaries, as the length of their periods are consistent with those expected for binary systems with solar-mass components. Several other models have been proposed to explain LSPs, including radial and non-radial pulsation, and dust effects. However, Wood et al. (2004) demonstrate that there are problems with all of these models. It is clear that further investigation is required in order to discover the cause of LSPs.

Here, we study a sample of sequence D stars in the Large Magellanic Cloud (LMC) in order to provide new constraints on models for the sequence D phenomenon. In particular, we present new spectra taken with the European Southern Observatory's Very Large Telescope (VLT). We derive velocity curves, radii and effective temperatures from these spectra and compare them in phase with light curves and associated quantities from the MACHO (Alcock et al. 1992) and OGLE (Udalski et al. 1992) data bases.

Our study follows the methods of Hinkle et al. (2002) and Wood et al. (2004), who also compared both light and radial velocity data for sequence D stars. However, we use a larger data set and have the advantage of high quality and simultaneous light and velocity data, which we hope will give a more complete and accurate picture of the behaviour of LSP variables. Since our sample comes from the LMC, we also have the advantage that the distances, and hence derived luminosities, are well determined.

2 OBSERVATIONS AND DATA REDUCTION

We chose a region of high star density in the middle of the LMC bar for this study, in order to get a large sample of sequence D stars in the 20 arcmin field of the FLAMES (Pasquini et al. 2002) multifibre system on the VLT. The field centre is at 05^h28′15″−69°45′43″ J2000. It contained a sample of 58 variable red giants exhibiting LSPs which were obtained from the MACHO Project data base. It was important to have simultaneous light and velocity data. Fortunately, OGLE light curves were taken at the times of our VLT observations. We have used both OGLE II and OGLE III I-band light curves in this study.

Spectra were obtained on 21 nights from 2003 November to 2006 March. The FLAMES/GIRAFFE spectrograph (Pasquini et al. 2002) with a grating setting of HR16 was used, giving spectra with a resolution of 23 900 and a spectral interval of 693.7–725.0 nm. This region includes the TiO bandhead at 705 nm, which can be used for spectral typing and to examine variation of T_eff, as well as derivation of radial velocities. An exposure time of 20 min was used for all spectra. The spectra were obtained in service mode. Unfortunately, no observations were taken between May and August when the LMC (RA =05^h23′) was only observable near twilight, leaving a long winter gap in 2004 and 2005.

The raw spectra of the 58 sequence D stars, most of which have 21 observations from different dates, were reduced using the FLAMES/GIRAFFE pipeline or the Swiss reduction pipeline. A number of sample spectra are shown in Fig. 1, where the TiO bandhead is clearly visible in the second and third panels. The majority of the program stars have spectra similar to those shown in the top three panels, i.e. they are oxygen-rich giants of spectral types M and K. However, the sample does contain a number of carbon-rich stars (C-stars), as shown in the fourth panel of Fig. 1.

Using the iraf software, each spectrum was visually examined and any obvious cosmic rays removed. There was a small percentage of bad data: in particular, spectra taken on two dates (2005 March 17 and 2005 September 10) had very few counts. This was due to the weather and seeing conditions at the time, and the spectra taken on these dates were discarded, leaving 19 good observations for most stars (some stars had fewer spectra due to the fibre allocation process in FLAMES). A small number of other spectra were discarded due to a low number of counts. These spectra did not fall on particular dates or particular target stars, and were distributed randomly throughout the data set. Presumably, the low counts were due to the fibres for these stars being misaligned.

Relative radial velocities were obtained through cross-correlation with the iraf task fxcor. From each star, a spectrum with a high number of counts and narrow lines was selected, and this acted as the template spectrum for the star's cross-correlation. The cross-correlation was performed in the wavelength region 6950–7160 Å, as this region is relatively free of telluric lines, as can be seen in the lower panel of Fig. 1. Errors in the radial velocities were also taken from fxcor. The mean velocity error of the sample is 0.42 km s⁻¹, and errors typically range between 0.1 and 0.6 km s⁻¹.

Because of uncertainty about the accuracy of the velocity calibration resulting from the pipeline reduction, the spectra were checked to see if any zero-point correction to the velocity was needed, before we calculated their absolute velocities. This procedure is now described.

2.1 Testing the relative telluric radial velocities

To check the velocity calibration, program stars were cross-correlated with a B star whose spectrum contained only telluric lines. Cross-correlation of a program spectrum to a telluric spectrum produces the relative telluric velocity of the two spectra, which should be zero if there are no calibration errors. This cross-correlation was performed in the wavelength region 7160–7220 Å, as it has a large concentration of telluric lines (see Fig. 1).

For most dates of observation, the relative telluric radial velocities were found not to be significantly different from zero, and so no systematic correction was applied to the spectra for those nights. However, observations from two nights, with Heliocentric Julian Dates (HJD) of 2453377 and 2453643, consistently showed a positive velocity offset from zero, so a zero-point correction for the velocities on these dates was included for all stars. The velocity corrections for these dates were −2.061 and −0.729 km s⁻¹, respectively.

2.2 The effect of telluric lines on relative radial velocities

The highest concentration of telluric lines was found in the wavelength region 7160–7220 Å, but there were a small number in the program object region (6950–7160 Å, see Fig. 1). There was a possibility that these lines would affect the calculated radial velocities, so this was examined. The telluric lines in the program object region of several program stars were removed using iraf's telluric command. Then the program spectra with telluric lines removed were cross-correlated against each other in the same way the original program spectra were. It was found that removing the telluric lines in this region had no significant effect on the velocities or errors, so telluric line removal was not performed.

2.3 Calculating the absolute radial velocity

Following the zero-point corrections, the absolute radial velocities of the program stars were calculated. First, the template of each program star was cross-correlated to a star of known radial velocity, giving the absolute radial velocity of each template. Then by adding a keyword ‘vhelio’ (which contained the absolute radial velocity) to the image header of the template, we could cross-correlate the spectra for each star with its template to obtain the absolute radial velocities.

For the O-rich program stars, the radial velocity reference used was α Cet. Its spectrum was taken using the echelle spectrograph (resolution 70 000) on the late 74-inch telescope at Mount Stromlo Observatory, Canberra, Australia. For the C-stars, the reference used was the C-star X Vel, for which a spectrum was also taken with the 74-inch. The radial velocity used for α Cet was −25.8 km s⁻¹ (taken from the Astronomical Almanac, Tettelbach & Holdaway 2004). The velocity for the spectrum of X Vel (−5.4 km s⁻¹) was obtained by cross-correlation with α Cet (even though the spectral types are different, there was a strong cross-correlation peak due to common metal lines).

Radial velocities as a function of heliocentric Julian date (HJD) for a sample of the sequence D stars are shown in Table 1. Stars are identified by their MACHO numbers. Dashes denote dates for which there is no spectrum. The full table can be viewed in the online version of this paper (see Supporting Information).

2.4 Correcting for the short-period variation

Because sequence D stars show two distinct modes of variability, it is necessary to correct the light and velocity data for the behaviour of the short (primary) period in order to accurately study the LSP. We need simultaneous velocity and light data for this purpose, so we use the OGLE I light curve in these procedures, even though it is missing the southern winter season.

First, the OGLE I light curve was interpolated and boxcar-smoothed over a time interval equal to the (primary) pulsation period. Once this LSP-only light curve was made, it was subtracted from the original light curve, leaving only the variations due to pulsation, di, or the ‘pulsation light curve’.

Similarly to the pulsation light curve, we obtained the pulsation velocity curve by making a binary fit to the velocity data using the period of the LSP. The deviations of the velocity data, dv, from this fit form the pulsation velocity curve. We emphasize that although a binary fit was used for the purpose of obtaining a smooth fit, binarity may not be the physical mechanism behind the LSP.

For radial pulsation, the light and velocity variations are shifted in phase relative to one another, so it is necessary to find the phase shift between the pulsation light and velocity curves in order to find the amplitude relation between the short-period light and velocity.

The phase correction was found by plotting the pulsation light versus pulsation velocity. By varying the phase shift of this plot, and judging by eye at which phase shift the pulsation light and velocity showed the clearest correlation, the best phase shift was found to be ∼0.25 of a cycle. This means that most often the I light minimum occurs close to the same phase as the mean rising radial velocity. Lebzelter, Kiss & Hinkle (2000) and Lebzelter & Hinkle (2002) find that the phase shift between the light and velocity variations of pulsating Semi-Regular Variables (SRVs) is 0.5. This means maximum radial velocity occurs at minimum light, and vice versa. The phase shift is discussed further below.

Next, the slope of the light–velocity amplitude correlation for pulsation was measured, to give the velocity correction. From examination of several stars with relatively large pulsation amplitudes, it was concluded that the velocity correction was ∼15 km s⁻¹ mag⁻¹ in I. Generally, the correction to the LSP velocity resulting from these processes is small, less than 0.6 km s⁻¹, as the median value of di is ∼0.04 mag.

In order to further investigate the phase relation between the velocity and light curves, we constructed a model semi-regular pulsator using the non-linear pulsation code described in Keller & Wood (2006). The model is a first overtone pulsator, appropriate for a sequence B variable. The model oscillations should be similar to the primary oscillation in the sequence D variables studied here, and the semi-regular variables studied by Lebzelter & Hinkle (2002). The model had M= 1.5 M_⊙, L= 3000 L_⊙, T_eff= 3740 K, helium abundance Y= 0.3 and metallicity Z= 0.004. M_V and M_I were computed from the M star model atmospheres of Houdashelt et al. (2000), assuming [Fe/H]=−0.5.

The model light and radial velocity curves are shown in Fig. 2. The graph also shows the time variation of T_eff and radius for comparison with observational estimates of these quantities in later sections.

The model shows that minimum light occurs between mean increasing radial velocity (as estimated here) and maximum radial velocity (as estimated by Lebzelter & Hinkle 2002). The model also predicts a ratio of velocity to I amplitude of 9.2 km s⁻¹ mag⁻¹, somewhat smaller than the estimate of 15 km s⁻¹ mag⁻¹ given here. Given the small size of the corrections to velocity that we have applied, the uncertainty in the relative phase and amplitude ratio of velocity and I amplitude will not significantly affect the resulting velocity curve of the sequence D variability.

2.5 Obtaining LSP parameters and narrowing the sample

Following the corrections described above, we were now in a position to analyse the light and velocity variations associated with the LSP. A binary fit was made to the corrected velocity data of the LSP using a fortran program, Fitall to obtain the parameters of the velocity curve. Fitall also made a Fourier series fit (with a frequency f and one harmonic 2f) to the MACHO M_B, M_R and OGLE I light curves, and to M_B−M_R. Some examples of the fits to the light and velocity data made by Fitall can be seen in Fig. 3. The complete set is available in the online version of this paper (see Supporting Information).

Due to the poor quality of the data for some stars, the original sample of 58 sequence D stars had to be reduced. Some stars had large velocity errors due to poor quality or few spectra or no OGLE light curve, so they were excluded from the sample. A few stars had LSPs near 365 d leading to a large unfilled gap in the velocity curve, and could not be analysed.

3 RESULTS

In this section, we discuss some general results that are not model specific. Model specific results are presented in Section 4.

3.1 The true distribution of the velocity amplitude

The distribution of velocity amplitude is an important parameter for most models of the sequence D phenomenon.

A histogram of the velocity amplitudes of our sample is shown in Fig. 4, with the values separated into bins of width 0.5 km s⁻¹. It is important to note that we use the full or peak-to-peak amplitude, as opposed to the semi-amplitude. These amplitudes are obtained from the binary orbit fits to the observed velocities. We find that the velocity amplitude is concentrated mainly between 2.5 and 5.0 km s⁻¹, with the highest number of stars falling in the bin centred at 3.5 km s⁻¹. The distribution has a median value of 3.53 km s⁻¹. Thus, our sequence D stars all have low, relatively similar velocity amplitudes. These results are consistent with the results of Hinkle et al. (2002)'s and Wood et al. (2004)'s analysis of Galactic LSP variables.

There are very few stars in our sample with velocity amplitudes <3.0 km s⁻¹. To find out if we could detect such stars we performed a Monte Carlo simulation. 50 stars were simulated for each of the velocity amplitudes of 3.0, 2.0, 1.0 and 0.5 km s⁻¹. The period, eccentricity and time of periastron were all selected randomly. The number of velocity points in each simulation was the same as the number of observations (19) and they were distributed in time identically to the observations. Random noise was added to each point, drawn from a normal distribution whose mean was the mean velocity error in the observations. In plots for visual examination, an error bar was added for each point, drawn from a uniform distribution between 0.2 and 0.5 km s⁻¹ (consistent with observed errors).

A binary fit was made to the simulated velocities in the same way as to the observed velocities (see Section 2.5). The simulated data were put through the same reduction steps as the observed data, to see what percentage would pass. The percentages can be seen in Fig. 5. This tells us that we should be able to detect velocity amplitudes <3.0 km s⁻¹. Using our simulated percentages, we can estimate the true number of stars at each velocity amplitude by dividing the observed histogram by the observable percentage. The resulting modified velocity amplitude histogram is shown overplotted on the observed histogram in Fig. 4.

The fact that we do not observe many sequence D stars with velocity amplitudes <3.0 km s⁻¹ suggests that there are not many to be observed. We conclude that the true distribution of sequence D velocity amplitudes, represented in Fig. 4, shows a strong peak near 3.5 km s⁻¹.

3.2 Correlation of light and velocity amplitudes

Light amplitude is plotted against velocity amplitude in Fig. 6. There is no obvious correlation between light and velocity amplitudes in either M_B, M_R or I. As noted previously the velocity amplitude has a small spread, and peaks at the small value of 3.5 km s⁻¹. In particular, stars with small light amplitudes have velocity amplitudes as large as those stars with the largest light amplitudes.

3.3 Correlation of velocity amplitude with magnitude and period

Velocity amplitude is plotted against K magnitude in Fig. 7, which shows a suggestion that velocity amplitude increases slowly with luminosity. A similar correlation is seen in Fig. 8 between velocity amplitude and LSP. However, these possible correlations are disrupted by the presence of a few stars with larger velocity amplitudes, that fall apart from the general relation.

4 TESTING AND CONSTRAINING MODELS

Here, we examine some of the more plausible and popular models of LSPs, in the light of our new data.

4.1 The binary model

The model of binary motion as an explanation for LSPs has been much debated since its proposal by Wood et al. (1999). In its favour, binary motion may be able to explain the sequence D velocity curves, and the light variation was originally thought by Wood et al. to resemble an eclipse of the star by a dusty cloud surrounding a small, orbiting companion. However, when this model was further examined by Hinkle et al. (2002) and Wood et al. (2004), two main problems were found: the companion mass was calculated to be always ∼0.1 M_⊙, and the characteristic shape of sequence D stars' velocity curves was shown to imply that the angle of periastron is not uniformly distributed between 0 and 2π. With our larger sample we can place tighter constraints on these properties.

4.1.1 The companion mass

In a binary system, the velocity amplitude and period give an estimate of the mass of the companion. For a typical sequence D star with an LSP of 500 d, a typical velocity amplitude of 3.5 km s⁻¹, and an assumed total system mass of M= 1.5 M_⊙, the orbital separation is of the order of 1.4 au and the companion has a mass of 0.09 M_⊙. This means it is a very low mass main-sequence star or a brown dwarf. There is an observed deficit of brown dwarf companions to main-sequence stars, compared to both more massive stellar companions and less massive planetary companions. This is known as the ‘brown dwarf desert’ (McCarthy & Zuckerman 2004; Grether & Lineweaver 2006). The range of likely companion masses can be found by calculating the companion mass for velocity amplitudes within 1σ of the mean of 3.5 km s⁻¹. We find a companion mass range of 0.06 to 0.12 M_⊙. Using the fits to the data in fig. 8 of Grether & Lineweaver (2006), we find that only 0.86 per cent of low-mass main-sequence stars have companions with masses in this range. Thus, binarity is an unlikely model, when we remember that ∼30 per cent of low-mass stars will exhibit LSPs when they pass through the AGB stage.

4.1.2 Correlation of velocity amplitude with light amplitude, magnitude and period

We would not expect the velocity and light amplitudes to be related in a binary model where the light variation is due to an eclipse phenomenon (by a cloud in the Wood et al. 1999 model). However, these amplitudes may be related if the light variations were caused by distortion of the red giant by its unseen companion (see Section 4.1.5).

As noted in Section 3.3, only a hint of a correlation is seen between velocity amplitude and K magnitude (Fig. 7), and with LSP (Fig. 8). We cannot draw any conclusions from this in the context of a binary model.

4.1.3 The distribution of the angle of periastron

One of the most conclusive ways to test the plausibility of the binary model as a cause for the LSPs is to examine the distribution of the angle of periastron, ω. As Hinkle et al. (2002) and Wood et al. (2004) note, the observed distribution of ω should be uniform over the whole range of angles, since one would expect binary orbits to be randomly aligned in space.

The angle of periastron was calculated by the program Fitall as part of the binary fit made to our radial velocity data. The distribution of the angle of periastron is plotted as a histogram in Fig. 9, with the data divided into bins 20° in width. Note that since the errors in the angle of periastron are typically around 5° (Table 2), the distribution of angle of periastron shown in Fig. 9 will not be significantly broadened by these errors. There is a clear bias in the distribution towards angles >180°, and the median of the sample is 227°. The distribution of the angles of periastron calculated by Hinkle et al. (2002) and Wood et al. (2004) for their small samples of Galactic stars has been overplotted. It is clear that their values of ω lie in the same region as the bulk of the present sample. This shows that the LSPs share a common cause, whether found in the LMC or in our own Galaxy.

Hinkle et al. (2002) and Wood et al. (2004) both found that the characteristic shape of the velocity curves suggests an eccentric orbit and a large angle of periastron, and that the distribution of the calculated angle of periastron was inconsistent with what would be expected if the sequence D stars were binaries, as no satisfactory explanation could be found for this bias towards higher angles. But given their small sample size, their conclusions were not highly significant. We can quantify and greatly improve these claims by using a Kolmogorov–Smirnov, or K–S Test, with our new data.

We have used the one-sample K–S Test, in the form of the fortran 77 subroutine ksone (provided in Press et al. 1986), to compare our distribution of the angle of periastron to the uniform distribution. We find that the probability that our distribution of ω is consistent with the uniform distribution is 1.4 × 10⁻³. In other words, the probability that the sequence D stars are binaries is extremely small. This is a major result from this study.

An angle of periastron between 180° and 360° means that at periastron, the red giant is closest to the observer, with the smaller companion further away. It is hard to see what selection effects could cause an LSP to be observed only in these binary orientations.

4.1.4 The distribution of the eccentricity

A histogram showing the distribution of the eccentricity is shown in Fig. 10. The eccentricity was calculated by Fitall as one of the orbital parameters of the binary fit to the velocity data.

The plot shows that if the sequence D stars are caused by binary motion, then they are in eccentric orbits. This confirms the suggestion of Hinkle et al. (2002) and Wood et al. (2004). The median eccentricity of the sample is 0.3.

4.1.5 Ellipsoidal variables

An ellipsoidal variable is a binary star which is tidally distorted by an orbiting companion into an ellipsoidal shape. Light variability is caused mostly by the variation of the apparent surface area of the star seen by the observer as the star rotates. Its velocity curve is dominated by orbital motion, with a small contribution from the rotation of the limb-darkened ellipsoid.

The light curve for an ellipsoidal variable is expected to have two maxima and minima per orbit, due to the star's shape. However, the velocity curve is expected to have only one maximum and minimum per orbit. We can use this characteristic to investigate the plausibility of ellipsoidal variability as an explanation for the LSP variation.

Soszyński et al. (2004b) show that the light curves of sequence E variables are satisfactorily explained by ellipsoidal variability. Their data show a partial overlap of sequences D and E in the period-luminosity diagram, and so they suggest that the LSPs may have a binary ellipsoidal origin. Soszyński (2007) explores this possibility further by searching for ellipsoidal or eclipsing shapes in residual LSP light curves. He finds these variations in ∼5 per cent of his sample, and adopts the binary model to explain LSPs.

However, with our new data set (which includes some sequence E stars), the difference between sequence D and sequence E mechanisms can now be vividly demonstrated by the difference between the light and velocity curves of sequence D and E stars. The phased-up light and velocity curves of sequence E stars show that the light curve completes two full cycles to every single cycle of the radial velocity, as expected for an ellipsoidal variable (Adams, Wood & Cioni 2006; Nicholls, Wood & Cioni, in preparation). However, Fig. 3 shows that sequence D stars do not display this behaviour: the phased-up light and velocity curves match each other cycle for cycle. Furthermore, sequence E stars typically have velocity amplitudes of 30 km s⁻¹ or more while the sequence D stars have significantly lower typical amplitudes of 3.5 km s⁻¹. Soszyński (2007) predicted that sequence D stars showing residual ellipsoidal or eclipsing-type variations should have large velocity amplitudes, similar to those observed for sequence E stars. One of our stars (77.7671.282) was noted by Soszyński as a double-humped LSPV and its velocity amplitude is only ∼4 km s⁻¹. Additionally, none of our LSP sample shows typical sequence E velocity amplitudes (see Fig. 4). Sequences D and E variables clearly have distinctly different velocity amplitudes.

There is an unambiguous distinction between the mechanisms of variability responsible for causing sequences D and E. We can now state with some certainty that binary star ellipsoidal variability is not the cause of the LSP.

Similar to the sequence E ellipsoidal variables are symbiotic binaries, which are red giants with an accreting white dwarf companion, often in an eccentric orbit. Radial velocity amplitudes of these systems are usually ≥8 km s⁻¹, and the red giant may be a Mira variable (Hinkle et al. 2006). However, symbiotic spectra show high temperature emission lines from accretion on to the white dwarf, something that has not been seen in either the sequence E or sequence D spectra. There does not seem to be any relation between the symbiotic stars and the sequence D stars.

4.2 Radial pulsation

There are several ways in which the present data set can be used to test radial pulsation models.

4.2.1 The phase difference between light and velocity curves

Minimum light and the mean of the velocity during its increase are relatively well defined for our sequence D stars. In order to examine their relative phases, a histogram of the phase difference between minimum light and mean rising velocity is shown in Fig. 11. The phases of minimum light and mean rising velocity were obtained from the fits made to the light and velocity curves. Each panel shows a strong peak of phase differences lying between 0 and 0.2 and there is little spread in phase.

Lebzelter et al. (2000) and Lebzelter & Hinkle (2002) find that for radially pulsating SRVs, the light lags the radial velocity by roughly 0.5 in phase, which means that minimum light occurs ∼0.25 in phase after the time of mean rising velocity. Our model SRV (Fig. 2) shows that minimum light occurs close to 0.23 in phase after mean rising velocity. Our sequence D stars are significantly different. Taking the main peak of light-velocity phase differences in Fig. 11 (values between −0.1 and 0.3), the mean phase shift is 0.1, and the standard deviation is 0.08. More simply, minimum light usually occurs ∼0.1 in phase after mean rising velocity. Radial pulsation therefore may not be consistent with this aspect of the LSP phenomenon, although this is not a strong conclusion.

4.2.2 The shape of the velocity curve

During large amplitude radial stellar pulsation (as in Mira variables), a shock wave from the interior of the star causes the surface layers to rapidly accelerate outwards, and then settle back more slowly. This is also shown in the shape of the velocity curve of our low amplitude model SRV, Fig. 2. The velocity curves of most of our sequence D stars show exactly the opposite of this behaviour: typically, we find a rapid increase in radial velocity with a correspondingly slower decline (see Fig. 3). In a radial pulsation model, this translates to a slow increase in stellar radius and a quick decrease. This has not been observed in any red giant stars known to be radially pulsating variables (e.g. Hinkle, Scharlach & Hall 1984; Lebzelter 1999, for Mira variables and SR variables, respectively).

4.2.3 Correlation of velocity amplitude with light amplitude, magnitude and period

The lack of correlation between velocity amplitude and light amplitude for the stars in our sample (see Fig. 6) is an unexpected result for pulsation mechanisms as they would be expected to show a positive correlation between light and velocity amplitude (Hinkle, Lebzelter & Scharlach 1997). Additionally, correlations of velocity amplitude with luminosity and with period are expected for radially pulsating variables. However, Figs 7 and 8 show that any correlations between these quantities in our sample are vague at the best. We are unable to draw any conclusion from this in the context of radial pulsation.

4.2.4 The variation of stellar radius

The stellar radius can be calculated from our data in several ways. Our first method was to calculate R using the Stefan–Boltzmann equation L= 4πR²σT⁴_eff, with L and T_eff calculated from MACHO photometry (see Appendix).

Examples of the variation of stellar radius calculated from the photometry, R_phot, can be seen in Fig. 12. There is generally a clear periodic curve at the period of the LSP as well as a scatter due to the primary pulsation. Plots for the entire sample are available in the online version of this paper (see Supporting Information).

Variation of stellar radius can also be calculated from the velocity, using where v is the fit to the velocity data made by Fitall. If we assume that the velocity changes are due to radial pulsation, then this gives the expected radius variation for radial pulsation, R_vel, where the resulting radius has been normalized to the median photometric radius. Here, the velocity fit data were scaled by a factor of 1.4 to convert from observed to true pulsation velocity in red giants (Scholz & Wood 2000; Wood et al. 2004). The radius change computed from the velocity fit is shown in Fig. 12. Again it shows a clear periodic curve.

It is also possible to calculate stellar radius from the Stefan–Boltzmann equation using an effective temperature derived from the spectra, as opposed to the photometry (spectral T_eff derivation is described in Section 4.2.5). The luminosities at the times of the spectra were in this case calculated from the fits to M_B, M_R and M_B−M_R as described in the Appendix. The radius calculated in this manner, R_spec, can also be seen in Fig. 12.

A histogram comparing the median values of R_phot and R_spec is shown in Fig. 13. For both radius estimates, most stars have median radii lying between 100 and 200 R_⊙. The median values for the entire sample for R_phot and R_spec are 135.4 and 132.1 R_⊙. R_vel is not shown, as it is normalized to R_phot.

The radial amplitudes of our stars mostly lie between 3 and 60 R_⊙ for all three variations. The median amplitudes of the sample for ΔR_phot, ΔR_spec and ΔR_vel are 5.92, 5.62 and 40.02 R_⊙, respectively.

Fig. 14 shows the distribution of radial amplitude ΔR divided by period, P. This quantity is a more sensitive indicator of differences between stars and models than the amplitude itself. In a radially pulsating model, one would expect the amplitude, ΔR, to be tightly correlated to the period if the velocity amplitudes are the same i.e. ΔR/P would be a constant. Given that our stars show very similar velocity amplitudes (see Fig. 4), we should expect tight, and similar, peaks for all three radius variations.

From Fig. 14, it is clear that the distributions of ΔR_phot/P and ΔR_spec/P, are similar, but they are not the same as for ΔR_vel/P. ΔR_vel is computed assuming radial pulsation. Given that it appears that the three radius variations do not agree, the assumption that ΔR is due to radial pulsation of the star is unlikely.

Fig. 15 shows the distribution of ΔR divided by R. For R_phot and R_spec, the fractional radius change is usually 10 per cent or less. For R_vel, the fractional radius change is much larger, usually between 10 and 50 per cent. Once again, the assumption that the measured radial velocity variations are due to radial pulsation does not seem to be correct.

A histogram comparing the phases of the three radius variations to the phase of the M_R light curve is shown in Fig. 16. We compared the phase of minimum radius with the phase of minimum light, as these are both well-defined parts of the curves. In order to calculate these phases, a Fourier series with four terms was fit to the radius curves. For the light curve, the previously used Fourier fit was used (see Fig. 3). Fig. 16 shows that R_phot and R_spec both vary almost in phase with the light variation (phase differences generally ≤0.1). On the other hand, R_vel is more likely to be out of phase with the light (and consequently with the other radius variation estimates). Again, the assumption that the observed radius variations can be attributed to radial pulsation leads to inconsistencies.

Using our model SRV (see Section 2.4), minimum radius occurs about 0.7 in phase before minimum light (Fig. 2). Of our variations, R_vel is the most similar to this. R_phot and R_spec, whose phases agree well with each other, do not show the radius–light phase relation expected for radial pulsation.

Overall, the discrepancy between the variations of R_vel and the other two radius estimates suggests that radial pulsation is unlikely to be the cause of the LSPs.

4.2.5 The variation of effective temperature

The variation of T_eff was calculated from the photometry, from the spectra and from the velocity fit.

The variation in the photometric T_eff(T_eff(phot)), whose derivation is described in the Appendix, is shown in Fig. 12.

T_eff was also calculated using the depth of the TiO bandhead at 7054 Å(T_eff(spec)). The depth was calculated as the ratio of the mean count level of a 4 Å band longward of the bandhead to the mean count level of a 4 Å band shortward of the bandhead, taking into account the fact that the LMC has an average velocity of 270 km s⁻¹, leading to a longwards wavelength shift of 6.3 Å in the spectra.

In order to calculate the effective temperature from the TiO depth, we calculated the depth of the TiO bandhead for nine Wing Spectral Standard stars. Then a quadratic fit was made to the TiO depth-spectral type data. This quadratic fit was used to calculate the spectral types of our sequence D stars from their TiO depths, and the effective temperatures were then calculated using the Ridgway et al. (1980)T_eff–spectral type relation. An example of the typical variation of T_eff(spec) in these stars can be seen in Fig. 12.

Some of the warmer sequence D stars in our sample showed little to no evidence of TiO in their spectra, and so we were unable to include these in the T_eff(spec) calculation. However, good-quality velocity curves were not necessary for calculating T_eff(spec) so a number of stars previously discounted from the velocity studies were added in where consideration of radius and effective temperature was required. Overall there was a net increase in sample size from 30 to 37 (the stars added in were excluded only from the velocity amplitude, angle of periastron, eccentricity and velocity phase parts of this work).

We also calculated the effective temperature variation expected for radial pulsation, T_eff(vel). This was calculated from the radius change computed for pulsation in Section 4.2.4, the luminosity calculated from the fits to the photometry, and the Stefan–Boltzmann equation. The variation of T_eff(vel) can be seen in Fig. 12.

The distributions of the median values for each star of T_eff(phot) and T_eff(spec) are shown in Fig. 17. Most stars have median effective temperatures lying between 3300 and 3800 K. The median values for the entire sample for T_eff(phot) and T_eff(spec) are 3662.9 and 3725.2 K. T_eff(vel) is not shown, as it is calculated from R_vel, which is normalized to R_phot.

Fig. 18 shows the distribution of T_eff amplitude for the three variations. The amplitude distributions of T_eff(phot) and T_eff(spec) are both strongly peaked at lower values, but ΔT_eff(vel) shows a much wider distribution, and tends towards much larger amplitudes. The amplitudes of T_eff(phot) and T_eff(spec) mostly lie between 20 and 200 K, but T_eff(vel) amplitudes tend to lie between 200 and 800 K. The median amplitudes of the sample are 42.3, 72.56 and 582.78 K for ΔT_eff(phot), ΔT_eff(spec) and ΔT_eff(vel), respectively. T_eff(vel), like T_eff(phot), is a predicted temperature, whereas T_eff(spec) is measured from our data. T_eff(vel) has very large amplitudes in the context of radial pulsation. Adding this to the fact that T_eff(vel) does not agree with T_eff(spec) suggests that radial pulsation cannot explain the LSPs.

Fig. 19 shows the relation of light amplitude to T_eff amplitude for both T_eff(phot) and T_eff(spec). The graph shows that T_eff amplitude increases with increasing M_R amplitude, for both T_eff(phot) and T_eff(spec).

5 DISCUSSION AND CONCLUSIONS

5.1 The movement of the visible surface of a sequence D star during its LSP

Possibly the most remarkable feature of the sequence D stars is the large change in the radial distance to the visible surface of a sequence D star as it passes through one LSP. By integrating the radial velocity with time, we calculated in Section 4.2 a radius change for radial pulsation between 3 and 60 R_⊙, with the typical radius full amplitude being ∼41 R_⊙. The typical median radius of a sequence D star is around 135 R_⊙.

We can use this property to demonstrate just how unlikely radial pulsation is in these stars. An amplitude of 41 R_⊙ in a 135 R_⊙ radially pulsating star corresponds to a fractional radius change of over 30 per cent from minimum to maximum radius. These radial amplitudes are very large in the context of radial pulsation: fundamental-mode Mira variables pulsate with radius changes of the order of 50 R_⊙ (Ireland, Scholz & Wood 2004), but they have visible light amplitudes of ∼6 mag, while typical sequence D stars have light amplitudes of ΔM_R≤ 0.8 mag. It seems unlikely that such a modest light change could be associated with so large a radius change, suggesting that radial pulsation is not the cause of LSPs.

This comparative analysis backs up a huge problem with the radial pulsation model which was raised in section 4.2: the changes in radius would lead to changes in T_eff that are vastly greater than the directly observed changes from spectra or photometric colour.

In addition to the problems with radial pulsation raised by our own data, there are also other problems raised by Wood et al. (2004): the length of the primary period does not vary during the LSP as expected; and to date there is no known mode of radial pulsation that has the required periods.

Wood et al. (2004) proposed non-radial pulsation as a possible cause for the LSPs. They showed that the best match to the observed periods is given by the g mode with n= 2 and l= 2. In this case, the net apparent radial velocity arises from some parts of the visible surface approaching while other parts are receding. Thus at a given position on the stellar surface, the velocity amplitudes and overall radial motions will need to be even larger than in the radial pulsation case. Given the large amplitudes calculated in this work for radial pulsation, non-radial pulsation would lead to distortions of the star that appear far too large to be credible, especially since g modes have significant amplitudes only in radiative regions (Wood et al. 2004).

It would therefore seem that the only way to generate the observed change in radial distance to the visible surface of a sequence D star is to have the star itself move, i.e. it must be in a binary system. Binary models, however, are inconsistent with various observational results described above. Most particularly, the non-uniform distribution of the angle of periastron and the low companion masses of these stars appear incompatible with a binary system model. The problem of the companion mass might be overcome if the companions were initially planets in eccentric orbits whose mass has been increased via accretion of material from the red giant. This would mean that although the companions currently have masses in the Brown Dwarf range, they were not born as Brown Dwarfs. To remove the problem of the non-uniform distribution of angle of periastron, we would need to have an observational selection effect that rejects angles of periastron in the range 0°–180°, i.e. we would need an effect that masks the sequence D variation when periastron occurs with the smaller companion in front of the red giant. An effect which could produce this bias is difficult to imagine.

5.2 Attributes of the Long Secondary Period

As with past studies, we are unable to find definitive evidence for any model which explains the LSPs in sequence D stars. To conclude, we therefore list all the currently known properties of LSPs. Any proposed model for LSPs should be able to explain all of these attributes.

• Stars exhibiting LSPs occupy a clearly defined period-luminosity sequence.

• LSPs are of length ∼250–1400 d.

• LSP light variation is not regular and minima, in particular, vary in depth from cycle to cycle (see fig. 2 of Wood et al. 1999).

• The primary pulsation is visible in the light curve at all times throughout the LSP and the primary period does not significantly change with LSP phase (Wood et al. 2004).

• The ratio of LSP to pulsation period is ∼8–10. The shorter period variation lies usually on sequence B and is thought to be the first or second overtone radial pulsation (Wood et al. 1999; Ita et al. 2004; Fraser et al. 2005).

• Lack of spectroscopic line broadening in observations of sequence D variables indicates that any rotational velocities are ≤3 km s⁻¹ (Olivier & Wood 2003).

• The radial velocity curves show a characteristic shape: observed radial velocity increases quicker than it decreases.

• The velocity amplitudes cluster tightly around 3.5 km s⁻¹.

• The light–velocity phase shift for both the short period and the LSP is ∼0.25, with the light lagging the velocity. Minimum light is roughly aligned with mean rising velocity.

• The ratio of colour to light variation of the LSP and primary oscillations are similar (Wood et al. 2004). Similarly, Derekas et al. (2006) showed that the ratio of blue to red amplitude of the LSP was similar to the ratio for stellar pulsation, and somewhat different from that due to ellipsoidal light variations.

• There is no correlation between LSP light and velocity amplitude.

• Stellar radius tends to be between 100 and 200 R_⊙, with large radius variations of 3–60 R_⊙ (if it is assumed that the LSP is caused by radial pulsation).

• The equivalent width of the Hα absorption line varies with the LSP, indicating chromospheric activity (Wood et al. 2004).

PRW has been partially supported in this work by the Australian Research Council's Discovery Projects funding scheme (project number DP0663447). We are grateful for the multiple allocations of VLT service time over several semesters for this extended series of observations (program identifiers 072.D-0387, 074.D-0098, 075.D-0090 and 076.D-0162). This paper utilizes public domain data obtained by the MACHO Project, jointly funded by the US Department of Energy through the University of California, Lawrence Livermore National Laboratory under contract No. W-7405-Eng-48, by the National Science Foundation through the Center for Particle Astrophysics of the University of California under cooperative agreement AST-8809616, and by the Mount Stromlo and Siding Spring Observatory, part of the Australian National University.

REFERENCES

Appendix

APPENDIX A: CALCULATING LUMINOSITY AND EFFECTIVE TEMPERATURE FROM PHOTOMETRY

The luminosity was calculated from the MACHO Red photometry. The MACHO red and blue magnitudes were first dereddened using a colour excess of E(B−V) = 0.08 for the LMC (Keller & Wood 2006) and absorption A(λ) calculated by the reddening law from Cardelli, Clayton & Mathis (1989). This gave A(M_B) = 0.2531 and A(M_R) = 0.1783.

Supporting Information

Table 1. Heliocentric radial velocities and Julian dates for MACHO Sequence D stars.

Figure 3. Light and velocity curves plotted against phase for sequence D stars.

Figure 12. Different properties of sequence D stars plotted against Julian Date.

Please note: Wiley-Blackwell are not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the article.