Free Access
Issue
A&A
Volume 559, November 2013
Article Number A130
Number of page(s) 17
Section Stellar structure and evolution
DOI https://doi.org/10.1051/0004-6361/201322390
Published online 28 November 2013

© ESO, 2013

1. Introduction

The most important property of massive, hot stars is their mass loss expelled via stellar winds. These winds can be extremely strong, which has a significant effect on their evolution and affects their surface abundances (for a review see, e.g., Meynet & Maeder 2007, and references therein).

The line-driven wind theory, first proposed by Lucy & Solomon (1970) and later elaborated on by Castor et al. (1975), can explain the physical mechanism by which massive stars lose mass. However, the accurate values of the wind parameters are still under debate. One of the most challenging problems is to determine reliable mass-loss rates, which are derived from observations with the aid of some physical models. However, discordances between different mass-loss rate diagnostics were found (for a review see Puls et al. 2008).

A major complication of mass-loss estimates comes from the fact that stellar winds are inhomogeneous, as indicated by several observational evidences (see, e.g., Eversberg et al. 1998; Lépine & Moffat 2008; Prinja & Massa 2010), but also predicted by numerical simulations (e.g., Feldmeier et al. 1997; Runacres & Owocki 2002; Dessart & Owocki 2005). These simulations show that the instability of wind line driving leads to the formation of shocks and spatial structures in both density and velocity, that is to say, instability forms “clumps”.

Table 1

Input stellar and wind parameters and element abundances by mass fraction.

Analyses of massive star spectra rely on their comparison with sophisticated model simulations. State-of-the-art model-atmosphere codes account for non-LTE radiative transfer in a spherically symmetric wind and incorporate detailed model atoms along with an approximate treatment of the line-blanketing effect from iron-group elements. A few such codes have been developed, e.g., CMFGEN (Hillier & Miller 1998), PoWR (Hamann & Gräfener 2004), and FASTWIND (Puls et al. 2005). In all these codes, clumping is only included in the approximation that the clumps are assumed to be optically thin at all frequencies (the so-called “microclumping” approximation). The clumps have a density that is enhanced by a “clumping factor” D compared to a smooth wind with the same mass-loss rate. The clumps move according to the adopted velocity law of the wind. In most simulations, the interclump space is assumed to be void.

The main effect of microclumping is that empirical mass-loss rates that are derived from recombination lines, which is a process that depends quadratically on density, are overestimated by a factor of when microclumping is neglected, and have to be reduced accordingly by a mild factor of about 2 or 3. However, Massa et al. (2003) and Fullerton et al. (2006) have studied O-star spectra in the far-ultraviolet (FUV), obtained with the Far Ultraviolet Spectroscopic Explorer (FUSE) satellite, which show the unsaturated resonance line doublet of P v at λλ 1118,1128 Å. The formation of resonance lines depends only linearly on density, and is therefore not sensitive to microclumping. They find mass-loss rates that are lower by factors of 10 to 100, compared with those obtained from recombination lines under the assumption of no clumping. Consequently, these authors conclude that the mass-loss rates have been greatly overestimated when derived from recombination lines, obviously because the clumping contrast D is in fact extremely high (and the volume-filling factor of the clumps correspondingly tiny). Such low mass-loss rates would have far-reaching consequences, say, for the evolution of massive stars.

The same strategy for resolving the mass-loss rate discrepancy was employed in papers by Bouret et al. (2003, 2005). They also reduced the mass-loss rates in order to weaken the UV resonance lines, while the Hα emission is kept at the observed strength by adopting extremely large clumping factors up to D = 100. As an additional means to achieve a consistent fit, Bouret et al. (2012) reduced the phosphorus abundances to values that are lower by a factor 1.4 to 5.1 than the solar abundance according to Asplund (2005). However, subsolar abundances for O stars are not expected and have no justification.

Alternatively, Oskinova et al. (2007) suggest that the mass-loss rate discrepancy can be explained by the effect of optically thick clumps, which was hitherto been neglected within the microclumping approximation. They promoted the “macroclumping” (porosity) approach, taking into account that clumps may become optically thick at certain frequencies. These authors show that, while the optically thin Hα line is not affected by wind porosity, the P v resonance doublet becomes significantly weaker when macroclumping is included. This leads to the conclusion that clumping must be included in the modeling to get different diagnostics to agree on the mass-loss rate and that the microclumping approximation is not adequate for modeling optically thick transitions.

To derive more reliable mass-loss rates from observation, and to resolve discrepancy between different diagnostics, a more detailed treatment of wind clumping is needed. Unfortunately, a general treatment of 3D clumps in full non-LTE radiative transfer simulations is not possible. However, the formation of resonance lines can be treated in the much simpler pure-scattering approximation. This allows the use of Monte-Carlo techniques, which can be adapted to complicated geometrical situations. Sundqvist et al. (2010, 2011) used a 2D and pseudo-3D stochastic wind model, and achieve a reasonable line fit for HD 210839, which is also in our sample (see below).

Šurlan et al. (2012a,b) have developed a full 3D description of clumping and investigated how the properties of clumping (e.g., the velocity dispersion inside the clumps, the radius where clumping sets on, and the density of the interclump medium) affect the resonance line profiles and, consequently, the derived mass-loss rates.

Table 2

Optical, FUV, and NUV observation logs.

