© ESO 2018

1 Introduction

Low- to intermediate-mass stars evolve to red giant and asymptotic giant branch (AGB) stars. Mass loss increases during the AGB evolution. The AGB mass loss is driven by an interplay between pulsations, which extend the atmosphere, dust formation in the extended atmosphere, and radiation pressure on the dust (e.g. Wachter et al. 2002; Mattsson et al. 2010; Höfner & Olofsson 2018). However, details of these interrelated processes are a matter of debate, in particular regarding the atmospheric levitation and wind acceleration for oxygen-rich stars (e.g. Bladh et al. 2015; Höfner et al. 2016; Bladh et al. 2017). These processes are better understood for carbon-rich AGB stars (e.g. Nowotny et al. 2010, 2011; Rau et al. 2015, 2017).

Dynamic 1D model atmospheres of oxygen-rich Mira stars based on self-excited radial pulsation models (CODEX models) by Ireland et al. (2008, 2011) predicted a regular sinusoidal variation of the photospheric radius and an irregular chaotic variability of the outer molecular layers. The latest 3D radiation hydrodynamic (RHD) simulations of AGB stars (CO5BOLD) by Freytag & Höfner (2008) and Freytag et al. (2017) show non-radial structures such as long-lasting giant convection cells and short-lived surface granules. These dynamical phenomena trigger large-scale atmospheric shock waves that expand roughly spherically and are similar to those of 1D models, except they do not cover the full surface at a given instance.

Thompson et al. (2002) conducted long-term narrow-band interferometric monitoring over 2–3 stellar cycles of the oxygen-rich Mira S Lac and the carbon-rich Mira RZ Peg. Their data shows the expected sinusoidal variation in the continuum angular radius and different phase lags in the continuum and molecular bands. For the oxygen-rich Mira S Lac, the continuum minimum size tracked the visual maximum brightness, i.e. the continuum size and the visual lightcurve were anti-correlated. For the carbon-rich Mira RZ Peg, the phase lag was 0.28. A few other interferometric studies included fewer epochs and confirmed the presence of variability in the continuum size and in molecular layers (e.g. Woodruff et al. 2009; Haubois et al. 2015; Wittkowski et al. 2016).

Here we report on measurements of the variability in the continuum radius and in extended molecular layers using the example of the oxygen-rich Mira star R Peg obtained during science verification and early science operations of the newly available near-infrared K-band beam combiner GRAVITY (Gravity Collaboration 2017) at the VLT Interferometer (VLTI).

2 Observations and data reduction

We used the GRAVITY instrument during science verification and early science operation periods to perform high-precision multi-epoch spectro-interferometry of the Mira variable R Peg.

R Peg is an oxygen-rich Mira variable of spectral type M6–M9 with a mean period of 378 d (Samus’ et al. 2017). We adopt the parallax to R Peg from the Gaia DR 2 (Gaia Collaboration 2016, 2018) of ϖ = 2.8300 ± 0.2544 mas. The relative error of f = 0.09 is quite small, so we can simply estimate the distance by inverting the parallax (Luri et al. 2018). This distance is $\rho ={353}_{-29}^{+35}$ pc. Figure A.1 shows the recent visual lightcurve of R Peg based on data obtained from the AAVSO (American Association of Variable Star Observers) and AFOEV (Association Francaise des Observateurs d’Etoiles Variables) databases. Table A.1 provides the log of our observations. We fit a sine curve to the ten most recent cycles of the visual lightcurve relative to our observations (JDs 2453800–2457600) to estimate a current period of 377.1 ± 0.2 d and to assign visual phases of 0.15, 0.36, 0.45, and 0.53 to our four epochs of observation in June, September, October, and November 2016. R Peg transitioned from post-maximum to minimum phases during our four epochs. We used the compact baseline configuration for all observations, where the 1.8 m auxiliary telescopes (ATs) were located at stations A0, B2, C1, and D0, giving ground baseline lengths between 11.3 m and 33.9 m, and projected baseline lengths between 9.8 m and 32.0 m. We bracketed observations of R Peg by observations of interferometric calibrators in sequences of calibrator-science-calibrator observations. Information on the calibrators and their adopted angular diameters is available in Tables A.1 and A.2. The spectral range of our GRAVITY observations included the full near-infrared K band between 1.99 μm and 2.45 μm at a spectral resolution of R ~ 4000 (high spectral resolution mode). The GRAVITY instrument allows us to choose between splitting the light into two polarization angles (split polarization mode), which increases the internal fringe contrast, or to use the combined light (combined polarization mode), which increases the sensitivity of the instrument. As our target is bright and sensitivity was not an issue, we chose to use the split polarization mode. The detector integration time (DIT) on the science spectrometer was 5 s for the observations in June and September, and 10 s for those in October and November. On each night the same value was used for R Peg and its calibrators. The fringe tracker was operated with a DIT of 0.0085 sec for all observations.

