© ESO, 2018

1 Introduction

The angular momentum evolution of young low-mass stars has been investigated for several decades (e.g. Weber & Davis 1967; Skumanich 1972; Kawaler 1988; Keppens et al. 1995; Bouvier 2008; Reiners & Mohanty 2012). The associated theoretical efforts led to a better understanding of the main physical mechanisms at work in this evolution (star–disk interaction, magnetic braking, and internal redistribution of angular momentum; see, e.g., Bouvier et al. 2014, and references therein). Strong theoretical constraints have been added to the processes that drive angular momentum transport in stellar interiors (e.g. Amard et al. 2016) and the extraction of angular momentum by magnetized stellar winds (Matt et al. 2015; Réville et al. 2015; See et al. 2017). The physics that controls the angular velocity evolution of low-mass stars from the pre-main sequence (hereafter PMS) up to the end of the main sequence (hereafter MS) is thus now relatively well understood (see, e.g., Matt & Pudritz 2005; Gallet & Bouvier 2015; Somers & Pinsonneault 2015; Lanzafame & Spada 2015; Johnstone et al. 2015; Amard et al. 2016; Sadeghi Ardestani et al. 2017). Despite the uncertainties about the correct physical description to use to describe the mechanisms involved in the internal transport of angular momentum between the radiative core and the convective envelope, current models grasp the main trends of stellar rotational evolution.

In parallel to these theoretical developments, numerous exoplanets have been detected since 1995 (Mayor & Queloz 1995) and now reach a number of confirmed objects between 2950 and 3786 (June 5 2018)¹. These exoplanets can be found in a wide range of star-planet configurations that encompass a large distribution of planetary masses ranging from 10⁻⁴ to 100 M_Jup, orbital periods from 10⁻¹ to 10⁵ days, and a (sub)stellar mass ranging from 2 × 10⁻² to 4 M_⊙. Nevertheless, most angular momentum evolution models mainly focus on isolated stars, thus neglecting the possible impact of a planetary companion on the rotational evolution of the central star. However, the presence of close-in planets ought to be included in such numerical codes as pointed out by Bolmont et al. (2012), Zhang & Penev (2014), Lanza & Mathis (2016), Privitera et al. (2016a), and Rao et al. (2018), who show the strong impact of the stellar rotational history on the orbital evolution of massive close-in planets, and vice versa.

During the last decade, thanks to the inauguration of the Kepler satellite and most recently because of Kepler’s second life mission K2, we have entered a new era of improved rotational period measurements that allows for advanced astrophysical quests. Indeed, the precision of measured stellar surface rotation periods (through photometric variation induced by the presence of magnetic stellar spots) is now good enough to detect specific features in the rotation period distribution of open clusters. This is, for instance, the case for the rotation period distribution of the Pleiades cluster (a 120-Myr-old MS open cluster located at about 140 pc from the Earth) that has been analysed by Rebull et al. (2016) and Stauffer et al. (2016) using Kepler-K2 (Howell et al. 2014). In this cluster, they found K-type starswith a faster rotation rate than expected from current theoretical angular velocity evolution tracks (hence producing a “kink” in the rotational distribution). These “classical” rotational tracks are produced by numerical models that only invoke star–disk interaction, angular momentum extraction by stellar wind, and internal redistribution of angular momentum within the stellar interior; but this anomaly could potentially result from the presence of an exoplanet that may affect, through tidal interaction, the surface rotation rate of its host star. Indeed, Mathis (2015b; see also Bolmont & Mathis 2016, Gallet et al. 2017a) showed that the tidal dissipation inside the star is maximum, during the early-MS phase, for early K-type stars due to a specific configuration of their internal structure.

In the theoretical framework, tides in stars can be described by two components: the equilibrium tide, which corresponds to the large-scale hydrostatic adjustment induced by the gravitational interaction between the star and its companion (of stellar or planetary nature; see Zahn 1966; Remus et al. 2012; Ogilvie 2013) and which is made up of a large-scale non-wavelike/equilibrium flow; and the dynamical tides, which correspond to the dissipation of tidal inertial waves (mechanical waves that are generated inside rotating fluid bodies)due to the turbulent friction in convective regions (Ogilvie & Lin 2007; Ivanov et al. 2013) and to thermal diffusion and breaking mechanisms acting on gravito-inertial waves (gravity waves influenced by the effect of rotation through the Coriolis acceleration) in radiative regions (e.g. Zahn 1975; Terquem et al. 1998; Barker & Ogilvie 2010). The dissipation of the dynamical tides thus strongly depends on the internal structure of the star (see, e.g., Chernov et al. 2013; Ogilvie 2014; Mathis 2015b; Gallet et al. 2017a).

Most of the studies dedicated to tidal star-planet interactions often assume solid body rotation for the whole star (e.g. Bolmont et al. 2012; Nordhaus & Spiegel 2013; Bolmont & Mathis 2016). However, recent works have allowed for stellar core-envelope decoupling so as to investigate the impact of the presence of a massive planet on the surface rotation of the star. From the literature, Zhang & Penev (2014) used constant tidal dissipation efficiencies along the stellar evolution; Penev et al. (2014, 2018) included core-envelope decoupling so as to add constraints on the evolution of the tidal dissipation; and Privitera et al. (2016a) focused on star-planet interaction during the red giant phase.

In this article, and in complement to the work of Bolmont & Mathis (2016), we investigate the impact of tidal dissipation evolution (controlled by the internal stellar structure during the PMS and by the surface rotation rate of the star during the MS) on the evolution of the rotation rate of the host star using a two-zone rotational model that allows for core-envelope decoupling. We explore the parameter space of star-planet systems, considering stellar mass, initial parameters (rotation, disk lifetime, and coupling timescale), planetary mass, and initial orbital distance to map the impact of star-planet interaction on the rotational evolution of low-mass stars. Rao et al. (2018) also recently studied the impact of the equilibrium and dynamical tides on the orbital evolution of massive close-in planets. In their work they focused on the initial conditions that affect the planetary survivability around stars more massive than 1.0 M_⊙, while we are more interested in how the surface rotation rate of the host star is modified by the star-planet tidal interaction. These two works are thus very complementary.

This paper is structured as follows. The numerical model used in this work is described in Sect. 2. In Sect. 3 we investigate the rotational evolution of low-mass stars in the presence of a close-in planet, and how it is impacted by the main star-planet parameters. We first study the case of a solar mass star in Sect. 3.1, and we study the evolution of a initial rotational distribution orbiting an early-K 0.8 M_⊙ type star in Sect. 3.2. Finally, we generalize these results to a broader mass range in Sect. 3.3. We finally compare the results of our simulations to the Pleiades data in Sect. 3.4 and conclude in Sect. 4.

2 Model

The numerical model used in this work combines the rotational evolution model of Gallet & Bouvier (2013, 2015) and the modified orbital evolution model of Bolmont et al. (2012; see Bolmont & Mathis 2016). The dissipation of the dynamical tide inside the star is treated in the convective envelope as in Gallet et al. (2017a), who followed Ogilvie (2013) and Mathis (2015b) and the stellar structure is from the stellar evolution code STAREVOL (see Amard et al. 2016, and references therein). We developed this code that combines these two numerical approaches so as to study, in a more realistic way through the addition of the decoupling between the radiative core and the convective envelope, the impact of the star-planet interaction on the surface rotation rate of low-mass stars. This code is specifically designed for stars between 0.3 and 1.2 M_⊙ and for binary star-planet systems consisting of one central object and one planet.

In the following, we recall the specificities of each part of our combined numerical model.

2.1 Stellar evolution

