© The Authors 2026

1. Introduction

Stellar models are crucial in the interpretation and analysis of astrophysical data (e.g., Brott et al. 2011a; Bressan et al. 2012; Ekström et al. 2012; Choi et al. 2016; Limongi & Chieffi 2018; Amard et al. 2019; Costa et al. 2025). One of the largest sources of uncertainty in their predictions comes from the internal transport mechanisms (e.g., Buldgen 2019; Fuller et al. 2019; Eggenberger 2024). Asteroseismology has proven successful in constraining the angular moment transport (AMT) efficiency in low- to intermediate-mass stars (e.g., Deheuvels et al. 2014; Van Reeth et al. 2016; Pápics et al. 2017; Aerts et al. 2019). These observations generally indicate that an efficient AMT is at play in radiative zones, which cannot be reproduced solely by hydrodynamic instabilities (e.g., Eggenberger et al. 2012, 2022; Fuller et al. 2019; Moyano et al. 2023). In contrast, early-B and O stars seem to show a diversity of internal rotation profiles, with asteroseismic observations of β-Cephei stars revealing both radial differential rotation and uniform rotation (see e.g., Aerts et al. 2003; Dupret et al. 2004; Briquet et al. 2007; Dziembowski & Pamyatnykh 2008; Marques et al. 2013; Ouazzani et al. 2020; Salmon et al. 2022; Burssens et al. 2023; Fritzewski et al. 2025).

It has long been identified that main-sequence (MS) massive stars often show observational evidence of carbon-nitrogen-oxygen (CNO) cycle signatures, suggesting that efficient mixing operates in their radiative envelope (e.g., Gies & Lambert 1992; Kilian 1992; Morel et al. 2008; Hunter et al. 2009; Brott et al. 2011b; Grin et al. 2017; Mahy et al. 2020). Various physical processes have been proposed to explain these signatures, notably rotational mixing (RM; e.g., Maeder & Meynet 2000; Brott et al. 2011b), internal gravity waves (IGWs; e.g., Rogers & McElwaine 2017; Mombarg et al. 2025), and binary interactions (e.g., de Mink et al. 2009; Langer 2012).

The majority of massive stars possess companions; many of these systems are expected to interact at some stage of their evolution (e.g., Sana et al. 2012; Moe & Di Stefano 2017; Offner et al. 2023; Bodensteiner et al. 2025; Patrick et al. 2025). de Mink et al. (2009) proposed using pre-interacting close binary systems to test RM. In these close configurations, tides are expected to synchronize the stars to the orbital motion, which should boost RM. This theoretical prediction was confirmed by the studies of Song et al. (2013, 2016), who found that tidal interactions always boost RM. However, observational studies challenge this picture: the majority of the short-period systems in the samples of Martins et al. (2017), Pavlovski et al. (2018, 2023), and Abdul-Masih et al. (2019, 2021) do not show any evidence of strong nitrogen enrichment. Moreover, they do not suggest a clear anticorrelation between period and nitrogen enrichment – a direct prediction of the models by de Mink et al. (2009).

For this study, we investigated the interplay between tidal torques and RM, with a particular focus on the impact of the AMT assumptions. To this aim, we computed several single and binary star grids with the GENeva Evolutionary Code (GENEC, Eggenberger et al. 2008) with different AMT treatments: purely hydrodynamic models with an advective-diffusive AMT and magneto-hydrodynamic models accounting for the calibrated Tayler-Spruit instabilities (Eggenberger et al. 2022). Tides are often seen as an additional, independent source of chemical mixing in radiative envelopes. However, in the way they are accounted for in most state-of-the-art binary codes (for instance in MESA, Paxton et al. 2015, 2018, 2019; Marchant et al. 2016, BONN, de Mink et al. 2009; Yoon et al. 2010; Szécsi et al. 2022, or GENEC, Song et al. 2013, 2016), they only alter mixing by modifying the stellar rotation, and thereby the RM efficiency (but this picture may be too simplistic; see discussion in Sect. 4.3). We aim to demonstrate that in this framework, model reactions to tidal torques strongly depend on the AMT treatments and that in certain configurations, stars in binary systems may experience less mixing –“tidally suppressed mixing”– than single stars with identical initial conditions. The approach and conclusions of the present paper are reminiscent of Zahn (1994), who demonstrated that tidal synchronization in late-type binaries reduces RM, thereby providing a theoretical explanation for their reduced lithium depletion relative to single-star counterparts.

The structure of this paper is as follows: in Sect. 2, we provide a description of our methods, and the single and binary star physics assumptions we made. In Sect. 3, we present our binary models with tidal mixing. In Sect. 4, we discuss the limitations of our results and how they compare to those of previous studies and observations. Finally, we summarize our findings in Sect. 5.

2. Single and binary star physics

2.1. Stellar physics ingredients

In this study we use GENEC to compute detailed stellar evolution 1D models. The adopted stellar physics is generally the same as in Ekström et al. (2012), although with some differences that are explicitly listed below. We recall here the main ingredients.

The adopted wind prescription is the standard combination used in the GENEC grids. It consists of the de Jager et al. (1988), Vink et al. (2001), and Nugis & Lamers (2000) prescriptions lowered by a scaling factor 0.85.

During the MS, convection in the core is assumed to be very efficient, i.e., angular momentum transport and chemical mixing are assumed to be instantaneous in this region. Small convective zones near the stellar surface are in contrast treated following the mixing length theory (MLT) framework, with α_MLT = 1.6. The size of the convective core is increased following a step overshoot treatment with α_ov = 0.1.

We computed models with two different AMT treatments. Purely hydrodynamic models¹ (hereafter hydro models) account for shear instabilities (Zahn 1992) and meridional circulation (Eddington 1925; Sweet 1950; Zahn 1992; Maeder & Zahn 1998). Magneto-hydrodynamic models (hereafter magnetic models) in addition account for the Tayler-Spruit instability (Tayler 1973; Spruit 2002; Maeder & Meynet 2003, 2004, 2005). For the magnetic models we used the asteroseismic-calibrated version of the Tayler-Spruit dynamo as formulated by Eggenberger et al. (2022).

In hydro models, the AMT by meridional circulation is treated as an advective process as in Ekström et al. (2012). In this case, the advective-diffusive equation reads as

$$\begin{array}{c}\hfill \rho \frac{\partial \left(r^2\mathrm{\Omega }\right)}{\partial t}=\frac{1}{5r^2}\frac{\partial }{\partial r}\left(\rho r^4\mathrm{\Omega }U(r)\right)+\frac{1}{r^2}\frac{\partial }{\partial r}\left(\rho r^4\mathrm{\Omega }{D}_{\mathrm{AM}}\frac{\partial \mathrm{\Omega }}{\partial r}\right),\end{array}$$(1)

where ρ is the mean isobar density, r the radius, U the radial dependence of the vertical component of the meridional circulation, Ω the average angular velocity on an isobar, and D_AM the diffusion coefficient for AMT. As demonstrated in Maeder & Meynet (2005), when the magnetic instabilities are accounted for, they are strong enough for achieving near solid-body rotation during the whole MS evolution, so that the AMT by meridional circulation can effectively be neglected and a fully diffusive treatment can be followed.

GENEC models also account for the RM imparted by the three mentioned processes. In the present work, shear instabilities were modeled using the expression of Maeder (1997),

$$\begin{array}{c}\hfill D_{\mathrm{shear}}=f_{\mathrm{energ}}\frac{H_P}{g\delta }\frac{K}{\frac{\varphi }{\delta }{\mathrm{\nabla }}_{\mu }+{\mathrm{\nabla }}_{\mathrm{ad}}-{\mathrm{\nabla }}_{\mathrm{rad}}}{\left(\frac{9\pi }{32}\mathrm{\Omega }\frac{\text{d}ln\mathrm{\Omega }}{\text{d}lnr}\right)}^2,\end{array}$$(2)

where f_energ is the fraction of the excess energy in the shear that contributes to mixing, H_P the pressure scale height, g the local gravity, δ and ϕ thermodynamic derivatives, K the thermal diffusivity, ∇_μ, ∇_ad and ∇_rad respectively the mean molecular weight, adiabatic and radiative gradients. It is worth noting that the shear diffusion coefficient depends on both Ω and its gradient.

As demonstrated by Chaboyer & Zahn (1992), the combined effect of the meridional circulation and the horizontal turbulence mitigates the advection of the chemical elements through homogenization of horizontal layers. The vertical RM then reduces to a diffusion process, which can be described by the effective diffusion coefficient,

$$\begin{array}{c}\hfill D_{\mathrm{eff}}=\frac{1}{30}\frac{{|rU(r)|}^2}{D_{\mathrm{h}}}·\end{array}$$(3)

The horizontal turbulence was modeled using the prescription of Zahn (1992):

$$\begin{array}{c}\hfill D_{\mathrm{h}}=\frac{1}{c_{\mathrm{h}}}r|2V(r)-\alpha U(r)|,\end{array}$$(4)

where c_h is a constant on the order of 1, V the radial dependence of the horizontal component of the meridional circulation, and $\alpha =\frac{1}{2}\frac{\text{d}ln(r^2\overline{\mathrm{\Omega }})}{\text{d}lnr}$. The transport equation for the chemicals therefore reads:

$$\begin{array}{c}\hfill \rho \frac{\partial X_i}{\partial t}=\frac{1}{r^2}\frac{\partial }{\partial r}\left(\rho r^2{D}_{\mathrm{chem}}\frac{\partial X_i}{\partial r}\right),\end{array}$$(5)