We reduced and calibrated the data with the latest release of the GRAVITY pipeline (version 1.10.11)¹ and its Reflex workflows gravity-wk and gravity-viscal, respectively. Each observation was executed as a sequence between object and sky positions (object-sky-object-sky), and we computed averages of the object and sky files. The pipeline gives results for each of the two polarization directions separately, and we averaged the two results during post-processing. There were no significant differences between them. Fringe tracker data, obtained at a low spectral resolution of five spectral channels across the K-band, were processed in the same way.

As a reference, all calibrated visibility spectra are shown in Figs. A.2–A.5, together with the model fits as described below. The visibility spectra show the typical shapes of oxygen-rich Mira stars as observed previously with the VLTI-AMBER instrument (e.g. Wittkowski et al. 2011, 2016). The visibility spectra show a maximum at the near-continuum bandpass around 2.25 μm, corresponding to a minimum angular diameter, where the bandpass is only slightly contaminated by molecular layers. Toward shorter and longer wavelengths the visibility drops due to the presence of extended molecular layers; in the K band most importantly it drops due to H₂O and CO. The GRAVITY data show an increased precision with respect to the AMBER data; provide the data for six baselines in one snapshot, while the AMBER data provide three; and provide a spectral resolution of ~4000 across the full K band compared to a resolution of ~1500 for AMBER across two separatehalves of the K band. The closure phase spectra are not an essential part of the analysis of this work which focuses on overall diameter variations based on visibility data in the first lobe. For reasons of completeness, we show the closure phase data in Figs. A.6–A.9.

The data obtained on the science spectrometer mostly agree well with those obtained on the fast fringe tracker, which confirms a high accuracy of the absolute visibility calibration. In a few visibility spectra from September and October, the fringe tracker results deviate from the science spectrometer results for some of the baselines, especially toward the red part of the spectrum for some baselines. Based on the shape of the visibility spectra and the consistency of the science camera results across different baselines, we attribute these differences to a problem with the absolute calibration of the fringe tracker data at these epochs, which does not affect the quality of the fringe tracking itself and whose origin is not yet known. Some visibility spectra show a dip centred at 2.00 μm, close to the blue edge of the wavelength range, which is an artefact and most likely caused by low flux at these wavelengths owing to absorption by water, or due to a contamination from the metrology laser operating at 1.91 μm.

3 Data analysis

We are interested in the changes of the photospheric radius and that of extended molecular layers across our four epochs transitioning from post-maximum to minimum visual phases. As a model-independent analysis, we first provide fits of a uniform disc (UD) model at a few selected bandpasses before discussing a comparison of our data with dynamic model atmosphere predictions. The intensity profile at the photospheric layer in a near-continuum bandpass is expected to be well described by a UD (e.g. Wittkowski et al. 2016, Fig. 3b), so that we chose this simple geometrical model. The intensity profile in bandpasses corresponding to molecular layers may be more complex, and may typically resemble geometrical profiles with two or more components. With our strategy of snapshot observations at different epochs, we are not able to constrain the true shape of the intensity profile at these wavelengths. This would require time series of image reconstructions, which is much more expensive in terms of observing time. Nevertheless, the deviation of a uniform disc diameter in a molecular band compared to that in a near-continuum band can be used as a simple measure of the overall contribution of the molecular layer to the intensity profile, which includes effects of both its radial extension and its relative intensity level. Followed over different epochs, the UD analysis gives information about the time variability of the molecular contribution. We chose observations in the first lobe of the visibility function which are sensitive to first-order structure, i.e. the overall radial extent, and much less sensitive to higher orderstructures such as the limb-darkening effect or even surface inhomogeneities caused by convection cells.

