EDP Sciences
Free Access
Issue
A&A
Volume 590, June 2016
Article Number A12
Number of page(s) 14
Section Stellar structure and evolution
DOI https://doi.org/10.1051/0004-6361/201527873
Published online 28 April 2016

© ESO, 2016

1. Introduction

Wolf-Rayet (WR) stars of spectral class WNE (nitrogen rich) are very hot, highly luminous stars that are thought to be the bare cores of evolved massive stars. Owing to strong stellar winds, these stars have lost almost all their hydrogen-rich envelope allowing us to model them as H-free helium stars (Chiosi & Maeder 1986; Langer 1989).

Such stars are expected to develop a convective core, due to the high central energy production, and a radiative envelope (Kippenhahn & Weigert 1990). In addition, the temperature in the outer layers decreases so much that another convective region arises because of a bump in opacity at log (T) ≅ 5.3, which is associated with iron and iron-group elements (FeCZ; Langer et al. 1994; Iglesias & Rogers 1996; Gräfener et al. 2012a; Langer 2012). As the opacity increases, the increased radiative acceleration brings the layers close to their local Eddington luminosity and, to avoid exceeding this luminosity, the star is forced to expand forming a low density inflated convective envelope (Petrovic et al. 2006; Gräfener et al. 2012a; Sanyal et al. 2015). This convective region is inefficient in transporting energy, but may lead to observable effects, such as turbulence in the surface layers, and may directly influence the wind at its base (Owocki & Rybicki 1986; Heger & Langer 1996; Glatzel 2005; Moffat 2008; Cantiello et al. 2009; Langer 2012; Gräfener & Vink 2013; Grassitelli et al. 2015b).

The mass-loss rate in massive stars is a crucial ingredient in determining correct stellar parameters, the evolutionary scenario, and consequently the fate of these stars (Langer 2012). In the WR phase the mass-loss rate and wind density are sufficiently high for the stellar winds to be optically thick, obscuring the hydrostatic layers from direct observations.

The spectra of WR stars of type WNE are dominated by broad emission lines of helium and, in part, nitrogen. Spectroscopic analysis of line profiles reveal the presence of systematic variability in the form of emission sub-peaks migrating from the line centre to the wings (Moffat et al. 1988; Lépine & Moffat 1999). Phenomenological investigations of these discrete wind emission elements involve the presence of a high number of randomly distributed stochastic density inhomogeneities, i.e. clumps in the wind of WR stars. The non-periodical nature of these radially propagating clumps and the migration of the sub-peaks have partially ruled out non-radial pulsation as the origin of these inhomogeneities (Cranmer & Owocki 1996; Lépine & Moffat 1999).

These structures in the outflow dynamically arise from instabilities that may be triggered at the base of the wind (Owocki 2015). Given the stochastic nature of this phenomenon, we investigate the conditions in which clumping can be triggered by sub-surface convection and can eventually be enhanced during the radial propagation in the wind by radiative instabilities (Owocki et al. 1988; Owocki 2015). In fact, observational evidence tends to show that clumps are formed at the base of the wind and do not develop when the wind is already significantly supersonic (Davies et al. 2005; Vink 2014; Owocki 2015). The phenomenon called clumping affects the parameters derived from model atmosphere analysis introducing further uncertainties in the observational analysis (Moffat & Robert 1994; Hillier 2008). Therefore an understanding of the origin of the observed variability and instabilities in the outer layers of these stars is fundamental for correctly estimating mass-loss rates and for understanding the physics of radiation pressure dominated winds (Hamann & Koesterke 1998).

As a result of the buoyancy force the convective elements move upwards and downwards and, in some cases, they may reach velocities of the order of the local sound speed. The interaction with the material in the upper radiative layer may consequently generate gravity waves that propagate to the surface, potentially leading to observable phenomena such as turbulence or small-scale velocity fields, as in the case of micro- and macroturbulence (Goldreich & Kumar 1990; Cantiello et al. 2009; Grassitelli et al. 2015a,b). The appearance of a velocity field at the surface is connected with the outermost region of the sub-surface convective zone and the physical conditions in the low density envelope. However, observing these perturbations at the stellar surface is hampered by the optically thick wind and only indirect evidence might be inferred (Langer 1989; Gräfener & Vink 2013).

Another type of spectral variability is present in the wind of some WR stars. This variability is related to the presence of the corotating interaction regions (CIRs; Mullan 1984, 1986; Cranmer & Owocki 1996; Dessart & Chesneau 2002; St-Louis et al. 2009; Chené & St-Louis 2011), connected with the rotation of WR stars of the order of 50 km s-1 (Meynet & Maeder 2003; Yoon et al. 2006; Gräfener et al. 2012b). This large-scale variability in the wind of hot massive stars has a periodic behaviour, although this behaviour is epoch dependent, and is thought to be a signature for spiral-like structures appearing in the wind and carried around by rotation (Cranmer & Owocki 1996; Dessart & Chesneau 2002; St-Louis et al. 2009). The triggering mechanism at the base of the wind is probably not related to stochastic convective motions given the periodical stream-like nature of CIRs. Instead, studies have suggested that CIRs are analogous to the spiral structures in the solar corona and solar wind (Hundhausen & Gosling 1976; Gosling & Pizzo 1999; St-Louis et al. 2009), which are linked with magnetic fields and most likely also connected with shocks observed in the wind of hot luminous stars (Berghoefer et al. 1997; Marchenko et al. 2006; Lépine & Moffat 2008). In this paper we neglect to include both CIRs and large-scale variability in our direct analysis, but rather try to enlighten some aspects related to the observed small-scale variability in WR stars.

Furthermore, theoretical works (e.g. Glatzel et al. 1993; Glatzel 1994; Saio et al. 1998) predict the appearance of periodic variabilities arising from pulsations in the envelopes of stars with high luminosity-to-mass ratios. However, up to now, this periodical variability has not been clearly observed and identified, challenging both observations and models. In the context of WR stars, both observation and pulsational analysis need to take into account the effects of mass loss via stellar wind in the outer layers, a key feature in this class of objects. Therefore we also investigate the instability to pulsations of our helium star models in a hydrodynamically consistent way, including the effects of mass loss for the first time.

