Exploring the influence of different velocity fields on Wolf–Rayet star spectra Free

ABSTRACT

Given their strong stellar winds, Wolf–Rayet (WR) stars exhibit emission line spectra that are predominantly formed in expanding atmospheric layers. The description of the wind velocity field |$\upsilon (r)$| is therefore a crucial ingredient in the spectral analysis of WR stars, possibly influencing the determination of stellar parameters. In view of this, we perform a systematic study by simulating a sequence of WR-star spectra for different temperatures and mass-loss rates using β-type laws with 0.5 ≤ β ≤ 20. We quantify the impact of varying |$\upsilon (r)$| by analysing diagnostic lines and spectral classifications of emergent model spectra computed with the Potsdam Wolf–Rayet (PoWR) code. We additionally cross-check these models with hydrodynamically consistent – hydro – model atmospheres. Our analysis confirms that the choice of the β exponent has a strong impact on WR-star spectra, affecting line widths, line strengths, and line profiles. In some parameter regimes, the entire range of WR subtypes could be covered. Comparison with observed WR stars and hydro models revealed that values of β ≳ 8 are unlikely to be realized in nature, but a range of β values needs to be considered in spectral analysis. UV spectroscopy is crucial here to avoid an underestimation of the terminal velocity |$\upsilon _\infty$|⁠. Neither single- nor double-β descriptions yield an acceptable approximation of the inner wind when compared to hydro models. Instead, we find temperature shifts to lower T_2/3 when employing a hydro model. Additionally, there are further hints that round-lined profiles seen in several early WN stars are an effect from non-β velocity laws.

1 INTRODUCTION

Massive stars (with |$M_\ast \gtrsim 8\, \mathrm{M}_\odot$|⁠) are known to have a large influence on the energetics of their host galaxies via ionizing radiation, strong stellar winds, and end-of-life supernovae. As a subset of massive stars, the massive Wolf–Rayet (WR) stars have exceptionally powerful stellar winds leading to typical mass-loss rates of |$\dot{M} = 10^{-5}\, \dots 10^{-4}\, \mathrm{M}_\odot \, \mathrm{yr}^{-1}$| (e.g. Crowther 2007), severely impacting their stellar evolution and enriching their circumstellar media. The emission lines in this stellar wind are typically optically thick, acting as a strong driver for the stellar wind and causing a large outflow of material. This material typically renders the hydrostatic surface of the WR star itself invisible and causes most, if not all, observed light to originate from the stellar wind of the star. Due to the effect of line-driving, a profound understanding of the stellar wind is required to correctly derive the properties of a WR star.

The WR stars were first discovered and recognized by their peculiar spectra, containing a multitude of strong and broad emission lines, causing them to be assigned to a separate spectral class. Based on the presence of certain elements, the WR star spectra are further divided into subclasses: the WN stars with strong nitrogen lines in the spectra, WC stars with prominent carbon lines and the WO stars with strong oxygen features. WR stars are in most cases hydrogen-depleted or hydrogen-free, having lost their outer layers through continuous winds or eruptions (e.g. Smith & Willis 1982; Conti et al. 1983; Smith 2014) or binary mass-transfer (e.g. Paczyński 1967; Vanbeveren, De Loore & Van Rensbergen 1998; Neugent & Massey 2014; Shenar et al. 2020).

About |$\sim 90~{{\ \rm per\ cent}}$| of the known WR star population consists of massive stars that evolved past the main sequence and are typically core-helium burning (Hainich et al. 2014; Shenar et al. 2019). These evolved WR stars, called the classical WR stars, are thought to be progenitors of stellar-mass black holes (e.g. Woosley, Heger & Weaver 2002). In contrast, there is a second group of very massive WR stars, displaying significant fractions of hydrogen in their spectra and still reside on, or close to, the main sequence. Thus, they are often considered to be core-hydrogen burning (de Koter, Heap & Hubeny 1997).

To infer the parameters of WR stars, sophisticated model atmospheres accounting for the expanding atmosphere are required, such as CMFGEN (e.g. Hillier 1990; Hillier & Miller 1998), FASTWIND (e.g. Santolaya-Rey, Puls & Herrero 1997), or PoWR (e.g. Gräfener, Koesterke & Hamann 2002; Hamann & Gräfener 2003; Sander et al. 2015b). These codes treat the wind in a non-LTE and chemically homogeneous framework, where the statistical rate equations and the full radiation transfer is computed instead of using the LTE approximations. Given the numerical demands of the models, the velocity structure in the stellar wind is commonly pre-specified in the form of a so-called β law. In the framework of line-driven winds, the β law can be obtained from a solution of the (modified) CAK-theory (Castor, Abbott & Klein 1975; Friend & Abbott 1986). It can be expressed as

where |$\upsilon _\infty$| denotes the terminal velocity of the wind, R_* the stellar radius, and r ∈ [R_*, +∞) the radius in the stellar wind. The name-giving parameter β is an exponent describing the steepness of velocity increase. In the case of O-type stars, Pauldrach, Puls & Kudritzki (1986) derived values of β ≈ 0.8. All model codes mentioned above allow their user to pick their own value of β. Hence, the value is often used as another free parameter to improve the overall fit, e.g. by trying to reproduce the H α-profile in B supergiants (e.g. Searle et al. 2008). For WR stars, the value of β = 1 is commonly adopted (e.g. Hamann, Schmutz & Wessolowski 1988; Hillier 1988; Hillier & Miller 1999), and is also standard in pre-calculated model grids (e.g. Hamann & Gräfener 2004; Sander, Hamann & Todt 2012; Todt et al. 2015). However, this assumption has also been challenged. Spectral modelling of observations carried out by Robert (1994) resulted in better fits of WR-star spectra for velocity laws with β > 3. Studies by Lépine & Moffat (1999) suggest values of β ∼ 5–20, using the temporal evolution of line profile variations interpreted as clumping in the stellar wind. The latter values probe mostly the outer regions of the stellar wind, where the continuum is optically thin. Thus, one could argue that these high values might only apply to this part of the wind and a lower β could be necessary at smaller radii. Interestingly, the analysis of the sonic point by Nugis & Lamers (2002), coming from a completely different angle, also favours β values of 5 and larger.

Instead of pre-specifying a velocity law in WR-star models (e.g. Gräfener & Hamann 2008), studies by Gräfener & Hamann (2005), Sander & Vink (2020), and Sander, Vink & Hamann (2020) have constructed wind models where the velocity field is computed from solving the hydrodynamic equation consistently with the radiative transfer. However, hydrodynamically consistent studies of classical WR stars are computationally expensive and their resulting velocity fields cannot be easily parametrized (Gräfener & Hamann 2005; Sander et al. 2015a; Sander et al. 2020). From the modelling of an early-type WC star, Gräfener & Hamann (2005) suggested a double-β law

(as proposed in Hillier & Miller 1999) with β₁ = 1 and β₂ = 50 as an approximation of the hydrodynamic result. The parameter q is a weighting factor between 0 and 1. This could qualitatively still be in line with the results of Lépine & Moffat (1999). This double-β law is used to apply two separate β laws to the inner wind and outer wind regime separately. As most of the observable line spectrum of a WR star is formed in its outer wind, a full parameter study of the double-β law is unnecessary as the emergent spectral imprint of the inner β law will be minimal.