Fig. 1 Variability in the V magnitude and in the UD angular diameters in the near-continuum band, and in bands dominated by H2O (2.05 μm) and CO 2–0 (2.29 μm). Also shown are sinusoidal fit results. The minimum continuum size tracks the maximum light, which can be understood by the increasing effective temperature when the star gets smaller. The minimum contribution of H2O tracks the maximum light as well, which can be understood as the destruction of water vapour at maximum light and formation at minimum light. The contribution from CO instead is greatest at maximum light, indicating a different and more stable behaviour of CO compared to water vapour.

3.1 Uniform disc angular diameters

We obtained at each epoch UD diameters at the near-continuum bandpass at 2.25 μm to characterize the photospheric continuum radius. The bandpass at 2.25 μm is expected to be almost free of molecular contamination for oxygen-rich Mira stars (Wittkowski et al. 2008, Fig. 4). We chose a bandpass at 2.05 μm to monitor the strength of water vapour layers, excluding the instrumental artificial dip in some visibility functions at 2.00 μm mentioned above. Finally, we chose narrow bandpasses at the lowest points of the visibility drop in the CO (2–0) bandhead at 2.29 μm and in the CO (3–1) bandhead at 2.32 μm.

Table 1 shows the UD fit results at these bandpasses for each of our epochs. We adopt errors of ~1% including calibration uncertainties, while the formal errors are much lower. The synthetic visibility values of the UD fits are represented by blue dots in Figs. A.2–A.9. They show a good consistency between the visibility spectra corresponding to the different baselines, which again demonstrates the high accuracy of the absolute visibility calibration. Figure 1 plots the resulting UD diameters compared to the contemporaneous visual lightcurve, together with sinusoidal fits, as described below.

Uniform disc angular diameters at the near-continuum bandpass steadily increase with values of 8.7 ± 0.1 mas, 9.4 ± 0.1 mas, 9.8 ± 0.1 mas, and 9.9 ± 0.1 mas from post-maximum to minimum visual phases of 0.15, 0.36, 0,45, and 0.53, respectively. This means that the continuum radius is anti-correlated to the visuallightcurve across our epochs. This behaviour is consistent with the multi-epoch narrow-band interferometric study of the oxygen-rich Mira S Lac by Thompson et al. (2002), which – to our knowledge – represents the only other such study in the literature.

Uniform disc diameters at the water vapour bandpass around 2.05 μm follow the variability of the near-continuum radii at larger UD diameters between 10.7 mas and 11.7 mas. UD diameters at the CO 2–0 bandhead, instead, are correlated with the visual lightcurve and anti-correlated with the near-continuum UD diameters with values between 12.3 mas and 11.7 mas.

We computed best-fit sinusoidal functions to the UD diameters in the different bandpasses. The results are listed in Table 2 and plotted in Fig. 1. Although the measured points lie well on the sinusoidal curves, the fit is not well constrained due to the low number of epochs. More epochs are needed to study the intra- and inter-cycle variability in more detail. Nevertheless, this fit provides us with a characterization of the phase lags with respect to the visual lightcurves, and with estimates of the mean values and amplitudes of the diameter and lightcurve variation. The continuum diameter lags behind the visual lightcurve by 0.57. The water vapour diameter similarly lags behind the visual lightcurve by 0.60. The CO (2–0) and CO (2–0) diameters have different phase lags of –0.05 and 0.10. Thompson et al. (2002) derived similar phase lags for S Lac of 0.48 in the continuum and 0.53 in a water vapour bandpass at 2.0 μm. Their bandpass at 2.4 μm had a phase lag of 0.34, but is mostly dominated by water vapour as well, and may not track the CO bandheads.

