Free Access
Issue
A&A
Volume 617, September 2018
Article Number A121
Number of page(s) 7
Section Stellar atmospheres
DOI https://doi.org/10.1051/0004-6361/201731614
Published online 28 September 2018

© ESO 2018

1. Introduction

The advent of photometers on board dedicated spacecrafts like Kepler, CoRoT, MOST, and BRITE revolutionized the field of astronomical photometry. The space-borne photometers provide not only light curves with unprecedented coverage, but they also lead to a significant improvement of the precision of measurements. As a result, stars that were deemed constant were found to be variable leading to discoveries of new types of variable stars.

Single O stars are a typical examples of stars for which no variability was expected, yet that still show low-amplitude light variations. For example, significant low-frequency variability was detected in Kepler photometry of a blue supergiant HD 188209 and was attributed to the gravity waves (Aerts et al. 2017). A dedicated CoRoT observing run focused on the young open cluster NGC 2244 detected stochastic variations in O stars, which were attributed to the red noise that is possibly due to the subsurface convection (Blomme et al. 2011) and multiple frequencies due to stellar oscillations (Degroote et al. 2010; Briquet et al. 2011; Mahy et al. 2011). Similar variations were also found after the subtraction of binary light curve in δ Ori Aa (Pablo et al. 2015), in early-B supergiant HD 2905 (Simón-Díaz et al. 2018), and in ζ Pup (Ramiaramanantsoa et al. 2018).

Low-amplitude light variability in O stars may be connected to the line-driven wind instability. As a result of the Doppler effect, wind driving by multiple lines in hot stars is unstable (Lucy & Solomon 1970; Owocki & Rybicki 1984). The numerical simulations have shown that the line-driven wind instability leads to the variability of the wind density, velocity, and mass-loss rate on a typical scale of hours (Owocki et al. 1988; Feldmeier et al. 1997). Due to supersonic nature of hot star winds, the wind instability leads also to appearance of wind shocks and the wind X-ray emission (Feldmeier et al. 1997). The line-driven wind instability is likely connected with clumping, which affects observed wind spectral features (e.g. Sundqvist et al. 2010, 2011; Šurlan et al. 2012, 2013; Shenar et al. 2015) and therefore affects the wind mass-loss rate estimates.

Wind instabilities are self-initiating (Owocki et al. 1988), but may be also initiated and modulated by photospheric motions (Feldmeier et al. 1997). Consequently, the instabilities may be connected with photospheric turbulence whose presence affects, for example, the widths of photospheric lines (Aerts et al. 2009; Cantiello et al. 2009; Jiang et al. 2015).

Part of the flux emerging from the photospheres of stars is absorbed by the wind and emitted back to the stellar photosphere. The backwards emission causes heating of the photosphere especially at low optical depths and also leads to redistribution of the flux mainly from the short-wavelength part of spectrum to longer wavelengths. This effect is termed wind blanketing (Abbott & Hummer 1985) and is important for obtaining precise effective temperatures of hot stars with strong winds (Bohannan et al. 1986; Crowther et al. 2002; Martins et al. 2005). The wind blanketing redistributes the emergent flux in dependence of wind opacity, and amount of redistribution depends on the wind mass-loss rate. Therefore, wind blanketing does not lead to any light variability for fixed mass-loss rate.

On the other hand, if the mass-loss rate is variable, then the wind blanketing is modulated by the mass-loss rate, resulting in the photometric variability. For example, in magnetic hot stars the wind mass flux depends on the tilt of the magnetic field (Owocki & ud-Doula 2004). Therefore, the strength of the wind blanketing varies across the stellar surface, which, due to stellar rotation, leads to the photometric variability of magnetic O stars. This effect can partly explain the photometric light curve of HD 191612 (Krtička 2016).

The line-driven wind instability is another process that leads to the mass-loss rate variation (Owocki et al. 1988; Feldmeier et al. 1997; Runacres & Owocki 2002). Therefore, also the line-driven wind instability causes photometric variability in hot stars. This can be connected with the low-amplitude variability observed in O stars. To study this possibility, we used hydrodynamical simulation of Feldmeier et al. (1997) to derive the mass-loss rate variation and the relation between the mass-loss rate and photometric flux (Krtička 2016) to predict the light variability.

2. Modelling of the light variability: spherically symmetric case

