EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 708 (April 2026) A&A, 708 (2026) A206 Full HTML
Open Access
Issue
A&A
Volume 708, April 2026
Article Number A206
Number of page(s) 19
Section Stellar atmospheres
DOI https://doi.org/10.1051/0004-6361/202558762
Published online 09 April 2026

© The Authors 2026

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. This email address is being protected from spambots. You need JavaScript enabled to view it. to support open access publication.

1 Introduction

Mass loss is a key ingredient that can shape the entire evolution of a massive star (M ≳ 8 M⊙) since its birth as an O or early-B star until their eventual core collapse (see, e.g., Langer 2012; Smith 2014; Vink 2022, for recent reviews). While the post-main-sequence stage represents only ∼ 10% of a star’s lifetime, dramatic changes in the stellar properties, such as the surface chemical composition, luminosity, and temperature, can be observed. These changes considerably affect the resulting stellar wind properties, such as the mass-loss rate Ṁ, which in turn is a major ingredient for the subsequent evolution of a star. Therefore, a quantitative understanding of the interplay between the stellar and wind properties is key to better connecting the observed stages of massive post-main-sequence stars and predict their evolutionary fates.

Toward the cool end of the line-driven wind regime reside late-B and A supergiants and hypergiants1. Many of these are luminous blue variables (LBVs) or LBV candidates. In its raw definition, the LBV phenomenon as such is not an evolutionary stage, but an observational umbrella term (Conti 1984) that encompasses several different objects and phenomena, namely eruptive variables (e.g., η Car), P-Cygni supergiants (e.g., P Cyg and HD 316285), and S-Dor variables (e.g., S Dor and MWC 930).

Members of this class of objects, which frequently produce nebulae around themselves, are expected to play a pivotal role in the evolution of massive stars as they are likely the gateway to stages such as Wolf-Rayet stars (WRs) or core-collapse supernovae (CCSNe), especially of type IIn (e.g. Kotak & Vink 2006; Trundle et al. 2008; Nyholm et al. 2020). Moreover, they are expected to be the second main source of metal enrichment of their host galaxies (e.g., Agliozzo et al. 2021). They are also producers of complex carbon molecules because their temperatures and mass-loss rates are sufficiently high and they are highly variable. This is extremely relevant for astrobiology (Arun et al. 2025).

This zoo of stars spans a wide range of effective temperatures and luminosities (Teff=6 … 40 kK; log L=5.0 … 6.5; e.g., Smith et al. 2019). This parameter range encompasses regions in which theoretical models traditionally predict a discontinuous behavior in the wind parameters, which is attributed to the also found discontinuities in a bistability of solutions found originally by Pauldrach & Puls (1990) when modeling P Cyg under the framework of the “modified CAK theory” (Castor et al. 1975; Pauldrach et al. 1986; Friend & Abbott 1986, hereafter mCAK). Lamers et al. (1995) later found discontinuities in the ratio of the terminal velocity v∞ and the effective escape velocity, v esc ,Γe=2GM⋆r(1-Γe),Mathematical equation: $v_{\mathrm{esc}, \Gamma_{\mathrm{e}}}=\sqrt{\frac{2 G M_{\star}}{r}\left(1-\Gamma_{\mathrm{e}}\right)},$(1) for a sample of OBA supergiants.

Incorporating this discontinuity in v∞/vesc,Γe, modeling efforts in a broader parameter regime using a Monte Carlo approach (Vink et al. 1999) yielded two particular sharp increases (jumps) in the behavior of Ṁ as a function of decreasing effective temperature (Teff): a first jump at Teff ∼ 25 kK, the so-called first bistability jump (1-BiSJ), which can be traced back to a recombination from Fe iv to Fe iii in the atmosphere, and a second jump at Teff ∼ 10 kK, the second bistability jump (2-BiSJ), which is attributed to a recombination from Fe iii to Fe ii. These are also included in the mass-loss description from Vink et al. (2001), which is widely employed in stellar evolution. Further, the jumps have been found independently of any assumption on v∞/vesc, Γe (e.g., Müller & Vink 2008; Vink & Sander 2021). Commonly, only the 1-BiSJ is explicitly incorporated, while the switch to the Ṁ treatment for cooler stars approximately incorporates the effect of the 2-BiSJ when a recipe with high mass loss is used, such as that by de Jager et al. (1988).

Recent multiwavelength quantitative spectroscopic studies of BSGs challenged the 1-BiSJ scenario (e.g., Bernini-Peron et al. 2023, 2024; Verhamme et al. 2024; Alkousa et al. 2025; de Burgos et al. 2024). They did not find signs of discontinuous increases and reported a rather constant trend at best that was better captured by the Ṁ prescriptions of Krtička et al. (2024). Moreover, the discontinuity in v∞/vesc,Γe is also contested (Crowther et al. 2006; Markova & Puls 2008; Alkousa et al. 2025). However, there is some evidence for the 1-BiSJ in LBVs based on analyses of AG Car throughout its variability (e.g., Vink & de Koter 2002; Groh et al. 2011; Groh & Vink 2011) and dynamical structure modeling (Grassitelli et al. 2021).

The scientific panorama for the 2-BiSJ remains very unexplored spectroscopically because the UV spectral coverage of late-B and A supergiants is far sparser, and modeling this parameter regime is challenging. Therefore, hardly any investigations were dedicated specifically to the wind behavior in this regime. In yellow supergiants and hypergiants (YSGs and YHGs), the 2-BiSJ is thought to trigger dust shells around some of these stars (e.g., Koumpia et al. 2020). Using a grid of CMFGEN models with global wind consistency, Petrov et al. (2016) found that the 2-BiSJ is associated with the change in ionization of Fe iii to Fe ii for models with a high luminosity-to-mass ratio (L/M). Sabhahit et al. (2026), using a grid of fully hydrodynamically consistent models, also found an increase in Ṁ that was driven by Fe ii for models with high L/M. However, their investigation was motivated by the enhanced Ṁ of very massive stars, and thus, their grid in this region is sparse.

We use hydrodynamically consistent atmosphere models to show that the behaviour of mass loss and other wind properties of stars closer to the Eddington limit in the Hertzsprung gap region (usually occupied by hypergiants and LBVs) is indeed complex. This is strongly tied to the interplay of changes in the ionization of Fe, H, and He, and to the action of the star’s atmospheric turbulent pressure.

In Sect. 2, we briefly describe our tools, approach, and assumptions. In Sect. 3, we briefly illustrate the resulting optical spectra and their resemblance to those of known LBVs and BHGs before we present the trends in mass-loss rate and wind structure we obtain. We then discuss the implications of our results in Sect. 4 before we draw our conclusions in Sect. 5.

2 Methods

We used the stellar atmosphere code PoWRhd (Sander et al. 2018) to produce a fine-spaced sequence of hydrodynamically consistent models of massive stars with moderate mass and luminosity that were still close to the Eddington limit. We selected the option to fix all stellar parameters to PoWRhd when we computed the velocity field (including v∞) and Ṁ.

We fixed the mass (20.6 M⊙), luminosity (log L/L⊙=5.58; representative values of that lead to an L/M typical for hypergiants and LBVs), and a mildly enriched He abundance (with a He mass fraction XHe=0.33) as well as the same CNO composition as in Bernini-Peron et al. (2025). Motivated by time-dependent hydrodynamical simulations of Driessen et al. (2019), who found that stars below Teff<20 kK have lower clumping, we fixed the volume filling factor fvol=0.33 at all depth points.

We adopted a depth-independent turbulence velocity (vturb= 20 km s−1), as is typically employed in atmospheres of LBVs and even some BSGs and BHGs (e.g., Groh et al. 2009). However, we note that this value might be an underestimation, especially for the hotter models of our sequence, as recent multidimensional simulations of atmospheres for stars with high Γe(Rcrit) (Moens et al. 2025) predict much higher values of vturb ∼ 150 km s−1. As we aim for a homogeneous sequence of models isolating the effect of a T* variation, we kept v∞ (and other parameters) fixed.

In the atomic physics, we considered the elements H, He, C, N, O, Si, S, Al, P, and Mg, and the iron-group elements (G), which, albeit largely dominated by Fe, conflate levels and transitions from Sc to Ni. In Appendix A, we disclose the considered levels and transitions for each element.

In our main model sequence, we only varied the temperature at the inner boundary T*=Teff(τ=50), where τ is the Rosseland continuum optical depth, from 13 kK to 38 kK. The photospheric temperature, that is, Teff(τ=2/3) (hereafter simply T2/3) is an output quantity that depends on the resulting Ṁ. As our model sequence has high Γe values, the models’ atmospheres are expanded, such that the difference between T* and T2/3 amounts to ∼ 7 kK. In Fig. 1 we show the location of our model sequence, indicating T* and T2/3, in the Hertzsprung–Russell diagram (HRD) together with a set of known objects in the hypergiant and LBV domain. All comparison data were taken from the literature, although the luminosities were updated when newer distance information from Gaia DR3 was available (see Appendix F).

A more direct illustration of the difference between T* and T2/3 at different optical depth is shown in Fig. 2. Despite some scatter on small scales, the atmospheres clearly become more extended with cooler T*, and the difference between T* and Teff(τ<50) increases. For the coolest models, the effective temperature of the critical point becomes lower than T2/3=Teff(τ= 2/3), meaning that the wind stars to launch beneath the layers forming the continuum. This is a signature of an optically thick wind and is known in particular from models describing classical WRs (e.g., Sander et al. 2020, 2023; Lefever et al. 2023, 2026).

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

Hertzsprung-Russell diagram depicting our model sequence against LBVs (cyan squares), B[e] SG (blue triangles), OBA hypergiants or LBV candidates (green circles), YSGs (small yellow pentagons), and YHGs (large gold pentagons) from the literature (see Appendix F). The blue lines connecting the LBV data indicate the same star in different epochs. Likewise, the dark orange lines connect YHGs in different epochs. The luminosities are updated to distances mostly based on Gaia DR3 (Bailer-Jones et al. 2021; Gaia Collaboration 2023). The black line indicates the luminosity of our sequence. The small diamonds below indicate T*, and the large open diamonds indicate T2/3. The thin dotted lines indicate the Galactic evolutionary tracks by Ekström et al. (2012). The thick cyan line indicates the LBV instability strip (Groh et al. 2011), and the thick violet line indicates the HD limit (Humphreys & Davidson 1979).

3 Results

3.1 Spectral appearance

Our PoWRhd models are not specifically designed to fit any particular star and have a fixed luminosity. However, as illustrated in Fig. 3, we found that the normalized synthetic spectra of our models generally resemble spectra of existing hypergiants and LBVs in their temperature domain. This indicates that many hypergiants, LBVs, and even some supergiants that spectroscopically appear to be very different from each other might belong to this L/M regime. Moreover, this is an important sanity check for the robustness of the atmosphere structures and indicates that our sequence reflects their physics in the hypergiant regime.

Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