This work is based on a grid of non-rotating stellar models we computed with the code STAREVOL (see, e.g., Amard et al. 2016) for a range of initial masses between 0.5 and 1.0 M_⊙ at solar metallicity (Z = 0.0134; Asplund et al. 2009). The references for the basic input microphysics (equation of state, nuclear reactions, and opacities) can be found in Amard et al. (2016) and in Lagarde et al. (2012). The initial helium abundance and mixing length parameter are calibrated without atomic diffusion to reproduce a non-rotating Sun with respect to the solar mixture of Asplund et al. (2009) with a 10⁻⁵ precision for the luminosity and the radius at the age of the Sun. The corresponding mixing length parameter and initial helium abundance are α_MLT = 1.6267 and Y = 0.2689. The stellar structure provided by the STAREVOL code is the foundation of our model. It is used, on the one hand, to follow the evolution of the stellar rotation rate and, on the other hand, to estimate the tidal dissipation inside the stellar interior.

2.2 Tidal dissipation

Following the studies of Bolmont & Mathis (2016), Gallet et al. (2017a), and Bolmont et al. (2017), we compute the tidal evolution of planets using a model improved with respect to the classical orbital evolution models, which only take into account the equilibrium tide (e.g. Mignard 1979; Hut 1981; Bolmont et al. 2011, 2012, 2015). The equilibrium tide consists of a large-scale flow driven by the hydrostatic adjustment of the body due to the perturbing gravitational potential of the planet. For rotating bodies, inertial waves, which are driven by the Coriolis acceleration, can be excited for a certain range of excitation frequencies (ω ∈ [−2Ω_⋆, 2Ω_⋆], where $\omega \equiv 2\left(n-{\Omega }_{\star }\right)$ is the tidal frequency in the case of circular coplanar systems, n is the mean orbital motion, and Ω_⋆ is the surface rotation rate of the star) in the convective envelope of low-mass stars. Assuming a two-layer star in solid body rotation, Ogilvie (2013) derived the frequency-averaged tidal dissipation induced by the dynamical tide, constituted by tidal inertial waves, in the convective envelope. In the case of a coplanar star-planet system in which the orbit of the planet is circular, this dissipation is given by $$\begin{array}{lll}{\langle \mathcal{D}\rangle }_{\omega }\hfill & =\hfill & {\displaystyle {\int }_{-\infty }^{+\infty }\text{Im}}\left[k_2^2(\omega )\right]\frac{\text{d}\omega }{\omega }=\frac{100\pi }{63}{ϵ}^2\left(\frac{{\alpha }^5}{1-{\alpha }^5}\right){\left(1-\gamma \right)}^2\hfill \\ \hfill & \hfill & \times {\left(1-\alpha \right)}^4{\left(1+2\alpha +3{\alpha }^3+\frac{3}{2}{\alpha }^3\right)}^2\left[1+\left(\frac{1-\gamma }{\gamma }\right){\alpha }^3\right]\hfill \\ \hfill & \hfill & \times {\left[1+\frac{3}{2}\gamma +\frac{5}{2\gamma }\left(1+\frac{1}{2}\gamma -\frac{3}{2}{\gamma }^2\right){\alpha }^3-\frac{9}{4}\left(1-\gamma \right){\alpha }^5\right]}^{-2},\hfill \end{array}$$(1)

with $$\alpha =\frac{R_{\text{rad}}}{R_{\star }}, \beta =\frac{M_{\text{rad}}}{M_{\star }}, \gamma =\frac{{\alpha }^3(1-\beta )}{\beta (1-{\alpha }^3)}<1,ϵ=\frac{{\Omega }_{\star }}{{\Omega }_{\text{c}}},$$(2)

where Ω_c is the critical angular velocity of the star, $k_2^2$ is the Love number of degree 2 (corresponding to the quadrupolar component, l = 2 and M = 2, of the time-dependent tidal potential proportional to the spherical harmonic $Y_l^m$), R_rad and M_rad are the radius and mass of the radiative core, respectively, and R_⋆ and M_⋆ are the radius and mass of the whole star, respectively. When present, the convective envelope surrounds the radiative core and both regions are assumed to be homogeneous with respective average densities ρ_c and ρ_e. The ratio between the perturbation of the gravitational potential induced by the presence of the planetary companion and the tidal potential evaluated at the stellar surface is given by $k_2^2$. Its imaginary component $\text{Im}\left[k_2^2(\omega )\right]$ is a direct estimation of the tidal dissipation. This formalism is very convenient because it allows us to take into account the dependence of tidal dissipation on stellar structure and rotation along their evolution, filtering out the complex frequency dependence of the dissipation of tidal inertial waves (Ogilvie & Lin 2007). This constitutes a first step for the study of the secular tidal evolution of star-planet systems (we refer the reader to Mathis 2015b; Bolmont & Mathis 2016; Gallet et al. 2017b, for detailed discussion).

In order to compute the orbital evolution of close-in planets, we use the model introduced in Bolmont & Mathis (2016). The evolution of the semi-major axis a of a planet on a circular orbit is given by (Hansen 2010; Leconte et al. 2010; Bolmont et al. 2011, 2012) $$\frac{1}{a}\frac{\text{d}a}{\text{d}t}=-\frac{1}{T_{\star }}[1-\frac{{\Omega }_{\star }}{n}],$$(3)

where T_⋆ is the evolution timescale given in Eq. (6). It depends on the stellar equivalent structural quality factor $\overline{Q_{{\text{s}}^{\prime }}}$, which can be expressed in terms of ${\langle \mathcal{D}\rangle }_{\omega }$ as $${\langle \mathcal{D}\rangle }_{\omega }=\frac{3}{2\overline{Q^{\prime }}}={ϵ}^2\frac{3}{2\overline{Q_{s^{\prime }}}},$$(4)

where $\overline{Q^{\prime }}$ is the equivalent modified tidal quality factor as defined in Ogilvie & Lin (2007) and Mathis (2015a). The tabulated values of $\overline{Q^{\prime }}$ and $\overline{Q_{{\text{s}}^{\prime }}}$ used in this study can be found here².

The tidal torque exerted by the planet on the convective envelope of the star is given by (e.g. Bolmont & Mathis 2016; Gallet et al. 2017a) $${\Gamma }_{\text{tide}}=\frac{h}{2T_{\star }}\left[1-\frac{{\Omega }_{\star }}{n}\right],$$(5)

where h is the orbital angular momentum of the planet, n is its orbital frequency, and T_⋆ is an evolution timescale given by $$T_{\star }=\frac{2}{9}\frac{M_{\star }}{M_{\text{p}}(M_{\text{p}}+M_{\star })}\frac{a^8}{R_{\star }^5}\frac{\overline{Q_{{\text{s}}^{\prime }}}}{{\widehat{ϵ}}^2}\frac{|n-{\Omega }_{\star }|}{\mathcal{G}},$$(6)

which relies on the semi-major axis a of the planet, the mass (M_⋆) and radius (R_⋆) of the star, the mass M_p of the planet, the stellar equivalent structural tidal quality factor $\overline{Q_{{\text{s}}^{\prime }}}$ (see Gallet et al. 2017a; Bolmont et al. 2017, and references therein), and $\widehat{ϵ}=\left({\Omega }_{\star }/\sqrt{\mathcal{G}{M}_{\odot }/R_{\odot }^3}\right)\equiv {\Omega }_{\star }/{\Omega }_{\odot ,\text{c}}$ with Ω_⊙,c the critical angular velocity of the Sun. The stellar equivalent structural tidal quality factor $\overline{Q_{{\text{s}}^{\prime }}}$ strongly and primarily depends on the evolution of the stellar internal structure (i.e. the relative mass and size of the radiative core compared to the whole star). It evolves over five orders of magnitude during the stellar evolution (see Gallet et al. 2017a). For the equilibrium tide, we use a constant value for the dimensionless dissipation factor ${\overline{\sigma }}_{\star }$ (which is linked to $1/\overline{Q_{{\text{s}}^{\prime }}}$; see Bolmont & Mathis 2016)³, which only depends on the stellar mass and is calibrated on the observation (see Hansen 2012).

