Zodiacal Exoplanets in Time. X. The Orbit and Atmosphere of the Young "Neptune Desert"-Dwelling Planet K2-100b

The Kepler mission showed that planets are common around other stars and revealed structure in the distribution of planets with size and irradiance (or orbital separation) from their host stars (Fig. 1): This is an important clue to how planets formed, migrated, and evolved. One such feature is the paucity of planets with radii of 2-4$R_{\oplus}$ on close-in (P $\lesssim$3 day) orbits. This region is known as the “Neptune desert" (Lundkvist et al., 2016; Berger et al., 2018), since measured masses and radii of planets in that radius range indicate a hydrogen-helium gas-rich composition analogous to Neptune (Weiss & Marcy, 2014; Rogers, 2015). It has been proposed that the "desert" is the product of the loss of H-He atmospheres by core accretion-powered escape (Ginzburg et al., 2016), X-ray/UV (XUV)-powered (photoevaporative) escape (e.g., Owen & Lai, 2018), or, if the planets accreted in situ (Chiang & Laughlin, 2013), inhibition of accretion of such envelopes by the proximity of the host star, or lack of gas in the inner protoplanetary disk.

Recently, interlopers in the desert have been identified (e.g., Berger et al., 2018; West et al., 2019; Hobson et al., 2019). There are several possible explanations for this population: Planets could move into the desert as a result of the central star evolving from the main sequence and increasing in luminosity; many of the exceptions found by Berger et al. (2018) fall in this category. Some planets may be much more massive, with high molecular-weight envelopes that resist escape to space. They could be the result of catastrophic escape of the H-He envelopes of close-in giant planets ("hot Jupiters") that have shrunk to Neptune size (Dong et al., 2018). The mass of these planets can be measured by precision radial velocity measurements (e.g., Espinoza et al., 2016). Finally, in young planetary systems, escape of the H-He envelope may be incomplete and still ongoing. Since an envelope of at least 5% by mass is required to produce a Neptune-like object (Ginzburg et al., 2016), and the time scale for escape is supposed to be similar to that of decreasing stellar magnetic activity (a few hundred Myr), one expects an escape rate of at least 1 $M_{\oplus}$ Gyr^-1 during that epoch.

The last two possibilities can be tested with observations that constrain the orbit and atmosphere of such a planet. Some hot Jupiter orbits are inclined with respect to the rotational axis of the star, and thus also the presumed plane of the natal disk. Such orbits can arise when the planet is originally driven onto a highly inclined and eccentric orbit by planet-planet scattering, or a Kozai resonance with an outer stellar companion before circularization by the planetary tide (Dawson & Johnson, 2018, and references therein). Sufficiently rapid escape of most of the mass of the giant planet would halt the process of circularization, leaving a Neptune-mass planet on an inclined, moderately eccentric orbit. Orbital inclination or, in the case of a transiting planet, stellar obliquity, can be derived from the Rossiter-McLaughlin effect – the change in a star’s apparent RV due to the partial occultation of the rotating stellar disk by a transiting planet (Triaud, 2017). In the case of rapidly rotating stars it is sometimes possible to perform a detailed analysis of the changes in the stellar line shape, a method often referred to as “Doppler tomography" (Albrecht et al., 2007; Collier Cameron et al., 2010). While large orbital eccentricities of giant planets are readily revealed by RV measurements, modest orbital eccentricities of smaller planets are more difficult to detect. One approach for transiting planets is to compare the duration of a transit to that expected from a stellar density established by other observations and an inferred transit impact parameter (e.g., Van Eylen & Albrecht, 2015; Xie et al., 2016; Mann et al., 2017b; Van Eylen et al., 2019).

Around young stars, ongoing loss of an atmosphere should manifest itself as an extended, escaping atmosphere that could be detected by spectroscopy during transit. Escaping neutral hydrogen (H I) has been detected by observations in the line of Lyman $\alpha$ (e.g., Vidal-Madjar et al., 2003; Lecavelier Des Etangs et al., 2010; Kulow et al., 2014; Ehrenreich et al., 2015) but this can be done only from space and is impeded by the absorption of the line core by interstellar H I, interfering geocoronal Lyman $\alpha$ emission, and variable stellar chromosphere emission. Neutral helium (He I) in its metastable 2s³S or “triplet" state can be detected via absorption in a near-infrared set of lines around 10830Å that are accessible from the ground (Seager & Sasselov, 2000; Oklopčić & Hirata, 2018a). The first detection was with HST from space (Spake et al., 2018) but the introduction of high-resolution infrared spectrographs at facilities on the ground has led to detections or claimed detection of He I among several transiting planets (Allart et al., 2018; Mansfield et al., 2018; Nortmann et al., 2018; Salz et al., 2018; Allart et al., 2019; Ninan et al., 2019) and some non-detections (Moutou et al., 2003; Nortmann et al., 2018; Kreidberg & Oklopčić, 2018; Crossfield et al., 2019).