Effective temperature at different optical depths τ, normalized by T*=Teff(τ=50) for each point of the sequence. The pink curve indicates T2/3, and the blue curve indicates Teff at the critical point.

3.2 The behavior of the mass-loss rate

Figs. 4 and 5 show a complex pattern for the derived mass-loss rate as a function of T*. The transformed radius, Rt=R⋆(v∞M˙f cl ·10-4M⊙ yr -12500 km s-1)2/3,Mathematical equation: $R_{\mathrm{t}}=R_{\star}\left(\frac{v_{\infty}}{\dot{M} \sqrt{f_{\mathrm{cl}}}} \cdot \frac{10^{-4} \mathrm{M}_{\odot} \mathrm{yr}^{-1}}{2500 \mathrm{~km} \mathrm{~s}^{-1}}\right)^{2 / 3},$(2) describes the ratio of the emission measure of the wind and the surface area at R*, and in the case of our T* sequence, it separates the geometrical effect from the wind behavior. Two main types of wind solutions can be identified: dense-wind solutions with Rt> 200 R⊙, and rarefied-wind solutions2 with smaller Rt (see Fig. 4). Significant jumps to higher Ṁ appear around the temperatures associated with the traditional bistability jumps, but these jumps are not the only abrupt transitions. Instead, Fig. 5 shows a more complex pattern with jumps to higher and also drops to lower Ṁ. In the following, we refer to these as valleys in Ṁ(T*).

In the background of Fig. 5, we show a stacked area plot indicating the contributions to the total outward pushing forces Γwind:=Γrad+ΓpressMathematical equation: $\Gamma_{\text {wind }}:=\Gamma_{\text {rad }}+\Gamma_{\text {press }}$(3) at the critical point Rcrit, with Γrad summarizing the radiative line and continuum accelerations, while Γpress includes the acceleration due to gas and turbulent pressure gradients, all normalized to gravity g. The temperatures associated with the two classical BiSJs are clearly associated with the change in the ionization of Fe. However, the other end of the valleys, that is, the decline of Ṁ to the low values before the jumps, is less abrupt and seems to be associated with the gradual increase in the contribution of other elements and the pressure term. In turn, this seems to be indirectly caused by a slow transition of the dominance of different Fe ions. A similar behavior was recently mapped for stars close to the Eddington limit (Lefever et al. 2025; Sabhahit et al. 2026), but for a much higher luminosity regime.

A further inspection of the relative contributions to Γwind in Fig. 5 shows that radiatively driven processes dominate the wind acceleration. However, the sheer proximity to the Eddington limit of the model sequence (Γe ∼ 0.4) leads to the effect that in all models, less than half of the acceleration at the critical point is provided by line opacities. In addition to a large contribution from free-electron (i.e., Thomson) scattering, the continuum opacity of hydrogen, and the acceleration due to gas and turbulent pressure gradients, need to be taken into account. In the following, we discuss the different features in the mass-loss trend. We also consider for this the work ratios Qk=∫R⋆∞a rad ,kdr/∫R⋆∞(vdvdr+GM⋆r2)dr,Mathematical equation: $Q_{k}=\int_{R_{\star}}^{\infty} a_{\mathrm{rad}, k} \mathrm{~d} r / \int_{R_{\star}}^{\infty}\left(v \frac{\mathrm{~d} v}{\mathrm{~d} r}+\frac{G M_{\star}}{r^{2}}\right) \mathrm{d} r,$(4) where the arad,k is the radiative acceleration due to the respective contribution k (e.g., electron scattering, gas pressure, Fe, H, and CNO). In Fig. 6, we plot Γe(rcrit) and Qe, which shows that the valleys in Ṁ match regions of low Qe, that is, a low global contribution of the electron scattering and thus a high amount of line driving.

The hot dense-wind plateau. The hottest models in the sequence have moderately high log Ṁ ∼−6 and probe the temperature range of OHGs. The winds in this region are largely driven by Γe – and Fe iv. Pressure and other elements contribute up to ∼ 25% to launch the wind.

The obtained terminal velocities match the v∞ determined by Prinja et al. (1990) for GS Mus (O9.7Ia+, Carretero-Castrillo et al. 2023), one of the few Galactic OHGs that has a UV spectrum. The T2/3 and log g(τ=2/3) also match the results obtained by Holgado et al. (2020), which results in a very similar spectral appearance. We also find a similar spectral appearance of HD 151804 (formally an O8Iaf) in the optical. In the UV, we find an apparently slow wind. It is possible that elements like Ar, which are shown to affect v∞ in earlier O stars and are not included in our sequence, might start to play a role. While we did not perform a dedicated spectral analysis, we scaled the flux of our models to the observed SED (see Appendix D) to determine the luminosities of GS Mus and HD 151804 (log L/L⊙ ∼ 5.2 and log L/L⊙ ∼ 6.0) and placed them in the HRD (see Fig. 1).

The hot rarefied-wind valley. As Fe iv slowly loses dominance, the available opacity is insufficient to maintain a deeper wind-launching associated with a higher mass-loss rate. With decreasing τcrit, CNO and other metal lines contribute quite significantly to the wind launching, as does the pressure term. In this regime, the winds reach ∼ 1000 km s−1, which resembles winds of BSGs (∼ B1 Ia) on the verge of the predicted classical 1-BiSJ (Vink et al. 2001). Spectroscopically, the models in this range and the B 1 Ia supergiant ρ Leo match as well.

The classical first jump and the warm dense-wind plateau. At T* ∼ 27 kK(T2/3 ∼ 22 kK), the sharp increase in Ṁ is directly associated with the increase in Fe iii that contributes to the wind driving. This is the 1-BiSJ already seen in the Monte Carlo models by Vink et al. (1999); Vink & Sander (2021), albeit with shifts in the temperature. In Appendix C we show that we found a bistability of solutions in this regime transition, but we only selected one of them for clarity of the overall plots, such as Fig. 5. Toward cooler temperatures, the regime of high M persists until T* ∼ 23 kK(T2/3 ∼ 15 kK). In contrast to the results by Vink et al. (1999), we found no gradual decrease in Ṁ after the jump, but a plateau. As for the hot plateau before, the contributions other than from Fe and Γe only amount to ∼ 20%.

The wind solution we obtained for this region is very similar to the solution obtained for the PoWRhd model of ζ1 Sco (B1.5Ia+ Bernini-Peron et al. 2025). The corresponding spectra for models in this regime resemble the spectrum of BP Cru (B1Ia+/HMXB), a hypergiant with similar properties to ζ1 Sco (Clark et al. 2012). We also note that it resembles the spectrum of HR Car in 2017, suggesting that the star is currently in this densewind regime. In this temperature range, we also find LBVs such as P Cyg (B1Ia+/LBV), whose spectra are much more extreme, which suggests either a considerably different clumping structure or a higher L/M ratio. In contrast, we also see stars like HD 190603 (B1.5Ia+) and ζ1 Sco, which have less prominent P Cygni features (see Sect. 6.3 in Bernini-Peron et al. (2025) for a discussion on the effects of clumping in the spectra of early-type BHGs).

The cool rarefied-wind valley. At T* ∼ 22 kK(T2/3 ∼ 15 kK), the mass-loss rates drop by more than an order of magnitude. This behavior is the opposite of what was predicted for the classical BiSJ and is not associated with any major changes in the ions (in particular Fe) driving the wind. Nonetheless, the wind solution changes and τcrit decreases again because the Fe iii opacities are no longer sufficient to launch the wind at similar radii as obtained for the solutions of the warm dense-wind plateau. Fe ii is not sufficiently populated to take over the wind launching, however. Instead, the pressure term supports (quasi-)hydrostatic equilibrium, until an outflow is launched farther out, at the expense of a considerably reduced Ṁ. Unlike the hot rarefiedwind valley, non-Fe metal lines do not significantly contribute to the wind launching (≲ 10%), while the contribution of the pressure term is even stronger.

The synthetic model spectra in this regime resemble the appearance of several late-type BHGs (B2...B9), for instance, HD 80077, HD 183143, and HD 168625. Notably, many of the known BHGs in this temperature regime (T2/3 ∼ 14 … 10 kK) also display the obtained weak hypergiant characteristics, that is, Hβ, Hγ, Hδ, with possibly He and metal lines in P Cygni or emission. The hypergiant status is questioned for some of those stars because Hβ is variable and often fully in absorption (see Negueruela et al. 2024). Our finding of the cool rarefied-wind valley indicates that these weaker characteristics might be a natural outcome of the opacity situation in the atmosphere and do not reflect a lower L/M. This would mean that earlier-type BHGs such as HR Car, BP Cru, or ζ1 Sco are not necessarily closer to the classical Eddington limit, but rather have more favorable wind conditions.

The classical second jump and the cold dense-wind wall. Toward even lower temperatures, Ṁ increases significantly at T* ∼ 17 kK (T2/3=10 kK), which is associated with the change from Fe iii to Fe ii in the wind driving. The higher wind density causes the model spectra in this regime to resemble AHGs and LBVs that have P Cygni and emission features, such as HD 160529 and the LBV MWC 930, crowded with Fe ii profiles in the case of the coolest models. Approximately similar in T2/3, this essentially represents the 2-BiSJ predicted by Vink et al. (1999) and Petrov et al. (2016). The jump also occurs in the LBV models by Vink & de Koter (2002), albeit at temperatures more similar to our T* value. The 1-BiSJ also has a bistability of solutions at the onset of the cold dense-wind wall. We illustrate this further in Appendix C.

Notably, the 2-BiSJ in Petrov et al. (2016) only appeared in models closer to the classical Eddington limit (Γe=0.39, log L/L⊙=5.5, M=20 M⊙, and Y=0.25) when they transitioned from T2/3=10.0 kK to T2/3=8.8 kK. The jump in their analysis appears in a very similar parameter space and between very similar Ṁ values (−7.04 to −5.70) as those that we obtain (−6.7 to −5.8). Interestingly, this similarity in outcome does not coincide with the same driving explanation. Petrov et al. (2016) reported that Fe ii causes ∼ 70% of the total work ratio, while in our case, Fe ii contributes only up to ∼ 10% (first panel of Fig. B.1) and most is contributed by electron scattering and pressure. Moreover, when we consider the total work ratio Qwind for our models in the cold dense-wind wall, all the Fe ions contribute similarly. In addition, the total Fe contribution to Qwind is equivalent to that of H and much smaller than that of electron scattering or pressure.

