KEPLER PLANETS: A TALE OF EVAPORATION

James E. Owen; Yanqin Wu

doi:10.1088/0004-637X/775/2/105

1. INTRODUCTION

The spectacular success of the Kepler mission has yielded thousands of planetary candidates (e.g., Borucki et al. 2011; Batalha et al. 2013). Most of these exoplanets are smaller than Neptune (4 R_⊕) and are likely to have masses between a few to tens of Earth masses. We can now hope to make substantial progress in understanding the origin and evolution of planetary systems by studying the interior composition and orbital structure of these objects.

However, well-known degeneracies prevent us from resolving the interior compositions of the planets. A typical Kepler planet is likely composed of a dense core and a tenuous envelope: the core can be volatile-rich, rock-rich, or iron-rich; the envelope can consist of steam or hydrogen/helium. Unlike the case of main-sequence stars, radius measurements of these systems (all that is possible for most Kepler candidates) cannot be used to constrain their structure. Even in the case where planetary masses are measured (via transit-timing variations, TTVs, or radial velocity), multiple solutions for the interior structure exist. For instance, a Neptune-mass planet, with a density of 2 g cm⁻³, can have either a thin hydrogen/helium envelope surrounding a dense rocky core or be entirely ice/water dominated (e.g., Adams et al. 2008; Rogers & Seager 2010).

This degeneracy can be broken, however, by studying a population of planets and investigating how planetary sizes and densities correlate with physical environment (Wu & Lithwick 2013). In particular, for those planets at close separations from their parent star (<0.2 AU), the total received high-energy irradiation over Gyr timescales can represent a considerable fraction of their gravitational binding energy (e.g., Lammer et al. 2003; Lecavelier des Etangs 2007; Davis & Wheatley 2009). If these planets have hydrogen-rich envelopes, evaporation of these envelopes could markedly reduce their sizes and increase their bulk densities, compared with planets further out.

Two recent studies have looked for correlations with separation. Ciardi et al. (2013) note that for pairs of Kepler planets bigger than Neptune (4 R_⊕; these planets almost certainly have hydrogen envelopes), the inner planets are more frequently the smaller ones. However, they do not observe such a hierarchy in pairs smaller than Neptune. Therefore, they concluded that the correlations may be a result of the planet formation process, as opposed to any post-formation process. This conclusion contrasts with that of Wu & Lithwick (2013); in this study, the authors measured planetary masses in a sample of TTV pairs and find that, within the pairs, the inner planets tend to be denser. Moreover, when considering all available mass measurements together, they find that planet densities increase with decreasing orbital periods. Therefore, they conclude that the low-mass planets observed by Kepler are composed of dense rocky cores overlaid with various amounts of hydrogen in their envelopes, which are then sculpted by evaporation. The same density and radius trends are encapsulated in many of the Kepler multi-planet systems, e.g., Kepler-11 (Lissauer et al. 2011), Kepler-18 (Cochran et al. 2011), Kepler-20 (Gautier et al. 2012), and Kepler-36 (Carter et al. 2012), as well as reported by recent radial velocity studies (Weiss et al. 2013). Thus, an important question has arisen. Are the observed correlations (radius/density versus distance) results of formation or evaporation?

On the theoretical side, there is a large body of literature discussing planet evaporation, ever since the first hot Jupiters were discovered (e.g., Lammer et al. 2003; Yelle 2004; Tian et al. 2005; Hubbard et al. 2007a, 2007b; Murray-Clay et al. 2009; Koskinen et al. 2010; Owen & Jackson 2012). Most of these studies deal with gas giants and only a few focus specifically on low-mass planets. We summarize here two works that are of direct relevance to our discussion.

Owen & Jackson (2012) develop realistic models of hydrodynamic evaporation for low-mass planets, including both X-ray and EUV irradiation. They perform the calculations on planet models of different densities and radii. They identify that only low mass planets (M_p < 60 M_⊕) could experience a significant effect from evaporation and suggest that evaporation may lead to a stability boundary below which all observed planets should lie (e.g., Koskinen et al. 2007). Lacking a thermal evolution model, they could not make predictions for the final planet radius and density. The detailed evaporation models allow Owen & Jackson (2012) to calibrate the evaporation efficiency and they argue that this parameter can vary by orders of magnitude over the stellar and planetary lifetimes. Their work also highlights the relative importance of X-ray versus EUV evaporation.

Lopez et al. (2012) provide the only significant attempt to study the combined effects of thermal evolution and evaporation, adopting a simplified energy-limited formalism for evaporation and a constant efficiency factor—the caveats of which are discussed in Section 2—and they also infer the presence of an evaporation threshold. They find that since the evaporation timescale scales with planet mass × planet density, planets with mass × density above a critical value and an age older than the threshold should have been evaporated. However, as is demonstrated by Owen & Jackson (2012), a simplistic evaporation model may severely overestimate or underestimate mass-loss rates, especially in the first stages of the stars' and planets' lives.

Now with the large Kepler data set, evaporation of low-mass planets deserves a better study and a direct comparison with the observations. Moreover, it is useful to be able to predict the final planet size and density for a given initial model. This procedure allows us to backtrack the original planetary structures at formation. To accomplish these goals, we need to improve on previous work; in particular, we need to couple realistic calculations of evaporation with thermal evolution models.

We describe our approach in Section 3, after briefly reviewing the evaporation theory in Section 2. We then compare our theoretical results directly with observations in Section 4. Our evaporation theory successfully explains a number of observed facts, establishing that evaporation is the driving process of evolution for close-in Kepler planets. In fact, some 50% of the currently known Kepler planets may have been significantly affected by evaporation. Finally, we discuss the caveats and implications of our work in Section 5 and we conclude in Section 6.

2. PHYSICS OF PLANETARY EVAPORATION

Planet evaporation can take place through a variety of mechanisms: non-thermal escape, thermal Jean's escape, and hydrodynamic escape. Hydrodynamic evaporation occurs when the density of the heated region is sufficiently high so that the gas is collisional even in the supersonic region of the flow. Thus, it is only hydrodynamic evaporation that can produce high enough mass loss to affect the structure and sculpt the planet, which occurs at high incident fluxes; we focus on this mechanism.

