X-ray irradiation of three planets around Hyades star K2-136

We study the X-ray irradiation and likely photoevaporation of the three planets around the star K2-136. These are the Earth-sized K2-136 b, the mini-Neptune K2-136 c, and the super-Earth K2-136 d. XMM-Newton observations of the star indicate an X-ray luminosity of $(1.18\pm0.1)\times10^{28}$ erg s$^{-1}$ in the range 0.15 to 2.4 keV, resulting in an activity of $L_X/L_{\rm bol}=(1.80\pm0.68)\times10^{-5}$. The evaporation past of the planets were modelled using the XUV stellar tracks by Johnstone et al. (2021), the energy limited mass loss formulation (Lecavelier des Etangs et al. 2007, Erkaev et al. 2007), and the thermal evolution formulation by Lopez&Fortney (2014). Our results suggest that planets b and d are most probably purely rocky and have been stripped of their envelopes. On the other hand, planet c likely still has an envelope consisting of 0.8% of the total mass, with its core being relatively large (2.3 $R_{\rm E}$).


INTRODUCTION
The Kepler mission has uncovered a surprisingly large population of planets between Earth and Neptune in size. Their radii have been observed to follow a bimodal distribution (Fulton et al., 2017), with twin peaks at 1.3 and 2.4 E , and a gap at 1.8 E . Whilst the first group is populated by dense rocky planets, the second one can be explained by the additional presence of H/He gaseous envelopes that can double a planet's radius and yet comprise less than a percent of its mass. This phenomenon can be observed as a valley on the radius-period parameter space in Figure 1 . There is evidence that points to stellar high-energy radiation being the main driver for this behaviour (e.g. Owen, 2019;Owen & Wu, 2013), although other channels such as core-powered mass loss have also been suggested (e.g. Gupta & Schlichting, 2020).
The X-ray and extreme ultraviolet radiation (together XUV) that stars emit is absorbed by the upper layers of the atmospheres of close-in exoplanets. This heats the upper atmosphere driving a hydrodynamic wind that overflows the Roche lobe and escapes the planet. Gradual mass loss sculpts the observed planet populations, particularly during the first 100 Myr when the XUV emission is the strongest, stripping the envelopes from planets that are too close or too light. There is also evidence that EUV radiation exerts significant influence on envelopes on much-longer Gyr timescales (King & Wheatley, 2021).
Observations, however, only provide a snapshot of a planet's current evolutionary state. One can attempt to model its current physical properties but its history can be hard to constrain due to uncertainty in the stellar rotation and hence X-ray luminosity history. A straightforward solution to break the degeneracy is to study multiplanetary systems: where all the planets in the system share the same the XUV irradiation history. Their current states can then be used to constrain each other's pasts, especially if the planets are found on both sides of the evaporation valley (Owen & Campos Estrada, 2020).

THE TARGET
K2-136 is a K-type star located 59 parsecs away (as determined by Gaia EDR3). It is a member of the Hyades cluster, making it around 600 Myr old (Perryman et al., 1998), although more recent estimates go as high as 800 Myr (Brandt & Huang, 2015;David & Hillenbrand, 2015). In this work, we adopt an age of 700±100 Myr. The star also hosts a three-planet system: the Earth-sized K2-136 b, the mini-Neptune K2-136 c, and the super-Earth K2-136 d, which are all situated within 0.2 AU of their host star (Ciardi et al., 2018;Livingston et al., 2018;Mann et al., 2018). Ciardi et al. (2018) also reports the existence of a candidate stellar-mass companion, a M7/8V dwarf, at a projected separation of 40 AU, although follow-up observations are still needed to confirm the stars are gravitationally bound to each other.

