Stellar mid-life crisis: subcritical magnetic dynamos of solar-like stars and the breakdown of gyrochronology

TL;DR

This paper proposes that subcritical magnetic dynamos in aging solar-like stars explain the breakdown of gyrochronology, the bimodal sunspot distribution, and the transition to magnetic inactivity, challenging traditional stellar spin-down models.

Contribution

It introduces the concept of subcritical magnetic dynamos as a unifying physical mechanism for various observed stellar phenomena and their evolution with age.

Findings

Reproduces bimodal sunspot distribution with subcritical dynamos

Suggests aging stars access weak magnetic field regimes

Proposes disrupted magnetic braking causes slower spin-down

Abstract

Recent observations reveal the surprising breakdown of stellar gyrochronology relations at about the Sun's age hinting that middle-aged, solar-like stars transition to a magnetically inactive future. We provide a theoretical basis for these intriguing observations inspired by simulations with a mathematical dynamo model that can explore long-term solar cycle fluctuations. We reproduce the observed bimodal distribution of sunspot numbers, but only for subcritical dynamos. Based on a bifurcation analysis, we argue that ageing of solar-like stars makes the magnetically-weak dynamo regime readily accessible. Weak magnetic field production in this regime compromises wind-driven angular momentum losses thus disrupting the hegemony of magnetic braking on stellar rotational spin-down. This hypothesis of {\emph{subcritical magnetic dynamos}} of solar-like stars provides a self-consistent, unifying physical basis for a diversity of solar-stellar phenomena such as why stars beyond their mid-life do not spin-down as fast as in their youth, the break-down of stellar gyrochronology relations, the observed bimodal distribution of long-term sunspot observations and recent findings suggesting that the Sun may be transitioning to a magnetically inactive future.

Equations

$$\frac{d{B}_{\varphi }(t)}{dt}=\frac{\omega }{L}A(t-T_0)-\frac{B_{\varphi }(t)}{\tau }$$

$$\frac{dA(t)}{dt}={\alpha }_{BL}\cdot f\left(B_{\varphi }\left(t-T_1\right)\right)\cdot B_{\varphi }(t-T_1)-\frac{A(t)}{\tau }$$

$$\frac{dA(t)}{dt}={\alpha }_{BL}\cdot f\left(B_{\varphi }\left(t-T_1\right)\right)\cdot B_{\varphi }(t-T_1)-\frac{A(t)}{\tau }+ϵ(t)$$

Full text

1 Introduction

The widely used diagnostic of stellar age rotation relationship that is the basis of gyrochronology has recently been found to be unreliable after stellar midlife (van Saders et al., 2016). The standard understanding of angular momentum loss mediated braking of rotation due to magnetically-coupled stellar wind cannot explain this. It has been reported that the Sun is currently going through a weakened braking phase at about its midlife of 4.5 billion years (van Saders et al., 2016). What causes this weakened magnetic braking is an intriguing question. Latest measurements confirm that the Sun is less magnetically active than similar stars (Reinhold et al., 2020b). This has spurred debates and controversy on whether the Sun is already transitioning to a magnetically inactive future (Metcalfe & van Saders, 2020; Reinhold et al., 2020a); studies exploring rotation period magnetic cycle period relationships support this hypothesis (Metcalfe & van Saders, 2017; Metcalfe et al., 2016). Since stellar magnetic activity is ultimately a product of the magnetohydrodynamic (MHD) dynamo mechanism, a causal basis could exist between the nature and evolution of stellar dynamo action and these observations. Inferring the nature of magnetic dynamos in Solar-like stars at about their midlife, therefore, could be the Rosetta stone to explain these apparent disparate observations.

With an aim to characterize solar-like dynamos at around their midlife – i.e., the current age of the Sun we look at long-term solar observations. These reveal intriguing phases of reduced sunspot activity that are termed as grand minima episodes, the most recent being the Maunder minimum between 1645 to 1715 AD (Eddy, 1976). Several long-term solar activity reconstructions based on cosmogenic isotopes like ${Be}^{10}$ in ice cores and $C^{14}$ in tree rings suggest that such grand minima have also occurred in the past (Usoskin et al., 2007); reconstructed data reveal 20 grand minima in the last 9000 years (Usoskin et al., 2016). Such grand minimum have recently been confirmed as a special mode of solar dynamo operation, distinct from the regular activity cycle mode (Usoskin et al., 2014).

Similar low activity phases have been observed in other Sun-like stars, which have been dubbed “Maunder-minima stars” (Baliunas et al., 1995); the nature of magnetic dynamos in these stars is hotly debated. A crucial observation is that the main sequence age of these low activity “Maunder-minima stars” are typically similar or older compared to the Sun. This is thought to be due to stellar spin down through steady angular momentum losses mediated via stellar winds and consequent reduction in the efficiency of the dynamo mechanism (Nandy, 2004; Brun et al., 2015). It has been speculated that this might be related to a drastic change in the nature of solar-stellar dynamos at the current solar age (Metcalfe & Egeland, 2019; Metcalfe & van Saders, 2017; Metcalfe et al., 2016).

