Magnetic braking saturates: evidence from the orbital period distribution of low-mass detached eclipsing binaries from ZTF

ABSTRACT

We constrain the orbital period (P_orb) distribution of low-mass detached main-sequence eclipsing binaries (EBs) with light-curves from the Zwicky Transient Facility (ZTF), which provides a well-understood selection function and sensitivity to faint stars. At short periods (P_orb ≲ 2 d), binaries are predicted to evolve significantly due to magnetic braking (MB), which shrinks orbits and ultimately brings detached binaries into contact. The period distribution is thus a sensitive probe of MB. We find that the intrinsic period distribution of low-mass (0.1 ≲ M₁/M_⊙ < 0.9) binaries is basically flat (⁠|${\rm d}N/{\rm d}P_{\rm orb} \propto P_{\rm orb}^0$|⁠) from P_orb = 10 d down to the contact limit. This is strongly inconsistent with predictions of classical MB models based on the Skumanich relation, which are widely used in binary evolution calculations and predict |${\rm d}N/{\rm d}P_{\rm orb} \propto P_{\rm orb}^{7/3}$| at short periods. The observed distributions are best reproduced by models in which the magnetic field saturates at short periods with a MB torque that scales roughly as |$\dot{J}\propto P_{\rm orb}^{-1}$|⁠, as opposed to |$\dot{J} \propto P_{\rm orb}^{-3}$| in the standard Skumanich law. We also find no significant difference between the period distributions of binaries containing fully and partially convective stars. Our results confirm that a saturated MB law, which was previously found to describe the spin-down of rapidly rotating isolated M dwarfs, also operates in tidally locked binaries. We advocate using saturated MB models in binary evolution calculations. Our work supports previous suggestions that MB in cataclysmic variables (CVs) is much weaker than assumed in the standard evolutionary model, unless mass transfer leads to significant additional angular momentum loss in CVs.

1 INTRODUCTION

Stellar magnetic fields cause mass lost in a stellar wind to co-rotate with the star, such that the terminal specific angular momentum of the wind is larger than that of the stellar surface. This causes single stars to spin-down over time, a process known as magnetic braking (MB; e.g. Schatzman 1962; Weber & Davis 1967; Mestel 1968; Mestel & Spruit 1987). In close binaries, tides lock rotation periods to the orbital period, preventing spin-down of the component stars. In this case, the angular momentum debt of the stellar wind is instead repaid by loss of orbital angular momentum, ultimately shrinking the binary orbit and increasing the spin rate of the component stars (e.g. Verbunt & Zwaan 1981).

MB provides a mechanism to bring detached binaries into contact and is thus an important ingredient in the evolution of close binaries. It determines the lifetime, mass transfer rate, and many of the basic observables of cataclysmic variables (CVs; Warner 2003; Knigge, Baraffe & Patterson 2011). Only in short-period, CVs with P_orb ≲ 2 h is MB thought to become subdominant to gravitational wave radiation in driving angular momentum loss. Even at short periods, there is evidence of additional angular momentum loss beyond gravitational waves; MB is one possible explanation (Knigge et al. 2011). MB plays a similarly important role in the formation and evolution of low-mass X-ray binaries (LMXBs) and main-sequence contact binaries (Podsiadlowski, Rappaport & Pfahl 2002; Li, Han & Zhang 2004).

Despite its importance for binary evolution, MB is poorly understood. The most widely used MB prescriptions are empirical, and were motivated by early observations of the rotation rate evolution of single, solar-type stars. These prescriptions are often used in models of close binaries hosting a compact object and a rapidly-rotating, low-mass star. When extrapolated to the low masses and rapid rotation rates relevant for CVs and LMXBs, prescriptions for the MB torque in the literature differ by several orders of magnitude (e.g. Knigge et al. 2011). The prescription that is by far the most widely used in binary evolution studies (Verbunt & Zwaan 1981; Rappaport, Verbunt & Joss 1983) was calibrated to pioneering measurements of the rotation rates of solar-type stars (Skumanich 1972; Smith 1979). More recent studies have found evidence that the rotation period evolution of low-mass main-sequence stars is more complicated than suggested by those early works, particularly at high-rotation rates (e.g. Stauffer et al. 1997; Barnes 2003; McQuillan, Aigrain & Mazeh 2013; Newton et al. 2016; van Saders et al. 2016).

Attempts have also been made to infer a MB law directly from the observed mass transfer rates and period distribution of CVs (e.g. Patterson 1984). However, such efforts are complicated by several factors: (a) mass transfer rates in CVs are difficult to measure and vary on both short and long time-scales, (b) a wide range of mass transfer rates are observed at fixed orbital period in the period range where MB is expected to drive mass transfer (e.g. Townsley & Gänsicke 2009; Pala et al. 2017, 2022), (c) CVs may experience additional angular momentum losses as a consequence of the mass transfer process, which are not easily decoupled from MB, and (d) observed CV samples are affected by complicated selection effects because their observable properties depend strongly on mass transfer rate, and thus on the strength of MB itself.

In this paper, we constrain the strength of MB in close binaries by studying the period distribution of short-period main-sequence binaries. This approach has two clear advantages over previous works: (1) the stars involved are rapidly-rotating, low-mass stars in tidally locked binaries; i.e, the same kind of stars found in CVs and LMXBs, and (2) the selection function of a search for eclipsing main-sequence binaries is straightforward to model. High-quality light-curves from ZTF and astrometry from Gaia allow us to probe lower masses and shorter periods than was possible until recently.

Many previous studies have assembled large samples of EBs and studied their period distribution (e.g. Farinella et al. 1979; Giuricin, Mardirossian & Mezzetti 1983, 1984; Maceroni & van’t Veer 1991; Rucinski 1998, 2007; Maceroni & Rucinski 1999; Derekas, Kiss & Bedding 2007; Norton et al. 2011). However, these studies have had poorly understood selection effects and in most cases included binaries with a wide range of masses and evolutionary states, making population modelling difficult. As expected for flux-limited samples, most studies have been dominated by binaries with components more luminous than the Sun, in which MB is expected to be less important. The primary innovations of our work are well-known distances and stellar masses, a reasonably well-understood selection function, and sensitivity to faint, low-mass sources down to the bottom of the main sequence.

This paper is organized as follows. Section 2 summarizes classical MB theory and the relation between the MB law and the orbital period distribution. Section 3 describes our search for detached eclipsing binaries (EBs) using Gaia and ZTF. Section 4 presents our constraints on the incompleteness-corrected period distribution and compares to predictions of several magnetic braking models. We discuss the implications of our results for CVs and LMXBs in Section 5, and we summarize the main results in Section 6. Additional details about our modelling and EB sample are provided in the Appendices.

2 MAGNETIC BRAKING THEORY

2.1 Classical models

For a single star with constant mass and radius, the angular momentum reservoir is rotation with J_rot ∝ MR²/P_rot. Equation (3) is then a differential equation that can be solved for P_rot(t): for a given initial rotation period, it implies a spin-down with asymptotic |$\dot{P}_{{\rm rot}}\propto t^{-1/2}$|⁠. If their spin-down evolves according to equation (3), the rotation period distribution of a sample of stars with a uniform age distribution will be skewed toward long periods, since individual rotation periods evolve more slowly at long P_rot.

The RVJ MB law has been extremely widely used in models of close binary evolution over the last four decades (e.g. Kolb 1993; Tauris & Savonije 1999; Howell, Nelson & Rappaport 2001; Podsiadlowski et al. 2002; Knigge et al. 2011; Istrate, Tauris & Langer 2014; Paxton et al. 2015; Kalomeni et al. 2016; Belloni et al. 2018; Schreiber et al. 2021, to list only some of the most influential studies). It should be emphasized, however, that equations (1)–(3) were inferred from a sample of G-type main sequence stars with rotation periods ranging from about 2–30 d. MB is expected to produce significant orbital evolution only at short periods (⁠|$P \lesssim \hbox{few days}$|⁠). Moreover, while equation (1) was motivated by observations of G dwarfs, most of the binaries in which the role of MB is of greatest interest (e.g. CVs and LMXBs) have K and M type donors. In short, applying equation (3) to the donors in CVs and LMXBs requires considerable extrapolation.

2.2 Application to binaries

In close binaries, tides lead to the synchronization of rotation with the orbital period (e.g. Zahn 1977). Observations show that synchronization is efficient in main-sequence binaries with orbital periods below about 10 d with weak dependence on mass (e.g. Lurie et al. 2017). In a tidally synchronized binary, the MB torque due to winds from either star cannot spin the stars down, since their rotation is locked to the orbital period. The result is that angular momentum is removed from the orbit instead, while the component stars remain tidally locked. This leads to a decrease of the orbital period and ultimately a spin-up of the component stars in contrast to the single-star case.