The intention of the present paper is to check for a small sample of stars whether our detailed treatment of clumping, together with solar phosphorus abundance and moderate D, may resolve the discordance between the mass-loss rates derived from Hα and P v diagnostics, and also to establish some global properties of wind clumping. We selected 5 O-type supergiants and analyzed their spectra first by means of the Potsdam Wolf-Rayet (PoWR) model atmosphere code, and then applied our 3D Monte-Carlo radiative transfer code for simulating the UV resonance lines.

In the following section we present our stellar sample, observations, and data reduction. The PoWR models and the 3D Monte-Carlo code for the resonance line formation in a clumped wind are introduced in Sects. 3 and 4, respectively. Our spectral fitting procedure is explained in Sect. 5. The results of the consistent analysis of the optical and UV spectra are presented in Sect. 6 and discussed in Sect. 7. Finally, a summary is given in Sect. 8. The complete spectral fits are available in Appendix A.

2. Stellar sample and observation

2.1. Stellar sample

We selected five O-type Galactic supergiants covering spectral types O4If to O7If (see Table 1). All these stars are very luminous and show evidence of a strong wind. Owing to the intrinsic instability of the line driving mechanism it is expected that these winds exhibit pronounced clumping.

These stars have been already analyzed in the optical, UV, infrared, and radio spectral regions (e.g., Markova et al. 2004; Repolust et al. 2005; Puls et al. 2006; Bouret et al. 2012). Stellar parameters of the sample were reliably determined, and mass-loss rates from different diagnostics are also available for comparison. All stars from our sample are presumably single, showing no indications of binarity (Bouret et al. 2012; Mason et al. 1998; De Becker et al. 2009; Turner et al. 2008). The FUV spectra that include the P v resonance doublet are available for all stars of the sample, which is prerequisite to studing the so-called “P v problem”.

2.2. Optical spectra

For the four northern stars in our sample, Hα and blue spectra were observed with a CCD SITe ST-005 800 × 2000 pix camera attached to the Coudé spectrograph of the 2-m telescope at the Ondřejov Observatory (Czech Republic), with the slit width set to 0.6″. Two grating angles were chosen, one centered on the Hα line and the second one covering the He ii 4686 Å and Hβ lines. The achieved spectral resolution is 13 600 and 19 400, respectively.

All spectra were wavelength-calibrated with a ThAr comparison arc spectra obtained shortly after each exposure. The telluric features in the spectra were removed using spectra of the fast-rotating stars 27 Vul and 116 Tau. The data reduction (including telluric and heliocentric velocity corrections) was performed with standard IRAF1 tasks. The program Cosmic Ray Removal2 (dcr, Pych 2004) was used to clean the spectra.

The optical spectrum of HD 66811 was taken at the Complejo Astronómico El Leoncito (CASLEO) in Argentina. The observation was carried out with the 2.15-m Jorge Sahade telescope using a REOSC echelle spectrograph in cross-dispersion mode with a Tek 1024 × 1024 pixel CCD as detector. The adopted instrumental configuration was a 316 l/mm grating at a tilt angle of 5  °40  ′ and a 350 μm slit width, resulting in a resolving power of 12 500. A ThAr lamp was used as a comparison source, with a reference exposure taken immediately after the stellar target at the same telescope position. The data reduction was performed with standard IRAF tasks. The spectroscopic observation logs are documented in Table 2.

2.3. Ultraviolet spectra

The spectral region of the P v resonance doublet was covered by high-resolution observations with FUSE, which we retrieved from the MAST3 archive (see Table 2). To increase the signal-to-noise ratio, the P v spectra of HD 15570 and HD 14947 were smoothed using the splot task in IRAF. Low-resolution near-ultraviolet (NUV) spectra (1200 to 2000 Å), taken with the International Ultraviolet Explorer (IUE), were downloaded from the INES Archive Data Server4 (see Table 2). The FUV spectrum of HD 66811 had been observed with the Copernicus satellite.

All observed spectra were corrected for the radial velocity of the individual star before the comparison with model simulations.

3. 1D spherically symmetric wind models

To analyze the observed spectra we calculated wind models using the PoWR unified model atmospheres code (see Hamann & Gräfener 2004, and references therein). The PoWR code is able to solve non-LTE radiative transfer in a spherically expanding atmosphere simultaneously with the statistical equilibrium equations, and it accounts for energy conservation. Detailed model atoms of the most relevant elements (H, He, C, N, O, Si, and P) are taken into account in the present paper. Line blanketing is also taken into account, with the iron-group elements treated in the super-level approach. Mass-loss rate and wind velocity are among the free parameters of the models.

3.1. Stellar parameters and chemical composition

Stellar and wind parameters of the stars as obtained by Puls et al. (2006), Oskinova et al. (2007), and Bouret et al. (2012) serve as input parameters for our PoWR models (see Table 1). The chemical abundances of H, He, C, N, and O are adopted from Bouret et al. (2012) (see Table 1). For the mass fractions of Si, P, and Fe-group elements, we take the solar values (6.649 × 10-4, 5.825 × 10-3, and 1.292 × 10-3, respectively) as determined by Asplund et al. (2009).

3.2. Velocity field

thumbnail Fig. 1

Dependence of wind velocity on radius. Comparison of the standard (single) β-law (dotted red line) and double-β law (solid black line).

Open with DEXTER