In the models from Petrov et al. (2016), the 2-BiSJ marks the last point of their model sequence. In their models, Fe ii is not fully recombined at T2/3=8.8 kK, raising the question of the further trend of Ṁ versus temperature. Our model sequence extends slightly farther, and Ṁ continues to increase together with the contribution of Fe ii. During the S Doradus cycle, these high values of Ṁ might explain nebulae around many LBVs and hypergiants, including those that currently do not present themselves as bona fide LBVs. Ahigh Ṁ past the 2-BiSJ like this might also explain the signatures of type IIn SNe (e.g., Kotak & Vink 2006; Trundle et al. 2008). Within the capability of our current simulations, we did not reach a stage in which the trend in Ṁ turns over into another plateau or valley. However, we found an indication of a possible stagnation as the increase in Ṁ becomes smaller for the coolest models in the series, which was also suggested by the Ṁ behavior found by Vink & de Koter (2002).

Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Nontailored comparison between the output spectra of some of the models (red lines) and observed spectra (black lines) of OB supergiant/hypergiants and LBVs.
Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Transformed radius (Rt) vs. T*. Models with Rt<200 R⊙ are in the dense-wind regime, i.e., with higher Ṁ and optical spectra crowded with emission and P Cygni features (shaded region). The colored texts in the lower part of the diagram name the wind-solution regions identified in the T* sequence. The 1-BiSJ and the 2-BiSJ are indicated by the dotted and dashed lines in the transitions from a valley to a plateau (i.e., increase in Ṁ) toward lower T*.
Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Contributions of different processes, elements, and ions to the wind launching at the critical point (stacked area plot) and corresponding mass-loss rates (red line) for different inner boundary temperatures (T*). From bottom to top, the light beige area indicates the electron-scattering acceleration, dark beige shows the pressure, and smoke white indicates H i continuum driving. The teal tones indicate the contribution of CNO and other metals, and the shades of blue indicate the contributions of different Fe ions.
Thumbnail: Fig. 6 Refer to the following caption and surrounding text. Fig. 6

Comparison between the trend of mass-loss rates (red curve) and the classical Eddington parameter at the critical point Γe(Rcrit) (thick dotted curve) and the total contribution of the electron scattering to the work ratio QewindMathematical equation: $Q_{\mathrm{e}}^{\text {wind }}$ (shaded area).

3.3 Velocity structure

The velocity structures in the T* sequence reflect the different dense and rarefied solutions. The upper panel of Fig. 7 shows that the trends in v(r) trace the five different solution regimes presented above, namely models in the valleys see an increase in v∞ compared to models in the plateaus. The differences become even clearer in the shapes of v(r), which show major structural differences for the different solutions (lower panel of Fig. 7). None of the solutions is sufficiently described by a β-velocity law, especially in the wind-launching region. However, when we fit the solutions of the hot dense-wind plateau (OHG regime, blue curves in Fig. 7) with a single β-law, we obtained solutions resembling the outer wind part with β-values between 1.8 and 2.0, explaining why higher β-values are often favored in traditional analyses of B super- and hypergiants.

The hot dense-wind plateau has v∞ ∼ 700 km s−1 and a relatively sharp velocity gradient. In the hot rarefied-wind valley (BSG-like regime, green curves), this increases to ∼ 1000 km s−1 and an even more extreme gradient at the critical point. When we fit this regime with β-laws, we obtained slightly higher values of β ∼ 2.3. From the 1-BiSJ into the warm dense-wind plateau (gray curves), v∞ sharply drops to ∼ 300 km s−1. In parallel, the shapes noticeably switch from a shallower velocity increase and a reduction of the maximum velocity gradient. The corresponding best β-values increase to ∼ 4 in this regime, although the representation of the shape in the outer wind clearly deviates from a β-law.

Within the cool, rarefied-wind valley (yellow curves), the solutions switch back to a steeper gradient (like in the hot valley previously), but commence farther out in the wind (i.e., larger Rcrit). The v∞, however, does not increase much relative to the warm plateau, while the best β increases further to values between 8 and 20. After the 2-BiSJ, in the cold dense-wind wall, v∞ is rather constant despite a rather sharp increase in the mass-loss rates. This is a consequence of the support by gas and turbulent pressure in the outer wind (see Fig. B.3). As the pressure term drops only with 1/r, it eventually becomes the dominant term in the outer wind and causes up to ∼ 50% of the obtained v∞ value. This regime can no longer be meaningfully described by a β-law. Like in the warm plateau, the differences in v(r) are strong, with the coolest models at ∼ 10 R* showing a velocity plateau due to the decreasing effect of Fe ii (lowest panel in Fig. B.2).

The different wind regimes and velocity structures also leave their imprint on the resulting synthetic spectra. In Fig. 8 we plot the normalized spectra around the first three Balmer lines as well as He i 5875 Å. In addition, we mark the locations of v∞ leftward and rightward of the central line wavelength. In the hottest regime as well as in the rarefied regimes, the lines do not reflect v∞, even when parts are in emission. However, for the warm plateau below the 1-BiSJ, at least H α traces v∞ when the full emission line width is considered. Below the 2-BiSJ, the PCygni troughs of all Balmer lines trace the terminal velocity well, while for He i 5875 Å this is only true for the coolest models, in which He is entirely recombined in the outer atmosphere and can produce the blue absorption.

Thumbnail: Fig. 7 Refer to the following caption and surrounding text. Fig. 7

Terminal velocities (upper panel) and velocity structures (lower panel) of the models along the T* sequence. Each wind-solution domain is color-coded following the scheme of Fig. 4.
Thumbnail: Fig. 8 Refer to the following caption and surrounding text. Fig. 8

Hydrogen and helium profiles in velocity space compared to the terminal velocity (black dots). The color of each spectrum marks the wind-solution regime to which the model belongs, following the same color scheme as Fig. 7. The temperatures rightward in the upper panels indicate T2/3, and rightward in the lower panel, they indicate T*.

4 Discussion

4.1 Terminal and escape velocity ratio

From conservation of energy, we expect that the terminal velocity of radiation-driven winds directly depends on vesc,Γe at the stellar surface. For line-driven winds within the mCAK formalism, v∞/v esc ,Γe∼[α/(1-α)]1/2,Mathematical equation: $v_{\infty} / v_{\mathrm{esc}, \Gamma_{\mathrm{e}}} \sim[\alpha /(1-\alpha)]^{1 / 2},$(5) where α is the parameter that condenses the blend of optically thin and thick sets of lines (see, e.g. Lamers et al. 1995). This means that changes in the total line opacities as a result of temperature changes, for example, will affect how v∞ and vesc,Γe scale. Using literature values and spectral type calibrations as well as Γe estimates based on evolutionary tracks, Lamers et al. (1995) found two jumps in the v∞/vesc,Γe ratio in the OBA supergiant regime. These were later interpreted as a consequence of changes in the dominant Fe ion (Vink et al. 1999). More recent studies contested the sharpness of these jumps (e.g., Crowther et al. 2006; Markova & Puls 2008).

With our PoWRhd model sequence, we tested the v∞/vesc,Γe ratio and compare it to that of Lamers et al. (1995) in Fig.9. For stars with extended atmospheres, the choice of r in Eq. (1) (where Γe has to be evaluated as well) is not trivial, and we thus plot results for the inner boundary radius (corresponding to τ=50), R2/3, and the critical radius Rcrit. As the latter essentially reflects the wind onset, we considered this as the most coherent choice. For models on the hot side of the classical 1-BiSJ, OBHGs have a lower v∞/vesc,Γe than OBSGs. This is mainly the result of the lower v∞ values of the hypergiants, which are lower by at least a factor of two than those of comparable supergiants (i.e., similar temperatures). Hence, although Γe is higher for hypergiants, the change to vesc,Γe resulting from the 1-ΓeMathematical equation: $\sqrt{1-\Gamma_{\mathrm{e}}}$ dependence in Eq. (1) is much smaller.

In the hot rarefied valley, the terminal velocities become significantly higher because the lower mass-loss rates are lower. In this regime, some models even exceed the ratio found by Lamers et al. (1995) for BSGs (but in agreement with some results reported in Markova & Puls 2008), which is a combined effect of the higher v∞ and a lower vesc,Γe caused by a larger Rcrit. In the 1-BiSJ regime, where Fe iv recombines to Fe iii, the v∞/vesc,Γe ratio drops sharply to values only slightly above unity, which agrees well with the results of Lamers et al. (1995). While the v∞ values are very similar to those in the hot valley, the drop in the v∞/vesc,Γe ratio is mainly caused by the deeper wind-launching, resulting in lower Rcrit values (i.e., higher vesc,Γe).

In the cool valley, before the 2-BiSJ, the ratios increase again. While a glimpse of this might be seen in the results from Lamers et al. (1995), as well as in some of the recent METUJE calculations by Krtička et al. (2025), and possibly also in the empirical determinations by Markova & Puls (2008), we obtained a much more pronounced effect, again for the same reasons as in the first valley. Crossing the 2-BiSJ, where the Fe iii → Fe ii switch occurs, then coincides with another drop in v∞/vesc,Γe, which is also reflected in the data from Lamers et al. (1995), albeit with an offset. However, since the turbulent pressure term affects v∞ considerably in this regime, as discussed above, we might overestimate the ratio by up to a factor of two, which would bring our data in line with those of Lamers et al. (1995), despite their rough temperature assumptions from spectral types for the B and A stars.

Our v∞/vesc,Γe in the A-star regime are also larger than the values obtained from the so-called δ-slow solutions by Curé et al. (2011) and empirical values listed therein. Again, we attribute this to the effect of the pressure term in the outer wind, which consists of a gas pressure and a turbulent pressure component. In the current framework, both components are treated together as an isotropic pressure. The turbulent pressure is further described with a constant vturb. While the value is clearly motivated by observations (e.g., Searle et al. 2008), simulations of the multidimensional time-dependent hot star atmosphere in the WR- and O-type star domain from Moens et al. (2022) and Debnath et al. (2024) suggested that radiatively driven turbulence is no longer isotropic in the supersonic wind. While a more detailed treatment of radiatively driven turbulence requires a major update of the 1D hydrodynamic framework as well as calibration with multidimensional models, it is clear that solutions with a leading 1/r-acceleration in the outer wind have to be critically reviewed. For the coolest models, Fig. B.3 clearly shows that the turbulent velocity exceeds the sound speed, meaning that the outermost wind solution is dominated by vturb, similar to the models in Sabhahit et al. (2026). These results, as well as the ongoing discussions about turbulent pressure as a driver for cool massive star outflows (e.g., Stothers 2003; Klochkova 2019; Kee et al. 2021), underline that turbulent pressure is an important ingredient at the cool end of the massive star regime, requiring the development of a more coherent framework addressing the interplay between radiative acceleration and radiatively driven turbulence in launching mass outflows.

