Core-envelope Decoupling Drives Radial Shear Dynamos in Cool Stars

Magnetic fields are ubiquitous in cool stars. Young stars have significant X-ray fluxes (Stauffer et al. 1994) and emission features (Wilson 1963) in the cores of strong lines that provide evidence for magnetically driven nonthermal chromospheric and coronal heating of their outer atmospheres. Stellar activity diagnostics are variable on moderate to long timescales for stars like our Sun (Radick et al. 1998), which supports a dynamo origin for the magnetic fields that generate them (Parker 1993). The standard interface dynamo model for the Sun relies on a radiatively stable core where differential rotation can amplify fields at the base of the turbulent convection zone (Parker 1993; Charbonneau & MacGregor 1997). Alternative dynamo modes have been postulated but not observed (Charbonneau & MacGregor 1997) as stars less massive than 0.35 solar masses (lacking an interface) appear to follow the same activity-scaling relations as more massive stars (Wright et al. 2011). Dynamos driven by surface convection launch magnetized winds, which are efficient engines for extracting angular momentum (Weber & Davis 1967), thus linking observed stellar rotation rates to dynamo action. Perturbations in our standard dynamo framework therefore also impact gyrochronology, the use of stellar rotation as an age diagnostic (Barnes 2007).

Standard magnetic braking models have been unsuccessful in explaining recent observations of low-mass stars. Surface rotation evolution is observed to stall but then resume, in a strongly mass-dependent fashion (Curtis et al. 2020). The epoch of stalling occurs at tens of megayears for solar analogs and hundreds of megayears for lower-mass stars (Spada & Lanzafame 2020). Apparent rotational stalls with ongoing angular momentum loss can be explained with a lag time between the spin down of the outer turbulent layers and the inner radiative core, which is referred to as core-envelope decoupling (MacGregor & Brenner 1991). Decoupling implies the temporary presence of strong radial shears below the surface convection zone, which could induce a distinct departure from standard braking laws.

In this paper we demonstrate that stars experiencing stalled surface spin down also have unusually strong surface magnetic fields. In Praesepe, these stars have spot levels comparable to saturation, which is a maximal state of surface magnetism. This provides evidence for an alternative dynamo mechanism driven by radial shears, supports core-envelope decoupling as an astrophysical phenomenon, and demonstrates that stars in the epoch of decoupling experience shear-enhanced magnetism for hundreds of megayears to gigayears on the main sequence.

Prior attempts to explain rotational stalling with core-envelope decoupling have not been definitive. Core-envelope decoupling gyrochrones appear to reproduce the morphology of open cluster rotation trends with mass and age (Spada & Lanzafame 2020), but a temporary weakening of stellar winds or surface magnetism as the source of the stalling cannot be ruled out (Curtis et al. 2020). Rotationally induced mixing can also explain main-sequence lithium depletion. Evolutionary models of lithium depletion favor larger internal shears for longer timescales in lower-mass stars, consistent with core-envelope decoupling predictions (Somers & Pinsonneault 2016); however, such models assume an invariant birth rotation distribution between clusters and deviate from observations at still older ages. Nevertheless, recent ground-based observations of field stars observe that a gap in the period distribution of more massive stars disappears in fully convective stars (Lu et al. 2022), which is expected if core-envelope decoupling is the source of the feature.

If core-envelope decoupling actually occurs, we expect strong internal shears during that epoch to have a direct observational signature. The large rotational shear between stellar surface and core predicted by decoupling is thought to drive dynamo action (Bouvier et al. 1997). This contrasts with the standard interface dynamo model, which predicts that stars with slow surface rotation would be inactive. However, we argue that a conclusive enhanced activity signal has not previously been seen in the literature because of poor sampling and precision. Starspot measurements have been documented as overactive in a handful of young stars in the Pleiades but are limited to an unconverged rotational sequence (Cao & Pinsonneault 2022). UV measurements suggested overactivity in a small sample of K dwarfs in the Hyades but with considerable scatter (Richey-Yowell et al. 2022).