The calculation of the light variability of hot stars due to the line-driven wind instability is a very complex problem. In general, 3D time-dependent hydrodynamical simulations of the stellar wind are required to obtain the density and velocity structure of the circumstellar environment of hot stars and its variability. Such simulations have to be coupled to a global (unified) 3D wind model that includes the influence of the wind on the stellar photosphere (accounting for wind blanketing) and that is able to solve the radiative transfer equation for moving media whose ionization and excitation state is out of equilibrium. There is no such code that would be able to treat the whole complex problem. Therefore, we used existing wind simulations and global wind models and combine their output to derive the expected light variations in our stars. Due to the simplifications involved, we have not focussed on a particular star or stellar type, but instead we aim to obtain general magnitude of the effect and its properties.

2.1. Hydrodynamical simulations

The hydrodynamical simulations we have employed here were described by Feldmeier et al. (1997) in detail. The simulations solved continuity, momentum, and energy equation for a spherically symmetric non-rotating wind flow. The simulations were performed using the smooth source function method (Owocki et al. 1991; Owocki & Puls 1996). The hydrodynamical part of the code is based on a standard van Leer solver. The simulations use staggered mesh, operator splitting of advection and source terms, advection terms in a conservative form using van Leer’s (1977) monotonic derivative as an optimised compromise between stability and accuracy, Richtmyer artificial viscosity, and non-reflecting boundary conditions (Hedstrom 1979; Thompson 1987, 1990). A turbulent variation of the velocity at a level of roughly one third of the sound speed was introduced as seed perturbation for unstable growth at the wind base.

The hydrodynamical simulations were performed for specific stellar parameters corresponding to ζ Ori A, that is the effective temperature Teff = 31 500 K, radius R* = 24 R, and mass M = 34 M, and yielded the mean mass-loss rate . The simulations are robust against changes in stellar parameters. Consequently, the relative variations of wind density derived for such specific parameters are also applicable to other O stars.

The calculation of the wind blanketing directly from the hydrodynamic simulations would require solution of radiative transfer and kinetic equilibrium equations (also called NLTE equations) for non-monotonic velocity law. To our knowledge, there is no such code that is able to treat this problem in its full complexity. Therefore, to make the problem tractable, we calculate a mean value of the mass-loss rate in the region close to the star, where the wind blanketing originates,(1)

where ρ and v are the wind density and radial velocity in particular time. Most of the flux in the optical domain is emitted from the region below the sonic point, which typically appears at radii 1.01 − 1.02 R* in the numerical simulations. Therefore, we selected r1 = 1.01 R* and r2 = 1.02 R*. However, the final light curve does not significantly depend on the choice of r1 and r2. As a result of this choice, in Eq. (1) we average over 50 grid points of hydrodynamical simulation. A plot of is given in Fig. 1.

thumbnail Fig. 1.

Mass-loss rate from the numerical simulations averaged over r ∈ [1.01 R*, 1.02 R*] (see Eq. (1)).

Open with DEXTER

2.2. Global wind models

The variations of the emergent flux with wind mass-loss rate were derived from spherically symmetric stationary METUJE global wind models (Krtička & Kubát 2017). The code calculates the wind model in a global (unified) approach that integrates the description of the hydrostatic photosphere and supersonic wind. The model solves the radiative transfer equation in the co-moving frame (CMF) with opacities and emissivities calculated using occupation numbers derived from the kinetic equilibrium (NLTE) equations. For given stellar parameters, the model enables us to predict the radial dependence of the density, velocity, and temperature from hydrodynamical equations and to derive the wind mass-loss rate and the stellar emergent flux that accounts for the wind blanketing.

For our modelling we assumed a typical parameters of O stars corresponding to HD 191612, Teff = 36 000 K, R* = 14.1 R, and M = 29.2 M. To model the dependence of the emergent flux on the mass-loss rate with fixed stellar parameters, we artificially scaled the radiative force in the wind. This yields a series of wind models and emergent fluxes parameterized by the wind mass-loss rate. Because the amount of the flux redistributed by the wind depends on the effective temperature, we expect slightly different dependence for stars with different temperature. However, because we are concerned with typical light variations, we limited our study to just one star. On the other hand, such selection affects only the amplitude of the light curve. The preliminary results derived for 300-1 model of Krtička & Kubát (2017, the model has parameters close to that used in hydrodynamical simulations) indicate that the amplitude of the variability due to the wind blanketing does not significantly depend on stellar parameters.

The light variability is derived from a relation between the Hp magnitudes and mass-loss rate (Krtička 2016, Eq. (2))(2)

where ΔF = 1.6 × 107 erg s−1 cm−2, F0 = 9.0 × 108 erg s−1 cm−2 is the emergent flux for a reference mass-loss rate, and is the mean mass-loss rate. The relation was derived for HD 191612 parameters using METUJE global wind models (Krtička & Kubát 2017).