Are the wandering of the Sun’s magnetic activity between the two distinct modes of low and regular activity cycles, the occurrence of “Maunder minima stars” and the observed break-down of the rotation rate–age relationship connected? Based on a reduced dynamo model that imbibes the spiritual essence of Babcock-Leighton and turbulent dynamos, and utilizing stochastically forced, non-linear delay differential equations, we show these diverse phenomena can be explained on the basis of bimodality of the underlying dynamo mechanism. Based on this analysis, we argue that a stellar evolution mediated shift in the dynamo efficiency enables a transition across supercritical to subcritical states that makes the full dynamical regime – including very weak activity states – accessible to an aging star.

The reconstructed sunspot number suggests that the occurrence of grand minima is not a cyclic phenomenon but an outcome of an underlying generation mechanism that is inherently stochastic and/or chaotic (Usoskin et al., 2007). In a chaotic dynamo, the long time-scale ‘supermodulation’ could be related to symmetry changes in magnetic fields (Weiss & Tobias, 2016). The current understanding of large-scale magnetic field generation in Sun-like stars is based on the recycling of toroidal and poloidal magnetic field components (Parker, 1955; Charbonneau, 2020). The toroidal field is generated from the poloidal field due to stretching by differential rotation. For the current Sun, observations and data-driven model reconstructions of solar activity (Dasi-Espuig et al., 2010; Cameron & Schüssler, 2015; Bhowmik & Nandy, 2018) demonstrate that the dominant mechanism for poloidal field creation is the so-called Babcock-Leighton (BL) mechanism, i.e., the decay and dispersal of bipolar sunspot pairs which contribute to a large-scale dipolar field (Babcock, 1961; Leighton, 1969). This mechanism relies on the existence of strong toroidal flux tubes in the Sun’s convection zone, which is not the case during grand minima episodes. It has been demonstrated recently that a mean-field dynamo $\alpha$-effect capable of working on weaker fields can “rescue” the solar cycle from a grand minimum (Passos et al., 2014; Hazra et al., 2014) and therefore, must be operational in stellar convection zones. Direct numerical simulations utilizing the full magnetohydrodynamic (MHD) system of equations or spatially extended mean-field dynamo models are too expensive computationally for the necessary long-term simulations to explore the range of possible dynamical solutions a solar-like dynamo mechanism may exhibit over long-term. For our numerical bifurcation analysis, we adapt a reduced BL dynamo model based on delay differential equations that has been found very useful in providing insights on long-term solar cycle behavior (Wilmot-Smith et al., 2006; Hazra et al., 2014).

In Babcock-Leighton dynamos, the sources of the toroidal and the poloidal fields are spatially segregated as they act across a distributed domain in the solar convection zone (SCZ). Such spatial segregation necessitates a finite time for flux transport across the SCZ, introducing time delays in the communication between the two source layers (Jouve et al., 2010b). These time delays have been shown to create a dynamical memory even in stochastically forced dynamo systems which allows a predictive window for the activity (Yeates et al., 2008; Karak & Nandy, 2012; Bhowmik & Nandy, 2018).

The evolution of the toroidal field ($B_{\phi}$) and the vector potential ($A$) for the poloidal field in the reduced BL dynamo model are governed by the following delay differential equations:

[TABLE]

where $T_{0}$ and $T_{1}$ are the time delays in the conversion of the poloidal flux to the toroidal one and vice-versa. The interplay amongst meridional circulation, turbulent diffusion and turbulent pumping mediated flux transport (Hazra & Nandy, 2016) determine $T_{0}$. While there have been suggestions that meridional flow variations can explain the observed stellar activity amplitude-period relationship (Nandy, 2004), numerical simulations exploring meridional flow dependence on stellar rotation rates indicate the necessity of more complex flow profiles or additional mechanisms to explain this (Jouve et al., 2010a). Thus we reiterate that $T_{0}$ must be interpreted as a generic flux transport timescale, not necessarily determined by stellar meridional flows alone. The turbulent dissipation timescale of the field is $\tau$, $\omega$ is the difference in rotation rates across the SCZ of length $L$, and $\alpha_{BL}$ represents the BL source term for the poloidal field. The quenching function $f$ constrains the poloidal source to work within a finite range of toroidal field strength with a non-zero lower bound. This is motivated from the perspective that weak toroidal flux tubes are shredded by turbulence and cannot form sunspots and very strong flux tubes are immune to the Coriolis force, which causes them to emerge at the solar surface without any tilt. This effectively quenches the BL source (for further details and motivation, see (Wilmot-Smith et al., 2006; Hazra et al., 2014)). The efficiency of the dynamo is measured by a non-dimensional number called the dynamo number, which is defined as $N_{D}=\omega\alpha_{BL}\tau^{2}/L$.