As high-energy photons from the star ionize gas in the upper envelope (either hydrogen or metals), the newly freed electrons heat up the local gas and the atmosphere expands. A flow may be initiated that eventually escapes from the gravitational well of the planet. Let the efficiency of converting received energy to PdV work be η, so the mass-loss rate is simply

$\begin{equation} \dot{M}=\eta \frac{L_{\rm HE} R_p^3}{4GM_pa^2}, \end{equation} \tag{ 1 }$

where L_HE is the high-energy luminosity (X-rays or EUV), M_p is the planet mass, R_p is the planet radius, and a is the separation from the parent star.

The usual so-called energy-limited approach is equivalent to taking η to be an order-unity constant (e.g., Watson et al. 1981; Lammer et al. 2003; Lecavelier des Etangs 2007; Erkaev et al. 2007). Under such an assumption, the evaporative timescale is $M_p/{\dot{M}} \propto M_p\rho /F_{{\rm HE}}$ , where F_HE is the high-energy flux. So, Equation (1) lends itself to an evaporation threshold in terms of mass × density (Lopez et al. 2012).

However, the evaporation efficiency is not always constant but may vary significantly with planet mass, radius, and ionizing flux (see Section 5). Murray-Clay et al. (2009) demonstrated that in the case of EUV evaporation of hot Jupiters, the "energy-limited" approximation is only valid at low fluxes. At high fluxes, the mass-loss process is controlled by the ionization/recombination balance, yielding mass-loss rates that scale as $L_{\rm EUV}^{1/2}/a$ or the "efficiency" (η) decreasing with flux. This limit (the "recombination-dominated" regime) is similar to EUV-driven evaporation of gas clumps (Bertoldi & McKee 1990) and protoplanetary disks (e.g., Johnstone et al. 1998; Hollenbach et al. 2000).

For X-ray ionization, Owen & Jackson (2012) similarly showed that the evaporation is not "energy-limited," but rather the "efficiency" is controlled by line cooling, and is a strong function of planet mass. They found that cooling is most important for high-mass Jovian planets, since the higher planet escape temperature and larger physical size mean the flow timescale is long. Then, there is sufficient time to radiate away the received heating when the flow is still subsonic. For lower mass, Neptune-like planets, the escape velocity is much lower and the physical sizes smaller, meaning the flow timescale is shorter, so much less energy is radiated away and the corresponding "efficiency" (η) is higher. While Owen & Jackson (2012) found that X-ray evaporation gives a similar flux scaling as given in Equation (1), the flow is not close to being "energy-limited" (where the PdV work is the dominant energy loss channel), and only begins to approach an energy-limited case at low masses < 3 M_⊕. Furthermore, Owen & Jackson (2012) noted that unlike the EUV case, since the X-ray driven mass loss scaled as L_X/a² at high and low fluxes, there is no transition from line-cooling limited to energy-limited at low X-ray fluxes.

The ionizing flux from a main sequence star varies by orders of magnitude throughout its life. In addition, the size of a planet also changes over time, under both thermal contraction and mass loss, along with the mass of the planet decreasing. With η being a function of planet mass and radius, as well as the ionizing flux (see Figure 12), one should not adopt a constant η (as in Lopez et al. 2012) in tracing out the evaporation history of a planet. It is important that the correct prescription of evaporation is used to follow the planetary evolution, particularly for low-mass planets where evaporation can significantly affect the evolution. This fact is essential if one wants to make inferences about the initial state of these planets.

**Figure 1.** Benchmark calculations of our modified version of mesa against the calculations of Fortney & Nettelmann (2010) for 95 (squares) and 32 (circles) M_⊕ planets at a separation of 0.045 AU from a Sun-like star. Our calculations are shown for three initial cooling times of 10⁵ yr (solid line), 10⁶ yr (dashed line), and 10⁷ yr (dot-dashed line). The disagreement between the models is small, at the <5% level.
Download figure:
Standard image High-resolution image

The last issue worth our attention is the nature of the ionizing flux. Which source of high-energy luminosity is driving the evaporation has only been recently tackled. Some authors have purely used EUV flux in their studies (e.g., Lecavelier des Etangs 2007), while others have used only X-rays (Jackson et al. 2012). Owen & Jackson (2012) solved the flow problem including both X-ray and EUV irradiation. They found that the position of the transition from subsonic to supersonic flow (the "sonic point"), relative to the respective ionization fronts, determines whether X-rays or EUV drive the mass loss. At a similar ionizing energy flux, they found that X-rays drive the mass loss when the flux is high, and EUV dominates the mass loss when the flux is low. For the first ∼100 Myr of their lives, main sequence stars have X-ray luminosities that reach up to ∼10⁻³ of their bolometric luminosities; this fraction falls off with time roughly as 1/t as the stars age (Güdel 2004; Ribas et al. 2005; Jackson et al. 2012).

3. EVOLUTIONARY MODELS OF PLANETARY EVAPORATION

In order to follow the evolution of a close-in planet, we need to include the effects of thermal cooling of the planet, irradiation of its upper atmosphere, and mass loss in the form of hydrodynamic evaporation. Previously, Lopez et al. (2012) addressed this problem by coupling the Fortney et al. (2007b) and Nettelmann et al. (2011) planetary structure models to a simplistic energy-limited estimate for the mass loss. We aim to achieve a similar goal here, but we will base our mass-loss rate on a realistic calculation that includes both the X-ray and EUV radiation (Owen & Jackson 2012). We make use of the mesa stellar evolution code (Paxton et al. 2011, 2013) to model the thermal evolution and couple it to the Owen & Jackson mass-loss calculations.

3.1. Method