Another caveat when discussing terminal velocities in hydrodynamically consistent modeling is the choice of wind clumping, which can affect the calculated v∞ by enhancing opacities. We fixed fcl=3 throughout the whole sequence to a value motivated for the BSG regime by the line-deshadowing instability simulations from Driessen et al. (2019). As discussed in their work and as also empirically supported, the clumping in the supersonic region in early- B and O supergiants (i.e., on the hot side of the 1-BiSJ) is expected to be much larger (fcl ∼ 20) than in later-type BSGs. Therefore, it is plausible that v∞ is systematically underestimated for the hottest models of our sequence (i.e., in the hot plateau). This, however, partially contradicts our finding of a very good match between the UV P Cygni profiles of the OHG GS Mus, but agrees with the mismatch for HD 151804, which has a much higher v∞.

Thumbnail: Fig. 9 Refer to the following caption and surrounding text. Fig. 9

Behavior of the v∞/vesc,Γe relation. The gray crosses indicate the data of Lamers et al. (1995) for OBA supergiants. The blue curve with diamonds indicates the values for vesc,Γe applied to the critical point. The pink curve indicates the ratio at τ=2/3, and the thin light blue curve indicates the ratio at τ=50, i.e., the inner boundary.

4.2 The v∞−Teff relation

From v∞ ∝ vesc,Γ and the Stephan-Boltzmann equation, we derive the relation v∞=C(κ[X He ,Z,ρ,f cl ])L1/4T eff 1-ΓeΓe,Mathematical equation: $v_{\infty}=C\left(\kappa\left[X_{\mathrm{He}}, Z, \rho, f_{\mathrm{cl}}\right]\right) \quad L^{1 / 4} T_{\mathrm{eff}} \sqrt{\frac{1-\Gamma_{\mathrm{e}}}{\Gamma_{\mathrm{e}}}},$(6) where C is a constant that depends on the opacity κ, which depends on the chemical composition (XHe, Z), density ρ, and clumping factor fcl. This relation explains the empirical scaling of v∞ with Teff. In Fig. 10 we plot our resulting values using Teff=T2/3. The different dense and rarefied solutions leave clear imprints as they show different slopes in the v∞−Teff relation with sharp breaks precisely where the wind becomes denser (i.e., in the dense-wind plateaus). While there is more scatter for the hotter models, the rarefied solutions generally match the empirically obtained slope from Hawcroft et al. (2024) (but shifted in Teff for the cooler valley), while in the dense-wind regimes, the relation becomes much flatter.

In Fig. 10, we also plot recent other temperature sequences of hydrodynamically consistent models (with otherwise homogeneous parameters) in the OB and WR regime. Other studies with Γe ∼ 0.3 tend to agree well with the slope of the relation reported by Hawcroft et al. (2024). Likewise, when their winds become thicker, the slopes become flatter. While the slopes in the dense-wind regime disagree with those of the WR models from Sander et al. (2023), from which we show the 12.9 M⊙ sequence with Γe ∼ 0.32, the thin-wind regime also matches the (shifted) empirical relation by Hawcroft et al. (2024), indicating that it might be generally valid throughout the whole line-driven wind regime down to ∼ 10 kK, independent of Γe, as long as the winds do not become dense. Even the near-Eddington models from Lefever et al. (2025) agree mostly with the Hawcroft et al. (2024) slope, similar to our rarefied-wind solutions with a departure toward the trend from the WR sequence for T2/3>28 kK.

Below ∼ 10 kK, that is, past the 2-BiSJ, the stars approach the edge of the line-driven wind regime, and their mass loss might be generated and sustained by other sources, such as thermal and turbulent pressures and pulsations (e.g, Glatzel & Kraus 2024; Glatzel & Kiriakidis 1993; Fadeyev 2010). As discussed above, most of the outer wind in our models in this temperature regime is sustained by gas and turbulent pressure. Turbulent pressure also plays a role in launching, although radiative acceleration remains the major contributor there, as shown in Figs. B.1 and B.2. For the coolest models, the dependence on the turbulence is strongest as a viable wind (i.e., Γwind>1) is not expected to occur without this term, and different values of vturb will change v∞. This is also shown in Fig. 10, where the model sequence of Sabhahit et al. (2026) with similar T2/3 has higher v∞ values precisely due to the much higher vturb=70 km s−1.

Thumbnail: Fig. 10 Refer to the following caption and surrounding text. Fig. 10

Relation of v∞−Teff for different studies that made hydrodynamically consistent sequences of models in which T* varied. The black lines indicate the empirical relation by Hawcroft et al. (2024) (the solid line represents the actual relation, and the dashed line indicates the relation shifted arbitrarily). The pink and violet squares indicate the models by Sander et al. (2023) and Lefever et al. (2025), respectively. The blue circle indicates the models of Sabhahit et al. (2026) with Γe ∼ 0.35, and the green circle indicates the models of Björklund et al. (2021) with similar Γe. The red and brown diamonds indicate our sequence of models and show the rarefied and dense solutions, respectively. The lower panel is a zoom-in of the region within the dotted rectangle in the upper panel.

4.3 Interpretation of the mass-loss trend and comparison with existing recipes

As shown in Sect. 3, the Ṁ trend we obtained with plateaus and valleys is benchmarked by the resemblance of the corresponding synthetic spectra to several hypergiants and LBVs in the OBA star domain, and it is thus expected to capture the atmospheric physics in this high- Γe regime. The jumps in Ṁ (1-BiSJ and 2-BiSJ) are clearly caused by the recombination of Fe, as predicted by Vink et al. (1999). In contrast, no clear jumps are observed in BSGs, despite changes in Fe ionization at the base of the wind (e.g., Bernini-Peron et al. 2023, 2024; de Burgos et al. 2024; Alkousa et al. 2025). This clearly indicates a difference that to first order originates in the proximity to the Eddington limit as the BSGs are in a lower Γe regime. Together with the recent calculations in Sabhahit et al. (2026) probing different Γe-sequences (although with fixed L), the complexity of the derived trend and the bistability of the solutions obtained near the jumps suggest a general bimodality of dense and rarefiedwind solutions. For sufficient conditions, a dense-wind solution with high mass-loss rates can be reached where Ṁ even seems to increase with decreasing temperature. However, when these conditions are not met, the solutions switch to the rarefied branch, where mass-loss rates tend to decline with decreasing temperature. When we adopt this generalized paradigm, the spectra of our findings reveal that hypergiants and some late-type LBVs reside in a transition domain in which the atmospheric conditions can switch between rarefied and dense solutions, depending on the temperature.

While the two solutions do not correspond to the transition to optically thick winds (i.e., τcrit>1), the transition is still triggered by the availability or absence of sufficient opacity (plus turbulent pressure), thus marking a striking similarity with the regime differences found for stripped He-burning stars in Sander & Vink (2020) and Sander et al. (2023). Comparing our findings with the different sequences from Sabhahit et al. (2026) and Björklund et al. (2023), we suggest that depending on Γe, there might be no switch between dense and rarefied regimes because for high Γe (and sufficient metallicity, see below), the whole sequence can be in the dense-wind regime where continuum and geometrical effects will dominate. In contrast, for low Γe, all solutions fall onto the rarefied branch, as likely occurred in the models from Björklund et al. (2023) and Krtička et al. (2024, 2025). However, to confirm this hypothesis, a systematic exploration of the precise conditions is required to isolate the contributions of supplementary effects, such as clumping and turbulent pressure.

Comparing our results with existing mass-loss recipes, we found that the recipe from Krtička et al. (2024) applied to our model properties (considering T2/3) matches Ṁ in our models in the rarefied-wind regime very well. Their recipe successfully reproduces the typical mass-loss rates determined for cooler BSGs (e.g., Bernini-Peron et al. 2023, 2024; Verhamme et al. 2024). The recipe by Krtička et al. also predicts a mild increase in the warm plateau region, albeit with a temperature offset. The recipe by Björklund et al. (2023) (applied to T2/3) in the hot plateau and relatively cool valley matches very well, but underestimates Ṁ for cooler temperatures. When we compared our results to those by Vink et al. (2001, applied to T2/3), we found a decent match to our sequence up to the cool valley, with a systematic offset of 0.8 dex. Interestingly, when we applied Vink et al. (2001) to our sequence, considering T* and shifting the results down by this offset, we obtained an excellent match with our results for the hot plateau, hot valley, and warm plateau (including the 1-BiSJ). All the recipes fail catastrophically to explain the behavior past the 2-BiSJ. At the coolest temperatures, the unscaled Vink et al. (2001) rates again match the Ṁ values we derived.

Thumbnail: Fig. 11 Refer to the following caption and surrounding text. Fig. 11

