Open Access
Issue
A&A
Volume 677, September 2023
Article Number L5
Number of page(s) 9
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/202347454
Published online 31 August 2023

© The Authors 2023

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. Subscribe to A&A to support open access publication.

1. Introduction

Over the next decades, stellar structure and evolution (SSE) theory will begin to advance from one-dimensional models towards solving the stellar structure equations in three spatial dimensions. So far, all SSE codes rely on the approximation of spherical symmetry as this greatly simplifies the numerics. But asteroseismic and interferometric studies of stars on the main sequence (MS) have revealed that a significant fraction of stars rotate with velocities where the approximation of spherical symmetry is not justified. Many SSE codes that account for rotational effects are based on the pioneering works by Zahn (1992) and Chaboyer & Zahn (1992), and sequential works (e.g., Talon & Charbonnel 1998; Palacios et al. 2003; Mathis et al. 2004), which provide expressions for the efficiency of chemical mixing due to both shear-driven turbulence and advective transport (cf. Ekström et al. 2012). As this paradigm is constructed for one dimension, several assumptions have been made regarding the relative strengths of the horizontal and vertical diffusion, and the validity domain of these assumptions is currently unknown, except that they should apply in the limit of slow rotation.

Attempts at solving the stellar structure of rotating stars in more than one dimension began in the sixties, but it only recently became possible to compute the stellar structure in a self-consistent manner in two dimensions, that is, with the ESTER1 code (Espinosa Lara & Rieutord 2013; Rieutord et al. 2016). These models have been restricted to a steady state, where the hydrogen mass fraction profile is modelled as a step function, namely a constant value in the convective core and a constant value in the radiative envelope (see e.g., Gagnier et al. 2019a,b; Bouchaud et al. 2020; Howarth et al. 2023). Furthermore, ESTER is designed for early-type stars and the applicability is currently limited to stars with a convective core and radiative envelope.

In this Letter, we present the results of a novel approach whereby the temporal dimension is added to the two-dimensional steady-state models by solving, in addition, the equation of chemical evolution. The overall aim of the present work is to demonstrate the capabilities of the state-of-the-art 2D evolution ESTER models of a rotating massive star, and present predictions for the rotation and chemical profiles. This includes fully solving the large-scale velocity field inside the star. We successfully performed 2D evolution of a slow and a moderate rotator with a mass of 12 M, from the zero-age main sequence (ZAMS) up to near the terminal-age main sequence (TAMS). For the predicted rotation profile, we also aim to make a comparison with observations. In Sect. 2, we describe the fundamentals of the ESTER code and the improvements made. In Sect. 3, we compare the predictions of the 2D models, especially the predicted rotation profile, with asteroseismic measurements of a massive star. We show that the predictions of the 2D model are consistent with the observations. In Sect. 4, we discuss the contribution of meridional circulation to nitrogen enrichment as predicted by 2D models. Finally, we conclude in Sect. 5.

2. ESTER models

In this work, we make use of new ESTER models where time evolution has been implemented. We recall that ESTER models are 2D models of rotating stars that include the centrifugal distortion of the star and the associated baroclinic flows. These models therefore predict the structure (pressure, density, and temperature distributions) and the associated large-scale flows, namely the differential rotation and meridional circulations. The first version of the ESTER code, which computes the steady states of early-type star models (Espinosa Lara & Rieutord 2013), has now been completed with a time evolutive version that includes the unsteady terms needed to follow the thermal and nuclear evolution of a star. The code solves the set of equations given in Appendix A, still assuming the axisymmetry of the star. Additionally, in this new version, the evolution of the mass fraction of hydrogen is solved,

(1)

where ρ is the density, v is the meridional velocity, and [D] is a tensor constructed from the horizontal and vertical chemical diffusion coefficients.

As for the steady-state version, we are using OPAL for opacities and equation of state (Rogers et al. 1996). Nuclear energy production is modelled via a simple law for the CNO cycle (see appendix) from which we derive nuc.