If both components of a binary remove orbital angular momentum according to equation (3) while the binary remains tidally locked, one can predict the evolution of P_orb(t) analytically. Details are given in Appendix  A. We are primarily interested in low-mass main sequence binaries in which both components have lifetimes longer than a Hubble time, so we assume that the component masses and radii do not evolve.

Fig. 1 shows predicted tracks of P_orb(t) for main-sequence binaries with a range of masses and initial periods, assuming the angular momentum loss from each component follows equation (3). We assign radii for each mass using solar-metallicity zero age-main sequence MIST models (Choi et al. 2016), and show two different mass ratios for each primary mass. The predicted inspiral is much more rapid at higher binary masses. For example, a binary containing two |$1{\rm M}_{\odot }$| stars in a 1-d orbit is predicted to come into contact within 0.6 Gyr, while the same inspiral in a binary containing two |$0.35{\rm M}_{\odot }$| stars is predicted to take 9 Gyr. This is primarily a consequence of the strong (∝R⁴) dependence in |$\dot{J}$| predicted by equation (3).

2.3 Other magnetic braking prescriptions

Various alternatives to equation (3) for the MB torque due to a single star have been proposed. We review several models below and compare them in Fig. 2.

2.3.1 ‘Saturated’ magnetic braking

For the saturation threshold, we adopt P_crit = 0.1P_⊙(τ_c/τ_{c, ⊙}). Here τ_c is the convective turnover time-scale, which we calculate from the empirical determination of Wright et al. (2011, their equation 11), and |$P_{\odot } = 28\,{\rm d}$|⁠. The saturation threshold can be more succinctly expressed in terms of the Rossby number, Ro = P_rot/τ_c: saturation occurs at |${\rm Ro} \lesssim 0.1\, {\rm Ro}_{\odot }$|⁠. The convective turnover time-scale – and the period below which magnetic braking is saturated – increases with decreasing stellar mass for main-sequence stars. For the saturation threshold adopted here, the saturation period increases from |$P_{\rm crit} = 2.7\, \rm d$| at |$M=1{\rm M}_{\odot }$|⁠, to |$P_{\rm crit} = 5.4\, \rm d$| at |$M=0.6\, {\rm M}_{\odot }$|⁠, to |$P_{\rm crit} = 11.9\, \rm d$| at |$M=0.3{\rm M}_{\odot }$|⁠. This means that essentially all binaries in the period range where MB sets the period distribution are expected to be in the saturated regime.

Another quasi-saturated MB prescription was presented by Ivanova & Taam (2003), who inferred a MB law from the observed scaling of X-ray luminosity with rotation rate in single stars. There is no explicit saturation of the magnetic field in their model, but their analysis implied that at high-rotation rates an increased fraction of magnetic field lines are closed, leading to less efficient MB. They found |$\dot{J}\propto P_{\rm rot}^{-1.3}$| in the saturated regime, only moderately steeper than the models from S00 and Matt et al. (2015). A stronger scaling with R in their model makes the predicted MB in low-mass stars significantly weaker than in the S00 and Matt et al. (2015) models.

2.3.2 Complexity-modulated magnetic braking

Another class of models, popularized by Garraffo, Drake & Cohen (2015), proposes that MB becomes inefficient as the complexity (i.e. multipole order) of the magnetic field increases. They find that for quadrupolar and higher-order magnetic fields, the mass loss rate and magnetic open flux decline, and the Alfven surface shrinks, compared to a simple dipolar field. All these factors diminish angular momentum loss. Observations indeed indicate that the field complexity of low-mass main-sequence stars increases with rotation rate (e.g. Donati & Landstreet 2009), potentially leading to less efficient MB for fast rotators.

The G16 MB law results in a torque with a local maximum at a period of order a few days (Fig. 2). The model can thus predict a bimodal period distribution for young clusters, as observed (e.g. Garraffo et al. 2018). Stars are born with a range of rotation rates after decoupling their rotation from their protostellar discs. Those with rotation periods below the local maximum are predicted to experience a relatively weak braking torque due to their complex magnetic fields, and thus can remain at high-rotation rates for a long time. However, the torque strengthens as the stars spin-down, eventually causing them to move rapidly from the fast-rotator to the slow-rotator sequence when they are near the maximum torque. Such a torque is similar to that required by the phenomenological ‘metastable dynamo’ spin-down model introduced by Brown (2014).

2.3.3 CARB magnetic braking

For low-mass main-sequence stars, equation (15) predicts a torque that is significantly stronger than the other models we consider (see Fig. 2). The scaling with period is similar to that in the RVJ model, transitioning from |$\dot{J}\propto P_{\rm rot}^{-11/3}$| at slow rotation rates to |$\dot{J}\propto P_{\rm rot}^{-7/3}$| at the fastest rotation rates. V19 found that this relation could explain the high-mass transfer rates inferred for persistent LMXBs with low-mass donors, for which the predicted rates for other MB models were too low (Van, Ivanova & Heinke 2019). The model was also tested by Chen et al. (2021) and Soethe & Kepler (2021), who found that it could better explain the population of millisecond pulsars with helium white dwarf companions than other prescriptions.

2.4 Comparison of different models

We compare the different MB laws discussed above in Fig. 2 for equal-mass main-sequence binaries with a range of masses and orbital periods. We use the same values of τ_c and the solar rotation period in all models for consistency, leading to minor differences between our calculated torques and those in the papers that introduced each model.

Although all models make similar predictions for solar-mass stars with periods longer than a few days – likely because they all considered the Sun as a calibration point – the models diverge at lower masses and higher rotation rates. Both the RVJ and V19 models scale roughly as |$\dot{J} \propto P_{\rm orb}^{-3}$| in the period range of interest, but the V19 model predicts a stronger torque at low masses, mainly because of the |$\dot{J}\propto \tau _c^{8/3}$| scaling. Both the S00 and Matt et al. (2015) models predict a shallower |$P_{\rm orb}^{-1}$| scaling at short periods, leading to weaker MB at the shortest periods, where the G16 model also predicts strong MB suppression. For the models considered here, the torque due to gravitational radiation is almost always negligible compared to MB at all periods where main-sequence binaries are detached.

It is worth noting that all the MB laws we consider, except the one from V19, are empirically inferred. An alternative approach would be to predict the angular momentum loss theoretically based on models for the stars’ winds and magnetic fields. This is in principle straightforward to do, and toy models developed elsewhere for stellar magnetospheres provide some conceptual guidance. But in practice, it is difficult to predict from first principles how the magnetic field strength and geometry, and the mass loss rate due to winds, change as a result of rotation (e.g. Taam & Spruit 1989; Garraffo et al. 2015).

2.5 Sensitivity of the main-sequence binary period distribution to MB

The sensitivity of the period distribution to MB is illustrated schematically in Fig. 3, where we compare expected orbital evolution and period distributions of a population of short-period binaries with different MB laws. We consider equal-mass binaries containing |$0.75{\rm M}_{\odot }$| stars with a uniform distribution of ages, approximately as expected for a stellar population with constant star formation rate.

For the RVJ MB law, the period distribution at short periods is predicted to exhibit a deficit of short-period binaries with |${\rm d}N/{\rm d}P_{\rm orb}\propto P_{\rm orb}^{7/3}$| (equation 8). Fig. 3 shows that this conclusion is not very sensitive to the adopted birth period distribution. Even if binaries are preferentially born at the closest periods, their predicted lifetime is sufficiently short that an equilibrium is reached and the period distribution is set by MB. On the other hand, for saturated MB, the orbital period evolves almost linearly with time, and so the observable period distribution is also sensitive to the birth period distribution. In practice, the ‘birth’ distribution for short-period binaries refers to the distribution after various dynamical effects (such as three-body interactions in birth environments; Fabrycky & Tremaine 2007) have occurred. The main point we wish to convey from Fig. 3 is that with a MB torque that increases strongly at short periods (e.g. RVJ or V19), a strong deficit of short-period detached binaries is expected. We now test this prediction.

3 DATA

3.1 Search for eclipsing binaries with ZTF and Gaia

Many large samples of EBs have been collected in the literature with more than half a million binaries catalogued to date (e.g. Paczyński et al. 2006; Prša et al. 2011; Soszyński et al. 2016; Jayasinghe et al. 2018; Chen et al. 2020; Gaia Collaboration et al. 2022; Rowan et al. 2022). Our goal here is not to produce a large sample, but a pure one, and one that has a reasonably well-understood selection function.

We use light-curves from ZTF (Bellm et al. 2019). ZTF provides densely sampled, high-quality photometry across the northern sky with improved sensitivity to faint sources compared to other wide-field surveys. The typical r-band uncertainty of individual photometric epochs is 0.03 mag for 18th mag sources and 0.06 mag for 19th mag sources.