In this work, the tidal dissipation (and thus the tidal torque) in the radiative core is not included (Zahn 1975; Goodman & Dickson 1998; Terquem et al. 1998; Barker & Ogilvie 2010; Barker 2011; Ivanov et al. 2013). This additional torque could have a strong influence on the orbital evolution of close-in planets (it might make them fall significantly faster into the stars) and thus changes the way the surface rotation is impacted by the tidal interaction. The internal differential rotation inside the star has also a clear effect on the rotational evolution of the convective envelope (especially during the MS phase; see Gallet & Bouvier 2013, 2015). A change of rotation rate due to the tidal torque applied in the radiative core could then affect the rotation of the envelope. We also do not include the dissipation inside the planet itself. In this work, we do not consider massive close-in planets that were formed with the high-eccentricity scenario. This hypothesis is consistent with the fact that we neglect planetary tides that arise if the orbit is eccentric or if the rotation is not synchronized with zero obliquity. We also consider that the planet is alone (or at least that it does not feel any gravitational pull from another planet), which is also consistent with the disk migration ofa single massive planet. Indeed, this latter should have destabilized its close environment along its evolution.

Finally, in the case of a planetary engulfment, we assume that the whole angular momentum of the planet is instantaneously transferred to the star which end up in the rapid increase of the stellar surface angular velocity.

Fig. 1 Sketch of the stellar structure studied here. We consider a two-layered star: a radiative core of rotation rate Ωrad and a convective envelope of rotation rate Ωconv. Both are assumed to rotate as a solid body with different rotation rates. In red we indicate the different processes included in the model. The convective envelope is subjected to external torques (stellar wind and tidal torques), as well as torques due to the radiative core development and coupling with the convective envelope.

2.3 Stellar rotational evolution

We model the evolution of the rotation of the central star using the formalism of Gallet & Bouvier (2015) in which the star is assumed to be composed of two parts: a radiative core surrounded by a convective envelope (see Fig. 1). As in the averaged-dissipation tidal model, both regions are assumed to rotate as a solid body with different rotation rates. The stellar angular momentum is described by J_⋆ = I_⋆Ω_⋆, where $I_{\star }\propto M_{\star }{R}_{\star }^2$ is the moment of inertia. The total stellar angular momentum is thus defined by J_⋆ =J_rad + J_conv, where J_rad is the angular momentum of the radiative region, with R_rad its radius and Ω_rad its rotation rate, and J_conv is the angular momentum of the convective region and Ω_conv its rotation rate. The total angular momentum evolution rate is given by $$\frac{\text{d}{J}_{\star }}{\text{d}t}=\frac{\text{d}{I}_{\star }}{\text{d}t}{\Omega }_{\star }+I_{\star }\frac{\text{d}{\Omega }_{\star }}{\text{d}t}={\Gamma }_{\text{all}},$$(7)

where Γ_all is the sum of the external torques; here we account for the stellar wind Γ_wind and the tidal torque Γ_tide; we neglect other external torques such as the accretion torque Γ_acc and the star–disk interaction torque Γ_disk since our simulations start after proto-planetary disk dissipation. The magnetic star-planet interactions Γ_mag are also not included in this work (see Strugarek et al. 2017 for detailed discussion about the range of application of the magnetic star-planet torque). While the wind always acts to remove angular momentum (Γ_wind < 0), the tidal interaction can either spin up or spin down the star depending on whether the planet is located inside or outside of the corotation radius (e.g. Bolmont & Mathis 2016) $$R_{\text{co}}={\left(\frac{\mathcal{G}{M}_{\star }}{{\Omega }_{\star }^2}\right)}^{1/3}={\left(\frac{\mathcal{G}{M}_{\star }{P}_{\text{rot},\star }^2}{{(2\pi )}^2}\right)}^{1/3},$$(8)

where P_rot,⋆ is the surface rotation period of the host star, $\mathcal{G}$ the gravitational constant, M_⋆ the stellar mass, and Ω_⋆ = 2π∕P_rot,⋆ the surface angular velocity of the host star.

Following Gallet & Bouvier (2015), the evolution of the angular momentum rate of the convective envelope is given by $$\frac{\text{d}{J}_{\text{conv}}}{\text{d}t}={\Gamma }_{\text{wind}}+{\Gamma }_{\text{tide}}+{\Gamma }_{\text{c-e}}+{\Gamma }_{{\text{rad}}_{\text{evol}}},$$(9)

where Γ_wind, Γ_tide, Γ_c-e, and Γ_rad,evol are the torques that act on the convective envelope; Γ_c-e and Γ_rad,evol are torques applied on the convective envelope by the radiative core. Figure 1 shows the different torques exerted on the convective envelope of a low-mass star. For each of these components, we use the following prescriptions.

2.3.1 Wind torque Γ_wind

The wind braking torque is given by (Schatzman 1962; Weber & Davis 1967) $${\Gamma }_{\text{wind}}\propto K_1^2{\Omega }_{\star }{\dot{M}}_{\text{wind}}{r}_{\text{A}}^2.$$(10)

It depends on the stellar rotation rate Ω_⋆, mass-loss rate Ṁ_wind, and Alfvén radius r_A. We use here the prescription of Matt et al. (2012) for the Alfvén radius and the revised prescription of Cranmer & Saar (2011) for the mass-loss rate (see Gallet & Bouvier 2015 for details). In our numerical code, K₁ is a free parameter and is associated to the wind braking efficiency. It is set so as to reproduce the observed rotation rate of the present Sun and rotational distribution of late-MS cluster.

The stellar wind braking Γ_wind relies on two parameters: the mass-loss rate and the value of the Alfvénic radius r_A, which both depend on the evolution of the mean stellar magnetic field B_⋆ f_⋆ (see Cranmer & Saar 2011; Matt et al. 2012), where B_⋆ is the stellar magnetic field strength and f_⋆ is the magnetic filling factor (i.e. the fraction of the surface of the star that is magnetized). As in Gallet & Bouvier (2015), we follow Cranmer & Saar (2011) by assuming that the magnetic field of the star is at the thermal equilibrium with the stellar photosphere and can thus be expressed as a function of the equipartition magnetic field strength, and f_⋆ is a unique function of the Rossby number Ro = P_rot,⋆∕τ_conv, that is, the ratio of the rotation period P_rot,⋆ to the convective turnover timescale τ_conv. Finally, the mean magnetic field can be expressed as $$B_{\star }{f}_{\star }=1.13\sqrt{\frac{8\pi {\rho }_{\star }{k}_{\text{B}}{T}_{\text{eff}}}{\mu m_{\text{H}}}}\frac{0.55}{{\left[1+{(x/0.16)}^{2.3}\right]}^{1.22}},$$(11)

where x = Ro/Ro_⊙ and Ro_⊙ = 1.96, ρ_⋆ is the photospheric mass density, k_B is the Boltzmann’s constant, T_eff is the effective temperature, μ is the mean atomic weight, and m_H is the mass of a hydrogen atom (see Eq. (3) from Gallet & Bouvier 2015).