The spatial discretization is based on a spectral element method where spectral elements are spheroidal shells bounded by isobars Rieutord et al. (2016). Typically, a model uses 12 spheroidal shells with 30 points in a Gauss–Lobatto grid radially and 24 points in a Gauss–Legendre grid in latitude. Time evolution is insured by a first-order backward Euler method where the time step is manually set to 0.5 Myr, and is decreased when the model fails to converge.

Following the theoretical work of Zahn (1992) on the predicted diffusion constant resulting from rotationally induced chemical mixing, and the 1D implementation by Mombarg et al. (2022), the vertical chemical diffusivity is taken as,

(2)

Here, K is the thermal diffusivity, n is the unit vector normal to the isobars, and Ω is local angular velocity. The term is the volume-averaged squared Brunt–Väisälä frequency of the initial steady-state model, where the average is taken over the volume where . We therefore smooth out the rapid radial variations of N2 that often give rise to numerical difficulties. Furthermore, η is a free non-dimensional parameter, which we adjust to avoid numerical instabilities arising if chemical diffusion is too low, while keeping it close to unity. For model M1, η needs to be increased to compensate for the smaller shear compared to model M2. The values of η are given in Table C.1. Furthermore, we take the angular average of the term between ⟨ ⋅ ⟩θ in Eq. (2) and assume a fixed horizontal diffusion coefficient Dh = 105 cm2 s−1. In the ESTER evolution models, there is no ad hoc enhanced mixing at the core boundary (overshooting) as is typically applied in 1D evolution models to account for the discrepancy between predicted and observationally inferred core masses.

We successfully ran 2D evolution models for a 12 M star: once for an initial rotation (at the equator) of 15% the critical angular velocity Ωc (model M1, veq = 112 km s−1) and once for 50% (model M2, veq = 358 km s−1)2. The evolution models start from a steady-state model at the ZAMS with a uniform chemical composition with X = 0.71 and Z = 0.012 (no assumption of solid-body rotation). The evolution in the Hertzsprung-Russell diagram (HRD) of the models computed with ESTER is shown in Fig. 1 for the effective temperature at the pole and at the equator. Figure 2 shows the angular velocity distribution at 1 Myr and 15.75 Myr for model M1. This stellar model starts with equatorial regions rotating more rapidly than polar ones, namely solar-like, and gradually evolve towards a more shellular rotation with a slightly anti-solar surface rotation. The same behaviour is also seen in M2 (see Appendix B). Contrary to an evolution modelled with a series of steady-state models, where one decreases the hydrogen mass fraction in the core (Gagnier et al. 2019b), our model evolution shows that Ω/Ωc decreases as the star evolves. We understand this behaviour to be a consequence of the slow redistribution of AM through baroclinic modes, which are damped on a timescale similar to the nuclear one. The timescale on which baroclinic modes are damped is given by

(3)

thumbnail Fig. 1.

HRD showing the evolution tracks of the 2D models for a 12 M star (Z = 0.012, Xi = 0.71) for (Ωeqc)i = 0.15 (dashed-dotted), and 0.5 (dashed). The black lines represent the properties at the pole, the grey lines the properties at the equator. The tracks are terminated slightly before the TAMS when the solver can no longer converge, which occurs at Xc/Xini = 0.141 and 0.131 for the respective aforementioned rotation rates. The dots mark the model of HD 192575 that is discussed in Sect. 3.

thumbnail Fig. 2.

Map of the angular velocity as a fraction of the critical angular velocity for a model of 1 Myr (top panel) and 15.75 Myr (bottom panel). The plots show cuts in the meridian plane. These plots are for a 12 M star with (Ωeqc)i = 0.15.

following the work of Busse (1981). The τbaro is usually of the order of the Eddington–Sweet timescale (Rieutord et al. 2006b). In its expression, denv is the thickness of the radiative envelope in the polar direction. The nuclear evolution timescale, τevol = Xc/c, is found to be roughly 3 and 30 times larger at the start of the MS for M1 and M2, respectively. Hence, and especially for M2, baroclinic modes that may be excited by initial conditions can be damped during nuclear evolution. We can therefore assume that initial conditions are of little importance and are forgotten during the first part of the MS. Yet, at the end of the MS, the ratio of timescales is reversed: the damping of baroclinic modes happens on a timescale that is 100 (M1) and 10 (M2) times longer than the nuclear one. This implies that the dynamical evolution, and especially the rotation rates, cannot be computed as a succession of stationary states as in Gagnier et al. (2019b). Table C.1 lists the values of these time scales for our two models. We note that the growth of the stellar radius is the main contributor to the growth of τbaro with MS evolution.