Due to its reference radius R_*, equation (1) is connected to the stellar radius of a WR star. From equation (1), we would expect the velocity to be zero there, meaning that R_* should correspond to a hydrostatic radius. Typically, the stellar surface radius is defined at the point where the Rosseland mean optical depth τ_Ross ≈ 2/3. Because of the optically thick winds of WR stars, this point – which we denote as R_2/3 with a corresponding effective temperature T_2/3 – is usually far out in the wind of the star, where we have a significant bulk movement away from the actual star and the hydrostatic equation is significantly violated. Instead, subsonic velocities are only reached much further inside at τ_Ross ≈ 20 (e.g. Nugis & Lamers 2002; Sander & Vink 2020) although some hydrodynamical models also yield values of τ_Ross ≈ 5 (e.g. Gräfener & Hamann 2005). To be safely in the hydrostatic regime, the inner boundary of PoWR models, e.g. in the published WR grids (Hamann & Gräfener 2004; Todt et al. 2015) is traditionally set to τ_{Ross, cont} = 20. The corresponding radius to this (continuum) optical depth is termed as R_* with a formally defined effective temperature T_* following from the Stefan–Boltzmann equation. Hence, for PoWR models, T_* is typically an input parameter, while T_2/3 is an output of a converged atmosphere model.

In stellar structure models, the effect of optically thick winds cannot be taken into account accurately as this would require a hydrodynamic treatment and the computation of the flux-weighted opacity, which deviates significantly from the Rosseland opacity in an expanding non-LTE environment (illustrated e.g. in Sander et al. 2020; Sander 2022). This makes it difficult to empirically compare derived temperatures T_* for WR stars in the Hertzsprung–Russell diagram (HRD). Some evolution codes provide an estimate of the wind correction, but these can be severely off compared to actual atmosphere calculations as shown by Groh et al. (2014). However, also the usage of T_* raises questions as indicated in Fig. 1. For He-burning classical WR (cWR) stars, the large spread in temperature is not expected, especially if the stars are almost or completely depleted in hydrogen (e.g. Ekström et al. 2012).

Rather, these H-stripped stars should be located close to the helium zero-age main sequence (He-ZAMS). This discrepancy between observations and structure models is referred to as the Wolf–Rayet radius problem (e.g. Grassitelli et al. 2018). As we will show later on in this work, the T_2/3-temperatures, which can be directly inferred from spectral analyses, are less affected by uncertainties in R_*. However, they do not reflect the hydrostatic layers of the stars. This is illustrated in the HRD from Hamann et al. (2019) now with T_2/3 instead of T_* in Fig. 2. In the T_2/3-HRD, almost no star of the WR population is located on the He-ZAMS, even for the completely hydrogen-free stars, emphasizing that these ‘observable’ T_2/3 temperatures only reflect the extended wind layers and thus are problematic for any general comparison with evolutionary tracks.

Two major solutions are suggested: either the actual hydrostatic radii from where the winds are launched are significantly smaller than derived from the empirical models using β = 1 velocity laws, such as obtained in the prototypical calculations of hydrodynamically consistent atmosphere models (Gräfener & Hamann 2005; Sander & Vink 2020; Sander et al. 2020), or rather the hydrostatic radii are inflated (e.g. Petrovic, Pols & Langer 2006; Gräfener, Owocki & Vink 2012) and the winds are actually launched at approximately the radii determined from the empirical models. In fact, both phenomena might exist in nature with a branching depending on the precise stellar parameters Grassitelli et al. (2018); Poniatowski et al. (2021).

In this work, we will investigate how much the empirically derived results are affected by the common choice of β = 1 for the velocity law. Moreover, we will investigate whether the use of different β-type velocity laws could bring the empirical values of T_* closer to the predictions from stellar structure theory. Following up on the early suggestion by Hillier (1991) that the β law may result in ambiguous parameters such as the stellar radius R_*, we probe the impact of different β values for a range of mass-loss rates and temperatures. Hillier (2003) concluded that the β parameter did not have a significant impact on the line spectrum when varying β between 1 and 3, but as we will show, this does not hold when considering larger β values. Moreover, we compare the β model efforts to trends from hydrodynamically consistent models and perform an exemplary approximation of a hydrodynamically consistent solution to gain first insights on the wind acceleration, T_2/3 and spectral appearance, preparing us for more tailored investigations of individual targets with hydrodynamically consistent models in the future.

The paper is organized as follows: the methodology and modelling used to conduct this study is explained in Section 2. In Section 3, we present the comparison and classification of the emergent spectra from our β-law models and further discuss the outcomes and implications in Section 4. Afterwards, we compare our findings with results from hydrodynamically consistent modelling in Section 5. A summary of our work, along with future prospects, is given in Section 6.

2 METHODS

The strong coupling between radiation field and population numbers cause the modelling of WR-star winds to be numerically demanding as the local thermodynamic equilibrium (LTE) approximation cannot be realistically assumed. The Potsdam Wolf–Rayet (PoWR) code (Gräfener et al. 2002; Hamann & Gräfener 2003; Sander et al. 2015b) is one of the few codes capable of modelling these stellar winds in a non-LTE setting in the comoving frame, assuming a spherically symmetric, stationary outflow, and chemical homogeneity. Unless the hydrodynamic branch is used (Sander et al. 2017), the wind is described by the mass-loss rate |$\dot{M}$| and a pre-specified velocity field |$\upsilon (r)$| with the help of a β law with parameters β and |$\upsilon _\infty$|⁠. The code computes the atmospheric stratification and emergent model spectra for a given set of stellar parameters |$(M_\ast , T_\ast , R_\ast , L_\ast , \dot{M}, \upsilon _\infty , \beta)$|⁠. The emergent spectra are computed in the observer’s frame such that they can be properly compared to observed spectra.

To study the impact of different values of β, we calculated several model sequences employing a range of T_* values that maps the bulk of the empirically analysed Galactic WR stars (Hamann et al. 2019; Sander et al. 2019). The corresponding area in the Hertzsprung–Russell diagram (HRD) is illustrated in Fig. 1.

For a fixed velocity law, WR wind models can be characterized by two main values: the stellar temperature T_* and the transformed radius R_t (Schmutz, Hamann & Wessolowski 1989):

with D denoting the clumping factor (Hamann & Koesterke 1998). While having the dimension of a radius, the absolute value of R_t has no direct meaning in the sense of a physically relevant radius. Instead, R_t represents an emission measure introduced by Schmutz et al. 1989, who discovered that line strengths in stellar models with the same T_* are approximately equal for similar values of R_t (Schmutz et al. 1989), regardless of the individual values entering equation (3). This allows for a scaling of the models with regard to the luminosity of the object.

Instead of R_t one can also define the transformed mass-loss rate,