with X_i the considered chemical species, D_chem = D_AM + D_eff. We note that the diffusive treatment only holds for RM. In magnetic models where the AMT is highly efficient, the contribution of meridional circulation to the AMT can be neglected; however, its impact on RM remains significant. The shear diffusion coefficient is negligible for flat Ω-profiles (see Eq. (2)), as a result RM is in fact dominated by meridional circulation, i.e., D_chem ≈ D_eff (Maeder & Meynet 2005; Song et al. 2016; Nandal et al. 2024; Asatiani et al. 2025, see also Fig. 2). In this case, the RM efficiency is primarily determined by the angular momentum (AM) content of the stars (or equivalently their angular velocity).

Following this treatment of AMT and RM, the presented hydro models are computed with the same physics as in Ekström et al. (2012). These models were calibrated to reproduce the observed abundances of nitrogen of a population of Galactic MS B-type stars. We performed an equivalent calibration for the magnetic models. We found that using c_h = 0.70 in Eq. (4) (instead of c_h = 1 for the hydro models) offers a satisfactory reproduction of the observed abundances of the same stars. The details of the performed calibration are given in Appendix A. As a result of this calibration, the predicted RM of magnetic models is reduced compared to those presented in Song et al. (2016), which should alter the size of the predicted parameter space where chemically homogeneous evolution occurs in close binaries.

2.2. Binary interactions

We simulate binary systems in circular orbits. We effectively evolve single stars but account for binary interactions assuming a non-evolving point mass for the companion, as in Detmers et al. (2008), de Mink et al. (2009), Song et al. (2013, 2016, 2018), and Sciarini et al. (2024), for example. The simulations are stopped either when the star fills its Roche lobe (RL), which size is computed following Eggleton (1983), or at terminal-age-main-sequence (TAMS) if the RL has not been filled during the MS evolution.

We propose a refined treatment of tides in massive GENEC models, following the recent findings of Fragos et al. (2023) and Sciarini et al. (2024). The dynamical tides are treated as in Song et al. (2013) and Sciarini et al. (2024), following the Zahn (1977) prescription. In case of circular orbits, the tidal torque is:

$$\begin{array}{c}\hfill {\frac{\text{d}}{\text{d}t}\left(I{\mathrm{\Omega }}_{\mathrm{spin}}\right)|}_{\mathrm{Dyn}}=\frac{3}{2}\frac{G{M}^2}{R}{E}_2\left(q^2{\left(\frac{R}{a}\right)}^6\right)s_{22}^{8/3}\text{sgn}(s_{22}),\end{array}$$(6)

where I is the moment of inertia of the star, Ω_spin, M and R respectively its surface angular velocity, mass, and radius, G the gravitational constant, q the mass ratio, and a the semi-major axis of the orbit. For E₂, we use the prescription by Qin et al. (2018). Following the approach of Song et al. (2013, 2016, 2018), the tidal torques are applied to an outer layer comprising 3% of the total mass, assumed to rotate rigidly. This simplification provides a numerically stable way to deposit the torque in the envelope while allowing differential rotation in the interior. In reality, the dynamical tides torque results from gravity waves traveling trough the envelope (Zahn 1975, 1977; Goldreich & Nicholson 1989) and dissipating either gradually or when reaching critical layers (e.g., Talon & Kumar 1998; Mathis 2009). Our treatment should thus be viewed as an approximation of AM deposition in the outer radiative envelope, rather than a detailed physical model of wave damping.

It is worth highlighting that s₂₂ scales as Ω_orb − Ω_spin, with Ω_orb and Ω_spin the orbital and spin angular velocities. Hence, the dynamical tides torque scales as (Ω_orb − Ω_spin)^8/3. As demonstrated in Sciarini et al. (2024), in the Hurley et al. (2002) adaptation of the Zahn (1977) formalism the dynamical tides torque in contrast scales linearly with Ω_orb − Ω_spin, which alters the torque efficiency depending on the departure from synchronization. When the dynamical tides are treated consistently with the original formulation by Zahn (1977), they are very efficient when the star is far from synchronization, but become extremely inefficient close to synchronization (see Figs. A.1 and A.2 in Sciarini et al. 2024). As a result, under this formalism the stars may deviate from synchronization at Roche lobe overflow (RLOF).

Fragos et al. (2023) proposed that equilibrium tides acting on small convective zones near the surface can become the dominant contribution to the evolution toward synchronization already during the MS due the fast decrease of the term E₂. In the present work we modeled equilibrium tides acting on small subsurface convective zones following a treatment similar to that in Fragos et al. (2023). The treatment of equilibrium tides in Fragos et al. (2023) is inspired from Hurley et al. (2002), which is adapted from Hut (1981) and Rasio et al. (1996). Under their formalism, the ratio of the apsidal motion constant (AMC) to the tidal dissipation timescale k₂/T in the Hut (1981) formulas is replaced by:

$$\begin{array}{c}\hfill {\left(\frac{k_2}{T}\right)}_{\mathrm{conv}}=\frac{2}{21}\frac{f_{\mathrm{conv}}}{{\tau }_{\mathrm{conv}}}\frac{M_{\mathrm{conv}.\mathrm{reg}.}}{M},f_{\mathrm{conv}}=min(1,{\left(\frac{P_{\mathrm{tid}}}{2{\tau }_{\mathrm{conv}}}\right)}^2),\end{array}$$(7)

where τ_conv is the convective turnover timescale, M_conv.reg. the mass of the convective region, and M the total mass of the star. f_conv reduces the strength of the tides when the turnover timescale is greater than the pumping timescale P_tid = 2π|Ω_orb−Ω_spin|⁻¹ (see e.g., Zahn 1966; Goldreich & Nicholson 1977; Goodman & Oh 1997; Zahn 2008; Souchay et al. 2013). Equation (7) is obtained comparing the circularization timescales of Hut (1981) and Rasio et al. (1996). By doing so, the AMC dependence is lost.

In this work, we implemented a formalism closer to the original formulation by Hut (1981), keeping the AMC dependence. We include the term f_conv in our formalism to account for fast tides. Furthermore, as tides involve AM exchanges, we scale the strength of the tides by the ratio of the moment of inertia of the convective region I_conv.reg to that of the star I instead of the mass ratio. As a result, our equilibrium tides prescriptions for small subsurface convective regions can be written as:

$$\begin{array}{c}\hfill {\left(\frac{k_2}{T}\right)}_{\mathrm{conv}}\to \frac{k_2}{{\tau }_{\mathrm{conv}}}{f}_{\mathrm{conv}}\frac{I_{\mathrm{conv}.\mathrm{reg}.}}{I}·\end{array}$$(8)

The term 2/21 is not present in our prescription as it is not present in the original formalism by Hut (1981). The tidal torque in case of circular orbits is thus:

$$\begin{array}{c}\hfill {\frac{\text{d}}{\text{d}t}\left(I{\mathrm{\Omega }}_{\mathrm{spin}}\right)|}_{\mathrm{Eq}}=3\frac{k_2}{{\tau }_{\mathrm{conv}}}{f}_{\mathrm{conv}}\frac{I_{\mathrm{conv}.\mathrm{reg}.}}{I}M{R}^2{q}^2{\left(\frac{R}{a}\right)}^6({\mathrm{\Omega }}_{\mathrm{orb}}-{\mathrm{\Omega }}_{\mathrm{spin}}).\end{array}$$(9)

We compute k₂ consistently with the stellar structure solving the Clairaut-Radau equation:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill r\frac{\text{d}{\eta }_2(r)}{\text{d}r}& ={\eta }_2(r)-{\eta }_2^2(r)+6(1-\frac{\rho (r)}{\overline{\rho }(r)}({\eta }_2(r)+1)),\hfill \\ \hfill k_2& =\frac{3-{\eta }_2(R)}{4+2{\eta }_2(R)},\hfill \end{array}\end{array}$$(10)

where ρ and $\overline{\rho }$ are respectively the density and average density at distance r from the center and R is the radius at optical depth τ  =  2/3. Equation (10) is solved using a 4-th order Runge-Kutta integrator. The AMC is directly linked to the internal density profiles of the stars and serves as a powerful observable constraint, being related to the apsidal motion of close binary systems (e.g., Sterne 1939; Claret 2004, 2024; Claret & Giménez 2010; Rosu et al. 2020a,b, 2022a,b). This quantity is included among the outputs of the stellar grids provided with this paper (Sect. 3). The tracks are available on Zenodo, Yareta, and on the Geneva stellar group database.

It is worth noting that the Hut (1981) formalism relies on the simplifying assumption of constant time lag. As demonstrated by Zahn (2008) and Remus et al. (2012), this approximation breaks down in the fast tides regime (τ_conv > P_tid). In this case it is in principle needed to expand the tidal potential in its multiple Fourier components. However, as illustrated in Appendix B, we find that the dynamical (equilibrium) tides largely dominate when the stars are far from (close to) synchronization. As such, τ_conv < P_tid is always verified when equilibrium tides dominate.

Finally, in Hurley et al. (2002) and Fragos et al. (2023), the convective turnover timescale is obtained using an estimate adapted from Rasio et al. (1996). In Fragos et al. (2023), they used:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{conv}}={\left(\frac{M_{\mathrm{conv}.\mathrm{reg}.}}{3L}\frac{R_{\mathrm{t},\mathrm{conv}.\mathrm{reg}.}+R_{\mathrm{b},\mathrm{conv}.\mathrm{reg}.}}{2}(R_{\mathrm{t},\mathrm{conv}.\mathrm{reg}.}-R_{\mathrm{b},\mathrm{conv}.\mathrm{reg}.})\right)}^{1/3}.\end{array}$$(11)

We note that this timescale directly depends on the mass M_conv.reg. and radial extent R_{t, conv.reg.} − R_{b, conv.reg.} of the convective region. As such, it can get very small depending on the size of the considered convective zone. In this work, we obtained the value of τ_conv within the MLT framework (see Sect. 2.1) as:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{conv}}=\frac{l_{\mathrm{conv}}}{{\upsilon }_{\mathrm{conv}}},\end{array}$$(12)