2 Results and Discussion

Long-term simulations with this dynamo model shows hysteresis (see Fig. 1 left panel). There is a hint of similar behavior in earlier studies with non-linear mean field dynamos and $3$-dimensional MHD simulations of turbulent dynamos (Kitchatinov & Olemskoy, 2010; Karak et al., 2015b; Oliveira et al., 2021). However, these computationally expensive numerical simulations precluded an adequate analysis of the distribution function of a large number of magnetic cycles. The results from our model (Fig. 1 left panel) substantiate the proposition and demonstrate that hysteresis arises from two distinct modes of dynamo activity. Note that the results are qualitatively robust with respect to choice of parameters within the regimes established in Hazra et al. (2014) and Wilmot-Smith et al. (2006).

The dynamo cannot self-consistently recover from a grand minimum once it is nudged into such a phase by stochastic fluctuations in the BL source term (Fig. 1 middle panel). Our extensive parameter space study shows that the BL dynamo shut-down threshold (dashed horizontal line in the inset of the middle panel of Fig. 1) shifts upward to engulf more cycles into grand minima episodes when the timescale $T_{0}$ is increased in the diffusion-dominated regime (or decreased in the advection-dominated regime) as defined in Wilmot-Smith et al. (2006).

We find that the addition of a “magnetic noise” allows the dynamo to recover from grand minima episodes as well as sustains hysteresis. Incorporation of additive noise, $\epsilon(t)$, modifies equation (2) to the form:

[TABLE]

where $\epsilon(t)$ is uniform and white in time (Charbonneau, 2001; Charbonneau et al., 2004, 2007) and has zero mean. Since the noise has zero mean, it sometimes builds up the poloidal field while at other times, reduces the existing dipolar field of the Sun. Magnetic noise definitely abounds in the solar convection zone, especially near the sub-surface layers of the Sun where the convection is much more turbulent. Small-scale dynamo action in this region (Meneguzzi & Pouquet, 1986; Nordlund et al., 1992) produces magnetic field that is highly uncorrelated in space and time. Such noise continuously regenerates and replenishes the magnetic field in the upper layers of the Sun (Cattaneo et al., 2003; Hagenaar et al., 2003; Charbonneau et al., 2004) and our results (Fig. 1 right panel) indicate this could be enough to recover solar-stellar dynamos from a grand minima over a growth time that is dependent on the level of imposed noise.

A time series of the smoothed cycle-averaged sunspot number from an additive-noise incorporated Babcock-Leighton dynamo model is shown in Fig. 2, upper panel. The smoothing has been done using the same procedure used by Usoskin et al. (2016), i.e., we first average the monthly sunspot number (SSN) over the usual cycle period and apply the Gleissberg low pass 1-2-2-2-1 filter. The shaded regions below the horizontal line signify grand minima, which are defined as epochs during which the smoothed SSN falls below 50% of its average value for at least 3 consecutive decades. We observe 13 to 18 grand minima (18 in Fig. 2) in different realizations of stochastic fluctuations in a typical 9000 yr simulation, which is in close agreement with the reconstructed data identifying 20 grand minima in the last 9 millennia; this is of course sensitive to model parameter choices (Choudhuri & Karak, 2012).

What is fascinating is that the assimilation of magnetic noise in our model generates magnetic cycles that qualitatively reproduce the observed bimodal distribution of toroidal field energy; see Fig. 3 and for a comparison with observations, Fig. 3 in Usoskin et al. (2014). This implies that grand minima episodes cannot be considered as events that correspond to a tail of a single regular mode of solar activity, but it is in itself a separate state of dynamo operation. Most previous studies of grand minima relied on fluctuations in a single regular mode of activity cycles (Olemskoy et al., 2013; Moss et al., 2008; Brandenburg & Spiegel, 2008) and thus could not explain the observed bimodal distribution. This bimodal nature in the probability density of cycle-averaged SSN is, in fact, due to the hysteresis present in the model (Kitchatinov & Nepomnyashchikh, 2015). Hysteresis in a certain parameter regime emerges with the existence of bistable solutions, which emanate in our model due to the non-zero lower operating threshold of the Babcock-Leighton process and the shorter poloidal-to-toroidal flux conversion timescale. Numerical experiments with our dynamo model show that if we remove the lower operating threshold and increase the time delays simultaneously, hysteresis and the bimodal nature of the probability density of cycle-averaged sunspot number is no longer manifested.