2.3. Calculated light curve

The light curve calculated using Eq. (2) for the mass-loss rate variations plotted in Fig. 1 is given in Fig. 2. The amplitude of the light curve is of the order of 0.01 mag and the typical timescale of the light variability (hours) is the same as the timescale of the variability of the mass-loss rates as derived from the simulations.

The amplitude of the light curve is proportional to the mass-loss rate variations, which depend on the surface perturbations. The light variations are therefore stronger for stronger base perturbations. The wind structure is triggered by the photospheric perturbations, which therefore determine the time scale of the variability.

thumbnail Fig. 2.

Light curve calculated using Eq. (2) for the mass-loss rate variations plotted in Fig. 1.

Open with DEXTER

3. The case without spherical symmetry

Within a spherical symmetry we assumed that mass-loss rate variations are coherent across the stellar surface. This assumption is highly unrealistic for real wind variations. However, the coherence length in the horizontal direction is unclear. A typical patch size of the instability generated structure is of the order of degrees, as inferred from observations (e.g. Dessart & Owocki 2002) and the typical number of clumps in the observed part of the wind is of the order of 103–105 (Davies et al. 2007; Oskinova et al. 2007; Nazé et al. 2013; Šurlan et al. 2013).

To account for the deviation from spherical symmetry, we followed the approach of Dessart & Owocki (2002) and Feldmeier et al. (2003) and assumed that the stellar wind consists of N concentric cones. We assumed that each of these cones is independent, and that the mass-loss rate from the ith cone with stellar surface cross-section at time t is(3)

where Δti is random for each i and is given by Eq. (1). Moreover, we assumed that varies periodically with the period given by the time extent of available simulations (which is about 1.5 d, see Fig. 1).

The cones were selected in such a way that Ωi is roughly the same for all cones. We divided the visible stellar surface into elements bounded by a specified number of concentric rings. The first element that directly faces the observer is assumed to specify just one cone. The number of cones specified by other rings is selected in such a way that solid angle set by the cones is roughly equal to solid angle of the first cone. The number of cones is therefore a function of the number of rings. Because we specify the number of latitudinal rings, we cannot calculate the light curves for an arbitrary number of cones.

The magnitude difference between the observed flux fHp(t) at a given time and the reference flux in passband Hp is defined as(4)

The reference flux is obtained under the condition that the mean magnitude difference over the simulated light curve is zero. The radiative flux in a passband Hp at the distance D from the spherical star is (Hubeny & Mihalas 2014)(5)

The intensity IHp(θ, Ω, t), at angle θ with respect to the normal to the surface was obtained at each surface point with spherical coordinates Ω from the emergent flux taking into account the limb darkening uHp(θ)(6)

with limb darkening coefficients from Reeve & Howarth (2016). Here . The emergent flux is obtained from Eq. (2) modified for fluxes, and the mass-loss rate is given by Eq. (3) for a cone that corresponds to given coordinates Ω.

The resulting light curves are given in Fig. 3. From this figure it follows that the amplitude of the light variability due to wind instability quickly decreases with the number of cones. On the other hand, even with a relatively large number of cones, the light variability should be detectable with satellite photometry.

thumbnail Fig. 3.

Light curve calculated using Eq. (4) for different number N of concentric cones.

Open with DEXTER

4. Comparison with observations

We used Period04 (Lenz & Breger 2005) to analyse the simulated light curves (see Fig. 4). The general shape of the Fourier spectrum does not significantly depend on the number of cones used to calculate the light curve. The Fourier spectrum shows maximum at the frequency of about 3 d−1 which corresponds to the characteristic times scale of the variability of the mass-loss rate. The location of the maximum peak is nearly independent of the number of assumed cones. There are further peaks at higher frequencies. An additional low-frequency peak at about 1 d−1 is most likely connected with the length of available simulations.

Due to a stochastic nature of the light curves, a direct comparison of the simulated and observed light curves is not meaningful. On the other hand, it is possible to compute the Fourier spectrum of the light curves, fit the spectrum by some phenomenological model, and compare the fit parameters derived for observed and simulated light curves. Similar techniques are used to study time series in various contexts, for example, light curves of X-ray binaries (Burderi et al. 1993), solar flare data (Threlfall et al. 2017), or to search for granulation in Cepheids (Derekas et al. 2017).

To compare the theoretical light curves with observations, we fitted the Fourier spectrum A(ν) in Fig. 4 by a polynomial(7)

The fit parameters given in Table 1 demonstrate that the shape of the Fourier spectrum (given by a parameter) does not depend on the assumed number of cones, while the amplitude of peaks (described by b) depends on N. The fit of results from Table 1 shows that the parameter b varies with N on average as(8)