Comparison between the output mass-loss rates of the T* sequence (red curve) and different recipes applied to the stellar properties of each model. The thick purple and fuchsia lines indicate the mass-loss rates by Vink et al. (2001), applied to T* and T2/3, respectively. The thin lines of the same colors indicate the Vink rates with an offset of 0.8 dex. The yellow curve indicates the rates of Björklund et al. (2023), and the teal line indicates the rates of Krtička et al. (2024), both applied to T2/3. The green line indicates the rates computed using the tool LIME (https://lime.ster.kuleuven.be/, Sundqvist et al. 2025).
Thumbnail: Fig. 12 Refer to the following caption and surrounding text. Fig. 12

Contribution to the radiative acceleration from different sources, elements, and ions (stacked colored areas) and mass-loss rates (red curve) according to metallicity for a fixed temperature of T*= 16.1 kK.

4.4 The effect of metallicty on the cool models

One of the most striking results of our study is the increase in Ṁ toward the coolest temperatures. As we discussed above, this is considerably different from what has been obtained in recent modeling efforts in the B supergiant regime. To investigate in particular the increase due to the ionization change from Fe iii to Fe ii, we also created auxiliary shorter sequences for a fixed T*=16.1 kK in this regime, where we changed the He composition and metallicity. In Fig. 12, we present the results, where we varied the metallicity of the model, and hence, also the Fe content. It is noticeable that between 0.7 and 0.6 Z⊙, Ṁ drops by 1.5 dex, causing the models to switch back to a rarefied solution.

This additional sequence confirms that the line opacities are important for the switch between the two regimes, boosted in particular by the recombination of Fe iii to Fe ii. The Fe ii opacities cause the mass-loss rates to increase, thus increasing the wind density. In turn, this causes H to recombine and further reduce the wind temperature. When Z is sufficiently low, the wind is thinner, which allows more photons to ionize Fe ii into Fe iii. As noticeable in Fig. 12, the relative role of the pressure term tends to increase while the shift in the Fe ion dominance occurs, which we also found in the main model sequence (Fig. 5). A similar Z dependence of a solution regime change seems to occur for the Fe iv to Fe iii transition in the models in Lefever et al. (2025), where a sharp transition in the Ṁ trend occurs between 0.8 Z⊙ and 0.7 Z⊙, while a smoother behavior is found at the highest and lowest Z. This would indicate that for the lowest metallicity, all solutions are on the rarefied branch, while at the highest Z, the solutions might all be on the dense branch.

Similar to the solar metallicity sequence, we also found that existing late-type BHGs in the LMC and SMC resemble the normalized spectra resulting from our lower- Z models. Namely, for Z=0.2 Z⊙, we found a good match between the late BHGs AzV 393 (B3Ia+, SMC) and the synthetic UV and optical spectra (see Fig. 13). For Z=0.5 Z⊙, we also found signatures comparable to the late BHG Sk-68 8 (B5Ia+, LMC). The L, M, and XHe values of these stars are in particular very similar to our metallicity model sequence (Verhamme et al. 2024; Alkousa et al. 2026), such that the spectral mismatches, especially a systematic excess in the electron scattering wings, might be an indication of a different clumping stratification, which in case of hydrodynamically consistent models affect the output atmosphere structure and spectrum.

Thumbnail: Fig. 13 Refer to the following caption and surrounding text. Fig. 13

Nontailored comparison between the 0.5 Z⊙ and 0.2 Z⊙ models (olive and violet lines, respectively). The observed spectra of the BHGs SK −68°08 (B5Ia+, in the LMC) and AzV 393 (B3Ia+, in the SMC) are represented by black lines.

4.5 Comparisons with yellow hypergiants

When comparing the position in the HRD and sHRD of our models against the positions of some known YHGs, we noted that their empirically obtained Teff (e.g., Nieuwenhuijzen & de Jager 2000; Nieuwenhuijzen et al. 2012; Kasikov et al. 2024) and our coolest models overlap. Taking V509 Cas as a representative example, the effective temperature after the year 2000 remains consistently between ∼ 7.5 and ∼ 8.2 kK (Kasikov et al. 2024). Considering the Teff and log g values from Nieuwenhuijzen et al. (2012), we found that our model sequence crosses the variability track of this star in the sHRD parameter space essentially orthogonally (see Fig. 14).

In the lower panel of Fig. 14, we compare the spectrum of the YHG V509 Cas to that of the LBV MWC 930 during its maximum, which is closely resembled by our model with T*=14.6 kK. It is evident that the LBV-type spectrum essentially resembles the YHG spectrum, but is enhanced by wind features. Notably, both stars have T2/3 ≪ 10 kK, but in the LBV case, this is due to the photosphere being formed beyond the critical point, thus not representing the Teff of the (quasi-)hydrostatic atmosphere. Our model resembling MWC 930 has T2/3 ∼ 7 kK, while V509 Cas is formally hotter.

The formal distinction between supergiants and hypergiants in the yellow regime is different from that in the BA-star regime and is not based on the same spectroscopic criterion as was used for OBA hypergiants (e.g., de Jager 1998; Kasikov et al. 2026). Therefore, LBVs such as MWC 930 at its maximum might very well be what more extreme BHGs with dense winds are to BSGs or even BHGs in the rarefied regime. While our model sequence does not extend to cooler temperatures where Fe i becomes more important, the trend of finding rarefied solutions when the Fe ionization stage at r ∼ Rcrit is shifting suggests that another rarefied-wind regime might appear at even cooler T*, which might then produce observed YHG spectra such as that of V509 Cas.

Thumbnail: Fig. 14 Refer to the following caption and surrounding text. Fig. 14

Upper panel: comparison between our model sequence (brown diamonds), MWC 930 at its maximum (Miroshnichenko et al. 2014, cyan square), and the YHG V509 Cas after 2000 (Nieuwenhuijzen et al. 2012, golden pentagons) in the Kiel diagram. The cyan diamond marks the model that resembles MWC 930 best. Lower panel: comparison between the respective spectra. The brown plot shows our model that is most similar to MWC 930, and the golden plot shows the spectrum of V509 Cas in 2012. The spectrum of the LBV is plotted underneath the two as a dark cyan line for comparison. The thin dotted line in the background indicates the tracks of Ekström et al. (2012).

5 Conclusions

We conducted a thorough analysis with hydrodynamically consistent atmosphere models of a parameter space representative of OBA hypergiants and LBVs at solar metallicity to investigate their mass-loss rates and wind velocity structure in this crucial snapshot in the evolution of high-mass stars. We highlight the following conclusions:

  • We found two solution branches, which we described as dense and rarefied. In the explored temperature sequence, the models mainly settled on the dense solution, but were unable to maintain this in certain parameter regions. These regions appear as valleys in the Ṁ(Teff) diagram. The dense solutions are associated with higher wind optical depths, but do not require values exceeding unity. Instead, the rarefied solutions in this work can be described by Rt>200 R⊙, although more work is needed to verify whether this value is universal;

  • The normalized spectra synthesized from our models resemble the main spectral appearance of various OBA hypergiants and LBVs. Spectral resemblance is found for the dense and rarefied solutions. Exemplary calculations indicate that these findings might even hold at lower metallicity;

  • We found that the contributions of the electron scattering and thermal pressure components cause ∼ 50% of the force at the wind onset (i.e., critical point);

  • Our model sequence results agree with recent near-Eddington PoWRhd calculations employing different luminosities and masses than in our study. Sabhahit et al. (2026) also reported valley-like features in some of their Ṁ(T*) sequences, in particular, those matching our Γe-regime and farther away from the Eddington limit. The higher- Γe models calculated by Lefever et al. (2025) at solar metallicity do not display the pronounced valleys, but seem to be in a transition regime between our solutions and even higher Γe sequences in Sabhahit et al. (2026). However, the Ṁ jumps occurring at lower metallicities in Lefever et al. (2025) show clear similarities to the jumps obtained in this work;

  • In contrast to modeling efforts for OB stars, we found a general increase in the mass-loss rate with decreasing temperature. We suggest that this behavior might be exclusive to the dense branch. The solutions in the rarefied branch are notably supported by turbulent pressure, and our absolute results thus differ from calculations such as those by Björklund et al. (2023) and Krtička et al. (2024), although each of the recipes matches partially in some temperature regimes. Down to T* ∼ 22 kK, our derived mass-loss pattern matches the prediction by Vink et al. (2001) well when it is shifted downward by 0.8 dex;

  • The solution switched from rarefied to dense in our sequence when Fe recombined, such that Fe iii or Fe ii can provide significant opacity at the critical point. This sharp increase in the mass-loss rate coincides with the traditional bistability jumps from Vink et al. (1999). We also found an actual bistability of solutions at these jumps;

  • For LBVs (P Cygni supergiants) with dense winds on the cool side of the second jump, the terminal velocity is well traced by the H lines, but not necessarily by the He I lines;

  • We found that the empirical linear relation by Hawcroft et al. (2024) of the wind terminal velocity with the stellar photospheric temperature seems to be valid for all rarefied-winded stars that are confirmed to be in the line-driven wind regime (Teff ≳ 10 kK), that is, it is valid from late BHGs to WRs. However, the different classes of stars seem to have a different offset;

  • We provide updated luminosities for galactic OBA hypergiants and LBVs based on corrected Gaia DR3 parallaxes when available or applicable. The updated luminosity distribution causes the population of LBVs and hypergiants (assuming properties publicly available in the literature) to agree better with the LBV instability strip. However, careful and thorough analysis is necessary to confirm whether this agreement between theory and empirical findings actually holds.

Acknowledgements

The authors thank A. Kasikov for inspiring discussions and sharing some forthcoming results before publication. The authors thank S. Simón-Diaz and A. de Burgos for providing spectroscopic data from the IACOB project as well as inspiring discussions. The authors also thank J. Josiek, R. R. Lefever, G. González-Torà, E. Schösser, and D. Pauli for helpful discussions during this manuscript’s preparation. Moreover, the authors thank the anonymous referee for the constructive comments and suggestions, which were fundamental to improve the quality of this manuscript. MBP and AACS are supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) in the form of an Emmy Noether Research Group – Project-ID 445674056 (SA4064/11, PI Sander). VR acknowledges financial support by the Federal Ministry for Economic Affairs and Energy (BMWE) via the Deutsches Zentrum für Luftund Raumfahrt (DLR) grant 50 OR 2509 (PI Sander). MBP is a member of the International Max Planck Research School for Astronomy and Cosmic Physics at the University of Heidelberg (IMPRS-HD). This project was co-funded by the European Union (Project 101183150 – OCEANS). Based on observations collected at the European Southern Observatory under ESO programme(s). This work presents results from the European Space Agency (ESA) space mission Gaia. Gaia data are being processed by the Gaia Data Processing and Analysis Consortium (DPAC). Funding for the DPAC is provided by national institutions, in particular the institutions participating in the Gaia MultiLateral Agreement (MLA). FN acknowledges support by PID2022-137779OB-C41 funded by MCIN/AEI/10.13039/501100011033 by “ERDF A way of making Europe” and grant MAD4SPACE, TEC-2024/TEC-182 from Comunidad de Madrid (Spain).