Our model also qualitatively reproduces the toroidal field’s phase space collapse; compare Fig. 4 upper panel here with fig. 10 in Lopes et al. (2014), which Passos & Lopes (2011) unsuccessfully tried to achieve using their low dimensional model with stochastic fluctuations in $\alpha$-effect alone. The collapse in the phase space implies that a grand minimum episode results in strong intermittency in magnetic cycles and that weak intermittency (or moderate amplitude cycles) cannot be attributed to a grand minimum epoch. It is interesting to note in Fig. 4 (upper panel) that the polarity reversal of magnetic cycles, albeit far less regular, still occurs during the grand minimum when the Babcock-Leighton source is dysfunctional. This is because the stochastically injected poloidal field due to magnetic noise still provides a seed for differential rotation to generate toroidal fields. A combination of statistical fluctuations in the right direction can then nudge the weak toroidal field such that it exceeds the lower operating threshold and revive the large-scale BL dynamo that produces sunspots. The lower panel of Fig. 4 portrays the corresponding temporal variation of toroidal field energy depicting this dynamics.

We observe the bimodal distribution of sunspot numbers only in the (significantly) subcritical regime of the dynamo. Earlier works in this context were inconclusive and conflicting (Karak et al., 2015b; Kitchatinov & Nepomnyashchikh, 2017). The dynamo hysteresis in our model sustains coexisting oscillating and decaying solutions for similar dynamo parameters, which is not possible for supercritical dynamos. Thus, our theoretical analysis establishes that the observed bimodal distribution of decadal sunspot numbers is a direct manifestation of the subcritical nature of the current solar dynamo (at midlife of its evolution as a star).

Over its life time, a Sun-like star spins down due to angular momentum loss through magnetically and thermally driven stellar winds. Aging and lower rotation rates (i.e., higher rotation period $P_{\Omega}$) increases the stellar Rossby number (${R_{o}}$ = $P_{\Omega}/\tau_{c}$, where $\tau_{c}$ is the convective turn-over time) – a parameter used widely in observational studies of stellar activity. The efficiency of the dynamo generation of magnetic fields represented by the dynamo number scales as ${N_{\mathrm{D}}}\sim 1/{R_{o}}^{2}$. Based on our earlier analysis we postulate that the aging of a star – with a concomitant increase in ${R_{o}}$ and decrease in ${N_{\mathrm{D}}}$ – makes the subcritical dynamo regime with bimodal activity distribution accessible around the midlife of its evolution on the Main Sequence.

3 Conclusion

Our analysis is relevant only for stars which continue to exhibit solar-like rotation and thus not applicable for anti-solar type differential rotation, see, e.g., Karak et al. (2015a); Karak et al. (2020); Brandenburg & Giampapa (2018). A systemic change in global coronal magnetic field configuration (e.g., from closed to open loops) may also impact wind properties and the nature of magnetic braking which we do not address here.

Our findings have several important consequences. We establish that Maunder-like grand minima are a distinct mode of dynamo operation possible only in sub-critical dynamo regimes. Stellar observations indicate the existence of the low activity “Maunder minima stars”, which are similar or older than the Sun (Baliunas et al., 1995) and have high ${R_{o}}$. Our framework suggests these stars are in the sub-critical regime, where accessing the low activity mode of the bimodal stellar dynamo becomes possible. Based on our analysis, we postulate that this is not due to a drastic change in the underlying dynamo mechanism itself, but merely the opening up of the possibility of existing in the low activity mode – wherein, very low dynamo output drastically reduces magnetically driven stellar winds and angular momentum loss rates. At about stellar mid-life, the magnetically-weak sub-critical states are readily attainable. These weak-magnetic epochs disrupt the hegemony of magnetic braking on stellar rotation. This argument is consistent with the suggestion by Metcalfe et al. (2020) – namely, magnetic cycles become weaker while still intermittently exhibiting large-scale dynamo action. This regime allows for the Rossby number (or equivalently, the dynamo number) to meander across the critical Rossby number for transition between different activity states or remain nearly constant, as hypothesized by Metcalfe et al. (2020) for 94 Aqr Aa.

To summarize, the subcritical solar magnetic dynamo hypothesis self-consistently explains the breakdown of stellar gyrochronology, the bimodal nature of the dynamo mechanism in the current Sun, solar cycle transition between grand minima and regular activity modes and the Sun’s possible transition to a less magnetically active future.

Acknowledgments

The authors thank Travis S. Metcalfe for insightful discussions. B.T. wishes to acknowledge the support extended to him during a Visiting Fellowship at the Center of Excellence in Space Sciences India (CESSI) at IISER Kolkata. CESSI is funded by the Ministry of Education, Government of India under the Frontier Areas of Science and Technology Scheme.

Data Availability

The data underlying this article will be shared on reasonable request to the corresponding author.