The simulated light curve can be compared, for example, with the light curve of ϵ Ori (HD 37128, B0Ia) derived using the BRITE satellites. The BRITE-Constellation consists of six satellites working in blue and red domains of visible spectra (Weiss et al. 2014). The satellites provide high-precision photometry of bright objects including chemically peculiar stars, pulsating stars, and Be stars (Baade et al. 2016; Weiss et al. 2016; Daszyńska-Daszkiewicz et al. 2017; Handler et al. 2017). The data were downloaded from the BRITE public data archive1. We used observations of ϵ Ori obtained by BAb (BRITE-Austria) satellite between JD 2456628 and 2456734. We used observations in a blue filter which covers the wavelength range of 390–460 nm (Kaiser et al. 2008). The observed differential light curve shows complex variability with amplitude comparable to the simulated light variations in Fig. 2 and a typical timescale of about 1 d (see Fig. 5 and also David-Uraz et al. 2017a). A corresponding peak appears also in the Fourier spectrum in Fig. 6. Consequently, while the time scale of the variability is comparable to that derived from simulations, the interpretation of the amplitude of the observed variability would require either a coherent wind variability or a very strong base perturbation.

The inspection of the BRITE archive revealed that a similar type of variability is also likely present in two additional stars ζ Ori A (O9.2Ib, HD 37742, HR 1948) and J Pup (B0.5Ib, HD 64760, HR 3090). The light variations of these stars have roughly the same amplitude and the time scale as light variations of ϵ Ori (see also Buysschaert et al. 2017). The star ζ Ori A has weak magnetic field and the rotational period of about 6.8 d (Blazère et al. 2015) and therefore it enables us to study the interaction of atmospheric turbulent motions and magnetic field. Both stars show the line profile variability due to corotating interacting regions (Howarth et al. 1998; Kaper et al. 1999), which are likely unrelated to large-scale magnetic fields (David-Uraz et al. 2014). Because the line-driven wind instability is possibly present in all hot star winds being one of the sources of their X-ray emission (e.g. Antokhin et al. 2008; Nazé 2009), the light variations due to to corotating interacting regions and wind instabilities coexist.

The open cluster NGC 2244 observational run of the CoRoT satellite provides even more promising results. The photometric study of the O9V star HD 46202 (Briquet et al. 2011) revealed multiple frequencies in the range 2−5 d−1 and amplitude of the order of 0.01 mmag. This can be explained as a result of wind blanketing due to multiple small-scale wind perturbations. O stars HD 46150, HD 46223, and HD 46966 show small amplitude stochastic variations (Blomme et al. 2011, see also Fig. 7) that can be connected to the subphotospheric convection (Aerts & Rogers 2015). Comparing with Fig. 3 follows that such variability can be caused by a relatively large number of wind perturbations and wind blanketing.

This conclusion is further supported by the comparison of the Fourier spectra (Fig. 8) and their fit parameters in Table 1. The slope parameters a of the spectra based on predicted light curves are close to the theoretical results. From Eq. (8) it follows that the observed light variability is likely caused by a very large number of cones N ∼ 102 − 104, which is consistent with results of numerical simulations (Sundqvist et al. 2018). Numerical tests using the uncertainties of observed data showed that the Fourier spectrum in Fig. 8 is at least one order of magnitude higher than that of observational noise. Therefore, the Fourier spectra in Fig. 8 should correspond to stellar variations.

thumbnail Fig. 4.

Fourier analysis of the light curves given in Fig. 3.

Open with DEXTER
Table 1.

Parameters of the fit Eq. (7) of the Fourier spectrum.

thumbnail Fig. 5.

Portion of the light curve of ϵ Ori derived using the BRITE satellite. Difference between the actual and mean magnitudes is plotted.

Open with DEXTER
thumbnail Fig. 6.

Fourier spectrum of the light curve of ϵ Ori.

Open with DEXTER
thumbnail Fig. 7.

Portion of the light curve of O stars derived using the CoRoT satellite. Difference between the actual and mean magnitudes is plotted.

Open with DEXTER
thumbnail Fig. 8.

Comparison of the Fourier spectra of light curves given in Fig. 7 (black curves) with Fourier spectrum of the theoretical light curve for N = 1 (red curves). The theoretical spectrum was scaled by a factor of 10b taken from Table 1 accounting for N > 1. For the analysis we used the CoRoT light curves secured between JD 2456203 and 2456234 for HD 46150 and HD 46223 and between JD 2454748 and 2454783 for HD 46966.