Planetary system
The three discovery papers of the planetary system (Ciardi et al., 2018;Livingston et al., 2018;Mann et al., 2018) agree within 1 in their estimates for the periods and radii of the three planets. These are = 1.05 ± 0.16 R E and = 7.975 days for planet b, = 3.14 ± 0.36 R E and = 17.307 days for planet c, and = 1.55 ± 0.24 R E and = 25.575 days for planet d. All periods have uncertainties of less than a percent. These are put in context with the evaporation valley on Figure 1 . It is interesting how the middle planet (c) is the one that remains above the valley whereas the one furthest away is not. This points to planet c having a much heavier core than d in order to avoid having its envelope stripped.
The planet masses have also been recently measured by Mayo, Dressing, Buchhave, & The Harps-N Collaboration (2021); with a direct estimate for the mini-Neptune (15.9 ± 2.4 M E ) and upper limits for the other two (< 2.67 and < 6.47 M E for b and d, respectively). We make use of the mass-radius relationships by Otegi, Bouchy, & Helled (2020) to obtain mass estimates for b and d, and to put planet c in context ( Figure  2 ). The small radii of planets b and d place them below the evaporation valley, probably with purely solid cores with no envelope. Their masses are estimated to be 1.07 +0.67 −0.46 M E and 4.09 +2.63 −1.80 M E respectively. In contrast, the mass of planet c indicates a large core, and its location above the valley points to the presence of an envelope. As its envelope likely contributes less than a percent of its mass, we fit its mass to the rocky relation and obtain an estimate for the core size of , = 2.30±0.13 R E , hence its envelope would comprise 25% of the total radius.

The Hyades
The Hyades is a well studied nearby open cluster, with several extensive X-ray studies (Collura et al., 1993;Micela et al., 1988;Stern, Zolcinski, Antiochos, & Underwood, 1981) as well as rotational surveys (Delorme et al., 2011;Douglas et al., 2019). Freund, Robrade, Schneider, & Schmitt (2020) compiled 1100 Hyades members from several studies, both from the core of the cluster and its recently discovered tidal tails. About 130 of these cluster members have measurements for both the rotation periods and the X-ray luminosity. The data presents a tight rotation-spectral type relation for F, G, and K-dwarfs ( Figure 3 , left panel), with M-dwarfs being much more scattered. According to Johnstone et al. (2020), M-dwarfs retain their initial distribution of rotation rates for longer before settling into a single mass-dependent value, joining the tight relation of F,G, and K-type stars. K2-136 fits well on the tight relation. The stray K-type stars at faster rotation rates below the tight relation are likely close binaries (and indeed we see they are offset from the main sequence in an H-R diagram).

OBSERVATIONS
We observed the target with the XMM-Newton telescope in November 2018 for 43 kiloseconds. Data reduction was performed on the EPIC-pn data, which has an operative energy range of 0.15 to 12 keV. The spectrum was modelled with a TBABS model (Wilms, Allen, & McCray, 2000) to account for interstellar absorption as well as an APEC model (Smith, Brickhouse, Liedahl, & Raymond, 2001) for the emission spectrum from a collisionally-ionized diffuse plasma. Solar abundances by Asplund, Grevesse, Sauval, & Scott (2009) were used. The hydrogen column density was set to = 10 18 cm −2 , as estimated by Redfield & Linsky (2001) for the Hyades. Furthermore, the ROSAT-range flux was estimated by

DATA ANALYSIS
Since the source is only detected up to 2.4 keV, only the range 0.15 to 2.4 keV was used for model fitting. Due to the low number of counts per bin, we use the C statistic as the fit parameter (Cash, 1979). We find that a three-temperature APEC model works best for our data, with temperatures 0.14 +0.05 −0.06 , 0.29 +0.17 −0.06 , and 0.97 +0.17 −0.14 keV. The spectrum is plotted in Figure 4 . Integrating the model yields a flux of = 2.83 ± 0.25 × 10 −14 erg cm −2 s −1 , resulting in a luminosity of 1.18 ± 0.1 × 10 28 erg s −1 in the range 0.15 to 2.4 keV. Extrapolating down to an energy of 0.1 keV yields a ROSAT-range luminosity of 1.32 ± 0.12 × 10 28 erg s −1 .
Some of this X-ray flux is likely to arise from the closeby M-dwarf companion. We estimate this from the measured 2MASS fluxes of the candidate companion, which Ciardi et al. (2018) find to be = 14.1 ± 0.1, = 13.47 ± 0.04, and = 13.03 ± 0.03. Interpolating with the stellar quantity sequences by Pecaut & Mamajek (2013), we estimate its bolometric luminosity to be bol = 2.00 ± 0.35 × 10 30 erg s −1 , about 5 × 10 −4 times lower than that of K2-136. Since Wright, Drake, Mamajek, & Henry (2011) finds that stars tend to be saturated at an X-ray activity of 10 = −3.13 ± 0.08, and we assume that such a low mass star remains X-ray saturated even after 800 Myr (Johnstone et al., 2020), we can estimate its X-ray luminosity to be ≈ 1.48 ± 0.38 × 10 27 erg s −1 , about 10% of the flux from K2-136 measured in this work, which is comparable to the uncertainty in our measurement.
Finally, we estimate the EUV flux using the relations by King et al. (2018) and thus obtain a value of EUV = 2.38 ± 0.17×10 28 erg s −1 in the range 0.0136 to 0.15 keV. This results in a combined XUV luminosity of XUV = 3.7 ± 0.2 × 10 28 erg s −1 , three times stronger than X-rays alone.
Plotting the star on the rotation-activity relation with its cluster siblings (Figure 3 , right panel) indicates good agreement with the other Hyades K-dwarfs. It is also worth noting the large scatter in activity around the fit, up to one order of magnitude. It is yet to be explained whether this originates from random variability or intrinsic tracks for each star.

