Magnetospheric Truncation, Tidal Inspiral, and the Creation of Short-period and Ultra-short-period Planets

We create ${N}_{\mathrm{disk}}={\rm{20,000}}$ disks, each truncated at its own innermost edge at period ${P}_{\mathrm{in}}$. The disks are assumed to be magnetospherically truncated at corotation with their host stars (i.e., we assume disk locking): each value of ${P}_{\mathrm{in}}$ is drawn directly from an empirical distribution of stellar rotation periods ${P}_{\star }$. A given MC model utilizes rotation periods measured for one of three stellar clusters: the Orion Nebula Cluster (ONC, age of 1 Myr), NGC 2362 (5 Myr), and NGC 2547 (40 Myr). For MC models of planets around FGK dwarfs, we draw ${P}_{\mathrm{in}}={P}_{\star }$ using stars with masses $0.5\,{M}_{\odot }\lt {M}_{\star }\leqslant 1.4\,{M}_{\odot }$ (upper panel of Figure 2). For MC models of planets around M stars, we do the same for stars with masses $0.2\,{M}_{\odot }\lt {M}_{\star }\leqslant 0.5\,{M}_{\odot }$ (lower panel of Figure 2).

Each disk is taken to contain ${N}_{\mathrm{planet}}=5$ planets. Their initial periods—initial as in before stellar tides have acted—are decided according to one of two schemes. "In situ" formation models are very simple: the pre-tide planet periods are distributed randomly and uniformly in log period from ${P}_{\mathrm{in}}$ to ${P}_{\mathrm{out}}=400$ days. The logarithmic spacing befits planetary systems that form by in situ, oligarchic accretion, whereby planets are separated by a multiple number of Hill radii (e.g., Kokubo & Ida 2012). Our specific choices of ${N}_{\mathrm{planet}}=5$ and ${P}_{\mathrm{out}}=400$ days affect only the normalization of the occurrence rate profile; this normalization is adjusted later to fit the observations (see the captions to Figures 4 through 8). The aim of this paper is to reproduce the breaks and slopes of the observed occurrence rates, not their absolute normalizations.

The alternative "disk migration" models start the same way as the in situ models, but transport planets to shorter pre-tide periods according to the Type I migration rate:

$$\begin{eqnarray}&&\dot{a}=\displaystyle \frac{2{T}_{{\rm{L}}}}{{M}_{{\rm{p}}}{\rm{\Omega }}a},\end{eqnarray} \tag{ 1 }$$

where a is a given planet's orbital radius, ${M}_{{\rm{p}}}$ is its mass, Ω is its orbital angular frequency, and T_L is the combined Lindblad and corotation torque exerted by the disk on the planet (e.g., Kley & Nelson 2012):

$$\begin{eqnarray}&&{T}_{{\rm{L}}}=-\xi {\left(\displaystyle \frac{{M}_{{\rm{p}}}}{{M}_{\star }}\right)}^{2}{\left(\displaystyle \frac{a}{h}\right)}^{2}{{\rm{\Sigma }}}_{\mathrm{gas}}{a}^{4}{{\rm{\Omega }}}^{2},\end{eqnarray} \tag{ 2 }$$