The adopted velocity field in the model consists of two parts. In the inner part, the hydrostatic equation is integrated with the stratification of temperature and mean particle mass, yielding the density stratification. The velocity in this part of the wind is then defined via the continuity equation. This hydrostatic part of the atmosphere is connected smoothly to the wind, with the so-called double-β law (Hillier & Miller 1999). This law consists of the sum of two beta-law terms with different exponents β1 and β2, each of them contributing a prespecified fraction to the total wind velocity. Compared to the standard “one-beta” law, this allows for a smaller velocity gradient in the lower part of the wind, while the second term, for which we always adopt β2 = 6 and a contribution of 35% to the final velocity, causes some noticeable acceleration even at relatively large distances from the star (Fig. 1). The values for ν and β1, as included in Table 1, were also adopted from the references.

In the PoWR models the lines are broadened by thermal and microturbulent motion with νD = 20 km s-1. In addition, radiation damping and pressure broadening are accounted for in the formal integral.

3.3. Clumping in the 1D wind model

In our PoWR models, the wind inhomogeneities are treated in the microclumping approximation (for more details see Hamann & Koesterke 1998). The matter density in the clumps is enhanced by a factor D = 1/fV, where fV is the fraction of volume filled by clumps. In the present study, we have allowed the clumping factor to depend on radius, thereby starting to deviate from the homogeneous wind (D = 1) at about the sonic point (5 km s-1) and quickly reaching D = 10 at ν(r) = 40 km s-1.

4. 3D clumped wind model

thumbnail Fig. 2

Snapshot from a 3D Monte-Carlo simulation, showing in 2D projection an example of a clump realization (red full circles) with respect to star (yellow dashed circle) and the path of a line-scattered photon (blue dotted line).

Open with DEXTER

To model the P v resonance line profiles, we use our Monte-Carlo code for the radiative transfer in a 3D clumped wind. In the spirit of a core-halo approximation, the lower boundary of the wind is placed just above the photosphere. Here we employ the photospheric line spectrum as obtained from the PoWR model as inner boundary.

In the wind, we create a random distribution of spherical clumps. These clumps move with the wind velocity law, but also have an additional internal velocity gradient (see below). The number density of clumps obeys the continuity equation. The density in the clumps and in the interclump medium is specified from the mass-loss rate and the clumping parameters.

The photons that are now released at the lower boundary travel through the wind, where they can be repeatedly scattered in the considered line doublet (continuum opacities are neglected). The line scattering is assumed to be isotropic in the comoving frame of reference, while the frequencies are completely redistributed over the Doppler-broadened opacity profile. The opacity is computed according to the mass-loss rate, element abundance, and ionization fraction. The last is retrieved from the PoWR model. Traveling photons experience Doppler shifts due to the wind expansion. Once a photon crosses the outer boundary of the wind, it is counted for the emergent profile. The principle of this formalism is illustrated in Fig. 2, while more details of the code are given in Šurlan et al. (2012b).

A set of parameters describes the inhomogeneous wind. The clumping factor D specifies the density inside clumps with respect to the smooth wind density. We used the same value as in the PoWR code for the microclumping. Other clumping properties are the clump separation parameter L0, the density of the interclump medium d (for the case of a two component medium), and the radius rcl where clumping sets in. The velocity range inside each clump is described by the velocity deviation parameter m = νdis(r)/ν(r). For a more detailed description of these parameters we refer to Šurlan et al. (2012b).

5. Model fitting

For each star from the sample, we perform the following procedure:

  • a)

    1D models are calculated with the PoWR code in order toestablish the mass-loss rate from fitting the Hα line;

  • b)

    then, the obtained mass-loss rate, the P v ionization fraction, and the photospheric spectrum are used as input for the 3D Monte Carlo simulations of the clumped wind. The clumping parameters are determined by optimizing the fit of the P v resonance doublet.

These steps are described in more detail in the next two subsections.

5.1. 1D PoWR model fitting

As input to the 1D PoWR model calculations, we used the stellar and wind parameters and element abundances as compiled in Table 1. The mass-loss rates were slightly adjusted to optimize the fit with the optical observations (Hα, Hβ, Hγ, and He ii lines). The finally adopted are listed in Table 3. All spectral fits are documented in Appendix A. The synthetic spectra were flux-convolved to simulate instrumental and rotational broadening, taking νsini from Bouret et al. (2012).

Table 3

Finally adjusted stellar and wind parameters.

The UBVJHK photometry of all stars is taken from the GOS catalog (Maíz-Apellániz et al. 2004), and the color excess EB  −  V is adopted from Bouret et al. (2012). We applied the reddening law from Cardelli et al. (1989) and adjusted the RV parameter to optimize the fit between the spectral energy distributions (SED) of the model and the flux-calibrated observations. Moreover, since we kept the luminosity at the literature value, we adjusted the stellar distances to achieve the SED fits. Our final values for RV and the stellar distance are listed in Table 3. The SED fits are documented in the upper panels of Figs. A.1–A.5.

5.2. 3D Monte-Carlo model fitting

thumbnail Fig. 3

The upper and the middle panels show the synthetic spectra (dotted black lines) of HD 66811 (left panels) and HD 15570 (right panels) obtained from the PoWR models, i.e. without macroclumping. Thin solid blue lines are observed spectra. The dashed green lines in the middle panels are from the same model, but only accounting for the photospheric lines while wind lines from P v and Si iv are suppressed. These photospheric spectra are used as input for the 3D Monte-Carlo calculations with macroclumping (lower panels). The parameters of these models are given in Tables 15.