where l_conv = α_MLTH_P and υ_conv are the local mixing length and velocity of the bubble in the MLT framework. In a convective region where the local value of τ_conv differs among several shells, we take the minimum value across the convective region. Models presented in this work typically have several subsurface convective regions. We account for the combined effect of all the regions by adding up the tidal torques. The evolution toward synchronization of an illustrative system –considering only dynamical tides or both dynamical and equilibrium tides, treated either following Fragos et al. (2023) or Eq. (8)– is shown in Appendix B.

We account for orbital evolution imparted by mass-loss and spin-orbit AM exchange according to:

$$\begin{array}{c}\hfill \frac{\dot{a}}{a}=2\frac{{\dot{J}}_{\mathrm{orb}}}{J_{\mathrm{orb}}}-2\frac{{\dot{M}}_1}{M_1}+\frac{{\dot{M}}_1}{M_1+M_2},\end{array}$$(13)

where M₁ and M₂ are the masses of the stars, a the semi-major axis of the orbit, and J_orb the orbital AM. Equation (13) is obtained following an AM balance approach under the assumption that Ṁ₂ = 0. ${\dot{J}}_{\mathrm{orb}}$ is expressed as:

$$\begin{array}{c}\hfill \begin{array}{cc}& {\dot{J}}_{\mathrm{orb}}={\dot{J}}_{\mathrm{tides}}+{\dot{J}}_{\mathrm{ml}},\hfill \\ \hfill & \text{where}{\dot{J}}_{\mathrm{tides}}=-({\frac{\text{d}}{\text{d}t}\left(I{\mathrm{\Omega }}_{\mathrm{spin}}\right)|}_{\mathrm{Dyn}}+{\frac{\text{d}}{\text{d}t}\left(I{\mathrm{\Omega }}_{\mathrm{spin}}\right)|}_{\mathrm{Eq}})\hfill \end{array}\end{array}$$(14)

and ${\dot{J}}_{\mathrm{ml}}$ is the AM carried away by mass loss assuming that the removed mass has the specific orbital AM of the star.

In case of twin (equal mass) systems, we assume the companion has the same mass loss and initial parameters so that its contribution to the orbital evolution is the same as that of the simulated star. In this case, the right hand side of Eq. (13) is simply multiplied by two.

3. Binary models with tidal mixing

We present several sets of binary models at solar metallicity (Z_⊙ = 0.014), incorporating the ingredients described in Sect. 2, to illustrate how tides can either enhance or suppress RM, and how closely this behavior is linked to the choice of AMT mechanism. To this end, we systematically compare the outputs of magnetic and hydro models, using single-star models as a reference.

3.1. Comparison with Song et al. (2013)

We start by reproducing the results of Song et al. (2013)² and comparing them with magnetic models. This allows us to illustrate in a simple way that magnetic and hydro models react differently to tidal torques. We provide illustrative examples of the evolution toward synchronization for both types of models in Appendix C.

3.1.1. Spin-down models

The spin down-models presented in Song et al. (2013) are primaries of initial mass M₁ = 15 M_⊙, with a companion of fixed mass M₂ = 10 M_⊙. The 15 M_⊙ models are initiated at 60% of the critical velocity, defined as in Ekström et al. (2012). The chosen periods are P = 1.1, 1.4, 1.6 and 1.8 days. Figure 1 shows the nitrogen enhancement, defined as:

Fig. 1. Nitrogen and helium enhancement evolution of the same spin-down systems as in Song et al. (2013). Left panel: Hydro models. Right panel: Magnetic models.

$$\begin{array}{c}\hfill \mathrm{\Delta }log(\text{N}/\text{H})=[\text{N}/\text{H}]-{[\text{N}/\text{H}]}_0,\end{array}$$(15)

and the surface helium mass fraction Y_surf evolution. Both the nitrogen and the helium enhancement are used as tracers of RM. Figure 1 allows to observe several key results:

1. In this spin-down configuration, tides increase mixing in hydro models and decrease mixing in magnetic models. Compared to single-star evolution, hydro models accounting for spin-down by tides are more enhanced in both nitrogen and helium (the picture is more complex for helium, see remark 4). In contrast, magnetic binary models show less enrichment in both nitrogen and helium than the single-star model.

2. The nitrogen enrichment of hydro models increases with increasing period, which can seem counterintuitive as a longer period corresponds to a lower synchronization angular velocity. Mixing is more efficient at lower synchronization angular velocities because they generate a higher degree of differential rotation (Ω-gradient), which increases shear mixing (see Eq. (2) and Fig. 2).

3. The nitrogen and helium enrichments of magnetic models decrease with increasing period, which can be explained by the fact that the strong AMT in magnetic models makes shear mixing inefficient. The Ω-profiles are flat, and hence mixing is not increased by the tidal torques, but decreases as the stars spin down (see Fig. 2).

4. The nitrogen and helium enhancement of hydro models do not follow exactly the same evolution. Nitrogen is efficiently enhanced at the beginning of the evolution of models with tides, while the simultaneous helium enrichment is less pronounced. Hydro models are subject to high degrees of differential rotation at the beginning of the evolution caused by the tidal torques, which lead to strong shear mixing. Nitrogen is more sensitive to this early shear as the nitrogen profile quickly changes due to CNO burning. Although helium is also impacted by CNO, the profiles are less steep than the nitrogen profiles and thus RM at the beginning of the evolution is less efficient in enhancing He. The single-star model ends up more enriched in He than the binary models with P = 1.4, 1.6 and 1.8 days as it has kept more AM and therefore He is more efficiently mixed later in the evolution, when the He profile becomes steeper (see also Sect. 3.1.3 and Appendix D).

Figure 2 shows the angular velocity and diffusion coefficient profiles at t = 500 kyrs of the hydro and magnetic models (the shear diffusion coefficient for the former, the effective diffusion coefficient for the latter). It clearly illustrates the points discussed above: magnetic models have flat rotation profiles and therefore mixing is dominated by meridional circulation (D_eff). Even when subject to strong tidal torques, their profiles remain flat and therefore they are not subject to shear mixing. Systems of longer periods are more braked and therefore have lower diffusion coefficient profiles. In hydro models, the AMT is less efficient and the stars rotate differentially. The degree of differential rotation increases with increasing period, and as a result D_shear is larger for longer periods.

Fig. 2. Upper panel:Ω profiles of the hydro and magnetic models at t = 500 kyrs of the same spin-down systems as in Song et al. (2013). Lower panel: Diffusion coefficient profiles.

3.1.2. Spin-up models

We here focus on spin-up models, using the same initial conditions as in Song et al. (2013), and compare the results with magnetic models. In this case, the initial velocity is υ_ini/υ_crit = 0.2, the periods P = 1, 1.1, 1.2 and 1.4 days, and the initial masses are the same as in Sect. 3.1.1. Figure 3 shows the nitrogen and helium enhancements evolution. In this case hydro and magnetic models both predict that tides increase mixing. This illustrates that the impact of tidal torques on RM is simpler in magnetic models: in spin-up cases, they enhance mixing, whereas in spin-down cases, they suppress mixing. The Ω-profile of magnetic models being flat, RM by shear is inefficient and mixing only depends on the value of Ω.

Fig. 3. Same as Fig. 1 but for spin-up models.

In the case of hydro models, the picture is more complex: the efficiency of RM depends on both the value and the gradient of Ω. Spin-up by tides momentarily decreases the differential rotation, but it also increases the value of Ω (see Fig. 12 in Song et al. 2013). The net effect in this case is an increase in mixing with increasing synchronization angular velocity (decreasing period). We show in Sect. 3.2 that with hydro models there also exists configurations where spun-up models are less enriched than their single-star counterparts, a result that may appear even more counterintuitive.

3.1.3. Spin-down massive star models

In order to test the effect of mass on the previous results, we compute spin-down models with the same initial velocity as the models of Sect. 3.1.1, but with larger initial masses (M₁ = 45 M_⊙ and M₂ = 30 M_⊙) and longer periods (P = 1.8, 3, 5 and 8 days). Figure 4 shows the nitrogen and helium enhancements evolution of the hydro and magnetic models. We note that although the same trends as in Fig. 1 are observed, the conclusions are different. The magnetic models still behave the same: spun-down stars are less enriched than the single star, and the longer the period (lower Ω_orb) the lower the enrichment.

Fig. 4. Same as Fig. 1 but for models with M1 = 45 M⊙, M2 = 30 M⊙, and P = 1.8, 3, 5, and 8 days.

Hydro models also show the same trend as in Sect. 3.1.1 (the longer the period, the higher the enrichment). However, while tidal spin-down initially boosts shear mixing –mostly enhancing nitrogen– after this phase of fast enrichment the nitrogen enhancement slows. At RLOF binary models are less enriched in both helium and nitrogen than the single-star model, for which the enrichment steadily increases throughout the evolution. This can be attributed to the fact that in the long term, mixing is quenched in the binary models as the tidal torques make them rotate at lower angular velocity and reduce the Ω–gradients (once the stars are synchronized, the tidal torque acts in the opposite direction, transferring AM from the orbit to the spin in order to maintain synchronization. By keeping the surface synchronized, this prevents the development of significant differential rotation).

The difference in enrichment between the single star and the binary systems is more pronounced for helium. For the single star, the helium enrichment accelerates at the end of the evolution, as the He profile gets steeper. Models including tides synchronize quickly and as a result show less helium enrichment, as the angular velocity and the differential rotation have significantly decreased when CNO burning establishes a gradient in the helium profile (see Fig. D.1 in Appendix D for more details).

These results highlight the complexity of the reaction of hydro models to tidal torques. Depending on the configuration, spin-down systems may experience more or less mixing than the single star with same initial conditions.

3.2. Grids of binary systems with tidal interactions