introduced by Gräfener & Vink (2013), which denotes the mass-loss rate the star would have if it had no clumping (D = 1), |$\upsilon _\infty = 1000$| km s⁻¹ and |$L_\ast = 10^6\, \mathrm{L}_\odot$|⁠. In a similar fashion as with R_t, also |$\dot{M}_\mathrm{t}$| can be seen as a measure for the general emission line strength (see e.g Sander & Vink 2020, for details). With the unit of |$\dot{M}_\mathrm{t}$| being slightly more intuitive, we will denote our WR-star models via their respective T_* and |$\dot{M}_\mathrm{t}$| values: letters A, B, and C denote T_*, while numbers 1, 2, and 3 denote |$\dot{M}_\mathrm{t}$|⁠, see e.g. Table 1, with higher numbers corresponding to a higher |$\dot{M}_\mathrm{t}$|⁠.

For each of our selected T_*-values, we pick the same three values of |$\dot{M}_\mathrm{t}$|⁠. In our selection, we deliberately avoid the regime of very dense winds where we already know that T_* becomes meaningless as the whole spectrum is generated in layers where |$\upsilon \approx v_\infty$| (Hamann & Gräfener 2004; Hamann, Gräfener & Liermann 2006).

We investigate three different types of WR models in this work: hydrogen-free WN stars, hydrogen-rich WNh stars, and carbon-rich WC stars. The first type of models are resembling cWR stars with strong nitrogen lines in their spectra and are for convenience referred to as of ‘WN type’ in the rest of this work. The WNh-star models have a significant amount of hydrogen at their surface (X_H = 0.5) and are denoted as WNh models. The carbon-rich WC models represent the hydrogen-free WC stars with carbon and oxygen surface mass fractions of X_C = 0.4 and X_O = 0.05. These are denoted as WC models.

For each of the three WR-types and each of the nine |$(T_\ast ,\dot{M}_\mathrm{t})$|-combinations, we calculate a sequence of models using β-velocity fields with β = 0.5, 1, 2, 4, 8 and 20, leading to a total of 27 model sequences or 262 models. An overview of their specific parameters is shown in Table 1. For the remainder of this paper, we will denote a model which uses the β-velocity law with β = x as a βx model. Hence, if β = 1, we denote it as a β1 model, for β = 4, we denote it as a β4 model, and so forth.

To quantify the differences between the resulting model spectra, we employ the spectral classification criteria and diagnostics from Smith, Michael M & Moffat (1996) for WN stars and from Crowther, De Marco & Barlow (1998) for WC stars. An overview of the diagnostic spectral lines used in this work is given in Table 2. These lines are especially sensitive to changes in stellar parameters. Both Smith et al. (1996) and Crowther et al. (1998) preferably use peak-to-continuum values – the spectral line peaks in continuum units subtracted by 1 – of spectral lines as indicators for line strengths and not the equivalent widths for emission lines. This is due to the often blended emission lines in WR stars, caused by heavy broadening of the lines due to large Doppler shifts in the stellar wind. Using the normalized spectra of our models, we automatize the classification scheme accounting for the broadening of spectral lines and the often accompanying line-blending. Based on the line strength ratios, the spectra are classified from the early to late subclasses WN2 to WN9 or WC4 to WC11. For the comparison with hydrodynamically consistent models, we make use of models calculated in the framework of Sander & Vink (2020).

3 RESULTS FROM THE β SEQUENCES

An overview of the different β-velocity laws from our model sequences is given in Fig. 3 with the β1-velocity law highlighted. The velocity fields in Fig. 3 are normalized to their (terminal) value at infinity. As clearly evident from the figure, the β parameter from equation (1) reflects the steepness of the velocity field, meaning that with higher values of β, a given velocity is only reached further out in the wind.

3.1 Line strengths

To convey the influence of the β parameter, we compare the emergent spectra of our atmosphere models along each β-sequence. With this term, we denote the sequences of models that are identical in their stellar parameters and only differ in the β value applied for their wind description. The influence of β is most obvious when studying the diagnostic lines. An example is shown in Fig. 4 where we plot the spectra of the WN B2 models (see Table 1 for the corresponding stellar parameters) in the N iii-v|$\lambda \, 4640 \text{-} 04$| Å line region. The change in line strengths (peak-to-continuum ratios) is immediately apparent. In the β0.5 and β1 cases, the N iii|$\lambda \, 4640$| Å line is too weak to be distinguished from the continuum and the N v|$\lambda \, 4604$| Å line dominates. While emergent spectra with β2 and β4 laws show the N iii and N v lines being of comparable strength, the N iii line dominates while the N v line diminishes for β8 and β20 laws. This implies a strong decrease in the |$\rm{N\,\,\small {\rm V}}\, \lambda \, 4604\, \big /\, \rm{N\,\,\small {\rm III}}\, \lambda \, 4640$| line-strength ratio with increasing β. While the changes in line strength for low β values of 0.5 to 2 are small – similarly to Hillier (1991) – and may fall below the detection limit in real observations, more extreme β values (4 to 20) lead to very significant changes. The remaining diagnostic lines and accompanying ratios are affected in a similar way. Analogous to the WN-star models, the C iii|$\lambda \, 5696$| and C iv|$\lambda \, 5808$| lines are affected in WC stars. We illustrate this in Fig. 5, where we present the WC B2 models with different β values, highlighting that the C iii line increases in strength while the C iv line diminishes when increasing β. This causes the |$\rm{C\,\,\small {\rm IV}}\, \lambda \, 5808\, \big /\, \rm{C\,\,\small {\rm III}}\, \lambda \, 5696$| line-strength ratio to drastically change as well, with similar trends visible for other diagnostic lines and line strength ratios.

3.2 Continuum

In addition to the line spectrum, also the continuum itself is affected by changes of the velocity law. As an example, we take the B2 WN models (similarly to Fig. 4) and compare the continua for different β, shown in the upper panel of Fig. 6. The figure shows the continuum (without the lines) at a distance of |$10\, \mathrm{pc}$| and is similar in shape to how it would be observed from Earth. For all models in Fig. 6, the total (bolometric) luminosity is the same. The initial impression seems to be that the continuum shifts to cooler temperatures with increasing β, an effect similar to the results above for the line strength ratios. Interestingly, the continuum of all models aligns at the He i ionization edge (⁠|$504.3\,$| Å) in the extreme UV (EUV). When looking towards longer wavelengths in the UV, optical, and infrared (IR) regimes, we further see that the flux is higher for higher values of β. In the optical regime, the shift is mostly wavelength-independent, corresponding to a simple shift in log–log space or a multiplication factor in linear units (see Section 4.2 for further discussion). These continuum differences start to diminish for longer wavelengths, eventually converging around 100 |$\mu$|m. When considering the full flux-calibrated spectrum (bottom panel of Fig. 6), the differences are less obvious in the UV, but the overall trends remain. In the EUV, the differences become even more pronounced compared to the raw continuum case. All of our exemplary models, independent of β, are optically thick below the He ii ionization edge, meaning that essentially no photons from this region escape from the stellar winds. However, the He i ionizing flux is notably affected by the choice of the velocity field. For β = 0.5, the number of He i ionizing photons is about an order of magnitude higher than in the β20 case. For the H i ionizing flux (⁠|$\lambda \lt 912\,$| Å), the trend is more moderate with a decrease of only ≈20 per cent. The trend is reversed for the H₂-dissociating Lyman–Werner photons, which are emitted between |$912\,$| and 1110 Å. For them, we obtain an increase of about a factor of 2 between β = 0.5 and 20. However, the total photon budget below 1110 Å remains approximately constant due to the different trends almost cancelling each other in our considered model sequence.