We present a set of helium star models computed with a state-of-the-art hydrodynamical stellar evolution code. We present the methods used to compute stellar models and to study the interaction between convective and radiative layers in Sect. 2, and the results of our analysis are shown in Sect. 3. In Sect. 4 we present the results obtained from the pulsating models, in Sect. 5 we discuss the possible observational signatures comparing our predictions to previously and newly analysed observations, and in Sect. 6 we give our conclusions.

Table 1

Properties of our set of helium zero age main-sequence models with Z = 0.02.

2. Method

We computed a set of chemically homogeneous helium star models with masses ranging from 2 M to 17 M with a one-dimensional (1D) hydrodynamical stellar evolution code described in Heger et al. (2000), Petrovic et al. (2006), Brott et al. (2011), and references therein. The approach to use chemically homogeneous models is justified since the properties of H-free WN stars are insensitive to the internal chemical structure or evolutionary history (Langer 1989). The code adopts the standard non-adiabatic mixing length theory (MLT) to model the convective regions and physical quantities associated with the convective layers such as the convective velocities (Böhm-Vitense 1958; Kippenhahn & Weigert 1990; Heger et al. 2000). The mixing length parameter α is assumed to be equal to 1.5 as in Abbett et al. (1997). The models are computed with outer boundary conditions derived from the assumption of a plane parallel grey atmosphere. A consequence of this assumption is that the feedback from an optically thick wind on the outer hydrostatic layers is neglected. The models are computed with OPAL opacity tables from Iglesias & Rogers (1996) and with a metallicity of Z = 0.02. We adopt the mass-loss rate prescription from Nugis & Lamers (2000) and neglect the effects of rotation and magnetic fields.

In order to determine whether a layer is unstable to convection, we adopt the Schwarzschild criterion, which can be expressed as (1)where β is the ratio of gas pressure to total pressure, M is the stellar mass, L the luminosity of the star, G the gravitational constant, κ the Rosseland mean opacity, and Γ the local Eddington factor (Joss et al. 1973; Langer 1997).

Below, we compare the convective velocities vc with the isothermal sound speed to define the Mach number Mc. The low density inflated envelope is defined as the layers above the point at which β is equal to 0.15 (Sanyal et al. 2015). The comparison with the isothermal sound speed in the envelope can be justified by analysing the timescales in the radiation pressure dominated envelopes. The thermal timescale is (2)where \begin{lxirformule}$ \Delta M_{\rm env}$\end{lxirformule} is the mass of the inflated envelope and R the stellar radius. This timescale is of the order of 100 s in our set of models (see Table 2), whereas the dynamical timescale τdyn, which is written as , is one order of magnitude larger. This suggests that a perturbation would rapidly thermalize, therefore the isothermal sound speed is thought to be a better approximation for the real sound speed in these kinds of models.

Cantiello et al. (2009) presented a study of the convective regions in the outer envelopes of massive main-sequence stars and of the induced velocity field at the surface via propagation of uncorrelated gravity waves. Following Cantiello et al. (2009), we introduce vc as the convective velocity averaged over one mixing length starting from the upper boundary of the sub-surface convective zone at r = Rc, defined as in Cantiello et al. (2009), (3)where the mixing length αHP is obtained by the product of the pressure scale height HP (with the pressure given by the sum of gas and radiation pressure) and α. Gravity waves propagate typically with a frequency ω ~ Mcωac, where ωac is the acoustic cut-off frequency, which is of the order of 10-2 Hz in our models, while acoustic waves are expected to have frequencies ω>ωac. Goldreich & Kumar (1990) showed that it is possible to estimate the fraction of convective energy flux going into acoustic and gravity waves from (4)where Fc is the total convective energy flux and Fa and Fg are the energy fluxes transported by acoustic and gravity waves, respectively. As a result, in the hydrostatic sub-sonic regime (where Mc< 1) we expect most of the convective energy to go into gravity waves.

In order to investigate the influence of the sub-surface convective zone in helium stars and relate the convective motion to observational phenomena at the surface without the ability to compute the energy dissipation by the gravity waves along the path exactly, we introduce an upper limit to the expected velocity amplitude at the surface based on the conservation of energy, namely (5)where \begin{lxirformule}$\rho_{\rm s}$\end{lxirformule} and ρc are the densities at the surface and at the upper border of the convective zone and where the Mach number is defined as (6)where cT is the isothermal sound speed at the top of the FeCZ (Goldreich & Kumar 1990; Cantiello et al. 2009).

3. Sub-surface convection

Our helium star models with masses ranging from 2 to 17 M are within the luminosity range log (L/L) ≈ 4−5.5 and the temperature range log (Teff/ K) ≈ 4.9−5.1 (Table 1), which is very close to the helium zero age main-sequence computed by Petrovic et al. (2006) and by Gräfener et al. (2012a). We encounter numerical difficulties in treating models more massive than 17 M, especially in connection to mass loss by stellar wind and, therefore, we do not investigate them. We find inflated envelopes in the outer layers of models with M ≥ 6 M, which become more extended with increasing luminosity-to-mass ratios. For nearly all inflated regions we find Γ(r) ≈ 1 (Gräfener et al. 2012a; Sanyal et al. 2015, see Fig. 1).

All computed models have a convective region close to the surface extending over a significant fraction of the stellar radius (10%) and that comprises a very small amount of mass, i.e. less than 10-8M (see Table 1). This convective region arises because of the opacity bump at log (T/ K) ≅ 5.3 (see Fig. 2). An estimate of the relative amount of flux carried by convection (Kippenhahn & Weigert 1990) in our models can be given via (7)where Fc and Ftot are the convective and total flux, respectively, cP is the specific heat at constant pressure, T is the temperature, ρ is the density, and g is the gravitational acceleration. The fraction of convective flux decreases when the convective layers are increasingly radiation pressure dominated, i.e. with increasing stellar mass, going from 10-6 in the low mass to 10-10 in the more massive helium star models of our sample. The partial ionization zones due to helium and hydrogen recombination are absent as the surface temperatures of our models are log (Teff/ K) ≥ 4.9 (see Fig. 2).