Open with DEXTER

5. Discussion and conclusions

We simulated the light variability of O stars due to variable wind blanketing modulated by line-driven wind instability. We used the output from hydrodynamical wind simulations and combined the derived time dependence of the mass-loss rate with the dependence of optical flux on the wind mass-loss rate derived from global wind models.

The resulting light curve has amplitude of the order of hundredths of magnitude with a typical time scale of several hours. The amplitude is of the order of millimagnitudes (and lower) when assuming that the wind consists of concentric cones in which the wind behaves independently. Such variability is still observable using high precision photometry.

We compared the derived light curve with light curve of ϵ Ori obtained using the BRITE satellite and with light curves of O stars from open cluster NGC 2244 obtained by CoRoT satellite. The observed light curves show stochastic variability with amplitude and time scale comparable to that derived from simulations. The light curves of NGC 2244 stars yield Fourier spectra that have the same shape as the Fourier spectra of the predicted light curves and the amplitudes that imply presence of large number (N ∼ 102 − 104) of independent surface patches. The shape of the Fourier spectrum is slightly different for ϵ Ori, but this may be connected with a narrow range of frequencies available for analysis.

The variability due to the wind blanketing is another source of the light variability in O stars, which contributes to the light variability due to the pulsation and convection (Aerts & Rogers 2015; Buysschaert et al. 2015). Consequently, it may be problematic to distinguish the wind blanketing variability from other effects. Ultraviolet observations may solve this problem and test the proposed variability mechanism, because the light variability in the ultraviolet (e.g. around 1300 Å) should be in the anti-phase with visual variability (Krtička 2016). Moreover, there may be stars that are variable only because of the wind blanketing, because only perturbations of surface density and velocity (and no variations of the frequency-integrated radiative flux) are required to trigger this type of variability.

The proposed mechanism requires relatively large mass-loss rate variations. Consequently, we do not expect strong light variations connected with corotating interaction regions (Cranmer & Owocki 1996; Krtička et al. 2004; Lobel & Blomme 2008). These regions have been invoked to explain the discrete absorption components commonly observed in the ultraviolet lines of hot stars. The discrete absorption components are connected with velocity plateaus in the wind and require order of magnitude lower mass-loss rate variations than those invoked by wind instabilities (David-Uraz et al. 2017b). Therefore, from Eq. (2) we expect up to millimagnitude variations due to corotating interaction regions.

The light variability due to the wind blanketing is another piece of puzzle to the general picture of variability in O stars. The subsurface convective motions trigger surface light variations connected with pulsations and variable wind blanketing leading to perturbations that disseminate in the wind initiating the line-driven wind instability.


Acknowledgments

JK acknowledges support by grant GA ČR 16-01116S. 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 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). The CoRoT space mission, launched on 2006 December 27, was developed and operated by the CNES, with participation of the Science Programs of ESA, ESA’s RSSD, Austria, Belgium, Brazil, Germany and Spain.

References

All Tables

Table 1.

Parameters of the fit Eq. (7) of the Fourier spectrum.

All Figures

thumbnail Fig. 1.

Mass-loss rate from the numerical simulations averaged over r ∈ [1.01 R*, 1.02 R*] (see Eq. (1)).

Open with DEXTER
In the text
thumbnail Fig. 2.

Light curve calculated using Eq. (2) for the mass-loss rate variations plotted in Fig. 1.

Open with DEXTER
In the text
thumbnail Fig. 3.

Light curve calculated using Eq. (4) for different number N of concentric cones.

Open with DEXTER
In the text
thumbnail Fig. 4.

Fourier analysis of the light curves given in Fig. 3.

Open with DEXTER
In the text
thumbnail Fig. 5.

Portion of the light curve of ϵ Ori derived using the BRITE satellite. Difference between the actual and mean magnitudes is plotted.

Open with DEXTER
In the text
thumbnail Fig. 6.

Fourier spectrum of the light curve of ϵ Ori.

Open with DEXTER
In the text
thumbnail Fig. 7.

Portion of the light curve of O stars derived using the CoRoT satellite. Difference between the actual and mean magnitudes is plotted.

Open with DEXTER
In the text
thumbnail Fig. 8.

Comparison of the Fourier spectra of light curves given in Fig. 7 (black curves) with Fourier spectrum of the theoretical light curve for N = 1 (red curves). The theoretical spectrum was scaled by a factor of 10b taken from Table 1 accounting for N > 1. For the analysis we used the CoRoT light curves secured between JD 2456203 and 2456234 for HD 46150 and HD 46223 and between JD 2454748 and 2454783 for HD 46966.

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.