3.3 Spectral classification

The spectral sub-classification for WR stars uses diagnostic line strengths and line strength ratios. With the strong influence on these lines and ratios by the β parameter, the spectral subclasses are affected as well. In Fig. 7, we show the classification – here based on Smith et al. 1996 – of nine WN-star models with different stellar parameters (see Table 1) as a function of the β parameter. For observed WR stars, the early and late-type spectral subclasses are commonly interpreted as an indicator for the stellar temperature T_* at the hydrostatic layers: early-type classes are associated with hotter WR stars and later-type classes with cooler WR stars. In our β study, we now find for all nine WN-star model sequences, the spectra are classified as a later subtype when increasing the β parameter compared to their classification when β = 1 is used.

As an example, the β1 B2 model in Fig. 7 is classified as an early-type WN3 star, while the β20 model is classified as a late-type WN8 star despite the stellar parameters being identical, save for β. Moreover, the β20 C1 and β1 A3 models are both classified as intermediate-type WN5 stars, despite a temperature difference of ΔT_* ≈ 40 kK between the two models. Given that large model grids are so far only available for β1 laws, one could associate average temperatures with accompanying spectral subclass. This has been done in the right-hand side axis in Fig. 7, where we average over the results from Hamann et al. 2006, employing β1 models. An overview of these average temperatures is also shown in Table 3. Our compiled T_*(β = 1) values essentially reflect which temperatures would be inferred for a specific spectral subclass when employing only β1 models. Assuming that β values are also distributed over a range of values in nature, this implies that there is no coherent relation between spectral subtypes and hydrostatic temperatures of WR stars.

In addition to the nine model sequences representing the bulk of observed WN stars, we further calculated a series of additional (denoted as high T_* in Fig. 7) models with different β laws representing a WN star close to the He-ZAMS (cf. Fig. 1). Motivated by the steady-state model from Poniatowski et al. (2021), the stellar parameters of this model sequence are T_* = 130 kK, |$\log (\dot{M}_\mathrm{t}\, [\mathrm{M}_\odot \,\mathrm{year^{-1}}]) = -4.434$|⁠, |$R_\mathrm{t} = 2.0893\, \mathrm{R}_\odot$|⁠, |$\log (L_\ast \, [\mathrm{L}_\odot ]) = 5.416$| and |$\log (\dot{M}_\ast \, [\mathrm{M}_\odot \,\mathrm{year^{-1}}]) = -4.831$|⁠. Interestingly, the spectral classification of the high T_* series is less affected by the different β values in the stellar wind than the other, lower-T_* models, as evident from Fig. 7. Here, the different models are classified as WN2 for all β values with the exception β = 20. However, the value of |$\dot{M}_\text{t}$| in the model motivated by Poniatowski et al. (2021) is lower than in all of our other models, so we cannot rule out that the effect might be stronger for higher mass-loss rates (i.e. even more dense winds).

The same analyses have been performed on the WNh-type models, of which the results are shown in Fig. 8. While there are differences in the individual classifications, we see the same effect of increasing β as in the WN sequences: models with higher β values in their stellar winds tend to be classified as later-type stars. This thus leads to the same ambiguity between the apparent T_* and the β value.

Unsurprisingly, the classification of the WC-star models is affected as well. Applying a similar analysis as for the WN-type models in this work, but now employing the classification criteria from Crowther et al. (1998), the resulting subtypes are shown in Fig. 9. As evident from the figure, the model spectra tend to be again classified as later-type stars for higher β values. The right-hand axis in Fig. 9 shows the stellar temperatures per subclass (compiled from Sander et al. 2012 and Leuenhagen, Hamann & Jeffery 1996, see also Table 3). This again yields the same ambiguity in stellar temperature and the β parameter as with the WN-star models.

3.4 Line profile changes

Beside affecting peak heights, the use of different β values for |$\upsilon (r)$| also leads to different line widths. An example is shown in Fig. 10, where the He ii λ 5411 spectral line is shown for the WN B2 models. The velocities associated with the Doppler-shift and accompanying line broadening from the centre of the line are shown on top. The figure shows a strong decrease in both line-strength (peak-to-continuum value) and line width, despite the |$\dot{M}$| and |$\upsilon _\infty$| values being the same for the different spectra.

Spectral classification of WN stars sometimes also incorporates line widths, discerning between strong and weak-lined spectra. As a measure for the line widths, Smith et al. (1996) use the full width at half maximum (FWHM) of the He ii|$\lambda \, 4686$| line, where thresholds of FWHM > 30 Å and FWHM < 30 Å are used for broad-line and weak-lined spectra, respectively. While the measure as such is purely morphological, our study highlights that any immediate conclusions about the overall wind strength or mass-loss of the objects based only on emission-line widths are ambiguous. Without additional knowledge on their formation regions, terminal wind velocities measured from the width of emission lines need to be considered as lower limits. This is due to the strong influence of the velocity field on the line-widths. In case of β-type laws, weaker emission lines can even be accompanied by relatively high mass-loss rates if β ≫ 1.

3.5 UV spectra

The main lines used for classification of WR stars are located in the optical spectrum. At Galactic metallicities, the optical spectrum of WR stars is usually dominated by emission lines. This is not necessarily the case for the UV spectrum. Given the previously discussed ambiguities for deriving the wind parameters only from (optical) emission lines, we now take a look at the effect of different βlaw assumptions on prominent UV lines. As a case study, we take the C iv λλ 1548–50 doublet – a P-Cygni profile – shown in Fig. 11 for the WN B2 models, applying different β values. The P-Cygni profile in the figure shows a varying emission part, similar to the optical emission lines, while the blue-most wavelength of the P-Cygni absorption part remains unaltered when changing the β parameter. Similar results are obtained for other WN and WC star models.

To illustrate the ambiguity caused by varying the β parameter and the potential of UV data in removing it, we mimic an observation in the visual domain by degrading the resolution of an emergent PoWR spectrum. Additionally, Gaussian noise is added to simulate a finite signal-to-noise (SNR) ratio. The simulated observation used here is a WN B2 (see Table 1) β20 model with |$\upsilon _\infty = 1600$| km s⁻¹, degraded to a resolution of R = 3600 in similarity to observational data used in Hamann, Gräfener & LiermannHamann et al. (e.g. 2006) – originally from Hamann, Koesterke & Wessolowski (1995) – and with a simulated SNR = 100. We then proceed by comparing PoWR β1 models to the mock observation considering only the optical regime and compare the C iv λλ 1548–50 doublets of both the mock observation and the best-fitting model. The best-fitting β1 model is a WN model with T_* = 44.7 kK, |$\dot{M}_\mathrm{t} = -4.06$| and |$\upsilon _\infty = 600$| km s⁻¹. Both are shown in Fig. 12, along with a zoom-in on diagnostic lines and the C iv λλ 1548–50 profile.