thumbnail Fig. 1

Density profiles (in g / cm3) in the outer regions of our helium stars models with stellar masses ranging from 8 M to 16 M. The extent of the inflated envelope increases while the density decreases with increasing mass.

Open with DEXTER

thumbnail Fig. 2

Opacity κ as a function of temperature for helium stars of different masses using OPAL opacity tables. The opacity bump from the iron group elements is clearly visible at log (T/ K) ≅ 5.3.

Open with DEXTER

Our models show that the convective velocities within the sub-surface convective zone are a function of the stellar mass and of the radial coordinate within the convective zone. Figure 3 shows the presence of a convective zone in all the investigated models, whose spatial extent gradually increase for higher masses.

The helium star models up to 10 M display a convective zone that extends over less than 1% of the stellar radius and with convective elements that reach on average velocities not higher than few km s-1. The models with M ≥ 10 M have convective zones that show convective velocities of the same order as the local sound speed, i.e. velocities up to 30 km s-1. Figure 3 shows how convection moves deeper inside in normalized radial coordinate as the envelope becomes larger, and, therefore, it occupies a larger fraction of the stellar radius. The convective region also moves closer to the surface, reaching it in the 17 M model.

thumbnail Fig. 3

Normalized radial coordinate as a function of the mass of our helium star models for the outer 6.5% of the stellar radius. The white region represents radiative layers, while the coloured region denotes convection. The colours indicate the convective velocities within the FeCZ (in km s-1); the black dashed line superposed lies one mixing-length away from the upper border of the convective region. The colours in the top bar indicate the expected velocity fluctuations at the surface (see also Fig. 4 for a direct comparison between vc and vsurf).

Open with DEXTER

The increase in size of the convective zone is connected to the larger inflated envelope when moving towards higher luminosity-to-mass ratios (Petrovic et al. 2006; Sanyal et al. 2015, see Fig. 1) The envelopes become more and more dominated by radiation pressure dominated for higher mass stellar models, density, and β decrease while Γ(r) stays approximately equal to unity and the Schwarzschild criterion for convection (see Eq. (1)) is fulfilled in a larger region. The convective velocities also significantly increase for higher mass models, with the maximum occurring in the lowest density regions (see Fig. 1), i.e. where β is at its minimum.

We compute the average velocity of the convective elements in the upper part of the sub-surface convective zone (Eq. (3)), assuming implicitly that only the elements in the last mixing length can interact with the upper radiative layer (Fig. 3). The models show that these average velocities are smaller than 2 km s-1 for stars with masses below 10 M. For models with M ≥ 10 M, vc becomes larger with values of the order of 20 km s-1 in the case of 15 M (Fig. 4). This result can be related again with the fact that when we approach higher stellar masses, the luminosity-to-mass ratios increase. This inflates further the stellar envelope which in turn increases the mixing length and, combined with the higher luminosities, increases the convective velocities (Kippenhahn & Weigert 1990; Sanyal et al. 2015, see Fig. 1). This is in agreement with Cantiello et al. (2009) and Grassitelli et al. (2015b), who find that convective velocities in the FeCZ are higher for higher luminosities, while here, conversely from main-sequence stellar models, the convective velocities are in general higher for higher effective temperatures.

Figure 4 shows the average convective velocities at the top of the iron convective zone and the expected velocity perturbation at the surface. We can see how the velocities are attenuated in the radiative zones, i.e. the expected velocity fields at the surface are smaller than vc. This is due to the factor in Eq. (5). However the expected vsurf roughly follow the increase of the convective velocities in the case of the most massive models. Moreover in the models computed with a stellar mass larger than 12 M, the relative attenuation becomes lower and the difference between the average convective velocities and surface velocity fields is reduced. In fact the convective layers approach the surface and the radiative layers are less extended. Both the Mach number and ratio ρc/ρs approach unity in the models more massive than 15 M with the convective zone reaching the surface in the 17 M model. In this configuration one can forsake the assumption made in Eq. (5) concerning the propagation of the waves, given the absence of a radiative layer separating the convective region from the surface.

We thus infer from Figs. 3 and 4 that low-mass helium stars are not expected to show strong velocity fields originating in the iron convective zone at the base of the wind. A steeper trend of increasing surface velocity fluctuations is expected for helium stars with higher luminosity-to-mass ratios due to the higher convective velocities in the region of the iron bump and the lower distance from the surface. An extrapolation of these results to higher masses (17 M) needs to be performed with caution. This is because the computed convective velocities are expected to become supersonic and the standard MLT may not apply anymore (Goldreich & Kumar 1990; Canuto & Mazzitelli 1991).

thumbnail Fig. 4

Average velocity (black solid line, vc) of the convective elements at the top of the convective zone and expected surface velocity fluctuations (blue dashed line, vsurf) as a function of stellar masses.

Open with DEXTER

The number of clumps expected to be triggered by these perturbations may be roughly estimated, assuming a transversal correlation length at the surface of the order of the local pressure scale height (Cantiello et al. 2009). This is usually of the order of 109 cm in our massive models, a scale that is very similar to the lateral coherence scale of few degrees that is able to reproduce observations by Dessart & Owocki (2002) while investigating the line profile variability of small-scale structures in the wind of WR stars. With this approach, we estimate a number of clumps NClumps, (8)with a decrease in number for the more massive models.

In summary, convective velocities of the order of the local sound speed have been found in the case of the highest mass models considered here. A steep increase in vc has been found starting from 10 M. Similarly, the expected velocity field at the surface sharply increases for M> 10 M, approaching the local sound speed for the most massive models.

4. Pulsations

Table 2

Physical parameters of the models that are unstable to pulsations.

thumbnail Fig. 5

Evolution of the radius as a function of time for a 10 M radially pulsating helium star model. The pulsation amplitude grows rapidly before saturating.

Open with DEXTER