K2-100b is a 3.5$R_{\oplus}$ planet that was detected by the K2 mission on a 1.67-day orbit transiting the G0-type star K2-100 (TYC 1398-142-1 or EPIC 211990866) (Mann et al., 2017a), and which occupies the Neptune desert (Fig. 1). The star is a member of the nearby (189 pc, Gaia Collaboration et al., 2018) young open cluster M-44, the Beehive or Praesepe. The age of the cluster is not known precisely; it is variously estimated as 625 Myr (Perryman et al., 1998), 850 Myr (Brandt & Huang, 2015), and 590 Myr (Gossage et al., 2018; Schröder et al., 2019), with an uncertainty of about 100 Myr. The planet’s orbital ephemeris were refined with ground-based detections by Stefansson et al. (2018), allowing for transits to be predicted with accuracy of a minute or less. Mann et al. (2017a) pointed out that the transit duration of “b" and the estimated density of the star are suggestive of an eccentric orbit ($e=0.24\pm^{0.19}_{0.12}$), depending on the impact parameter. Stefansson et al. (2018) also derived a transit-based density that was $\sim 3\times$ higher than the one based on multi-band photometry and a Hipparcos parallax of the cluster (Mann et al., 2017a), consistent with an eccentric orbit. However, this is based on transit photometry obtained with a 30-min cadence, only one-third of the transit duration. K2-100b’s host star is highly rotationally variable, challenging both transit and RV measurements. Barragán et al. (2019) used a Gaussian process-regression of RV measurements with spectroscopic indicators of stellar activity to infer a mass of $22\pm 6$$M_{\oplus}$ , i.e. at 3$\sigma$ significance. The mass and radius combination point to a Neptune-like composition and thus a significant H-He envelope which is normally not associated with planets this close to their host star. To constrain the orbit and atmosphere of K2-100b, we analyzed short-cadence photometry obtained by K2 during Campaign 18 (2018), and high-resolution infrared spectra obtained with the IRD spectrograph on the Subaru telescope during a transit in 2019.

In Sec. 5, we translate our non-detection of transit-associated 10830Å absorption into a limit on the amount of triplet-state He I in a planetary wind and thence to a constrain on atmospheric mass loss rate. The metastable triplet (2³S) state of the He I transition is primarily populated by recombination of He ionized by extreme ultraviolet (EUV) photons with energies $>26.4$ eV ($\lambda<504$Å) and primarily depopulated by ionizing near ultraviolet (NUV) photons with energies $>4.8$ eV Å ($\lambda<2583$, Oklopčić & Hirata, 2018b). For this reason, interpretation of our observations requires knowledge of stellar emission and irradiation of its planet in the different UV bands.

To translate our upper limit on 10830Å He I line absorption associated with the transit of K2-100b into statements about mass loss, we modeled the escaping atmosphere of the planet and its photochemistry as a spherically symmetric, isothermal Parker wind (Parker, 1958) with temperature $T_{w}$. Our model largely recapitulates the model of Oklopčić & Hirata (2018a), and we refer the reader to that work for the details, however, we calculate the wavelength-dependence optical depths and photoionization rates explicitly, without averages as approximations. Temperature-dependent recombination coefficients for single and triplet He are from Osterbrock (1989), collision induced transition strengths and energies are from Bray et al. (2000), and neutral collisional de-excitation coefficients are from Roberge & Dalgarno (1982). We assumed a solar-like composition for the wind of 91.3% H and 8.7% He by number. Line profiles were calculated using oscillator strengths from Wiese & Fuhr (2009), assuming Voigt profiles, and integrating over the accessible wind out to the Roche radius, while accounting for the Doppler shift due to the projected wind velocity. These intrinsic profiles do not include the effect of the finite resolution of the spectrograph ($\approx 4$ km s^-1), or the orbital motion of the planet during the transit, which will blur profiles by $2\pi R_{*}t/(P\tau)$ or about 1.5 km s^-1, since these are small due to the expected thermal broadening ($>6$ km s^-1).

Figure 10 plots the densities of electrons (equivalent to the density of ionized H), and He I in the singlet and triplet electronic configurations for the case of $\dot{M}=1$$M_{\oplus}$ Gyr^-1 and $T_{w}=10000$ K. Escaping neutral hydrogen from K2-100b is rapidly and completely ionized by EUV photons from K2-100. We can understand the model prediction that H I is completely ionized (except very close to the planet) by performing a simple Strömgren-sphere-like comparison which competes the rate of impinging photoionizing photons against the rate of escaping H I atoms:

We estimate the rate of impinging photoionizing photons (the left hand side) to be $3\times 10^{35}$ sec^-1, while the escape rate of H I atoms for a mass loss of 1 $M_{\oplus}$  Gyr^-1 is $8\times 10^{34}$ sec^-1. Although this simple calculation (but not the model) ignores recombination, it shows that the EUV irradiation field at the orbit of K2-100b is capable of completely ionizing the flow of H. (We verified that at substantially lower EUV irradiance, the model predicts significant H I in the wind.)

Figure 11 plots the rates of processes that populate (solid lines) or deplete (dashed lines) the 2s³ ground state of the 10830 Å line, normalized by the wind timescale $R_{s}/c_{s}$, where $R_{s}$ is the sonic radius and $c_{s}$ is the isothermal sound speed. “Electron collisions" is the net depopulation of the 2s³S level by electron collisions. The 2s³S state is almost entirely populated by recombination of ionized He (blue solid line) and depopulated by a combination of electron collisional de-excitation (dashed brown line), dominant in the inner part of the wind, and photoionization (dashed green line), dominant in the outer part. Oklopčić & Hirata (2018b) find that photoionization is pervasively dominant, but used a mass loss rate and hence densities that are five times lower than that assumed here.

Figure 12 shows the calculated profile of the triplet He I line complex produced by a planetary wind with $\dot{M}=1$$M_{\oplus}$ Gyr^-1 and $T_{w}=10000$ K. The total EW is 18 mÅ and the two strongest lines of the the three are blended by Doppler broadening ($\approx$15 km s^-1 or about 0.5Å).

These model calculations suggest that for a solar-like composition, an upper limit of total EW $<5.8$ mÅ rules out mass loss rates $>$1 $M_{\oplus}$  Gyr^-1 for all wind temperatures $T_{w}$ considered here, and loss rates $>$0.3 $M_{\oplus}$  Gyr^-1 for $T_{w}$$<10000$K (Figure 13). He I line observations are maximally sensitive to winds with $T_{w}$=4000-5000K, and least sensitive to those at 20000K. Wind temperature affects the sound speed, the sonic point, recombination rates, collisional de-excitation and excitation rates. For a given $\dot{M}$, increasing $T_{w}$ decreases density and hence the EW of the line. The recombination rate coefficients of Osterbrock & Ferland (2006) also decrease with temperature and the density is lower so the actual rate is lower still. Photoionization, the primary mechanism that depopulates the triplet state, is independent of $T_{w}$. The “a" pathway of collisional de-excitation ($q_{31}a$) peaks at 18000K, and this may explain some of the high-temperature structure in Fig. 13, but we caution that both the low- and high-$T_{w}$ calculations depend on extrapolated values for the recombination rates of ionized He to singlet and triplet He I and may not be reliable.

The EUV and Ly $\alpha$ emission of K2-100 have not been measured and presumably will never be directly measured due to the neutral H I along the intervening 186 pc to this star. Both the X-ray and FUV emission used to estimate those fluxes were only marginally ($3\sigma$) detected and formal uncertainties in estimates of the EUV flux are of order unity. It is therefore prudent to consider the sensitivity of the modeling results to changes in the EUV emission at different wavelengths. We calculated the EW of the 10830 Å He I line for a mass loss rate of 1 $M_{\oplus}$  Gyr^-1 and wind temperature of 10000K, and repeated the calculation while doubling the UV irradiance in 10 Å-wide windows over the entire spectral region, recording the fractional change in the calculated EW. We did not vary the mass loss rate, which of course, would also be expected to change with varying EUV flux. Rather, this numerical experiment examines our ability to robustly translate an observed EW (or an upper limit) into a given mass-loss rate and indicates at wavelengths it is most important to establish the stellar emission.