Observations (e.g. X-ray luminosity; Pizzolato et al. 2003, and magnetic flux density; Reiners & Basri 2010) show that the evolution of the magnetic field of low-mass stars reaches a plateau (hereafter the saturation or saturated regime) at a maximum value (see Fig. 6 of Reiners et al. 2009) when the rotation of the star exceeds a certain threshold ( ≈ 16 Ω_⊙ during the MS for a 1.0 M_⊙ star).

Above this limit the mass-loss rate and magnetic field strength of fast rotating stars become constant regardless of the evolution of the surface rotation rate. In the saturated regime, the extraction of angular momentum through the stellar wind scales as Γ_wind ∝Ω_⋆ (since B_⋆ and; see Eq. (10) and Fig. 4 from Gallet & Bouvier 2015) instead of scaling as ${\Omega }_{\star }^3$ (i.e. the empirical Skumanich 1972 relationship).

2.3.2 Core-envelope coupling torque Γ_c-e

As in MacGregor & Brenner (1991), we assume that both the core and the envelope are in solid body rotation with two different rotation rates and that a quantity ΔJ of angular momentum is transferred between the core and the envelope over a timescale τ_c-e. The torque Γ_c-e associated with the angular momentum transfer rate between the radiative core and the convective envelope is expressed as $${\Gamma }_{\text{c-e}}=\frac{\Delta J}{{\tau }_{\text{c-e}}},$$(12)

where ΔJ is the quantityof angular momentum to be transferred, instantaneously, between the two regions to obtain a uniform rotation (MacGregor & Brenner 1991): $$\begin{array}{lll}\Delta J\hfill & =\hfill & \frac{I_{\text{conv}}{J}_{\text{rad}}-I_{\text{rad}}{J}_{\text{conv}}}{I_{\text{rad}}+I_{\text{conv}}},\hfill \\ \hfill & =\hfill & \frac{I_{\text{conv}}{I}_{\text{rad}}}{I_{\text{rad}}+I_{\text{conv}}}\left({\Omega }_{\text{rad}}-{\Omega }_{\text{conv}}\right).\hfill \end{array}$$

The coupling timescale τ_c-e is a free parameter of our numerical code that is determined using the observed rotational period distribution of early-MS stellar clusters (see Gallet & Bouvier 2015). It is related to internal transport mechanisms (e.g. Maeder 2009; Mathis 2013, and references therein)

2.3.3 Radiative core torque Γ_rad,evol

Along the growth of the radiative core during the PMS phase, part of the convective envelope becomes radiative. Thus, we can consider that angular momentum is extracted from the envelope and transferred to the expanding radiative core. The associated torque Γ_rad,evol is (Allain 1998) $${\Gamma }_{\text{rad},\text{evol}}=\frac{2}{3}{R}_{\text{rad}}^2{\Omega }_{\text{conv}}\frac{\text{d}{M}_{\text{rad}}}{\text{d}t},$$(15)

with M_rad the mass of the radiative core.

2.4 Rotational evolution

Finally, the evolution of the rotation rate of the convective envelope and that of the radiative core are given by $$\begin{array}{lll}\hfill & \hfill & \frac{\text{d}{\Omega }_{\text{conv}}}{\text{d}t}=\frac{1}{I_{\text{conv}}}\left({\Gamma }_{\text{wind}}+{\Gamma }_{\text{tide}}+{\Gamma }_{\text{c-e}}+{\Gamma }_{\text{rad},\text{evol}}\right)-\frac{{\Omega }_{\text{conv}}}{I_{\text{conv}}}\frac{\text{d}{I}_{\text{conv}}}{\text{d}t},\hfill \\ \hfill & \hfill & \hfill \\ \hfill & \hfill & \frac{\text{d}{\Omega }_{\text{rad}}}{\text{d}t}=\frac{1}{I_{\text{rad}}}\left({\Gamma }_{\text{rad},\text{evol}}-{\Gamma }_{\text{c-e}}\right)-\frac{{\Omega }_{\text{rad}}}{I_{\text{rad}}}\frac{\text{d}{I}_{\text{rad}}}{\text{d}t}.\hfill \end{array}$$

The rotation rate of the convective envelope can be decomposed into five contributions: the contraction, the stellar wind, the tides, the core-envelope coupling, and the development of the radiative core. The rotation rate of the radiative core is controlled by only three terms: the contraction, the development of the radiative core, and the transport of angular momentum between the core and the envelope. We therefore expect the rotational evolution of the envelope to be quite complex.

Each run for a given star-planet configuration is thus characterized by a number of initial parameters. For the star, they are the initial rotation rate of the star P_{rot, init}, the coupling timescale between the radiative core and the convective envelope τ_c-e, the disk lifetime τ_disk, the wind braking efficiency K₁, and the stellar mass M_⋆. For the planet they are the initial semi-major axis a (hereafter SMA), the planetary mass M_p, and the initial time at which we consider the presence of the planet t_init.

Since we start our simulations at the end of the disk lifetime, t_init ≡ τ_disk and corresponds to the age at which the tidal interaction starts. The initial SMA is also taken as a fraction of R_co and is thus a function of the initial stellar rotation rate P_{rot, init}. There is a degeneracy between t_init ≡ τ_disk and P_{rot, init} since using a longer t_init for a given P_{rot, init} can be similar (in terms of rotational evolution) to a short t_init but with a small P_{rot, init}. However, each t_init is associated with only one stellar internal structure and tidal dissipation efficiency couple.

In the rest of this work, we adopt the parametrization extracted from Gallet & Bouvier (2015) that links the mass and initial rotation rate to the free parameters of our numerical code: $$\begin{array}{lll}\hfill & \hfill & {\tau }_{\text{c-e}}=P_{\text{rot},\text{init}}^{0.66}(\text{days})\times {10}^{2.06}\times 0.06{\left(\frac{M_{\star }}{M_{\odot }}\right)}^{-3.97}\text{Myr},\hfill \\ \hfill & \hfill & {\tau }_{\text{disk}}=P_{\text{rot},\text{init}}^{0.86}(\text{days})\times {10}^{0.17}\times {\left(\frac{M_{\star }}{M_{\odot }}\right)}^{1.55}\text{Myr},\hfill \\ \hfill & \hfill & K_{\text{1}}=26.6+22.6{\left(\frac{M_{\star }}{M_{\odot }}\right)}^2-47.5\frac{M_{\star }}{M_{\odot }}.\hfill \end{array}$$

3 Rotational evolution and star-planet systems

3.1 Case of a solar mass star

The reference case used for comparison is given in Table 1. We consider a one solar-mass star, with a rotation period of three days at 5 Myr (between a moderate and a fast rotation rate regarding the observed rotation period distribution of early PMS clusters; see Rodríguez-Ledesma et al. 2009), and a core-envelope coupling timescale of 0.5 Myr (this short timescale is chosen to mimic a solid body rotation where both convective and radiative regions havethe same rotation rate). We consider a close-in Jupiter-mass planet and investigate the angular velocity of the star under the action of planetary tidal migration. In this work we focus on initial SMA close to and within the corotation radius (see Eq. (8)), since type I and type II migrations should bring the planets up to the corotation in the accretion disk on a short timescale (0.01–0.1 Myr; see Masset et al. 2006).

3.1.1 Generalities about planetary orbital evolution