In low density extended envelopes of stars with high luminosity-to-mass ratio, such as the envelopes of our massive helium main-sequence stars (Fig. 1), the appearance of so-called strange mode pulsations is expected (Glatzel et al. 1993; Saio et al. 1998; Owocki 2015). The main characteristic of these pulsations is their occurrence in regions where the ratio of the local thermal-to-dynamical timescale is small (τth/τdyn ≪ 1, see also Sect. 2 and Table 2). This is the case for WR stars and makes them good candidates for the appearence of pulsational instabilities in their radiation pressure dominated envelopes. The short thermal timescales exclude a thermal origin for these kind of pulsations, unlike κ and ϵ mechanisms (Glatzel 1994; Saio et al. 1998; Blaes & Socrates 2003).

We use the same hydrodynamical stellar evolution code introduced in Sect. 2 to analyse the pulsational properties of our models. This code has already been used to investigate pulsations in the red supergiant phase by Yoon & Cantiello (2010) and is tested against linear pulsation analysis for basic pulsational properties in Heger et al. (1997). The code adopted here is fully implicit and therefore numerical damping may be present (Appenzeller 1970).

Radial pulsations in our models were found to be excited in the mass range 914 M. They show a short growth time of the order of few dynamical timescales, which is a characteristic feature of strange mode instabilities (see Fig. 5 and Glatzel 1994). After the growth phase, all these models reach saturation. Stability tests with respect to radial perturbations for helium main-sequence models have been investigated in Glatzel et al. (1993), who identified the set of unstable modes as strange modes. We find that our pulsation periods are in good agreement with those of the lowest order unstable modes in Glatzel et al. (1993). The periods obtained with our models are plotted in Fig. 6 and listed in Table 2. The pulsation periods almost linearly increases with stellar mass (Fig. 6). This is related to the larger spatial extent of the envelope and increased dynamical timescale in the high mass models. No pulsations are found for models below 9 M in contrast to the stability analysis performed by Glatzel et al. (1993), where modes are excited for stellar masses as small as 5 M.

From Fig. 7 we can see how the extent and position of the convective zone within the stellar model evolves during a pulsation cycle. The convective zone in Fig. 7 is no longer strictly shaped by the iron opacity bump, but follows the perturbed layers where the density, and therefore the opacity, is higher. It moves periodically from deeper inside the star where the convective zone forms during the contraction of the envelope, up to the surface when the compression wave reaches the surface and forces the envelope to expand. During a pulsation cycle part of the envelope shows significantly high values of the local Eddington factor, up to Γ ≈ 1.5−2 in the contraction phase.

thumbnail Fig. 6

Period of a saturated pulsational cycle for different pulsating model masses.

Open with DEXTER

thumbnail Fig. 7

Radial coordinate and surface temperature as a function of time for a 10 M pulsating stellar model. The black thick line is the surface of the stellar model and the green striped region shows the layers unstable to convection during a pulsation cycle. The blue dot-dashed line tracks the location as a function of time of the temperature at which the iron opacity peak is maximum in a static configuration, i.e. log (T/ K) ≈ 5.3 (see also Fig. 2).

Open with DEXTER

thumbnail Fig. 8

Saturated amplitude of the pulsations appearing in the helium star models. The decreasing amplitude above 13 M is believed to be an effect of the increased mass-loss rate.

Open with DEXTER

thumbnail Fig. 9

Maximum radial velocity reached at the surface during a saturated pulsation cycle in the helium star models. The sound speed in these models is of the order of 30 km s-1.

Open with DEXTER

The saturated amplitude of the pulsations ΔR is plotted in Fig. 8. The amplitude is defined as the difference between the minimum and maximum radius in a saturated pulsation cycle. In Fig. 9 we show the maximum surface radial velocity amplitude due to the radial pulsations (). Figure 8 shows an increase of the maximum pulsation amplitude from 9 to 13 M, which is connected to the higher radial velocity (see Fig. 9). In contrast, the 14 M model has a lower maximum amplitude and maximum velocity than lower mass models. No other pulsating model has been found for masses M> 14 M, when mass loss was applied. Conversely, we find that higher mass models without mass loss are also unstable to more than one mode. We therefore interpret the decrease of the amplitude for the models with higher masses as the result of the higher mass-loss rate. For the pulsating 14 M model, which shows ΔR/R of the order of 0.03 with no mass loss applied, the pulsation amplitude is smaller when mass loss is applied (ΔR/R = 0.0085). In these cases the amount of mass depleted from the inflated envelope during a pulsation cycle becomes larger than 5%, which is equivalent to 40% of the envelope radial extent in a hydrostatic non-pulsating configuration.

Although the WR mass-loss rates proposed by Nugis & Lamers (2000) are approximately 2 orders of magnitude smaller than the critical mass-loss rate necessary to remove the inflated envelope (Petrovic et al. 2006), the absence of radial pulsations for masses M> 14 M can be explained as a stabilizing influence of mass loss on the strange modes instability. These models are in fact expected to be unstable in the context of the analysis carried out by Glatzel et al. (1993). However, neither Glatzel et al. (1993) nor any other author (at the best of our knowledge) included mass loss in their stability analysis.

From Fig. 9 we notice that pulsations occur with supersonic velocities with values up to a factor 4 higher than the isothermal sound speed. This is another peculiar feature of the strange mode instability (Glatzel 1994, 2005). Instead, the luminosity variations during a saturated pulsation cycle are very small (see Table 2), corresponding to a variability of the order of 10-7 bolometric magnitudes. The results obtained here show that strange modes in mass-losing, H-free WR stars occur in the mass range 914 M with periods of the order of minutes, amplitude of the order of 1% of the stellar radius, and supersonic radial velocities. The amplitude of the pulsations increases as a function of mass in the range 913 M, while pulsations are damped or inhibited for the more massive models considered here, i.e. M ≳ 13 M and with ≳ −10-5M/ yr. Consequently, the most massive H-free WR stars may not be unstable against pulsations. Therefore, the role of pulsations in driving the wind or in inducing instabilities in WR stars may not be crucial, or their importance can be restricted to the mass range 914 M.

5. Comparison with observations