We now extend our analysis to broadly assess the role of tidal mixing in the evolution of massive close binary systems. We investigate whether tides enhance or suppress mixing by computing two grids of stellar models, and comparing each time the nitrogen and helium enrichment of the stars in binary systems to that obtained for the single star. We compute twin systems models (see Sect. 2.2) with initial mass M = 10, 15, 20, 30, 45 and 60 M_⊙, periods P = 1, 1.5, 2, 3, 4, 5 and 7 days and compare predictions of hydro and magnetic models. For each simulated binary system, we compute a single-star model with identical initial conditions (including the initial velocity).

In the first grid, we assume that tides during the pre-MS phase (not computed) are efficient enough for the stars to be already synchronized at zero-age-main-sequence (ZAMS). This assumption is used in most close binary systems simulations (e.g., de Mink et al. 2009; Paxton et al. 2015; Marchant et al. 2016; Fragos et al. 2023). Depending on the initial mass and period, this leads to different initial values of υ_ini/υ_crit that we report in Table 1. This first grid consists of a total of 7 (periods/velocities) × 6 (masses) × 2 (AMT assumptions) × 2 (binary/single) = 168 models. In the second grid, we make the opposite assumption that tides are inefficient during the pre-MS phase. In this case the models (whether single or binaries) are initialized with the typical velocities used in the single-star grids, i.e. υ_ini/υ_crit = 0.4. We recall that in the grids by Ekström et al. (2012) and Choi et al. (2016), this value was chosen to reproduce the observed average equatorial velocities of Galactic massive stars (Dufton et al. 2006; Huang & Gies 2006; Ekström et al. 2012) (see also Appendix A). The second grid consists of a total of (7 + 1) (periods  +  one single) × 6 (masses) × 2 (AMT assumptions) = 96 models.

Figure 5 shows the time evolution of the surface angular velocity of the 30 M_⊙ magnetic models³. With this choice of initial conditions, the binary models of the first grid are accelerated as compared to the single stars (more precisely, the stars are maintained at synchronization, which means that they receive AM from the orbit). Single-star models experience a decrease in Ω_surf with time due to the stellar winds and the MS expansion of the stars. In contrast, when the stars are initialized at υ_ini/υ_crit = 0.4 (second grid), binary models are braked until they reach synchronization. Tidal interactions transfer AM from the star to the orbit and binary models generally rotate slower than single-star models. Indeed, with our choice of initial periods and masses, the synchronization velocities at ZAMS are always smaller than υ_ini/υ_crit = 0.4 (see Table 1).

Fig. 5. Surface angular velocity evolution of the 30 M⊙ magnetic models. Binary models are represented in solid lines, single-star models in dashed lines. Upper panel: Models initialized at synchronization (first grid). Lower panel: Models initialized at υini/υcrit = 0.4 (second grid).

In order to quantify whether tides enhance or suppress mixing compared to what is found for the single star, we introduce the quantities $\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}$, defined as:

$$\begin{array}{cc}\hfill \mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}& =\frac{{[\text{N}/\text{H}]}_{\mathrm{binary}}-{[\text{N}/\text{H}]}_{\mathrm{single}}}{X_0-X_{\mathrm{c},\mathrm{final}}};\hfill \\ \hfill \mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}& =\frac{Y_{\mathrm{surf},\mathrm{binary}}-Y_{\mathrm{surf},\mathrm{single}}}{X_0-X_{\mathrm{c},\mathrm{final}}}·\hfill \end{array}$$(16)

The signs of $\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}$ indicate whether the star experiences more chemical mixing when found in a binary system (positive) or single (negative). The stopping condition is either the reaching of the TAMS (hence X_c, final = 0), or the onset of mass transfer (hence X_c, final = X_c, RLOF). For the binaries, we use the nitrogen and helium enrichment of the last computed model. For the single-star models, the values are taken at the same X_c, final as the corresponding binary model. The term in the denominator computes how much H has been burnt before the model is stopped, with X₀ the initial hydrogen mass fraction. This normalization aims at better quantifying the strength of mixing. As the different systems are initialized with different orbital periods, they do not fill their RL at the same moment of their evolution. Longer period systems are wider and thus evolve longer before overfilling the RL. One may conclude that mixing is more efficient in these longer period systems simply because they had more time to evolve, even though at a fixed evolutionary stage (characterized by the central hydrogen mass fraction X_c), their enrichment are weaker than those of shorter period systems.

3.2.1. First grid: Models initialized at synchronization

$\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}$ of the models initialized at synchronization are shown in Fig. 6. Models with P = 1 day, M = 20 − 60 M_⊙ and the 60 M_⊙ model at P = 1.5 days already fill their RL at ZAMS (gray crosses). In this configuration, binary models are accelerated compared to single-star models (Fig. 5). One may thus expect them to be systematically more enriched than single-star models ($\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}>0$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{{}_{\mathrm{surf}}}>0$). While this is systematically the case for magnetic models, the picture is more complex for hydro models. As illustrated in Sect. 3.1, this is attributed to the fact that mixing in magnetic models only depends on the value of Ω, whereas in hydro models it also depends on its gradient. By maintaining the stars at synchronization, tides reduce the steepening of the Ω-gradients, which occurs during the MS of the single stars (see Fig. E.1 in Appendix E, upper panel). As can be observed in the left panel of Fig. 6, the competition of these two effects of tidal torques (increase in Ω, decrease in Ω-gradients) splits the parameter space of (P, M) into two different regions. Binary models with low mass and short period are generally more enriched than their single-star counterpart, whereas at high mass and long period, they are generally less enriched than their single-star counterpart. Finally, all magnetic and hydro models with M < 30 M_⊙ and P > 3 days have $\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}\sim 0$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{{}_{\mathrm{surf}}}\sim 0$, which is explained by their low initial velocities (Table 1). RM is not very efficient, resulting in only minor differences in enrichment between single and binary models.

Fig. 6. Period-mass diagram comparing the enrichment of binary models initialized at synchronization to those of single-star models with identical initial conditions. The green (magenta) color is used for models with Δ > 0 (Δ < 0). A symmetric colorbar is used to illustrate the magnitude of the difference in enrichment. Models with Δ ∼ 0 appear in white. Left panels: Hydro models. Right panels: Magnetic models. Upper panels: Nitrogen enrichment. Lower panels: Helium enrichment. Systems overfilling their RL at ZAMS are represented by gray crosses.

In order to more clearly illustrate the period dependence of the nitrogen enrichment, we show in Fig. 7 the time-averaged nitrogen enrichment of the binary models, defined as

$$\begin{array}{c}\hfill \langle log(\text{N}/\text{H})\rangle =\frac{1}{t_{\mathrm{final}}}{\displaystyle {\int }_0^{t_{\mathrm{final}}}}log(\text{N}/\text{H})\text{d}t,\end{array}$$(17)

Fig. 7. Period-average nitrogen enrichment diagrams across the masses of the grid. Models are color-coded according to their average equatorial velocity. Left panel: Hydro models. Right panel: Magnetic models.

as a function of period for the considered masses. t_final can be t_final = t_RLOF or t_final = t_TAMS depending on whether the stars fill their RL during the MS. This average enrichment represents a suitable quantity for comparison with observations. Models are color-coded based on their average equatorial velocity, defined analogously to the average enrichment, highlighting that in synchronized systems, short orbital periods are associated with high rotational velocities.

The magnetic models show a clear period-nitrogen enrichment pattern. The enrichment increases with decreasing period, which can be explained by the fact that at short period models rotate faster. This is consistent with the findings of previous studies (e.g., Yoon et al. 2006; de Mink et al. 2009; Mandel & de Mink 2016; Marchant et al. 2016; Song et al. 2016). The low nitrogen enrichment of the 45 M_⊙ model at P = 1.5 day is due to the fact that it fills its RL shortly after the beginning of the evolution. We note that RM is expected to continue during a mass transfer phase (or in a contact phase, which is expected to occur at short periods), which would increase the expected average nitrogen enrichment of these close systems. Thus, the predictions of Fig. 7 are relevant for pre-interaction systems. The nitrogen enrichment of magnetic models also shows a strong mass dependence, with massive models being more enriched than lower mass ones.

The average nitrogen enrichment of hydro models shows a more complex period and mass dependence. This is mainly due to the reduction of the Ω-gradients by tides, which mitigates or even suppresses mixing. In the mass range 10–30 M_⊙, the overall period dependence is similar, but less pronounced than in magnetic models. Similarly to what is found for the 45 M_⊙ magnetic model, we notice a drop at P = 1.5 day, which is explained by the fact that models at P = 1.5 day experience RLOF early, before significant enrichment has time to occur. This effect is more pronounced in hydro models as their enrichment is slower. We also observe that models with M = 45, 60 M_⊙ are more enriched at P = 7 days than at P = 5 days. This can explained by two reasons. First, models at P = 7 days evolve longer before filling their RL, which increases their time-averaged enrichment. Second, in this high mass regime shear mixing is strong and very sensitive to the degree of differential rotation. At P = 7 days, tides are slightly less efficient in maintaining the stars at synchronization (see Fig. 5), and therefore M = 45, 60 M_⊙ models have steeper Ω-gradients than at P = 5 days, increasing the strength of shear mixing.

Figure 7 illustrates that depending on the chosen AMT assumption, models predict a clear trend in the period-enrichment diagram (magnetic models) of short-period pre-interaction systems, or a more complex period-enrichment dependence (hydro models). Although both types of models predict a mass dependence in the nitrogen enrichment, it is much more pronounced in the magnetic models. Moreover, hydro models generally predict lower nitrogen enrichment than magnetic models, providing further evidence that tides often mitigate mixing in the former. These contrasting results represent key model predictions to be tested against observations.

Finally, we note that no matter the AMT treatment, the predicted enrichment of binary models with P > 4 days are very low. Consistently with the findings of de Mink et al. (2013, their Fig. 3), in this regime tides are efficient and the synchronization velocities relatively low compared to the average velocities used for the calibration of the models (see Table 1 and Appendix A), suppressing RM.