Open with DEXTER

thumbnail Fig. 4

Same as Fig. 3 but for HD 14947 (left panels) and HD 210839 (right panels).

Open with DEXTER

Once the mass-loss rate is established from the Hα fitting, the stratification of the P v ionization fraction and photospheric spectra is extracted from the final PoWR model and used as input for calculating the 3D Monte-Carlo radiative transfer in the clumped wind as described in Sect. 4. To be consistent with the PoWR models, the velocity of the wind is described by a double-β law with the same parameters. The clumping factor D and the Doppler-broadening velocity is taken consistent to the PoWR models. While some of the clumping parameters are fixed (see Table 4), the interclump medium density factor d and the velocity deviation parameter m are varied in order to find the best fit to the observed P v doublet. Their final values are given in Table 5.

Table 4

Fixed model parameters used in the 3D Monte-Carlo code.

Table 5

Clumping parameters that give the best fit to the observed P v line profiles.

6. Results

We first review the global results of the modeling. As seen in the upper panels of Figs. 35, the Hα line fits reasonably well, although not perfectly (see also Hβ and He ii profiles in Appendix A). Comparing our fits with the fits obtained in other investigations, such as the one by Bouret et al. (2012) for the star HD 210839, our fits are not worse. Both our PoWR models and the models of Bouret et al. assumed microclumping throughout the wind. The overall shape of the Hα, Hβ, and He ii 4686 Å lines is fitted well.

For an additional check of the fit, we compared the synthetic UV spectrum with low-resolution IUE spectra (third panels in Figs. A.1A.5). All lines fit well except for the N v resonance doublet. This line is partly formed in the tenuous interclump medium (Zsargó et al. 2008), which is not included in the PoWR model calculations.

The P v resonance doublet, however, is predicted as being far too strong by the PoWR models in all cases (see Figs. 35, middle panels). In contrast, excellent agreement with observations is achieved with the 3D Monte-Carlo simulations (see Figs. 35, lower panels). The disagreement in the red component of the P vλλ 1118,1128 Å line for the star HD 192639 (Fig. 5) is caused by blending with the Si iv 1128 Å line, which is not included in our 3D Monte-Carlo simulations.

In the following sections we describe how we chose the clumping parameters for achieving the best agreement between calculated and observed P v profiles. We demonstrate the effect of these parameters by taking HD 14947 as an example.

thumbnail Fig. 5

Same as Fig. 3 but for HD 192639.

Open with DEXTER

6.1. Number of clumps

The clump separation parameter L0 controls the number of clumps, Ncl, in the wind. Decreasing this parameter causes more clumps. For very small L0 (L0 → 0), the smooth wind is approached (upper panel of Fig. 6). It can be seen that with the smooth wind approximation the absorption is deeper and the emission peaks are higher. For a better fit with observation, we adopt the clumped-wind model and set L0 = 0.5. This implies 1.13 × 104 clumps within 100 R (cf. Eq. (24) in Šurlan et al. 2012b). The calculated line profile (upper panel of Fig. 6) is drastically reduced.

Now we increase the number of clumps, setting L0 = 0.2, which implies 1.75 × 105 clumps within 100 R. Again, neither the strength of the emissions nor the depth of absorptions can be reproduced (upper panel of Fig. 6). Even when we create as many as 1.4 × 106 clumps in the wind by setting L0 = 0.1 (upper panel of Fig. 6), the observed P v line profile is not reproduced.

One may compare these numbers with independent estimates for the number of clumps. For example, Nazé et al. (2013) have recently found that more than 105 clumps are required in the wind of HD 66811 to explain its very low level of stochastic X-ray variability.

6.2. Interclump medium density

For a satisfactory fit of the observed P v profile, additional matter must be located between the clumps. The interclump medium density parameter d (see Sect. 2.1.2. in Šurlan et al. 2012b) defines its density. We find that a reasonable fit to the observation can be achieved this way, even with fewer of clumps. We set L0 = 0.5 and then increase d until satisfactory agreement is reached with about d = 0.2 (upper panel of Fig. 6).

If we decrease L0 (i.e., more clumps), this can be compensated for lower values of d to reproduce the observations. This means that different combinations of L0 and d may give equally good agreement with observations. From our clumped wind model, it is not possible to tell with certainty which combination of L0 and d corresponds to reality. We can only say that for winds that consist of less than about 106 clumps, interclump medium density is a necessary ingredient of the wind in order to satisfactorily reproduce the P v resonance doublet. But in any case, the interclump space cannot be void.

6.3. Onset of clumping

The parameter rcl controls the radius where clumping sets in. Since Sundqvist & Owocki (2013) have shown that structures in the wind may already develop very close to the wind base at rcl ≲ 1.1 R, we check which effect the onset of clumping may have on the calculated line profile.

First, we assume a one-component wind (D = 10, d = 0) and adopt the rcl = 1.1 R. As a result, absorption dips appear close to the laboratory wavelength of both P v doublet components (lower panel of Fig. 6). To get rid of these sharp absorptions, we set the interclump density to d = 0.1. The result (lower panel of Fig. 6) shows that the absorption dips almost disappear, but still the level of absorption is not fitted well. However, after increasing the interclump density to d = 0.2, the absorption totally disappears, and the level of absorption is reproduced (lower panel of Fig. 6).