A technique for measuring starspot filling factors from APOGEE high-resolution near-infrared spectra was recently published (Cao & Pinsonneault 2022), introducing a population-scale technique for studying stellar magnetism. The ensemble of spectral features in the H band permits precision measurements of starspots with typical dispersions ranging from ∼6.6% in active stars to ∼1.1% in solar-type stars, which enabled the construction of a mean relation between the Rossby number (the ratio of rotation period to convective overturn timescale) and starspots (Cao & Pinsonneault 2022). The recent 17th data release of the APOGEE survey (Beaton et al. 2021) includes a deep sampling of the ∼670 Myr (Douglas et al. 2019) anomalously rotating open cluster Praesepe, enabling a detailed study of stellar magnetism during the epoch of stalled spin down.

2.1. Sample Selection and Stellar Parameters

2.2. Gyrochronology Models

We measure starspot filling fractions and temperature contrasts jointly with spectroscopic stellar parameters (T_eff, [M/H], logg, , v_dop) in 179 Praesepe members, obtaining a final sample of 88 single stars with rotation periods (Appendix A.1). Figure 1 shows the dependence of starspot filling factor with stellar effective temperature. Stars highlighted in red are over-spotted for their Rossby number by >3σ after accounting for total cycle and inclination variation, found empirically to be 6.0% in Praesepe. Assuming starspot filling fractions are normally distributed about the mean Rossby relation, the probability of finding ≥15 stars to be 3σ discrepant out of 88 single stars is 0.00003%. Combined long-term and short-term starspot cycle variations estimated using scaling relations inferred from the Mt. Wilson solar analog (Radick et al. 1998) and starspot cluster samples (Cao & Pinsonneault 2022) yield f_spot,rms ∼1.2% for solar-type stars and f_spot,rms ∼0.9% for active stars. As overactive stars here are discrepant by at least 18%–20% in f_spot, the total cycle variation is an order of magnitude too low to explain these departures. Overactive stars are confined to a narrow mass range, making starspot cycles an unlikely explanation. These stars are confirmed members and pass checks against binary contamination.

Zoom In Zoom Out Reset image size

In Figure 2 we show the 670 Myr Praesepe sequence as points scaled by f_spot and colored by the ratio between observed and predicted f_spot using an empirical starspot relation (Cao & Pinsonneault 2022). We overlay these points on the 1 Gyr old cluster NGC 6811 in blue, showing that stars which have not spun down during the ∼300 Myr of evolution between Praesepe and NGC 6811are significantly spot enhanced.

Zoom In Zoom Out Reset image size

Stellar activity is related to the predicted spin-down rate in wind models, so the observed starspot enhancements disfavor models where stalled spin down is associated with a reduction in starspot activity. Our data represent a striking counterexample to standard models of stellar dynamo and angular momentum evolution. No other proposed mechanism naturally produces the simultaneous magnetic and rotational anomalies in Praesepe without core-envelope decoupling. These stars are on the fast branch of the McQuillan gap (McQuillan et al. 2014), and their effective temperatures are affected by starspots (Cao & Pinsonneault 2022), which suggests that the morphology of the gap may be influenced by a transition in surface spot coverage.

Models initialized with rotation rates spanning the 7.5–92.5th percentiles in Upper Sco (see Appendix A.2) are plotted in Figure 3. In the top panel, we plot S_ce, the core-envelope radial shear , where τ_CZ is the convective overturn timescale, normalized to the solar latitudinal shear. The observed overactivity in Praesepe coincides with transitions to high internal shears in two-zone rotational models, with the strongest signal seen at intermediate shears.

Zoom In Zoom Out Reset image size

In the bottom panel of Figure 3 we demonstrate the impact of recoupling on core angular momentum content. We introduce a coupling diagnostic, , where denotes the angular momentum of the stellar core and denotes the same in a corresponding decoupled model with recoupling switched off. We observe that stars with strong shears (S_ce ≳ 200) but limited A.M. transport (F_couple ≲ 0.05) from the core are not overactive, while stars with intermediate shears (200 ≳ S_ce ≳ 50) and ongoing A.M. transport (F_couple ≳ 0.05) are over-spotted.