References

  1. Agliozzo, C., Noriega-Crespo, A., Umana, G., et al. 2014, MNRAS, 440, 1391 [NASA ADS] [CrossRef] [Google Scholar]
  2. Agliozzo, C., Phillips, N., Mehner, A., et al. 2021, A&A, 655, A98 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Alkousa, T., Crowther, P. A., Bestenlehner, J. M., et al. 2025, A&A, 699, A314 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Alkousa, T., Crowther, P. A., Bestenlehner, J. M., et al. 2026, A&A, 707, A253 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Arun, R., Prasoon, S. A., Mathew, B., et al. 2025, MNRAS, 543, 3214 [Google Scholar]
  6. Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., Demleitner, M., & Andrae, R., 2021, AJ, 161, 147 [Google Scholar]
  7. Bakhytkyzy, A., Miroshnichenko, A. S., Klochkova, V. G., et al. 2025, Galaxies, 13, 17 [Google Scholar]
  8. Barniske, A., Oskinova, L. M., & Hamann, W.-R., 2008, A&A, 486, 971 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Berlanas, S. R., Herrero, A., Comerón, F., et al. 2020, A&A, 642, A168 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Bernini-Peron, M., Marcolino, W. L. F., Sander, A. A. C., et al. 2023, A&A, 677, A50 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Bernini-Peron, M., Sander, A. A. C., Ramachandran, V., et al. 2024, A&A, 692, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Bernini-Peron, M., Sander, A. A. C., Najarro, F., et al. 2025, A&A, 697, A41 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Björklund, R., Sundqvist, J. O., Puls, J., & Najarro, F., 2021, A&A, 648, A36 [EDP Sciences] [Google Scholar]
  14. Björklund, R., Sundqvist, J. O., Singh, S. M., Puls, J., & Najarro, F., 2023, A&A, 676, A109 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Cardelli, J. A., Clayton, G. C., & Mathis, J. S., 1989, ApJ, 345, 245 [Google Scholar]
  16. Carretero-Castrillo, M., Ribó, M., & Paredes, J. M., 2023, A&A, 679, A109 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Castor, J. I., Abbott, D. C., & Klein, R. I., 1975, ApJ, 195, 157 [Google Scholar]
  18. Chentsov, E. L., & Luud, L., 1989, Astrophysics, 31, 415 [Google Scholar]
  19. Clark, J. S., & Negueruela, I., 2004, A&A, 413, L15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Clark, J. S., Egan, M. P., Crowther, P. A., et al. 2003, A&A, 412, 185 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Clark, J. S., Crowther, P. A., Larionov, V. M., et al. 2009a, A&A, 507, 1555 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Clark, J. S., Davies, B., Najarro, F., et al. 2009b, A&A, 504, 429 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Clark, J. S., Najarro, F., Negueruela, I., et al. 2012, A&A, 541, A145 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Conti, P. S., 1984, in IAU Symposium, 105, Observational Tests of the Stellar Evolution Theory, eds. A. Maeder, & A. Renzini, 233 [Google Scholar]
  25. Crowther, P. A., Lennon, D. J., & Walborn, N. R., 2006, A&A, 446, 279 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Curé, M., Cidale, L., & Granada, A., 2011, ApJ, 737, 18 [CrossRef] [Google Scholar]
  27. Cutri, R. M., Wright, E. L., Conrow, T., et al. 2013, Explanatory Supplement to the AllWISE Data Release Products, by R. M. Cutri, et al. [Google Scholar]
  28. de Burgos, A., Keszthelyi, Z., Simón-Díaz, S., & Urbaneja, M. A., 2024, A&A, 687, L16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. de Jager, C., 1998, A&A Rev., 8, 145 [Google Scholar]
  30. de Jager, C., Nieuwenhuijzen, H., & van der Hucht, K. A., 1988, A&AS, 72, 259 [Google Scholar]
  31. Debnath, D., Sundqvist, J. O., Moens, N., et al. 2024, A&A, 684, A177 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Driessen, F. A., Sundqvist, J. O., & Kee, N. D., 2019, A&A, 631, A172 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Egan, M. P., Clark, J. S., Mizuno, D. R., et al. 2002, ApJ, 572, 288 [Google Scholar]
  34. Ekström, S., Georgy, C., Eggenberger, P., et al. 2012, A&A, 537, A146 [Google Scholar]
  35. Fadeyev, Y. A., 2010, Astron. Lett., 36, 362 [Google Scholar]
  36. Friend, D. B., & Abbott, D. C., 1986, ApJ, 311, 701 [Google Scholar]
  37. Gaia Collaboration (Vallenari, A., et al.,) 2023, A&A, 674, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Glatzel, W., & Kiriakidis, M., 1993, MNRAS, 263, 375 [NASA ADS] [Google Scholar]
  39. Glatzel, W., & Kraus, M., 2024, MNRAS, 529, 4947 [Google Scholar]
  40. Gräfener, G., Koesterke, L., & Hamann, W. R., 2002, A&A, 387, 244 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Grassitelli, L., Langer, N., Mackey, J., et al. 2021, A&A, 647, A99 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Groh, J. H., & Vink, J. S., 2011, A&A, 531, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Groh, J. H., Hillier, D. J., Damineli, A., et al. 2009, ApJ, 698, 1698 [NASA ADS] [CrossRef] [Google Scholar]
  44. Groh, J. H., Hillier, D. J., & Damineli, A., 2011, ApJ, 736, 46 [CrossRef] [Google Scholar]
  45. Groh, J. H., Hillier, D. J., Madura, T. I., & Weigelt, G., 2012, MNRAS, 423, 1623 [NASA ADS] [CrossRef] [Google Scholar]
  46. Gvaramadze, V. V., & Menten, K. M., 2012, A&A, 541, A7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Haucke, M., Cidale, L. S., Venero, R. O. J., et al. 2018, A&A, 614, A91 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Hawcroft, C., Sana, H., Mahy, L., et al. 2024, A&A, 688, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Hillier, D. J., Crowther, P. A., Najarro, F., & Fullerton, A. W., 1998, A&A, 340, 483 [NASA ADS] [Google Scholar]
  50. Holgado, G., Simón-Díaz, S., Haemmerlé, L., et al. 2020, A&A, 638, A157 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Humphreys, R. M., & Davidson, K., 1979, ApJ, 232, 409 [Google Scholar]
  52. Hur, H., Lim, B., & Chun, M.-Y., 2023, J. Korean Astron. Soc., 56, 97 [NASA ADS] [Google Scholar]
  53. Kasikov, A., Kolka, I., Aret, A., Eenmäe, T., & Checha, V., 2024, A&A, 686, A270 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Kasikov, A., Mehner, A., Kolka, I., & Aret, A., 2026, A&A, in press, https://doi.org/10.1051/0004-6361/202558527 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Kee, N. D., Sundqvist, J. O., Decin, L., de Koter, A., & Sana, H., 2021, A&A, 646, A180 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Klochkova, V. G., 2019, Astrophys. Bull., 74, 475 [Google Scholar]
  57. Kniazev, A. Y., Gvaramadze, V. V., & Berdnikov, L. N., 2016, MNRAS, 459, 3068 [NASA ADS] [CrossRef] [Google Scholar]
  58. Kotak, R., & Vink, J. S., 2006, A&A, 460, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Koumpia, E., Oudmaijer, R. D., Graham, V., et al. 2020, A&A, 635, A183 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Krtička, J., Kubát, J., & Krtičková, I., 2024, A&A, 681, A29 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  61. Krtička, J., Kubát, J., & Krtičková, I., 2025, A&A, 702, A9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Lamers, H. J. G. L. M., Snow, T. P., & Lindholm, D. M., 1995, ApJ, 455, 269 [Google Scholar]
  63. Langer, N., 2012, ARA&A, 50, 107 [Google Scholar]
  64. Larkin, C. J. K., Hawcroft, C., Mackey, J., et al. 2025a, A&A, 703, L14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Larkin, C. J. K., Mackey, J., Haworth, T. J., & Sander, A. A. C. 2025b, A&A, 700, A60 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Lefever, K., Puls, J., & Aerts, C., 2007, A&A, 463, 1093 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Lefever, R. R., Sander, A. A. C., Shenar, T., et al. 2023, MNRAS, 521, 1374 [NASA ADS] [CrossRef] [Google Scholar]
  68. Lefever, R. R., Sander, A. A. C., Bernini-Peron, M., et al. 2025, A&A, 700, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  69. Lefever, R. R., Sander, A. A. C., Bernini-Peron, M., et al. 2026, A&A, 707, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Leung, H. W., Bovy, J., Mackereth, J. T., et al. 2023, MNRAS, 519, 948 [Google Scholar]
  71. Lobel, A., Groh, J. H., Martayan, C., et al. 2013, A&A, 559, A16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Mahy, L., Hutsemékers, D., Royer, P., & Waelkens, C., 2016, A&A, 594, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Markova, N., & Puls, J., 2008, A&A, 478, 823 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Martins, F., Genzel, R., Hillier, D. J., et al. 2007, A&A, 468, 233 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Maryeva, O. V., Karpov, S. V., Kniazev, A. Y., & Gvaramadze, V. V., 2022, MNRAS, 513, 5752 [Google Scholar]
  76. Miroshnichenko, A. S., Manset, N., Zharikov, S. V., et al. 2014, Adv. Astron., 2014, E7 [Google Scholar]
  77. Miroshnichenko, A. S., Zharikov, S. V., Vaidman, N. L., & Khokhlov, S. A., 2025, arXiv e-prints [arXiv:2511.13058] [Google Scholar]
  78. Moens, N., Poniatowski, L. G., Hennicker, L., et al. 2022, A&A, 665, A42 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  79. Moens, N., Debnath, D., Verhamme, O., et al. 2025, A&A, 704, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Müller, P. E., & Vink, J. S., 2008, A&A, 492, 493 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Najarro, F., 1995, PhD thesis [Google Scholar]
  82. Najarro, F., Krabbe, A., Genzel, R., et al. 1997, A&A, 325, 700 [NASA ADS] [Google Scholar]
  83. Najarro, F., Figer, D. F., Hillier, D. J., Geballe, T. R., & Kudritzki, R. P., 2009, ApJ, 691, 1816 [NASA ADS] [CrossRef] [Google Scholar]
  84. Nazé, Y., Rauw, G., & Hutsemékers, D., 2012, A&A, 538, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  85. Negueruela, I., Alfaro, E. J., Dorda, R., et al. 2022, A&A, 664, A146 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  86. Negueruela, I., Simón-Díaz, S., de Burgos, A., Casasbuenas, A., & Beck, P. G., 2024, A&A, 690, A176 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  87. Nieuwenhuijzen, H., & de Jager, C., 2000, A&A, 353, 163 [NASA ADS] [Google Scholar]
  88. Nieuwenhuijzen, H., De Jager, C., Kolka, I., et al. 2012, A&A, 546, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  89. Nyholm, A., Sollerman, J., Tartaglia, L., et al. 2020, A&A, 637, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  90. Pasquali, A., Comerón, F., & Nota, A., 2006, A&A, 448, 589 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. Pauldrach, A. W. A., & Puls, J., 1990, A&A, 237, 409 [NASA ADS] [Google Scholar]
  92. Pauldrach, A., Puls, J., & Kudritzki, R. P., 1986, A&A, 164, 86 [NASA ADS] [Google Scholar]
  93. Paunzen, E., 2022, A&A, 661, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  94. Petriella, A., Paron, S. A., & Giacani, E. B., 2012, A&A, 538, A14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  95. Petrov, B., Vink, J. S., & Gräfener, G., 2016, MNRAS, 458, 1999 [NASA ADS] [CrossRef] [Google Scholar]
  96. Prinja, R. K., Barlow, M. J., & Howarth, I. D., 1990, ApJ, 361, 607 [NASA ADS] [CrossRef] [Google Scholar]
  97. Ritchie, B. W., Clark, J. S., Negueruela, I., & Najarro, F., 2009, A&A, 507, 1597 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  98. Rivet, J.-P., Siciak, A., de Almeida, E. S. G., et al. 2020, MNRAS, 494, 218 [NASA ADS] [CrossRef] [Google Scholar]
  99. Rubio-Díez, M. M., Sundqvist, J. O., Najarro, F., et al. 2022, A&A, 658, A61 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  100. Sabhahit, G. N., Vink, J. S., & Sander, A. A. C., 2026, A&A, 706, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  101. Sander, A. A. C., & Vink, J. S., 2020, MNRAS, 499, 873 [Google Scholar]
  102. Sander, A. A. C., Fürst, F., Kretschmar, P., et al. 2018, A&A, 610, A60 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  103. Sander, A. A. C., Vink, J. S., & Hamann, W.-R., 2020, MNRAS, 491, 4406 [Google Scholar]
  104. Sander, A. A. C., Lefever, R. R., Poniatowski, L. G., et al. 2023, A&A, 670, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  105. Searle, S. C., Prinja, R. K., Massa, D., & Ryans, R., 2008, A&A, 481, 777 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  106. Smith, N., 2006, ApJ, 644, 1151 [Google Scholar]
  107. Smith, N., 2014, ARA&A, 52, 487 [Google Scholar]
  108. Smith, N., & Stassun, K. G., 2017, AJ, 153, 125 [Google Scholar]
  109. Smith, N., Aghakhanloo, M., Murphy, J. W., et al. 2019, MNRAS, 488, 1760 [NASA ADS] [CrossRef] [Google Scholar]
  110. Sterken, C., Gosset, E., Juttner, A., et al. 1991, A&A, 247, 383 [NASA ADS] [Google Scholar]
  111. Stothers, R. B., 2003, ApJ, 589, 960 [NASA ADS] [CrossRef] [Google Scholar]
  112. Sundqvist, J. O., Debnath, D., Backs, F., et al. 2025, A&A, 703, A284 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  113. Trundle, C., Kotak, R., Vink, J. S., & Meikle, W. P. S., 2008, A&A, 483, L47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  114. Ueta, T., Meixner, M., Dayal, A., et al. 2001, ApJ, 548, 1020 [NASA ADS] [CrossRef] [Google Scholar]
  115. van Genderen, A. M., 2001, A&A, 366, 508 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  116. Verhamme, O., Sundqvist, J., de Koter, A., et al. 2024, A&A, 692, A91 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  117. Vieu, T., Larkin, C. J. K., Härer, L., et al. 2024, MNRAS, 532, 2174 [NASA ADS] [CrossRef] [Google Scholar]
  118. Vink, J. S., 2022, ARA&A, 60, 203 [Google Scholar]
  119. Vink, J. S., & de Koter, A., 2002, A&A, 393, 543 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  120. Vink, J. S., & Sander, A. A. C., 2021, MNRAS, 504, 2051 [NASA ADS] [CrossRef] [Google Scholar]
  121. Vink, J. S., de Koter, A., & Lamers, H. J. G. L. M., 1999, A&A, 350, 181 [Google Scholar]
  122. Vink, J. S., de Koter, A., & Lamers, H. J. G. L. M., 2001, A&A, 369, 574 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  123. Weßmayer, D., Przybilla, N., & Butler, K., 2022, A&A, 668, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  124. Weßmayer, D., Przybilla, N., Ebenbichler, A., Aschenbrenner, P., & Butler, K., 2023, A&A, 677, A175 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  125. Zacharias, N., Monet, D. G., Levine, S. E., et al. 2004, in American Astronomical Society Meeting Abstracts, 205, 48.15 [Google Scholar]