Figure 3 shows the evolution of the vertical diffusion coefficient as described by Eq. (2) for model M1. At the start of the MS, the latitudinally average angular velocity decreases from the core towards the surface until roughly 60% of the fractional radius, after which increases towards the surface (see right panel of Fig. 4), creating a layer where the shear is weak and thus Dv reaches a local minimum. Further along the MS evolution, keeps decreasing all the way towards the surface. Figure 3 also shows the difference between taking an average value for N2 as we do and using the full profile. This latter case shows the rapid variations of D that occur near the core and give rise to numerical problems as mentioned above.

thumbnail Fig. 3.

Profile of the vertical diffusion coefficient according to Eq. (2) for a 12 M model with an initial rotation frequency of 15% the critical angular velocity. The grey dotted lines show at a few ages (from left to right: Xc/Xini = 0.14, 0.30, 0.58, 0.82, 0.98), that is, if we were to take the local value of the Brunt-Väisälä frequency instead of the volume averaged value.

3. Calibration to HD 192575

To calibrate our rotating models, we use the results of Burssens et al. (2023, hereafter B23), who performed an asteroseismic modelling of the β Cephei pulsator HD 192575 based on 1D non-rotating stellar evolution models. This star is a unique case with which to test the theory of angular momentum transport as it has a relatively precise age estimate, and an inferred rotation profile from the rotational splittings of the observed multiplets. The stellar parameters of HD 192575 derived by B23 are summarised in the middle column of Table 1. We computed ESTER evolution tracks with a mass (M), metallicity (Z), and initial hydrogen-mass fraction (X) fixed to the values found by B23. We then picked the model for which, after evolution, the hydrogen-mass fraction in the convective core is closest to that inferred by these authors. As the shape of the internal rotation profiles in massive stars is currently unknown, these authors modelled the rotational splittings by assuming a rotation profile that is described as

(4)

Table 1.

Comparison of ESTER with the observations.

where ΔΩ = Ωcore − Ωsurf, rcore is the radius of the convective core, and rshear is the outer radius of the shear zone where the rotation frequency decreases linearly until it reaches the surface value. The values of Ωcore and Ωsurf are then optimised to best reproduce the observations for a given assumption for rshear. However, in the ESTER models, the rotation profile is computed self-consistently. Figure 4 shows the predicted evolution of the rotation profile for model M1. As can be seen from this figure (left), the latitudinal differential rotation changes from solar-like to anti-solar. This change may be interpreted as follows. The solar-like latitudinal profile is indeed the relaxed steady baroclinic state of a rotating radiative envelope (Espinosa Lara & Rieutord 2013). As time evolution proceeds, the core shrinks and spins up due to angular momentum conservation. Thanks to the Taylor–Proudman theorem, which states that the velocity field cannot vary in the direction of the rotation axis3, we understand that polar regions tend to follow the core and therefore rotate more rapidly than equatorial regions. This is a consequence of the rapid nuclear evolution, which prevents the star from relaxing to a quasi-steady state.

thumbnail Fig. 4.

Evolution of the rotation profile of model M1. Left panel: predicted surface rotation as a function of colatitude throughout the evolution (indicated by colour) normalised by the rotation frequency at the pole. Right panel: profile of the angular velocity (averaged over θ) throughout the evolution. The thick coloured line indicates the location of the convective core boundary. The black lines in both panels correspond to the best-matching model for HD 192575 (last column in Table 1).