As described in Sect. 2.2, there are two components of the stellar tides:equilibrium and dynamical. While the equilibrium tide is always present during the whole stellar life, the dynamical tide is only triggered when P_orb > 1∕2 P_rot,⋆ (Bolmont & Mathis 2016). The orbital evolution of a given planet then depends on the value of the ratio P_orb∕P_rot,⋆, and in substanceon which dominant tide the planet will be subject to during its evolution. Moreover, for a given initial P_orb∕P_rot,⋆ configuration, a change of t_init leads to different stellar angular velocity and planetary orbital evolutions (since each t_init is associated with a single stellar internal structure–tidal dissipation properties couple). In this work, the impact of the planetary orbital evolution on the surface rotation rate of the host star is estimated using the quantity δP_rot = 1 −Ω_conv,⋆∕Ω_conv,isol. = 1 − P_rot,isol.∕P_rot,⋆; δP_rot corresponds to the departure of the rotational evolution of the planet host star (P_rot,⋆) compared to an isolated (i.e. without planet) star (P_rot,isol.).

Figure 2 shows the evolution of the SMA, δP_rot, and stellar rotation P_rot,⋆ for the reference model with different initial SMAs ranging from 0.8 R_co to 1.2 R_co, with R_co = 4.06 × 10⁻² au at t_init. During the first Myr after t_init, there is no planetary migration and the planet marginally affects the surface rotation rate of its host star. Outward and inward planetary migrations start to be observed at about 1 Myr after t_init, which corresponds to the age (i.e. around 6 Myr) at which the dissipation of the dynamical tide inside the stellar convective envelope is maximum in regard to the stellar internal structure (see Fig. 4 of Gallet et al. 2017a). The corresponding effect on the surface rotation rate of the star is shown in the middle panel of Fig. 2. Inward and outward migrations globally lead respectively to the acceleration (δP_rot < 0) and deceleration (δP_rot > 0) of the star compared to an isolated star.

During the PMS phase, when the star is contracting, the surface rotation rate of the star is marginally impacted (δP_rot < –0.1) either by the planetary inward migration or by the planetary engulfment (since J_orb∕J_⋆ ≈ 10⁻², with J_orb the orbital angular momentum; see Hut 1981).

On the MS, and because of the extraction of angular momentum by the stellar wind, the δP_rot of each configuration tends towards zero (i.e. the rotational convergence), which erases the knowledge of the presence of a massive planet from the rotational history of the star.

Finally, during the sub-giant phase and red-giant branch (hereafter RGB), all low-mass stars engulf their planet (because the planet either reaches the Roche lobe⁴ or the stellar radius). The surface rotation rate of the host stars is then strongly affected by this planetary engulfment leading to the rapid decrease of δP_rot (see Privitera et al. 2016a,b, for more details about the planetary engulfment during the RGB phase).

Fig. 2 Evolution of the orbital distance of a Jupiter-mass planet with different initial SMA (top panel) ranging from 0.8 Rco (black), 1.0 Rco (dark grey), and 1.2 Rco (light grey), δProt (middle panel), and the stellar rotation period (bottom panel) during the evolution of a 1.0 M⊙ star with an initial rotation period of three days and a coupling timescale τc-e of 0.5 Myr. Upper panel: the orbital distance of the planet is represented in full coloured lines. The upper dotted lines correspond to Porb = Prot,⋆ (the corotation limit) and the lower dotted lines correspond to Porb = 1/2 Prot,⋆ limit, which marks the region above which the dynamical tide operates. The red long-dashed line is the evolution of the Roche lobe radius. Middle panel: departure of the rotation rate of the host star (Prot,⋆ ) from an isolated star (Prot,isol.). The horizontal dashed line indicates no difference with an isolated star, and the vertical dashed line represents the localisation of the ZAMS (around 50 Myr for a 1.0 M⊙ star). Lower panel: rotation period (in days) of the host star. The horizontal dashed line correspond to the transition between saturated and unsaturated magnetic regimes (1.6 days ≈ 16 Ω⊙), and the circle represents the rotation rate of the present Sun. The time on the x-axis is given from an initial time tinit of 5 Myr.

Fig. 3 Effect of planetary companion mass on the rotational evolution of a 1.0 M⊙ mass star. The initial rotation period of the star is three days, the star is assumed to rotate as a solid body with a core-envelope coupling timescale of 0.5 Myr; SMAinit = 0.8 Rco ≈ 3.25 × 10−2 au. The masses of the planets are 10 M⊕, 1 MJup, 5 MJup, and 13 MJup, displayed as darkish to greyish coloured lines, respectively. In the upper panel, the Roche lobe limits are displayed as long-dashed lines and the radius of the star as red dash-dotted line. The horizontal dashed line in the lower panel corresponds to the transition between the saturated and unsaturated magnetic regimes, and the circle represents the surface rotation rate of the present Sun. The time on the x-axis is given from an initial time tinit of 5 Myr.

3.1.2 Exploration of the star-planet system parameters

The path followed by the surface rotation rate of the star impacts both the orbital evolution of the planet (through its impact on the corotation radius and on the dynamical and equilibrium tide limit) while it also has an impact on the difference δP_rot of the planet host’s star, compared to an isolated star, through the impact of orbital migration on the rotation of the star. Rotational evolution has thus a major role in star-planet interaction processes. In the following, we describe the impact of the star-planet system’s physical quantities on the rotational evolution of the star and its feedback effect on the orbital evolution of massive planets.

Impact of core-envelope decoupling timescale and initial time

In the literature (e.g. Bolmont et al. 2012; Bolmont & Mathis 2016), the study of the orbital evolution of massive planet and its effect on thestellar surface rotation rate is often treated using solid body rotation for the whole star (see the cases in Fig. 2 with a coupling timescale of 0.5 Myr). Increasing the decoupling between the core and the envelope of the star allows the convective envelope to be braked earlier during the PMS phase. As a consequence, the corotation radius increases (since $R_{\text{co}}\propto P_{\text{rot},\star }^{2/3}$; see Eq. (8)) and the tidal torque is reduced (since Γ_tide ∝|1 −Ω_⋆∕n|; see Eq. (5)). This allows the planet to survive longer, by 10–15 Myr compared to the solid body rotation case. Indeed, in the case of the stellar reference star and with a 1 M_jup planet initially located at 0.8 R_co ( ≈0.03 au), increasing the coupling timescale between the core and the envelope from 0.5 Myr (solid body) to 10 and 30 Myr (that correspond to the parametrization of the fast and slow rotators; see Gallet & Bouvier 2015) induces a shift in the planetary engulfment of ≈5 and 15 Myr, respectively.

Similarly, decreasing t_init ≡ τ_disk induces an earlier spin-up of the star because of the contraction. This allows the planet to survive the stellar evolution due to a decrease in R_co. In the case of the stellar reference star and with a 1 M_jup planet initially located at 0.8 R_co ( ≈0.03 au), increasing the disk lifetime from 2 to 5 Myr leads to the inward (engulfment) migration of the planet.

Impact of planetary mass

Figure 3 shows the implication of a change of the mass of the planet, ranging from 10 M_⊕ to 13 M_Jup, on the planetary orbital evolution and stellar rotational evolution. Here we assume the same configuration as in the reference case, except for the planetary mass, and we consider an initial SMA of 0.8 R_co≈ 3.25 × 10⁻² au. The effect of the planetary mass is indeed quite strong on both SMA and stellar surface rotation rate evolution.

For a given star-planet configuration, when the mass of the planet increases, the migration timescale decreases causing the planet to migrate earlier during the stellar evolution (see Eq. (6)). The mass of the planet also affects the strengthof the impact of the planetary migration on the surface rotation evolution of the host star. Higher mass means higher orbital angular momentum. At the Zero age main-sequence (hereafter ZAMS), this leads to δP_rot = –1.5 to –3 for the 13 M_Jup and to δP_rot = –0.5 for the 5 M_Jup. For planets less massive than 1 M_jup, no significant effect on the surface rotation of the star is seen up to the end of the MS phase, consistent with previous theoretical works (Bolmont & Mathis 2016) and observations (Ceillier et al. 2016).