One of the dominant features characterizing WR stars is their strong, partly optically thick wind. Such a feature shrouds the hydrostatic surface of the star and we are only able to have indirect evidence of the physical conditions at the base of the wind.

5.1. Spectral variability

Seeking to test the influence of the velocity field generated by the sub-surface convective zone and to relate it to the formation of clumps, we compare our results with the recent studies by Michaux et al. (2014), Chené & St-Louis (2011), and St-Louis et al. (2009). These works focus on WR stars and the variability of WR spectral lines and are based on atmosphere models from Sander et al. (2012) and Hamann et al. (2006). Additionally, we analyse the spectra of eight WR stars of the WN spectral type as described in Appendix A.

Our goal is to verify whether an increasing small-scale variability level as a function of the estimated mass can be found in some spectral lines, as expected from Fig. 4. Following St-Louis et al. (2009), Chené & St-Louis (2011), Michaux et al. (2014), and the formalism of Fullerton et al. (1996), we consider the rms variation across the Heii wind spectral lines relative to the local continuum strength identified as variable σ. We did this to investigate the spectroscopic variability in the sub-sample of H-free WR stars of the WN subtype (see Chené & St-Louis 2011, and reference therein for further details). The degrees of spectral variability σ were derived in general from only four to five spectra per star, which does not allow a direct distinction between stochastic and periodic variability. This in turn does not always allow for a direct distinction between WR stars showing CIRs and WR stars with only clumping, which is achieved by characterizing the kinematics of the excess emission on top of the emission lines that either move from one line edge to the other in the case of CIRs, or from the line centre to the line edges in the case of clumps.

From Michaux et al. (2014), which is mostly a compilation of the data of St-Louis et al. (2009) and Chené & St-Louis (2011) for the WN class, we exclude the objects showing to some extent hydrogen lines. In fact, H-rich WN stars can be structurally different from the H-free models that have been computed in this work. Additionally, we exclude stars that have been reported to most likely have CIR-type variability (St-Louis et al. 2009; Chené & St-Louis 2011; St-Louis 2013). The CIR-type variable spectra are associated with large-scale variability, which most likely does not find its origin in the small-scale inhomogeneity investigated here. Therefore these objects are not included in the linear fit in Fig. 10 (empty circles) and in the following analysis. In other words, only the sample of single1 H-free WN stars has been plotted in Fig. 10 as a function of the estimated mass.

thumbnail Fig. 10

Variability of rms relative to the line strength σ as a function of mass of the Galactic H-free WN stars. The blue dots correspond to the small-scale variability in WN stars from St-Louis et al. (2009), Chené & St-Louis (2011), Michaux et al. (2014), while the pink dots correspond to the variability for the stars analysed in this work (see Appendix A). Empty black circles correspond to CIR-type variability (St-Louis et al. 2009; Chené & St-Louis 2011; St-Louis 2013). The stellar masses are derived by Hamann et al. (2006). The linear fit in blue dashed is computed considering only the blue and pink dots. The Spearman’s rank correlation coefficient is 0.7.

Open with DEXTER

Most of the spectral variability studies on which the σ values are extracted are presented in the works of St-Louis et al. (2009) and Chené & St-Louis (2011). An additional eight WN stars, however, partially presented in Chené et al. (2012) were also included (pink dots in Fig. 10). They are the Galactic stars WR 7 (WN4b), WR 20 (WN5O), WR 34 (WN5o), WR 37 (WN4b), WR 51 (WN4o), WR 62 (WN6b), WR 84 (WN7o), and WR 91 (WN7b). Each of these stars was observed four to eight times with the Gemini Multi-Object Spectrograph (GMOS) at Gemini South, under programme numbers GS-2008B-Q-87, GS-2010B-Q-58, GS-2014A-Q-42, and GS-2014A-Q-73. In Appendix A, we present the spectra and a detailed variability analysis comparable to what is available in the literature for other WN stars.

The observed WR stars only partially fall in the mass range presented in Fig. 4. However, the linear fit in Fig. 10 shows larger variability at higher stellar masses. This trend is in agreement with the results shown in Fig. 4. In addition the zero point of the linear fit at 10 M matches our prediction. Moreover, we do not find any correlation when the degrees of spectral variability are compared to the terminal wind velocities tabulated by Hamann et al. (2006), suggesting that the spectral variability is not related to the specific wind structure but rather by the conditions in the deeper hydrostatic structure of the WR stars.

While these results have to be taken with caution given the small sample and the associated uncertainties, we conclude that our results are supporting the idea that sub-surface convection may trigger the formation of clumps in the wind of Wolf-Rayet stars.

5.2. Characteristic timescale and number of clumps

A further comparison with observational data could involve the characteristic timescale associated with convection, τconv, defined as (9)This timescale is of the order of 102−103 s in our models with M ≥ 10 M (Table 1). Clumping could be detected by means of linear polarimetry (see Davies et al. 2005, 2007; Vink 2012; St-Louis 2013) and a Fourier analysis of its variability could lead to the direct detection of the sub-surface convection in WR stars where the optical depth of the wind is sufficiently low.

The derived timescale appears to be in good agreement with spectroscopic variations associated with small-scale structures in the wind of WR stars. Lépine et al. (1999) obtained a time series of high-resolution and high signal-to-noise ratio spectra of the WC8 star γ2 Velorum with a relatively high sampling rate of one spectrum every 300 s. From these observations, they concluded that the temporal variability of the sub-peaks associated with clumps in the wind took place on a timescale of the order of minutes to hours2.

The timescale might also be observable by detecting short time variations in different parameters of linear polarimetry such as the total linear polarization or the angle of maximum polarization. The detection of the stochastic variability on such a short timescale requires the use of extremely short exposures and readout times, therefore limiting both the number of usable instruments and observable targets. In this respect it might be best to use instruments capable of obtaining both linear polarimetry Stokes parameters Q and U simultaneously to halve the observing time.