1

We use the definition for hypergiants by Clark et al. (2012) based exclusively on spectral morphology (cf. Negueruela et al. 2024).

2

We did not label them “thick” and “thin” as they do not coincide with the wind-optical depth τF crossing unity. The thick solutions have a higher τF, but there is no uniform transition value among the sequence.

Appendix A Atomic data

In Table A.1, we list the elements and ions with the number of levels and line transitions considered in our atmosphere models. For Fe – which in PoWR refers to all iron-group elements conflated-, the levels and line numbers refer to superlevels and “superlines” (see Gräfener et al. 2002, for more details). The vdop parameter, which controls the resolution of the frequency-dependent cross section of all the elements (Gräfener et al. 2002) is set to the same value of the turbulent motion – i.e., vdop=20 km s−1.

Table A.1

Ions, levels, and transitions (lines) considered in the models

Appendix B Wind structure in different regimes

In Fig. 5 we compare Ṁ with the contribution of the main processes, elements, and ions to the wind acceleration at the critical point for each value of T*. In Fig. B.1 we provide an extension of this figure showing the respective main contributions to the wind driving at different points in the atmosphere. Namely, at τ=20, near the inner boundary (τ=50), probing the forces in the inner atmosphere, (ii) at the critical point, where the wind is effectively launched, and (iii) at τ=0.01, which probes conditions far out in the wind.

A more detailed view on the acceleration structure for representative models at each regime described in Sect. 3 is given in Fig. B.2. In Fig. B.3 we illustrate the velocity structure of models in different regimes. One can see that the turbulent velocity is equal to or larger than the sound speed in the wind, indicating that the turbulent pressure dominates over the gas thermal pressure. This is particularly important for models in the cold wall.

Thumbnail: Fig. B.1 Refer to the following caption and surrounding text. Fig. B.1

Upper panel: Contribution of different processes, elements, and ions to the wind work ratio Qwind. The other panels indicate each contribution to the wind driving at specific depths – namely, at τ=30 (near the inner boundary), τ=τcrit (where the wind is launched, like in Fig. 5), and τ=0.01 (far out in the wind). The color code is the same of Fig. 5.

Appendix C Bistability of the solutions at the jumps

When computing hydrodynamically consistent models, the complex interplay between radiation field, species populations, and velocity field does not necessarily converge to unique solutions.

Thumbnail: Fig. B.2 Refer to the following caption and surrounding text. Fig. B.2

Acceleration structure and main contributions to the driving of the winds for models in different wind-solutions domains.

Such complexity is due to the nonlinearity of the equations being solved, coupled in space and frequency. Therefore, the radiative acceleration structure, obtained via the integration of the momentum equation inward and outward from the critical point, depends on conditions elsewhere in the atmosphere.

This means that, for instance, different initial physical conditions can yield different results (with them all satisfying the criteria for numerical convergence) or even fail to reach convergence in case the initial “guess” is far from a valid solution. Given the high computational cost of each model, exhaustively mapping these degeneracies is beyond the scope of this work. However, as a sanity check of the obtained trend, we compute models with different initial conditions in distinct regions of our Teff sequence.

Thumbnail: Fig. B.3 Refer to the following caption and surrounding text. Fig. B.3

Wind velocity structure (black line) for models in the cold wall, hot valley, and hot plateau wind-solution regions. The beige line indicates the sonic speed structure (cs) while the dotted beige line indicate the turbulent velocity vturb=20 km s−1. The thin vertical lines indicate the sonic radius (Rsonic, beige) and critical radius (Rcrit, black).

This is particularly important when operating in the transition between the rarefied and dense solutions. In Fig. C.1, we show that at the edge of the 1-BiSJ and 2-BiSJ, the same T* can produce different values of Ṁ when the initial conditions are a model on the rarefied or on the dense branch. This also manifests in the spectral appearance of the models’ synthetic spectra, which also match different observed stars.

The strong dependence on initial conditions prohibits us from pinpointing one-to-one relations between temperatures (or other parameters) and wind properties. This exercise of computing different models with different initial conditions can also provide a certain interval of confidence in our Ṁ(Teff) behavior and the associated transition between regimes. However, it indicates that stars with similar properties may exist in different regimes – which recycles the original discussions on bistability by Pauldrach & Puls (1990).

Appendix C.1 Mechanical wind luminosity

To investigate the impact of hypergiants on their environment, we calculate the mechanical wind luminosity Lmech:=0.5M˙v∞2Mathematical equation: $L_{\text {mech }}:=0.5 \dot{M} v_{\infty}^{2}$. In general, up to the cool rarefied valley, we find Lmech ∼ 1036 erg s−1, comparable to late-O and B stars (see, e.g., Vieu et al. 2024). In the valley, Lmech plunges due to the combination of very low Ṁ and v∞. Past the 2-BiSJ, the high mass-loss rates reached in the cold wall region invert the trend and cause Lmech to increase up to ∼ 1036 erg s−1. However, the low terminal velocities lead to generally lower Lmech values, making hotter stars and in particular Wolf-Rayet stars significantly more important with respect to kinetic feedback unless much higher luminosities and thus mass-loss rates would be reached than mapped in our model sequence. Nonetheless, significant mass loss from cooler massive stars can affect cluster properties via mass loading (e.g., Larkin et al. 2025a,b).

Thumbnail: Fig. C.1 Refer to the following caption and surrounding text. Fig. C.1

Demonstration of the bistability of solutions in the 1-BiSJ (upper panel) and 2-BiSJ (lower panel). In either case, we find models with the same T* yielding different results based on their initial conditions.
Thumbnail: Fig. C.2 Refer to the following caption and surrounding text. Fig. C.2

Manifestation of the bistability of solutions in the model spectrum for the 1-BiSJ in our sequence.

Appendix D Luminosity of GS Mus (O9.7la+) and HD151804 (O8laf)

In the literature, we were able to find spectroscopic properties of GS Mus and HD151804 obtained via analysis of the normalized spectra (e.g., Prinja et al. 1990; Holgado et al. 2020). However, we could not find luminosities based on SED fitting.

Thumbnail: Fig. C.3 Refer to the following caption and surrounding text. Fig. C.3

Analogous to Fig. C.2 but for the 2-BiSJ.
Thumbnail: Fig. C.4 Refer to the following caption and surrounding text. Fig. C.4

Mechanical luminosity of the model sequence compared to literature data from Crowther et al. (2006), Berlanas et al. (2020) – light blue circles–, and WR Lmech estimates in Vieu et al. (2024) – pink squares and arrows, which indicate WRs outside the plotting area.

For GS Mus Lefever et al. (2007) derived log L/L⊙=5.75 based on the V magnitude, spectral type and luminosity class. As our models yielded similar T2/3 and log g, we estimated their luminosities via fitting the optical and IR photometric fluxes and flux-calibrated far-UV – see Fig. D.1 and Fig. D.2. The reddening was estimated via fitting the near-UV 2175 Å absorption bump applying the extinction law by Cardelli et al. (1989). As we do not perform a dedicated fitting, we consider a conservative uncertainty of 0.2 dex to the luminosity of each star.

Thumbnail: Fig. D.1 Refer to the following caption and surrounding text. Fig. D.1

Spectral energy distribution of GS Mus. The photometry used in were sourced from Gaia Collaboration (2023), Paunzen (2022), Zacharias et al. (2004), and Cutri et al. (2013).
Thumbnail: Fig. D.2 Refer to the following caption and surrounding text. Fig. D.2

Spectral energy distribution of HD151804. The photometry is sourced from the same references as in Fig. D.1.

Appendix E Model fundamental wind parameters

In Tables E.1 and E.2, we disclose the fundamental atmospheric properties of the sequences of models discussed in this paper.

Appendix F Luminosities of LBVs and cLBVs with updated distances