where $\xi =(1.36+0.62\beta +0.43\gamma )$ is an order-unity constant calibrated from simulations (D'Angelo & Lubow 2010), and β and γ characterize an assumed power-law gas disk of surface density

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{\mathrm{gas}}={{\rm{\Sigma }}}_{\mathrm{gas},0}{(a/{a}_{0})}^{-\beta }\exp (-t/{t}_{\mathrm{disk}})\end{eqnarray} \tag{ 3 }$$

and temperature

$$\begin{eqnarray}&&T={T}_{0}{(a/{a}_{0})}^{-\gamma }.\end{eqnarray} \tag{ 4 }$$

Other variables include the disk scale height $h=\sqrt{{kT}/\mu {m}_{{\rm{H}}}}/{\rm{\Omega }}$, the gas mean molecular weight μ, the mass of the hydrogen atom ${m}_{{\rm{H}}}$, and Boltzmann's constant k.

For simplicity, regardless of whether the host stars modeled are FGK dwarfs or M dwarfs, we fix ${a}_{0}=0.1\,\mathrm{au}$, ${T}_{0}=1000$ K, $\gamma =1/2$, $\mu =2$, and ${M}_{{\rm{p}}}=5\,{M}_{\oplus }$. We fix ${M}_{\star }=1\,{M}_{\odot }$ for FGK dwarfs and ${M}_{\star }=0.5\,{M}_{\odot }$ for M dwarfs. The errors introduced by fixing these variables are subsumed to some extent by large uncertainties in the disk parameters β, ${{\rm{\Sigma }}}_{\mathrm{gas},0}$, and the dissipation timescale ${t}_{\mathrm{disk}}$. We now describe and comment on our choices for these three parameters:

• 1.  
  For a given MC model, we set β to one of three values: $\beta \in \{1,2,3\}$. If we assume that a disk's gas traces its solids (i.e., a spatially constant gas-to-solids abundance ratio), then these β-values encompass the range of quoted slopes in the literature, as inferred from observations of exoplanets comprising mostly solids: ${{\rm{\Sigma }}}_{\mathrm{solid}}\propto {a}^{-1.5}$ to ${a}^{-2.4}$ (e.g., Hansen & Murray 2012; Chiang & Laughlin 2013; Schlichting 2014; Hansen 2015). Note that ${{\rm{\Sigma }}}_{\mathrm{solid}}\sim {dN}/d\mathrm{log}P\times {M}_{{\rm{p}}}/{a}^{2}\propto {a}^{-2}$ follows immediately from ${dN}/d\mathrm{log}P\propto {P}^{0}$ (this flat slope is seen at long periods in Figure 1) if the characteristic planet mass ${M}_{{\rm{p}}}$ does not vary much with distance.
• 2.  
  We experiment with two normalizations for the gas surface density at ${a}_{0}=0.1\,\mathrm{au}$: ${{\rm{\Sigma }}}_{\mathrm{gas},0}\in \{40,400\}$ g cm⁻². These choices yield rather gas-poor disks, depleted by 3–4 orders of magnitude relative to solar-composition "minimum-mass nebulae." Such highly gas-depleted environments are motivated by the need to defeat gas dynamical friction before sub-Neptune cores can form by the mergers of protocores (Lee & Chiang 2016, their Figure 5). We assume that Type I migration proceeds normally in these gas-poor disks and discuss this assumption in the Appendix.
• 3.  
  For the gas exponential decay timescale, we try ${t}_{\mathrm{disk}}\in \{1,10\}$ Myr, which overlaps with the range of typical disk lifetimes reported from observations (Mamajek 2009; Alexander et al. 2014; Pfalzner et al. 2014).

Having assigned all parameters for our disk migration model, we then migrate every planet inward by integrating Equation (1) over a time interval of $10\times {t}_{\mathrm{disk}}$ (10 e-foldings of the disk gas density).

2.2.1. What Happens when Planets Migrate into the Innermost Edge

The first planet to migrate by disk torques toward the disk truncation radius will park in its vicinity—either interior to the disk edge where the last of the principal Lindblad torques peters out (Goldreich & Tremaine 1980; see also Lin et al. 1996); at the disk edge because of principal corotation torques (Tanaka et al. 2002; Masset et al. 2006; Terquem & Papaloizou 2007); or exterior to the disk edge following the reflection of planet-driven waves off the edge (Tsang 2011). In our MC model, we assume for simplicity that the first planet that can migrate by disk torques to arrive at the inner edge does so, and stays there (until step 3, after which it and all other planets move further in by stellar tides). Experiments that do not park this first planet at $P={P}_{\mathrm{in}}$ but instead at $P={P}_{\mathrm{in}}/2$ (so that the planet is located at the interior 2:1 resonance with the disk edge) change none of our qualitative conclusions.

We have two options for the planets that subsequently follow this first parked planet. In one class of MC model ("resonance lock"; see Figure 3 for an illustration), we stop a planet when it migrates convergently into mean-motion resonance with its interior neighbor. For simplicity, we consider only 3:2 resonances (the choice of 3:2 is motivated by the observation that, of the small subset of Kepler planets located near resonances, most are situated near 3:2; Fabrycky et al. 2014). In a second class of MC model ("merged"), planets are allowed to migrate unimpeded to the innermost disk edge and assumed to merge there with the first planet. We count the merger product as a single planet, equivalent to all other planets in our final tally of planet occurrence rates (this simplification should be acceptable insofar as observed occurrence rates based on transit data are typically tallied by planet radius, which at fixed bulk density does not change appreciably when a planet doubles or triples its mass).

Zoom In Zoom Out Reset image size

The intention behind these two flavors of MC model is to bracket the range of possible planet–planet interactions near the disk edge. On the one hand, convergently migrating planets can become trapped in a resonant chain (e.g., Terquem & Papaloizou 2007) or linger near resonances (Deck & Batygin 2015). Although observed short-period planets are generally not found in resonance (e.g., Sanchis-Ojeda et al. 2014, their Figure 8), orbital instabilities over the Gyr ages of observed systems can lead to mergers of resonant planets, transforming them into (lower multiplicity) non-resonant systems. Numerical N-body integrations of high-multiplicity systems show that resonances are preferentially unstable; non-resonant neighbors appear to destabilize otherwise stable resonant pairs (Mahajan & Wu 2014; Pu & Wu 2015, their Figure 5). Our "resonance lock" model captures mergers of destabilized resonant pairs insofar as merger products should have the same period distribution, statistically, as their progenitors.

On the other hand, planets that collect near the inner disk edge may never lock into resonances and instead collide into each other—hence our "merged" model. Given the potential complexity of planet–planet and planet–disk interactions near the disk edge (e.g., Tanaka et al. 2002; Masset et al. 2006; Tsang 2011; Deck & Batygin 2015; Liu et al. 2017)—interactions that may be modulated by the looming stellar magnetosphere (e.g., Romanova & Owocki 2016) and/or disk magnetic fields (e.g., Terquem 2003)—it would be unsurprising if planets fail in general to be captured into resonance in this messy environment.

The asynchronous equilibrium tide raised by a planet on its host star causes its orbital semimajor axis to decay at a rate

$$\begin{eqnarray}&&\dot{a}=-\displaystyle \frac{9}{2}{a}^{-11/2}{G}^{1/2}\displaystyle \frac{{M}_{{\rm{p}}}}{{M}_{\star }^{1/2}}\displaystyle \frac{{R}_{\star }^{5}}{{Q}_{\star }^{\prime }}\end{eqnarray} \tag{ 5 }$$

where ${R}_{\star }$ is the stellar radius and ${Q}_{\star }^{\prime }$ is the effective tidal quality factor (Goldreich & Soter 1966; see also Ibgui & Burrows 2009). For every planet, we integrate Equation (5) over a time interval of 5 Gyr to calculate its final "post-tide" orbital period.⁷ Equation (5) assumes that planets remain on circular orbits, and that planet orbital periods are shorter than the stellar rotation period (otherwise planets would migrate outward). The latter assumption becomes safer as stars age beyond the zero-age main sequence and slow their spins. Although the assumption is violated for the planets with the longest periods (say ≳ 10 days), such planets hardly migrate by tides anyway because of the steep dependence of tidal friction on a. We do not consider planet–planet interactions during tidal migration, as tides act to separate orbits: divergently migrating planets can cross resonances but do not lock into them, and the eccentricities that sub-Neptunes impart to one another when crossing resonances are small (see, e.g., Dermott et al. 1988).

For FGK dwarfs, we fix ${R}_{\star }=1\,{R}_{\odot }$ and vary ${Q}_{\star }^{\prime }\in \{{10}^{6},{10}^{7},{10}^{8}\}$, a range that overlaps with estimates based on observations and modeling of hot Jupiters (e.g., Matsumura et al. 2008; Hansen 2010; Schlaufman et al. 2010; Penev et al. 2012). For M dwarfs, we fix ${R}_{\star }=0.5\,{R}_{\odot }$ and vary ${Q}_{\star }^{\prime }\in \{{10}^{4},{10}^{5},{10}^{6},{10}^{7},{10}^{8}\}$. The smaller values of ${Q}_{\star }^{\prime }$ for M dwarfs are motivated by the idea that stellar tidal dissipation mostly occurs in convective zones, which are larger for lower-mass stars (see, e.g., Hansen 2012, their Figure 3).

Having computed the final periods of all planets in all disks, we tally them up to calculate model occurrence rates versus orbital period, using the same period bins as reported by the observations (Fressin et al. 2013 and Sanchis-Ojeda et al. 2014 for FGK dwarfs; Dressing & Charbonneau 2015 for M dwarfs). We construct both probability distribution functions (PDFs) and corresponding cumulative distribution functions (CDFs) to compare directly with their observational counterparts. This comparison is done by eye—we are looking for broad trends and aiming to identify, only in an order-of-magnitude sense, those regions of parameter space that are compatible with the observations. Our choices on how to normalize our PDFs and CDFs are given in the figure captions.

We begin by reviewing the results of MC models for disk migration, first omitting the effects of tidal orbital decay in Section 2.4.1 and then incorporating them in Section 2.4.2. In situ + tidal models are presented in Section 2.4.3.

2.4.1. Disk Migration without Tidal Friction

Figure 4 shows that despite the many free parameters, our disk-migration-only models give practically the same results: overprediction of the number of short-period planets at 0.5 day ≲ P ≲ 10 days, a tendency to underpredict longer-period planets at P ≳ 10 days, and a failure to capture ultra-short-period planets at P ≲ 0.5 day. Disk migration scoops out too many planets at large P to overpopulate small P.

Zoom In Zoom Out Reset image size

That bins at P < 0.5 day are empty follows from $\min {P}_{\star }=0.5$ day (for NGC 2362, the stellar cluster used to make Figure 4), and our assumption that planets do not migrate past the disk truncation period ${P}_{\mathrm{in}}={P}_{\star }$. We have tried relaxing this assumption, allowing planets to migrate to ${P}_{\mathrm{in}}/2$ (i.e., parking them at 2:1 resonance with the disk edge; see Section 2.2.1). This succeeds in creating USPs, but merely extends the problem of overproducing short-period planets (relative to the observations) to P < 0.5 day (data not shown).

The occurrence rate profiles at P ≲ 10 days appear almost identical over a wide range of disk parameters (β, ${{\rm{\Sigma }}}_{\mathrm{gas},0}$, ${t}_{\mathrm{disk}}$) because "all roads lead to Rome": various disk migration histories all end with the same outcome of planets converging onto disk inner edges, whose period distribution is fixed by the distribution of stellar rotation periods. The details of the planet distribution in the vicinity of the edge depend on what assumption we adopt: "merged" models pile up planets at P ∼ 0.5–1 day, while "resonance lock" models spread planets out more evenly in $\mathrm{log}P$ from ∼1–3 days (under our assumption that planets lock into 3:2 resonance; in general, for a $j:(j-1)$ resonance, the pile-up sharpens with increasing j). At longer P ≳ 10 days, there is more variation between disk models. In particular, β = 1 yields such high disk densities at large orbital radius that all planets at long periods are flushed inward by disk migration.

2.4.2. Disk Migration with Tidal Friction

Incorporating tidal inspiral into our disk migration models can yield improved fits to the occurrence rate of USPs at P ≲ 1 day. In Figure 5 we apply tidal decay to one of the disk-migration-only models shown in Figure 4 (${{\rm{\Sigma }}}_{\mathrm{gas},0}=40\,{\rm{g}}\,{\mathrm{cm}}^{-2}$, β = 2, and ${t}_{\mathrm{disk}}=1$ Myr). From the grid of models in Figure 5, it appears that using stellar rotation periods from NGC 2362 and ${Q}_{\star }^{\prime }={10}^{7}$ works best for reproducing USPs, either in merged or resonance lock models (middle column, blue curves).

Zoom In Zoom Out Reset image size

Drawing disk truncation periods from older clusters, which contain more rapidly rotating stars (Figure 2), produces more USPs. In the oldest cluster considered (NGC 2547 at 40 Myr), typical disk truncation periods are so short that tidal inspiral is required to reduce the number of USPs and better match observations. In the youngest cluster considered (ONC at 1 Myr), tidal inspiral is also required, but for the opposite reason: disk truncation periods here are too long and no USPs are generated without tidal migration. The case of NGC 2362 is intermediate and exhibits both behaviors.

All disk+tidal models, however, overproduce planets with periods between ∼1 and 10 days. This is the same problem as seen in disk-migration-only models (Section 2.4.1). Tidal orbital decay has too short a reach to eliminate this excess.

2.4.3. In Situ Models with Tidal Friction

Figure 6 showcases in situ MC models, with and without tidal friction. At P ≳ 3 days, these models fit the observations remarkably well, performing better than the disk+tidal migration models in this period regime. Simply assuming that sub-Neptunes are randomly distributed in their disks (in a uniform $\mathrm{log}P$ sense) avoids the overproduction-at-short-P/underproduction-at-long-P problem introduced by disk migration.

Zoom In Zoom Out Reset image size

Ultra-short-period planets (USPs) can be reproduced, approximately, by tidal orbital decay of planets located near the inner edges of their disks. Drawing the shorter disk truncation periods from older clusters like NGC 2362 and NGC 2547 works best by giving planets a "head start" toward the shortest periods.

Figure 7 zooms in to directly compare one of the better fitting in situ+tide models with its disk+tide counterpart. With respect to the observations, the in situ+tide model does better than the disk+tide model at P ∼ 1–10 days by not overproducing the number of planets. But for that same reason, the disk+tide model does better than the in situ+tide model at P ≲ 1 day: disk migration generates more USPs by supplying more planets for tidal inspiral, and matches the observations at ultra-short periods more closely. Both classes of model overpredict the number of planets at P ∼ 1–3 days.

Zoom In Zoom Out Reset image size

Figure 8 compares model occurrence rates for planets orbiting M dwarfs against observations. The best agreement is obtained for in situ+tidal models calculated by drawing ${P}_{\mathrm{in}}$ from the ONC and by adopting ${Q}_{\star }^{\prime }\sim {10}^{5}$. Using these same assumptions in disk-induced migration models gives worse fits characterized by an excess of planets inside ∼10 days. Appealing to stronger tidal dissipation can remove much of this excess but simultaneously destroys USPs by shuttling them to the stellar surface.

Zoom In Zoom Out Reset image size

We require a significantly smaller ${Q}_{\star }^{\prime }$ to produce a best-fitting in situ+tidal model for M dwarfs (${Q}_{\star }^{\prime }\sim {10}^{5}$) as compared to FGK dwarfs (${Q}_{\star }^{\prime }\sim {10}^{7}$). M stars tend to rotate faster than FGK stars, so stronger tides are needed to remove the excess of planets near ∼1–2 days. The smaller ${Q}_{\star }^{\prime }$ for M dwarfs supports the idea that the bulk of stellar tidal dissipation occurs in convective zones, which are larger for lower-mass stars (e.g., Hansen 2012, and references therein).

Drawing ${P}_{\mathrm{in}}$ from clusters older than the ONC overestimates the number of planets inside 10 days, as stars spin up with age and truncate their disks at shorter periods. This effect is amplified for M stars over FGK stars because M dwarfs spin up more rapidly as they get older; see in Figure 2 how the peak of the period distribution shifts more dramatically with age for lower-mass stars. The excess population is especially difficult to remove when planets enter into and remain in a "resonance lock" near the disk edge. Resonance lock models distribute planets broadly over this period range, beyond the reach of tidal friction in M stars even for ${Q}_{\star }^{\prime }={10}^{4}$.

We have shown in this paper the essential role tides play in creating USPs. Another line of evidence pointing to the influence of tides can be found by examining planetary orbital spacings. Because tidal orbital decay generally proceeds more quickly at smaller periods, interplanetary spacings should be stretched out more (in a fractional sense) at small P than at large P. Just such a signature is observed: USPs are more widely spaced than their longer-period counterparts (Steffen & Farr 2013). In Figure 10 we demonstrate how orbital spacings widen preferentially at the shortest periods because of tides. Comparison of our model with observations shows only qualitative agreement; the empirical data (gray circles and orange line) suggest that tidal friction in reality is more effective at drawing planets apart than our constant-${Q}_{\star }^{\prime }$ theory allows (blue curves). A better theory for tidal friction, one incorporating dynamical tides and the evolution of stellar spin and structure, is needed to better reproduce the data (see, e.g., Bolmont & Mathis 2016, and references therein).

Zoom In Zoom Out Reset image size

Mulders et al. (2015) point out how the sub-Neptune occurrence rate breaks at ${P}_{\mathrm{break}}\sim 10$ days regardless of whether the host star is of spectral type F, G, K, or M. In our theory, the break and dropoff at shorter periods reflect the distribution of stellar spin periods ${P}_{\star }$. We can reproduce the invariance of ${P}_{\mathrm{break}}$ by drawing ${P}_{\star }$ from a distribution that does not vary with stellar mass. The Orion Nebula Cluster (ONC) provides such a spin distribution, which peaks at ∼10 days irrespective of whether the stars are ∼1 M_⊙ or ∼0.3 M_⊙ (compare the top and bottom panels of Figure 2). Drawing ${P}_{\star }$ from the ONC implies that planets emerge from disks that dissipate over the ONC age of ∼1 Myr. One problem with our ONC-based models is that they underestimate the number of USPs, particularly around FGK stars (left column of Figure 6); not enough stars rotate with periods faster than a day. USPs are produced more easily around older, more rapidly rotating stars (middle and right columns of Figure 6), but these have ${P}_{\star }$ distributions that do vary with stellar mass (Figure 2) and therefore have trouble reproducing the invariance of ${P}_{\mathrm{break}}$ (note how the model break period in the middle and the right columns of Figure 8 misses the observed break period by a factor of 5). Simultaneously fitting both USPs and their longer-period counterparts near the break is a challenge, requiring improvements in our understanding of disk truncation and tides.

There is also room for improvement for observations: the USP occurrence rates from Sanchis-Ojeda et al. (2014) appear systematically offset from the Fressin et al. (2013) data at longer periods; see how the blue points in Figure 9 appear to be shifted higher than the black points. Whether this offset is real or an artifact of combining two different data sets is unclear; removal of the offset would ameliorate if not solve the problem noted above with our ONC-based models. A uniform analysis is needed to measure planet occurrence rates from ∼0.1 to >10 days, not just for FGK stars but also for M stars.

Detecting and characterizing planets around main-sequence A stars is the next frontier. The invariance of ${P}_{\mathrm{break}}$ with respect to stellar spectral type may not extend to A stars. Compared to T Tauri stars, Herbig Ae stars are observed to be faster rotators, with typical spin periods of ∼1 day (Hubrig et al. 2009; Alecian et al. 2013). It follows that their gas disks should truncate at shorter periods.⁸ We therefore predict the occurrence rate of sub-Neptunes to break at a correspondingly short period ${P}_{\mathrm{break}}\sim 1$ day around A stars, as compared to ${P}_{\mathrm{break}}\sim 10$ days for FGKM stars. This shift in ${P}_{\mathrm{break}}$ between high-mass and low-mass stars is illustrated in Figure 11. The break around A stars might occur at even shorter sub-day periods, as our spin period data for A stars are taken from stellar magnetic field studies whose samples are biased to include more slowly rotating stars (Hubrig et al. 2009; Alecian et al. 2013). Figure 11 also shows that tidal inspiral does not much alter the period break for A stars for ${Q}_{\star }^{\prime }={10}^{7}\mbox{--}{10}^{9}$. These larger values of ${Q}_{\star }^{\prime }$ are motivated by main-sequence A stars lacking outer convective envelopes (see, e.g., Mathis 2015, their Figure 4). A caveat to our prediction that sub-Neptunes around A stars exhibit a shorter ${P}_{\mathrm{break}}$ is that we have so far discussed only disks of gas, not disks of solids—and sub-Neptunes are primarily composed of solids. Dust disks around A stars may not extend down to periods of ∼1 day because of sublimation (e.g., Dullemond & Monnier 2010, their Figure 7). Nevertheless, larger solids (planetesimals orders of magnitude larger than dust grains) are more resistant to vaporization and can drift, aerodynamically or otherwise, through gas disks toward their inner edges, assembling into planets there.

Zoom In Zoom Out Reset image size