The reason for this effect is that the interclump medium above rcl shields the lower, smooth part of the wind. If the former is dense enough (as for d = 0.2), both absorption and emission are strong there to hide the layers below.

If we compare the solid black lines in the upper and lower panels of Fig. 6, which only differ by rcl, we cannot say which line fits the observations better, since both agrees with observations. Even if we set rcl = 1.3 R and a slightly different value of d, the agreement with the observed line remains similar. With appropriate values of d, the P v line profile can thus be fitted equally well regardless of whether clumping starts from the surface of the star or a bit above. The interclump medium hides the spectral signature of the onset of clumping. A nonvoid interclump medium has to be assumed always to fit the overall shape of the P v profile.

For our final models presented in this paper we assume that clumping starts at the base of the wind, for both the PoWR and the Monte-Carlo simulations.

6.4. Dependence on the clumping factor D

The parameter D defines the density inside clumps with respect to the smooth wind density. To check how different values of D influence formation of the P v resonance line profile, we varied this parameter for one selected model (HD 14947). The value D = 10 gives a good fit to the observation (see left panel of Fig. 4). In Fig. 7 we now compare the profiles that results for different values of D between 5 to 400. While a slight preference exists for the fit with our standard value D = 10, the results depend only very little on that parameter. In our previous paper we demonstrated that enhancing the clump density parameter D leads to a more pronounced porosity effect (cf. Fig. 6 in Šurlan et al. 2012b). Since in this work we include the interclump medium, this dependence is apparently reduced.

6.5. Velocity dispersion inside clumps

The porosity effect in the line radiation transfer depends on the gaps in frequency between the line absorption of the individual clumps. This kind of porosity in the velocity coordinate is sometimes called “vorosity” (Owocki 2008). To investigate this effect, we study the impact of the velocity deviation parameter m (see Eq. (20) and Fig. 3 in Šurlan et al. 2012b), which allows the velocity inside the clumps to deviate from the monotonic dependence on radius. As indicated by hydrodynamic simulations, the velocity gradient is assumed to be negative there, while the ambient interclump medium moves monotonically according to the wind velocity law.

The effects of m on the P v line profile are shown in Fig. 8. From our modeling it follows that vorosity mainly affects the outer part of the wind, extending the absorption beyond ν and leading to a softening of the blue edge. This improves the fit in all stars in the sample. However, we do not find any significant reductions of the overall line strength. On the other hand, Sundqvist et al. (2010) find that vorosity is just as important as porosity. We cannot confirm their conclusion.

6.6. Clumping parameter degeneracies

In our models, we use several free parameters that describe the properties of clumping (L0, D, d, rcl, m). However, it happens that for some different combinations of clumping parameters we obtain similar or exactly the same P v line profiles. There is a question, of whether it is possible to break these degeneracies of parameters. This is quite a complicated problem, and if using a single line it is probably impossible. However, adding more lines to the analysis with different line opacities and different depths of formation may put some additional constraints on the adopted set of parameters. For example, Zsargó et al. (2008) show that the O viλλ 1032,1038 resonance doublet (created by Auger ionization by X-rays from O iv) could be used to characterize the inter-clump medium and break the degeneracy in the interpretation of fits to the P v doublet. However, their conclusion is based on a 1D wind model that assumes microclumping (consequently with fewer free parameters) and requires verification for the more general case of 3D geometry. It is an interesting problem, that goes beyond the scope of the present paper, but we intend to study it in future by simultaneous fitting of P v and O vi doublets using their ionization structure obtained from PoWR model, provided we find appropriate observational data. Even though it is difficult to disentangle the parameters of our model, our key conclusion is the accordance of different mass-loss rate diagnostics when macroclumping is included in the models.

6.7. Velocity law and P v ionization stratification

In our models we prefer the so called double-β law (cf. Sect. 3.2). The main motivation for adopting the double-β law is to get rid of the absorption dip close to the the blue edge of the profile, which appears notoriously in models with the standard β-law. With the double-β law, the persistent acceleration in the outer wind enhances the velocity gradient there, and thus reduces the line optical depth at the highest blueshifts. A physical reason for such dynamic behavior has been already suggested by Lucy & Abbott (1993), who speculated that such a persistent acceleration could arise from changes in the ionization structure. The mentioned blue-edge absorption is clearly seen in upper panel of Fig. 9.

In the lower panel of Fig. 9 we demonstrate the influence of the ionization stratification assuming constant qv = 1 and the ionization stratification that has been taken from the corresponding PoWR model. For the parametric range of the stars in our sample, the PoWR models predict that more than 80% of phosphorus is in ionization stage of P v, except close to the photosphere. The double-β law was employed in both cases. As can be seen, the decrease in the P v ionization fraction results in a shallower absorption in the outer wind.

thumbnail Fig. 6

Comparison between calculated P v and observed (thin solid-blue lines) line profiles of HD 14947. Upper panel: effect of the number of clumps. The purple line with crosses is calculated with d = 0, rcl = 1, m = 0.1, and L0 = 0.5. The green line with triangles and the red line with squares differ only by L0 = 0.2 and L0 = 0.1, respectively. The thick solid black line is calculated with rcl = 1, m = 0.1, L0 = 0.5, and d = 0.2. The dotted orange line is from the PoWR model and corresponds to a smooth wind. Lower panel: effect of the onset of clumping and the interclump medium density. The green line with triangles is calculated with d = 0, rcl = 1.1, m = 0.1, and L0 = 0.5. The red line with asterisks differs only by d = 0.1, while the thick solid black line is for d = 0.2 so it differs against the thick solid black line in the upper panel only by the different rcl = 1.1. The remaining clump parameters are fixed as given in Table 4.