3.2.2. Second grid: Models initialized at υ_ini/υ_crit = 0.4

$\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}$ of the models of the second grid are shown in Fig. 8. In the chosen mass and period ranges, the synchronization velocities are lower than the value υ_ini/υ_crit = 0.4 generally used in the single-star grids (Table 1). As a result, stars in binaries are braked by tides (Fig. 5), which reduces the enrichment of magnetic models compared to the single-star counterparts with identical initial conditions. The only exception is the system with P = 2 days, M = 60 M_⊙. In this case, although the initial velocity is higher than the synchronization velocity of the binary components, the single star slows down due to winds and MS expansion, so that its time-averaged rotational velocity is lower than that of the binary model, which explains the positive values found for $\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}$ and $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}$.

Fig. 8. Same as Fig. 6 but for models initialized at υini/υcrit = 0.4.

The enrichments of binary hydro models show a more complex mass and period dependence. The strong tidal torques at the beginning of the evolution increase the Ω–gradients, boosting the nitrogen transport, which explains the positive sign of $\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}$ for the 10–30 M_⊙ models (see also Fig. E.1 in Appendix E, lower panel). Consistently with the results of Sect. 3.1.1, the tides are less efficient in enhancing helium, which explains why some models have $\mathrm{\Delta }{\stackrel{\sim }{\text{N}}}_{\mathrm{surf}}>0$ but $\mathrm{\Delta }{\stackrel{\sim }{Y}}_{\mathrm{surf}}<0$. In the mass range 45–60 M_⊙, single-star models are more enriched than binary models, as in Sect. 3.1.3.

The period-average nitrogen enrichment diagrams of the second grid are shown in Fig. 9. Opposite trends between magnetic and hydro models are observed in this case. In hydro models, the average nitrogen enrichment generally increases with the period, consistently with the results of Sects. 3.1.1 and 3.1.3. This can be explained by the dominating effect of shear mixing, which is boosted by the larger degree of differential rotation imposed by the tidal torques at longer periods and the fact that longer period systems evolve longer before interacting. The model predictions significantly differ from those of Sect. 3.2.1, highlighting their sensitivity to the initial conditions.

Fig. 9. Same as Fig. 7 but for models initialized at υini/υcrit = 0.4.

In contrast, binary magnetic models’ predictions are insensitive to the initial velocity: the period-enrichment diagram is almost identical to that of the first grid. This can be explained by the fact in magnetic models, tidal torques quickly synchronize the whole star (see Figs. 2 and C.1) and mixing only depends on the value of Ω. The general trend in the period-enrichment diagram is the exact opposite than that predicted by the hydro models in this case: the enrichment decreases with the period. Period-enrichment diagrams can be interpreted as velocity-enrichment diagrams (originally introduced by Hunter et al. 2009; Brott et al. 2011b), where a long period corresponds to a low velocity, as also illustrated in Figs. 7 and 9. Interestingly, spin-down systems simulated with a hydrodynamic AMT predict higher enrichments for longer periods, in contrast with standard single-star models, where faster rotators typically exhibit stronger enrichment. We discuss in Sect. 4.2 how confronting these predictions to observations of close pre-interacting systems may lend support to or argue against particular types of AMT.

3.3. Structural reaction to chemical mixing

In the most massive stars of the grids presented in Sect. 3.2, RM is strong enough for significantly altering the surface abundance of helium, which induces substantial structural changes. The helium enhancement lowers the opacity of the surface layers, which prevents the stars from expanding, or even makes them contract (e.g., Meynet & Maeder 2000; de Mink et al. 2009; Song et al. 2016). In extreme cases this effect can lead to unusual evolutionary pathways, such as a mass transfer initiated by the secondary star (de Mink et al. 2009) or the so-called chemically homogeneous evolution (e.g., Mandel & de Mink 2016; Marchant et al. 2016), where both components stay within their RL and mass transfer is avoided.

In Fig. 10, we show the evolution of the Roche lobe filling factors R/R_L of 60 + 45 M_⊙ systems at solar metallicity in the period range 1.8–5 days. Systems are initialized at synchronization, therefore tides transfer AM from the orbit to the spin, as in Sect. 3.2.1. In this case, we successively simulated the evolution of the primary and the secondary, and stopped the simulations at RLOF, as was done in de Mink et al. (2009) for SMC models. We compare the evolution predicted by hydro and magnetic models. We also include model variations where the stars are initiated with identical initial conditions but that do no account for the tidal torques. This allows us to assess whether tides increase or suppress mixing as compared to the single-star evolution. For these model variations, we also stop the simulations at RLOF.

Fig. 10. Roche lobe filling factor evolution of 60 + 45 M⊙ systems at different orbital periods. Models accounting (not accounting) for tidal interactions are represented in solid (dashed) lines. Left panels: Hydro models: Right panels: Magnetic models. Upper panels: Primary star evolution. Lower panels: Secondary star evolution. For each simulated system a dot indicates whether the primary or the secondary fills its RL first (solid dot for the models with tides, empty dot for the models without tides).

When the tidal torques are neglected and the models behave like single stars, the primaries of the hydro models with P = 2.5 − 4 days experience strong mixing, consistently with the results of Sect. 3.2.1, which makes them contract after reaching a maximum radius smaller than the RL radius. As a result, the secondaries fill their RL first. Magnetic models neglecting tides experience little mixing and in this case the primaries systematically fill their RL first. In contrast, when the tidal torques are accounted for, in hydro models the primaries systematically fill their RL first. In magnetic models, the primaries of the P = 1.8 − 2 days systems experience strong mixing, which slows their expansion and makes the secondaries fill their RL first. This highlights opposite models’ chemical and structural reactions to tidal torques depending on the AMT assumptions.

4. Discussion

4.1. Comparison with previous studies

The primary conclusion of the studies by Song et al. (2013, 2016) is that tidal interactions in binary models consistently enhance mixing as compared to single-star models with identical initial conditions. In this paper, we expand on this finding by testing its robustness across various system configurations and AMT treatments. By reproducing the systems presented in Song et al. (2013) with different AMT assumptions, we show that in spin-down configurations, mixing is generally suppressed in binary magnetic models. The picture is different in hydro models where the spin-down boosts shear mixing. This may lead the stars in binary systems to experience more mixing than the single-star counterpart, although they synchronize at low angular velocities (Sect. 3.1.1). We demonstrate that this is not a general result, as there also exists spin-down configurations where binary hydro models predict less enrichment than what is found for the single star (Sects. 3.1.3 and 3.2.2). Additionally, in Sect. 3.2.1 we highlight that when stars are initiated at synchronization, binary hydro models may experience less mixing than the single star with identical initial velocity, although tides in this case increase their AM content. Finally, the conclusions may change depending whether helium or nitrogen is considered as tracer of RM. Although both are enhanced by the CNO process, the mixing of the processed element not only depends on the diffusion coefficients –which are the same for all species– but also on the abundance profiles themselves. The fact that nitrogen is initially less abundant than helium implies that it is more efficiently mixed at early stages of evolution, as a strong chemical gradient is quickly established.

In the study by de Mink et al. (2009), models were initiated at synchronization, assuming that tides are efficient during the pre-MS phase. The AMT treatment is similar to the magnetic models of the present study and imposes a strong core-envelope coupling. In this case, shear mixing is inefficient, and the authors found a clear anticorrelation between period and nitrogen enrichment, which we also obtain in magnetic models. In contrast, hydro models do not show such a clear anticorrelation (Fig. 7). Moreover the trend depends on the assumed initial velocity, unlike in magnetic models (Fig. 9). In Sect. 5.2 of their paper, de Mink et al. (2009) study the evolution of 50 + 25 M_⊙ short-period systems at SMC metallicity. At very short periods (P < 1.8 days), the primary experiences enough mixing for undergoing substantial structural changes, which makes it stay compact long enough for the secondary to fill its RL first. We find similar results with magnetic models at solar metallicity for 60 + 45 M_⊙ systems. However, in this mass range hydro models maintained at synchronization experience little mixing and as a result the primary fills its RL first across the whole period range, highlighting the importance of the AMT efficiency in determining the evolutionary pathways of close massive binary systems.

4.2. Comparison with observations

The identification of the main process(es) responsible for chemical mixing is crucial and still heavily debated. RM is one of the most promising explanations, but it faces several challenges. Hunter et al. (2009) and Brott et al. (2011b) proposed a direct test to the predictions of RM by inspecting correlations between the observed nitrogen enrichment of massive stars and their surface projected rotational velocity. Although the majority of stars in their sample have velocities and nitrogen enrichment compatible with models accounting for RM, they identified two categories of stars challenging the predictions: slow rotators with strong enrichments at relatively high surface gravity (typically reached during the MS phase) and fast rotators with low enrichments. de Mink et al. (2009) expanded on their idea by proposing an equivalent test in massive close binaries, where tidal synchronization is expected to boost RM, so that any system with P ≤ 3 days is expected to show significant nitrogen enrichment. Since the synchronization velocity inversely scales with the period, one expects to find an anticorrelation between period and nitrogen enrichment. Several spectroscopic studies of close binary systems have challenged these predictions. The majority of short-period systems with orbital periods smaller than 10 days in the studies by Martins et al. (2017), Pavlovski et al. (2018, 2023), and Abdul-Masih et al. (2019, 2021) do not show any evidence of strong enrichment. None of these studies report any clear anticorrelation between period and nitrogen enrichment; however, this may be attributable to insufficient sample sizes.

Mahy et al. (2020) determined the surface abundances of the 31 double-lined spectroscopic binaries in the Tarantula Massive Binary Monitoring sample, in which at least one of the components is classified as an O-type star. Their subsample 3 is composed of 18 stars in short-period (P < 10 days) circular systems. The authors found a trend between the nitrogen enrichments and the projected rotational velocities, the three stars with the highest nitrogen enrichment rotating faster than 150 km s⁻¹ at relatively short periods (P < 2 days). On the other hand, the remaining stars in their sample do not exhibit any evidence of nitrogen enrichment, despite some having rotational velocities up to 200 km s⁻¹, thereby diluting the overall trend.