EVOLUTION MODELLING
In order to explore the pasts of the planetary envelopes, three ingredients are necessary: the XUV emission track of the host star as it spins down, a formulation for the mass loss that updates the envelope mass on each time step, and a description of the planetary interior that relates the envelope mass fraction to its size, vital for recalculating the envelope radius after mass is lost on each time step.

Stellar tracks
The XUV emission of a star is tied to its Rossby number, the ratio between rotation period and convective turnover time (e.g. Wright et al., 2011). Stars spin down with age, which results in a decrease of the star's emission with time, eventually diminishing its capability to evaporate atmospheres past a few Gyr.
We make use of the physically motivated X-ray tracks proposed by Johnstone et al. (2020), which model the radiative core and convective envelope of F, G, and K-type stars as two independently rotating solid components. Each component has certain angular momentum affected by torques of different origins (stellar wind, momentum exchange, and PMS spin-up), which can be used to predict the star's rotation and XUV dependence on age. These tracks are plotted in Figure 5 (left panel) for a star of the mass of K2-136 (0.7 M ⊙ ) and a distribution of initial rotation rates. K2-136 seems to fit well with the average rotation track (50th percentile).

Thermal evolution
We adopt the envelope structure model described by Lopez, Fortney, & Miller (2012), which models the cooling and contraction of adiabatic envelopes assuming a hot-start model.

FIGURE 3
Left: rotation period of Hyades members versus colour index (Gaia blue and red magnitudes, BP-RP). Right: Xray luminosity against Rossby number for Hyades members. The rotation-activity relation by Wright et al. (2018) is plotted in green and the relation for Hyades members by Freund et al. (2020) in magenta. Data have been taken from Freund et al. (2020). K2-136 is plotted as a black star in both panels. The planet's energy budget is tracked by taking into account heat from different sources, such as stellar radiation, radioactive decay from heavy elements, and thermal inertia from the isothermal rocky core. Although they begin their models at 1 Myr, the hot-start assumption might result in over-inflated radii (>10 R E ) for low-mass planets at early ages. Despite that, these envelopes cool down and contract quickly and should be insensitive to the choice of initial entropy by the age of 10 Myr. Computed model tables and polynomial fits are provided by Lopez & Fortney (2014), which only run from 100 Myr to 10 Gyr, however. Regardless, we choose a lower age limit of 10 Myr in our simulations as most of the mass loss occurs before 100 Myr.

Mass loss
The photoevaporation is modelled using the energy-limited approach (Erkaev et al., 2007;Lecavelier Des Etangs, 2007;Watson, Donahue, & Walker, 1981). This balances the incident energy flux from the star against the energy necessary to displace matter from the planet's upper atmosphere to its Roche lobe. The mass loss rate is given bẏ where and are the mass and radius of the planet, respectively, XUV is the incident combined X-ray and EUV flux into the planet, is the ratio between the radii at which the atmosphere becomes optically thick for XUV wavelengths and for visible wavelengths ( XUV ∕ ), is the energy efficiency of the process, and accounts for the potential energy difference between the surface and the height of the Roche-lobe (see Erkaev et al. (2007) for more details). This approach is attractively simple, although it does mask more complex physics behind the and parameters.
Generally, constant canonical values for these two parameters are used, but there more sophisticated formulations as well. For instance, Salz, Schneider, Czesla, & Schmitt (2016) perform hydrodynamic simulations to derive scaling laws for the efficiency and the XUV radius based on both the gravitational potential of a planet and its irradiation level. In this work, we adopt a canonical efficiency of = 15% and a XUV radius ratio = 1 (King et al., 2018), and we compare with the formulation by Salz et al. (2016) for .