Initial SMA inside the corotation radius

Finally, we fix the initial SMA of the planet at 0.03 au and explore the impact of the initial rotation rate of the star using one, three, six, and eight days. Figure 4 shows the evolution of the planetary SMA and stellar surface rotation rate for a Jupiter-mass planet orbiting a 1.0 M_⊙ mass star. For each initial rotation period, the initial SMA corresponds to a fraction of the corotation radius. With an initial SMA of 0.03 au, we have

• one day: SMA_init = 1.53 R_co,

• three days: SMA_init = 0.73 R_co,

• six days: SMA_init = 0.46 R_co,

• eight days: SMA_init = 0.38 R_co.

In Fig. 4, three different cases can be described. The first one is the case where the planet migrates outward already during the PMS. This happens because the initial rotation period of the star is sufficiently fast (one day) so that the initial value of the SMA (0.03 au) is located outside of the corotation radius. Because the planet is initially outside of the corotation radius and its evolution is determined by the dynamical tide, it is submitted to an efficient outward torque.

The second case is the one where the planet falls onto the stellar surface during the PMS phase. It appears when the initial rotation rate of the star is three days. In that configuration, the planet is initially located between 1∕2 P_rot,⋆ < P_orb < P_rot,⋆ and is thus dominated by the dynamical tides. In that case, the Jupiter-mass planet falls into the star during the PMS phase and does not significantly affect its rotation rate.

The lastcase is the one where the planet falls into the stellar surface during the MS phase. This happens when the initial rotation rate of the star is low (six to eight days). With such low initial rotation rate, the planet is initially inside of the 1∕2 P_rot,⋆ limit and is thus only subject to the equilibrium tide. Because the equilibrium tide is less efficient than the dynamical tide, the planet migrates inward on much longer timescales (of the order of 10⁹ years) and thus falls onto the stellar surface much later. For this case, the impact of the planet engulfment on the surface rotation rate of the star is quite important as it reaches a spin rate four times larger than that of an isolated star.

Fig. 4 Evolution of the orbital distance of a Jupiter-mass planet with an initial SMA of 0.03 au (top panel), δProt (middle panel), and the stellar rotation period (bottom panel) during the evolution of a 1.0 M⊙ star with four different initial rotation rates: one, three, six, and eight days from darkish to greyish coloured lines, respectively. In the upper panel, the dotted lines correspond to the Porb = 1∕2 Prot,⋆ limit, the red dashed line is the Roche lobe limit, and the red dashed-dotted line is the stellar radius.

3.2 Evolution of star-planet initial distributions

The most important parameters that are involved in the star-planet rotational–orbital evolution are the initial location SMA_init of the planet around the host star (expressed in R_co) and the mass of the planet M_p. To analyse the impact of the SMA_init–M_p space, we first evolve in time an initially uniform (in terms of surface rotation rate) distribution consisting of 0.8 M_⊙ solar metallicity stars with an initial rotation period ranging from 1 to 11 days (with a 0.2 days step).

3.2.1 Corotation radius

The greatest impacts of tidal interactions on surface rotation rate are located, for a 0.7–0.8 M_⊙ mass star, in a narrow range of initial rotation period around five to eight days. We consider here that around each star orbits a 2 M_jup mass planetwhose initial SMA is taken to be 10, 30, 40, 45, 50, 55, 70, and 100% R_co. Figure 5a shows each of thesedistributions for a 0.8 M_⊙ star at 120 Myr (the age of the Pleiades) as a function of the initial localization of the planet. During the MS phase, the rotational period distributions are the most impacted by the planet when it is located around 50 ± 5% of the corotation radius.

In the case of an initial SMA of 50% R_co, the planets are initially inside the equilibrium tide regime. Hence, they remain on the same orbit during the PMS phase since equilibrium tide is less efficient than the dynamical tide. However, while the rotation rate of the star increases (due to the stellar contraction), the corotation radius slowly moves close to the stellar surface. At some point during the end of the PMS phase, the planet crosses the limit 1/2 P_rot,⋆ and becomes sensitive to the dynamical tide and thus starts to migrate faster. The planet then falls into the stellar surface during the early-MS phase, and strongly impacts its surface rotation rate.

At 30 or70% R_co, the planet falls during the PMS which, as pointed out above, has almost no effect on the MS stellar rotational evolution. For the first case, it is becausethe planets are initially too close to the star and thus either fall into the stellar surface during the PMS phase, or during the MS phase at older ages than the age of the Pleiades, hence not creating any kink in the rotation period distribution at 120 Myr. For 70% R_co, the most close-in planets fall into the stellar surface during the PMS phase, the others migrate outward. The limit P_orb = 1∕2 P_rot,⋆ corresponds to 2^−2∕3 R_co ≈ 63% R_co.

3.2.2 Massive planet

We also explored the impact of the planetary mass on the opening of the rotational kink. Figure 6 shows the temporal evolution of the rotational distribution related to the 0.8 M_⊙ mass star in the case of a 1 M_jup (Fig. 6a) and 2 M_jup (Fig. 6b). In these figures, we can see a rotational kink that starts to be opened between 20 and 30 Myr up to the early MS phase in each case. The size of the rotational kink is more important when the planetary mass increases. The amplitude of the rotational kink is at its maximum during the early MS phase and is slowly reduced by the rotational convergence due to the magnetic braking (see Gallet & Bouvier 2015).

The tidal star-planet interaction is at its maximum for SMA_init = 0.5 R_co. In the case of a 0.8 M_⊙ star, only the intermediate rotators (P_rot,init = 5–8 days) are impacted by the tidal interactions, which implies that not all planets in our initial rotational distribution are involved in the production of the rotational kink.

Fig. 5 Distribution at 120 Myr in log Prot,⋆ obtained for a 0.8 M⊙ from an initial distribution of rotation period ranging from 2 to 11 days. Effect of the initial SMA (10, 30, 45, 50, 55, 70, and 100% Rco) in the case ofa 2 Mjup (panel a). Impact of the mass of a planet (between 1 and 5 Mjup) initially located at 50% Rco (panel b). These distributions are extracted at 120 Myr that is about the age of the Pleiades cluster.

Fig. 6 Temporal evolution of the rotation distribution related to the 0.8 M⊙ mass star. The mass of the planet is taken to be 1 Mjup (panel a) and 2 Mjup (panel b) and is initially located at 0.5 Rco. Each line correspond to one stellar configuration. The bottom grey line corresponds to an initial period of two days, the upper black line corresponds to an initial period of 11 days. The consecutive difference between each line is 0.2 day.

3.3 Exploration of stellar mass

We explore a larger range of stellar mass, namely 0.5-1.0 M_⊙, and consider a 2 M_jup planet so as to increase the rotational kink in log P_rot observed in Fig. 5 around log P_rot = 0.5 − 0.6. We then evolve the same initial distribution as described above for initial planetary SMA around 50% R_co. Table 2 shows part of the parametrization used for each of the stars from the star-planet configurations explored in this article. This table summarizes the value of the coupling timescale τ_c-e (Myr) and initial time t_init (Myr), as well as the value of the wind braking constant K₁ for each stellar mass and non-decimal rotation period from one to 11 days. These values are estimated using Eqs. (18)–(20) and allow the Gallet & Bouvier (2015) models to reproduce the observed rotational distribution of low-mass stars. In our simulations, the rotation period step is 0.2 days.