While we do not aim at a fine-tuned analysis here – a few optical lines do not agree perfectly – the optical spectra in the figure resemble each other reasonably well in the context of a grid-based approach. Such an approach would be used when studying a larger sample or using a grid model to account for the WR imprint in an unresolved population. On the contrary, the P-Cygni profile shows a large discrepancy in the blue-most wavelength of the absorption trough. Moreover, there are significant differences in other stellar parameters between both models: a difference in T_* of 26.1 kK and a |$\dot{M}_\mathrm{t}$| difference of 0.3 dex (the true |$\dot{M}$| difference is 0.12 dex). Thus, beside underestimating |$\upsilon _\infty$| by 1000 km s⁻¹, the β1 model further severely underestimates T_* and underestimates |$\log \dot{M}$| of the simulated observation. We also note that the difference in T_2/3 between these models is 6.8 kK, noticeably smaller than the T_* difference and more in line with the similarity of the line spectra. Also the continuum differences are fairly small. While there is some discrepancy in the optical and infrared regimes, in a similar way as in Fig. 6, the effect is much smaller and amounts only to an apparent luminosity difference of ∼0.004 dex in the optical regime (while the real model luminosities are the same). Lastly, we want to highlight that the differences in the C iv λλ 1548–50 absorption troughs are not caused by the different velocity laws, in contrast to Fig. 11, but due to different |$\upsilon _\infty$| in the input of the models.

Implications for stellar feedback

Underestimations of the terminal velocity can also have implications for the derived stellar feedback. Regarding the ionizing feedback, the β1 model has ∼10 per cent more hydrogen ionizing flux than the mock observation. When considering the mechanical luminosities

the consequences of underestimating |$\upsilon _\infty$|⁠, e.g. due to considering only β1-velocity law, can be significant as well. In our example, the ∼35 per cent difference in |$\dot{M}$| and the 2.7 times higher |$\upsilon _\infty$| both need to be taken into account, resulting in an underestimation of L_mech in our mock observation in Fig. 12 by nearly an order of magnitude if the values from the β1 ‘fit’ are adopted.

3.6 Comparison to observations

To get a first look on whether our findings from the model calculations are realized in nature, we take three objects with available optical and UV observations as case studies: WR128, SMC AB2, and M31 WR096 (LGGS J004412.44+412941.7). These targets were chosen due to the availability of UV spectra as well as their difference in metallicities and galactic backgrounds (Milky Way, SMC, and M31, respectively). Their derived properties are shown in Table 4, where the data of the analysis of WR128 from Hamann et al. (2019) and of SMC AB2 from Hainich et al. (2015) are also given for comparison. In both cases, a β1-velocity law was used to determine the stellar parameters in the original work. For SMC AB2, we use the stellar parameters found by Hainich et al. (2015) as the best β1 representation. For WR128, the terminal velocity was not derived by custom modelling in Hamann et al. (2006) and we thus calculated a dedicated β1 model with an adjusted |$\upsilon _\infty$| for WR128. The same was done for M31 WR096, where no previous analysis has been performed. Without aiming for fine-tuned spectral fits for each of the targets, which would require a separate dedicated study, we calculate models with β ≠ 1 and compare them to the best β1 solution. The detailed results for three targets are discussed in Appendix  A. We give a short overview here as well.

For all three targets – WR128, SMC AB2 and M31 WR096 – β4 models seemed to give the best agreement with the observed spectra. Higher β velocity laws were increasingly difficult to reconcile with the observed spectra. Given the UV lines, the spectrum of WR128 is reproduced better with a β4 velocity law, while the optical regime showed similar agreement between the β1, β4 models and the observation. The β4 model for SMC AB2 shows a comparable quality to the β1 model. However, an unfeasibly large clumping factor of D = 130 was needed to attain this. Additionally, the observed UV spectrum is also recreated well with a β1-velocity law, suggesting the original β1 model should be used in this case. In a similar way to WR128, also for M31 WR096, the β4 model can be favoured based on the slightly better quality of the UV spectrum.

With these three targets we investigated whether models with β ≠ 1 could generally reproduce the spectrum in a similar quality as a ‘traditional’ β1 model. We found that for two of the three targets one could justify a β4-velocity law, while one target (SMC AB2) was best reproduced with the β1-model. As models with even higher β values yielded less good results, we thus conclude that stellar models with β > 20 will only rarely, if at all, provide a suitable representation of existing WR winds. This, fortunately, limits the impact of the aforementioned degeneracy with the β law. None the less, a range of β values around 1 should be explored when analysing WR stars with prescribed β law. To determine the ‘best’ β value, optical data are only of limited help. Instead, the saturated P-Cygni profiles in UV, which are still available for most WR stars, are important to not underestimate the value |$\upsilon _\infty$|⁠.

Considering the derived parameters in Table 4, we see that for stars analysed with different β models, the most striking differences are seen between the inferred T_* and |$\upsilon _\infty$| values, with noticeable but less severe discrepancies between the |$\dot{M}$| values. On the other hand, while the T_2/3 values are hardly affected by the differences in wind parameters, there can be significant changes in the luminosity, when other parameters are adjusted simultaneously with β. This is reflected most clearly in our analysis of SMC AB2 (see Section  A). As the other two studied targets are hardly affected, we conclude that different choices of β, within the acceptable values for observations, only have a mild effect on the resulting spectral energy distribution (SED). The similarity of the T_2/3 temperatures could be attributed to the motivation of getting the same spectral appearance. This would require that the shape of the velocity field hardly affects T_2/3, which seems to be valid for β-type velocity laws, but will not hold in general as we will show later in Section 5.

4 DISCUSSION: THE β APPROACH FOR WR STARS

The differences in the emergent spectra seen in Section 3 can physically be attributed to the changes in the different β-velocity fields, shown in Fig. 3, and to their coupling with densities and temperatures (see Section 4.2 and further) in the stellar wind.

4.1 Determination of terminal velocities from line widths

The change in diagnostic line widths, shown in e.g. Figs 4, 5, and 10, has implications in determining the correct |$\upsilon _\infty$|⁠. Physically, the different velocity fields can cause the line formation regions to shift closer to or further away from the star. In the case of the β law, velocity descriptions with higher β values have a shallower velocity increase. Hence, higher velocities are reached at larger radii from the star. This is demonstrated in Fig. 13 for the example of the He ii λ 5411 line, depicted in Fig. 10 where a measure of the line-formation strength ξ (from Hillier 1987b) is used. While the value of ξ does not reflect the intrinsic line strength, its distribution reflects where an emission line forms in the wind. The plot in Fig. 13 clearly shows that the maximum of xi is shifting to larger radii but lower wind velocities when increasing the β value. When the line emission happens at lower velocities, the corresponding lower Doppler shifts lead to less broadened spectral lines, as also evident from Fig. 10. Consequently, the derived |$\upsilon _\infty$| could be ambiguous if only inferred from the optical spectrum. For each of the β models in Fig. 13, also the radii and velocities where τ_Ross = 2/3 are indicated. We see the R(τ_Ross = 2/3) = R_2/3 consistently moving outwards for increasing β, while |$\upsilon (\tau _\mathrm{Ross}=2/3) = \upsilon _{2/3}$| peak at β = 1 with then an overall decreasing trend for higher β.