Open with DEXTER

thumbnail Fig. 7

Dependence of the P v profile on the clumping factor D. All four calculations are for rcl = 1, L0 = 0.5, d = 0.2, and m = 0.1. The thick solid black line is for our standard value D = 10. The red line with asterisks shows the profile for D = 5, the green line with triangles for D = 50, and the dashed orange line for D = 400. The thin solid blue line is the observed spectrum.

Open with DEXTER

thumbnail Fig. 8

Upper panel: effects of the “vorosity” on the P v profile. All three calculations are for L0 = 0.5 and d = 0.2. The green line with triangles shows the profile for m = 0.01, the red line with asterisk for m = 0.3, and the thick solid black line for m = 0.1. Lower panel: P v line profiles calculated with standard β-law and constant ionization fraction qv = 1 (dotted-red line), compared to the simulation with double-β law and the ionization stratification from the corresponding PoWR model (thick solid black line). Both profiles are calculated for L0 = 0.5, d = 0.2, rcl = 1, and m = 0.1. The thin solid blue lines in the panels are the observed spectrum. λ1∞ and λ2∞ represent the wavelength associated with the assumed ν.

Open with DEXTER

thumbnail Fig. 9

Influence of velocity law and ionization fraction of P v on profiles of the P v resonance doublet. Profiles were calculated with d = 0.2, rcl = 1, m = 0.1, and L0 = 0.5 for the case of HD 14947. Upper panel: comparison of profiles of the P v doublet for standard β-law (the dotted red line) and double-β law (solid black line). Lower panel: comparison of profiles of the P v doublet for constant ionization fraction qv = 1 (the dotted red line) and ionization fraction following from the double-β law (the thick solid black line). The thin solid blue line is the observed spectrum.

Open with DEXTER

7. Discussion

7.1. Clumping in the inner wind

Time-dependent 1D hydrodynamic models of radiation-driven winds has always predicted that the line deshadowing instability grows only in the acceleration zone. Thus, strongly inhomogeneous structures were expected to develop only at radii rcl ≳ 1.3 R (Feldmeier et al. 1997; Runacres & Owocki 2002; Dessart & Owocki 2005). Recently, however, Sundqvist & Owocki (2013) have obtained structures that are already close to the wind base (rcl ≲ 1.1 R) when they accounted for the effect of limb darkening.

There are observational arguments to suggest that clumping may start close to the stellar surface (e.g., Puls et al. 2006). The early onset of clumping close to the stellar photosphere may be a consequence of the subsurface convection as investigated by Cantiello & Braithwaite (2011). The presence of X-ray emission very close to the photosphere (Waldron & Cassinelli 2007) is an argument to support these considerations.

In our last paper (Šurlan et al. 2012b), we also concluded that clumping sets in close to the wind base, because otherwise deep, unshifted absorption features should be seen in the P v resonance doublet, which are not observed. In the light of our present results, this conclusion is not all that firm anymore. As mentioned in Sect. 6.3, the observations can also be reproduced if we assume rcl = 1.3 R while setting the interclump medium parameter to d = 0.2. Obviously, the interclump medium reprocesses the radiation coming from the base of the wind, so it is not possible to tell if the lower wind is already clumped or not, at least not for the relatively dense supergiant winds investigated here.

7.2. Clumping in the outer wind

The hydrodynamic simulations of the line deshadowing instability mentioned above predict that the clumps persist to large distances from the star. This justifies our assumption that clumping extends to large radii, even to rmax ≳ 100 R. Smooth-wind models predict much deeper absorption in the blue part of the line profile than observed (see Fig. 8). Obviously, the effective opacity in the outer part of the wind that has been overestimated. Most importantly, we showed that this opacity can be effectively reduced by the porosity effect. Additionally, the correct ionization fraction of P v lowers the absorption a bit, and the double-β law distributes the line opacity more uniformly across the line profile.

Our 3D Monte-Carlo radiative transfer calculations were performed with constant clumping parameters throughout the wind, which implies that the clumps become larger and more separated from each other with growing distance from the star. We should note here that the clumping parameter D may in fact dependent on depth (see Puls et al. 2006).

7.3. P v ionization fraction

The importance of the P v ionization fraction has been pointed out by Crowther et al. (2002). The studies that state the mass-loss rate discrepancy (Massa et al. 2003; Fullerton et al. 2006) have used simplified methods for simulating the P-Cygni profiles of the P v resonance doublet, in particular the “Sobolev with exact integration” (SEI) method (Lamers et al. 1987). For this approach, the ionization fraction of the P v ion must be adopted, and it has been assumed to be constant throughout the wind and close to unity (i.e., all phosphorus is in the P v ground state). Our detailed non-LTE models show that the ionization fraction of P v is actually somewhat lower and can vary with radius.

Additionally, the ionization fraction of P v may be affected by X-rays. This influence was examined by Krtička & Kubát (2009), who show that the X-rays of the observed intensity cannot deplete the P v ionization fraction significantly. Still, Waldron & Cassinelli (2010) suggest that strong emission line radiation in the XUV energy band can significantly reduce the abundance of P v and thus explain the discrepant low mass-loss rates estimates. However, Krtička & Kubát (2012) show that if the XUV radiation were strong enough to reduce the ionization fraction of P v, it would also change the ionization balance of other elements and significantly reduce the wind driving force, chence also stellar mass-loss rates.

