A test of the high-eccentricity migration scenario for close-in planets
Steven Giacalone Titos Matsakos Arieh Königl
 


1Department of Physics, The University of Chicago, Chicago, IL 60637, USA
2Department of Astronomy & Astrophysics and The Enrico Fermi Institute,
The University of Chicago, Chicago, IL 60637, USA
Giacalone, Matsakos, & Königl
A test of the high-eccentricity migration scenario

In the high-eccentricity migration (HEM) scenario, close-in planets reach the vicinity of the central star on high-eccentricity orbits that become circularized—with a concomitant decrease in the semimajor axis—through a tidal interaction with the star. Giant planets that arrive with periastron distances that are smaller than the Roche limit a_R lose their gaseous envelopes, resulting in an inner edge to the surviving planets' distribution. The observational evidence for this effect, while extensive, is nonetheless somewhat ambiguous because of the effect of tidal orbital decay. Here we consider another key prediction of the HEM scenario—the existence of a spatial eccentricity gradient near the location where the circularization time becomes comparable to the planet's age for typical parameters. Previous studies already found evidence for this gradient and demonstrated that its properties are consistent with the circularization process being dominated by tidal dissipation in the planet (encapsulated by the tidal quality factor Q^'_p). Our work extends these treatments by constructing explicit model distributions for comparison with the data and by carrying out backward-in-time integrations using observed system parameters. We show that circularization generally occurs outside the distribution's inner edge (which defines the boundary of the so-called sub-Jovian desert) and that typically Q^'_p≈10^6 in the circularization zone (to within a factor of 3). We also find tentative evidence for an eccentricity gradient in lower-mass planets, indicating that formation through HEM may be relevant down to Neptune scales.


A test of the high-eccentricity migration scenario for close-in planets
    Steven Giacalone1, Titos Matsakos2, and Arieh Königl2
    December 30, 2023
=======================================================================



§ INTRODUCTION


The growing number of observed close-in exoplanets (planets with orbital periods P_orb≲ 10 days) has motivated researchers to look for trends in the distribution of their physical and orbital parameters that might help clarify the origin of these planets and the nature of their interaction with the host star. In an early study of this type, <cit.> drew attention to two such trends for giant planets: the prevalence of circular orbits for very short periods, and the inverse correlation between the planet's mass M_p and P_orb for the closest planets. Under the prevailing view, wherein giant planets form beyond the water–ice line, the first of these trends has two distinct potential interpretations: in the disk migration scenario, the observed planets reach the central star by drifting inward through the protoplanetary disk on nearly circular orbits <cit.>, whereas in the high-eccentricity migration (HEM) picture they arrive on high-eccentricity orbits that become tidally circularized when the planets approach the star <cit.>. In this connection, PHMF11 identified an unambiguous transition from eccentric to circular orbits on going from long to short orbital periods (for a given value of M_p) and from high-mass to low-mass  planets (for a given value of P_orb), and pointed out that this behavior is consistent with the expected outcome of a circularization process that is dominated by tidal dissipation in the planet <cit.>. These results were confirmed in a recent study by <cit.>.

In trying to interpret the second identified trend—the pile-up of the shortest-period planets in such a way that those with higher values of M_p have lower values of P_orb—PHMF11 speculated that tidal circularization and the stopping mechanism of close-in planets might be related. However, this mass–period relation was subsequently recognized to be part of a more general feature in the P_orb–M_p plane: a nearly empty area, outlined roughly by two oppositely sloped lines, in the region of sub-Jupiter-mass planets on short-period orbits <cit.>. <cit.> showed that this feature (dubbed the sub-Jovian desert) can be interpreted in terms of HEM, with the two distinct segments of the desert's boundary reflecting the different slopes of the empirical mass–radius relation for small and large planets <cit.>.  A plausible physical origin for the boundary is the Roche limit a_R, the distance from the star where the planet starts to be tidally disrupted. If the planet arrives on an orbit with initial (subscript 0) semimajor axis a_0 and eccentricity e_0, its distance of closest approach will be a_per,0=(1-e_0) a_0 and it will circularize (assuming conservation of orbital angular momentum) at a_cir=(1+e_0) a_per,0 (i.e., at ≃ 2 a_per,0 for a highly eccentric orbit; ). The upper boundary of the sub-Jovian desert was already interpreted in this way by <cit.>, although the data available at the time was insufficient for a definitive model fit.[As was pointed out by MK16, the observed shape of the desert's upper boundary is not adequately reproduced unless one also takes into account the tidal dissipation in the star, which, on a timescale much longer than the circularization time, causes the planet's orbit to decay. In a previous study, <cit.> interpreted the finding of giant planets with semimajor axes < 2 a_R in terms of orbital decay of this type.] The high-e_0 orbit could originate in a sudden planet–planet scattering event or in a slower interaction such as Kozai migration (involving either a stellar or a planetary companion) or secular chaos. In the latter case <cit.> suggested that the planet's inward drift might be arrested at the location where the rate at which its longitude of pericenter precesses due to a secular interaction with a more distant planet is equal to the orbit-averaged precession rate associated with the tidal quadrupole induced on the planet by the star. The locus of circularization radii in the P_orb–M_p plane is, however, similar in this case to the one obtained by setting a_per,0=a_R. Note that the planet's stopping mechanism is not directly related to the circularization process in either of these two explanations of the pile-up. However, in order for the data to be compatible with the prediction a_cir≈ 2  a_per,0, the circularization radius r_cir (obtained by equating the planet's circularization time to the time that has elapsed since its arrival at the stellar vicinity) must exceed a_cir. As we demonstrate in Section <ref>, this condition is typically satisfied for planets in the pile-up zone.