To construct the HR diagram in Fig. 1, we first collected literature values of effective temperatures and luminosities of OBA and yellow hypergiants, LBVs, B[e]SG, and yellow supergiants (YSGs). The list of OBA hypergiants and LBVs is largely based on previous literature compilations from Nazé et al. (2012) and Smith et al. (2019). The B[e]SG stars are sourced from Miroshnichenko et al. (2025) and the YHGs and YSGs values are taken from Kasikov et al. (2026). We updated the luminosities based on updated distance measurements – mostly from the catalog of Bailer-Jones et al. (2021) based on Gaia DR3 parallaxes. We considered the geometric distances listed in their catalog. In Table F.1 we disclose the literature stellar temperatures, luminosities (with the respective distances), and updated luminosities adopting the corresponding new distances.

In general, we do not find a systematic increase in the obtained luminosities compared to previous measurements, a result that differs from what Smith et al. (2019) found when updating previous literature values to Gaia DR2-based distances. Namely, they found an overall reduction in luminosity. However, we find that the population of LBVs with updated distances tends to agree much more with the LBV instability strip (Groh et al. 2009). Moreover it seems to connect smoothly with the brightest hot YSGs (the edge of the so-called Yellow Void de Jager 1998), whose luminosities were recently updated by Kasikov et al. (2026).

Table E.1

Atmosphere fundamental properties for each model in the T* sequence.

Table E.2

Atmosphere fundamental properties for each model in the metallicity sequence for a fixed T*=16.1 kK.

Table F.1

Literature properties and updated distances and luminosities of LBVs and OBA hypergiants.

All Tables

Table A.1

Ions, levels, and transitions (lines) considered in the models

In the text
Table E.1

Atmosphere fundamental properties for each model in the T* sequence.

In the text
Table E.2

Atmosphere fundamental properties for each model in the metallicity sequence for a fixed T*=16.1 kK.

In the text
Table F.1

Literature properties and updated distances and luminosities of LBVs and OBA hypergiants.

In the text

All Figures

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

Hertzsprung-Russell diagram depicting our model sequence against LBVs (cyan squares), B[e] SG (blue triangles), OBA hypergiants or LBV candidates (green circles), YSGs (small yellow pentagons), and YHGs (large gold pentagons) from the literature (see Appendix F). The blue lines connecting the LBV data indicate the same star in different epochs. Likewise, the dark orange lines connect YHGs in different epochs. The luminosities are updated to distances mostly based on Gaia DR3 (Bailer-Jones et al. 2021; Gaia Collaboration 2023). The black line indicates the luminosity of our sequence. The small diamonds below indicate T*, and the large open diamonds indicate T2/3. The thin dotted lines indicate the Galactic evolutionary tracks by Ekström et al. (2012). The thick cyan line indicates the LBV instability strip (Groh et al. 2011), and the thick violet line indicates the HD limit (Humphreys & Davidson 1979).
In the text
Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

Effective temperature at different optical depths τ, normalized by T*=Teff(τ=50) for each point of the sequence. The pink curve indicates T2/3, and the blue curve indicates Teff at the critical point.
In the text
Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Nontailored comparison between the output spectra of some of the models (red lines) and observed spectra (black lines) of OB supergiant/hypergiants and LBVs.
In the text
Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Transformed radius (Rt) vs. T*. Models with Rt<200 R⊙ are in the dense-wind regime, i.e., with higher Ṁ and optical spectra crowded with emission and P Cygni features (shaded region). The colored texts in the lower part of the diagram name the wind-solution regions identified in the T* sequence. The 1-BiSJ and the 2-BiSJ are indicated by the dotted and dashed lines in the transitions from a valley to a plateau (i.e., increase in Ṁ) toward lower T*.
In the text
Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Contributions of different processes, elements, and ions to the wind launching at the critical point (stacked area plot) and corresponding mass-loss rates (red line) for different inner boundary temperatures (T*). From bottom to top, the light beige area indicates the electron-scattering acceleration, dark beige shows the pressure, and smoke white indicates H i continuum driving. The teal tones indicate the contribution of CNO and other metals, and the shades of blue indicate the contributions of different Fe ions.
In the text
Thumbnail: Fig. 6 Refer to the following caption and surrounding text. Fig. 6

Comparison between the trend of mass-loss rates (red curve) and the classical Eddington parameter at the critical point Γe(Rcrit) (thick dotted curve) and the total contribution of the electron scattering to the work ratio QewindMathematical equation: $Q_{\mathrm{e}}^{\text {wind }}$ (shaded area).
In the text
Thumbnail: Fig. 7 Refer to the following caption and surrounding text. Fig. 7

Terminal velocities (upper panel) and velocity structures (lower panel) of the models along the T* sequence. Each wind-solution domain is color-coded following the scheme of Fig. 4.
In the text
Thumbnail: Fig. 8 Refer to the following caption and surrounding text. Fig. 8

Hydrogen and helium profiles in velocity space compared to the terminal velocity (black dots). The color of each spectrum marks the wind-solution regime to which the model belongs, following the same color scheme as Fig. 7. The temperatures rightward in the upper panels indicate T2/3, and rightward in the lower panel, they indicate T*.
In the text
Thumbnail: Fig. 9 Refer to the following caption and surrounding text. Fig. 9

Behavior of the v∞/vesc,Γe relation. The gray crosses indicate the data of Lamers et al. (1995) for OBA supergiants. The blue curve with diamonds indicates the values for vesc,Γe applied to the critical point. The pink curve indicates the ratio at τ=2/3, and the thin light blue curve indicates the ratio at τ=50, i.e., the inner boundary.
In the text
Thumbnail: Fig. 10 Refer to the following caption and surrounding text. Fig. 10

Relation of v∞−Teff for different studies that made hydrodynamically consistent sequences of models in which T* varied. The black lines indicate the empirical relation by Hawcroft et al. (2024) (the solid line represents the actual relation, and the dashed line indicates the relation shifted arbitrarily). The pink and violet squares indicate the models by Sander et al. (2023) and Lefever et al. (2025), respectively. The blue circle indicates the models of Sabhahit et al. (2026) with Γe ∼ 0.35, and the green circle indicates the models of Björklund et al. (2021) with similar Γe. The red and brown diamonds indicate our sequence of models and show the rarefied and dense solutions, respectively. The lower panel is a zoom-in of the region within the dotted rectangle in the upper panel.
In the text
Thumbnail: Fig. 11 Refer to the following caption and surrounding text. Fig. 11

Comparison between the output mass-loss rates of the T* sequence (red curve) and different recipes applied to the stellar properties of each model. The thick purple and fuchsia lines indicate the mass-loss rates by Vink et al. (2001), applied to T* and T2/3, respectively. The thin lines of the same colors indicate the Vink rates with an offset of 0.8 dex. The yellow curve indicates the rates of Björklund et al. (2023), and the teal line indicates the rates of Krtička et al. (2024), both applied to T2/3. The green line indicates the rates computed using the tool LIME (https://lime.ster.kuleuven.be/, Sundqvist et al. 2025).
In the text
Thumbnail: Fig. 12 Refer to the following caption and surrounding text. Fig. 12

Contribution to the radiative acceleration from different sources, elements, and ions (stacked colored areas) and mass-loss rates (red curve) according to metallicity for a fixed temperature of T*= 16.1 kK.
In the text
Thumbnail: Fig. 13 Refer to the following caption and surrounding text. Fig. 13

Nontailored comparison between the 0.5 Z⊙ and 0.2 Z⊙ models (olive and violet lines, respectively). The observed spectra of the BHGs SK −68°08 (B5Ia+, in the LMC) and AzV 393 (B3Ia+, in the SMC) are represented by black lines.
In the text
Thumbnail: Fig. 14 Refer to the following caption and surrounding text. Fig. 14

Upper panel: comparison between our model sequence (brown diamonds), MWC 930 at its maximum (Miroshnichenko et al. 2014, cyan square), and the YHG V509 Cas after 2000 (Nieuwenhuijzen et al. 2012, golden pentagons) in the Kiel diagram. The cyan diamond marks the model that resembles MWC 930 best. Lower panel: comparison between the respective spectra. The brown plot shows our model that is most similar to MWC 930, and the golden plot shows the spectrum of V509 Cas in 2012. The spectrum of the LBV is plotted underneath the two as a dark cyan line for comparison. The thin dotted line in the background indicates the tracks of Ekström et al. (2012).
In the text
Thumbnail: Fig. B.1 Refer to the following caption and surrounding text. Fig. B.1

Upper panel: Contribution of different processes, elements, and ions to the wind work ratio Qwind. The other panels indicate each contribution to the wind driving at specific depths – namely, at τ=30 (near the inner boundary), τ=τcrit (where the wind is launched, like in Fig. 5), and τ=0.01 (far out in the wind). The color code is the same of Fig. 5.
In the text
Thumbnail: Fig. B.2 Refer to the following caption and surrounding text. Fig. B.2

Acceleration structure and main contributions to the driving of the winds for models in different wind-solutions domains.
In the text
Thumbnail: Fig. B.3 Refer to the following caption and surrounding text. Fig. B.3

Wind velocity structure (black line) for models in the cold wall, hot valley, and hot plateau wind-solution regions. The beige line indicates the sonic speed structure (cs) while the dotted beige line indicate the turbulent velocity vturb=20 km s−1. The thin vertical lines indicate the sonic radius (Rsonic, beige) and critical radius (Rcrit, black).
In the text
Thumbnail: Fig. C.1 Refer to the following caption and surrounding text. Fig. C.1

Demonstration of the bistability of solutions in the 1-BiSJ (upper panel) and 2-BiSJ (lower panel). In either case, we find models with the same T* yielding different results based on their initial conditions.
In the text
Thumbnail: Fig. C.2 Refer to the following caption and surrounding text. Fig. C.2

Manifestation of the bistability of solutions in the model spectrum for the 1-BiSJ in our sequence.
In the text
Thumbnail: Fig. C.3 Refer to the following caption and surrounding text. Fig. C.3

Analogous to Fig. C.2 but for the 2-BiSJ.
In the text
Thumbnail: Fig. C.4 Refer to the following caption and surrounding text. Fig. C.4

Mechanical luminosity of the model sequence compared to literature data from Crowther et al. (2006), Berlanas et al. (2020) – light blue circles–, and WR Lmech estimates in Vieu et al. (2024) – pink squares and arrows, which indicate WRs outside the plotting area.
In the text
Thumbnail: Fig. D.1 Refer to the following caption and surrounding text. Fig. D.1

Spectral energy distribution of GS Mus. The photometry used in were sourced from Gaia Collaboration (2023), Paunzen (2022), Zacharias et al. (2004), and Cutri et al. (2013).
In the text
Thumbnail: Fig. D.2 Refer to the following caption and surrounding text. Fig. D.2

Spectral energy distribution of HD151804. The photometry is sourced from the same references as in Fig. D.1.
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.