The general purpose mesa code provides the framework to simulate planetary structure and evolution. Our planet models have solid inner cores, experience irradiation (under the two-stream approximation; e.g., Guillot 2010), and evaporation. The mesa equations of state are discussed in detail in Paxton et al. (2011), but in the planetary regime the equations of state are typically based on the SCVH equation of state (Saumon et al. 1995). The opacity tables used for the irradiation of the atmosphere are based on an updated version of the Freedman et al. (2008) tables, as detailed in Paxton et al. (2013), where we adopt an opacity κ_* = 4 × 10⁻³ cm² g⁻¹ for the incoming stellar irradiation suitable for a solar-type star (Guillot 2010).

To include evaporation, we tabulate the mass-loss rates as functions of ionizing flux, planet mass, and planet radius, using the results of Owen & Jackson (2012). The tabulated results span a range in planet mass from 1 M_⊕ to 3 M_J and in planet radius from 1 R_⊕ to ∼0.01 AU (roughly the Hill radius at 0.1 AU for our most massive planet). For the ionizing X-ray flux, we adopt the observed relation between this flux and stellar spectral type and stellar age (Jackson et al. 2012). In general, the X-ray luminosity saturates at ∼10⁻³ of the stellar bolometric luminosity during the first ∼100 Myr of a star's life and decays approximately as 1/t afterward. Furthermore, following Owen & Jackson (2012), we take the EUV luminosity to follow the same time evolution as the X-ray luminosity.

Since convective transport is the dominant heat transport mechanism in a planet's envelope, the thermal structure of the envelope is set to be initially adiabatic, in the absence of irradiation. The radius of the planet is defined to be the radius where the optical depth to the incoming bolometric radiation is 2/3, and is typically around millibar (for younger and lower mass planets) to bar (for older and more massive planets) pressures for the planets considered here. This planet radius is also taken as the input radius in the evaporation model. Since the atmosphere's underlying scale height is typically small compared to the planet radius, such an approximation is accurate. For the same reason, we have ignored the difference between the above radius and the radius a planet exhibits at transit (Hubbard et al. 2001). However, such a simplification may break down for the closest in planets at the earliest times, where there may be a ∼10% difference in a planet's optical photosphere and the base of the evaporative flow. However, given that planets quickly contract through such a phase, the effect will be negligible when integrated over the Gyr of planetary evolution.

Our planets are composed of two separate regions: an envelope of hydrogen/helium with a solar abundance metallicity (as used for the evaporation models) and a solid core. Given the large range of possible core compositions, as a starting point we focus on a pure rock core, using the mass–radius profile provided by Fortney et al. (2007a, 2007b), and adopting core masses (M_c) of 6.5, 7.5, 10, 12.5, and 15 M_⊕. These cores are assumed to be of fixed radius and do not evolve during the planet's evolution. However, the cores do have a thermal content from both radioactive decay and thermal heat capacity; we adopt an Earth-like value (see Nettelmann et al. 2011; Lopez et al. 2012). Both these thermal sources are implemented in mesa using the core luminosity function (Paxton et al. 2013). As current planet formation models are unable to provide a good handle on the initial thermal properties of formed planets (i.e., a "hot" or "cold" start), we consider models with a wide range of initial radii. These initial radii (or, more correctly, entropy) are parameterized in terms of an initial cooling time (t_cool), which we define as the ratio of a planet's initial internal energy to its initial luminosity. We could of course parameterize the initial entropy in terms of some other variables; however, the initial cooling time is perhaps the most interesting to compare with protoplanetary disk lifetimes of ∼3 Myr (e.g., Haisch et al. 2001; Hernández et al. 2007; Owen et al. 2011; Armitage 2011).

In order to make sure our model is accurately following the evolution of low-mass planets, we benchmark our calculations against the models presented in Fortney & Nettelmann (2010, their Figure 10); the calculations were performed using the Fortney et al. (2007a) code for a 95 and 32 M_⊕ planet with a 25 M_⊕ core³ that is a 50/50 mix (by mass) of ice and rock. We follow the evolution of these planets under the influence of irradiation by a Sun-like star, but no evaporation. We find excellent agreement with these calculations at the ⩽5% level, giving us confidence that our modified version of mesa is performing as expected at low masses. In Figure 1, we show the radius evolution of our calculations for cooling times of 10⁵ yr (solid line), 10⁶ yr (dashed line), and 10⁷ yr (dot-dashed line). These calculations are compared against the results from Fortney & Nettelmann (2010) at a separation of 0.045 AU, shown as points for both the 95 (squares) and 32 (circles) M_⊕ planets.

For our calculations, we choose values of t_cool (computed in the absence of irradiation) in the range 3 × 10⁶–10⁸ yr to span a range of "hot" start and "cold" start models, similar to those values chosen by Lopez et al. (2012). We note that with cooling times <3 × 10⁶ yr, most low-mass planets at separations <0.1 AU have initial radii larger than their Hill radii, and they cannot be considered hydrostatically bound objects.

For each cooling time and core mass, we construct 40 models with envelope masses ranging from 3 M_J to a few percent of the core mass, logarithmically spaced in envelope mass. The planets are then evolved forward in time, under the influence of evaporation and irradiation for 10 Gyr or until the entire envelope is lost.⁴ Our integrations begin at 3 Myr, a time at which the protoplanetary disk clears and the planet is fully exposed to X-ray and EUV irradiation. In general, a planet's evaporation begins in the X-ray driven phase and may switch to EUV driven at some late time, where the planet's final mass is set by a few hundreds of Myr.

3.2. Jupiter-like Planets

Jupiter mass planets close to their central stars represent the case where there is direct observational evidence (e.g., Vidal-Madjar et al. 2003, 2004) of evaporation occurring. Thus, it is worth investigating whether there are any evolutionary consequences for the evaporation of high-mass planets. For example, Baraffe et al. (2004) noted that if the evaporation rate was high enough, so that the evaporation time ( ${\sim } M_p/\dot{M}$ ) became comparable to the thermal time of the planet's envelope, then a Jupiter-mass planet could lose its entire envelope rapidly. Such an inference was based upon the rather unrealistic assumption of 100% mass-loss efficiency. In reality, at high masses, the mass-loss rates never reach such high values (see Owen & Jackson 2012 for a more detailed discussion). To illustrate this fact, we show the evolution of a Jupiter-like planet at a very close separation of 0.025 AU around a solar-type star in Figure 2. We plot the radius and mass evolution of planets with a 15 M_⊕ rock core and cooling times ranging from the very short (10⁶ yr, solid line) to the very long (10⁸ yr, dotted line). Figure 2 clearly shows that evaporation is unable to affect the planet's evolution and planets with very different initial cooling times end up on almost identical evolutionary tracks at late times, with the mass change in the planets at the <1% level, as argued for by Hubbard et al. (2007a).