The observational support for the role of the HEM mechanism in shaping the spatial distribution of close-in planets has so far been based primarily on the apparent paucity of planets with semimajor axes a≲ 2 a_R <cit.> and on the corresponding dearth (the sub-Jovian desert) in the period–mass plane (; MK16). While the observational evidence for this effect is strong, it is not entirely unambiguous on account of the (already noted) additional orbital evolution induced by tidal dissipation in the star. In this paper we consider a complementary observational test of this scenario, the expected gradient in planet eccentricities in the vicinity of the locus of the circularization radii in the period–mass plane. The existence of such a gradient was already demonstrated in PHMF11, but here, in addition to updating the database, we compare it explicitly with the predictions of the HEM model. We describe our modeling approach in Section <ref>, present results in Section <ref>, discuss the main implications in Section <ref>, and summarize in Section <ref>. In a separate paper <cit.> we use the model employed in this work to study the fate of high-mass planets that arrive by HEM and end up crossing the Roche limit, which results in the loss of their gaseous envelopes: we argue that the remnant rocky cores of these planets can plausibly account for the recently identified population of dynamically isolated hot Earths <cit.>.


§ MODELING APPROACH


Our treatment is based on the formulation presented in MK16, whose work was similarly concerned with the properties of close-in planets that arrive by HEM and undergo orbital circularization through internal tidal dissipation. That paper examined the shape of the boundary of the sub-Jovian desert in the P_orb–M_p plane under the assumption that planets reach the vicinity of the Roche limit a_R (or, alternatively, the point of closest approach in the secular-chaos model of ) with e_0 ≈ 1 and that their orbits then undergo an effectively instantaneous circularization. In contrast with that work, in which only the post-circularization orbital evolution of planets due to tidal dissipation in the star was calculated, in this paper we also account explicitly for tidal dissipation in the planets and we consider its effect on orbital circularization for a range of initial eccentricities and without assuming a priori that this process always runs to completion. We carry out Monte Carlo simulations using the same distributions of P_orb,0, R_p, M_p, and t_arr (the planet arrival time at the stellar vicinity) as in MK16, but we update their choices for the distributions of radii and masses of large planets as well as of  t_age (the system's age) using data downloaded from the Extrasolar Planets Encyclopedia database at exoplanet.eu.

Any given system is specified by six parameters: P_orb,0, e_0, R_p, M_p, t_age, and t_arr. We assume that the planets that reach the stellar vicinity by HEM originate in the P_orb,0 range [10,100] days with e_0 in the range [0.5,0.9]. Though consistent with our current understanding of how the HEM process operates in real systems <cit.>, these ranges are not meant to represent any particular physical model and are chosen for illustrative purposes only. For P_orb,0 we adopt the empirical distribution ∂f/∂logP_orb,0∝ P_orb,0^0.47 given in <cit.>, whereas for the e_0 distribution we adopt the form ∂f/∂ e_0 = constant (which corresponds to the steady-state distribution obtained by  and  for planets that undergo eccentricity oscillations due to secular gravitational interactions with an outer companion, as well as to the “coplanar” distribution obtained in the secular-chaos model and shown in figure 4 of ).[The results presented in this paper are not sensitive to the details of the P_orb,0 and e_0 distributions. However, in KGM17 we examine the dependence of the fraction of planets that end up crossing the Roche limit on the form of the e_0 distribution and on the maximum value of e_0.] The values of R_p are also sampled from an empirical distribution given in <cit.>, ∂f/∂logR_p∝ R_p^-0.66. As in MK16, we distinguish between small and large planets, separated at R_p = 12  R_. For the small planets we adopt 

    M_p = (R_p/R_)^2 M_  ,     R_p<12 R_

<cit.>. In the case of the large planets—for which R_p is nearly independent of M_p—we resample R_p from the interval [9,20] R_ and independently sample M_p from the interval [0.3,10] M_J using the exoplanet.eu database (see Appendix <ref>; note that this was done in MK16 using the data originally compiled by ).[The radii of giant planets evidently depend on the incident flux from the host star <cit.> and could thus vary systematically with orbital period. We do not explicitly account for this dependence since it was found to be fairly weak <cit.> and because the orbital circularization zone that we investigate corresponds to a rather narrow range of P_orb values.] In another modification of the MK16 implementation, we replace the age distribution given in <cit.>—which was determined by applying gyrochronology relationships to comparatively rapidly rotating stars and is therefore biased toward young systems—by an empirical distribution obtained from exoplanet.eu. We restrict attention to the P_orb interval [2.5,7] days (where the upper bound defines the regime of hot Jupiters and the lower limit roughly corresponds to the distance from the star below which the effect of tidal orbital decay becomes significant) and consider separately the age distributions for small and large planets (see Appendix <ref>). Finally, we assume a uniform distribution in logt_arr for the arrival times, with t_arr∈[0.01,10] Gyr (see MK16).

After sampling for P_orb,0 and e_0, one can determine a_per,0(a_0,e_0) and a_cir(a_0,e_0) (setting a_0 =(GM_*P_orb,0^2/4π^2)^1/3, where G is the gravitational constant and M_* is the stellar mass). These values can, in turn, be used to select the systems that are relevant to the present calculation. In this work we only consider the Roche-limit interpretation of the sub-Jovian desert boundary, taking the lower bound on a planet's initial periastron distance to be given by 

    a_R   =   q(M_*/M_p)^1/3 R_p
       =    0.016 (q/3.46) (M_*/M_)^1/3(M_p/M_J)^-1/3(R_p/R_J) au

(where the normalization of the coefficient q is based on the results of MK16). The systems whose evolution we follow are defined by the requirements a_per,0>a_R and P_orb(a_cir)≲ 7 days. We simplify the evolution equations by assuming that the orbital angular momentum vector is aligned with the spin vector of the star as well as with that of the planet, and that the planet's rotation period does not change with time (corresponding to pseudosynchronicity). We also neglect the time variation of the stellar rotation period P_* (which we assume to be distributed uniformly in the interval [5,10] days). The evolution equations are then given by

    da/dt=9/Q^'_p (GM_*/a^3 )^1/2M_*/M_pR_p^5/a^4(1-e^2)^-15/2
    ×[(f_2(e^2))^2/f_5(e^2)-f_1(e^2) ]
    
    + 9/Q^'_* (GM_*/a^3 )^1/2M_p/M_*R_*^5/a^4(1-e^2)^-15/2
    ×[f_2(e^2)P_orb/P_*(1-e^2)^3/2-f_1(e^2) ]

and

    de/dt=81/2 Q^'_p (GM_*/a^3 )^1/2M_*/M_pR_p^5/a^5e(1-e^2)^-13/2
    ×[11/18f_4(e^2)f_2(e^2)/f_5(e^2)-f_3(e^2)]
    
    + 81/2 Q^'_* (GM_*/a^3 )^1/2M_p/M_*R_*^5/a^5e(1-e^2)^-13/2
    ×[11/18f_4(e^2)P_orb/P_*(1-e^2)^3/2-f_3(e^2) ]

<cit.>, where the eccentricity functions f_1, . . ., f_5 (each of which equals 1 at e=0) are given in <cit.>, Q^'_p and Q^'_* are, respectively, the (modified) planetary and stellar tidal quality factors, and R_* is the stellar radius. We express both Q^'_p and Q^'_*  in the form Q^'=Q^'_1(P_orb/P_1), which was employed in previous studies as a representation of equilibrium tides in the weak-friction approximation <cit.>.[Q^'_* is probably better modeled in terms of a dynamical tide, but such models also infer a power-law dependence on P_orb with a positive (albeit >1) index in both the weakly and the strongly nonlinear regimes <cit.>. In this paper we are primarily interested in the behavior of Q^'_p, which underlies the orbital circularization process.] We set P_1=4 days and adopt Q^'_*1=10^6.[Note in this connection that MK16 treated Q^'_*  as a spatial constant equal to 10^6.] We use the results presented in Section <ref> to constrain the value of  Q^'_p1.