3.1.1 Parent samples

We assembled four different parent samples from Gaia DR3 within which to search for EBs. The different samples correspond to different extinction-corrected absolute magnitude (and thus, mass) ranges and are summarized in Fig. 4. The selection criteria were as follows:

• Magnitude limit: 12.5 ≤ G ≤ G_max, where G is the phot_g_mean_mag reported in the Gaia archive. We set G_max = 19.5 for the lowest-mass sample and G_max = 19 for the three higher-mass samples, since the lowest-mass stars are faint. The bright limit is to avoid saturated sources, and only removes a small fraction of stars in the highest-mass bin.

• Distance limit and uncertainty: We required Gaia DR3 parallax_over_error > 5, so all sources have well-constrained distances and luminosities. We also require parallax > ϖ_min, where |$\varpi _{\rm min}=2\, \rm mas$| (500 pc) for the three highest-mass bins, and 1 mas (1 kpc) for the lowest-mass bin.

• Colour: To exclude most white dwarf binaries, cataclysmic variables, and related objects, we remove objects below the main sequence, requiring M_G < 3.25 (G_BP − G_RP) + 9.625. This does not remove binaries with cool and faint white dwarfs; these are identified later based on light-curve shape.

• Astrometric fidelity: We require fidelity_v2 > 0.75, where fidelity_v2 is the Gaia astrometric reliability diagnostic calculated by Rybizki et al. (2022) to filter out sources with spurious parallaxes.

• Uncontaminated photometry: To minimize spurious variability in the ZTF light-curves due to blended sources, we require norm_dG < −4, where norm_dG is the blending diagnostic calculated by Rybizki et al. (2022). This removes sources that have bright neighbors: e.g. a companion at least 0 magnitudes brighter at a distance of 4 arcsec; a companion less than 2 magnitudes fainter at a distance of 2 arcsec, or a companion at least 2 magnitude brighter at a distance of 6 arcsec.

• Light-curve sampling: We required each object to have at least 100 epochs of clean photometry in the r-band as of ZTF DR8, where for ‘clean’ photometry we required catflags < 32768.

• Sky position: We excluded sources with declination δ < −28°, corresponding roughly to the southern pointing limits for ZTF and PanSTARRS (which the extinction map relies on).

We retrieved the public ZTF light-curves in r-, g-, and i-bands of all Gaia sources passing the cuts above, excluding data with catflags ≥ 32768. To increase the number of usable photometric epochs, we combined photometry from the primary and secondary ZTF grids (see Bellm et al. 2019). In cases where there are multiple ZTF object IDs within <0.5 arcsec of a Gaia source, we concatenate photometry from both IDs into a single light-curve.

We estimated extinctions for each source using the 3D dust map from Green et al. (2019), assuming A_G = 2.67E (g − r) and E (G_BP − G_RP) = 1.33E (g − r). We neglect the weak colour-dependence of these transformations and assume a universal extinction curve. Because the sources in our sample are nearby, extinction is generally modest.

We divided the lower main sequence into four bins of extinction-corrected absolute magnitude, M_{G, 0}. Basic properties of the resulting parent samples are summarized in Fig. 4. The four bins correspond roughly to primary spectral types M3.5–M8 (red, M_{G, 0} > 10, |$M_1 \lesssim 0.3{\rm M}_{\odot }$|⁠); M3.5–M2.5 (yellow, 9 < M_{G, 0} < 10, |$0.3 \lesssim M_1/{\rm M}_{\odot } \lesssim 0.4{\rm M}_{\odot }$|⁠); M2.5–K7 (green, 7 < M_{G, 0} < 9, |$0.4 \lesssim M_1/{\rm M}_{\odot } \lesssim 0.65{\rm M}_{\odot }$|⁠), and K7–K0 (blue, 5 < M_{G, 0} < 7, |$0.65 \lesssim M_1/{\rm M}_{\odot } \lesssim 0.9\, {\rm M}_{\odot }$|⁠). We estimated these ranges from solar-type MIST isochrones (Choi et al. 2016), assuming the secondary brightens the unresolved source by 0.5 mag.

The divisions between bins are somewhat arbitrary, but were chosen such that the lowest-mass bin is expected to contain only fully-convective stars,³ and comparable numbers of EBs are found in each subsample. The mapping from absolute magnitude to mass is not unique, because a range of mass ratios (and ages and metallicities) are found in the sample. In Fig. 4, we show primary mass estimates calculated from M_{G, 0} assuming q = 1, solar metallicity, and an age of 1 Gyr. In the mass range of interest, a lower mass ratio of q ≲ 0.5 would imply a primary mass typically 10 per cent larger. By design, all stars in all samples are expected to have main-sequence lifetimes exceeding the age of the Universe and are expected to be near the zero-age main sequence.

3.1.2 Search for eclipsing binaries

Having obtained ZTF light-curves for ≈6 million Gaia sources, we used a multistep process to identify EBs. We began by computing a box least squares (BLS; Kovács, Zucker & Mazeh 2002) periodogram for each r-band light-curve, using a frequency grid running from 45 min to 10 d. The 45 min lower limit was chosen to be half of the shortest plausible orbit inside which a main-sequence star would fit inside its Roche lobe. We used the GPU-accelerated implementation of BLS in cuvarbase⁴ and a linear frequency grid with an oversample factor of 30; this typically yields about 10⁶ trial frequencies per source. We used a logarithmic grid of 20 eclipse phase durations between 0.005 and 0.5, selecting the ‘best’ duration as the one that yielded the highest likelihood. In order to mitigate the effects of flares and spuriously bright photometric epochs, we masked points that exceeded the median flux by more than 1.5 times the inter-quartile range.

Because low-mass main-sequence stars are red, the r-band light-curves have smaller uncertainties and more photometric points in most cases than the g-band light-curves. The i-band is generally sparsely sampled. We thus only calculated periodogram metrics for the r-band, though we also inspected the phased g- and i-band data during the visual inspection phase.

In addition to the BLS periodogram, we calculated the conditional entropy periodogram (CE; Graham et al. 2013) on the same frequency grid. Unlike the BLS periodgram, which assumes a specific light-curve shape, CE quantifies how much the level of disorder improves in phased versus unphased light-curves. In particular, CE is often more robust than traditional periodograms against spurious periods associated with the observing cadence. We used 20 phase bins and 10 magnitude bins in calculating CE.

We then selected candidate EBs as sources satisfying the following:

• Best-fit BLS eclipse depth of at least 10 per cent in flux.

• BLS power of at least 50.

• At least six distinct eclipses covered by the r-band data. Here the best-fit BLS model determines whether a data point is ‘in eclipse’. The number of ‘in eclipse’ datapoints is often larger than the number of eclipses covered by the data, because multiple epochs can fall into a single eclipse.

• Period not close to 1 sidereal day: |$|P- 1\,{\rm d}| > 0.005$|⁠. This eliminates a large fraction of spurious detections, including sources displaying aperiodic long-term variability. It causes the search to miss real EBs with periods between 0.995 and 1.005 d, but this is of little consequence for our analysis.

• For sources with best-fit periods less than one day, we additionally require that the CE associated with the best-fit BLS period is at least 10 per cent smaller than the median of the CE periodogram. We do not impose this requirement for longer periods, because there the eclipse phase duration becomes short enough that the eclipse is sometimes not well-resolved with the 20 phase bins used to calculate CE.

These criteria were chosen with the goal of minimizing both the number of spurious detections and the number of genuine EBs missed. We chose the detection thresholds based on injection and recovery tests with simulated eclipse light-curves injected into the data (see Section 3.1.7).

We also experimented with using Gaia variability statistics to identify variable stars, as described by Guidry et al. (2021). We found that cuts on variability statistics successfully identify most EBs with P_orb < 1 d. However, at longer periods, many sources that are unambiguous EBs have unremarkable Gaia variability statistics and would be missed with cuts on such statistics. This reflects the fact that the eclipse duty cycle falls with increasing period, such that the finite sampling of Gaia light-curves might never catch longer-period binaries in eclipse.

3.1.3 Visual inspection

The final stage of our EB search is visual inspection of candidates identified in the automated search described above. In short, we looked at the phased ZTF light-curves of each source in all available bands to decide whether the light-curve shape could best be attributed to a detached main-sequence EB or to another source of variability, either astrophysical or instrumental.