Figure 4 (right panel) also shows the ‘radial’ rotation profile averaged over the latitude. This latter is similar in appearance to the linear piece-wise profile given in Eq. (4) when rshear is taken equal to the outer radius of the region with a non-zero gradient in the mean molecular weight (μ). This behaviour appears to be independent of the initial rotation, because model M2 shows the same profile. Therefore, we compare our results with the core and surface rotation frequencies derived with this assumption. The value of Ωcore in the ESTER models is defined as the value of , where an average is taken over all points in latitude. The last column in Table 1 shows the parameters of the ESTER M1 – model. The predicted core- and surface rotation of this model are consistent with the observationally derived values for HD 192575. Moreover, this model is able to reproduce, at the inferred age, the core mass, the asteroseismic radius from the 1D modelling, the astrometric luminosity from Gaia (Gaia Collaboration 2023), and the spectroscopically derived effective temperature within the uncertainties quoted in B23. It should be noted that with the implemented rotational mixing, the age is consistent with the age predicted from 1D (MESA) models with core-boundary mixing, although B23 were not able to precisely constrain its efficiency. In summary, the 2D evolution model M1 is able to explain the measured core and surface rotation of HD 192575 at the asteroseismically inferred age.

Additionally, we computed 1D SSE models with MESA (r22.11.1; Paxton et al. 2011, 2013, 2015, 2018, 2019; Jermyn et al. 2023) using the same physics as B23, except that we also include rotation in our models (more details in Appendix D). The rotation profiles predicted by the 1D MESA models are similar to taking rshear = R (Eq. (4)). Therefore, we compare the MESA profiles with the core and surface rotation frequencies derived by B23 with this assumption. When the same uniform viscosity is assumed as the one used in the ESTER models, we find that the 1D MESA model for (Ωeqc)i = 0.15 does not match the core and surface rotation rates. To reproduce the measured core and surface rotation frequencies at the current age with the 1D MESA model, the viscosity has to be increased to 109 cm2 s−1 compared to the value of 107 cm2 s−1 used in the ESTER models (vertical and horizontal components). While both the 1D MESA and 2D ESTER models can reproduce the core and surface frequencies of HD 192575, the shape of the rotation profile is significantly different between the two (see Fig. 5). Future studies on asteroseismic rotation inversion might be able to rule out one of these profiles (Vanlaer et al. 2023).

thumbnail Fig. 5.

Predicted rotation profile of HD 192575 with ESTER (red solid line) and MESA (black lines). The profiles indicated by a solid line and dotted line are for a constant viscosity of 107 cm2 s−1, while the dashed line is for 109 cm2 s−1. The dots correspond to the location of the core boundary. The red and grey shaded areas indicate the measured core (top one) and surface (bottom one) rotation frequencies by Burssens et al. (2023), assuming rshear equal to the outer boundary of the μ-gradient zone, and rshear = R, respectively.

4. Nitrogen enhancement at the surface

The observed abundance of 14N at the surface is often used to probe the efficiency of (rotational) chemical mixing in massive stars (e.g., Brott et al. 2011). The nuclear reactions in the ESTER models are described by analytical formulae for the energy generation rates of the pp-chain and CNO cycle (Kippenhahn & Weigert 1990; Rieutord et al. 2016). To track the abundance of 14N as a result of hydrogen-burning via the CNO cycle, we make the following assumptions. First, at the start of the evolution, all 12C present in the core is converted to 14N. The reaction that controls the rate at which 14N is generated is therefore the proton capture of 16O to create 17F, which then rapidly decays into 17O, itself reacting on protons to yield 14N and 4He. Therefore, the evolution of the mass fraction of 14N due to nuclear reactions is described by

(5)

where n(⋅) denotes the number densities, m the atomic mass, and λ(16O→17F) is the Maxwellian-average reaction rate ⟨σV⟩ taken from Angulo et al. (1999),

(6)

where 𝒩A is Avogadro’s number. With this simple modelling, we can reproduce the evolution of 14N core abundance as calculated through a more realistic network of nuclear reactions, such as in the MESA code.

In 1D stellar evolution codes, the effect of chemical transport via meridional circulation is typically treated in a diffusive way by adding an extra term to the vertical diffusion coefficient. This additional term takes on the form Deff = |ruv|2/(30Dh) (Chaboyer & Zahn 1992), where uv is the vertical component of the meridional flow, and Dh the horizontal diffusion coefficient (see Eq. (18) in Palacios et al. 2003).

