{\gamma} Doradus stars as test of angular momentum transport models

Helioseismology and asteroseismology of red giant stars have shown that the distribution of angular momentum in stellar interiors, and its evolution with time remains an open issue in stellar physics. Owing to the unprecedented quality of Kepler photometry, we are able to seismically infer internal rotation rates in \gamma Doradus stars, which provide the MS counterpart to the red-giants puzzle. We confront these internal rotation rates to stellar evolution models with rotationally induced transport of angular momentum, in order to test angular momentum transport mechanisms. We used a stellar model-independent method developed by Christophe et al. in order to obtain seismically inferred, buoyancy radii and near-core rotation for 37 \gamma Doradus stars observed by Kepler. We show that the buoyancy radius can be used as a reliable evolution indicator for field stars on the MS. We computed rotating evolutionary models including transport of angular momentum in radiative zones, following Zahn and Maeder, with the CESTAM code. This code calculates the rotational history of stars from the birth line to the tip of the RGB. The initial angular momentum content has to be set initially, which is done by fitting rotation periods in young stellar clusters. We show a clear disagreement between the near-core rotation rates measured in the sample and the rotation rates obtained from evolutionary models including rotationally induced transport following Zahn (1992). These results show a disagreement similar to that of the Sun and red giant stars. This suggests the existence of missing mechanisms responsible for the braking of the core before and along the MS. The efficiency of the missing mechanisms is investigated. The transport of angular momentum as formalized by Zahn and Maeder cannot explain the measurements of near-core rotation in main-sequence intermediate-mass stars we have at hand.


Introduction
Distribution of angular momentum in stellar interiors, and its evolution with time remains an open issue in stellar physics. The first reason is that internal angular momentum distribution has been poorly constrained by observations until recently. The tightest constraint on internal angular momentum is provided by seismic measurements of rotation profiles, thanks to helio-and asteroseismology.
One recent success of this approach was obtained for subgiant and red giant stars. Being in the late stages of stellar evolution, these stars have a highly condensed core, and hence present non-radial mixed modes of oscillation. These modes are of particular interest for the determination of the rotation profile throughout the star, as they carry information on the star's innermost layers and are detectable at the surface. The NASA Kepler spacecraft has allowed a significant leap forward Evolutionnary sequences (Fig. 7) and evolution of stellar global parameters are only available at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsarc. u-strasbg.fr/viz-bin/qcat?J/A+A/626/A121. providing exquisite seismic observations of thousands of such stars. Their rotationally split multiplets give access to nearcore rotation rates (Beck et al. 2012;Deheuvels et al. 2012;Mosser et al. 2012), which were found to be much slower than predicted by the current models including physically motivated angular momentum transport mechanisms (Eggenberger et al. 2012;Marques et al. 2013;Cantiello et al. 2014;Fuller et al. 2014;Belkacem et al. 2015). The disagreement with observations points toward the existence of missing mechanisms which would extract angular momentum from the core. Currently, the community effort is into two complementary approaches. On the one hand increasing the number of measurements and inferring correlations with other stellar properties (Gehan et al. 2016), and on the other hand developing parametrized models of core-toenvelope coupling processes (in terms of e.g., time-scales and efficiency, see for instance Eggenberger et al. 2017).
A crucial path forward is to get an insight on the angular momentum distribution at earlier stages of evolution, that is on the main sequence. Unfortunately, solar-like stars on the main sequence only exhibit pressure modes to the level of detection, and these modes probe their superficial layers. The only solarlike star on the main sequence for which seismology provides information on the deeper layers is the Sun. The long-lived satellite SoHO has successfully provided seismic observations of the Sun's interior, which allowed the inversion of the solar internal rotation profile (e.g. Schou et al. 1998;García et al. 2007;Fossat et al. 2017). As with red giant stars, these seismic measurements first questioned the understanding of angular momentum transport mechanisms.
In order to obtain inner rotation profiles for a number of stars on the main sequence, and for stars that are progenitors of red giants stars, one has to consider stars that are slightly more massive than the Sun. Among these stars, a quite interesting sample is that of γ Doradus (γ Dor) stars which are late Ato early F-type stars with masses ranging from 1.3 to 1.9 M . They burn hydrogen in their convective cores, surrounded by a radiative region where a shallow convective layer subsists. Such shallow convective layers give birth to gravity oscillation modes (g-modes) excited by the convective blocking mechanism (Dupret et al. 2005). Due to the structure of these stars, the detected oscillations are able to probe the deep internal region of these stars. With periods typically of the order of one day, these oscillations were extremely difficult to detect from ground. One had to wait for the four years of nearly continuous observations from Kepler to obtain data of sufficient quality in order to perform seismic studies of these stars. Moreover, unlike solar-type stars, the convective envelopes in γ Dor stars are too shallow to generate a magnetic field able to act as a magnetic torque and slow down their surfaces (Schatzman 1962). For that reason these stars have projected rotation velocities of around 100 km s −1 in average (Abt & Morrell 1995). These highrotation velocities have long hampered the interpretation of γ Dor stars seismology.
The first Kepler γ Dor stars analysed were slow rotators, for which it was still possible to retrieve rotationally split gmodes. Hence, the first internal rotation periods, by choice of method, were found to be very slow (of the order of a few tens of days, see Kurtz et al. 2014;Saio et al. 2015;Schmid et al. 2015;Keen et al. 2015;Murphy et al. 2016). New methods had to be developed before the rapid rotators could be analysed. Based on patterns in the oscillation periods of their seismic spectra, it was finally possible to infer the near-core rotation rates, for rapid rotators (Van Reeth et al. 2016;Ouazzani et al. 2017). In particular, Ouazzani et al. (2017) established a one-to-one relation between an observable of the periodogram and the inner rotation rate, which is valid for the whole γ Dor instability strip. These studies allowed us to find near-core rotation frequencies for these stars ranging from 5 to 25 µHz -periods of 0.5 to 2.5 days-for dozens of these stars.
The present study is a first attempt to confront these findings with the models of angular momentum (AM) evolution in intermediate-mass stars. In particular, we have made use of the stellar evolution code cestam (Marques et al. 2013), which models AM evolution from the stellar birth-line to the tip of the red-giant branch. These models are presented in Sect. 2, where the choices for angular momentum transport processes are explained. Such calculations require the definition of appropriate initial conditions, which set the initial AM content. These initial conditions are chosen using observations of young stellar forming regions as a constraint. This constitutes Sect. 3. In order to follow AM distribution along evolution, we defined a core-averaged property, the buoyancy radius, and show that it can be readily used as an age indicator on the main sequence of these stars (Sect. 4). Finally, we report on the results in Sect. 5, before discussing them, and drawing conclusions in Sect. 6.