In Figs 5 and 6, we show examples of periodic variables that are (Fig. 5) and are not (Fig. 6) likely detached main-sequence EBs. All the examples we show are from the lowest-mass bin. Objects in the higher-mass bins are generally brighter and thus have higher-quality light-curves. Common false-positives that were removed during our visual inspection include detached white dwarf  + main sequence binaries, cataclysmic variables, non-eclipsing ellipsoidal variables, contact binaries, rotating stars with spots, and young stellar objects. In addition to these astrophysically periodic variables, we also removed sources in which variability was likely due to spurious instrumental effects, such as stars falling on a bad column in the ZTF detector. Many of the short-period variables we reject as rotators are likely still binaries, in many cases with the rotation period tracking the orbital period (e.g. Simonian, Pinsonneault & Terndrup 2019). However, we limit our sample to detached EBs in the interests of maintaining a selection function that is straightforward to model.

We also experimented with a range of automated tools for distinguishing EBs from other periodic variables, including the random forest classifier Upsilon (e.g. Kim & Bailer-Jones 2016) and a variety of bespoke cuts on light-curve summary statistics. Such approaches can achieve reasonably good purity and completeness in a small fraction of the time required by human classifiers. However, all approaches we tested performed worse than human classifiers. Because our main goal, constraining the intrinsic period distribution at short periods, can be achieved with a modest sample size (and because we find visual inspection of light-curves a pleasantly meditative activity), we opted for the more labour-intensive approach of inspecting candidates visually.

The false-positive rate of our initial automated selection varies significantly between absolute magnitude bins. In the highest-mass bin, about half of all candidates that pass our automated selection criteria also pass visual inspection. In the lowest-mass bin, the success rate of the automated selection is lower than 20 per cent. This owes primarily to the poorer signal-to-noise ratio in the fainter light-curves, which results in many spurious periodic detections due to correlated noise.

3.1.4 Removing contact binaries

At the shortest periods in each absolute magnitude subsample (e.g. |$P_{\rm orb}\lesssim 0.1\,{\rm d}$| for M_{G, 0} > 10; |$P_{\rm orb}\lesssim 0.3\,{\rm d}$| for 5 < M_{G, 0} < 7), we find candidate EBs with smoothly varying light-curves (i.e. no clear ingress or egress and no flat region between eclipses). In most cases, this is due to tidal deformation of component stars that either nearly or completely fill their Roche lobes. In at least some such binaries, mass transfer leads to a long-lived, stable common envelope phase during which binaries are observed as ‘W UMa’ contact systems (e.g. Lucy 1968). Since mass transfer complicates the angular momentum evolution, we exclude such binaries from the sample. We find that they are common in the highest-mass bin, where there is a clear pileup in the period distribution at |$P_{\rm orb} \lesssim 0.3\,{\rm d}$|⁠, but are rare (and possibly completely absent) in the three lower-mass bins. Jiang et al. (2012) have suggested that contact binaries are unstable, or at least relatively short lived, at primary masses below |${\approx}0.6{\rm M}_{\odot }$|⁠. Our period distributions are consistent with such a scenario, which predicts a pileup of binaries at the period where main-sequence stars overflow their Roche lobes. Although there are some systems with light-curves resembling contact binaries in the lower mass samples (e.g. the rejected source in the centre panel of Fig. 6), follow-up by Drake et al. (2014) found that at least a large fraction of such systems are M dwarfs with white dwarf companions.

Although we attempt to exclude contact binaries, some cases are ambiguous: detached binaries with components that are nearly Roche lobe filling can have light-curves quite similar to contact binaries. Fortunately, the orbital period at which a main-sequence star will overflow its Roche lobe is primarily a function of mass. To prevent unrecognized contact binaries from biasing our results, we only consider binaries in subsequent sections that have periods well above the longest Roche lobe overflow period for their mass.

3.1.5 Period aliases

As is common in analysis of EBs, it is sometimes unclear whether the peak of the BLS periodogram corresponds to P_orb or P_orb/2. To decide, we phase each light-curve to both periods during our visual inspection. In most cases, a secondary minimum is evident with different depth from the primary minimum, removing the ambiguity.

For cases where no secondary eclipse is evident, we performed simulations of eclipsing binaries using the ellc code (see Section 3.1.7). We found that in all but the highest-mass bin, a stellar companion large enough to cause an eclipse with >10 per cent depth would also produce a detectable secondary eclipse; this reflects the fact that temperature varies slowly with mass along the lower main sequence. For the highest mass bin, there is some degeneracy between an edge-on orbit with q ≈ 0.3 and a grazing orbit with q = 1 and twice the orbital period. Given the precise Gaia distances to objects in our sample, we can distinguish between these two scenarios based on CMD position: a q ≈ 0.3 companion will increase the luminosity of the unresolved binary by only a few per cent compared to the primary in isolation, so such binaries are expected to be found near the main sequence. On the other hand, a q ≈ 1 companion will double the luminosity of the unresolved source, shifting it ≈0.75 mag above the main sequence. This allows us to determine the orbital period unambiguously in most cases, though we expect that some period aliases still exist in the catalogue.

In the course of visual inspection, we identified one short-period EB in which there is no secondary eclipse, but the orbital period is still well-determined due to ellipsoidal variation. This object, Gaia DR3 2085085710293765120, has the 2nd shortest period in our sample, with P_orb = 1.90 h. We suspect that the companion is a brown dwarf and will analyse this object in more detail in future work.

3.1.6 Eccentricities

To assess whether EBs have circular or eccentric orbits, we checked the relative spacing of the primary and secondary eclipses. In EBs with circular orbits, these will always be separated by half an orbit, whereas in eccentric EBs, the spacing will in general be unequal. Searching the full catalogue, we identified only a few short-period systems with clearly eccentric orbits. The shortest-period eccentric system is Gaia DR3 518184643170677376 with P_orb = 2.15 d and absolute magnitude M_G, 0 = 8.2. Most importantly for our purposes, we did not find any eccentric EBs in the period range where we expect MB to be efficient (P_orb ≲ 2 d). This validates our assumption that binaries at these periods are tidally synchronized, a requirement for MB to remove orbital angular momentum.

3.1.7 Completeness

To quantify the selection function of our search, we performed injection-and-recovery simulations using the ZTF light-curves. For each absolute magnitude subsample, we simulated a population of close binaries, injected them into a random subset of the observed light-curves (chosen to have apparent magnitudes and cadences that are representative of the parent sample), and then repeated the search procedure we used with the real data, including both the automated search and visual inspection. We assumed a |$p(i)\, {\rm d}i=\sin (i)\, {\rm d}i$| distribution of orbital inclinations and a uniform mass ratio distribution |$\mathcal {U}(q_{\rm min},1)$|⁠, where |$q_{\rm min}=0.1\, {\rm M}_{\odot }/M_{1}$|⁠. We assumed both components are on the zero-age main sequence with radii calculated from solar-metallicity MIST models (Choi et al. 2016). We calculated model light-curves for each binary using the ellc code (Maxted 2016), assuming a uniform distribution of orbital phases and gravity and limb darkening coefficients appropriate for the component stars (Claret 2004). For each bin of absolute magnitude and orbital period shown in Fig. 7, we proceeded until at least 20 simulated binaries had been recovered.

Fig. 7 presents the results of these experiments. The top panels show the overall recovery rate of all binaries, including those which do not eclipse because their inclinations are too low. The middle and bottom panels show the median and middle 68 per cent of the inclinations and mass ratios of the detected binaries. As expected, the range of inclinations that result in a detectable eclipse is largest at the shortest periods. The range of detectable mass ratios varies only weakly with primary mass and period. At short periods, the detection efficiency is set primarily by the eclipse probability, which scales as |$P_{\rm orb}^{-2/3}$| at fixed mass and radius. The eclipse duty cycle (i.e. the fraction of the orbit spent in eclipse) also scales as |$P_{\rm orb}^{-2/3}$|⁠, and the combination of these factors leads to a roughly |$P_{\rm orb}^{-4/3}$| scaling at longer periods.

3.1.8 Wide tertiaries

It is well-established (e.g. Tokovinin et al. 2006) that most close binaries have distant tertiary companions. Indeed, outer tertiaries likely play a decisive role in the formation of almost all close binaries (e.g. Mazeh & Shaham 1979; Fabrycky & Tremaine 2007). These tertiaries are expected to have a separation distribution similar to that of normal binaries (a lognormal distribution peaking at log (P_orb/days) ≈ 5 with logarithmic dispersion of 2.3), which is truncated at the closest separations due to dynamical stability limits (e.g. Tokovinin 2014).

Given a typical distance of 300 pc and the Gaia angular resolution of ≈1 arcsec, most tertiaries to EBs in our sample will be unresolved at separations less than ≈300 au. This means that, absent other selection effects, we expect roughly 30 per cent of tertiaries to be resolved, and 70 per cent to be unresolved. In the three lower mass bins, the unresolved fraction may be even higher, since the separation distribution is skewed toward closer separations in low-mass binaries compared to solar-type binaries (e.g. Moe & Di Stefano 2017).