**Figure 2.** Top panel: radius evolution of Jupiter-like planets as a function of time since gas disk clearing, at a very close separation to the parent star (∼0.025 AU); the bottom panel shows the mass evolution of this planet. The solid line represents a planet with an initially high entropy with an initial cooling time of 10⁶ yr, the dashed line represents an initial cooling time of 10⁷ yr, and the dotted line represents an initially low entropy with an initial cooling time of 10⁸ yr.
Download figure:
Standard image High-resolution image

Therefore, while Jupiter-mass planets currently provide the best opportunity for actually studying the hydrodynamics of evaporation by directly probing the flow, they do not provide a good laboratory for studying the evolutionary consequences of evaporation and we must go toward lower mass planets where the evolutionary effect will be more pronounced.

3.3. Low-mass Planets

Unlike the Jupiter-type planets discussed above, several authors have suggested that the effects of evaporation will be more prominent for lower mass planets (e.g., Hubbard et al. 2007a; Baraffe et al. 2008; Jackson et al. 2012; Owen & Jackson 2012; Lopez et al. 2012). At lower masses, evaporation can begin to sculpt the planet population, removing significant amounts of a planet's envelope during its lifetime. To investigate this effect, we show the evolution of a low-mass planet in Figure 3, where we consider the evolution of a "standard" model that is an initially 20 M_⊕ planet with a 12.5 M_⊕ core. Figure 3 shows the evolution of the "standard" model for the full range of initial cooling times considered. The planets are at a close separation of ∼0.05 AU and have equilibrium temperatures of 1300 K, where we define the equilibrium temperature as the blackbody temperature at a given separation, i.e.,

$\begin{equation} T_{\rm eq}=T_*\sqrt{\frac{R_*}{2a}}, \end{equation} \tag{ 2 }$

where T_* and R_* are the parent star's temperature and radius, respectively.

The panels in Figure 3 show the qualitative features of the evolution of low-mass planets. In particular, planets with higher initial entropies have initially larger radii and therefore lose more mass. The result is that planets with shorter initial cooling times end up with smaller radii and higher densities than planets with longer initial cooling times. At low masses, evaporation plays a strong role in planetary evolution, which in the case of our "standard" model makes the planets a factor of approximately two smaller and a factor of approximately four denser compared with the same planet that is not undergoing evaporation. In this case, an initial envelope containing ∼40% of the original mass is evaporated down to an envelope containing only few percent of the total mass; in the most extreme case (for the planet with an initial cooling time of 3 × 10⁶ yr), evaporation leaves a planet with 0.5% of the total mass in a hydrogen/helium envelope.

The bottom right-hand panel of Figure 3 shows the evolution of $\dot{M}\times {\rm age}$ , which indicates when mass loss is most significant. This plot shows that in all cases the mass loss is most important at roughly the point where the X-ray luminosity begins to decline. This fact is easy to understand since at early times planets are large and fluxes are high, so a planet can absorb a significant amount of high-energy radiation. Once the X-ray fluxes begin to decline, the planet evolution is less affected by evaporation. Once the evaporation switches to EUV-driven, the final planet properties have already been "frozen" in, similar to the results from the previous models calculated by Owen & Jackson (2012) and Lopez et al. (2012). This feature is rather generic to all evolutionary models, where the saturation timescale for the X-rays sets the length of time over which evaporation is important in driving planetary evolution; there is very little change in planetary masses from 100 Myr to 10 Gyr.

4. POPULATION STUDY AND COMPARISON WITH OBSERVATIONS

We have shown that for low-mass planets close to their parent star, a significant amount of a planet's gaseous envelope, or even its entire atmosphere, can be removed over Gyr timescales.⁵ Given the prominent role evaporation can play, we can use this evolution to make inferences about the initial structure of observed close-in planets that may provide clues as to their formation. In the following, we provide comparisons between our theoretical models and observations, mostly using data from the Kepler mission. The effects of evaporation are clearly visible in the data. In fact, the evaporation model naturally explains a number of remarkable correlations seen in the Kepler data.

In Figure 4, we plot the final planet mass and radius as a function of separation for planets with initial properties that vary from our "standard" model, where the central star is 1 M_☉ in all cases. At a given planet mass, a lower core mass results in initially larger planets, which means that the planet can absorb a larger fraction of the X-ray flux and drive stronger evaporative flows. So, planets with smaller cores may lose their entire envelope at larger separations than planets with larger cores (see Figure 8). The initial cooling time shows the same effect as discussed above, with shorter cooling times resulting in more mass loss due to the initially larger radii; however, this effect is much less important than the effect of initial core mass. Finally, evaporation can drive convergent evolution. For instance, planets at <0.05 AU that initially had 15 M_⊕ cores, but envelopes of 2 and 5 M_⊕, respectively, end with somewhat identical structures. This result means that, at least in some cases, it may be difficult to retrieve a planet's initial structure, particularly at low envelope mass fractions.

4.1. An Upper Envelope in Planet Radius

The Kepler transiting-planet catalog is the most extensive catalog of planets in close orbits around their parent star; for better statistics, we use the Kepler object of interest (KOI) catalog, which lists the radii (not masses) of planet candidates. Most of the KOIs have not been confirmed as planets and there is a certain, but low, percentage of false positives (Morton & Johnson 2011). To minimize contamination, we choose to consider only KOIs that have been identified as being in multiple systems, where the significance of a planetary nature is considerably higher (≳95%; Lissauer et al. 2012). The use of only multi-planet systems, while the most robust, may introduce some implicit biases. In particular, planets that may have been dynamically moved to small orbital periods at late times will have undergone a different evolutionary path. Since we are interested in the long-term evolution of planets due to evaporation, the multi-planet KOIs represent the cleanest sample to begin with.