Identifying significant internal shears and core-envelope coupling as parameters correlated with shear-enhanced activity in Praesepe, we model the time evolution of shear-enhanced magnetism at fixed rotation and metallicity in Figure 4. Core-envelope shears enhance activity in K dwarfs in Praesepe, but at older ages the magnetic anomaly moves to lower-mass stars (Figure 4)—a consequence of a coupling timescale varying strongly with stellar mass (MacGregor & Brenner 1991; Denissenkov et al. 2010; Gallet & Bouvier 2015; Lanzafame & Spada 2015; Somers & Pinsonneault 2016). Only unusually rich young clusters deviate from standard rotational distributions (Moraux et al. 2013), and the structural impact of metallicity appears to dominate observed trends with rotation and activity (Amard et al. 2020; See et al. 2021), suggesting that this result is robust. We predict starspot enhancements in Ruprecht 147 (∼2.7 Gyr, Curtis et al. 2020) between ∼3800 and 4300 K for a [M/H] of +0.1, and in M67 (∼4 Gyr) between ∼3700 and 4300 K at solar metallicity, noting that a metallicity offset of +0.1 dex roughly decreases the T_eff at fixed mass by ∼100 K. This work is evidence for an evolving shear front affecting cool stars, differentially enhancing stellar surface magnetism and impacting their high-energy radiation environments for hundreds of megayears to gigayears during the recoupling phase.

Zoom In Zoom Out Reset image size

The large starspot filling fractions measured in Praesepe overactive stars are not merely a general feature of K-dwarf spectra; in the Pleiades, there is no activity enhancement in K dwarfs (Cao & Pinsonneault 2022). In the ∼4 Gyr old open cluster M67, early K dwarfs appear unspotted, with filling fractions consistent with zero (Cao & Pinsonneault 2022). Though sparsely sampled by DR17, there also appears to be a spot enhancement signal for K dwarfs in the ∼727 Myr (Douglas et al. 2019) Hyades cluster. Overactive stars have similar spectroscopic chi-squared values to other members and do not preferentially show convergence failures in spectroscopic fits. The contamination rate of unresolved binaries is minimal due to our comprehensive binary rejection (Appendix A.1). APOGEE targets members to lower masses than our sample (Beaton et al. 2021), the rotation catalog (Rebull et al. 2017) is highly complete (∼67%), and starspot measurements are precise for stars both more and less active than our sample (Cao & Pinsonneault 2022). In the range we identify for overactivity in Praesepe (1.5 < B_P − R_P < 1.9, or 3950 < T_eff < 4410 K), all 21 single stars, and all nine members flagged as binaries, have spot fractions >1σ above the mean relation.

The mechanism mediating angular momentum transport between the core and envelope is actively debated. Hydrodynamic mechanisms cannot explain the internal solar rotation profile alone, so waves or magnetic fields are probable causes (Pinsonneault 1997). Our observations disfavor models involving strong magnetic coupling as they predict rigid internal rotation on short timescales (Fuller et al. 2019). Instead of direct coupling, models that exchange angular momentum nonmagnetically between a magnetically coupled core and convective envelope across thin shear layers permit long coupling timescales and persistent core-envelope shears (Oglethorpe & Garaud 2013). Simple estimates of the characteristic angular momentum transport timescale via internal gravity waves within 0.5 pressure scale heights of the convective zone boundary (following Section 4.2 of Fuller et al. 2014) produce sufficiently slow timescales that scale strongly (∝M^−7.6) with stellar mass. Alternatives to the decoupling model exist, such as the increase in complexity of magnetic field morphologies (Garraffo et al. 2018). Although interesting in the saturated domain, such mechanisms do not predict significant departures from the mean Rossby scaling and cannot explain the overactivity of stars experiencing stalled spin down. Magnetic field measurements in field populations do not appear to show a decline in the dipole field strength (See et al. 2019), but future work will more conclusively discriminate between these models.