There is, however, an additional selection effect: we remove sources with bright, close companions from the parent sample in order to minimize light-curve issues related to blending (Section 3.1.1). This means that sources with a bright companion separated by 1 ≲ ρ/arcsec ≲ 5 are missing from the sample with the exact limit depending on luminosity ratio. As a result, we expect almost half of all systems with a resolved companion (corresponding to |${\approx}15\ \hbox{per cent}$| of all systems) to be removed due to the norm_dG (blending) cut.

To compare the wide binary fraction of sources in the EB sample to that of the Galactic field, we cross-matched the EB sample with the wide binary catalogue constructed by El-Badry, Rix & Heintz (2021a). Considering the subsample with M_{G, 0} > 10, we found that of the 162 EBs with P_orb < 0.5 d, 24 have a wide companion with projected physical separation |$s > 2000\,{\rm au}$|⁠. We focus on companions with |$s > 2000\,{\rm au}$|⁠, because our blending cut will preferentially remove resolved companions with closer separations. This implies a detectable wide binary fraction of ≈ 15 per cent for the short-period EBs. This is consistent with a scenario in which most of the short-period EBs have an outer tertiary.

On the other hand, among the 207 436 sources with M_G > 10 in the Gaia Catalog of Nearby Stars (GCNS; Gaia Collaboration et al. 2021), 2694 are in a wide binary with |$s > 2000\,{\rm au}$| – a wide binary fraction of only 1.3 per cent, compared to 15 per cent for the EBs. Completeness to wide companions is likely higher in the GCNS sample due to its |${\lesssim}100\,{\rm pc}$| distance limit. It is thus clear that tertiaries play an important role in the formation of the short-period EBs in our sample. Most of the tertiaries that are spatially resolved are sufficiently wide that the time-scale for Kozai–Lidov oscillations in the current orbital configurations likely exceed the Hubble time, but this may not have been the case in the past (see also Hwang 2022). Since only companions with separations of 300 ≲ s/au ≲ 1500 are removed by the norm_dG cut, we expect the population properties of EBs in the sample to be reasonably representative of the population.

4 RESULTS

4.1 Summary of the eclipsing binary sample

The results of our search are summarized in Fig. 8. The left-hand panel compares our EB candidates (after visual inspection) to the parent sample, which is shown in black. The apparent discontinuity at M_{G, 0} = 10 is a consequence of the looser distance and apparent magnitude limits we use in the faintest sample. As expected, most of the EB candidates are overluminous compared to the main sequence. Those which are not overluminous are primarily systems with weak secondary eclipses, indicative of an unequal mass ratio system in which the secondary contributes little additional light.

The righ-hand panels show the period distributions of visually vetted EBs. In all samples, these rise steeply toward shorter periods. This is primarily an observational selection effect: the dashed lines in each panel show the detection efficiency for that mass bin (Section 3.1.7) with arbitrary rescaling to match the normalization of the data. Hints of a discontinuity at P_orb = 1 d can be seen in both the observed period distributions and detection efficiencies; this is a result of the conditional entropy cut we impose only for systems with inferred periods below one day.

Table 1 lists all the binaries in all absolute magnitude subsamples. The full catalogue contains 3879 EBs with P_orb < 10 d. The number in each absolute magnitude subsample is shown in Fig. 8. The reported P_orb values represent the peak of the BLS periodogram, or twice this value. We have not attempted to model light-curves in detail or account for apparent period variations due to light traveltime effects in triples, so more accurate periods can likely be inferred with careful modelling. Several of the binaries in the lowest-mass subsample have shorter orbital periods than any other main-sequence binaries discovered to date (e.g. Maceroni & Montalbán 2004; Nefs et al. 2012; Soszyński et al. 2015; Koen 2022; Kurtenkov 2022). The shortest-period system has an orbital period of 109 min. Follow-up analysis of some of the shortest-period systems will be presented in future work.

The three lower-mass subsamples contain many binaries well below 0.22 d ‘cut-off’ widely reported for main sequence binaries (e.g. Rucinski 1992; Paczyński et al. 2006). Indeed, in the lowest-mass bin, the peak of the observed distribution is well below this limit. The cut-off is, however, observed in the period distribution of the highest-mass subsample. We suspect that the cut-off found in previous works thus primarily reflects the fact that these surveys had poor sensitivity to faint, low-mass stars. Some previous works (e.g. Stepien 2006) interpreted the lack of observed binaries with periods below 0.22 d as a result of inefficient MB in M dwarfs.⁵ Our results suggest that this is unlikely to be the full story, since short-period binaries are indeed observed at all periods down to the contact limit at all masses.

In the highest-mass sample, the shortest-period bin in Fig. 8 with P_orb < 0.3 d extends well beyond the plot limits, containing almost 600 EBs. This is very likely due to a pileup of unrecognized long-lived contact binaries, which are absent in the lower-mass bins. This population does not affect our inference in later sections, where for this absolute magnitude bin we only consider EBs with P_orb > 0.35 d.

4.2 Intrinsic period distribution

The top panel shows P_orb < 2 d, roughly the regime in which we expect the effects of MB to be important. As discussed in Section 3.1.4, we only show results for orbital periods at which main-sequence binaries must be detached in order to avoid complications from mass transfer and a potentially long-lived contact phase. This is why the minimum period for which we calculate dN/dP_orb increases from lower to higher mass samples.

The dashed magenta line shows |${\rm d}N/{\rm d}P_{\rm orb}\propto P_{\rm orb}^{7/3}$|⁠, the predicted equilibrium period distribution expected if all binaries evolve according the RVJ MB law (equation 3). In stark disagreement with this prediction, the observed period distributions in all absolute magnitude bins are rather flat. Fitting a power-law model |${\rm d}N/{\rm d}P_{\rm orb}\propto P_{\rm orb}^{\alpha }$|⁠, we find |α| < 0.4 in for all absolute magnitude bins. This remains true if we limit our analysis to P_orb < 1 d. The steepest trend is found in the 9 < M_{G, 0} < 10 bin, where we find |${\rm d}N/{\rm d}P_{\rm orb}\propto P_{\rm orb}^{0.35}$| at P_orb < 2 d, but this trend flattens again at longer periods.

In the middle and lower panels, we zoom out further to show the separation distribution at longer periods. Here the dashed black line shows the canonical lognormal period distribution observed for solar-type stars in the solar neighbourhood (Duquennoy & Mayor 1991; Raghavan et al. 2010; Duchêne & Kraus 2013; Moe & Di Stefano 2017), which has a peak at 10⁵ d, logarithmic dispersion of 2.3 dex, and a total binary fraction of ≈ 0.5. The red dashed line shows the analogous distribution for M dwarfs, which we model with a peak at 10^3.75 d, logarithmic dispersion of 1.3 dex, and binary fraction 0.25 (Fischer & Marcy 1992; Duchêne & Kraus 2013). We note that while dN/dlog P_orb (bottom panel) decreases at short periods, dN/dP_orb (top and middle panels) increases. In the period range probed by our data, the observed period distributions in the two higher-mass bins (early M and K dwarfs) are basically consistent with the solar-type lognormal distribution. In the lowest-mass bin, the binary fraction at all periods probed by the data is lower than for solar-type stars. This is basically consistent with trends found in previous work at longer periods (e.g. Duchêne & Kraus 2013; Moe & Di Stefano 2017).

4.3 Implications for magnetic braking

We now compare the observed incompleteness-corrected period distributions to Monte Carlo simulations, which we carry out using the same approach used in Fig. 3. We assume a uniform age distribution between 0 and 10 Gyr, and experiment with a range of birth period distributions and several of the MB laws summarized in Section 2.3.

4.3.1 Birth period distribution

The birth period distribution – i.e. the period distribution before MB-driven orbital evolution – is uncertain. At the periods where MB is expected to be important, observational constraints are scarce (e.g. Kounkel et al. 2019). Pre-main sequence stars are larger and less dense than main-sequence stars, and thus would overflow their Roche lobes at the orbital periods of greatest interest. This is illustrated in Fig. 10, where we show as a function of age the orbital period at which main-sequence stars with a range of masses would overflow their Roche lobes. We calculate Roche lobe radii using the fit of Eggleton (1983). At an age of 3 Myr – comparable to the lifetime of a typical young stellar association in the Galactic disc – there can be no binaries with |$P_{\rm orb}\lesssim 0.75\,{\rm d}$|⁠. This has led some studies to assume a minimum birth period set by protostar radii at a fixed age. For example, Stepien (2006) modelled the evolution of short-period binaries under the assumption that the minimum birth period is two days.