In the present work, we tested the contribution of the meridional circulation to the transport of chemicals using 2D models. Typical values for the maximum of uv obtained from the 2D models range from 10−4 to 3 × 10−3 cm s−1, which turn out to be too small to have any noticeable influence on element transport. As a result, surface nitrogen abundance for the ESTER model of HD 192575 is weak: Δ[N/H]< 0.05 dex compared to the initial value.

In Appendix E, we show the stream function of the meridional velocity field for a star at the start of the MS, and for a star near the end. As the rotation profile gradually evolves to a more shellular configuration, the angular dependency of the meridional circulation is mostly described by a spherical harmonic of degree  = 2. Therefore, the 2D models show that the assumption that higher order spherical harmonics can be neglected in the expansion of the meridional velocity field is not true in stars at the beginning of the MS, even at the low rotation rate of model M1.

5. Conclusions

In this Letter, we present the first results of a 2D modelling of rotating stars including time evolution. These new models are an update of the 2D ESTER models designed by Espinosa Lara & Rieutord (2013), which compute steady 2D models of fast rotating stars. As in the steady models, the time-dependent ones take into account the centrifugal flattening of the star as well as the large-scale flows (differential rotation and meridional circulation) driven by the baroclinicity of the star.

We ran two 12 M models of a massive star with initial angular velocities of 15 and 50% of the critical one, respectively. The first model has a rather mild rotation rate but can be compared to the recent observation of a massive star, while the second model allowed us to test the performance of the code and revealed some new features of the internal dynamics of a massive star (see below).

The first model was actually designed to reproduce the observations of the β Cephei pulsator HD 192575 as derived by Burssens et al. (2023). We found that our model gives a luminosity, effective temperature, and core mass that are consistent with the observationally derived values. Moreover, the rotation profile derived from the ESTER model is also in accordance with the measurements of the core and envelope rotation rates of Burssens et al. (2023). We note that these authors used a simplified rotation profile (e.g., Eq. (4)) to derive the rotation rates. Fortunately, this profile turns out to be similar to the actual one predicted by the models, allowing a consistent comparison.

Another result provided by the above 2D models is the weakness of the meridional circulation. In the dynamics of baroclinic flows, this is controlled by viscosity to ensure the balance of angular momentum flux. The transport of chemicals by this flow seems to be negligible, but this needs to be confirmed by a more detailed analysis, because our models do not include the jump in viscosity expected at the core–envelope boundary, which is expected to drive a Stewartson layer along the tangential cylinder of the core (Rieutord 2006a; Gagnier & Rieutord 2020).

Our 2D models raise new questions as to the dynamics of rotating stars. In particular, the possibility of using a succession of steady-state models to monitor the rotational evolution of early-type stars, as done in Gagnier et al. (2019b), is now questionable and requires fresh investigation. Our results indeed show that for the 12 M model we computed, nuclear evolution is only slow enough to relax the star to a quasi-steady state at the beginning of the MS. When the star nears the end of the MS, the nuclear evolution becomes faster than the damping time of baroclinic modes (Busse 1981). The question then arises as to when a succession of steady models is liable to represent the evolution of a rotating star. This is a complex question related to the unsolved general problem of angular momentum transport in stars known since the first measurements of differential rotation in red giant stars (e.g., Beck et al. 2012; Deheuvels et al. 2012, 2015; Mosser et al. 2012). New investigations with the present 2D time-dependent ESTER models will be presented in forthcoming articles.


2

The models are available on https://doi.org/10.5281/zenodo.8228904.

3

This is exactly true for a steady solution of an incompressible inviscid rotating fluid when the Coriolis term dominates all other terms (Rieutord 2015), but it can be extended to fluids of varying density if momentum ρv is used instead of v.

4

am_D_mix_factor in MESA.

Acknowledgments

The authors are grateful to the referee prof. Georges Meynet for his comments and suggestions. The research leading to these results has received funding from the French Agence Nationale de la Recherche (ANR), under grant MASSIF (ANR-21-CE31-0018-02). The authors thank Siemen Burssens for providing the MESA setup. Computations of ESTER 2D-models have been possible thanks to HPC resources from CALMIP supercomputing center (Grant 2023-P0107).