Figure 14 shows the results for a restricted range in the EUV where the sensitivity was the non-negligible. Outside this range, the sensitivity was $<1$%, except for a range in the NUV where increasing irradiance produced a slightly lower EW due to elevated photoionization. In essence, this is a convolution of the EUV spectrum of the star (in this case, a scaled version of the Solar spectrum) and the photoionization cross-section of He I. This shows that it is most important to accurately measure the EUV emission in the 15-500Å range, and in particular at the 305Å He II resonant line for which there is already strong emission from Sun-like stars (Golding et al., 2017). Our planetary wind models show that the 2s³S state of He I is primarily produced by recombination of ionized helium, and the range of non-negligible sensitivity corresponds to the range over which He I ionization occurs ($\lambda<504$Å). This result should serve as a guide for where future efforts to estimate the EUV emission from stars for the purpose of He I line observation should concentrate. Oklopčić (2019) discusses how the He I line strength is governed by the balance between EUV radiation, which produces 2s³S He I via ionization of the ground state and subsequent recombination, and mid-UV radiation, which ionizes tripleet He I. For a fixed mass loss rate, wind temperature, and separation from the host star, there is a spectral type for which the triplet-state He I signal is maximum; hotter stars have more mid-UV emission, and cooler stars do not ionize enough ground-state He I. For $\dot{M}$ of 0.2 $M_{\oplus}$  Gyr^-1 and a wind temperature of 7300K this optimal spectral type is K-type, consistent with the observation that all four successful detections of He I are for planets orbiting K-type stars, and there are only non-detections for planets transiting A-, G-, and M-type stars. Two planets which were not detected, WASP-12b of Kreidberg & Oklopčić (2018) and HD 209485b of Nortmann et al. (2018), orbit G0 stars like K2-100, and are at separations which are well within a factor of two of K2-100b. However, those stars are older, less magnetically active, and presumably would have lower EUV luminosities. The upper limits for WASP-12b and HD209458b, converted into EW, are 24 mÅ (90%) and 8 mÅ, respectively (both $2\sigma)$, comparable to our formal 2.5$\sigma$ limit of 5.7 mÅ.

We can use the UV irradiances calculated and reported in Table 2 to estimate the energy-limited escape rate $\dot{M}_{\rm EL}$ of the envelope of K2-100b, an upper limit on the escape rate in the absence of other energy sources and ignoring processes such as erosion by the stellar wind. The relation (Watson et al., 1981; Erkaev et al., 2007) is:

where $\eta$ is an efficiency factor that accounts for energy loss by radiative and conductive cooling, ionization, as well as the remaining heat content of the escaping gas, $F_{\rm XUV}$ is the combination of X-ray and UV irradiances, $R_{\rm XUV}>R_{p}$ is the effective radius at the level fo the atmosphere which the XUV radiation is absorbed and from which escape occurs, $G$ is the gravitational constant, and $K$ is a correction factor $\approx 1$ that accounts for the fact that the escaping gas only needs to reach the Roche radius (Erkaev et al., 2007). Adopting $\eta=0.1$ (Shematovich et al., 2014) and $R_{\rm XUV}=R_{p}$, and the mass estimate of Barragán et al. (2019), we obtain an estimate for energy-limited escape rate of 1.1 $M_{\oplus}$  Gyr^-1 driven by the combination of X-ray, EUV, and Ly $\alpha$ irradiation. This scenario is ruled out by our observations. X-rays, EUV, and Lyman $\alpha$ may drive drive atmospheric escape with different efficiencies depending on proximity to the host star (Owen & Jackson, 2012), but X-rays are a small fraction of the total energy budget. Asymmetry of the escaping atmosphere may also affect our ability to observe triplet He I and thus limit mass loss. Otherwise, the efficiency of escape could be low or the irradiance significantly lower than estimated. Our analysis of K2-100b’s orbit based on transit observations limit the eccentricity to 0.15 (90% confidence) or 0.28 (99% confidence) and do not support an origin for this object as a more massive planet on a scattered orbit. Perhaps photoevaporative escape, expected to be important during the first few hundred Myr when the star was rapidly rotating and magnetically active (Owen, 2019), has largely halted, leaving the planet stranding in the desert. Alternatively, K2-100b’s atmosphere is escaping at a lower rate but over a longer (Gyr) interval, consistent with the scenario of core accretion-powered escape (Ginzburg et al., 2016).

We thank Evgenya Shkolnik for a helpful discussion on stellar UV emission. EG acknowledges support as a visiting professor from the German Science Foundation through the DFG Research Unit FOR2544 “Blue Planets around Red Stars". A.W.M. and D.O. were supported by the K2 Guest Observer Program (NASA grant 80NSSC19K0097) to the University of North Carolina at Chapel Hill. This work is supported by JSPS KAKENHI Grant Numbers 16K17660, 19K14783, 18H05442, 15H02063, and 22000005, and by the Astrobiology Center Program of the National Institutes of Natural Sciences (NINS) (Grant Number AB311017). This research has made use of the SVO Filter Profile Service (http://svo2.cab.inta-csic.es/theory/fps/) supported from the Spanish MINECO through grant AYA2017-84089, the LASP interactive Solar Irradiance Datacenter, NASA’s Astrophysics Data System Bibliographic Services, Centre de Données astronomiques de Strasbourg, NIST’s atomic line database, Astropy (Astropy Collaboration et al., 2013), and Scipy (Virtanen et al., 2019).