The tension between observations and predictions of stellar models led (Pavlovski et al. 2018) to wonder about “the role of binarity in terms of tides on damping internal mixing in stars residing in binary systems.” The magnetic models presented in this study incorporate physics similar to that of the models of de Mink et al. (2009), and also predict a clear period-nitrogen enrichment anticorrelation. In this sense, they are also challenged by the observations. In contrast, hydro models offer a larger variety of predictions, as the period-enrichment trends also depend on the choice of initial velocities. When tides are assumed to be efficient during the pre-MS, the models do not predict any clear period-enrichment anticorrelation. Moreover, the time-averaged enrichments are lower than in magnetic models⁴, especially at high masses and short periods. The lower enrichment in binary hydro models is due to the competing effects of tidal torques on shear mixing: in this configuration they increase the AM content of the stars, but also decrease their Ω-gradients. When tides are assumed to be inefficient during the pre-MS, hydro models even predict a correlation between period and enrichment– exactly opposite to the predictions of magnetic models. This prediction would be compatible with observed systems at intermediate periods (5 days ≲P ≲ 20 days)⁵ or apparently single stars with an undetected companion, rotating at low velocities due to the spin-down by tides and showing large nitrogen enrichments.

We stress that the efficiency of tides during the pre-MS of massive binaries is extremely uncertain. Observational studies that robustly assess tidal synchronization during the pre-MS are lacking so far (Lennon et al. 2024). On the theoretical side, GENEC pre-MS single-star models with M ≥ 3 M_⊙ by Haemmerlé et al. (2019) all undergo a fully radiative phase before the onset of H-burning. The dynamical tides considered in this study, as formulated by Zahn (1977), are expected to be completely inefficient in fully radiative stars. Consequently, one may argue that the ZAMS rotation rates of massive close binaries are more likely set by their pre-MS contraction phase and other processes at play, i.e., accretion disks, magnetic fields, and jets (e.g., Seifried et al. 2011; Pudritz & Ray 2019; Commerçon et al. 2022; Oliva & Kuiper 2023) than by tidal interactions, unless we are missing an additional, efficient dissipation mechanism.

4.3. Main caveats

The main limitations of the models presented in this work are the following:

1. The value of α_ov used in Ekström et al. (2012) grid and the present works may be underestimated for stars above 9 M_⊙ (e.g., Castro et al. 2014; Tkachenko et al. 2020; Martinet et al. 2021; Baraffe et al. 2023). For the sake of the comparison with the results of Song et al. (2013) we considered necessary to adopt consistent stellar physics inputs as we mainly focused on the interplay between tides and RM in this study. We note that the results of the magnetic models are overall in agreement with those of de Mink et al. (2009), who used α_ov = 0.355.

2. The models do not directly account for mixing nor AMT by IGWs. Zahn (1977) formalism of dynamical tides is based on the assumption of IGWs traveling from the convective-core boundary to a region close to the surface, where they dissipate energy. Theoretically IGWs can also transport AM (see e.g., Talon & Charbonnel 2003; Rogers et al. 2013; Rogers & Ratnasingam 2025) and chemicals (e.g., Montalban 1994; Montalban & Schatzman 1996; Herwig et al. 2023; Brinkman et al. 2024, 2025; Varghese et al. 2024; Mombarg et al. 2025), which was neglected in this work. However, a proper treatment of tidally excited IGWs is challenging to implement, and to date, thorough investigations of the interplay between RM and AMT accounting for the transport by tidally excited IGWs are lacking. Talon & Kumar (1998) studied the joint action of meridional circulation, shear, and tidally excited IGWs in the AMT transport of 9 M_⊙ models, but did not discuss the resulting impact on RM. Recent studies by Varghese et al. (2024), Mombarg et al. (2025), Brinkman et al. (2025) compare the efficiency of mixing by stochastically excited IGWs and rotation. They conclude that mixing by IGWs offers a promising explanation to the “outliers” of the Hunter diagram, as this process is not incompatible with high enrichment at low rotations. Varghese et al. (2024) even find the mixing efficiency to decrease with rotation, providing an explanation for fast rotators with low enrichment. Interestingly, the models presented in this study display a large variety of trends, offering alternative (and not incompatible) explanations, while relying on the physics of RM. It should be noted that in binaries, the competition between stochastically and tidally excited IGWs further complicates the picture.

3. Resonant locking is not considered in this study, which may alter the synchronization velocities of the models (e.g., Witte & Savonije 2002; Burkart et al. 2012; Ma & Fuller 2021; Fellay & Dupret 2025), and in turn the efficiency of RM.

4. The models are stopped at RLOF, as we can not follow mass transfer and contact phases in GENEC at this stage. We therefore caution that the period-enrichment diagrams (Figs. 7 and 9) are relevant for pre-interacting systems.

5. We restricted our analysis to solar-metallicity models. Extending the study to lower metallicities would be of particular interest, as they influence both mass-loss rates and RM. We plan to pursue this in future work.

5. Conclusion

In this study we computed several grids of single stars and close, pre-interacting binary systems with GENEC in the mass and period range 10–60 M_⊙, 1–7 days. We took benefit of the versatility in rotation treatment of the GENEC code to study the interplay between tidal interaction and rotational mixing with various angular momentum transport assumptions. We investigated whether tidal interactions enhance or suppress mixing by computing single-star model variations with identical initial conditions. We also proposed a refined treatment of tidal interactions, inspired from Fragos et al. (2023).

We find that when a strong AMT is accounted for (magnetic models with calibrated Tayler-Spruit dynamo), stars rotate nearly as solid-bodies, making shear mixing inefficient. In this case, the efficiency of rotational mixing is only determined by the angular velocity of the stars. As a result, the impact of tidal interactions on chemical mixing comes down to the question of whether tides spin up or spin down the stars. The corollary of this result is that mixing in binary magnetic models is rather insensitive to the assumed initial velocity and is only determined by the system configuration, as long as the synchronization timescale is short with respect to the nuclear timescale.

The picture is more complex in hydro models, which are subject to shear mixing due to their moderate core-envelope coupling. In this case, the efficiency of rotational mixing also depends on the angular velocity gradients. This intrinsic difference of the models has major consequences on their reaction to tidal torques. Spin-down models may experience more mixing than single stars with identical initial conditions, although tides make them rotate at lower angular velocities. Conversely, binary models maintained at synchronization may experience less mixing than their single-star counterparts. Contrary to magnetic models, in hydro models the efficiency of rotational mixing not only depends on the system configuration, but also on the assumed rotational velocity at ZAMS.

These disparities have multiple consequences. First, the predicted evolutionary pathways of massive close binary systems may differ. The results of Sect. 3.3 indicate that the chemically homogeneous evolution (de Mink et al. 2009; Mandel & de Mink 2016; Marchant et al. 2016) is likely more difficult to achieve with hydro models. Thus, we expect the predicted merger rates to sensibly depend on the considered AMT, even though the RM efficiency of both types of models are calibrated with the same observables (see Appendix A). Second, the predicted period-enrichment diagrams –suitable for direct comparison with observations– significantly differ. While magnetic models predict a clear anticorrelation between period and nitrogen enrichment, the picture is less clear in hydro models, and the predictions strongly depend on the ZAMS velocity. Notably, when tides are assumed to be inefficient during the pre-MS, the latter even predict a correlation between period and enrichment.

Recent observations challenge the predictions of models accounting for rotational mixing, as numerous massive close binaries rotating at high velocities do not show any evidence of nitrogen enrichment (e.g., Martins et al. 2017; Pavlovski et al. 2018, 2023; Abdul-Masih et al. 2019, 2021). Conversely, some longer period systems in the sample by Mahy et al. (2020) show larger nitrogen enrichment than systems in tighter orbits. No clear anticorrelation between period and enrichment is reported in these studies. Consistently with the findings of de Mink et al. (2009), our magnetic models show a clear period-nitrogen anticorrelation, and are thus challenged by these observations. Without claiming that these observations strongly favor hydro models, we note that they better stand the comparison, owing to the greater versatility of their predictions –stemming from their higher sensitivity to initial conditions. The need for statistically significant observational constraints on RM remains critical– a gap that current and upcoming large-scale surveys, such as IACOB (Simón-Díaz & Herrero 2014; de Burgos et al. 2024; Martínez-Sebastián et al. 2025) and BLOeM (Shenar et al. 2024), will help to narrow.

Data availability

The tracks are available on Zenodo, Yareta, and on the Geneva stellar group database.

Acknowledgments

We sincerely thank the anonymous referee for their constructive feedback and valuable suggestions, which enabled us to improve the manuscript. LS, SR and SE have received support from the SNF project No 212143. MM and PE have received support from the SNF project No 219745.

References

Appendix A: Calibration of the magnetic models

In this section, we present the calibration of RM performed for the magnetic models. For this calibration to be consistent with that of the hydro models, the rotational velocities must also be calibrated to reproduce the same observations.