References

  1. Angulo, C., Arnould, M., Rayet, M., et al. 1999, Nucl. Phys. A, 656, 3 [Google Scholar]
  2. Beck, P. G., Montalban, J., Kallinger, T., et al. 2012, Nature, 481, 55 [Google Scholar]
  3. Bouchaud, K., Domiciano de Souza, A., Rieutord, M., Reese, D. R., & Kervella, P. 2020, A&A, 633, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Brott, I., Evans, C. J., Hunter, I., et al. 2011, A&A, 530, A116 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Burssens, S., Bowman, D. M., Michielsen, M., et al. 2023, Nat. Astron., 3, 760 [NASA ADS] [Google Scholar]
  6. Busse, F. 1981, Geophys. Astrophys. Fluid Dyn., 17, 215 [Google Scholar]
  7. Chaboyer, B., & Zahn, J. P. 1992, A&A, 253, 173 [NASA ADS] [Google Scholar]
  8. Deheuvels, S., García, R. A., Chaplin, W. J., et al. 2012, ApJ, 756, 19 [Google Scholar]
  9. Deheuvels, S., Ballot, J., Beck, P. G., et al. 2015, A&A, 580, A96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Ekström, S., Georgy, C., Eggenberger, P., et al. 2012, A&A, 537, A146 [Google Scholar]
  11. Espinosa Lara, F., & Rieutord, M. 2013, A&A, 552, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Gagnier, D., & Rieutord, M. 2020, J. Fluid Mech., 904, A35 [NASA ADS] [CrossRef] [Google Scholar]
  13. Gagnier, D., Rieutord, M., Charbonnel, C., Putigny, B., & Espinosa Lara, F. 2019a, A&A, 625, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  14. Gagnier, D., Rieutord, M., Charbonnel, C., Putigny, B., & Espinosa Lara, F. 2019b, A&A, 625, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Gaia Collaboration (Vallenari, A., et al.) 2023, A&A, 674, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. Heger, A., Langer, N., & Woosley, S. E. 2000, ApJ, 528, 368 [NASA ADS] [CrossRef] [Google Scholar]
  17. Howarth, I. D., Bailey, J., Cotton, D. V., & Kedziora-Chudczer, L. 2023, MNRAS, 520, 1193 [NASA ADS] [CrossRef] [Google Scholar]
  18. Jermyn, A. S., Bauer, E. B., Schwab, J., et al. 2023, ApJS, 265, 15 [NASA ADS] [CrossRef] [Google Scholar]
  19. Kippenhahn, R., & Weigert, A. 1990, Stellar Structure and Evolution (New York: Springer) [Google Scholar]
  20. Mathis, S., Palacios, A., & Zahn, J. P. 2004, A&A, 425, 243 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Mombarg, J. S. G. 2023, A&A, in press, https://doi.org/10.1051/0004-6361/202345956 [Google Scholar]
  22. Mombarg, J. S. G., Dotter, A., Rieutord, M., et al. 2022, ApJ, 925, 154 [NASA ADS] [CrossRef] [Google Scholar]
  23. Mosser, B., Goupil, M. J., Belkacem, K., et al. 2012, A&A, 548, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Palacios, A., Talon, S., Charbonnel, C., & Forestini, M. 2003, A&A, 399, 603 [CrossRef] [EDP Sciences] [Google Scholar]
  25. Paxton, B., Bildsten, L., Dotter, A., et al. 2011, ApJS, 192, 3 [Google Scholar]
  26. Paxton, B., Cantiello, M., Arras, P., et al. 2013, ApJS, 208, 4 [Google Scholar]
  27. Paxton, B., Marchant, P., Schwab, J., et al. 2015, ApJS, 220, 15 [Google Scholar]
  28. Paxton, B., Schwab, J., Bauer, E. B., et al. 2018, ApJS, 234, 34 [NASA ADS] [CrossRef] [Google Scholar]
  29. Paxton, B., Smolec, R., Schwab, J., et al. 2019, ApJS, 243, 10 [Google Scholar]
  30. Rieutord, M. 2006a, A&A, 451, 1025 [CrossRef] [EDP Sciences] [Google Scholar]
  31. Rieutord, M. 2006b, in Stellar Fluid Dynamics and Numerical Simulations: From the Sun to Neutron Stars, eds. M. Rieutord, & B. Dubrulle, EAS Publ., 21, 275 [NASA ADS] [Google Scholar]
  32. Rieutord, M. 2015, Fluid Dynamics: An Introduction (New York: Springer), 508 [Google Scholar]
  33. Rieutord, M., Espinosa Lara, F., & Putigny, B. 2016, J. Comput. Phys., 318, 277 [NASA ADS] [CrossRef] [Google Scholar]
  34. Rogers, T. M., & McElwaine, J. N. 2017, ApJ, 848, L1 [CrossRef] [Google Scholar]
  35. Rogers, F. J., Swenson, F. J., & Iglesias, C. A. 1996, ApJ, 456, 902 [Google Scholar]
  36. Talon, S., & Charbonnel, C. 1998, A&A, 335, 959 [NASA ADS] [Google Scholar]
  37. Vanlaer, V., Aerts, C., Bellinger, E. P., & Christensen-Dalsgaard, J. 2023, A&A, 675, A17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Varghese, A., Ratnasingam, R. P., Vanon, R., Edelmann, P. V. F., & Rogers, T. M. 2023, ApJ, 942, 53 [NASA ADS] [CrossRef] [Google Scholar]
  39. Zahn, J. P. 1992, A&A, 265, 115 [NASA ADS] [Google Scholar]