Phenomenologically, our observed phase lags of sizes or contributions in continuum, water vapour, and CO bands with respect to the visual lightcurve can be interpreted as follows. When the continuum radius decreases, the effective temperature rises, and this increase in effective temperature dominates the lightcurve ($L\propto R{T}_{\text{eff}}^4$). At the same time, water vapour is an unstable molecule that can get destroyed at maximum effective temperature and luminosity, and can be formed at minimum effective temperature and luminosity (e.g. Fig. 9 in Bladh et al. 2013); its contribution is thus anti-correlated with the lightcurve and correlated with the continuum size. CO instead is a more stable molecule whose abundance may be more stable across the stellar cycle, and which may become more strongly excited at maximum light than at minimum light. Time series of images in water and CO bands may help to further investigate the different behaviour of these molecules as a function of visual phase, and to disentangle variations in geometrical size and in flux contribution.

3.2 Fundamental stellar parameters

We used the mean continuum angular diameter of 9.24 ± 0.09 mas together with our adopted distance and the bolometric flux to derive the mean radius, effective temperature, and luminosity of R Peg. We estimated a bolometric flux of 1.08 ± 0.05 10⁻⁹ W m⁻² using broad-band UBV RIJHK (Ducati 2002) and IRAS (Helou & Walker 1988) photometry, and using an A_V value of 0.20 mag from Whitelock et al. (2008). We obtained a mean continuum radius of 353${}_{-29}^{+35}$ R_⊙, a mean effective temperature of 2480 ± 40 K, and a mean luminosity of 4200${}_{-960}^{+1320}$ L_⊙. Assuming a mass of 1 ± 0.3M_⊙, this corresponds to a surface gravity of −0.65 ± 0.25. Table 3 provides an overview of these stellar parameters.

3.3 Comparison to one-dimensional CODEX dynamic model atmospheres

Ireland et al. (2008, 2011) presented dynamic model atmospheres for Mira stars that are based on self-excited pulsation models, and use the opacity sampling method to represent gas opacities. CODEX models are available in four different series that correspond to different stellar parameters of the underlying hypothetical non-pulsating parent star, covering radii between 209 R_⊙ and 201 R_⊙, luminosities between 5050 L_⊙ and 8160 L_⊙, effective temperatures between 2860 K and 3400 K, and periods between 307 d and 427 d, meant to represent the Miras o Cet, R Leo, and R Cas. R Peg has a larger radius and lower effective temperature compared to the parameters of the CODEX model series. However, general properties of the variability of the CODEX models may be compared to our observations.

The CODEX models predicted in general a regular sinusoidal variation in the photospheric radius and a more irregular variability in the outer molecular layers. We inspected the variabilities of the luminosity and the continuum radii as a function of model phase based on the tables provided by Ireland et al. (2011). The CODEX model phases match the visual lightcurves within 0.1. The CODEX model series show the highest luminosities at model phases between –0.1 and 0.2, and the lowest luminosities at model phases between 0.54 and 0.70, depending on the model series. They show the smallest continuum radii at model phases –0.1 and the largest continuum radii at model phases 0.3 to 0.6.

As a result, it is a typical characteristic of the CODEX models that the continuum radius increases between maximum and minimum light, when the luminosity decreases. This is consistent with our observations of R Peg and with the observations of S Lac by Thompson et al. (2002), which both showed an anti-correlation between the visual lightcurve and the continuum radius. A closer inspection showed that the amplitudes of the radius variations of the CODEX models lie between 45% and 67%, while we estimated a clearly lower observed amplitude of 14% (Tab. 2). For comparison, the S Lac observations by Thompson et al. (2002) showed an amplitude of 19%. Moreover, the shape of the radius variation in the CODEX models is steeper than the observations. The model amplitude varies from pulsation cycle to pulsation cycle; most of this variation occurs near minimum light when the models are most extended and coolest, and hence most uncertain (e.g. Figs. 11–13 of Ireland et al. 2011). However, even in the more regular part of the radius variation (phases about 0–0.25), the model amplitude is still larger than is the observed variation in R Peg. It is not clear whether these differences between models and observations are caused by the different stellar parameters, or whether they indicate intrinsic differences in the models compared to observations.