Here we have performed tests by including an X-ray field in the PoWR calculations, using the parameters of X-ray radiation as obtained from observations (Oskinova et al. 2006). We confirm the result found by Krtička & Kubát (2009) that the X-rays really have no effect on the ionization balance of phosphorus, especially on the abundance of P v. Similar conclusions are also drawn by Bouret et al. (2012).

7.4. Mass-loss rates

In principle, mass-loss rates through radiatively driven stellar winds can be predicted from adequate hydrodynamical models. Ideally, such models would yield and ν from a given set of stellar parameters (stellar luminosity, mass, radius, and chemical composition). Such codes have to accept severe approximations. In most of them, the radiative force is parameterized using the so-called force multipliers (see Castor et al. 1975; Abbott 1982). Some of these codes calculate the radiative force in detail from a list of spectral lines (Krtička & Kubát 2004, 2009, 2010).

Hydrodynamical stellar wind models that account for clumping are still missing. In a few test calculations, Krtička et al. (2008) and Muijres et al. (2011) studied the effect of clumping on the radiative force. Vink et al. (2000) performed Monte-Carlo calculations of the radiative force, also using detailed line lists. However, they assumed the velocity law, instead of a fully consistent hydrodynamical solution. Conveniently, they condensed their results into a fit formula, which is widely used as reference for mass-loss rates.

Therefore we also compare our mass-loss rates, which are consistent with the Hα emission and the unsaturated UV resonance doublet of P v, with the Vink formula (Table 3, last column). On average, our mass-loss rates are lower by a factor of two. However, one should keep in mind that our relies on the assumption that the clumping contrast is D = 10. Within some range of D, a simultaneous fit of Hα and P v may be possible as well, with somewhat different mass-loss rates  ∝ D−1/2.

8. Summary

For a sample of five O-type supergiants, we studied the effects of wind clumping on the mass-loss rate determination, simultaneously considering the Hα emission (and other Balmer and He ii lines) and the unsaturated resonance doublet of P v in the farultraviolet.

  • When accounting for macroclumping, it is possible tosimultaneously fit the Hα and the P v lines withthe same mass-loss rates.

  • The consistent fit is achieved when we simulate the P v resonance profile with our 3D Monte-Carlo code for the line radiation transfer in clumpy stellar winds. Obviously, the reported discrepancies between Hα and P v mass-loss rates were due to the inadequate treatment of clumping.

  • The mass-loss rates for our consistent fits are lower by a factor of 1.3 to 2.6, compared to the mass-loss formula by Vink et al. (2000).

  • In contrast to other studies, it was necessary neither to reduce the mass-loss rate by adopting an extremely high degree of clumping nor to assume a subsolar phosphorus abundance for our consistent fits.

  • The porosity that is needed to fit the P v lines implies that ~104 clumps populate the wind within 100 R at any given moment.

  • The velocity dispersion inside the clumps has a moderate effect on the porosity of the wind, hence on the P v profile. The lower this dispersion, the smaller is the effective line opacity.

  • Compared to the standard β-velocity law, the double-β law improves the detailed fit of the P v line profile. It smooths the blue absorption edge and removes the absorption dip close to that edge.

  • With the detailed ionization stratification of P v from the PoWR code, better agreement with the observed P v line profile can be achieved than with qv ≡ 1.

Our results emphasize that an adequate treatment of the line formation in inhomogeneous winds is a prerequisite for interpreting O-star spectra and determining mass-loss rates.


1

IRAF is distributed by the National Optical Astronomy Observatories, which are operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the National Science Foundation.

Acknowledgments

Some of the data presented in this paper were obtained from the Mikulski Archive for Space Telescopes (MAST). STScI is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. Support for MAST for non-HST data is provided by the NASA Office of Space Science via grant NNX09AF08G and by other grants and contracts. The authors would like to thank Dr. Petr Škoda for securing a spectrum of HD 192639 in the Hα region, and to night assistants Mr. Jan Sloup and Mrs. Lenka Kotková for their help with obtaining spectra used in this paper. The authors thank the anonymous referee for his comments which helped to improve the paper. This work was supported by grant GA ČR P209/11/1198. B.Š. acknowledges Ministry of Education and Science of Republic of Serbia, which supported this work through project 176002 “Influence of collisions on astrophysical plasma spectra”. A.A. acknowledges financial support from the research project SF0060030s08 of the Estonian Ministry of Education and Research. L.M.O. acknowledges support from DLR grant 50 OR 1302. A.F.T. acknowledges the financial support from the Agencia de Promoción Científica y Tecnológica (Préstamo BID PICT 2011/0885), PIP 0300 CONICET, and the Programa de Incentivos G11/109 of the Universidad Nacional de La Plata, Argentina. Financial support from International Cooperation of the Czech Republic (MŠMT, 7AMB12AR021) and Argentina (Mincyt-Meys, ARC/11/10) is acknowledged. The Astronomical Institute Ondřejov is supported by the project RVO:67985815.

References

Online material

Appendix A: Spectral fits

thumbnail Fig. A.1

Best fit from PoWR modeling (thick red solid lines) to the observed HD 66811 spectra (thin blue solid lines), together with the calculated P v line profile from 3D Monte-Carlo code (black dotted line). Blue labels with numbers in the upper panels are UBVJHK magnitudes.