Angular momentum transport
The stellar models used in this study have been computed with the cestam evolution code which originates from the cesam code Morel (1997), Morel & Lebreton (2008), where rotationally induced transport has been implemented (Marques et al. 2013). In convective zones, although there is differential rotation in latitude, the mean rotation at a given radius weakly depends on the radius. Therefore, in cestam models, that are one-dimensional, we assumed that convective zones rotate as solid bodies. In radiative zones, the transport of angular momentum and chemical elements is modelled following the formalism of Zahn (1992, hereafter Z92), refined in Maeder & Zahn (1998). According to these studies, because of the stable stratification in radiative zones, turbulence is much stronger in the horizontal than in the vertical direction. Thus we make the hypothesis of shellular rotation -that is the rotation rate is almost constant on isobars.
The transport of angular momentum obeys an advectiondiffusion equation: where ρ is the stellar density, r the radial coordinate, Ω the rotation angular velocity, U 2 is the vertical component of the meridional circulation, and ν v is the vertical component of the turbulent viscosity. Meridional circulation components are calculated following Maeder & Zahn (1998). The shear induced turbulence is considered to be a highly anisotropic diffusive process. For the diffusion coefficients, we chose the prescriptions for horizontal turbulent diffusion coefficient from , and for the vertical one Talon et al. (1997). According to Amard et al. (2016), this combination, together with Matt et al. (2012Matt et al. ( , 2015 for the loss of AM by magnetized winds, give the best fit to rotation periods of solar-like stars in clusters. In the current study, we did not consider loss of AM at the surface. The external convective zone in γ Dor stars being much shallower than that of solar-like stars, it has been assumed that generation of magnetic field by a dynamolike process would be inefficient. Turbulence in the vertical direction also mixes chemical elements. This mixing is further enhanced by the large-scale meridional circulation coupled to a strong horizontal turbulence (Chaboyer & Zahn 1992). As a result, the equation of the chemical composition evolution follows: where X i is the abundance by mass of the i-th nuclear species, D ν = ν v is the vertical component of turbulent diffusion, and D eff the vertical diffusivity and the diffusion coefficient associated with meridional circulation. It should be noted that in the present study, we neglect atomic diffusion, whose effects should be small compared to turbulent diffusion induced by differential rotation. In the following, for shortness, we refer to the formalism and prescriptions mentioned in this section as Z92.

Stellar models
These prescriptions have been used in order to compute stellar models for masses between 1.4 and 1.8 M . At each metallicity, we have derived the helium mass fraction using a primordial A121, page 2 of 11 Y p = 0.248 (Peimbert et al. 2007), and a helium-to-metal enrichment ratio of ∆Y ∆Z = 2.1 ± 0.9 (Casagrande et al. 2007). We adopted the agss09 solar metal mixture (Asplund et al. 2009) and corresponding opacity tables obtained with opal opacities , completed at low temperature (log T < 4.1) with Alexander & Ferguson (1994) opacity tables. We used opal equation of state . We used nacre nuclear reaction rates of Angulo et al. (1999) except for the 14 N + p reaction, for which we used the reaction rates given in Formicola et al. (2004). The Schwarzschild criterion was used to determine convective instability. Convection was treated using the mixing-length theory (MLT) formalism (Böhm-Vitense 1958) with a parameter α MLT = 1.70. The centrifugal acceleration is taken into account by adding the average centrifugal acceleration 2Ω 2 r/3 to gravity in the hydrostatic equilibrium equation. The atmosphere is matched to a T (τ) law at an optical depth of τ = 20.

Initial conditions
As mentioned in introduction, stellar evolution calculations require the definition of appropriate initial conditions, which set the initial AM. Stars with masses lower than about 2 M are born fully convective on the Hayashi track, and then develop a radiative zone on the Henyey track before they reach the zero-agemain-sequence and ignite nuclear reactions through the CNO cycle. There, they develop a convective core and the outer convective zone shrinks drastically in a mainly radiative envelope.
This whole pre-main-sequence (PMS) stage for a stellar mass typical of γ Dor stars lasts around 8-10 Myrs depending on the metallicity. The PMS is too short for internal transport of angular momentum by rotationally induced processes or stellar winds to slow down the star significantly. While contracting, these stars are prevented to spin up due to tight interaction with their residual accretion disc, before it fully dissipates after a few Myrs (up to around 5 Myrs for the low-mass end of γ Dor stars).
The star-disc interaction in this early stage involves a series of complex mechanisms which include accretion and magnetic interaction between the star and the disc, and the issue of AM exchange between the star and its environment still remains controversial. The choice here is to rely on an empirical description of that interaction, which simply assumes that the stellar angular velocity remains constant as long as the star interacts with its disc (see Bouvier et al. 1997). The problem is then reduced to two free parameters: the accretion disc lifetime τ disc , that is the time during which the star is forced to co-rotate with its disc, and the period of rotation of the disc P disc .

Rotation distributions in young clusters
In a similar fashion as in Amard et al. (2016), or before them Gallet & Bouvier (2013), we aim at anchoring the evolution of AM at the pre-main-sequence stage by reproducing the rotational distributions found in very young stellar clusters. The specificity here is that we narrow down the mass range of stars in these clusters to the corresponding masses of γ Dor stars, that is 1.3-1.9 M .
In order to set up the free parameters of the disc locking model (see Sect. 2.3), we have selected stellar clusters younger than 20 Myrs, which have surface rotation measurements available in the literature, and which would contain a significant number of stars in our mass range (1.3-1.9 M ). Three clusters fulfil these requirements: NGC 2264, NGC 2362, and hPer. The youngest, NGC 2264 has an age around 3 Myrs (Affer et al. 2013, and references therein), its stellar surface rotation periods have been measured by rotational signature in the stellar lightcurves by Venuti et al. (2017). NGC 2362 is about 5 Myrs old, and the rotational data have been taken from Irwin et al. (2008). Finally hPer (NGC 869) is around 13 Myr old, and its stellar surface rotations have been measured by Moraux et al. (2013). We have computed the median, 20th and 80th percentiles of these distributions. These have been taken respectively as reference points for initial conditions of typical models of 1.4, 1.6 and 1.8 M (computed as described in Sect. 2). In other words, τ disc and P disc have been tuned in order to fit these reference points in NGC 2264 and NGC 2362, and we have ensured that the models for the lowest mass also agree with the rotational distribution in hPer, when the disc has dissipated.
The result is shown in Fig. 1. The three sets of disc-locking parameters which allow us to reproduce the distributions satisfactorily are: -P disc = 2.4 days and τ disc = 3 Myrs, -P disc = 3.9 days and τ disc = 3 Myrs, -P disc = 7.2 days and τ disc = 5 Myrs.
The aim here is not to reproduce precisely the observed rotational distributions in these clusters, but to rather get the overall orders of magnitude correct. One should bear in mind that the uncertainties on age can be relatively important for these clusters ages. Indeed, the observed ages (data points: filled circles in Fig. 1) are obtained by isochrone fitting, and therefore are contaminated by inaccuracies in the stellar evolution models used to generate the isochrones. Moreover, dispersion in age can come from different generation of stars in the stellar forming region. Finally, concerning the stellar models which are fitted to these data points (solid lines in Fig. 1), they are also altered by the uncertainties on the birth lines location compared to what has been taken as age zero in the stellar evolution.
While the distributions in NGC 2264 and NGC 2362 are correctly reproduced by the three sets of models, only the 1.4 M models manage to fit the rotation distribution in hPer.
Unfortunately, given that hPer only contains stars with masses around 1.3-1.5 M , it does not allow us to constrain the higher mass models. At these very young age, we were not able to find clusters containing a significant number of stars of higher masses.

Buoyancy radius as an age indicator
For field stars on the main sequence, it is even less straightforward to determine the age with accuracy. Saio et al. (2015) mention the possibility to use the period spacing as a constraint of the evolution stage. Given the dependency of the period spacing on the rotation rate in the cavity probed by the g-modes, we would rather use the buoyancy radius as an age indicator for γ Dor stars on the main sequence.

Behaviour of the buoyancy radius P 0
The reason for this is that the buoyancy radius, P 0 , of a star is a monotonously decreasing function of age. P 0 depends on the internal structure only, via N, the Brunt-Väisälä frequency, and can be expressed as: where where the Brunt-Väisälä frequency is expressed by the help of g, the gravitational acceleration, P, the local pressure and Γ 1 the adiabatic exponent. ∇ = d ln T/d ln P is the thermal gradient, ∇ ad the adiabatic gradient, and ∇ µ = d ln µ/d ln P, where µ is the mean molecular weight. This last term shows the explicit contribution of the mean molecular weight gradient to the Brunt-Vaïssälä frequency, and hence to the buoyancy radius. The integral in Eq. (3) is computed over the g-modes cavity. We get back to this in Sect. 4.2.
In the convective core, the efficient convection induces an almost adiabatic stratification and mixes the material, N 2 is therefore equal to 0. Above the convective core, the N 2 profile depends on whether the star is on a core-contracting or core-expanding main sequence. For the physics we have considered, above around 1.3 M , the convective core retreats along the main sequence, and leaves behind, in the radiative zone, a region of strong molecular weight gradient ∇µ. The buoyancy radius strongly depends on both the extent of the convective core, and the chemical stratification above the core. As a result, it is a good indicator of evolution on the main sequence of γ Dor stars. The evolution of a star's P 0 is then straightforward to understand: along the main sequence evolution, the chemical abundance contrast increases at the edge of the core, thereby increasing the mean molecular weight gradient ∇µ, and the convective core contracts, therefore increasing the interval where N 2 is non-zero, which induces a monotonously increasing peak in the Brunt-Väissälä frequency.
The behaviour of P 0 is illustrated in Fig. 2, where it is plotted against the central hydrogen abundance, for models with varying masses (left) or varying metalicities (right). As expected, both a higher metallicity and a higher mass result in larger convective cores, and therefore a higher buoyancy radius at a given age. This degeneracy of mass and metallicity has to be carefully accounted for when comparing models with observations (see Sect. 5).
The wiggles in the buoyancy radius in Fig. 2 are likely due to the numerical implementations of the criterium for convective stability, and its sensitivity to numerical parameters (e.g. radial resolution). This is a numerical issue often encountered in stellar evolution codes, and which does not change significantly the buoyancy radius P 0 .

Seismic measurement of the buoyancy radius and the near-core rotation
Conveniently, by making reasonable approximations, P 0 can be determined from γ Dor stars g-modes spectra. As they pulsate in high radial order g-modes, their pulsations are located in the asymptotic regime of g-modes. Without rotation, their spectra can be well approximated by the first order asymptotic theory of Tassoul (1980), which predicts that the periods of oscillation can be approximated at first order as: where n is the radial order, the angular degree, m the azimutal order, and is nearly constant. For γ Dor stars, which are moderate to fast rotators, rotation can be accounted for through the traditional approximation of rotation (TAR) to a certain extent. The TAR assumes that a star is spherically symmetric and that the latitudinal component of the rotation vector in the Coriolis force can be neglected. We refer to Unno et al. (1989) for a complete mathematical derivation. The first authors to apply the TAR in the stellar case were Berthomieu et al. (1978). The principle is that under such assumptions, and assuming solid body rotation, the equation system of stellar oscillations including rotation is separable in terms of a radial component and the Hough functions. The Eigenvalues of the Hough functions are the λ ,m (s) functions which depend on the angular degree, the azimuthal order and the spin factor s = 2P co n, ,m /P rot . In the asymptotic regime, the TAR allows us to express the period of modes in the co-rotating frame, including the effect of rotation, at first order, as: In Christophe et al. (2018), the authors derived a new method to determine the rotation period P rot and the buoyancy radius P 0 from Kepler observations of g-modes series in γ Dor stars. The principle is to use a set of trial rotation periods, that allows us to change from the inertial frame to the co-rotating frame of reference, and then stretch the periodogram by multiplying it by the corresponding λ ,m (s) functions, and search for periodicities in the stretched periodogram by the mean of a direct Fourier transform. The systematic errors due to uncertainties in modes frequency determination are calculated by propagating them by means of a Monte-Carlo simulation. The authors also assess the biases induced by the use of the traditional approximation, by testing it against complete calculations performed with the acor code (Ouazzani et al. 2012(Ouazzani et al. , 2015. These inaccuracies are difficult to assess directly for each measurement, because they require the knowledge of the true stellar structure. They are therefore determined with the help of hare and hounds exercises on synthetic oscillation spectra -calculated with the non-perturbative method -and amounts to maximum 8.2% for the buoyancy radius and 16.8% for the rotation frequency in the worst case. This method was used in order to determine the buoyancy radii and the rotation periods for the four stars previously presented in Ouazzani et al. (2017). We complete the sample with A121, page 4 of 11 The X c /X c,init axis has been inverted to emphasize the behaviour with evolution. The models were computed with disc locking during 5 Myrs, to a disc rotating with a period of 7.2 days.  Kurtz et al. (2014) and Saio et al. (2015) are included in the sample. The results are given in a Kiel diagram Fig. 3, where we show that the stars in the sample cover the whole γ Dor stars instability strip. The distribution in terms of near-core rotation frequency (ν rot ) seem to indicate that slow core-rotators lie all over the instability strip, while the faster ones occur more in the young-and-less-massive area of the diagram.
These measurements are then compared with evolutionary models, for which we have computed P 0 using the integral expression (Eq. (3)). The rotation periods, on the one hand, are derived using the method of Christophe et al. (2018, hereafter SC18) based on the TAR that considers the star as a rigid rotator. In principle, we assume that it is equivalent to considering a rotation period averaged over the cavity of the observed modes, given that among a series of high-order g-modes the cavity extent varies little.
On the other hand, the models result from evolutionary calculations including transport of AM by the processes defined in Sect. 2. Transport by shear induced turbulence and meridional circulation Z92 give a shellular rotation profile. Rather than a central rotation period, one should consider the near-core rotation period for comparison with observations. In order to calculate these near-core rotation periods, or angular velocities, we choose to take the average velocity over the g-modes cavity, weighted by the Brunt-Väissälä frequency: In order to determine the integral boundaries in Eq.
(3) and Eq. (7), the extent of g-modes cavity has to be determined in each of the numerous computed stellar models. Our concern was to do so, without having to calculate the oscillation spectra at each evolutionary time step. Non-adiabatic calculations in γ Dor stars (Dupret et al. 2005;Bouabid et al. 2013) showed that the excited modes had their radial orders comprised between around −15 and −50. The lower boundary (−15) is set by the necessity of having the oscillation time-scale close to the convective turnover time-scale at the bottom of the convective envelope, while the higher boundary (−50) is set by the fact that high radial order modes are radiatively damped. Considering that the modes cavity vary little over the range of excited radial order, the g-modes cavity is determined as being the cavity for the n = −25 mode, which frequency is computed by the use of the Tassoul (1980) asymptotic formula, and the effect of rotation being accounted for by a Ledoux (1951) splitting. Then, the g-modes cavity is defined as the interval where the thereby determined n = −25 mode frequency is lower than the Brunt-Väissälä frequency and the Lamb frequency. We consider that this is close enough to an average on g-modes cavity, for a first attempt to compare with observations.

Sample of observed stars
We have gathered a total of 37 γ Dor stars, which have been observed with Kepler during the nominal mission. Because we want global parameters determined as precisely as possible, we have built up from the sample of spectroscopically observed stars in Van Reeth et al. (2015), to which we added the four stars mentionned in Ouazzani et al. (2017). For these last four stars, we have relied on the Kepler input catalogue for the atmospheric

Comparison with Van Reeth et al. (2016)
The Christophe et al. (2018) method, used to derive the values of buoyancy radius and internal rotation of this work, consists of recovering the equidistance of period spacings by stretching the pulsation periods following the asymptotic formulation of the traditional approximation (Eq. (6)) while simultaneously inferring the rotation frequency and buoyancy radius. While suffering of the limitations inherent to the use of the asymptotic formulation of the traditionnal approximation, extensively reported in SC18, one of the main advantages of this method is that it is stellar-model-independent.
The first comparison with the method developed by VR16 was made in SC18 for the star KIC 12066947. Here we provide a comparison for 29 stars in our sample, showing dipolar prograde g-modes series in their seismic spectra, stars which where part of the VR16 sample originally. The comparison is presented in Fig. 5. Both methods agree very well on the determination of the internal rotation frequency (Fig. 5, right). The differences between the two methods arise more drastically on the buoyancy radius (Fig. 5, left). In particular, for some stars (KIC 2710594, KIC 3448365, KIC 5708550, KIC 6468987, KIC 110099031 and KIC 110080103) VR16 seem to find much higher P 0 than SC18. Inspecting these stars g-modes in the ∆P versus P diagram, we find that they all present wiggly ridges. In that case, peak selection has a significant impact in the final P 0 value. The method used here for modes extraction (explained in SC18, Sect. 4.1) is more conservative than the one used for VR16. As a consequence, in our analysis, we include systematically fewer modes than VR16, but with higher confidence level. We note that in their mode selection process, after having selected a first set of modes in the spectrum, VR16 use models of period spacings including the effect of rotation in order to select additionnal modes. The additionnal peaks, selected in a very dense area of the spectrum, follow the trend imposed by the period spacing model. As a consequence, the addition of a handful of modes modifies the result noticeably for the buoyancy radius, without impacting much the rotation determination.
One should also mention that, in VR16, the buoyancy radii are taken on a grid of model, whereas our determinations are stellar-model-independent. In any case, it is worth noticing that the internal errors from SC18 are almost systematically lower than those from VR16.

Rotational transport of angular momentum
We compare the observations to the evolutionary models including angular momentum transport as described in Sect. 2 in a near-core rotation versus buoyancy radius plane. The results are given in Fig. 8. Stars with detectable binary companions were excluded from the sample. We compared rotational evolution on the main sequence for three typical masses of γ Dor stars: 1.  Fig. 8 (left) the evolutionary models account for transport of angular momentum via meridional circulation and shear-induced turbulence according to the formalism of Z92, as explained in Sect. 2.
The measured rotation frequencies are spread over a large range: from 0.8 µHz for KIC 9751996, to 26.1 µHz for KIC 7365537. But in terms of buoyancy radii, they only cover the higher half of the interval covered by the models. However, this can be explained by exploring the evolutionary tracks: the lower half of the buoyancy radius interval corresponds to the second contraction evolutionary stage, when the stars have exhausted their hydrogen reservoir, and gravitational contraction takes over. In other words, this stage lasts a very short amount of time at the end of the main sequence. For instance for a 1.8 M model, this represents around 30 Myr, that is 2.5 % of the total main-sequence duration. Statistically, there is much less chance to observe a star in that stage, than in the main sequence.
On each panel of Fig. 8, for each mass and initial condition, we have plotted two curves, corresponding to the ±1σ of the metallicity distribution. The lower the metallicity, the higher the buoyancy radius is. Indeed, a lower metallicity induces a lower opacity, which impacts the radiative gradient, and therefore the limit of the convective core. The resulting increase of the extent of the convective core, generates a decrease of the integral in Eq. (3), leading to an increase of the buoyancy radius.  Figure 6 shows the area in the P 0 − ν rot plane covered by the models in Fig. 8 (left) with parameters given above, and for the initial conditions set in Sect. 3. We assume that they should be representative of the population of γ Doradus stars.
Compared to the models, the measured rotation rates seem to be globally shifted towards lower rotation. If the models agreed perfectly with the observations, we would have 40% of the data points within the shaded area, and 60 % of them evenly distributed above and below the shaded area in Fig. 6. This is certainly not the case.
In the sample studied here, there is no observed γ Dor star with rotation frequency higher than about ∼26 µHz, while all the models corresponding to the fast initial conditions (P disc = 2.4 days and τ disc = 3 Myr) give such rotation values or higher. The stars which PMS progenitors were locked to their disc at 2.4 days are either not observed, or seem to be braked by a missing mechanism. Even for the low-mass stellar models, with a distribution of rotation periods constrained by hPer observations at 13 Myr, do not fit the observations later, on the main sequence. On the low-rotation end, even with the slowest initial conditions, the models cannot fit the observed slow rotators: 15 stars out of 37, that is 41%, fall out of the rotation interval covered by the models. 68% of the targets have metallicities within the ±1σ interval. Therefore the spread in metallicity cannot explain the discrepancy between the rotation interval covered by the models and the observations.

Effect of convective overshooting
We tested the impact of core overshooting on the results. The overshooting is formalized by simple instantaneous mixing over a distance of α ov H P above the fully convective core, where H P is the pressure scale height. The choice of α ov is motivated by the callibration performed by Deheuvels et al. (2016) on eight solarlike stars with convective cores, using seismic inferences from Kepler observations. They found these stars seismic properties to be compatible with values of the overshooting parameter ranging from 0.07 to 0.18. Hence, we have tested the impact of convective overshooting using an intermediate value of α ov = 0.15. We calculated one evolutionnary sequence for a 1.6 M star cor- responding to the fast initial conditions (P disc = 2.4 days, and τ disc = 3 Myr, see Sect. 3) and low metallicity ([M/H] = −0.06), and one corresponding to the slow initial conditions (P disc = 7.2 days, and τ disc = 5 Myr) and high metallicity ([M/H] = +0.028). Both these sequences have been computed with a core overshooting value of α ov = 0.15 (see Fig. 7, in purple). For comparison, the same sequences without overshooting have also been plotted in Fig. 7 (blue).
The main structural effect of adding overshooting above the convective core is to increase the size of the core. This has a twofold effect in the P 0 − ν rot plane. It reduces the interval on which the Brunt-Vaissala frequency N is non-zero, hence decreasing the integral which intervenes in the expression (Eq. (3)), hence resulting in an increase of the buoyancy radius, as seen on Fig. 7. Similarly, it reduces the extent of the g-modes cavity in radius, and therefore, reduces the part of rotation profile probed by the modes (see Eq. (7)).
If the models with overshooting reproduced perfectly the observations, we should have a data set evenly distributed above and below the set of two yellow lines in Fig. 7. This is clearly not the case. Even if overshooting improves the agreement, we still are unable to explain the fast rotators desert seen in Fig. 6. On the low rotation end, while the effect on the buoyancy radius is noticeable, the impact on the rotation frequency is almost negligible, and certainly does not allow us to explain the slow rotators accumulation.

Enforcing rigid rotation
We have tested a different angular momentum transport hypothesis: solid body rotation all along evolution. This is not physically motivated, but it can be considered as a limiting case of instantaneous and highly efficient transport inside the star. The corresponding evolutionary models are confronted to observations in Fig. 8 (right). Here again, very few observed stars (and only for the models with M = 1.8 M ) fall in the area corresponding to the fastest initial conditions. However, the comparison is more favourable than in the case of hydrodynamical angular momentum transport (Fig. 8, left). The real difference between the results obtained for the two transport prescriptions is at the low-rotation frequency end. The slow rotators (around ∼5 µHz) are much better reproduced by the models relying on the solid body assumption.
This suggests the existence of a mechanism transporting angular momentum from the core to the outer part of the star. In this study, we assumed that these stars cannot sustain generation of magnetic field through a solar-like dynamo process, therefore models do not loose AM at the surface. As a result, in the presented models, the AM transport mechanism rigidifies the rotation profile. However, the absence of constraint on the envelope rotation in these stars prevents any strong conclusions on the profile per say.
In the next subsections, we estimate the efficiency of the missing mechanism in a similar spirit as in Eggenberger et al. (2017) for the red-giant phase, in order to characterize the missing AM transport processes.

Uncertainties on the horizontal turbulent viscosity
We first attempt to extract AM from the core while remaining within the Z92 framework, by exploring the uncertainties in turbulent viscosity coefficients. It has been long known that large uncertainties remain concerning the turbulent coefficients in the formalism (see for instance , and references therein). In particular, the prescriptions for the horizontal coefficient of viscosity, ν h , can differ by as much as two orders of magnitude. The results given in Sect. 5.1 are obtained with the  prescription. Here we choose to multiply ν h by a factor of 10 2 . Such an increase of ν h results in an enhancement of the meridional circulation by ensuring attenuation of the horizontal variation of the mean molecular weight (Marques et al. 2013). While it indirectly impacts the transport of chemical elements, the effect on the chemical stratification is marginal. The main effect is the enhancement of advective transport of AM through the meridional circulation.
The results of such calculations are given in Fig. 7 (yellow solid lines), by comparison with the unchanged Z92 models (purple solid lines), and the solidly rotating models (blue solid lines). This figure shows that amplifying ν h by two orders of Fig. 7. Near core rotation rates as a function of the buoyancy radius.
The observed values are given by the black circles. The solid lines stand for models of angular momentum transport with different prescription: Z92 without overshooting in blue (cf. Sect. 5.1), with overshooting of α ov = 0.15 (purple), with enhanced horizontal viscosity in yellow (cf. Sect. 5.4), with enhanced vertical viscosity in green (cf. Sect. 5.5), and with solid body enforced along evolution in red (cf. Sect. 5.3). These models have been computed for two models of 1.6 M , one (top) corresponding to the fast rotation initial conditions given in Sect. 3, and for the lowest metallicity, and the other one (bottom) to slow initial conditions and high metallicity. magnitude radically increases the core-to-envelope coupling to such an extent that it is equivalent to simply enforcing solid body rotation all along evolution.
The two sets of lines in Fig. 7 correspond to two combinations of metallicity and initial conditions for rotation: [M/H] = −0.06, τ disc = 3 Myr, P disc = 2.4 days for the upper one, and [M/H] = +0.028, τ disc = 5 Myr, P disc = 7.2 days for the lower one. They delimit the area where we expect to find at least 68% of the observations if the prescriptions for transport of AM are correct. Figure 7 shows that there is a clear improvement at the lower rotation end of the diagram, where an increased ν h allows a much better agreement with observations. Indeed, 80% of the measured rotation values lie in the possible rotationversus-buoyancy radius area. However, at the high rotation end, such an increased transport does not allow us to slow down the models with the fast initial conditions.

Additional vertical diffusion of AM
The rotational transport of angular momentum, as formalized by Z92 does not seem to satisfactorily reproduce the observed core rotation in γ Doradus stars, particularly in the fast rotation regime. Hence, additional AM transport is suspected to operate. Because the physical nature of the missing mechanism remains unclear, we choose to model its impact as an additional vertically diffusive process, which efficiency, function of the vertical viscosity coefficient ν ν,add , has to be higher than provided by the standard Z92 processes. In the transport equation, it consists of adding ν ν,add to ν ν in the second term of the right hand side of Eq. (1), without modifying D ν in Eq. (2).
Different values for ν ν,add are tested. For this additional viscosity to have an impact on the internal transport, it has to be stronger than ν ν . Hence the process starts to be efficient for approximately 1 × 10 6 cm 2 s −1 . Higher values of ν ν,add are then tested as high In each panel, the three sets of models relate to the different initial conditions set up in Sect. 3: from top to bottom: disc locking at a period of 2.4 days during 3 Myrs, 3.9 days during 3 Myrs, and 7.2 days during 5 Myrs. Left: evolution of the angular momentum distribution has been calculated with transport by meridional circulation and shears (Z92). Right: solid-body rotation is enforced all along evolution. The x-axis has been inverted so that the behaviour with evolution goes from left to right. as 1 × 10 8 cm 2 s −1 , the viscosity for which the evolution of the near core rotation overlaps with the one obtained assuming rigid rotation. This is illustrated Fig. 7 in green solid lines.
The value of ν ν,add which allows us to reproduce more than 64% of the observed rotations is found to be approximately 5 × 10 6 cm 2 s −1 . In other words, in this case as in the previous sub-section, this additional process allows us to spin down the core rotation for the models corresponding to slow initial conditions sufficiently to reproduce the observed slowest rotators satisfactorily. However, as with horizontal viscosity (Sect. 5.4), the models with fast initial conditions are not slowed down enough to be shifted into the observed regime.

Discussion and conclusions
We have measured near-core rotation rates of a sample of 37 γ Doradus stars on the main sequence, in order to test angular momentum transport models for intermediate-mass stars. We calculated rotational evolution of stellar models with masses typical of these stars from the PMS to the TAMS. To set the initial conditions for such evolution, we constrained simple disc-locking models to rotation distribution in young stellar clusters.
After being released from the disc, the angular momentum evolution of the stellar model is then dictated by the prescription for transport. The aim of this study is to use measurements of near core-rotation in γ Dor stars in order to constrain such prescriptions.
Before confrontating evolutionary calculations to observations, it is necessary to define an evolution indicator. To that end, we have explored the dependence of the buoyancy radius P 0 (Eq. (3)) on stellar parameters, and showed that P 0 could play that role, provided the degeneracy in mass and metallicity is carefully handled.
The angular momentum transport by meridional circulation and shear-induced turbulence (as described by Z92) was first explored. Figure 6 summarizes the results obtained: such models do not allow us to reproduce the near-core rotation measured by g-modes period spacing in our sample of γ Dor stars. Not only does it fail to reach the particularly slow cores of 15 stars of the sample (slow rotators accumulation Fig. 6), but it also fails at slowing down models of stars which progenitors are the 20% of fast rotators observed in young stellar clusters (fast rotators desert Fig. 6). Considering these arguments, and as for the Sun and red giant stars, it is safe to assert that the transport of angular momentum as formalized by Z92 cannot explain the measurements of near-core rotation in main-sequence intermediate-mass stars we have at hand.
Before investigating other possible angular momentum transport processes, we evaluated the impact of core overshooting on the results obtained with the Zahn (1992)'s prescription for transport. By increasing the extent of the mixed region above the convective core, the effect of overshooting consists mainly in an increase of the buoyancy radius. Hence, while it only slightly reduces the disagreement on the fast rotators desert, the addintion of core overshooting certainly does not explain the slow rotators accumulation. The dependence on the mixing length parameter value has also been explored (taking α = 1.6−1.8, as suggested in Trampedach et al. 2014, for the relevant area of the HR diagram), only to reach the same conclusions.
Furthermore, we compared the rotation measurements to models where solid-body rotation is enforced accross evolution, mimicking the effect of an efficient, yet to be identified mechansism. But this hypothesis allows us to slow down the models corresponding to slow initial conditions so that more than 80% of the observed rotation velocites are recovered. The good agreement between these models and the observations, in the low rotation end, show that, as for red giant stars, in intermediate-mass stars on the main sequence, we also need a mechanism which efficiently transports angular momentum from the core of the star outward. However, the impact on the models evolved from fast initial conditions is marginal and does not allow us to explain the lack of observations in the high rotation regime, the so-called fast rotators desert (Fig. 6). Whether the thereby extracted angular momentum is preserved in the envelope, and results in rigid rotation, or is ejected outside the star through winds is not clear, and the data at hand cannot allow us to solve this issue.
Trying to characterize further such AM extraction mechanism, we have shown that a rotational evolution similar to the one with solid body rotation could be obtained in two different ways. First, remaining in the framework of Z92, such a transport can be mimicked by advective AM transport through meridional circulation enhanced by additional horizontal viscosity (ν h × 100). Second, extraction of AM to that extent can be reached by adding vertical viscous transport with a viscosity of ν ν,add = 1 × 10 8 cm 2 s −1 . Although both transport processes allow us to reach the slow core rotations observed in γ Dor stars, they cannot explain the lack of measurement in the high rotation regime.
These results must be taken with a word of caution. The observed sample used for comparison with models here is mainly composed of Kepler γ Dor stars in the VR16 sample. The distribution of projected rotation velocities in this sample can reach up to about 200 km s −1 with an average around ν sin i = 60 km s −1 . Royer (2009) compiled a large sample of stars with projected rotational velocities from the literature, and sorted them by spectral type. They found that for the spectral types covered by the know γ Dor stars (A3 to F3), the average velocities per spectral type are between 120 km s −1 for the stars of spectral type A3, to 50 km s −1 for the stars of spectral type F3. This can be explained by the fact that the sample was first gathered by Tkachenko et al. (2013) for spectroscopic analysis, avoiding too large line profile broadenning caused for high values of ν sin i. This shows that the sample of stars used in the current study might suffer selection biases, and should be further investigated for completeness.
However, there seem to be a lack of AM transport process in the evolutionary models presented here, which is at work during the PMS and the main-sequence phase, and which can be efficient enough to slow down the core of intermediate-mass stars.
In the evolutionnary models calculated in this work, braking through coupling to a magnetized wind has been neglected. The magnetic braking models available in 1D stellar evolution codes are much better suited for solar like stars. When stars are slightly more massive, and their external convective zones become much thinner than that of the Sun, braking laws callibrated on the solar case (Kawaler 1988) cease to be applicable. Even more sophisticated formalisms such that of Matt et al. (2015) give results, in the considered mass range, which are barely physical. We expect the stars of our sample to have masses between 1.3 and 1.9 solar masses. Therefore, some might undergo magnetic braking while others might not. The criteria for this distinction, as well as a well suited formalism which would allow us to account for this extraction of angular momentum will be investigated in future works.
However, we expect a significant part of the sample to be too massive to have a significant external convective zone. Hence, in absence of a dynamo-generated magnetic field, transport by IGWs, particularly generated by the convective core, should be explored as a possible explanation for the observed fast rotators desert (Fig. 6). Such transport could operate either during the PMS, or on the main sequence. Refining the developments of Talon & Charbonnel (2005), Charbonnel et al. (2013) investigated the effect of the transport by IGWs generated from both convective zones: from the core (for models with masses higher than about 1.6 M ) and from the envelope of PMS low and intermediate-mass stars. They have shown that IGWs could drastically reshape the internal rotation profile of these stars, in particular causing a spin down of the central rotation, when arriving to the ZAMS. Encouraged by the first rotation measurements for γ Doradus stars and slowly pulsating B stars in CoRoT and Kepler data, Rogers (2015) and Townsend et al. (2018) have A121, page 10 of 11 investigated the role of IGWs and transport by g-modes during the main-sequence for intermediate-mass stars, and have shown that those might play an important role in shaping their internal rotation profiles. However, those preliminary works rely on several assumptions which degrades their realism, and therefore do not allow us to conclude quantitatively on the matter of interest here, and the issue of a missing mechanism remains open.
While enforcing the core spin down, such a mechanism would not necessarily induce rigid rotation. Hence, on the asteroseismic side, additional constraints on the rotation radial differential rotation would allow us to put stronger contraints on the missing mechanism. Recent attempts to measure internal differential rotation (Van Reeth et al. 2018, and references therein) have resulted in finding a dozen of γ Doradus stars with rigid rotation profiles. Further efforts in that direction are needed in order to provide statistical significance and allow us to provide robust constraints for angular momentum transport processes.