Even after the dissolution of their birth environments, binary orbits can evolve due to three-body interactions. In particular, excitation of eccentricity through Kozai cycles with subsequent dissipation of orbital energy through tides can shrink orbits after stars have reached the main sequence (‘Kozai cylces with tidal friction’ (KCTF); e.g. Mazeh & Shaham 1979; Fabrycky & Tremaine 2007). Since most of the EBs in our sample likely have distant tertiaries (Section 3.1.8), interactions with them are likely important in setting the birth period distribution. Fabrycky & Tremaine (2007) performed simulations to assess the effects of KCTF on the period distribution of main-sequence binaries, assuming an initial period distribution similar to the canonical lognormal distribution for solar-type stars. Fitting a power law to their predicted final period distribution, we find that at short periods it is well-described by |${\rm d}N/{\rm d}\log P_{\rm orb} \propto P_{\rm orb}^{1/2}$|⁠; i.e. |${\rm d}N/{\rm d}P_{\rm orb} \propto P_{\rm orb}^{-1/2}.$|

Given the uncertainties discussed above, we show results for three different birth distributions:

• Uniform distribution between the contact limit for a main-sequence star (bottom panel of Fig. 10), and five days. The five day upper limit is chosen to be large enough that increasing it has no effect on the distribution at P_orb < 2 d.

• Uniform distribution between 0.75 and 5 d. The lower limit corresponds to the minimum possible period for a pre-main sequence star with age of a few Myr (Fig. 10).

• Power-law distribution, dN/dP_orb ∝ P^−1/2 between three times the contact limit and five days, motivated by the results of simulations by Fabrycky & Tremaine (2007). Here the lower limit of three times the contact period is chosen to represent the period, where tides become very strong, suppressing Kozai cycles. The typical period at which Kozai cycles are suppressed is likely longer than this. However, we do not find any dearth of short-period EBs among kinematically young stars (Appendix  B), suggesting that either MB or KCTF produces a population of short-period EBs on a time-scale that is short compared to the Hubble time.

4.3.2 Simulations

For each choice of birth period distribution and MB law, we simulate a population of 10⁶ binaries with a uniform age distribution between 0 and 10 Gyr, representative of the solar neighborhood (e.g. Cukanovaite et al. 2022). We assume a primary mass distribution matching the observed samples in each bin, and a uniform mass ratio distribution between |$0.1{\rm M}_{\odot }/M_1$| and 1, consistent with observations for close binaries (e.g. Duchêne & Kraus 2013). The lower limit of |$0.1{\rm M}_{\odot }$| is chosen to exclude companions below the hydrogen burning limit, which are rare or absent in the observed sample. We remove binaries from the simulation when their primaries overflow their Roche lobes, after which they would appear as contact binaries or merge. The results are shown in Figs 11–14. Each figure shows a different primary mass bin and compares the results of five MB prescriptions (including no MB) for three different birth period distributions. We normalize the final period distributions for all simulations and for the data, such that the slope is informative, but not the normalization. For the observational constraints, we again exclude the shortest-period bin, only showing systems that are unambiguosly detached.

The clearest message from these experiments is that for any plausible birth distribution, the RVJ and V19 MB prescriptions produce a period distribution that, at the short-period end, increases steeply toward longer periods. The predicted slope of the distribution is |${\sim}P_{\rm orb}^{7/3}$| for the RVJ law, and similar for the V19 law (Fig. 2). The data do not display any such trend in any of the mass bins we analyse and are thus inconsistent with these MB laws.

For the S00 and G16 laws, satisfactory agreement the data can be achieved if the birth period distribution was roughly uniform at short periods and extended to close to the contact limit. The birth period distribution with a sharp truncation at 0.75 d does not reproduce the data for the lowest-mass EBs with these MB laws (Figs 13 and 14), because MB is predicted to be weak in the S00 and G16 models, such that binaries born with P_orb > 0.75 d would not reach contact within a Hubble time. This suggests that either MB is stronger at low masses than assumed in these models (but still has a saturated dependence on orbital period), or another process, such as KCTF, is responsible for shrinking of binary orbits to near the contact limit.

Our comparison of simulations and data in Figs 11–14 is focused on the slope of the period distribution; i.e. the relative number of EBs at shorter versus longer periods. We have not attempted to model the absolute number of binaries in a given period bin, because we do not precisely know the initial binary fraction at the relevant periods (⁠|$P_{\rm orb} < \hbox{few days}$|⁠). The absolute number of short-period EBs can, however, be used as a constraint on some MB models even with very conservative assumptions. In particular, the CARB model introduced by V19 predicts that the pre-contact lifetime of an equal-mass binary with |$M_1= 0.4{\rm M}_{\odot }$|⁠, and initial period of 0.5 d is only 0.002 Gyr – this is a consequence of the very large MB torques it predicts at low masses (see Fig. 2).

The observed period distribution (Fig. 9) implies that about 1 in 2000 stars with |$M\approx 0.4{\rm M}_{\odot }$| is a close binary with P_orb < 0.5 d. Even in the extreme scenario where every 0.4M_⊙ star is born with a companion with P_orb < 0.5 d (a scenario clearly ruled out by measurements of the binary populations in young stellar associations), one would only expect |${\approx}(0.002\, \rm Gyr)/(10\, \rm Gyr)\approx$| 1 in 5000 binaries to be observed with P_orb < 0.5 d if their angular momentum loss were governed by CARB MB. Put simply, the model predicts that short-period EBs should almost instantly inspiral and become contact binaries or merge. It thus cannot explain the longer-term survival of a large population of short-period, low-mass main-sequence EBs. We stress that this is a statement about the typical MB law experienced by most short-period EBs. The data do not rule out a scenario in which a small fraction of binaries experience much stronger-than-average MB.

4.3.3 Other MB models

We briefly comment on the implications of our analysis for two other MB laws: the quasi-saturated model from Ivanova & Taam (2003), and the flexible model from Kawaler (1988).

The Ivanova & Taam (2003) law predicts a MB torque whose scaling with period, |$\dot{J}\propto P^{-1.3}_{\rm orb}$|⁠, is only slightly steeper than in the S00 law. The predicted MB at |$M=1{\rm M}_{\odot }$| is thus similar to that law, and the model can reasonably reproduce the observed period distributions in our highest-mass EB sample. However, at |$M\lesssim 0.5{\rm M}_{\odot }$|⁠, the model predicts MB that is too weak to significantly shrink the periods of most EBs within a Hubble time, even at periods near the contact limit. The model can thus only reproduce the data if the birth period distribution is already flat.

The model from Kawaler (1988, their equation 10) has a flexible parametrization that allow for a wide variety of different scalings with P_orb. We find that with a = 0 and n = 0.5, their parametrization produces a MB law similar to the S00 and Matt et al. (2015) parametrizations in the saturated regime, and thus predicts period distributions similar to those shown for the S00 law. Models with a = 0 have a magnetic field strength that is independent of rotation period, and thus are saturated; they also predict |$\dot{J}\propto P_{\rm orb}^{-1}$|⁠.

5 DISCUSSION

This work is not the first to suggest that a Skumanich-like MB law overestimates the strength of MB in rapidly-rotating low-mass stars. Indeed, many studies of cluster M dwarfs over the past three decades have found that rapid-rotators are more abundant than predicted by the Skumanich relation, which would predict them to spin-down in a few Myr (e.g. Stauffer & Hartmann 1987; Stauffer et al. 1997; Barnes 2003; Newton et al. 2016; Medina et al. 2020). These data are better explained by ‘saturated’ MB laws, in which the magnetic field strength ceases to increase with rotation rate above a critical rotation rate, leading to a weaker MB torque (e.g. Sills et al. 2000). The simplest such models predict a torque that scales as |$\dot{J}\propto P_{\rm orb}^{-1}$|⁠, not far from the |$\dot{J}\propto P_{\rm orb}^{-2/3}$| we infer. Our results thus imply that the same saturation occurs in tidally locked, close binaries.

5.1 Implications for CVs and LMXBs

Two decades ago, Andronov, Pinsonneault & Sills (2003) argued that the observed rotation periods of isolated M dwarfs in clusters are inconsistent with classical MB models assumed in the CV literature, which predict |$\dot{J}\propto P_{\rm orb}^{-3}$|⁠. They advocated a saturated MB torque that scales as |$\dot{J}\propto P_{\rm orb}^{-1}$|⁠, basically consistent with our findings. Andronov et al. found that such a MB torque would require significant revision of the standard evolutionary model for CVs⁶ for two reasons: (a) the evolutionary time-scales and inferred ages of CVs would increase by an order of magnitude and, (b) the donor stars in CVs above the period gap would be in thermal equilibrium and thus would not detach even if MB became weaker at the fully convective boundary. This would render untenable the disrupted MB model for the CV period gap. This tension remains unresolved, and most evolutionary models for CVs in the last two decades have continued to use a |$\dot{J}\propto P_{\rm orb}^{-3}$| prescription above the period gap.