This problem can be resolved by using saturated P-Cygni lines in the UV-spectrum, most notably the C iv λλ 1548–50 doublet. As the emission component of these lines form relatively close to the stellar surface, the absorption troughs are expected to be strong enough to overcome the emission in the P-Cygni profile, although this emission can still vary strongly with different velocity laws. For example, the P-Cygni profile in Fig. 11 is saturated and thus gives a suitable estimate for |$\upsilon _\infty$| (Prinja, Barlow & Howarth 1990) independent of the velocity law used in modelling the stellar wind. Moreover, fixing |$\upsilon _\infty$| via the UV P-Cygni absorption trough allows to exploit the β dependence of the optical diagnostic lines to estimate the best β value. We briefly tested this in our observational examples, although a broader range of model calculations would be required to give a quantitative estimate of the best β value and its error margin.

For O-type stars, UV P-Cygni profiles have been used extensively to determine the terminal wind velocity as typically most distinctive wind signatures are found in the UV (e.g. Prinja et al. 1990). As WR-star optical spectra usually contain plenty of emission lines formed in the wind of the star, these have typically been used over UV P-Cygni profiles to constrain |$\upsilon _\infty$|⁠. The above results showing the optical line width ambiguity for the β parameter now stress the importance of UV P-Cygni lines also for WR stars.

4.2 Degeneracies in spectral modelling

In Section 3.3, we demonstrated that the use of different β laws can affect the diagnostic line and also the spectral classification. To better understand this, we show the emission line formation strength ξ (defined by Hillier 1987a) in Fig. 13. We see that for different β values, the line formation regions occur at different radii and velocities, also corresponding to different local electron temperatures. For increasing β values, the line formation regions will in general occur at lower temperatures due to the increasing radii of the line formation regions. Consequently, stellar models with the same stellar temperature T_* but a higher β law will appear cooler as the diagnostic lines originate from cooler line formation regions. This causes the spectral types of WR-star models to shift in general to a later-type class when the β value is increased while keeping the stellar parameters constant. Naturally, this effect enhances when more extreme values of β are adopted (e.g. β = 20 or higher). As a result, the entire range of spectral subclasses can be covered by changing β while the model stellar temperature is left unchanged. This effect is so extreme for WR stars due to the line spectrum being formed largely or even completely in the stellar wind.

The model continuum is affected by different velocity laws as well (cf. Section 3.2). The continuum peaks shift to longer wavelengths, implying a cooler T_2/3. Due to the different fluxes of the continua, an additional ambiguity can arise in estimating the luminosity. When for example reproducing optical photometry with the B2 WN models, a factor of 4.3 is required to match the β1 flux to the β20 flux. This is illustrated in Fig. 14, where we look at the optical part of the flux-calibrated spectrum of both models. Taking the factor of 4.3 at face value, this would correspond to a significant difference of 0.6 dex in the derived luminosity. As we did not obtain good spectral fits with extreme values of β, the actual impact will usually be smaller, but can still be on the order of 0.3 dex (cf. Section 3.6).

The obtained differences in the continuum further raise the question whether our derived changes in the normalized line spectrum are intrinsic to the lines or arise from different continuum fluxes. The inset in Fig. 14 focuses on a few diagnostic lines in the optical spectrum for two different flux-calibrated models. For comparison, we further show the β1 model scaled to the level of the β20 flux. It is evident, that the lines are not diminished by the stronger continuum, but intrinsically change their strength and shape. Given that we see similar effects as in the normalized spectra (see e.g. Fig. 10), we conclude that most of the line changes originate due to the different velocity fields and cannot solely be attributed to normalization effects. A study by Ignace (2009) was able to show that continuum changes originating from wind structure differences can influence line widths in hot and dense stellar winds. Due to the different β-law assumptions, it is hard to exclude the possibility of additional continuum effects in our work.

We demonstrated previously that – at least to some degree – the changes to the spectrum due to a higher β value can be mitigated by changing other stellar and wind parameters of the model, causing a degeneracy in the solutions to an observed spectrum. This degeneracy acts on top of the already known degeneracy in the regime of very dense winds (see e.g. Hamann & Gräfener 2004). There, the same optical spectra (from models with fixed β = 1) can be obtained for different choices of R_* (or correspondingly T_*) if |$\log (R_\text{t}\, [\mathrm{R}_\odot ]) \le 0.4$|⁠. Our findings for the different β values now yields an additional parameter degeneracy on top of this effect.

4.3 Implications for the WR radius problem

Via the mass continuity equation |$\dot{M} = 4\pi r^2 \rho (r)\upsilon (r)$|⁠, we can directly connect the local wind density ρ(r) and velocity |$\upsilon (r)$|⁠. For the standard assumption of a constant mass-loss rate |$\dot{M}$|⁠, lower velocities directly imply higher wind densities. This causes wind models with β > 1 velocity laws to have a more extended higher density wind region. As a result, the optically thick region of a WR-star wind grows larger and hence the apparent ‘photospheric’ radius, i.e. R_2/3, does as well.

This radius dependency raises the question whether the application of higher β values for the velocity fields in the empirical spectral analysis of WR stars may provide a solution to the WR radius problem, see also Section 1. Assuming that in general WR stars have winds that could be modelled better using β > 1, the stellar temperatures derived with β = 1 would be underestimated. This would shift the observations in Fig. 1 closer to the He-ZAMS. However, to completely resolve the large temperature discrepancies between the He-ZAMS and the stellar temperatures derived in the literature, rather extreme values of β are required. For example, the temperature discrepancy between the β20 mock observation and the best-fitting β1-model in Fig. 12 adds up to 26.1 kK, causing a considerable temperature shift in the HR-diagram. Still, even with a β20 law, the model does not reach the He-ZAMS. Hence, to shift the coolest classical WR stars to stellar temperatures close to the He-ZAMS, either more extreme β values would be required or the β law description is generally insufficient to achieve this goal. In Section 3.6 and Appendix  A, we show that extreme values of β > 20 do not seem to be consistent with observations. This point will be discussed further in Section 5.

5 COMPARISON WITH HYDRODYNAMICALLY CONSISTENT MODELS

In the previous sections, we demonstrated that β-type velocity laws are able to yield a reasonable reproduction of the observed spectrum, but do not yield parameters for the deeper, hydrostatic layers that are in line with evolutionary expectations. Dynamically consistent atmosphere calculations (Gräfener & Hamann 2005; Sander et al. 2020) overcome this discrepancy by obtaining |$\upsilon (r)$| directly from solving the hydrodynamic equation of motion