We show in this work that core-envelope decoupling and stalled spin down are linked, producing observables not only in angular momentum evolution but also in activity. The existence of shear-enhanced magnetism suggests that stellar mixing processes can persist for hundreds of millions to billions of years—consistent with models where waves are the primary mechanism linking the cores and envelopes of stars. Evolutionary models will need to include core-envelope decoupling and angular momentum losses from dynamos with large core-envelope shears to obtain accurate rotation-based ages on the main sequence. Habitability models will need to be informed by the enhanced magnetism associated with core-envelope decoupling as shear-enhanced activity can last for billions of years in low-mass stars. Decoupling signals in activity and rotation are strongly mass- and age-dependent and are potential age diagnostics for low-mass main-sequence stars, where there are other age indicators are scarce.

We thank the SDSS-IV collaboration for making this research possible. Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy Office of Science, and the Participating Institutions. SDSS-IV acknowledges support and resources from the Center for High Performance Computing at the University of Utah. The SDSS website is www.sdss4.org.

L.C and M.H.P acknowledge computing resources and support made available by the Ohio State University College of Arts and Sciences.

L.C. and M.H.P. acknowledge support from TESS Cycle 5 GI program G05113 and NASA grant 80NSSC19K0597. J.v.S. acknowledges support from the National Science Foundation grant AST-1908723. This research was supported in part by the National Science Foundation under grant No. NSF PHY-1748958.

A.1. Binary Rejection

Binarity confounds measurements of starspot filling fraction for dwarfs in the case where the companion contributes significant flux (Cao & Pinsonneault 2022). In star clusters, targets with significant light from a companion appear unusually bright relative to single stars, and we identify such photometric binaries as lying more than 0.25 mag above the empirical 75th percentile main sequence. We use radial velocity and astrometric variations to detect the gravitational influence of a companion; we also identify sources that may overlap in the line of sight of the APOGEE fiber (Cao & Pinsonneault 2022). The results of this binary rejection technique are shown in Figure 5. We find 88/136 (65%) of our rotating sample are single stars in Praesepe.

Zoom In Zoom Out Reset image size

All single stars, including the overactive stars, pass each of the binary rejection stages where data were available; all overactive stars have Gaia astrometric RUWE measurements below 1.1, only two of them do not have an APOGEE radial velocity variability measurement, and all are on the single star photometric sequence.

The starspot measurements have been demonstrated to be sensitive even in inactive stars (Cao & Pinsonneault 2022), meaning that the underlying spot measurements are unlikely to be a source of selection effects. Figure 5 shows a Praesepe color–magnitude diagram with a histogram indicating the distributions of binaries, stars with no detected rotation periods, and single stars in our rotation sample. There do not appear to be any significant biases in the single stars. In particular, the overactive region of stars in Praesepe is not marked by an increase in rotation nondetections, nor are there any apparently unspotted members in that temperature window. We provide spectroscopic starspot parameters and flags for this sample in Table 1.

A.2. Gyrochronal Calibration

We tune the free parameters in our model to simultaneously reproduce the observed rotation of the Sun, Pleiades, and Praesepe open clusters. From Upper Sco, we characterize the range of possible initial angular momenta as a function of mass by adopting the 7.5th, 50th, and 92.5th percentile periods in mass bins spanning 0.2–1.3 M_⊙ in steps of 0.1 M_⊙. This choice captures the range but avoids the extrema, preventing rare tidally locked binaries and background sources from dominating our fits. We draw underlying stellar structure profiles from SPOTS model grids (Somers et al. 2020) at appropriate metallicities and spot filling fractions. These models are then subjected to a rotational braking calculation (van Saders & Pinsonneault 2013) with a two-zone angular momentum transport model (Denissenkov et al. 2010). Both metallicity and starspot filling fraction affect the mapping from mass to radius and convection zone properties, which enter at significant powers in the braking law.

We adopt the age of 112 ± 5 Myr for the Pleiades, based on the minimum mass for the destruction of lithium in brown dwarfs (Dahm 2015). Its metallicity is near-solar, with a literature estimate of 0.03 ± 0.02 (random) ± 0.05 (systematic) dex (Soderblom et al. 2009); in this work, we assume a solar metallicity for the Pleiades.