For the values of e_0 that we consider (1-e_0^2) is not ≪ 1 and one can define a characteristic orbital circularization time by τ_cir≡ |(1/e)de/dt|^-1. We estimate τ_cir  from the first term on the right-hand side of Equation (<ref>) by taking the limit e→ 0 and identifying a with r_cir. By equating τ_cir to t_age, we obtain an expression for the locus of the circularization radii of planets with the given age in the period–mass plane:[More precisely, one should consider the locus of the circularization radii of planets with a given (nonnegative) value of (t_age - t_arr). However, for the chosen distributions of ages and arrival times, the distribution of the nonnegative values of (t_age - t_arr) closely approximates that of t_age.] 

    P_orb,cir = 3.76 (P_1/4 days)^3/16(Q'_p1/10^6)^-3/16( t_age/1  Gyr)^3/16
    ×( M_*/M_)^-1/8( R_p/R_J)^15/16( M_p/M_J)^-3/16 days .

We express R_p in Equation (<ref>) as a function of M_p by using the relationship given in Equation (<ref>) for M_p<150 M_ and by approximating its behavior for more massive planets by R_p=constant. This implies that the P_orb,cir curve in the period–mass plane changes from having a positive slope (∝ M_p^9/32) for M_p<150 M_ to having a negative slope (∝ M_p^-3/16) for larger masses.[If Q^'_p were instead a spatial constant, as is sometimes assumed, then the listed scalings of P_orb,cir would change to ∝ M_p^9/26 and ∝ M_p^-3/13 for (respectively) small and large masses.] As was pointed out by MK16, an analogous behavior is found for the immediate-post-circularization boundary of the sub-Jovian desert P_orb,RL (identified as the orbital period that corresponds to a_cir(a_R)=(1+e_0)a_R) by using the above functional form of R_p(M_p) in Equation (<ref>): P_orb,RL∝ M_p^1/4 and ∝ M_p^-1/2 for M_p<150 M_ and >150 M_, respectively. 

One can similarly obtain the orbital decay isochrones by considering the dominant (second) term on the right-hand side of Equation (<ref>). For the assumed dependence of Q^'_* on P_orb, that term is ∝ a^-7, so we define the orbital decay time as τ_d≡ |(8/a)da/dt|^-1 <cit.>. This yields

    P_orb,d = 3.18 (P_1/4 days)^3/16(Q'_*1/10^6)^-3/16( t_age/1  Gyr)^3/16
    ×( M_*/M_)^-1/2( R_*/R_)^15/16( M_p/M_J)^3/16[P_orb/P_*-1]^3/16 days .

In all numerical evaluations of this expression we assume, for simplicity, that [(P_orb/P_*)-1]^3/16≈ 1.


§ RESULTS


Our Monte Carlo simulations each involve 30,000 samplings of planetary systems with a solar-type host (R_*=R_, M_*=M_). As the R_p distribution that we employ was corrected for observational selection effects and is dominated by small planets, we randomly reduce the number of small (R_p<12 R_) model planets that we exhibit by 90% to improve the presentation. The top panel of Figure <ref> shows the calculated planet distribution in the P_orb–M_p plane, color coded according to the value of the orbital eccentricity at the end of the modeled evolution. For this panel we adopt Q^'_p1=10^6. We also plot the circularization isochrones (Equation (<ref>)) for two values of t_age—1 and 5 Gyr (solid and dashed blue lines, respectively)—which mark off the range from which most of the system ages are drawn (see Figure <ref>). It is seen that these curves capture well the numerical results in that the region between them corresponds to the transition zone in the period–mass plane that separates mostly eccentric orbits (to the right of the dashed curve) from mostly circular orbits (to the left of the solid curve). In addition, we plot the immediate-post-circularization desert boundary curves (P_orb,RL(M_p)) for the two values of e_0 (0.5 and 0.9) that bracket our adopted range of initial eccentricities. These lines are seen to lie to the left of the circularization isochrones and well within the region of mostly circularized orbits, corroborating the assumption that the planets near the desert boundary satisfy r_cir>(1+e_0)a_R (see Section <ref>). Finally, we plot the orbital decay isochrones (Equation (<ref>)) for the same two values of t_age (solid and dashed red lines, respectively). These lines pass near the vertices of the P_orb,RL(M_p) curves, indicating that orbital decay is likely to affect the shape of the upper desert boundary (see MK16). However, the orbital decay isochrones intersect the circularization isochrones at a sufficiently large value of M_p (≃ 2 M_J) to insure that most of the modeled planets are not measurably affected by orbital decay during the circularization process.[It should, however, be possible for the orbits of sufficiently massive planets to decay before they are fully circularized, which is consistent with the comparatively high inferred frequency of eccentric orbits among transiting planets with M_p>3 M_J <cit.>.]

The top panel of Figure <ref> also exhibits observational data points from the exoplanet.eu compilation. They are shown as either triangles or stars and are color coded in the same way as the model dots. We only include planets with measured eccentricities that are listed with error values (generally both an upper and a lower one). For a data point to be considered reliable, we require that both of the associated values of δ e be < 0.5 max(0.1,e);[The form of this criterion is motivated by the separation of eccentric orbits in PHMF11, <cit.>, and <cit.> into those with e<0.1 and those with e≥ 0.1, and by the 1σ uncertainty limit δ e < 0.05 that <cit.> adopted as a reliability criterion for circular (e=0) orbits.] we represent such a data point by a triangle. If either one of the values of δ e does not satisfy the above inequality, we consider the associated data point to be questionable and display it as a star. These data points can be used to check a key prediction of the “HEM + circularization" scenario—the presence of an eccentricity gradient in the vicinity of the plotted circularization isochrones. Although the number of reliable data points in this region of the period–mass plane is relatively small, the predicted gradient is uniquely specified to point along the normal to these sloping lines, which should facilitate the test. A visual inspection of the top panel does indeed indicate consistency with this prediction, not just for the upper portions of the isochrone curves that were considered in PHMF11 but possibly also for the differently oriented lower branches of these curves. The additional data accumulated since the PHMF11 work was carried out also make it possible to resolve the gradient on smaller scales in the period–mass plane and therefore to localize it better in relation to the isochrone curves.


The location of the circularization radius depends on the magnitude of the planetary tidal quality factor: it shifts to lower values of P_orb as Q^'_p is increased. In an attempt to constrain the value of Q^'_p1, we consider the cases where it is changed to 3× 10^5 and 3× 10^6 (left and right panels, respectively, at the bottom of Figure <ref>). To simplify this exercise, we use the two selected circularization isochrones to demarcate the circularization zone in the period–mass plane, and we only display data points that have low associated errors; in addition, we only consider data points that lie to the right of the 1 Gyr P_orb,d curve to minimize the effect of orbital decay. It is seen that the low-e data points are concentrated too far to the left of the circularization isochrones for Q^'_p1=3× 10^5 and not far enough to the left for Q^'_p1=3× 10^6, pointing to Q^'_p1≈ 10^6 as the preferred value.


To check the extent to which our inferences from comparing model calculations with observational data depend on the error tolerance criterion used in selecting the data points, we modify the coefficient α in the condition δ e < α max(0.1,e) (where α = 1/2 corresponds to the fiducial case shown in Figure <ref>). The basic results from the top panel of Figure <ref> are shown in the top right panel of Figure <ref>, where we retain the various model curves (the circularization and orbital decay isochrones as well as the immediate-post-circularization boundaries of the sub-Jovian desert) and the reliable data points (triangles) but do not reproduce the simulation results (dots) and the high-δ e data points (stars). For comparison, we show the corresponding results using α = 1 and 1/3 in the top left and bottom left panels, respectively, of Figure <ref>. It is seen that, while the number of reliable data points decreases as α is decreased, the qualitative behavior—and in particular the appearance of a spatial eccentricity gradient in the vicinity of both the upper and the lower branches  of the circularization isochrones—is unchanged. We also confirmed that Q^'_p1≈ 10^6 remains the preferred value when either the more stringent selection criterion (α = 1/3) or the looser one (α =1) is used.

Figure <ref> presents the same results as Figure <ref> but in the period–eccentricity plane, with the planetary mass now being the color-coded variable. The three panels correspond to the same values of Q^'_p1 as in Figure <ref>. To bring out the effect of tidal dissipation in the planet, we only display model and data points that lie to the right of the 1 Gyr P_orb,d curve in Figure <ref>. The model planets shown in the top panel exhibit a clear transition from being dominated by comparatively high eccentricities for P_orb 6 days to acquiring low values of e closer to the star—the signature of the “HEM + circularization” process. The observed systems appear to have a similar distribution and thus to be compatible with this scenario. The model dots track the observational data points best in the top panel: they appear to lie too far to the right in the bottom left panel and too far to the left in the bottom right panel, reconfirming the choice of Q^'_p1≈ 10^6 as the best-fitting value.


A few of the data points displayed in Figure <ref> have very low (0.01) values of e, and yet their eccentricities are all distinct from zero. This raises two questions: (1) why are there no genuinely circular orbits in the dataset that we exhibit, and (2) could the very-low-e data points actually correspond to circular orbits? The answer to the first question is that most of the e=0 entries in the exoplanet.eu database are not listed with error values and we therefore do not consider them. (Several e=0 points that are listed with errors are included in the top panel of Figure <ref>; however, they all lie to the left of the 1 Gyr P_orb,d curve in that panel and have therefore been filtered out of Figure <ref>.) It is, however, evident from an inspection of Figure <ref> that this omission has little effect on the properties of the spatial eccentricity gradient—the focus of this work—since those are defined by data points with higher values of e. The second issue is related to a well-known intrinsic bias in the determination of eccentricity from radial velocity data <cit.>. This bias is a consequence of the fact that the value of e cannot be negative, which implies that observational uncertainties tend to yield positive (even if small) values of e for genuinely circular orbits. Mindful of this fact, PHMF11 and <cit.> carried out a homogeneous Bayesian analysis that explicitly addressed this issue. <cit.> extended these results by using a larger (by a factor of 3) sample of transiting planets with improved eccentricity and mass determinations. The latter authors considered planets in the mass range (0.1,25) M_J and classified their orbits as being circular (e=0 with 1σ uncertainty < 0.05), eccentric, or unconstrained (having e compatible with zero but δ e >  0.05 or else a slightly eccentric orbit that is not strongly supported by the Bayesian model). To check on the effect of this bias on our conclusions, we repeat the test presented in Figure <ref> using the data from this sample (and again including planet masses only up to 10 M_J). The result, shown in the bottom right panel of Figure <ref>, demonstrates that a spatial eccentricity gradient can be discerned also in this case in the vicinity of the upper branches of the circularization isochrones in the P_orb–M_p plane. This dataset does not, however, contain reliable eccentric orbits for masses that lie below the break in the circularization isochrones and thus provides no information on the possible presence of an eccentricity gradient also in lower-mass planets.

As a final check, we carry out backward-in-time integrations for observed systems that possess reliable eccentricity determinations. To isolate the effect of orbital circularization, we only consider planets for which exoplanet.eu lists an age estimate that satisfies t_age<τ_d (or, equivalently, P_orb>P_orb,d, where τ_d is evaluated using the listed values of R_*, M_*, M_p, and a), so that orbital decay is not important. Thus we only retain the planet-dissipation terms in Equations (<ref>) and (<ref>). Each system is integrated backward from its current location over a time interval equal to its age (see Footnote <ref>), with Q^'_p1 set equal to 10^6. Figure <ref> shows the calculated evolutionary tracks in the P_orb–e plane, with the three panels corresponding to the three error tolerance limits (specified by the parameter α) that were employed in Figure <ref>. The results indicate that a significant fraction of these planets have undergone orbital circularization and had initial orbital periods that lay outside the close-in range. Although the total number of systems that could be tested in this way is not large, the fact that this outcome was not inevitable—as evidenced by the presence among the systems considered in Figure <ref> of planets that exhibit little change in P_orb over their lifetimes—strengthens the conclusion that HEM is indeed relevant to the origin of many close-in planets.[<cit.> employed similar backward-in-time integrations in an attempt to estimate the values of Q^'_p and Q^'_* by matching the implied e_0 distribution to the observed eccentricity distribution for a>0.2 au. Although their derived best-fit values,  Q^'_p∼ 3× 10^6 and Q^'_*∼ 3× 10^5, are close to those obtained by other methods, this approach is subject to a number of caveats. For example, the calculated values of e_0 depend on the durations of the backward integrations, so their inferred distribution is affected by the (sometimes considerable) uncertainty in t_age. Furthermore, the general eccentricity distribution at a>0.2 au may not be representative of the initial distribution for planets that are transported by HEM to the center. Other uncertainties associated with this approach were noted by <cit.>. Such integrations are, however, useful for demonstrating that the data are consistent with the “HEM + circularization” scenario for a broad choice of values for Q^'_p and Q^'_* <cit.>.] It is noteworthy that a fraction of the systems that evolve back to P_orb> 10days in each of the panels (6 out of 14, 5 out of 11, and 2 out of 8 for α = 1, 1/2, and 1/3, respectively) have M_p<150 M_ and thus correspond to planets that lie below the vertices of the circularization isochrone curves in the P_orb–M_p plane. This supports our tentative inference from the results presented in Figure <ref> (and in the α = 1, 1/2, and 1/3 panels of Figure <ref>) that HEM may be implicated in the arrival of close-in planets that span a broad range of masses (from Jupiter to Neptune scales).



§ DISCUSSION


PHMF11 were the first to draw attention to the existence of an eccentricity gradient in the period–mass plane and to consider its implications. They pointed out that its orientation agreed with that of a circularization isochrone based on tidal dissipation in the planet, for which τ_cir∝ (a/R_p)^5M_p/M_*. (By contrast, τ_cir∝ (a/R_*)^5M_*/M_p if dissipation in the star is dominant.)  They further demonstrated that the location of the transition from e>0.1 to e<0.1 roughly corresponds to a 1 Gyr isochrone characterized by Q^'_p=10^6 <cit.>. These findings were confirmed in the more extensive recent study by <cit.>.  Our work extends these results by confronting the data with explicit predictions of the HEM model. This approach has made it possible to constrain the value of the planetary tidal quality factor in the circularization zone (Q^'_p≈ 10^6 for P_orb 4 days). Although a similar value is often adopted in the literature, based on theoretical calculations and solar-system observations <cit.>, the constraints derived directly from the exoplanet data have been much less restrictive. Specifically, by requiring that the orbits of planets with e≈ 0 circularize on timescales shorter than their ages and that those with finite eccentricities do not, one can obtain upper and lower limits, respectively, on Q^'_p. The range of values inferred in this way spans several orders of magnitude (∼ 10^5-10^9; e.g., ). By comparison, we were able to estimate a preferred characteristic value for Q^'_p in the circularization region to within a factor of three.

Although the database that we employed is not as reliable as the one assembled, for example, by <cit.> using a homogeneous statistical analysis, we have considered lower masses (down to M_p≈ 0.03 M_J) than in previous studies of this topic. This has led us to a tentative conclusion that the HEM process may also play a role in the formation of sub-Jovian-mass planets (down to Neptune size). The shape of the circularization isochrones in the period–mass plane resembles that of a bird's beak, and the presence of an eccentricity gradient for planets with M_p<150 M_ would establish the reality of the lower portion of that beak. The potential added significance of such a determination is that the HEM scenario implies the existence in the P_orb–M_p plane of a similar “bird's beak” structure at lower values of P_orb—the boundary of the sub-Jovian desert (see top panel in figure <ref> and MK16). Given that the shape and origin of the lower boundary of the desert are still being debated <cit.>, the presence of an eccentricity gradient in association with the lower branches of the circularization isochrones would support the HEM interpretation of the sub-Jovian desert even as it broadens the range of planetary masses in which this mechanism is found to operate. Clearly, more data are needed to validate the existence of this gradient: future space missions such as TESS <cit.> and PLATO <cit.> hold promise in this regard.

Despite its apparent success in explaining a variety of observational findings, the extent of the contribution of the HEM mechanism to the formation of close-in giant planets is still being investigated. In one test of this scenario, <cit.> looked for observational evidence for highly eccentric Jupiter-mass planets, which were predicted to exist if hot Jupiters originate outside the ice line and their HEM is induced by a distant (stellar) companion <cit.>. <cit.> did not find such evidence, but they suggested that this does not rule out the possibility that hot Jupiters originate interior to the ice line and that their orbits are perturbed by a planetary companion. In a subsequent study, <cit.> inferred that the probability of a giant planet having a Jupiter-mass companion capable of inducing HEM does not depend on whether the planet lies within the hot-Jupiter orbital range or on the location of the companion with respect to the ice line. These results do not support HEM models in which close-in giant planets originate beyond the ice line, but they are again compatible with the possibility (which  also recognized) that  such planets can travel at least part of the distance from their formation sites by means other than HEM (e.g., classical disk migration or a secular dynamical interaction with a companion).[As was noted in Section <ref>, the choice of initial conditions for our model is consistent with this emerging understanding of how the HEM mechanism operates in real systems.] Under these circumstances, and given that not all systems that harbor a close-in giant planet show evidence for an outer planetary companion <cit.>, it is natural to expect that some fraction of the observed giant planets—including hot Jupiters—have not experienced HEM. There have already been attempts to quantify this fraction based on the difference in the orbital characteristics of the two planet arrival modes, and it appears that it could be appreciable <cit.>. It is also worth keeping in mind that a large number of giant planets likely form in the protoplanetary disk and reach the star through classical disk migration <cit.>. Some of these planets may have been stranded near the host star for up to ∼ 1 Gyr before being tidally ingested and—notwithstanding the fact that their contribution to the observed number count of planets is small—could have left a lasting imprint on the obliquity and metallicity properties of their hosts (; KGM17).


§ CONCLUSION


We tested a key prediction of the HEM scenario for the origin of hot Jupiters and other close-in planets. In this interpretation, planets arrive in the vicinity of the host star on high-eccentricity orbits that become circularized through tidal interaction with the star if they reach orbital periods that are less than P_orb,cir (Equation (<ref>)). This picture implies that a spatial eccentricity gradient should be present in the period–mass phase space near the locus of circularization radii that correspond to typical system ages. The existence of such a gradient for close-in giant planets had been first pointed out by PHMF11 <cit.>, and this finding was recently confirmed by <cit.>. Our treatment is distinct from previous work in that it explicitly tests the HEM scenario by comparing the model predictions (obtained by integrating the evolution equations using observationally or theoretically constrained initial conditions) with the data. This approach has enabled us to extract valuable information about the circularization process. It was already deduced by PHMF11 (and, through alternative methods, by other workers) that this process is dominated by tidal dissipation in the planet. We verified that the effect of orbital decay (dominated by tidal dissipation in the star) for planets of mass M_p 2 M_J is not important in the circularization zone (which roughly spans the P_orb range ∼ 4–6 days) and inferred that the characteristic value of the planetary tidal quality factor Q^'_p in this region is ∼ 10^6. We reached this conclusion through a qualitative comparison between the model results and the data in the period–mass and mass–eccentricity planes. Although we checked that these results are not sensitive to the details of the error tolerance criteria used in selecting the data points, we did not perform a formal statistical test since the number of systems that can be used to demonstrate the existence of the gradient is not yet large enough to justify carrying out such an analysis over the three-dimensional (P_orb, M_p, e) parameter space. However, we have found that our procedure is reliable enough to pin down the value of Q^'_p in the circularization zone to within a factor of three, which can be compared with the ∼ 10^5-10^9 range obtained previously through an application of generic constraints. We stress, however, that the inferred value of Q^'_p pertains only to the narrow range of orbital periods where the transition from mostly eccentric to mostly circular orbits occurs in the model; thus, we cannot reliably test the possible spatial variation of Q^'_p. Furthermore, this parameter only provides a very basic description of tidal dissipation, and our procedure does not test realistic models of this process.

Planets that reach the Roche limit a_R (Equation (<ref>)) become tidally disrupted: this limit thus provides a natural edge to the observed distribution (corresponding to the boundary of the sub-Jovian desert). We demonstrated that this edge generally lies interior to the predicted location of the circularization zone, consistent with the data. We noted the possible observational indication of an eccentricity gradient for sub-Jovian-mass planets: if confirmed by additional data, such a gradient would attest to the relevance of the HEM mechanism also to non-giant close-in planets and would support the interpretation of the lower boundary of the desert  in terms of this mechanism.

Giant planets that cross the Roche limit—either on their initial high-eccentricity trajectories or at a later stage (after their orbits are circularized) due to orbital decay—can be expected to lose their gaseous envelopes and be converted into remnant cores. In KGM17 we extend the calculations presented in this paper by continuing to follow the evolution of these cores, and argue that such remnants are natural candidates for dynamically isolated hot Earths. If this interpretation is correct, it will provide an additional—and independent—argument in favor of the “HEM + circularization” scenario.
We acknowledge fruitful discussions with Dan Fabrycky and thank him and the referee for helpful suggestions. This work was supported in part by NASA ATP grant NNX13AH56G and by a University of Chicago College Research Fellows Fund award to S.G.


§ EMPIRICAL DISTRIBUTIONS FROM THE EXTRASOLAR PLANETS ENCYCLOPEDIA


As in MK16, we follow the approach of <cit.> and divide the planet population into two sets, “small” and “large,” with R_p=12 R_ and M_p=150 M_ serving as the rough dividing values of radius and mass, respectively. For the smaller planets we adopt the M_p(R_p) relation given by Equation (<ref>), whereas for the larger ones we employ the empirical distributions shown in Figure <ref> (which are sampled independently). It is seen that the variation in the value of R_p for this set is much smaller than that in M_p, justifying the adoption of the ansatz R_p = constant for the analytic approximations in Section <ref>.

We maintained the separation into “small” and “large” planets in constructing the age distribution. The left panel of Figure <ref> shows the results for the larger planets. It is seen that the bulk of the systems have t_age in the range ∼1–5  Gyr: this is considerably broader than the distribution used in MK16, which was biased toward younger systems. We employed the values that delineate this range (1 and 5 Gyr) in plotting the isochrone curves in Figures <ref> and <ref>. The age distribution for the smaller planets is shown in the right panel of Figure <ref>. Although the number of systems returned by the search in this case was small, there is a strong indication that the observed distribution is qualitatively different from the one in the left panel, justifying the separate catalog queries that we made.

The above empirical distributions were obtained without filtering on the basis of the listed errors. While these are typically not large for the planets' radii and masses, they can be significant for the systems' ages, with a factor of 2 uncertainty in the value of t_age being fairly common.



natexlab#1#1[Antonini et al.(2016)Antonini, Hamers, &
  Lithwick]Antonini+16
Antonini, F., Hamers, A. S., & Lithwick, Y. 2016, , 152, 174
[Barker(2011)]Barker11
Barker, A. J. 2011, , 414, 1365
[Barker & Ogilvie(2009)]BarkerOgilvie09
Barker, A. J., & Ogilvie, G. I. 2009, , 395, 2268
[Barker & Ogilvie(2010)]BarkerOgilvie10
—. 2010, , 404, 1849
[Beaugé & Nesvorný(2013)]BeaugeNesvorny13
Beaugé, C., & Nesvorný, D. 2013, , 763, 12
[Bonomo et al.(2017)Bonomo, Desidera, Benatti, Borsa,
  Crespi, Damasso, Lanza, Sozzetti, Lodato, Marzari, Boccato,
  Claudi, Cosentino, Covino, Gratton, Maggio, Micela, Molinari,
  Pagano, Piotto, Poretti, Smareglia, Affer, Biazzo, Bignamini,
  Esposito, Giacobbe, Hébrard, Malavolta, Maldonado, Mancini,
  Martinez Fiorenzano, Masiero, Nascimbeni, Pedani, Rainer, &
  Scandariato]Bonomo+17
Bonomo, A. S., Desidera, S., Benatti, S., et al. 2017, , 602, A107
[Bryan et al.(2016)Bryan, Knutson, Howard, Ngo, Batygin,
  Crepp, Fulton, Hinkley, Isaacson, Johnson, Marcy, &
  Wright]Bryan+16
Bryan, M. L., Knutson, H. A., Howard, A. W., et al. 2016, , 821, 89
[Dawson et al.(2015)Dawson, Murray-Clay, &
  Johnson]Dawson+15
Dawson, R. I., Murray-Clay, R. A., & Johnson, J. A. 2015, , 798, 66
[Dobbs-Dixon et al.(2004)Dobbs-Dixon, Lin, &
  Mardling]Dobbs-Dixon+04
Dobbs-Dixon, I., Lin, D. N. C., & Mardling, R. A. 2004, , 610, 464
[Eggleton et al.(1998)Eggleton, Kiseleva, &
  Hut]Eggleton+98
Eggleton, P. P., Kiseleva, L. G., & Hut, P. 1998, , 499, 853
[Essick & Weinberg(2016)]EssickWeinberg16
Essick, R., & Weinberg, N. N. 2016, , 816, 18
[Fabrycky et al.(2007)Fabrycky, Johnson, &
  Goodman]Fabrycky+07
Fabrycky, D. C., Johnson, E. T., & Goodman, J. 2007, , 665, 754
[Ford & Rasio(2006)]FordRasio06
Ford, E. B., & Rasio, F. A. 2006, , 638, L45
[Husnoo et al.(2012)Husnoo, Pont, Mazeh, Fabrycky,
  Hébrard, Bouchy, & Shporer]Husnoo+12
Husnoo, N., Pont, F., Mazeh, T., et al. 2012, , 422, 3151
[Hut(1981)]Hut81
Hut, P. 1981, , 99, 126
[Jackson et al.(2008)Jackson, Greenberg, &
  Barnes]Jackson+08
Jackson, B., Greenberg, R., & Barnes, R. 2008, , 678, 1396
[Königl et al.(2017)Königl, Giacalone, &
  Matsakos]Konigl+17
Königl, A., Giacalone, S., & Matsakos, T. 2017, , in press
[Laughlin et al.(2011)Laughlin, Crismani, &
  Adams]Laughlin+11
Laughlin, G., Crismani, M., & Adams, F. C. 2011, , 729, L7
[Lin et al.(1996)Lin, Bodenheimer, & Richardson]Lin+96
Lin, D. N. C., Bodenheimer, P., & Richardson, D. C. 1996, Natur, 380,
  606
[Lucy & Sweeney(1971)]LucySweeney71
Lucy, L. B., & Sweeney, M. A. 1971, , 76, 544
[Matsakos & Königl(2015)]MatsakosKonigl15
Matsakos, T., & Königl, A. 2015, , 809, L20
[Matsakos & Königl(2016)]MatsakosKonigl16
—. 2016, , 820, L8
[Matsumura et al.(2010)Matsumura, Peale, &
  Rasio]Matsumura+10
Matsumura, S., Peale, S. J., & Rasio, F. A. 2010, , 725, 1995
[Matsumura et al.(2008)Matsumura, Takeda, &
  Rasio]Matsumura+08
Matsumura, S., Takeda, G., & Rasio, F. A. 2008, , 686, L29
[Mazeh et al.(2016)Mazeh, Holczer, & Faigler]Mazeh+16
Mazeh, T., Holczer, T., & Faigler, S. 2016, , 589, A75
[Nelson et al.(2017)Nelson, Ford, & Rasio]Nelson+17
Nelson, B. E., Ford, E. B., & Rasio, F. A. 2017, ArXiv e-prints,
  arXiv:1703.09711
[Petrovich & Tremaine(2016)]PetrovichTremaine16
Petrovich, C., & Tremaine, S. 2016, , 829, 132
[Pont et al.(2011)Pont, Husnoo, Mazeh, &
  Fabrycky]Pont+11
Pont, F., Husnoo, N., Mazeh, T., & Fabrycky, D. 2011, , 414,
  1278
[Rasio & Ford(1996)]RasioFord96
Rasio, F. A., & Ford, E. B. 1996, Sci, 274, 954
[Rauer et al.(2014)Rauer, Catala, Aerts, Appourchaux,
  Benz, Brandeker, Christensen-Dalsgaard, Deleuil, Gizon, Goupil,
  Güdel, Janot-Pacheco, Mas-Hesse, Pagano, Piotto, Pollacco,
  Santos, Smith, Suárez, Szabó, Udry, Adibekyan, Alibert,
  Almenara, Amaro-Seoane, Eiff, Asplund, Antonello, Barnes,
  Baudin, Belkacem, Bergemann, Bihain, Birch, Bonfils, Boisse,
  Bonomo, Borsa, Brandão, Brocato, Brun, Burleigh, Burston,
  Cabrera, Cassisi, Chaplin, Charpinet, Chiappini, Church,
  Csizmadia, Cunha, Damasso, Davies, Deeg, Díaz, Dreizler,
  Dreyer, Eggenberger, Ehrenreich, Eigmüller, Erikson, Farmer,
  Feltzing, de Oliveira Fialho, Figueira, Forveille, Fridlund,
  García, Giommi, Giuffrida, Godolt, Gomes da Silva,
  Granzer, Grenfell, Grotsch-Noels, Günther, Haswell, Hatzes,
  Hébrard, Hekker, Helled, Heng, Jenkins, Johansen,
  Khodachenko, Kislyakova, Kley, Kolb, Krivova, Kupka, Lammer,
  Lanza, Lebreton, Magrin, Marcos-Arenal, Marrese, Marques,
  Martins, Mathis, Mathur, Messina, Miglio, Montalban, Montalto,
  Monteiro, Moradi, Moravveji, Mordasini, Morel, Mortier,
  Nascimbeni, Nelson, Nielsen, Noack, Norton, Ofir, Oshagh,
  Ouazzani, Pápics, Parro, Petit, Plez, Poretti, Quirrenbach,
  Ragazzoni, Raimondo, Rainer, Reese, Redmer, Reffert,
  Rojas-Ayala, Roxburgh, Salmon, Santerne, Schneider, Schou,
  Schuh, Schunker, Silva-Valio, Silvotti, Skillen, Snellen, Sohl,
  Sousa, Sozzetti, Stello, Strassmeier, Švanda, Szabó,
  Tkachenko, Valencia, Van Grootel, Vauclair, Ventura, Wagner,
  Walton, Weingrill, Werner, Wheatley, & Zwintz]Rauer+14
Rauer, H., Catala, C., Aerts, C., et al. 2014, ExA, 38, 249
[Ricker et al.(2014)Ricker, Winn, Vanderspek, Latham,
  Bakos, Bean, Berta-Thompson, Brown, Buchhave, Butler, Butler,
  Chaplin, Charbonneau, Christensen-Dalsgaard, Clampin, Deming,
  Doty, De Lee, Dressing, Dunham, Endl, Fressin, Ge, Henning,
  Holman, Howard, Ida, Jenkins, Jernigan, Johnson, Kaltenegger,
  Kawai, Kjeldsen, Laughlin, Levine, Lin, Lissauer, MacQueen,
  Marcy, McCullough, Morton, Narita, Paegert, Palle, Pepe,
  Pepper, Quirrenbach, Rinehart, Sasselov, Sato, Seager,
  Sozzetti, Stassun, Sullivan, Szentgyorgyi, Torres, Udry, &
  Villasenor]Ricker+14
Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2014, Proc. SPIE,
  9143, 914320
[Schlaufman & Winn(2016)]SchlaufmanWinn16
Schlaufman, K. C., & Winn, J. N. 2016, , 825, 62
[Socrates et al.(2012)Socrates, Katz, Dong, &
  Tremaine]Socrates+12
Socrates, A., Katz, B., Dong, S., & Tremaine, S. 2012, , 750, 106
[Southworth et al.(2009)Southworth, Hinse, Dominik,
  Glitrup, Jørgensen, Liebig, Mathiasen, Anderson, Bozza,
  Browne, Burgdorf, Calchi Novati, Dreizler, Finet, Harpsøe,
  Hessman, Hundertmark, Maier, Mancini, Maxted, Rahvar, Ricci,
  Scarpetta, Skottfelt, Snodgrass, Surdej, & Zimmer]Southworth+09
Southworth, J., Hinse, T. C., Dominik, M., et al. 2009, , 707, 167
[Steffen & Coughlin(2016)]SteffenCoughlin16
Steffen, J. H., & Coughlin, J. L. 2016, PNAS, 113, 12023
[Szabó & Kiss(2011)]SzaboKiss11
Szabó, G. M., & Kiss, L. L. 2011, , 727, L44
[Thommes et al.(2008)Thommes, Matsumura, &
  Rasio]Thommes+08
Thommes, E. W., Matsumura, S., & Rasio, F. A. 2008, Science, 321, 814
[Valsecchi & Rasio(2014)]ValsecchiRasio14
Valsecchi, F., & Rasio, F. A. 2014, , 787, L9
[Walkowicz & Basri(2013)]WalkowiczBasri13
Walkowicz, L. M., & Basri, G. S. 2013, , 436, 1883
[Weiss et al.(2013)Weiss, Marcy, Rowe, Howard, Isaacson,
  Fortney, Miller, Demory, Fischer, Adams, Dupree, Howell,
  Kolbl, Johnson, Horch, Everett, Fabrycky, & Seager]Weiss+13
Weiss, L. M., Marcy, G. W., Rowe, J. F., et al. 2013, , 768, 14
[Wu & Lithwick(2011)]WuLithwick11
Wu, Y., & Lithwick, Y. 2011, , 735, 109
[Youdin(2011)]Youdin11
Youdin, A. N. 2011, , 742, 38