Furthermore, we can compare the estimated order of magnitude for the number of clumps as from Eq. (8), i.e. about 103−104 clumps, with observations. The wind crossing time, usually of the order of hours (Lepine et al. 1996; Chené et al. 2008), is comparable to the excitation timescale above, therefore even if our estimate shall be considered as a lower limit, we do not expect this to differ significantly compared to the integrated number of clumps throughout the stellar outflow. We find that indeed our estimate is in agreement with the number of clumps observed by Lépine & Moffat (1999) for WR stars (104) as well as with Davies et al. (2007) for O and LBV stars.

5.3. Magnetic fields

The formation of magnetic flux tubes could also be connected to sub-surface convection in hot rotating stars (Prinja et al. 2002; MacGregor & Cassinelli 2003; Cantiello et al. 2009, 2011; Cantiello & Braithwaite 2011). A magnetic field can be generated from a dynamo mechanism taking place in the FeCZ and propagate through the outer layers of WR stars, affecting their winds and eventually the structures in these winds. The presence of extended FeCZs in the most massive WR models could be associated with different phenomena as was successfully accomplished in the case of O-stars (e.g. Prinja et al. 1995; Prinja & Crowther 1998; Cantiello & Braithwaite 2011). These phenomena are discrete absorption components (e.g. Cranmer & Owocki 1996), i.e. X-ray emission (e.g. Kramer et al. 2003; Oskinova et al. 2009) and large-scale variability induced by corotating interaction regions (e.g. Cranmer & Owocki 1996; Chené & St-Louis 2011; Chené et al. 2011b), which appear in several WR star winds.

Based on the model for the rise of toroidal magnetic flux tubes generated by a dynamo mechanism in the convective zone proposed by MacGregor & Cassinelli (2003), who investigated the buoyant transport of magnetic flux through the interior of hot stars, we estimate the field strength of the magnetic flux tubes in the sub-surface convective zone. Assuming equipartition between the magnetic and kinetic turbulent energy as in Cantiello et al. (2009), we can write the magnetic field strength B as (10)From Eq. (10) we derive a magnetic field at the interface between the convective and overlying radiative zone of the order of tens of gauss for the low-mass (10 M) helium stars, while for the more massive models we find magnetic fields of up to 500 G (cf. Table 1). If we assume conservation of the magnetic flux for the magnetic flux tubes (MacGregor & Cassinelli 2003), we find that the strength of the magnetic field arising at the surface is roughly the same as that at the interface between the convective and radiative layers, i.e. of the order of hundreds of gauss in the models with M ≳ 10 M. However the magnetic pressure at the surface is still only a fraction of the surface gas pressure, i.e. at maximum a third of it in the most massive model. This implies that the gas and radiation pressure most likely dominate and define the global structure of the outflow at the base of the wind (Gary 2001), but local twisting of the magnetic field lines might lead to locally enhanced magnetic field strength such that the magnetic pressure may exceed the local gas pressure with potential observational effects, such as CIRs.

Despite the simplified method and assumptions, and the observational difficulties in detecting and estimating the magnetic fields of WRs, we find that the predicted values are of the same order of magnitude as the most probable values determined for a sub-sample of WR stars by de la Chevrotière et al. (2013, 2014). These are marginal detections of the field in the observable region of the wind and that the field configuration was assumed to have a split monopole configuration.

5.4. Pulsations

It has been suggested that violent mass-loss events and large-scale variabilities may be connected with strange mode instabilities in hot massive stars (Glatzel et al. 1993; Glatzel & Kaltschmidt 2002). A clear observational signature for fast pulsating WR stars has not yet been found. After the non-confirmation of the presence of fast pulsations in WR40 of spectral sub-type WN8 (Blecha et al. 1992; Martinez et al. 1994), the only candidate at the moment is the H-depleted WR123 of the spectral sub-type WN8, which shows a 9.8 h period. This period however does not match the period range of the lowest order unstable modes in Fig. 6 (Balona 2010; Lefèvre et al. 2005; Chené et al. 2011a). However, WR123 is expected to be considerably more massive then the mass range considered here. Previous attempts to reproduce this 9.8 h period have been unsuccessful (see Chené et al. 2011a, and reference therein).

The predicted bolometric magnitude change of our models during a pulsation cycle is of the order of Δm ∝ log (Lmin/Lmax) ≈ 10-7 (see Table 2). However, the associated changes in effective temperature imply a variability as high as 0.1 mag in the X-ray and visible bands. While this may be above the level of non-periodical variability of WRs at the mmag level (Gosset & Rauw 2009; Balona 2010; Chené & Moffat 2011; David-Uraz et al. 2012), a periodic variability in the frequency range we predict has not been detected in the bright Galactic WC5 star WR111 by Moffat et al. (2008, see also Lefèvre et al. 2005; Chené et al. 2011b). The pulsations may be generally hidden by other types of variability appearing in WR outflows and by the optically thick winds, given the high optical depth of the hydrostatic layers. Whether any variability induced by the pulsations is detectable depends on how pulsations affect the stellar outflow and whether coherence is preserved throughout the wind. From Fig. 8 and the discussion in Sect. 4, it is suggested that the wind has a stabilizing effect on strange modes instability and inhibits them in most cases. This might explain why, even if expected, there has been no match between the theoretical predictions and observations of strange modes pulsations in H-free WR stars, as previous models neglected the effect of mass loss to the stability of strange modes.

A different explanation for small- or large-scale variability might instead involve the effects of turbulent pressure fluctuations and finite amplitude high order non-radial pulsations of non-thermal origin. These are expected to appear in the case of transonic convective velocities and inefficient convection, as has been suggested to connect the strength of turbulent pressure in the partial ionization zones to the appearence of the macroturbulent broadening in OB and late-type stars; this is a broadening of spectral lines on a scale larger than the line forming region (Grassitelli et al. 2015a,b). In our models the turbulent pressure contributes less then 1% to the equation of state, but the Mach numbers are close to unity in the more massive models.

6. Conclusions