Fig. A.1 shows the MS evolution of the equatorial velocity of 15 M_⊙ single-star models. The evolution of two magnetic models and one hydro model are shown. One of the magnetic model and the hydro model are initiated at initial velocity υ_ini/υ_crit = 0.4. The other magnetic model is initiated at υ_ini/υ_crit = 0.33. Soon after ZAMS, the velocity of the hydro model drops from ∼260 km s⁻¹ to ∼220 km s⁻¹. This is a consequence of the advective AMT by meridional circulation, which transports AM from the surface to the core (models are initiated with solid body rotation at ZAMS). Magnetic models with a purely diffusive AMT do not show the same behavior and hence the average equatorial velocity along the MS evolution is higher than that of the hydro model with the same initial value of υ_ini/υ_crit. The value υ_ini/υ_crit = 0.4 was chosen so that hydro models reproduce the observed average equatorial velocities of Galactic massive stars (Dufton et al. 2006; Huang & Gies 2006; Ekström et al. 2012). In order to ensure that the magnetic models also reproduce this constraint, we consider that the initial value of υ_ini/υ_crit required for performing the calibration should be slightly reduced. We found that υ_ini/υ_crit = 0.33 is a satisfactory choice as the average MS equatorial velocity of the 15 M_⊙ magnetic model in this case matches that of the hydro model with υ_ini/υ_crit = 0.4.

Fig. A.1. MS evolution of 15 M⊙ stars with different initial velocities and AMT treatments.

Then, the RM efficiency of the magnetic models is calibrated by tuning the value of c_h, which is a free parameter in the expression of D_h. We use the nitrogen enrichment as tracer of RM. We follow (Ekström et al. 2012) and use observed nitrogen abundances of Galactic MS B-type stars with initial masses up to 20 M_⊙ from (Gies & Lambert 1992; Kilian 1992; Morel et al. 2008; Hunter et al. 2009). Combining these observations, Ekström et al. (2012) provided two boxes of nitrogen enrichment Δlog(N/H) (at mid and end-MS) that models should reproduce. Fig. A.2 shows the MS evolution of Δlog(N/H) for GENEC models with initial masses M = 9, 12, 15 and 20 M_⊙. Three sets of models are shown: hydro models identical to those in Ekström et al. (2012), magnetic models without calibration, i.e., with c_h = 1 and calibrated magnetic models with c_h = 0.70. Uncalibrated magnetic models fail at reproducing the observations, as the mid-MS enrichment of the 9 M_⊙ model and all the enrichments of the 12-20 M_⊙ models are outside the boxes. The calibrated magnetic models show a better agreement with the observations, yet less satisfactory than the hydro models. The first reason is that the enrichment of the magnetic models occurs mainly at the beginning of the MS evolution, whereas hydro models show a more prolonged enrichment. This is explained by the fact that in the latter models mixing depends on both Ω and its gradient. The MS expansion and AM losses reduce the surface angular velocity of the models, but also steepens the Ω–gradients due to the moderate AMT, which sustains shear mixing. In magnetic models, shear mixing is inefficient and the enrichment is suppressed when the angular velocity decreases. As a result, it is harder to match both the mid- and the end-MS enrichments with magnetic models. The value of c_h was chosen so that a maximum number of models reproduce the observations across the mass range. With this calibration, the 9 M_⊙ model just enters the box at mid-MS. At end-MS, it is however lower than the box. In contrast, the 20 M_⊙ model is too enriched and lies outside the box both at mid- and end-MS. The 15 M_⊙ model enters the box at end-MS but not at mid-MS. This highlights that RM in magnetic models exhibits a stronger mass dependence than in hydro models, which is the second reason why magnetic models overall reproduce the observed nitrogen abundances less satisfactorily.

Fig. A.2. MS nitrogen enrichment of stars in the mass range 9 − 20 M⊙ with different AMT assumptions and RM efficiencies. Hydro models are shown in dashed lines, calibrated magnetic models in solid lines, uncalibrated magnetic models in dot-dashed lines. The gray boxes are compiled observations from Ekström et al. (2012).

Appendix B: Evolution toward synchronization with and without equilibrium tides

We show in Fig. B.1 the evolution toward synchronization of the same system as in Sciarini et al. (2024) with an initial period of P = 3 days. The evolution obtained with only dynamical tides is compared to that including equilibrium tides, treated either with the prescription adopted in this work (Eq. (8)) and that of (Fragos et al. 2023, i.e, Eq. (7)). This illustrative example highlights the importance of equilibrium tides in the tidal evolution and allows us to assess the relative strength of the present prescription against that of Fragos et al. (2023)⁶. Models are stopped at RLOF.

Fig. B.1. Evolution toward synchronization of the same system as in Sciarini et al. (2024) with P = 3 days including model variations accounting for the equilibrium tides. Upper panel: Angular velocity. Middle panel: Synchronization timescales (the different lines correspond to the dynamical tides and equilibrium tides synchronization timescales, not the "global" synchronization timescales of each computed model). Lower panel: Convective turnover timescales.

At the beginning of the simulation, the dynamical tides dominate the evolution as the three lines overlap until t ∼ 6 Myrs. At the end of the evolution, the equilibrium tides dominate, as shown by the evolution of the synchronization timescales, defined as ${\tau }_{\mathrm{sync}}=\left|\frac{{\mathrm{\Omega }}_{\mathrm{spin}}-{\mathrm{\Omega }}_{\mathrm{orb}}}{{\dot{\mathrm{\Omega }}}_{\mathrm{spin}}}\right|$. As a result, the models accounting for equilibrium tides are closer to synchronization at RLOF. The relative departure from synchronization $\left|\frac{{\mathrm{\Omega }}_{\mathrm{spin}}-{\mathrm{\Omega }}_{\mathrm{orb}}}{{\mathrm{\Omega }}_{\mathrm{orb}}}\right|$ is 0.66% for the model computed with the Fragos et al. (2023) prescription, 3% for the model computed with the updated prescription of this work, and 31% for the model accounting only for the dynamical tides. The large departure from synchronization at RLOF of the model without equilibrium tides can be imparted to the decrease of E₂ along the MS and the power 8/3 in Ω_orb − Ω_spin of the dynamical tides torque (Eq. (6)), which makes it very inefficient close to synchronization. Equivalently, the dynamical tides synchronization timescale diverges at synchronization, which can be seen at t ∼ 3.4 Myrs, when the star "crosses" the synchronization. Unlike the equilibrium tides, the dynamical tides synchronization timescale depends on Ω_orb − Ω_spin (Eq. (A.1) in Sciarini et al. 2024), which explains the significantly different evolution the timescales follow. The evolution of the equilibrium tides synchronization timescales follow the same trend as (a/R)⁶, which is shown as a reference.

Under the formulation by Fragos et al. (2023), the synchronization timescale is smaller by a factor ∼10 than the prescription used in this work. As shown in the lower panel of Fig. B.1, one of the reasons is that the convective turnover timescales obtained with Eq. (11) do not match the values predicted by the MLT approach. Two convective zones are present throughout the evolution of the models presented in Fig. B.1, consistently with the findings of Maeder et al. (2008) and He et al. (2025), for instance. The first one, the closest to the surface and the smallest in size and mass is due to the partial helium ionization zone, the second one, deeper in the envelope and larger in size, to the partial iron ionization zone. The helium zone being of very small size, Eq. (11) predicts it to have τ_conv ∼ 10² − 10³ s, whereas the MLT approach predicts τ_conv ∼ 5 ⋅ 10⁶ s throughout the evolution. The mass and radius extent of the zone being very small, it does not contribute much to the total torque (see Eqs. (7) and (8)). The equilibrium tides are dominated by the other zone for which the two timescales differ by about one order of magnitude.

Given that several terms differ between the present prescription (Eq. (8)) and that (Fragos et al. 2023, i.e, Eq. (7)), it is not trivial to understand why the resulting synchronization timescales differ by one order of magnitude in Fig. B.1. The ratio of the Fragos et al. (2023) timescale (τ_sync, F23) to the one obtained with Eq. (8) (τ_sync, S25) is:

$$\begin{array}{cc}\hfill \frac{{\tau }_{\mathrm{sync},\mathrm{F}23}}{{\tau }_{\mathrm{sync},\mathrm{S}25}}& =\frac{{\left(\frac{k_2}{T}\right)}_{\mathrm{S}25}}{{\left(\frac{k_2}{T}\right)}_{\mathrm{F}23}}\hfill \\ \hfill & =\underset{\sim 10}{\underbrace{\frac{21}{2}}}\underset{\sim {10}^{-2}}{\underbrace{k_2}}\underset{\sim {10}^{-1}}{\underbrace{\frac{{\tau }_{\mathrm{conv},\mathrm{F}23}}{{\tau }_{\mathrm{conv},\mathrm{MLT}}}}}\underset{=1\text{close to sync.}}{\underbrace{\frac{f_{\mathrm{conv},\mathrm{MLT}}}{f_{\mathrm{conv},\mathrm{F}23}}}}\underset{\sim 10}{\underbrace{\frac{\frac{I_{\mathrm{conv}.\mathrm{reg}.}}{I}}{\frac{M_{\mathrm{conv}.\mathrm{reg}.}}{M}}}}\hfill \\ \hfill & \sim {10}^{-1}.\hfill \end{array}$$(B.1)

Eq. (B.1) only gives an order of magnitude estimate for the system considered in Fig. B.1. The values of the different terms in the equations were empirically obtained from the models. Deviations between the two prescriptions larger than a factor 10 may arise depending on the initial parameters of the models. Eq. (B.1) aims at illustrating that the two prescriptions contain factors divergent by up to several order of magnitudes. As a result the predicted synchronization timescales may well vary by more than an order of magnitude. For instance, over a full MS evolution, k₂ may vary in the range 10⁻² − 10⁻⁴ (e.g., Claret 2004; Rosu et al. 2020a, 2022a; Sciarini et al. 2025). Although in our simulations the mass of the convective shells are quite stable over the evolution (they slightly increase), their moment of inertia generally increase more due to the radial expansion of the star and the r² dependence of I_conv.reg., thus the last term in Eq. (B.1) may vary accordingly.

Appendix C: Evolution toward synchronization of the P = 1.4 days models of Sect. 3.1.1

Fig. C.1. Ω–profiles of the P = 1.4 days models of Sect. 3.1.1 as they evolve toward synchronization. Each line represents a different epoch, indicated by its color. Dashed line: orbital angular velocity. Left panel: Hydro models. Right panel: Magnetic models.