Appendix A: Stellar equations

In addition to Eq. (1), the ESTER code solves the following set of fundamental stellar equations.

  • Mass conservation

    (A.1)

    where ρ is the density and v the meridional velocity.

  • Meridional momentum equation

    (A.2)

    where P is the pressure, ϕ the gravitational potential, s the radial distance to the rotation axis, Ω the local angular velocity, and the meridional components of the viscous force. In this equation, the time derivative of the meridional circulation has been neglected because of its extremely small value compared to other terms.

  • Angular momentum equation

    (A.3)

    where ν is the kinematic viscosity. The kinematic viscosities in both the horizontal and vertical directions are free parameters in ESTER. For both, we take 107 cm2 s−1, as this was found to be the order of magnitude needed to explain the rotation profiles of F-type stars by Mombarg (2023).

  • Entropy equation

    (A.4)

    where S is the entropy, χ the thermal conductivity, and ε the energy generation rate per unit mass.

  • Nuclear energy generation is computed with a simple law given by Kippenhahn & Weigert (1990), namely

    (A.5)

    from which we deduce the hydrogen consumption nuc. Here, T9 = T/109 K, A = 15.228, and C(T9) is a correction term.

Appendix B: Rotation profiles

In this Appendix, we also show the rotation profiles of the M2 model ((Ωeqc)i = 0.5).

thumbnail Fig. B.1.

Profiles of the angular velocity as a function of the radial coordinate in the polar direction (top panel), equatorial direction (middle panel), and averaged over latitude (bottom panel). The solid lines correspond to model M1, the dashed lines to model M2. The locations of the outer edge of the convective core are shown as vertical lines (dotted lines for model M1, dashed-dotted for model M2).

Appendix C: Timescales

In this Appendix, we show the nuclear and baroclinic timescales (see Sect. 2), once close to the ZAMS, and once close to the TAMS.

Table C.1.

Diffusivity factor and characteristic timescales of the models.

Appendix D: MESA setup

In this Appendix, we provide a short summary of the physics used for the 1D MESA model of HD192575 discussed in Section 3. The transport of AM in MESA is treated as a diffusive process (Eq. (B4) in Paxton et al. (2013)) and shellular rotation is imposed. As our aim here is to compare the 2D ESTER models with the 1D physics used by Burssens et al. (2023), we use the same description for the chemical mixing as these authors. This description is based on predictions of simulations of internal gravity waves (e.g. Rogers & McElwaine 2017; Varghese et al. 2023), where the chemical diffusion coefficient takes the form of