The purpose of this work was to study the condition in the envelopes of massive helium stars and analyse the instabilities appearing in these radiation pressure dominated regions. These models correspond to H-free Wolf-Rayet stars, which are located close to the Eddington limit. These instabilities can possibly provoke or influence the formation of structures observed at the surface or in the wind. We find that sub-surface convective zones are present in all the investigated models. Their spatial extent increases for higher stellar masses. The convective zones arise as a result of the iron opacity peak and they are inefficient in transporting energy. Convection reaches the surface in the 17 M model. In our models the iron opacity bump is also responsible for the envelope inflation (Sanyal et al. 2015) starting from the 6 M model. These envelopes are characterized by their proximity to the local Eddingont limit (i.e. Γ ≈ 1).

The convective motion may trigger velocity and density fluctuations at the stellar surface. An upper border for the surface turbulent velocity is estimated by considering the energy transported via gravity waves from the upper limit of the sub-surface convective zone to the surface. We find that the expected surface velocity amplitudes are small (1 km s-1) for stellar models with masses below M ≲ 10 M, but are more than 20 km s-1 for the more massive stellar models (M ≳ 10 M). Therefore if clumping is triggered at the base of the wind by these fluctuations, we expect that the formation of structures should be inhibited or spatially delayed (Gayley & Owocki 1995; Owocki 2015) in the case of low-mass WNE stars. However, our models predict strong small-scale variability for higher WR masses. A trend of stronger variability for higher stellar masses was identified from the observational data of St-Louis et al. (2009), Chené & St-Louis (2011), and from this work (see Sect. 5.1 and Appendix A) for single H-free WNE stars. We also conclude that if wind clumping is triggered by sub-surface convection, the assumption of a constant clumping factor as a function of mass may not be appropriate. We find also that the timescales associated with convection, of the order of hours, are supported by observations by Lépine et al. (1999). Moreover, the number of clumps that we estimate by considering the size of the biggest convective eddies and the correlation lengths associated with the gravity waves, i.e. 103−104, is also supported by observations by Lépine & Moffat (1999).

We find that the envelopes of the helium star models are unstable to pulsations in the range 9−14 M. The periods of the pulsations are in agreement with a low order mode in Glatzel et al. (1993). The variations in bolometric magnitude connected with these pulsations are very small (less than micromagnitudes), i.e. below the observable photometric variability and the observed variabilities connected with structures in the wind. However the supersonic radial velocities could give rise to observable periodical effects apparently not yet observed in the context of WR stars. Moreover, the strong mass loss applied to our stellar models has a stabilizing effect on the pulsations. It appears that mass loss has an inhibiting effect on the pulsational instability in the most massive cases. This may help to explain the lack of observational evidence for periodical variability in WR stars of the WNE type, which is in contrast to previous direct theoretical predictions that did not include mass loss.

We also find that the strength of the magnetic fields rising from the convective zone to the surface of our most massive models, of the order of hundreds of gauss, is comparable to the most probable values determined in some WR stars by de la Chevrotière et al. (2013, 2014). The pressure induced by the buoyant magnetic flux tubes is, however, globally not enough to define the structure of the wind, but local enhancement of the magnetic field might lead to or trigger the formation of CIRs.


1

Some WR stars might have an O-B companion, which in any case does not significantly affect the WR spectrum and derived physical quantities, according to Hamann et al. (2006).

2

The variability timescale does not have to be confused with the time necessary for a sub-peak to migrate from the line centre to the line edges.

3

iraf is distributed by the National Optical Astronomy Observatories (NOAO), which is operated by the Association of Universities for Research in Astronomy, Inc. (AURA) under cooperative agreement with the National Science Foundation (NSF).

Acknowledgments

L.G. is part of the International Max Planck Research School (IMPRS), Max-Planck-Institut für Radioastronomie and University of Bonn and Cologne. L.F. acknowledges financial support from the Alexander von Humboldt foundation. Based on observations obtained at the Gemini Observatory, processed using the Gemini IRAF package, which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), the National Research Council (Canada), CONICYT (Chile), the Australian Research Council (Australia), Ministério da Ciência, Tecnologia e Inovação (Brazil) and Ministerio de Ciencia, Tecnología e Innovación Productiva (Argentina). Further, L.G. thanks T. Moffat, Y. Michaux, G. Gräfener, J.-C. Passy, and the referee S. P. Owocki for precious discussions and comments on this manuscript.

References

Appendix A: Variability study of eight additional WN stars

Appendix A.1: Data extraction and analysis

The bias subtraction, flat fielding, spectrum extraction, sky subtraction, and wavelength calibration of all spectra were executed in the usual way using the gemini packages of the iraf3 software. Calibration lamp spectra were taken during the night after each spectrum. Special care was taken for the normalization of the spectra. First, a mean was made for each run. Then, each spectrum of a given star was divided by the mean spectrum and the ratio fitted with a low-order Legendre polynomial (between fourth and eighth order). The original individual spectrum was divided by this fit, and was therefore at the same level as the mean spectrum. This allowed us to put all individual spectra at the same level. The mean spectrum was then fitted in selected pseudo-continuum regions, i.e. wavelength regions where large emission lines do not dominate. Finally, the fitted continuum function is applied to each individual spectrum. The error on the continuum normalization measured as the standard deviation of individual spectra around the continuum function is typically of 0.5.

Appendix A.2: Determining the σ-value

The final spectra are shown in the top panels of Figs. A.1, and A.2 for the three strongest lines present in our wavelength interval, i.e. the Heiiλ4686, Heiiλ5411, and Heiλ5876.

To determine the level of line-profile variability, we first identify at which wavelengths the spectrum may be significantly variable, using the temporal variance spectrum (TVS) of each data set. The TVS was calculated using the formalism of Fullerton et al. (1996) and the quantity was obtained, where σ0 is the reciprocal of the rms of the noise level in the continuum in a time series of N spectra. The value of Σj(99%) quantifies the level of variabilityat each wavelength: a spectrum that reaches a value of n varies with an amplitude n times higher than the variability measured in the continuum (which is assumed to be pure noise) with a confidence level of 99%. The spectrum of a given star is considered significantly variable at a given wavelength j if the value of Σj(99%) is significantly greater than 1. When a line is identified as significantly variable, it is possible to calculate the amplitude of its variability relative to its intensity. To accomplish this, we calculated for each wavelength j a modified as defined in Chené (2007) and divided it by the line flux (), where is the weighted mean flux at wavelength j. This ratio is named the σ spectrum. One should note that the calculation of σ does not take into account instrumental variations due to the noise level. Thus, when the variation level of a given line is too close to the noise, which is the case for weaker lines, σ is artificially high. That is why we manually set σj = 0 when the variability at a given wavelength j is not clearly significant according to the value of Σj(99%) or when the line has a relative intensity that is lower than 1.5, which is the intensity limit of a spectral line from which we can investigate the variability for the present data set.