When we plot planet radius versus separation for planets in multi-planet systems around solar-type stars (T_* = 5200–6200 K), we notice an upper envelope in planet size that rises with separation. This effect has also been noted by Ciardi et al. (2013) and Wu & Lithwick (2013). The lack of large planets at small separations is statistically significant: we cut the sample into four regions (large and small planets with a division at R_p = 5 R_⊕ and hot and cold planets with a division at T_eq = 1000 K) and then compare the ratio of large, cold planets to the ratio of small, cold planets. Where one would expect to observe ∼20 hot, large planets if these are the same populations, we find only three.

We also argue that selection effect in the Kepler pipeline would not produce such an upper envelope. KOI candidates have to have a certain signal-to-noise ratio to be identified. This requirement favors detecting smaller planets closer to stars, as they have produced a larger number of transits during the mission lifetime. The completeness radius (sizes above which the KOI catalog is complete) will rise with planet period, similar to the upper envelope discussed here. We draw such a curve for the sub-sample considered by Petigura et al. (2013) in Figure 5, which lies far below the upper envelope. Though more noisy stars than those considered by Petigura et al. (2013) will move this curve upward, the movement will not explain the upper envelope.

**Figure 5.** Upper envelope of planet sizes as a function of equilibrium temperature. The open circles show the KOIs that are around solar-type stars (T_* = 5200–6200 K) and are in multiple transiting systems. The black curves are the theoretical final radii for planet models with initial masses of 20 M_⊕ and 10 M_⊕ (dashed line), 12.5 M_⊕ (solid line), and 15 M_⊕ (dot-dashed line) rocky cores. The dotted line corresponds to ∼30 M_⊕ planets with 12.5 M_⊕ cores. All models here have an initial cooling time of 10⁷ yr. The thin line shows the 50% completeness limit calculated by Petigura et al. (2013), and extrapolated to larger and smaller separations (solid line and dashed line, respectively).
Download figure:
Standard image High-resolution image

This upper envelope is naturally explained by evaporation: low density planets cannot survive in an environment of high ionizing flux. Quantitatively, evaporation of 20 M_⊕ planets with rocky cores of masses 10–15 M_⊕ provides a good fit to the upper envelope, as shown in Figure 5.

There is another implication to this agreement. If there were a significant population of planets with initial masses higher than, say, 30 M_⊕ (the dotted line in Figure 5), we would not expect to observe the same upper envelope. These more massive planets can hold on to their atmospheres, much like the hot Jupiters can (see Section 3.2), and would populate the upper left region in Figure 5. The absence of these more massive planets is interesting and perhaps not coincidental, as 20–30 M_⊕ roughly corresponds to the gap-opening mass in this region, suggesting that there is a ceiling to how much gas the planets can accrete in this neighborhood. The core mass also somewhat affects the upper envelope. We find that models with roughly half of their total mass in rocky cores best reproduce the data, consistent with the conclusion in Wu & Lithwick (2013). In contrast, the planet's initial cooling time makes little difference in the results.

4.1.1. Dependence on Stellar Spectral Types

The mass of the parent star is an important consideration. The parent star directly influences the evaporation, although its gravity is small. However, the total received X-ray flux at a fixed bolometric flux (a fixed equilibrium temperature) can vary greatly, with late-type stars being significantly more X-ray luminous when integrated over Gyr timescales (Güdel 2004; Güdel et al. 2007; Jackson et al. 2012). We demonstrate this fact by considering the evolution of the "standard" planet discussed above (M_p = 20 M_⊕, M_c = 12.5 M_⊕) at a fixed equilibrium temperature (1300 K) around a later-type (0.5 M_☉) and earlier-type (1.5 M_☉) star compared with a planet around a solar-type star, all with initial cooling times of 10⁷ yr.

The evolution of these planets is shown in Figure 6, where the radial evolution is shown in the top panel and the mass evolution is shown in the bottom panel.

Naively, one would expect a similar evolution as the bolometric flux received is identical in all cases. However, the variation of X-ray luminosity with stellar mass results in qualitatively different evolutionary paths for the planets, with order-unity differences in both the final planet mass and radius. The planet around the 0.5 M_☉ star has had its envelope completely removed, whereas the planet around the 1.5 M_☉ star still has a ∼3 M_⊕ envelope remaining after 10 Gyr.

We can go further and compare our evaporation threshold found above for solar-type stars with KOIs around other types of host stars. Lower mass stars (e.g., M dwarfs) have higher X-ray fluxes compared to their bolometric luminosities. They also remain choromospherically active for longer periods. Figure 6 shows that, indeed, at the same equilibrium temperature (measuring the bolometric flux), the upper envelope around lower mass stars appear to be associated with smaller planet sizes.

Currently, the numbers of candidates around A/F and M stars are not as large as those around G/K stars. Thus, a fully quantitative comparison is not possible at this stage. Moreover, planet radii determination around M stars suffer from large uncertainties (Muirhead et al. 2012; Mann et al. 2012; Morton & Swift 2013) and the radius determination for hot stars can also be polluted by the presence of sub-giants (Brown et al. 2011); it is important to bear these effects in mind when drawing inferences. We determine the theoretical evaporation threshold by extracting a linear relation between the maximum radii and equilibrium temperatures for planets with an initial mass of 20 M_⊕ and a core mass in the range 10–15 M_⊕, that are orbiting around a solar-type star and have equilibrium temperatures in the range of 500–2000 K. Since the total mass loss roughly scales linearly with the integrated X-ray flux, we expect the same radius threshold to apply to planets around all spectral types when we arrange them by their X-ray exposure. This result is shown in Figure 7, where the planets are separated according to the spectral type of their parent star.⁶ In contrast, we show the same planet radii plotted against their bolometric exposure. Planets satisfy the same evaporation threshold only when one considers X-ray exposure. This result argues that the ionizing flux, not the bolometric flux, is what determines the upper envelope. Furthermore, we also show the KOIs in single-planet systems in Figure 7 as small filled circles, which show the same behavior as the multi-planet KOIs, although with slightly more scatter.