with a denoting the isothermal sound speed, g the gravitational acceleration, a_rad the radiative acceleration, and Γ_rad ≔ a_rad(r)/g(r). From obtaining |$\upsilon (r)$|⁠, |$\upsilon _\infty$| is automatically known. Moreover, |$\dot{M}$|⁠, or alternatively M_*, has to be iteratively adjusted to obtain a consistent solution. We note that the dynamical consistency is reached in the stationary case.

The velocity structure of such a hydrodynamically consistent model may look significantly different from the typically assumed single β or double-β laws. To illustrate that, we take a consistent WN-star model from Sander & Vink (2020), having main parameters T_* = 141 kK, |$\dot{M}_\mathrm{t} = -4.20$|⁠, |$\upsilon _\infty = 1718$| km s⁻¹. A fit, via χ²-minimization, is performed to the resulting |$\upsilon (r)$| with both a single- and double-β description and then we run atmosphere models incorporating the β-type laws with |$\upsilon _\infty$| and all other parameters being identical to the hydrodynamic model. The velocity fields and resulting wind acceleration stratifications are shown in Fig. 15.

From the panel showing |$\upsilon (r)$|⁠, it is immediately evident that the hydrodynamically consistent solution is not well approximated by single- or double-β laws. The wind accelerations show significant differences as well. As expected, the hydrodynamically consistent model correctly balances the outward pushing forces (radiation a_rad and pressure gradient a_press) with the repelling terms (gravity g and inertia a_mech). The models using the single and double β laws do not only miss the local dynamical consistency, but their overall acceleration terms are also significantly different. In particular, the inner acceleration bump caused by the iron M-shell opacities cannot be reproduced. This is likely the main reason for the general failure of β laws to resolve the radius problem of Wolf–Rayet stars.

Focusing on a few diagnostic lines, the spectral imprint of the different velocity laws are shown in Fig. 16. While the single- and double-β law models show only minor differences, the hydrodynamically consistent model spectrum deviates strongly from the other two. Interestingly, the He ii 4686 Å line of the hydrodynamically consistent model in Fig. 16 has a more ‘round’ shape in comparison with the β model spectra. Round-shaped lines – also termed ‘bowler-hat’-shaped lines – have been observed in a sub-sample of WR stars (e.g. Hamann et al. 2006; Hainich et al. 2014) including the Galactic target WR2 and various LMC WN stars. In order to reproduced their shape, high rotational velocities have been adopted in spectral fitting (e.g. Hainich et al. 2014; Shenar, Hamann & Todt 2014). However, these would require a strong co-rotation of the wind, which would require a dedicated mechanism such as powerful magnetic fields (Shenar et al. 2014). A detailed study of WR2 star by Chené et al. (2019) demonstrated that fast rotation and strong magnetic fields to be unlikely, leaving the origin of the line shapes unclear. Given the obtained shape of the He ii 4686 Å line in our hydrodynamic model, it might just be the structure of |$\upsilon (r)$| itself yielding this peculiar emission line shape. However, a much more tailored study for these kinds of objects would be required in order to determine why some, but not all early-type WN stars show this phenomenon.

Beside the line shapes, the dynamically consistent model further has a notably cooler appearance than the β-law spectra, reflected in a T_2/3 difference of about 10 kK. This is quite significant given the former result of T_2/3 being hardly affected when changing β. Since the stellar parameters of the three models are the same, this effect is purely due to the different velocity and density structure of the models.

After verifying that the values of T_2/3 are not significantly affected when reducing the set of considered elements from the hydro models down to those in the typical WR grid models and our above β sequences, we now present a more general overview about the T_2/3-differences between β law and hydrodynamically consistent models. For this, we consider the transformed mass-loss rates |$\dot{M}_\mathrm{t}$| and use the recent result from Sander et al. (2023) that hydrodynamically consistent models predict a trend in the |$\dot{M}_\mathrm{t}$|-T_2/3 plane that is almost independent of chemical composition and only shifts with different clumping assumptions. In Fig. 17, we compare the T_2/3 values from our β-model sequences (see Table 1) with the results from the hydrodynamically consistent models. Given the similarities of the different hydrodynamic sequences obtained in Sander et al. (2023), we only plot the curves for the |$20\, \mathrm{M}_\odot$| sequences at Z_⊙ and X_H = 0.2. The three shown curves differ only in their clumping description, which is depth-dependent with the maximum clumping factor D being denoted in the figure. In contrast, the β models (A1 to C3) all have a constant clumping factor of D = 10.

As expected from our results in Section 3, the β-model sequences show consistently decreasing values of T_2/3 for increasing β. Moreover, the β models cover a parameter regime with cooler T_2/3 temperatures than predicted by hydrodynamically consistent models with winds launched by the hot iron bump. For higher values of |$\dot{M}_\text{t}$| – i.e. more dense winds – which are denoted by darker shaded colours, the hotter (blue coloured labels in Fig. 17) models with the lower β values are located near or even on the curve from the hydrodynamic models with a similar (D = 10) clumping. Especially the lower β C3 models agree well with the hydrodynamic model curve. These results indicate that the traditional WR model grids might partially probe a different parameter space than what we obtain with hydrodynamically consistent models, in particular for less dense WR winds with lower |$\dot{M}_\mathrm{t}$| (lighter shaded coloured labels in the figure). It is currently an open question whether some of the objects that are well analysed with β models can be explained in a dynamically consistent way without invoking additional effects such as an inflated envelope.

6 CONCLUSIONS

In this work, we have investigated how the velocity law impacts the emergent WR-star spectrum by computing model sequences with different β-velocity laws using the PoWR code. We analysed the resulting differences in the diagnostic lines and spectral-type classifications. Detailed modelling revealed that the emergent optical spectra are strongly dependent on the velocity field used in the atmosphere model. This effect is so profound that theoretically the entire range of spectral subtypes could be covered with a single set of stellar parameters and only varying the β parameter. Moreover, the continua are considerably affected by the choice of β with higher values yielding higher flux levels in the optical and near-IR regime. Below the He i ionization edge, the trend is reversed with lower β models providing higher numbers of He i ionizing photons. There are also effects on the |$H\, \small {\rm I}$| ionizing flux, which decreases more moderately with increasing β, and the produced Lyman–Werner photons, which tend to increase with increasing β.

To estimate the actual impact of different β laws in spectral analysis, we performed an exemplary fitting of observations in three different environments with β ≠ 1 models and compared to the β1-model fits. Our findings limit the range of reasonable β values to β ≈ 4 when restricting the analysis to single-β descriptions. While this limits the overall impact, we still find a notable ambiguity in determining stellar parameters when relying on β-type velocity laws, especially concerning the stellar temperature inferred for the hydrostatic layers T_*. The T_2/3-temperature is less affected within the framework of β-type laws, but does shift notably when considering a dynamically consistent velocity law. Within the constrained range of β values, the derived luminosity is only mildly affected by changing the velocity law. Tests with dynamically consistent models beyond Gräfener & Hamann (2005) will be required to check if this also holds for velocity descriptions beyond the β law.