As a result of these differences between observations and models, we cannot find a set of CODEX models that reproduces well the details of the observed variability of the continuum radius. We made an attempt to find a set of CODEX models that fits our observations by treating the model continuum diameter as a free parameter for each observation individually, effectively scaling the amplitude of the model variations to the observed ones. We selected a set of models that (1) matches as closely as possible the post-maximum to minimum visual phases of the observed epochs, and (2) reproduces the observed wavelength dependence of the visibility spectra and its variability. Table 4 lists the parameters of the identified CODEX model set that gives the best compromise between these criteria. Figures A.2–A.9 show the synthetic visibility spectra compared to the observed spectra. This set of models consists of adjacent post-maximum to minimum visual model phases as observed, reproduces the increasing continuum diameter at these phases, and provides a satisfactory fit to the wavelength dependent shape of the visibility spectra at the different epochs. The largest discrepancy between observations and model is seen for epoch 1, where the model phase (0.3) differs most with respect to the observed phase (0.15), and where the model still predicts deeper CO bandheads than observed.

3.4 Comparison to three-dimensional CO5BOLD model atmospheres

Recently, Freytag et al. (2017) presented 3D radiation-hydrodynamics models of AGB stars. They used different model parameters covering radii between 294 R_⊙ and 531 R_⊙, luminosities between 4978 L_⊙ and 10028 L_⊙, effective temperature between 2506 K and 2827 K, and periods between 338 d and 820 d. Models st28gm05n001 or st29gm06n001 from Table 1 of Freytag et al. (2017) are closest to the estimated stellar parameters of R Peg. We inspected variations in luminosity and radius of these 3D models, and found that the variability of both quantities is significantly more irregular than for the 1D CODEX models. This behaviour is expected because of the presence of additional non-radial structures such as long-lasting giant convection cells and short-lived surface granules. The model includes episodes where luminosity and radius are correlated, and episodes where these quantities are anti-correlated. The maximum full amplitude of the radius variations in these models is 14% and 8% for models st28gm05n001 or st29gm06n001, respectively, while it is lower during individual cycles. These numbers are comparable to our observed radius variation of R Peg of 14%. We chose not to attempt a direct comparison between our observations and predictions by 3D models because our present observations cover only four epochs and the 3D models show this strong irregularity. Observations at more phases are needed to attempt a meaningful comparison of observations with 3D model predictions of the stellar variability.

4 Summary and conclusions

We have reported on measurements of the variability in the continuum radius and in extended molecular layers at four epochs between post-maximum and minimum visual phases for the example of the oxygen-rich Mira star RPeg. To the best of our knowledge, our study represents the first characterization of the phase lags between Mira angular sizes in continuum and molecular bands and the lightcurve since the pioneering observations in this field by Thompson et al. (2002). We show that the continuum size and the size in a bandpass that is dominated by water vapour are anti-correlated with the visual lightcurve. This behaviour is consistent with the result by Thompson et al. (2002) for S Lac and with predictions by CODEX dynamic model atmospheres. The size in the CO (2–0) instead follows the visual lightcurve more closely, indicating a different and more stable – behaviour of CO compared to that of water vapour. The full amplitude of the near-continuum diameter variation of R Peg of 14% is comparable to the 19% for S Lac measured by Thompson et al. (2002). This amplitude is smaller than predicted by CODEX model atmospheres (45–67%), and closer to those predicted by 3D RHD simulations by Freytag et al. (2017) of up to 14%. The wavelength-dependent visibility variations at our epochs can be reproduced by a set of CODEX models at post-maximum to minimum visual model phases between 0.3 and 0.6. We show that the newly available GRAVITY instrument is well suited to conducting such high-precision time series of stellar sizes. Longer time series of these observations, and eventually of interferometric imaging, are needed to provide strong constraints on the model-predicted intra- and inter-cycle variability, and the interplay between pulsation, convection, and mass loss.

Acknowledgements