**Figure 7.** Observed radii of KOIs plotted against their bolometric exposure (bolometric flux received at the surface, integrated over 10 Gyr; left panels) and X-ray exposure (same but for X-ray exposure, using Jackson et al. 2012 values; right panels). The objects are roughly separated into those around A/F stars (top panels), G/K stars (middle panels), and M dwarfs (lower panels; the red squares represent the system KOI 961). Multi-planet KOIs are shown as open circles while single-point KOIs are shown as small dots. As predicted by our models of X-ray evaporation, planets around stars of different spectral types, while having very different bolometric exposures, have roughly the same maximum sizes at a given X-ray exposure. The solid lines show the theoretical evaporation threshold derived from a linear fit to the black evaporation curves shown in Figure 5. The right panels demonstrate that planets around different types of stars satisfy the same maximum radius/X-ray exposure relation as that found around solar-type stars. One cannot perform the same exercise in the left panels.
Download figure:
Standard image High-resolution image

4.2. Distribution of Radius

Our analysis in the previous sections suggests that the observed radius cut-off in KOIs is related to the fact that the most massive KOIs⁷ have masses not much exceeding 20 M_☉ and that their core masses are roughly half of their total masses. Here, we investigate the nature of all KOIs by studying the overall radius distribution.

In Figure 8, we present the final radii of multiple sequences of planet models with initial cooling times in the range 3 × 10⁶–10⁸ yr. These have core masses from 6.5 M_⊕ to 15 M_⊕, and atmosphere masses from approximately 1% of the core mass to much larger values. The total planet masses are restricted to <20 M_⊕. We do not consider atmospheres with masses below 1%, motivated by the discussion below. The radii are evaluated after 10 Gyr of orbiting around a Sun-like star, although the values differ little if we instead evaluate at 1 Gyr (e.g., Lopez et al. 2012; Lammer et al. 2013). Figure 10 shows the corresponding planet densities as a function of separation and this figure is discussed in Section 4.3.

The overall population shows the general feature noted previously; radius decreases with decreasing separation. One particularly interesting feature that appears is a gap in radius between planets that have gaseous envelopes and those that are bare cores. For instance, for the 6.5 M_⊕ core models, inward of ∼0.1 AU, all planets have their atmospheres stripped away with their final sizes reflecting that of their naked rocky cores; outside ∼0.1 AU, planets can retain atmospheres that are at least a fraction of a percent in mass, consequently they have sizes that are markedly larger. This bifurcation generates a gap in planet radius. The orbital separation at which this gap appears is smaller for planets that have bigger cores and stronger surface gravities. However, inside ∼0.03 AU, there are no surviving planets with gaseous envelopes.

The origin of this gap is easy to understand. As planets lose their hydrogen atmospheres, they become increasingly compact and dense, which reduces the mass-loss rate. However, there is also less hydrogen to lose. So, for planets inside the critical separation, the loss is total (e.g., Baraffe et al. 2006; Lopez et al. 2012), while for planets just outside the critical separation, there is a minimum atmosphere mass the planets need to avoid complete stripping. Any thinner atmosphere will be easily eroded. This mass is roughly 1% of the total mass and corresponds to roughly an order-unity modification to the planet radius. So, we expect to see a gap in planet size. Such a gap may become less pronounced when different core compositions are considered, but the basic property that small atmospheres are unstable to complete evaporation will always result in a region where planets are unlikely to end up. Observationally determining the presence of such a gap will place strong constraints on the model and further characterizing the gap will enable useful inferences about the dominant core mass/composition.

In Figure 9, we demonstrate that the radius distribution of KOI multi-planet systems is suggestive of a bimodal distribution. Most planets have sizes ∼1.5 R_⊕ or 2.5 R_⊕, with a deficit of objects at radius ∼2 R_⊕, indicative of the presence of a gap. To further test this suspicion, we divide the objects by their X-ray exposure into a high X-ray group (corresponding to <0.1 AU around a solar-type star) and a low X-ray group. Strikingly, objects with a high X-ray exposure mostly have sizes below the gap, at ∼1.5 R_⊕, while objects with a low X-ray exposure typically have sizes above the gap, at ∼2.5 R_⊕. This result argues that the deficit at 2 R_⊕ could be physical and is associated with X-ray exposure.

The same bimodal behavior appears when we consider only KOI singles, both single and multi-planet KOIs, or when we include only bright KOI targets or only dim KOI targets. In addition, the same behavior is seen when we restrict ourselves to planets that have periods shorter than 50 days and sizes above 1.3 R_⊕, a group of KOIs that suffer relatively little incompleteness (Petigura et al. 2013; Fressin et al. 2013). Although most careful studies to date have yet to have sufficient radius resolution to confirm this feature (Howard et al. 2012; Petigura et al. 2013; Fressin et al. 2013), it is hinted at by Morton & Swift (2013), who construct a probability distribution function rather than using histograms with large ranges. Currently, evaporation is the only process that can naturally explain such a bimodality. Any other processes (planet gas accretion, migration, orbital instability, planetary mergers) may lead to a correlation between planet size and location, but will not produce two separate peaks in radius.

The presence of this gap provides strong evidence that evaporation not only sculpts the upper envelope of planet sizes, but that it also drives the evolution of the majority of KOIs. If all the planets in the higher X-ray exposure peak originally had significant H/He envelopes comparable with the planets with lower X-ray exposures, then ∼50% of Kepler planet candidates have experienced significant mass loss during their lifetimes. Comparing the gap location (0.1 AU around a Sun-like star) against our theoretical calculations, we suggest that the planet population in the current KOI list⁸ have predominately low mass cores (∼6M_⊕), and that most started out their lives with hydrogen/helium envelopes of at least a few percent in mass. Certainly, the radius distribution of close-in planets requires further work—along the lines of Morton & Swift (2013)—before definitive conclusions can be drawn. An accurately determined bimodal distribution encodes valuable information about the initial and final planet mass and composition. Figure 8 illustrates that the gap radius, as well as the separation at which this gap appears, are direct probes of the core properties. An improved investigation on core composition and mass should be conducted when planet radii are better determined.