Open with DEXTER

thumbnail Fig. A.2

Same as Fig. A.1, but for HD 15570.

Open with DEXTER

thumbnail Fig. A.3

Same as Fig. A.1, but for HD 14947.

Open with DEXTER

thumbnail Fig. A.4

Same as Fig. A.1, but for HD 210839.

Open with DEXTER

thumbnail Fig. A.5

Same as Fig. A.1, but for HD 192639.

Open with DEXTER

All Tables

Table 1

Input stellar and wind parameters and element abundances by mass fraction.

Table 2

Optical, FUV, and NUV observation logs.

Table 3

Finally adjusted stellar and wind parameters.

Table 4

Fixed model parameters used in the 3D Monte-Carlo code.

Table 5

Clumping parameters that give the best fit to the observed P v line profiles.

All Figures

thumbnail Fig. 1

Dependence of wind velocity on radius. Comparison of the standard (single) β-law (dotted red line) and double-β law (solid black line).

Open with DEXTER
In the text
thumbnail Fig. 2

Snapshot from a 3D Monte-Carlo simulation, showing in 2D projection an example of a clump realization (red full circles) with respect to star (yellow dashed circle) and the path of a line-scattered photon (blue dotted line).

Open with DEXTER
In the text
thumbnail Fig. 3

The upper and the middle panels show the synthetic spectra (dotted black lines) of HD 66811 (left panels) and HD 15570 (right panels) obtained from the PoWR models, i.e. without macroclumping. Thin solid blue lines are observed spectra. The dashed green lines in the middle panels are from the same model, but only accounting for the photospheric lines while wind lines from P v and Si iv are suppressed. These photospheric spectra are used as input for the 3D Monte-Carlo calculations with macroclumping (lower panels). The parameters of these models are given in Tables 15.

Open with DEXTER
In the text
thumbnail Fig. 4

Same as Fig. 3 but for HD 14947 (left panels) and HD 210839 (right panels).

Open with DEXTER
In the text
thumbnail Fig. 5

Same as Fig. 3 but for HD 192639.

Open with DEXTER
In the text
thumbnail Fig. 6

Comparison between calculated P v and observed (thin solid-blue lines) line profiles of HD 14947. Upper panel: effect of the number of clumps. The purple line with crosses is calculated with d = 0, rcl = 1, m = 0.1, and L0 = 0.5. The green line with triangles and the red line with squares differ only by L0 = 0.2 and L0 = 0.1, respectively. The thick solid black line is calculated with rcl = 1, m = 0.1, L0 = 0.5, and d = 0.2. The dotted orange line is from the PoWR model and corresponds to a smooth wind. Lower panel: effect of the onset of clumping and the interclump medium density. The green line with triangles is calculated with d = 0, rcl = 1.1, m = 0.1, and L0 = 0.5. The red line with asterisks differs only by d = 0.1, while the thick solid black line is for d = 0.2 so it differs against the thick solid black line in the upper panel only by the different rcl = 1.1. The remaining clump parameters are fixed as given in Table 4.

Open with DEXTER
In the text
thumbnail Fig. 7

Dependence of the P v profile on the clumping factor D. All four calculations are for rcl = 1, L0 = 0.5, d = 0.2, and m = 0.1. The thick solid black line is for our standard value D = 10. The red line with asterisks shows the profile for D = 5, the green line with triangles for D = 50, and the dashed orange line for D = 400. The thin solid blue line is the observed spectrum.

Open with DEXTER
In the text
thumbnail Fig. 8

Upper panel: effects of the “vorosity” on the P v profile. All three calculations are for L0 = 0.5 and d = 0.2. The green line with triangles shows the profile for m = 0.01, the red line with asterisk for m = 0.3, and the thick solid black line for m = 0.1. Lower panel: P v line profiles calculated with standard β-law and constant ionization fraction qv = 1 (dotted-red line), compared to the simulation with double-β law and the ionization stratification from the corresponding PoWR model (thick solid black line). Both profiles are calculated for L0 = 0.5, d = 0.2, rcl = 1, and m = 0.1. The thin solid blue lines in the panels are the observed spectrum. λ1∞ and λ2∞ represent the wavelength associated with the assumed ν.

Open with DEXTER
In the text
thumbnail Fig. 9

Influence of velocity law and ionization fraction of P v on profiles of the P v resonance doublet. Profiles were calculated with d = 0.2, rcl = 1, m = 0.1, and L0 = 0.5 for the case of HD 14947. Upper panel: comparison of profiles of the P v doublet for standard β-law (the dotted red line) and double-β law (solid black line). Lower panel: comparison of profiles of the P v doublet for constant ionization fraction qv = 1 (the dotted red line) and ionization fraction following from the double-β law (the thick solid black line). The thin solid blue line is the observed spectrum.

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

Best fit from PoWR modeling (thick red solid lines) to the observed HD 66811 spectra (thin blue solid lines), together with the calculated P v line profile from 3D Monte-Carlo code (black dotted line). Blue labels with numbers in the upper panels are UBVJHK magnitudes.

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

Same as Fig. A.1, but for HD 15570.

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

Same as Fig. A.1, but for HD 14947.

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

Same as Fig. A.1, but for HD 210839.

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

Same as Fig. A.1, but for HD 192639.

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.