Michael Scholz was a co-I of this study. He unfortunately passed away before seeing the final results. We are grateful for many years of his inspiring work, dedicated collaboration, and friendship. GR thanks ESO for supporting her visit at ESO during her work on the GRAVITY data. This research has made use of the Jean-Marie Mariotti Center SearchCal service (http://www.jmmc.fr/searchcal) co-developed by LAGRANGE and IPAG; of CDS Astronomical Databases AFOEV, SIMBAD, and VIZIER (http://cdsweb.u-strasbg.fr) and of NASA’s Astrophysics Data System. We acknowledge with thanks the variable star observations from the AAVSO International Database contributed by observers worldwide and used in this research. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

Appendix A Additional material

Here we provide additional material. Figure A.1 shows the visual lightcurve of R Peg based on the AAVSO and AFOEV databases for the last few cycles together with a sine fit to the ten most recent cycles. Table A.1 details the log of our R Peg observations and Table A.2 provides the adopted properties of the calibrators. Figures A.2–A.9 show all the obtained visibility and closure phase data. The plots of the visibility also include predictions by a series of CODEX dynamic model atmospheres that best reproduce the observations.

Fig. A.1 Visual lightcurve of R Peg based on the AAVSO and AFOEV databases. The solid line represents a sine fit to the ten most recent cycles. The dashed vertical lines denote the epochs of our observations.

Fig. A.2 R Peg visibility data obtained in June 2016. The black lines denote the squared visibility amplitudes measured on the science camera. The errors are indicated by mean squared visibility amplitudes and their errors on each side of the plot for the first and second half of the wavelength interval. The red points with error bars denote the squared visibility amplitudes measured on the fringe tracker camera with low spectral resolution. The blue dots denote synthetic values based on UD fits at wavelengths of 2.05 μm (water), 2.25μm (near-continuum), 2.29 μm (CO 2-0), and 2.32 μm (CO 3-1). The green lines denote a fit of a CODEX model atmosphere with the parameters listed in Table 4.

Fig. A.3 As in Fig. A.2, but for the data obtained in September.

Fig. A.4 As in Fig. A.2, but for the data obtained in October.

Fig. A.5 As in Fig. A.2, but for the data obtained in November.

Fig. A.6 R Peg closure phase obtained in June 2016. The black lines denote the closure phases measured on the science camera. The errors are indicated by mean values and their errors on each side of the plot for the first and second half of the wavelength interval. The red points with error bars denote the closure phases measured on the fringe tracker camera with low spectral resolution.

Fig. A.7 As in Fig. A.6, but for the data obtained in September.

Fig. A.8 As in Fig. A.6, but for the data obtained in October.

Fig. A.9 As Fig. A.6, but for the data obtained in November.

References

All Tables

All Figures

Fig. 1 Variability in the V magnitude and in the UD angular diameters in the near-continuum band, and in bands dominated by H2O (2.05 μm) and CO 2–0 (2.29 μm). Also shown are sinusoidal fit results. The minimum continuum size tracks the maximum light, which can be understood by the increasing effective temperature when the star gets smaller. The minimum contribution of H2O tracks the maximum light as well, which can be understood as the destruction of water vapour at maximum light and formation at minimum light. The contribution from CO instead is greatest at maximum light, indicating a different and more stable behaviour of CO compared to water vapour.

In the text

	Fig. A.1 Visual lightcurve of R Peg based on the AAVSO and AFOEV databases. The solid line represents a sine fit to the ten most recent cycles. The dashed vertical lines denote the epochs of our observations.
In the text

Fig. A.2 R Peg visibility data obtained in June 2016. The black lines denote the squared visibility amplitudes measured on the science camera. The errors are indicated by mean squared visibility amplitudes and their errors on each side of the plot for the first and second half of the wavelength interval. The red points with error bars denote the squared visibility amplitudes measured on the fringe tracker camera with low spectral resolution. The blue dots denote synthetic values based on UD fits at wavelengths of 2.05 μm (water), 2.25μm (near-continuum), 2.29 μm (CO 2-0), and 2.32 μm (CO 3-1). The green lines denote a fit of a CODEX model atmosphere with the parameters listed in Table 4.

In the text

	Fig. A.3 As in Fig. A.2, but for the data obtained in September.
In the text

	Fig. A.4 As in Fig. A.2, but for the data obtained in October.
In the text

	Fig. A.5 As in Fig. A.2, but for the data obtained in November.
In the text

	Fig. A.6 R Peg closure phase obtained in June 2016. The black lines denote the closure phases measured on the science camera. The errors are indicated by mean values and their errors on each side of the plot for the first and second half of the wavelength interval. The red points with error bars denote the closure phases measured on the fringe tracker camera with low spectral resolution.
In the text

	Fig. A.7 As in Fig. A.6, but for the data obtained in September.
In the text

	Fig. A.8 As in Fig. A.6, but for the data obtained in October.
In the text

	Fig. A.9 As Fig. A.6, but for the data obtained in November.
In the text