Our constraints on the MB law are consistent with the model assumed by Andronov et al. (2003), at least as concerns the scaling of |$\dot{J}$| with P_orb: the MB law they assumed (which was taken from S00) can reproduce the observed EB period distributions for a suitable choice of birth period distribution, while the RVJ law and other laws that scale as |$\dot{J}\propto P_{\rm orb}^{-3}$| cannot. It thus seems worth revisiting whether other uncertain aspects of the CV standard model can change, such that the observed period distribution and mass transfer rates can be matched with a saturated MB law.

We also do not detect any significant difference between the period distributions of binaries containing fully and partially convective M dwarfs. Various works have proposed that MB ‘turns off’ below the fully convective boundary (e.g. Rappaport et al. 1983; Spruit & Ritter 1983; Howell et al. 2001; Knigge et al. 2011); this is the reason for the ‘period gap’ in the CV standard model. The shortest periods we observe in the fully convective EB sample (down to P_orb < 0.1 d) cannot be the birth periods: pre-main sequence M dwarfs are large and would not fit in orbits with periods less than ≈1 d (e.g. Fig. 10). The abundance of observed binaries at these short periods, including systems in which the presence of a stellar tertiary can be excluded, would appear to suggest that MB continues to operate efficiently below the fully-convective limit.

Previous works (e.g. Schreiber et al. 2010; Zorotovic et al. 2016) have found evidence for a reduction in MB at the fully convective boundary in the period and mass distributions of detached white dwarf + M dwarf binaries. Although the binaries studied in those works are detached, it is not easy to distinguish detached CVs (those which have temporarily paused mass transfer) from pre-CVs (those which have not yet come into contact). Inference about the MB law from detached white dwarf  + M dwarf binaries thus depends on modelling of the CV population and on assumptions about the post-common envelope period distribution.

The rotation period distributions of single stars (both in clusters and in the field) display no discontinuity at the fully convective boundary in tension with expectations in the disrupted MB paradigm. On the other hand, there are multiple lines of evidence that CV donors above the period gap are out of thermal equilibrium (e.g. Beuermann et al. 1998; Knigge 2006), which Andronov et al. (2003) showed is not expected if the only source of angular momentum loss is saturated MB following the prescription from S00.

The long CV lifetimes predicted for a saturated MB law could help explain the lack of period-bouncers in observed CV samples (e.g. Belloni et al. 2020; Pala et al. 2020), but would also increase the fraction of CVs observed to have evolved donors (e.g. El-Badry et al. 2021b).

Clearly, further work is needed to reconcile the apparently discrepant MB laws inferred from CVs and related objects, and those inferred from single stars and main-sequence binaries. One possibility is that MB operates differently in binaries with ongoing mass transfer than in detached binaries. Some phenomenological models have indeed been introduced to add an additional source of angular momentum loss associated with mass transfer (e.g. CAML; King & Kolb 1995), and these models have made significant progress in explaining several aspects of the observed CV population (e.g. Schreiber, Zorotovic & Wijnen 2016). However, among the small number of CVs with well-measured mass transfer rates at periods where MB is expected to dominate, mass transfer rates are generally lower than predicted by the standard model (e.g. Pala et al. 2022), and an additional source of angular momentum loss would exacerbate this. It is possible that the observed mass transfer rates are consistent with the saturated MB laws we examined with the addition of consequential angular momentum loss. For example, much of the total angular momentum loss could occur during or shortly after novae (e.g. Nelemans et al. 2016; Ginzburg & Quataert 2021), or be driven by circumbinary discs (e.g. Spruit & Taam 2001). In this case, the total angular momentum loss in CVs could even be similar to the predictions of RVJ, but this would be a coincidence, and there would still be little reason to expect that law’s predictions to hold over a range of periods.

It should in principle be possible to infer the strength of MB in CVs from measurements of their mass transfer rates. The most reliable estimates of CV mass transfer rates are likely those inferred from temperatures of the accreting white dwarfs, which are assumed to be set by compressional heating and trace the long-term average accretion rate (e.g. Townsley & Bildsten 2003; Pala et al. 2017). The inferred accretion rates from such studies are, however, difficult to interpret. In the period range of 3–5 hr – where MB is expected to be the dominant driver of CV evolution – there is a enormous dispersion in the inferred mass transfer rates. A few systems – mostly nova-like variables just above the period gap – have mass transfer rates even larger than predicted even in the standard model with the RVJ MB law. The remaining systems, mostly dwarf novae with periods above 4 h have mass transfer rates more than an order of magnitude below the predictions of the standard model. A full explanation of this phenomenology is still lacking (Belloni et al. 2020). Some of the systems with unexpectedly high-inferred mass transfer rates could be young systems that recently underwent thermal-time-scale mass transfer or a recent novae, but the data may simply indicate that the strength of MB varies significantly from CV to CV at fixed period.

The situation is similar for LMXBs. At orbital periods |$4\lesssim P_{\rm orb}/\rm hr \lesssim 24$|⁠, there are several persistent LMXBs with near-Eddington accretion rates (e.g. Van et al. 2019). It is the existence of these systems that motivated the ‘boosted’ CARB MB prescription (Van & Ivanova 2019), which is even stronger than the RVJ model and which we found inconsistent with the observed period distribution of main-sequence binaries. However, alongside these high-accretion rate systems, there is a larger population of LMXBs with mass transfer rates that are 2–3 orders of magnitude lower at fixed P_orb. The mass transfer rates of these systems are too low to explain with an unsaturated MB law, except in rather fine-tuned scenarios with an evolved donor just below the bifurcation period (Podsiadlowski et al. 2002; Podsiadlowski, Han & Rappaport 2003). Unlike the CVs with high-mass transfer rates, most of the known persistent LMXBs have sufficiently long periods that the donors must be somewhat evolved, and probably above the bifurcation period. In this case, their high-mass transfer rates could be driven mainly by the donor’s nuclear evolution and expansion, rather than by MB. As with CVs, it is possible that high mass transfer rates in LMXBs are due to angular momentum loss mechanisms that are a consequence of the mass transfer process (e.g. Ma & Li 2009; Chen & Podsiadlowski 2016).

Several millisecond pulsars have been observed with helium white dwarf companion in close orbits, and these systems presumably formed from LMXBs with evolved donors. It has been suggested (e.g. Chen et al. 2021; Soethe & Kepler 2021) that the CARB MB prescription can alleviate the ‘fine-tuning problem’ for the formation of these system (Istrate et al. 2014): CARB produces sufficiently strong MB that neutron star  + main-sequence binaries with a wide range of post-common envelope periods can evolve to short periods. Our results imply that if the MB law is universal, CARB MB is unlikely to be the right solution to that problem, since it dramatically overpredicts the inspiral rate of close main-sequence binaries. Consequential angular momentum loss is also unlikely to help in this case, because stronger MB appears to be required before the onset of mass transfer. Models in which weak MB below the bifucation period (due, e.g. to loss of the donor’s convective envelope) leads to detachment may be more promising.

In summary, there are multiple lines of evidence that MB in close detached binaries does not follow a Skumanich-like law with |$\dot{J}\propto P_{\rm orb}^{-3}$|⁠, and a saturated MB law with |$\dot{J}\propto P_{\rm orb}^{-1}$| is more capable of matching the data. In binaries with ongoing mass transfer, such as CVs and LMXBs, the situation may be more complicated. Observationally-inferred mass transfer rates in these systems are in most cases consistent with theory at short periods (where MB is expected to be subdominant), but not at long periods (e.g. Belloni et al. 2020; Pala et al. 2022). We conclude that the role of MB in driving angular momentum loss in LMXBs and CVs above the period gap is ripe for reassessment. Models in which much of the angular momentum loss in these systems is a consequence of mass transfer provide a potential solution to some, but not all, of these puzzles.

6 SUMMARY

Using light-curves from ZTF, we have constrained the period distribution of detached low-mass main-sequence eclipsing binaries (EBs) with primary masses 0.1 ≲ M/M_⊙ ≲ 0.9 and periods 0.1 ≲ P_orb ≲ 10 d. Although many previous works have assembled large samples of eclipsing binaries using photometric surveys, this is to our knowledge, the first study to model the survey selection function and thus to infer the intrinsic period distribution in this range of masses and periods. We then compared the observed period distribution to simple population models with different assumptions about how magnetic braking (MB) removes angular momentum from binary orbits. Our main results are as follows:

• The period distribution extends down to the contact limit at all masses. Many previous studies have reported a minimum period of ≈0.22 d for main-sequence binaries. This is roughly the period at which a 0.7M_⊙ main-sequence star will overflow its Roche lobe and form a contact binary (Fig. 10). However, lower-mass main sequence binaries can reach significantly shorter periods, down to ≈0.06 d. When we split the data into bins of primary mass, EBs are abundant down to the contact limit in each mass bin (Table 1; Fig. 8). Several of the binaries in the lowest-mass subsample have shorter orbital periods than any other main-sequence binaries discovered to date. This suggests that the 0.22 d cut-off found in previous studies was primarily a result of the dearth of low-mass main-sequence binaries in magnitude-limited samples.

• The incompleteness-corrected period distribution of detached binaries is flat at |$P_{\rm orb} \lesssim 10\,{\rm d}$|. We characterize the search’s sensitivity using injection and recovery simulations (Fig. 7). Unsurprisingly, the search is most sensitive to short-period EBs, which have higher eclipse probabilities and duty cycles. In all mass bins, we consider the incompleteness-corrected period distribution at short periods is nearly flat; i.e. dN/dP ∝ P⁰ (Fig. 9; Fitting a power law model |${\rm d}N/{\rm d}P_{\rm orb}\propto P_{\rm orb}^{\alpha }$|⁠, we find |α| < 0.4 in all subsamples at P_orb < 2 d.). This remains true for detached binaries at all periods down to the contact limit. In contrast, contact binaries pile up at the contact limit. This is consistent with distributions hinted at by previous studies of EBs (Maceroni & Rucinski 1999), but is here unambiguous over a range of masses and robust to selection effects.

  Although the absolute number of close binaries per unresolved source declines with stellar mass (Fig. 9), the shape of the period distribution at close periods varies little between mass bins. There is no obvious difference between the period distributions of binaries containing fully and partially convective stars.

• Observed period distributions are inconsistent with the MB models most widely used in binary evolution calculations: MB removes angular momentum from binary orbits, leading to inspiral. How rapidly this inspiral occurs depends on how the MB torque varies with a star’s rotation period. The most widely used MB laws in the binary evolution literature, which are based on the Skumanich (1972) spin-down relation, predict a torque that scales as |$\dot{J}\propto P_{\rm orb}^{-3}$| (e.g. Verbunt & Zwaan 1981, RVJ). This leads to rapidly accelerating inpiral (Fig. 1) and a dearth of binaries at the shortest periods. At periods where MB is predicted to be efficient, this prediction is insensitive to the initial period distribution (Fig. 3).

  A steeply falling period distribution at short periods is not found in the data in any mass bin (Figs 11–14). We conclude that the data cannot be understood as a result of a Skumanich-like MB law, and the scaling of |$\dot{J}$| with orbital period must be shallower than |$P_{\rm orb}^{-3}$|⁠. A flat observed period distribution is most naturally produced by a MB torque that scales as |$\dot{J}\propto P_{\rm orb}^{-2/3}$|⁠, which would lead to linear period decay with time. For such a MB law, the observed period distribution is also sensitive to the birth period distribution. Given uncertainties in initial short-period binary fraction, the magnitude of the MB torque is more difficult to measure with our approach than how it scales with period. The scalings we adopt (e.g. equation (9)) are calibrated to the rotation rates of single stars in clusters.

  The observed period distributions are more consistent with models in which the MB torque has a shallower period dependence. Saturated MB laws (e.g. S00, Matt et al. 2015) can reproduce the observed period distributions for a range of birth period distributions at the higher-mass end (i.e. K dwarfs). At the lowest masses, these relations predict that MB is too weak to change binary orbits significantly except at the shortest periods (P_orb < 0.5 d), so they can explain the data only if the birth period distribution is populated down to near the contact limit. This cannot occur in an isolated binary because pre-main sequence stars cannot fit in such tight orbits, but may occur with the help of three-body interactions. Similar considerations apply for the G16 model. Finally, the V19 suffers similar challenges in matching the data to the RVJ model with the additional complication that it predicts an unrealistically short lifetime for low-mass detached EBs. While it cannot apply to a majority of EBs, such a model might still explain rare systems with unusually efficient MB.

Our analysis is focused on low-mass main sequence stars with |$M\lesssim 0.9{\rm M}_{\odot }$|⁠. While such stars are representative of the donors in typical CVs and LMXBs, there are also several classes of interacting close binaries that are thought to have formed from donors that initially had higher masses. The naive expectation from models is that MB weakens in main-sequence stars with |$M\gtrsim 1.2{\rm M}_{\odot }$|⁠, which have thin or absent convective envelopes. There is some evidence for such weakening in the rotation period distribution of single stars (e.g. Kraft 1967) and in the mass transfer rates observed in evolved CVs (El-Badry et al. 2021b), but it would be useful to search for evidence of it in the orbital period distribution of higher-mass main-sequence binaries. Unlike the analysis in this work, such an analysis will need to account for significant evolution of the component stars during their main-sequence lifetimes.

SUPPORTING INFORMATION

Table 1. Eclipsing binary sample from ZTF, ordered by orbital period.

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the article.

ACKNOWLEDGEMENTS

We thank Saul Rappaport, Boris Gaensicke, Tom Marsh, Brian Metzger, Ken Shen, Dave Charbonneau, Mercedes López-Morales, Eliot Quataert, Hans-Walter Rix, Selma de Mink, and E. Sterl Phinney, for helpful discussions.

This research was supported in part by the National Science Foundation under Grant No. NSF PHY-1748958.

Based on observations obtained with the Samuel Oschin 48-inch Telescope at the Palomar Observatory as part of the Zwicky Transient Facility project. ZTF is supported by the National Science Foundation under Grant No. AST-1440341 and a collaboration including Caltech, IPAC, the Weizmann Institute of Science, the Oskar Klein Center at Stockholm University, the University of Maryland, the University of Washington, Deutsches Elektronen-Synchrotron, and Humboldt-Universität, Los Alamos National Laboratory, the TANGO Consortium of Taiwan, the University of Wisconsin at Milwaukee, and Lawrence Berkeley National Laboratory. Operations are conducted by COO, IPAC, and UW.

DATA AVAILABILITY

Data used in this study are available upon request from the corresponding author.

Footnotes

REFERENCES

APPENDIX A: ANALYTIC PERIOD EVOLUTION

Here, we solve for the expected orbital period evolution of a binary losing angular momentum via different MB prescriptions.

For the Garraffo et al. (2016) and Van & Ivanova (2019) MB laws, the analytic equations are more difficult to solve, so we calculate the period evolution numerically.

APPENDIX B: KINEMATIC AGES

To roughly estimate the ages of the EBs in our sample, we calculated their plane-of-the-sky tangential velocities, |$v_{\perp }=4.74\, {\rm km\, s^{-1}}\left(\frac{\mu }{{\rm mas\, yr^{-1}}}\right)\left(\frac{\varpi }{{\rm mas}}\right)^{-1}$|⁠. Older stars both have experienced more dynamical heating since their formation and were born from kinematically hotter gas (e.g. Seabroke & Gilmore 2007), and thus have larger typical velocities with respect to the local standard of rest.

Fig. B1 compares the v_⊥ distributions of EBs in each absolute magnitude bin to the v_⊥ distributions of all stars whose light-curves were searched. Overall, the distributions are similar in each absolute magnitude bin. This implies that the EBs in our sample have a similar age distribution to the local field population.

We also searched for potential difference between the period distributions of kinematically youthful and older binaries. Fig. B2 shows the incompleteness-corrected period distributions of binaries with tangential velocities above (red) and below (black) a threshold of |$v_{\perp } = 30\, \rm km\, s^{-1}$|⁠. Any differences between the distributions for young and old binaries (in either normalization or slope) are marginal. Most importantly, we can rule out a scenario in which the shortest-period EBs are missing in the young population (not having had sufficient time to reach short periods) or in the older population (having all been depleted by MB).

Previous work by Hwang & Zakamska (2020) has found an excess of short-period binaries with tangential velocities of 20-40 km s⁻¹. They interpreted this as evidence that typical main-sequence close binaries have lifetimes of a few Gyr, and then merge. In our sample, an excess of EBs with intermediate kinematic ages is weakly evident in the two higher-mass bins, and not in the lower-mass bins. We suspect this is primarily because the Hwang & Zakamska (2020) analysis was (a) focused on solar-type binaries, which are more massive than most of our sample, and (b) dominated by contact binaries, whereas we focus on detached binaries and attempt to explicitly remove contact systems. In contrast to the situation at lower masses, short-period solar-type binaries are dominated by contact systems. This suggests that solar-type contact binaries survive in a stable configuration for long time-scales – likely longer than the detached inspiral time-scale.