Figure 7 shows the rotation period distribution of stars between 0.5 and 1.0 M_⊙ at 120 Myr. We consider here that around each star orbits a 2 M_jup mass planet initially located at 50% R_co. The colour gradient corresponds to the value of the initial rotational period of the star (between two and 11 days). The presence of a planet more massive than 1 M_jup creates a kink in the rotational period distribution of stars with mass smaller than 0.9 M_⊙. The location of the rotational kink depends on the stellar mass (it is around log P_rot,⋆ = 0.6, corresponding to P_rot,⋆ = 4 days, for 0.8 M_⊙ stars, and around log P_rot,⋆ = 0, corresponding to P_rot,⋆ = 1 days, for the0.6 M_⊙ stars). However, at the age of the Pleiades cluster, the tidal star-planet interactions only affect stars with initial rotational period of around seven days (the greenish plain circle). Moreover, this rotational kink appears by the end of the PMS phase and persists (depending on the planetary mass, see Fig. 6) up to the early MS phase.

We recall that this article does not aim to reproduce the observed rotational distribution as this has been previously investigated (e.g. Gallet & Bouvier 2013, 2015; Sadeghi Ardestani et al. 2017) but rather to explore a large star-planet parameter space, in particular the initial rotational period. In this work, we do not consider the very external slow rotation part of the distributions but more likely the first and third percentiles (see Gallet & Bouvier 2013, 2015) of these distributions. As a consequence, we do not reproduce the outliers (at least the very slow rotating stars) of these distributions at any ages. Given our uniform initial rotational distribution that ranges between one and 11 days, it was expected that our model should miss parts of the observations’ properties.

There is also a deviation in the fit of τ_disk that produces a smaller disk lifetime than inferred from Gallet & Bouvier (2015) for the low-mass stars of our sample. As a consequence, the whole rotational distribution of 0.6 and 0.5 M_⊙ stars is shifted towards a small rotational period.

Fig. 7 Upper panel: rotation period distribution of stars between 0.5 and 1.0 M⊙ at 120 Myr. The mass of the planet is taken to be 2 Mjup. Each coloured plain circle corresponds to one stellar configuration. The initial SMA of the planet is 50% Rco. Lower panel:corresponding values of δProt for each of the star-planet systems as a function of the stellar mass.

Fig. 8 Rotation period distribution of stars between 0.1 and 1.2 M⊙ from the 120 Myr Pleiades cluster. The red circles correspond to the stars that display the rotational kink. Data are from Rebull et al. (2016) and Stauffer et al. (2016). The stellar masses are estimated using the absolute Ks magnitudes and the 125 Myr isochrones from Baraffe et al. (2015). The blue and red lines correspond to the least-squares fit to the slowest rotating stars inferred by Stauffer et al. (2016, see their Fig. 2). The rotational distributions are the same as in Fig. 7.

3.4 Comparison with the Pleiades observations

In the previous sections, we pointed out that the only way for the surface rotation rate of the star to be significantly affected by the star-planet interaction is 1) when the planet is initially located inside the corotation radius (around 50% of R_co ), and 2) when the planet is more massive than 1 M_jup. We thus decided to explore the parameter space of these two main parameters to try to reproduce the Rebull et al. (2016) and Stauffer et al. (2016) observations and the so-called rotational anomalies described below. We would like to point out that in this work we do not aim to reproduce all of the features observed in the Pleiades cluster (e.g. the wide rotational distribution of the less massive stars) but only the “kink” around 0.8 M_⊙ stars.

Figure 8 shows the rotation period distribution of stars between 0.1 and 1.05 M_⊙ from the 120 Myr Pleiades cluster (also known as M45 and Melotte 22). In this figure we can see the “anomaly” detected by Rebull et al. (2016) and Stauffer et al. (2016) in the Pleiades cluster obtained using Kepler-K2 (Howell et al. 2014). It is depicted by a rotational kink located around slow rotating (five to eight days) 0.8 M_⊙ stars (early K-type stars). In Fig. 8 it is highlighted by the red plain circles in which part of the slowly rotating stars appearto rotate with a faster rotation than the “classical” slow sequence. In this figure, there are two populations: low- (M_⋆∕M_⊙≲ 0.5) and intermediate- (0.5 ≲ M_⋆∕M_⊙≲ 1.1) mass stars. The low-mass stars display a quite large rotational distribution mainly located at small rotational period (P_rot,⋆ ≈ 0.3 days, corresponding to log   P_rot,⋆ = −0.52). Since these stars are largely (or totally if M_⋆∕M_⊙ < 0.3) convective, they are less efficiently braked by the stellar wind and thus are still on the fast rotating regime at the age of the Pleiades cluster. In contrast, the intermediate-mass stars harbour a more packed distribution centred around a rotational period that depends on the stellar mass (for instance the 0.8 M_⊙ stars peak around eight days, log P_rot,⋆ = 0.9). Because these stars are partly convective, they are more efficiently braked by the stellar wind and thus possess a more compact distribution even at as early an age as 120 Myr. Among the slowest rotating stars a scatter can be seen. This anomalous scatter could be created by a physical process that should accelerate part of the slow rotators at the age of the Pleiades cluster.

According to Stauffer et al. (2016), this anomaly is statistically robust and suggests two distinct regimes, one producing the slowest rotating stars that follows the classical rotational evolution, which only invokes magnetic braking and internal redistribution of angular momentum and the other that requires another external physical effect to produce faster rotators at the same mass and age. It is observed around the early K-type 0.8 M_⊙ (V –K_s ~ 2.6) stars for which the dissipation of the dynamical tides is maximum at ZAMS (Mathis 2015b; Gallet et al. 2017a). Moreover, it seems to be time dependent as it appears to evolve toward lower masses with increasing age (see Stauffer et al. 2016). This kink should then be time and mass dependent and should be present during the early MS phase (around 120–220 Myr) and for early K-type stars (around 0.7–0.8 M_⊙ ). Can tidal star-planet interactions reproduce these observations? If it is the case, then what could be the insights about the initial distribution of massive close-in planets around low-mass PMS stars, such as the initial location and the mass of the planets so as to reproduce the observations in the Pleiades?

Our model predicts that this rotational kink should also be present for slightly less massive stars, namely 0.7–0.6 M_⊙, corresponding to V − K_s ≈ 3.3. However, this behaviour is not observed in the data from Stauffer et al. (2016). It could not be explained by a shift in the mass estimation between the STAREVOL code used in this work and the Baraffe et al. (2015) code used in Stauffer et al. (2016) since both models are quite similar and isochrones from these two models produce the same theoretical stellar masses (Louis Amard’s PhD thesis and priv. comm.). However, given the quite large scatter in the rotational distribution of 0.5-0.6 M_⊙ stars, it is arduous to distinguish the presence or not of a given rotational kink.

In this work, we arbitrarily assumed an occurrence of a massive close-in planet of 100% so as to illustrate the impact of the effect of a tidal star-planet interaction on the rotational evolution of low-mass stars. The net result of this interaction on the stellar rotational distributions naturally depends on the initial occurrence of massive planets around young low-mass stars. To investigate this, more realistic (i.e. non-uniform) initial stellar rotational and planetary distributions are first of all required. Such a complete study is, however, beyond the scope of the present paper.

We demonstrate that the presence of a massive planet can have a notable impact on the rotation period distribution of ZAMS’s stellar clusters, especially for the K-type stars. While the modelled rotational kink is produced in the correct stellar mass range, the observed rotational kink in the Pleiades is located around slower rotating stars (five to eight days) than predicted by our model. It means that a physical process impacting preferentially very slow initial rotating stars is required to only impact the rotational period distribution of this stellar population. In this work, we only consider the tidal dissipation inside the convective envelope of the host star. Including the additional dissipation by gravity waves inside the radiative core could move the rotational kink towards a longer rotational period. Indeed, the planets that are engulfed during the early-MS phase would be translated outwards compared to their current location (with the present model) and would then modify the rotation rate of slower rotating MS stars (because of the magnetic braking). This additional dissipation should then be included in future works, so as to produce more realistic star-planet evolutions.