4.3. Comparison of Planet Density

While the KOI catalog only allows comparison of planet radius, the small but growing sample of low-mass planets with measured masses also provides another important comparison: planet density. In Figure 10, we illustrate the planet densities resulting from our integration. In the density-separation plane, the upper envelope in planet size is now translated into a lower envelope in planet density. There is a gap in planet density, similar to that in radius. However, if there is a planet population with low-mass cores (1–3 M_⊕), these low density cores may partially fill in the gap.

**Figure 10.** Same as Figure 8, but showing final planet density.
Download figure:
Standard image High-resolution image

To compare, in Figure 11, we plot the measured densities of a list of low-mass planets against their scaled separations. We scale the actual planet separations by their respective X-ray exposures, although this correction is typically small (<10%) since most host stars are solar type. The theoretically disallowed low-density region is shown. Most of the observed densities avoid this forbidden region, providing strong evidence for sculpting by evaporation.

Current density measurements (especially those using transit-timing fits) contain large uncertainties. This fact prevents us from making more quantitative comparisons at the moment. In particular, we could not discern the density gap as predicted by model calculations.

The particular case of Kepler-36bc (Carter et al. 2012) is worth commenting on. The two planets have nearly identical separations (0.115 AU and 0.128 AU) but largely discrepant densities, 6.8 and 0.86 g cm⁻³, respectively. At their present orbits, the minimum core masses for the two planets to retain their envelopes is ∼6.5M_⊕ (see Figure 10). The measured masses of the two planets are 4.45 M_⊕ and 8.08 M_⊕, respectively, naturally explaining the diversity in their structure. More systems like Kepler-36 will be able to provide strong constraints on the nature and strength of evaporation in close-in planetary systems.

5. DISCUSSION

We discuss some of the caveats and limitations of the presented calculations and how they bear on our inferences about the observed planet population. The two most important assumptions concern the use of a two-layer planetary model (a rock core plus a hydrogen/helium envelope) and the adoption of the Owen & Jackson (2012) evaporation model.

5.1. Variations in Planetary Structure

Low-mass planets could contain volatile-rich atmospheres, and a significant amount of iron and/or ice in their cores (Adams et al. 2008; Rogers & Seager 2010). Our simple two-layer model of rocky cores plus a hydrogen envelope has been successful in explaining the observations. But what about these other possibilities? Can we exclude them based on current observations?

We can exclude the possibility that the dominant primordial atmospheres of these low-mass planets are very rich in volatiles. Our arguments below run similarly to those given in Wu & Lithwick (2013) but are more informed by our detailed modeling and by the physics of evaporation. If one considers a water-rich envelope, for instance, evaporation will not proceed as is described here. First, water molecules have to be photodissociated, then oxygen has to settle out to produce a nearly pure hydrogen upper atmosphere (e.g., Kasting & Pollack 1983). If oxygen is present in the evaporative flow at a high enough concentration, it will produce strong cooling and increase the opacity to X-rays. This fact will severely suppress the gas temperature, leading to a much lower evaporation rate, similar to what occurs in high metallicity protoplanetary disks (Ercolano & Clarke 2010), but more extreme. Now suppose that all these conditions are satisfied and all hydrogen in the water atmosphere is lost, thus removing 1/8 of the atmosphere mass. However, since the original water atmosphere has a low scale height, such a removal can hardly change either the bulk size or the bulk density of the planet. One would therefore not be able to explain the upper envelope in the observed planet radius or the bimodal planet size distribution or the correlation between planet density and X-ray exposure in terms of planetary evaporation.

In contrast with atmospheric composition, we can be less certain about the core composition. For planets with hydrogen envelopes that are more than a percent in mass, the planet sizes are not sensitive to the core sizes (or equivalently, to the core composition). Planets that have been evaporated down to bare cores may be able to inform us about the core composition, if the core masses are known. A number of these objects show densities that are compatible with rocky or iron/rock compositions (e.g., Hatzes et al. 2011; Batalha et al. 2011). More detailed investigations are necessary to ascertain the core compositions.

5.2. Improving the Evaporation Model

As we have discussed previously, evaporation is key to the evolution of close-in planets. Thus, it is important to assess the role the assumed evaporation model plays in our conclusions. Most previous attempts to model the evolution of evaporating planets use a constant evaporation efficiency (η), which is typically taken to be ∼10% (e.g., Lopez et al. 2012). We have argued that this efficiency depends on planet mass and radius, as well as on the X-ray flux. We further demonstrate this point in Figure 12 by showing how the efficiency changes as a planet evolves. Following four planets with the same initial mass of 20 M_⊕ (but different core masses), we find that η can decrease by a factor of four as the planets evolve from the early puffy stage to the later denser stage, though the variations are not strictly monotonic in time. The overall decrease can be understood that as the planet evolves due to mass-loss and thermal contraction, the planet's density and the surface escape velocities increase with time. Therefore, it takes longer to accelerate the flow to the escape velocity. This fact leads to enhanced cooling and lower efficiencies. The non-uniform evolution of some of the models at early times is due to planets straddling the peak efficiency line (roughly when T_gas ∼ T_escape; see Owen & Jackson 2012) during their evolution and moving above and below it at early times.

**Figure 12.** Evaporation "efficiency" (η) as a function of time. Here, we trace the evolution of four 20 M_⊕ planets, with core masses of 7.5 M_⊕ (solid line), 10 M_⊕ (dashed line), 12.5 M_⊕ (dotted line), and 15 M_⊕ (dot-dashed line), at a separation of ∼0.03 AU from a Sun-like star. All models lose their entire envelopes at large times.
Download figure:
Standard image High-resolution image