When only relying on optical (or infrared) spectra, the terminal velocity |$\upsilon _\infty$| can be severely underestimated when considering only models with β = 1. This problem does not exist if UV spectroscopy with saturated P-Cygni profiles are available, underlining that despite their emission-line spectra, UV spectroscopy is crucial to fully determine the wind (and feedback) parameters of a WR star. Alternatively, WR-star terminal wind velocities can also be derived from forbidden emission lines (see e.g. Dessart et al. 2000; Ignace et al. 2007). However, these lines are located in the mid-IR and thus IR space telescopes (e.g. JWST MIRI) are required.

By comparing our data from the β-model sequences to results from hydrodynamically consistent atmosphere modelling, we show that β-type velocity laws are insufficient to describe the inner wind regime. This insufficiency also remains when using a double-β approach tailored to the hydrodynamic velocity solution. Given the fact that empirical analyses with models using β ≲ 4 provide good spectral reproductions, the outer regions of WR winds seem to be reasonably well described by a single β-law, despite the shortcomings in the inner wind. From a first test case where we compare a dynamically consistent solution to the ‘best’ β approximations, we conclude that the spectral line shapes are notably affected when trying to approximate with a β or double-β law. Interestingly, the emission lines are more roundish in the hydrodynamically consistent model, hinting at a possible origin for the special ‘round’/‘bowler-hat’ line shapes found for WR2 and various LMC stars.

In a follow-up study, we will investigate how to improve the description of WR-star velocity fields by performing detailed analyses of individual WR stars in different environments with hydrodynamically consistent atmosphere models.

Acknowledgement

This research is based on observations made with the NASA/ESA Hubble Space Telescope obtained from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5–26555. These observations are associated with program 15822 (PI A.A.C. Sander). This work is funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – Project-ID 138713538 – SFB 881 (‘The Milky Way System’, subproject P04). AACS and RRL are supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) in the form of an Emmy Noether Research Group – Project-ID 445674056 (SA4064/1-1, PI Sander) TS acknowledges support from the European Union’s Horizon 2020 under the Marie Skłodowska-Curie grant agreement no 101024605. LP acknowledges support from the KU Leuven C1 grant MAESTRO C16/17/007. KD acknowledges funding from the European Research Council under European Unions Horizon 2020 research programme (grant agreement no. 772225). We also want to thank the referee for their helpful comments and suggestions.

DATA AVAILABILITY

The data underlying this article will be shared on reasonable request to the corresponding author.

References

APPENDIX A: SPECTRAL FITS

The spectral fitting results from Section 3.6 are detailed here, giving a special emphasis on optical diagnostic lines and the C iv λλ 1548-50 P-Cygni profile in the UV. See also Table 4 for the resulting parameters.

A1 WR128

The spectral data from Hamann et al. (1995) are used for this star, which is classified as WN4(h)-w (Hamann et al. 2019). The resulting fit is shown in Fig. A1, where β1-and β4 models are compared to the observed spectrum. Higher β values do not longer provide a reasonable representation of the observed spectrum. There are obvious deficiencies between the two models and the observed spectrum in Fig. A1, but on a level which would e.g. be acceptable in a larger sample analysis with pre-calculated grids. In such a context, both models provide a reasonable representation of the optical spectrum (cf. the bottom two panels).

The corresponding model parameters are shown in Table 4, illustrating significant differences between the models. We see a T_* difference of ΔT_* = 18.4 kK, |$\Delta \upsilon _\infty = 800$| km s⁻¹, and |$\Delta \log \dot{M} = 0.2$| (amounting to about 60 per cent) between the new β1 and the β4 model. This example underlines our previously obtained ambiguity when using different β-velocity laws: a good agreement to the optical spectrum can be found with significantly different model parameters. When considering also the UV (second panel), none of the models fully reproduce the shape of the C iv λλ 1548–50 line. Still, it is clear that the absorption trough is broader than reproduced by the β1 model. The β4 model, which also requires a higher |$\upsilon _\infty$| to reproduce the optical part of the spectrum, provides a better representation of the wind velocity range here. The continua of both models agree fairly well, where there are some differences in the reddening and derived luminosities (see Table 4).

A2 SMC AB2

This star is classified as WN5ha in Hainich et al. (2015) and has an even weaker wind than WR 128. The spectral fitting to the observed spectrum is shown in Fig. A2, displaying the β1 model from Hainich et al. (2015) along with the observation and a new β4 model. Both models fit the optical spectrum (bottom two panels) quite well, with the β4-model overestimating some emission line strengths when compared to the β1 model. The UV spectrum (second panel) is also decently reproduced by both models. The continua (top panel) agree less well however, with significant differences arising in the mid and far IR regime for the β4 model. There is again a significant difference in the model parameters (see Table 4) with especially |$\Delta \log \dot{M} = 1.0$|⁠. It should be noted that the β4 model requires an unlikely high clumping factor of D = 130, making the β1 model more plausible here, especially when considering the deviating continuum of the β4 model. Additionally, the P-Cygni trough of e.g. the C iv λλ 1548–50 doublet is reproduced well by the β1 model. We thus conclude that in this case β = 1 remains the best approximation when focusing on the reproduction of the emergent spectrum. In a previous analysis of SMC AB2, Martins (2011) also used β = 1 and argued that it was best suited to reproduce the shape of the of optical He ii λ 4686 emission line. When considering both WR128 and SMC AB2, the line seems to indeed have some indicative power, but as we will see in the next example, this might be limited to this particular regime of weaker-wind WN stars if they are of intermediate subtype.

A3 M31 WR096 (LGGS J004412.44+412941.7)

Having a completely different spectral appearance – classified as WC6 in Neugent, Massey & Georgy (2012) – than the former two, we pick a WC-type target from M31, applying the naming scheme introduced in Sander et al. (2014). Novel UV spectra of this target were observed with the HST, programme 16170, and are shown in Fig. A3 together with a reasonable fit with both a β1 and a β4 model. Despite some discrepancies between the observed and the model spectra, most diagnostic lines (Table 2) are reproduced relatively well. Similar to the two former observed targets, there is a stellar parameter discrepancy between both models (see Table 4), most notably in T_* where ΔT_* = 17 kK. When only considering the optical spectrum in the bottom two panels, the two model spectra are mostly identical with some difference arising from slightly different mismatches for a few particular lines. The model ambiguity also largely continues in the UV-spectrum (second panel), although one could argue that the C iv λλ 1548–50 absorption trough is slightly better reproduced by the β4 model. Both continua of the two models are in good agreement as well, with only minor differences in reddening and luminosity.

Given that M31 WR096 has never been analysed with atmosphere models before, we also include silicon, phosphorus, sulfur, and neon in our models assuming solar abundances. The latter two (S and Ne) do not have noticeable spectral lines in our available observations, but affect the strength of other lines and improve the overall fit quality. The first two elements (Si, P) have clear wind lines in the UV, to be seen in Fig. A4. However, a further detailed abundance study is beyond the scope of this paper. The model spectra in Fig. A3 both contain these extra elements.