(D.1)

where D0 is a free parameter we set to 103 cm2 s−1 (same as Burssens et al. (2023)), and ρ0 the density at the core boundary. This means that we set the factor fC —which accounts for the different efficiencies between the transport of AM and chemical elements4 (Heger et al. 2000)— equal to zero.

Appendix E: Meridional circulation

In this Appendix, we show the stream lines of the meridional flow for models M1 ((Ωeqc)i = 0.15) and model M2 ((Ωeqc)i = 0.5).

thumbnail Fig. E.1.

Isocontours of the stream function for model M1 of 1 Myr (top panel) and 15.75 Myr (bottom panel). The different colours indicate counter-rotating cells, where the red cells rotate clockwise in the first quadrant of the plot. These models correspond to those shown in Fig. 2. The isobars at the edge of the convective core and at the surface are shown in grey.

thumbnail Fig. E.2.

Same as Fig E.1, but for model M2.

All Tables

Table 1.

Comparison of ESTER with the observations.

Table C.1.

Diffusivity factor and characteristic timescales of the models.

All Figures

thumbnail Fig. 1.

HRD showing the evolution tracks of the 2D models for a 12 M star (Z = 0.012, Xi = 0.71) for (Ωeqc)i = 0.15 (dashed-dotted), and 0.5 (dashed). The black lines represent the properties at the pole, the grey lines the properties at the equator. The tracks are terminated slightly before the TAMS when the solver can no longer converge, which occurs at Xc/Xini = 0.141 and 0.131 for the respective aforementioned rotation rates. The dots mark the model of HD 192575 that is discussed in Sect. 3.

In the text
thumbnail Fig. 2.

Map of the angular velocity as a fraction of the critical angular velocity for a model of 1 Myr (top panel) and 15.75 Myr (bottom panel). The plots show cuts in the meridian plane. These plots are for a 12 M star with (Ωeqc)i = 0.15.

In the text
thumbnail Fig. 3.

Profile of the vertical diffusion coefficient according to Eq. (2) for a 12 M model with an initial rotation frequency of 15% the critical angular velocity. The grey dotted lines show at a few ages (from left to right: Xc/Xini = 0.14, 0.30, 0.58, 0.82, 0.98), that is, if we were to take the local value of the Brunt-Väisälä frequency instead of the volume averaged value.

In the text
thumbnail Fig. 4.

Evolution of the rotation profile of model M1. Left panel: predicted surface rotation as a function of colatitude throughout the evolution (indicated by colour) normalised by the rotation frequency at the pole. Right panel: profile of the angular velocity (averaged over θ) throughout the evolution. The thick coloured line indicates the location of the convective core boundary. The black lines in both panels correspond to the best-matching model for HD 192575 (last column in Table 1).

In the text
thumbnail Fig. 5.

Predicted rotation profile of HD 192575 with ESTER (red solid line) and MESA (black lines). The profiles indicated by a solid line and dotted line are for a constant viscosity of 107 cm2 s−1, while the dashed line is for 109 cm2 s−1. The dots correspond to the location of the core boundary. The red and grey shaded areas indicate the measured core (top one) and surface (bottom one) rotation frequencies by Burssens et al. (2023), assuming rshear equal to the outer boundary of the μ-gradient zone, and rshear = R, respectively.

In the text
thumbnail Fig. B.1.

Profiles of the angular velocity as a function of the radial coordinate in the polar direction (top panel), equatorial direction (middle panel), and averaged over latitude (bottom panel). The solid lines correspond to model M1, the dashed lines to model M2. The locations of the outer edge of the convective core are shown as vertical lines (dotted lines for model M1, dashed-dotted for model M2).

In the text
thumbnail Fig. E.1.

Isocontours of the stream function for model M1 of 1 Myr (top panel) and 15.75 Myr (bottom panel). The different colours indicate counter-rotating cells, where the red cells rotate clockwise in the first quadrant of the plot. These models correspond to those shown in Fig. 2. The isobars at the edge of the convective core and at the surface are shown in grey.

In the text
thumbnail Fig. E.2.

Same as Fig E.1, but for model M2.

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.