The age of Praesepe is constrained by isochrones; for this work, we assume the value of 670 ± 67 Myr (Douglas et al. 2019). The metallicity of Praesepe is significantly nonsolar, with recent estimates from the literature varying from +0.11 to +0.27 (D’Orazi et al. 2020). We jointly obtain a spectroscopic metallicity estimate for the cluster of +0.188 ± 0.006 dex for 95 stars after removing two >3σ outliers. We therefore adopt a fiducial cluster metallicity of +0.188 dex for the Praesepe gyrochrone. Changes of 0.01 dex in metallicity and 10% in age, comparable in size to our ∼1σ confidence intervals, induce ∼5.0% and ∼13% offsets in predicted shear, which do not substantively alter our conclusions.

We divide the observed Pleiades rotation periods into three period bins (0–15th, 15–85th, and 85–100th percentiles) and 12 mass bins, determining the average f_spot in each bin. We initialize 12 model tracks (0.2–1.3 M_⊙ in steps of 0.1 M_⊙) using Pleiades-determined f_spot values, evolving each gyrochrone in the three Upper Sco period launch bins to Pleiades age (112 Myr adopted). For each star in each of the respective Pleiades period bins, we quantify agreement by computing with the corresponding gyrochrone. This procedure in its entirety is repeated for [M/H] = +0.188 models launched from Upper Sco to Praesepe age (670 Myr adopted), assuming the same mass dependence of the initial rotation periods. The 50th percentile track is evaluated at solar mass, composition, age, and f_spot = 0 compared to the solar rotation period (26.870 day). From empirical relations in the Pleiades (Cao & Pinsonneault 2022), we adopt a Gaussian prior on ω_crit with a central value of 3.341 × 10⁻⁵ rad s⁻¹ and 1σ range of 3.248 × 10⁻⁵–3.433 × 10⁻⁵ rad s⁻¹.

The calibration of the gyrochronology models involved the evaluation of a maximum likelihood estimator at each chain step of the form . This was evaluated as a sum including each star in both clusters at appropriate metallicities, spot filling fractions, and across each of three gyrochrones representing each of three rotation percentile bins; model expectations were produced through linear interpolation. The Sun was also included, with an assumed error of σ_P = 0.010 day. This produces gyrochrones that produce the solar rotation rate at solar age with 50th percentile rotators to within 0.006 day. We use emcee (Foreman-Mackey et al. 2013) with 32 walkers running for 2100 steps each. A burn-in of 300 steps was discarded, and the subsequent chains were found to be converged using an integrated autocorrelation time diagnostic (Foreman-Mackey et al. 2013).

The results from the Markov chain Monte Carlo fit are represented in Figure 6, with and α = −11.8 ± 1.0. We initially find a steeper α than the literature, which we attribute to the inclusion of metallicities and starspots in the stellar mass determination for our fit. After accounting for these differences, our slope is consistent with previous analyses (Somers & Pinsonneault 2016). The resultant model cluster sequences for the Pleiades and Praesepe are shown in Figure 7. These solar-calibrated gyrochrones are also in good agreement with the observed rotation sequence of M67 (4 Gyr adopted).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

A.3. Starspots and Color Anomalies

Color anomalies and flux excesses are known to correlate with chromospheric activity proxies and have been observed in the literature to be a result of starspots (Campbell 1984; Stauffer et al. 2003). With standard color relationships at the appropriate metallicities for the Pleiades and Praesepe (Casagrande & VandenBerg 2014), it is possible to check whether the overactive stars appear chromospherically discrepant in photometry. In Figure 8 we show sequences for both open clusters, showing that non-spotted stars have minimal color anomalies, while all spotted stars have strongly anomalous colors. We interpret this as a photometric confirmation of the large spot fraction on these stars. The divergence between the predicted and actual color sequences at V − K > 4 is attributed to poor model atmosphere fits for cool stars.

Zoom In Zoom Out Reset image size