Figures A.1 and A.2 show our results for the eight additional WN stars. The bottom panels show a superposition of all the observed spectra with the average spectrum overplotted in red. The line intensity is normalized to ease the comparison with the different lines. The normalizing factor of the lines is written on the top left. The second panels from the bottom show the Σ(99%) spectrum. The dotted line marks the threshold of significant variability. The third panels from the bottom show the σ spectrum. The most reliable values are extracted from the central region of the line, as the edges tend to be more sensitive to instrumental effects. Finally, the top panels show the montage of the residuals (individual spectra average). The residuals of the Heii λ 4686 line are divided by a factor 2 for a better visibility. The residuals are in chronological order, with the bottom to the top. The heliocentric Julian date corrected to start on 1 Jan. 2000 is written on the far left of the plot. All the spectra and associated values (Σ(99%)) are masked at the interstellar medium Na I D lines at 5885 and 5890 Å.

thumbnail Fig. A.1

Top: montage of the Heiiλ4686 (left), Heiiλ5411 (centre), and Heiλ5876 (right) residuals (individual spectra mean) for WR7, WR20, WR34, and WR37 (as indicated in the top left section of the top panel). For all cases, the scale factor of the ordinate is indicated in the top right-hand corner of the left column. HJD HJD(2000) is indicated in the y-axis. Second from top: σ spectrum. Second from bottom: Σ (99%) spectrum. Bottom: mean spectrum.

Open with DEXTER

thumbnail Fig. A.2

Same as Fig. A.1, but for WR51, WR62, WR84, and WR91.

Open with DEXTER

All Tables

Table 1

Properties of our set of helium zero age main-sequence models with Z = 0.02.

Table 2

Physical parameters of the models that are unstable to pulsations.

All Figures

thumbnail Fig. 1

Density profiles (in g / cm3) in the outer regions of our helium stars models with stellar masses ranging from 8 M to 16 M. The extent of the inflated envelope increases while the density decreases with increasing mass.

Open with DEXTER
In the text
thumbnail Fig. 2

Opacity κ as a function of temperature for helium stars of different masses using OPAL opacity tables. The opacity bump from the iron group elements is clearly visible at log (T/ K) ≅ 5.3.

Open with DEXTER
In the text
thumbnail Fig. 3

Normalized radial coordinate as a function of the mass of our helium star models for the outer 6.5% of the stellar radius. The white region represents radiative layers, while the coloured region denotes convection. The colours indicate the convective velocities within the FeCZ (in km s-1); the black dashed line superposed lies one mixing-length away from the upper border of the convective region. The colours in the top bar indicate the expected velocity fluctuations at the surface (see also Fig. 4 for a direct comparison between vc and vsurf).

Open with DEXTER
In the text
thumbnail Fig. 4

Average velocity (black solid line, vc) of the convective elements at the top of the convective zone and expected surface velocity fluctuations (blue dashed line, vsurf) as a function of stellar masses.

Open with DEXTER
In the text
thumbnail Fig. 5

Evolution of the radius as a function of time for a 10 M radially pulsating helium star model. The pulsation amplitude grows rapidly before saturating.

Open with DEXTER
In the text
thumbnail Fig. 6

Period of a saturated pulsational cycle for different pulsating model masses.

Open with DEXTER
In the text
thumbnail Fig. 7

Radial coordinate and surface temperature as a function of time for a 10 M pulsating stellar model. The black thick line is the surface of the stellar model and the green striped region shows the layers unstable to convection during a pulsation cycle. The blue dot-dashed line tracks the location as a function of time of the temperature at which the iron opacity peak is maximum in a static configuration, i.e. log (T/ K) ≈ 5.3 (see also Fig. 2).

Open with DEXTER
In the text
thumbnail Fig. 8

Saturated amplitude of the pulsations appearing in the helium star models. The decreasing amplitude above 13 M is believed to be an effect of the increased mass-loss rate.

Open with DEXTER
In the text
thumbnail Fig. 9

Maximum radial velocity reached at the surface during a saturated pulsation cycle in the helium star models. The sound speed in these models is of the order of 30 km s-1.

Open with DEXTER
In the text
thumbnail Fig. 10

Variability of rms relative to the line strength σ as a function of mass of the Galactic H-free WN stars. The blue dots correspond to the small-scale variability in WN stars from St-Louis et al. (2009), Chené & St-Louis (2011), Michaux et al. (2014), while the pink dots correspond to the variability for the stars analysed in this work (see Appendix A). Empty black circles correspond to CIR-type variability (St-Louis et al. 2009; Chené & St-Louis 2011; St-Louis 2013). The stellar masses are derived by Hamann et al. (2006). The linear fit in blue dashed is computed considering only the blue and pink dots. The Spearman’s rank correlation coefficient is 0.7.

Open with DEXTER
In the text
thumbnail Fig. A.1

Top: montage of the Heiiλ4686 (left), Heiiλ5411 (centre), and Heiλ5876 (right) residuals (individual spectra mean) for WR7, WR20, WR34, and WR37 (as indicated in the top left section of the top panel). For all cases, the scale factor of the ordinate is indicated in the top right-hand corner of the left column. HJD HJD(2000) is indicated in the y-axis. Second from top: σ spectrum. Second from bottom: Σ (99%) spectrum. Bottom: mean spectrum.

Open with DEXTER
In the text
thumbnail Fig. A.2

Same as Fig. A.1, but for WR51, WR62, WR84, and WR91.

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.