Fig. C.1 shows the Ω–profiles of the models at P = 1.4 days of Sect. 3.1.1, as they evolve toward synchronization. In the hydro model, the tidal torque efficiently slows down the layers close to the surface. After ∼0.2 Myrs, they are already close to synchronization. In contrast, at the same age the core is still rotating super-synchronously, due to the mild core-envelope coupling. This corresponds to the moment when shear mixing is the most efficient. Interestingly, during the early evolution, the core slightly accelerates instead of slowing down, which can only be explained by the advective nature of the meridional currents. The surface being braked, diffusive processes can only transport AM outward. Shortly after, the differential rotation becomes so important that AMT switches direction, slowing down the core. During this period, the surface layers are already close to synchronization, thus the tidal torque is not efficient anymore, (see Eq. 6 and remark below). The decay of the tidal torque, along with the AMT from core to surface delays the synchronization of the surface, which is achieved after ∼2.25 Myrs.

In the magnetic model, the AMT is very efficient and the star rotates almost rigidly throughout the whole evolution. A slight degree of differential rotation is observed at the beginning of the evolution, when the tidal torque is maximal. We note that the star reaches synchronization after ∼0.7 Myrs, i.e., faster than the hydro model, which can seem counterintuitive. One might expect that the stronger the AMT, the more AM needs to be transferred to the star to reach synchronization, which should delay the moment when it is achieved. This picture is reasonably valid during the early evolution. After 5 kyrs the surface velocity of the hydro model is lower than that of the magnetic model (Ω_{spin, hydro} ∼ 6.8 × 10⁻⁵ rad s⁻¹, Ω_spin, magn ∼ 7.5 × 10⁻⁵ rad s⁻¹). The dynamical tides formalism by Zahn (1977) predicts that the tidal torque scales as (Ω_orb − Ω_spin)^8/3, the synchronization timescale as |Ω_orb − Ω_spin|^−5/3. The tidal torque gets extremely inefficient close to synchronization (the synchronization timescale diverges), and in the hydro model, unlike the magnetic model, the core is still rotating super-synchronously when this occurs. As a result, it experiences a long phase where AMT is transported from the core to the surface, delaying synchronization.

Appendix D: Profiles of the hydro models of section 3.1.3

Fig. D.1. Same as Fig. 2 but for the hydro models of Sect. 3.1.3. Additionally, the helium mass fraction and nitrogen enhancement profiles are displayed. The profiles are shown at t = 40 kyrs and t = 3.4 Myrs.

The angular velocity, shear diffusion coefficient, helium mass fraction and nitrogen enhancement profiles of the hydro models of section 3.1.3 are displayed in Fig. D.1 at t = 40 kyrs and t = 3.4 Myrs. At t = 40 kyrs the binary models are evolving toward synchronization and therefore have sharp Ω-gradients. As a result, they have much larger shear coefficients than the single-star model. The nitrogen profiles are already significantly altered by CNO burning, and as a result, RM is efficient in diffusing nitrogen. Hence, the models show surface nitrogen enhancements. In contrast, the helium profiles are still almost flat, which is why RM is not efficient in enhancing helium at the surface at this point of the evolution, even though the shear coefficients are large. Later in the evolution, at t = 3.4 Myrs (which corresponds to the moment when the primary of the system with P = 3 days overfills its RL), the binary models have been synchronized for a long time and rotate with lower angular velocities and Ω-gradients than the single-star model. As a result, the single-star model has the highest shear coefficient. At this point, the model already shows larger helium enrichment than the binary models, and a larger nitrogen enrichment than the model with P = 3 days. Later in the evolution and before the binary models with P = 5 and 8 days fill their RL (see Fig. 4), even its nitrogen enhancement gets larger than those of the binary models.

Appendix E: Differential rotation in the hydro models of Sect. 3.2

Fig. E.1 shows the evolution of the core-to-surface rotation ratio of the hydro models of the two grids of Sect. 3.2 at period P = 3 days. The single-star evolutions are shown as a comparison. In the whole mass range, binary models initialized at synchronization evolve with a lower core-to-surface rotation ratio than the single stars initialized with identical velocities. Depending on the initial mass, in single-star models Ω_core/Ω_surf reaches values from ∼1.6 (60 M_⊙) to ∼4 (10 M_⊙). Smaller degrees of differential rotation are reached in higher mass models because the AMT is more efficient and the simulations are stopped earlier⁷.

Fig. E.1. Time evolution of the core-to-surface rotation ratio of the hydro models of the grids of Sect. 3.2 at P = 3 days. Upper panel: First grid (models initialized at synchronization). Lower panel: Second grid (models initialized at υini/υcrit = 0.4).

In contrast, in binary models Ω_core/Ω_surf does not exceed∼1.2 (the fast increase at the end of the evolution in the 10 M_⊙ model is due to its end-MS contraction), which is why tidal torques may suppress rotational mixing in this case. By maintaining synchronization, tidal interactions reduce the Ω–gradients, thereby diminishing the efficiency of shear mixing.

Binary models of the second grid reach high degrees of differential rotation (Ω_core/Ω_surf ∼ 2.3 for the 10 M_⊙ model, Ω_core/Ω_surf ∼ 1.6 for the 60 M_⊙ model) in the early evolution, due to the braking of the surface by tidal torques. Subsequently, the core-to-surface rotation ratio decreases as AM is transported from the core to the surface. The AMT is more efficient at higher mass, which is the reason why Ω_core/Ω_surf drops faster for these models. At the end of the evolution, models have Ω_core/Ω_surf ∼ 1.1 − 1.2, consistently with those of the first grid. In single-star models, the degree of differential rotation steadily increases along the MS evolution, reaching values up to ∼2. This illustrates that in the hydro binary models of the second grid, the spin-down by tides boosts rotational mixing in the early evolution, by increasing the Ω–gradients.

All Tables

All Figures

	Fig. 1. Nitrogen and helium enhancement evolution of the same spin-down systems as in Song et al. (2013). Left panel: Hydro models. Right panel: Magnetic models.
In the text

	Fig. 2. Upper panel:Ω profiles of the hydro and magnetic models at t = 500 kyrs of the same spin-down systems as in Song et al. (2013). Lower panel: Diffusion coefficient profiles.
In the text

	Fig. 3. Same as Fig. 1 but for spin-up models.
In the text

	Fig. 4. Same as Fig. 1 but for models with M1 = 45 M⊙, M2 = 30 M⊙, and P = 1.8, 3, 5, and 8 days.
In the text

	Fig. 5. Surface angular velocity evolution of the 30 M⊙ magnetic models. Binary models are represented in solid lines, single-star models in dashed lines. Upper panel: Models initialized at synchronization (first grid). Lower panel: Models initialized at υini/υcrit = 0.4 (second grid).
In the text

Fig. 6. Period-mass diagram comparing the enrichment of binary models initialized at synchronization to those of single-star models with identical initial conditions. The green (magenta) color is used for models with Δ > 0 (Δ < 0). A symmetric colorbar is used to illustrate the magnitude of the difference in enrichment. Models with Δ ∼ 0 appear in white. Left panels: Hydro models. Right panels: Magnetic models. Upper panels: Nitrogen enrichment. Lower panels: Helium enrichment. Systems overfilling their RL at ZAMS are represented by gray crosses.

In the text

	Fig. 7. Period-average nitrogen enrichment diagrams across the masses of the grid. Models are color-coded according to their average equatorial velocity. Left panel: Hydro models. Right panel: Magnetic models.
In the text

	Fig. 8. Same as Fig. 6 but for models initialized at υini/υcrit = 0.4.
In the text

	Fig. 9. Same as Fig. 7 but for models initialized at υini/υcrit = 0.4.
In the text

Fig. 10. Roche lobe filling factor evolution of 60 + 45 M⊙ systems at different orbital periods. Models accounting (not accounting) for tidal interactions are represented in solid (dashed) lines. Left panels: Hydro models: Right panels: Magnetic models. Upper panels: Primary star evolution. Lower panels: Secondary star evolution. For each simulated system a dot indicates whether the primary or the secondary fills its RL first (solid dot for the models with tides, empty dot for the models without tides).

In the text

	Fig. A.1. MS evolution of 15 M⊙ stars with different initial velocities and AMT treatments.
In the text

	Fig. A.2. MS nitrogen enrichment of stars in the mass range 9 − 20 M⊙ with different AMT assumptions and RM efficiencies. Hydro models are shown in dashed lines, calibrated magnetic models in solid lines, uncalibrated magnetic models in dot-dashed lines. The gray boxes are compiled observations from Ekström et al. (2012).
In the text

Fig. B.1. Evolution toward synchronization of the same system as in Sciarini et al. (2024) with P = 3 days including model variations accounting for the equilibrium tides. Upper panel: Angular velocity. Middle panel: Synchronization timescales (the different lines correspond to the dynamical tides and equilibrium tides synchronization timescales, not the "global" synchronization timescales of each computed model). Lower panel: Convective turnover timescales.

In the text

	Fig. C.1. Ω–profiles of the P = 1.4 days models of Sect. 3.1.1 as they evolve toward synchronization. Each line represents a different epoch, indicated by its color. Dashed line: orbital angular velocity. Left panel: Hydro models. Right panel: Magnetic models.
In the text

	Fig. D.1. Same as Fig. 2 but for the hydro models of Sect. 3.1.3. Additionally, the helium mass fraction and nitrogen enhancement profiles are displayed. The profiles are shown at t = 40 kyrs and t = 3.4 Myrs.
In the text

	Fig. E.1. Time evolution of the core-to-surface rotation ratio of the hydro models of the grids of Sect. 3.2 at P = 3 days. Upper panel: First grid (models initialized at synchronization). Lower panel: Second grid (models initialized at υini/υcrit = 0.4).
In the text