While the fixed efficiency of ∼10% adopted by, e.g., Jackson et al. (2012) and Lopez et al. (2012) does represent the approximate median value during the evolution, it can result in order-unity inaccuracies in the integrated mass loss. Therefore, any inferences about the initial planet structure should be made with caution.

5.2.1. Accuracy of the X-Ray Model

The Owen & Jackson (2012) model for X-ray evaporation contains a number of assumptions that may impact our conclusions here. First, only the soft <1–2 keV photons are responsible for heating, while the X-ray flux refers to the entire observed X-ray spectrum (0.1–10 KeV). If the adopted X-ray spectrum—while based on observed spectra (see Owen et al. 2010)—is overly soft or overly hard, then the X-ray heating will be overestimated or underestimated, respectively. Second, these calculations assumed that in the X-ray region, molecules are photodissociated or thermally dissociated by UV photons. Therefore, Owen & Jackson (2012) neglected cooling associated with molecular species. This fact is an important assumption that has not been vigorously tested. If instead the molecular species provide significant cooling (in the temperature range of 2000–5000 K), the X-ray driven flow will not reach T > 5000 K and may remain subsonic all the way out to the EUV ionization front. This fact will markedly reduce the mass-loss rate. More detailed calculations are necessary to address this issue in the future.

5.2.2. The Role of EUV Evaporation

We find that EUV evaporation contributes <10% of the total mass loss during a planet's evolution. The EUV portion of the stellar flux, although comparable in energy to X-rays, is much less efficient at driving a high integrated mass loss. This fact is because in the EUV flow region, the recombination of hydrogen and the subsequent cooling bleeds much of the energy to space. The radiation hydrodynamics are similar to those of an H ii region where ionization is balanced by recombination and the temperature is controlled by the cooling thermostat to ∼10⁴ K (Murray-Clay et al. 2009).

We re-calculate the evolution of our "standard" model (see Section 3.3) undergoing EUV evaporation only, the results of which are shown in Figure 13. At high EUV fluxes (≳ 10⁴ erg s⁻¹), the mass loss scales as $L_{\rm EUV}^{1/2}/a$ (Murray-Clay et al. 2009; Owen & Jackson 2012) and is significantly lower than the values that apply to the X-ray flow. This fact suppresses the mass loss at early times by a factor of ∼10, relative to the X-ray model. The integrated mass loss is ∼1M_⊕, as opposed to ∼7 M_⊕ in the X-ray model. As such, EUV evaporation alone is unable to sculpt the close-in planet population and cannot explain the observed features discussed in Section 4.

6. CONCLUSIONS

In this work, we have coupled the hydrodynamic evaporation models of Owen & Jackson (2012) to the stellar evolution code mesa in order to follow the mass and radius evolution of low-mass planets orbiting close to their stars. Evaporation, while having little effect on massive planets, can remove the entire hydrogen envelopes of the hottest low-mass planets. In all cases, we find that X-rays are the dominant sculpting force and that most of the mass loss occurs in the first few hundreds of Myr when the stars are most chromospherically active and the planets are still contracting thermally.

Our main conclusions are as follows.

1.
Evaporation produces an upper envelope in planet radius as a function of separation. The location of this upper envelope depends on planet mass and the X-ray luminosity of the host star. In particular, M dwarfs, having proportionally larger X-ray fluxes, should have stronger evaporative powers than those indicated by their low bolometric luminosities.
2.
To closely reproduce the observed upper envelope in Kepler candidates, both around Sun-like stars and around stars of other spectral types, we require that the most massive hot Neptunes should have total masses not much exceeding ∼20 M_⊕ and core masses roughly half of that.
3.
Very close-in planets can be stripped of their entire atmospheres. The boundary between complete loss and planets that can retain at least ∼1% of their atmosphere lies at ∼0.1 AU if the core masses are ∼6 M_⊕. At this distance (where most of the Kepler planets lie), a thinner envelope cannot survive. So, we expect a gap in planet size distribution—and this result is suggestively seen in the Kepler catalog, where there is a deficit of planets at ∼2 R_⊕ straddling objects with high X-ray exposures and those with low X-ray exposures. Comparison with our models suggests that most of the Kepler planets should have core masses ∼6 M_⊕ and should have primordial hydrogen/helium envelopes of at least a percent in mass. Moreover, about half of the Kepler planets belong to the high X-ray group and have been stripped down to naked cores.
4.
Evaporation naturally explains the observed correlation between planet density and separation. At closer separations, planets in general have higher densities.
5.
Combining all evidence that supports evaporation, we argue that Kepler planets were born with hydrogen envelopes, not volatile-rich atmospheres (also see Wu & Lithwick 2013).

Looking ahead, we expect that the approach we adopt here—coupling thermal evolution and evaporation—will perhaps be the only hope we have for recovering the initial structure of low-mass planets. With better determined stellar radii, and hence more accurate planetary radii, we may be able to retrieve the initial distribution of planet total masses and core masses. With more measurements of planet densities, it may be possible to reconstruct the histories of individual planets (as has been done for Kepler-11 by Lopez et al. 2012). A larger sample of planets with measured masses and radii will allow us to place constraints on the core compositions, as well as the initial planet entropies. It is at this point that we will begin to learn valuable information about the planet population at birth and make inferences about the planet formation process (e.g., Ida & Lin 2005; Mordasini et al. 2012a, 2012b; Hansen & Murray 2012, 2013).

We are grateful to the anonymous referee for comments that helped improve the paper. We thank Norman Murray, Chris Thompson, Jason Rowe, Alan Jackson, Eric Lopez, Adam Burrows, and Yoram Lithwick for interesting discussions. The calculations were performed on the Sunnyvale cluster at CITA, which is funded by the Canada Foundation for Innovation. Y.W. acknowledges support from NSERC and the Province of Ontario.

KEPLER PLANETS: A TALE OF EVAPORATION

Article metrics

Permissions

Author e-mails

Author affiliations

Dates

ABSTRACT

1. INTRODUCTION

2. PHYSICS OF PLANETARY EVAPORATION