FIGURE 5
Stellar tracks from the model by Johnstone et al. (2020). Left: possible rotation tracks for a 0.7 ⊙ star given a range of initial rotation rates; 5th percentile (blue), 50th percentile (red), and 95th percentile (green). The track that best fits K2-136 is the solid black line, with the uncertainties set by the two dotted black lines. Right: X-ray activity = ∕ bol (black line). The X-ray (red) and EUV (blue) tracks are plotted separately as well. K2-136 is plotted as a magenta star.

Simulation results
Planets c and d were de-evolved from their current states back to 10 Myr. This assumes that planet d just finished losing its now-missing envelope, which sets an upper limit for the total mass lost. Using an efficiency of = 15% and a XUV radius ratio of = 1, at 10 Myr we obtain envelope mass fraction of 0.83% for planet c, expanding from 3.14 to ∼3.8 R E , and 0.12% for planet d, inflating its radius from 1.55 to ∼2.15 R E .
We also evolved forwards two synthetic tracks for a planet with a core identical to planet b and envelope mass fractions equal to those of planets c and d at 10 Myr (as computed above). These tracks are evolved from the age of 10 Myr up to 700 Myr. The results show that, whichever envelope planet b had, it was likely lost within 5 to 10 Myr. This is expected, as b is the smallest and closest planet, thus more susceptible to photoevaporation. Planet b is therefore less useful in constraining the evaporation pasts of the system. Our evolutionary tracks are shown in Figure 6 (top panel).
We also compare these results with tracks generated using an XUV radius that follows the relation by Salz et al. (2016) (Figure 6 , bottom panel). In both cases, an efficiency of 15% is used. We note a mass loss enhancement 2 to 5 times greater. This results in planets c & d having a remarkably similar envelope mass fraction of about ∼0.9% at 10 Myr.

DISCUSSION & CONCLUSION
Using XMM-Newton measurements, we estimate the X-ray luminosity of K-type star K2-136 to be 1.18 ± 0.1 × 10 28 erg s −1 in the range 0.15 to 2.4 keV. The star, with an age of 700 ± 100 Myr, is about 50% X-ray fainter than the luminosity FIGURE 6 Top: planet radius against age. Planets c (green line) and d (blue) were de-evolved from their current age of 700 Myr. Using their envelope mass fractions at 10 Myr, planet b was evolved forward to 700 Myr, with initial fractions of ∕ = 0.83% (red dotted line) and ∕ = 0.12% (magenta dashed line). The envelope loss timescales for planet b were estimated to be ∼2 and ∼4 Myr. Bottom: comparison between an extended radius of = 1 and the formulation by Salz et al. (2016), ran using the tracks for planets c (green) and d (blue). An efficiency of = 15% was used.
predicted by Johnstone et al. (2020) in their stellar XUV tracks (see Figure 5 ). The discrepancy is not unexpected, however, given the relatively large scatter of X-ray luminosities in relation to age (Jackson, Davis, & Wheatley, 2012;Johnstone et al., 2020).
We modelled the evaporation pasts of the three planets in the system using the stellar tracks by Johnstone et al. (2020), the energy-limited mass loss formulation (Erkaev et al., 2007;Lecavelier Des Etangs, 2007), and the thermal evolution model by Lopez & Fortney (2014). We find that planets b and d are most likely rocky, given their current density and size, and planet c possesses a gaseous envelope that consist of 0.8% of its total mass. Using the extended radius formulation by Salz et al. (2016), we obtain similar initial envelope mass fractions for planets c & d of 0.9% at 10 Myr. Rogers & Owen (2021) finds the distribution of envelope mass fractions to be strongly peaked at 4% for their sample of exoplanets, which might be attained by evolving our planets further back in time to 1 Myr with a alternate thermal evolution formulation.
Further mass loss might be attained by taking into account additional mass loss sources, such as stellar wind and Coronal Mass Ejections (CMEs). The degree to which these phenomena influence the atmospheres of mini-Neptunes is still under study. Kislyakova et al. (2014) found it to be several times lower than thermal mass loss rates. Fast CMEs, however, could have a greater influence on planets that are extremely closein (<0.1 AU), or lack a significant magnetosphere (Lammer et al., 2009).