4 Conclusions

We investigated the theoretical effect of tidal interaction on the rotation rate of stars by coupling the parametric model described in Gallet & Bouvier (2015), which includes the decoupling between the core and the envelope, to the planetary orbital evolution code of Bolmont & Mathis (2016). We find that the stellar rotation is primarily impacted during planetary engulfment events.

With this work we showed that if the planets fall into the stellar surface during the early-MS phase, then planetary accretion has a strongeffect on the surface rotation rate of the star. This is especially true in the case of a massive close-in planet (one to two Jupiter-mass planet) orbiting a low-mass star. The surface rotation rate of the star can then be accelerated up to a factor of two to three. Within the right configuration, a planetary population between 1 and 5 M_jup initially located around 50% R_co can open a rotational kink during the early MS phase, at the age of the Pleiades cluster. This rotational kink is only present for stars with initial surface rotation rate around seven days and located around log   P_rot,⋆ = 0.6 (corresponding to P_rot,⋆ = 4 days) and stellar mass smaller than0.8 M_⊙. In Stauffer et al. (2016) this rotational kink is however observed around log   P_rot,⋆= 0.8 (corresponding to P_rot,⋆= 6.3 days).

The proposed model here thus does not allow us to exactly reproduce the location of the observed anomaly, possibly because we neglected the tidal dissipation in the stellar radiative core. However, it provides a promising direction to further investigate the influence of close-in planets on the rotational period distribution of young stars, beyond the well-known period-mass relationships, especially in the framework of the Gaia DR2 (Lanzafame et al. 2018) and future PLAnetary Transits and Oscillations of stars (PLATO)/Transiting Exoplanet Survey Satellite (TESS) missions (Rauer et al. 2014; Ricker et al. 2015).

Acknowledgements

We thank the anonymous referee and the editor for their fruitful comments and suggestions. We thank our colleague Philippe Delorme for the very useful discussions of some aspects of this paper. F.G acknowledges financial support from the CNES fellowship. The authors acknowledge financial support from the Swiss National Science Foundation (FNS) and from the French Programme National de Physique Stellaire PNPS of CNRS/INSU. F.G. and C.C. acknowledge financial support from the Swiss National Science Foundation and 294 SEFRI for project C.140049 under COST Action TD 1308 Origins (PI C.C.). This work results from the collaboration of the COST Action TD 1308. S.M. and E.B. acknowledge funding by the European Research Council through ERC grant SPIRE 647383 and support from the CNES-PLATO grant at CEA Saclay. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No 742095; SPIDI: Star-Planets-Inner Disk-Interactions).

References

All Tables

All Figures

Fig. 1 Sketch of the stellar structure studied here. We consider a two-layered star: a radiative core of rotation rate Ωrad and a convective envelope of rotation rate Ωconv. Both are assumed to rotate as a solid body with different rotation rates. In red we indicate the different processes included in the model. The convective envelope is subjected to external torques (stellar wind and tidal torques), as well as torques due to the radiative core development and coupling with the convective envelope.

In the text

Fig. 2 Evolution of the orbital distance of a Jupiter-mass planet with different initial SMA (top panel) ranging from 0.8 Rco (black), 1.0 Rco (dark grey), and 1.2 Rco (light grey), δProt (middle panel), and the stellar rotation period (bottom panel) during the evolution of a 1.0 M⊙ star with an initial rotation period of three days and a coupling timescale τc-e of 0.5 Myr. Upper panel: the orbital distance of the planet is represented in full coloured lines. The upper dotted lines correspond to Porb = Prot,⋆ (the corotation limit) and the lower dotted lines correspond to Porb = 1/2 Prot,⋆ limit, which marks the region above which the dynamical tide operates. The red long-dashed line is the evolution of the Roche lobe radius. Middle panel: departure of the rotation rate of the host star (Prot,⋆ ) from an isolated star (Prot,isol.). The horizontal dashed line indicates no difference with an isolated star, and the vertical dashed line represents the localisation of the ZAMS (around 50 Myr for a 1.0 M⊙ star). Lower panel: rotation period (in days) of the host star. The horizontal dashed line correspond to the transition between saturated and unsaturated magnetic regimes (1.6 days ≈ 16 Ω⊙), and the circle represents the rotation rate of the present Sun. The time on the x-axis is given from an initial time tinit of 5 Myr.

In the text

Fig. 3 Effect of planetary companion mass on the rotational evolution of a 1.0 M⊙ mass star. The initial rotation period of the star is three days, the star is assumed to rotate as a solid body with a core-envelope coupling timescale of 0.5 Myr; SMAinit = 0.8 Rco ≈ 3.25 × 10−2 au. The masses of the planets are 10 M⊕, 1 MJup, 5 MJup, and 13 MJup, displayed as darkish to greyish coloured lines, respectively. In the upper panel, the Roche lobe limits are displayed as long-dashed lines and the radius of the star as red dash-dotted line. The horizontal dashed line in the lower panel corresponds to the transition between the saturated and unsaturated magnetic regimes, and the circle represents the surface rotation rate of the present Sun. The time on the x-axis is given from an initial time tinit of 5 Myr.

In the text

Fig. 4 Evolution of the orbital distance of a Jupiter-mass planet with an initial SMA of 0.03 au (top panel), δProt (middle panel), and the stellar rotation period (bottom panel) during the evolution of a 1.0 M⊙ star with four different initial rotation rates: one, three, six, and eight days from darkish to greyish coloured lines, respectively. In the upper panel, the dotted lines correspond to the Porb = 1∕2 Prot,⋆ limit, the red dashed line is the Roche lobe limit, and the red dashed-dotted line is the stellar radius.

In the text

Fig. 5 Distribution at 120 Myr in log Prot,⋆ obtained for a 0.8 M⊙ from an initial distribution of rotation period ranging from 2 to 11 days. Effect of the initial SMA (10, 30, 45, 50, 55, 70, and 100% Rco) in the case ofa 2 Mjup (panel a). Impact of the mass of a planet (between 1 and 5 Mjup) initially located at 50% Rco (panel b). These distributions are extracted at 120 Myr that is about the age of the Pleiades cluster.

In the text

Fig. 6 Temporal evolution of the rotation distribution related to the 0.8 M⊙ mass star. The mass of the planet is taken to be 1 Mjup (panel a) and 2 Mjup (panel b) and is initially located at 0.5 Rco. Each line correspond to one stellar configuration. The bottom grey line corresponds to an initial period of two days, the upper black line corresponds to an initial period of 11 days. The consecutive difference between each line is 0.2 day.

In the text

	Fig. 7 Upper panel: rotation period distribution of stars between 0.5 and 1.0 M⊙ at 120 Myr. The mass of the planet is taken to be 2 Mjup. Each coloured plain circle corresponds to one stellar configuration. The initial SMA of the planet is 50% Rco. Lower panel:corresponding values of δProt for each of the star-planet systems as a function of the stellar mass.
In the text

Fig. 8 Rotation period distribution of stars between 0.1 and 1.2 M⊙ from the 120 Myr Pleiades cluster. The red circles correspond to the stars that display the rotational kink. Data are from Rebull et al. (2016) and Stauffer et al. (2016). The stellar masses are estimated using the absolute Ks magnitudes and the 125 Myr isochrones from Baraffe et al. (2015). The blue and red lines correspond to the least-squares fit to the slowest rotating stars inferred by Stauffer et al. (2016, see their Fig. 2). The rotational distributions are the same as in Fig. 7.

In the text