EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 646 (February 2021) A&A, 646 (2021) A15 Full HTML
Open Access
Issue
A&A
Volume 646, February 2021
Article Number A15
Number of page(s) 20
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/202039072
Published online 02 February 2021

© D. Blain et al. 2021

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

While fairly common among the thousands of exoplanets discovered to date1 (see Fig. 1), the atmospheres of super-Earths and sub-Neptunes – with masses between 2 and 10 M⊕ – are poorly understood. By essentially using the mass and radius of the planets, it has been well-established that planets with low masses (<2 M⊕) must primarily be made up of iron and silicates, and generally with a thin atmosphere. On the other hand, planets with high masses (>10 M⊕) must retain a thick atmosphere composed mainly of H2 and He, representing a significant portion of the planet mass (Chen & Kipping 2016). Within the transition between the two populations, however, studies with models struggle to give a clear answer (Valencia et al. 2013; Fulton et al. 2017; Zeng et al. 2019; Otegi et al. 2020), exhibiting a degeneracy between massive rocky planets, ocean planets, and small gaseous planets. Thus, spectral observations are crucial for characterising the atmosphere of these objects and to better constrain their internal composition.

The recently discovered, transiting exoplanet K2-18b (see Table 1) offers not only the opportunity to retrieve spectroscopic data on the atmosphere of a sub-Neptune, but also to study such atmospheres under nearly Earth-like conditions. Indeed, the stellar irradiance received by K2-18b ( W m−2) is very close to that of the Earth. Nine transits were acquired using the Wide Field Camera 3 on the Hubble Space Telescope (HST/WFC3), as well as data from Kepler K2 and Spitzer IRAC channels 1 and 2 (Benneke et al. 2019). These data have already been analysed by several teams (Benneke et al. 2019; Tsiaras et al. 2019; Madhusudhan et al. 2020; Scheucher et al. 2020; Bézard et al. 2020). Also, Scheucher et al. (2020) and Bézard et al. (2020) used self-consistent models (respectively, “1D-TERRA” and “Exo-REM”) to analyse the data, while the other used free retrieval algorithms. The first five investigations concluded that there is a presence of H2 O as well as a significant amount of H2 –He. They also derived upper limits for the abundance of CH4, which in all cases was lower than ≈ 3.5% at a 99% confidence level. In contrast, Bézard et al. (2020) found that the HST/WFC3 spectrum is dominated by CH4 absorption and found abundances of CH4 and H2O, respectively, between 3% and 10% and 5% and11% at a 1 − σ confidence level. They assumed an H2-dominated atmosphere and varied the atmospheric metallicity, but their simulations did not include H2 O clouds and they did not simulate non-solar C/O ratios.

Given the irradiance of the planet and the presence of H2 O in the atmosphere, the question of the existence of liquid H2O is naturally posed. Benneke et al. (2019) found that a cloud layer was needed to reproduce the data. They then used a self-consistent model and found that liquid H2O could condense at the right pressure to explain this cloud layer. Contrary to Benneke et al. (2019), Madhusudhan et al. (2020) did not find compelling evidence for clouds or hazes in the atmosphere. Using an interior model, Madhusudhan et al. (2020) indicated that if the planet had a small rocky core and a thin H2 /He atmosphere, an ocean of liquid H2O could exist. However, Scheucher et al. (2020) ruled out this possibility, arguing that an H2 O ocean wouldpartially evaporate in the atmosphere, giving a spectrum that would be incompatible with the data.

In the present work, we use Exo-REM, our self-consistent one-dimensional (1D) atmospheric model adapted for transiting exoplanets to study K2-18b atmospheric composition assuming a thick H2 –He atmosphere. We first present the extension of Exo-REM (Baudino et al. 2015, 2017; Charnay et al. 2018) to irradiated planets and detail the calculations of the absorption cross-sections. We then expand on the work of Bézard et al. (2020) by including clouds and studying their formation, and by investigating the effects of irradiation, the eddy diffusion coefficient, internal temperature, and metallicity on the atmosphere properties.

thumbnail Fig. 1

Mass-radius distribution of the 42662 confirmed exoplanets to date (17 June 2020). The scatter plot shows only the 454 planets for which the mass and the radius are known within ± 15%. The mass-radius curves for an Earth-like interior (orange) and a pure water ice planet (blue) are taken from Seager et al. (2007).
Table 1

Parameters of K2-18b and its star.

Table 2

Collision-induced absorption references.

2 Model

Exo-REM is a 1D radiative-equilibrium model first developed for the simulation of young gas giants far from their star and brown dwarfs (Baudino et al. 2015, 2017; Charnay et al. 2018). Fluxes are calculated using the two-stream approximation assuming hemispheric closure. The radiative-convective equilibrium is solved assuming that the netflux (radiative + convective) is conservative. The conservation of flux over the pressure grid is solved iteratively using a constrained linear inversion method. We take into account Rayleigh scattering from H2, He, and H2 O, as well as absorption and scattering by clouds – calculated from the extinction coefficient, single scattering albedo, and asymmetry factor interpolated from pre-computed tables for a set of wavelengths and particle radii (see Charnay et al. 2018).

We made several upgrades to this model in order to extend its simulation capabilities to high-metallicity (up to ≈ 1000 times the solar metallicity) and irradiated planets, namely sub-Neptunes and moderately hot Jupiters (equilibrium temperature lower than 2000 K). These upgrades are detailed in the following sections. The Exo-REM source code and the k-coefficients we used (see Sect. 2.1.3) are available online3.

2.1 Spectroscopic data

2.1.1 Collision-induced absorptions

The collision-induced absorptions (CIA) included in our simulations and their references are given in Table 2. We added the contributions of the H2 O–H2O CIA and, since spectroscopic data for the H2 O–H2 CIA were not available, of the H2O–air CIA. The reasons for this addition are detailed in Appendix A.

2.1.2 Absorption cross-sections

We calculated the absorption cross-sections of CH4, CO, CO2, FeH, H2O, H2 S, HCN, K, Na, NH3, PH3, TiO, and VO at 25 pressure levels equally spaced in the log-space between 0.1 and 107 Pa. At each of these pressure levels, we calculated the cross-sections at temperatures 100, 150, 200, 250, 300, 400, 500, 600, 800, 1000, 1200, 1500, 2000, 2500, and 3000 K – except for NH3, the reasons for which are detailed below. We calculated line absorption up to a given distance (Δν) from the line center, using the same procedure as described in Baudino et al. (2015). We used a sub-Lorentzian line profile with a χ factor, based on Burch et al. (1969) and Hartmann et al. (2002), for, respectively, CO2 and then all the other species. The cross-sections were calculated over the wavenumber range displayed in Table 3, at a resolution that is similar to the line half width.

When line lists of individual isotopes were available, we merged them by multiplying line intensities in order to reproduce the isotopic ratio found for Jupiter. Otherwise, we used the default isotopic ratio given by the database. We used H2 and He pressure-broadened halfwidths (γ) and temperatureexponents (n) whenever they were available. When both were available, we used a mixture of 90% H2 and 10% He, which roughly corresponds to the standard Solar System He/H ratio (Lodders 2019). More details are given in Tables 3 and 4. The specificities of some species are listed below.

Table 3

Species cross-section parameters.

CH4

We included the contribution from 12 CH4, 13 CH4, and CH3D, all taken from the TheoReTS database (Rey et al. 2017), which is more accurate than the ExoMol database. We used a 12 C/13C ratio of 89 (Niemann et al. 1998) and a D/H ratio of 2 × 10−5 (Lellouch et al. 2001). We recall that the TheoReTS line list of CH3D stops at 6500 cm−1. The line list provided by the TheoReTS database separates the “strong lines”, which have an intensity of >10−26 cm molecule−1, from the very-high-spectral-density “weak lines”, which are regrouped in “super lines”. This separation allows for much faster cross-section calculations, but at the price of losing information on the quantum numbers of individual lines. Hence, we had to use the same γ and n for all lines.

FeH, TiO and VO

H2 or He γ and n were not available for these molecules and we used the parameters from CO.

K and Na

Given that these species have very high intensity lines (up to ≈ 10−13 cm molecule−1), we needed to extend our Δν up to 9000 cm−1 in order to correctly account for far wing absorption. Following Burrows et al. (2000), we used a Voigt profile up to a detuning frequency of cm−1 for K and cm−1 for Na. Beyond that detuning frequency, we used the profile described by Baudino et al. (2015).

NH3

Cross-sections were calculated up to 1500 K because the line list we used lacks completeness above this temperature. We took a 14 N/15N ratio of 500 (Furi & Marty 2015).

Table 4

Line broadening references.

PH3

The line list we used is not complete above 1000 K. However, Sousa-Silva et al. (2014) provide the PH3 partition function up to 3000 K. This allows us to use the percentual loss of completeness to estimate the proportion of missing opacity, as suggested in the cited work.

2.1.3 k-coefficients

We used these high-resolution absorption coefficients to calculate k-coefficients according tothe method described in Baudino et al. (2015), using 16 Gauss–Legendre quadrature points (eight for values of the cumulative distribution of the absorption coefficients between 0 and 0.95, and eight between 0.95 and 1). Four sets of k-coefficients were calculated, at resolutions of 0.5, 20, 200, and 2000 cm−1, although only the 20 cm−1 step one was used in the main part of this study.

Table 5

Species included in thermochemical calculations.

2.2 Radiative-convective equilibrium model: stellar irradiance

We added the planetary averaged stellar irradiance E↓, e,ν reaching the planet to Exo-REM. We used stellar spectra from the BT-Settl model (Allard et al. 2012). We chose a spectrum modelled at an effective temperatureTS = 3500 K, log10(g[cm s−2]) = 5, a null metallicity and no alpha enhancement. We neglected the spectral dependency on surface gravity and metallicity. Then we interpolated the BT-Settl spectrum on a spectral grid with a 0.1 cm−1 wavenumber step. This interpolated spectrum was then convolved at the resolution of Exo-REM flux calculations (see Sect. 3) to obtain the modelled radiosity JS, e,ν. Then we obtained E↓, e,ν from: (1)

where JB, e,ν(T) is the radiosity of a black body at temperature, T, and T*, eff is the effective temperature of the star, R* is the radius of the star, ap is the distance between the star and the planet, and the geometric factor 1∕4 is used to represent planet-averaged conditions. The ratio of the JB, e,ν terms is used to obtain a stellar spectrum at the effective temperature of K2-18. This was done instead of interpolating on the BT-Settl grid as TS, eff is close to T*, eff.

2.3 Atmospheric model: thermochemistry

We used the recommended present atomic solar system abundances given by Lodders (2019) to define our standard abundances and our reference metallicity ((Z/H)⊙). The atomic abundances at a metallicity Z/H are obtained by keeping the H, He, Ne, Ar, Kr, and Xe standard abundances constant while multiplying the standard abundances of other elements by Z/H. The list of all the species included in our thermochemical calculations is displayed in Table 5. It is possible to adjust the abundances of all individual elements listed in this table, although it is always assumed that the atmosphere is dominated by H2, limiting Exo-REM capabilities to metallicities ≾1000 (Z/H)⊙. We included species that either: significantly affect the volume mixing ratio (VMR) profiles of our absorbing species (see Sect. 2.1.2), are susceptible to the formation of a major cloud layer, or substantially impact the molar mass of the atmosphere (according to Lodders 2010).

Most of the equilibrium abundances are derived from the equation of conservation of each species, using the standard Gibbs free energy of formation listed in Chase (1998). For PH3(g), we used the from Lodders (1999), while for MnS(cr,l) and ZnS(cr,l), we used the values from Robie & Hemingway (1995), and for CaTiO3(cr) we used values from Woodfield et al. (1999). We used a grid of temperatures between 200 and 4000 K, with a 100 K step, to map our . When the temperaturerange of the references used were smaller than our target temperature range, we used a linear extrapolation to fill up our grid. The at a given temperature are linearly interpolated from this grid. We calculated the saturation pressure of the following species directly: H2 O (Wagner & Pruss 1993; Wagner et al. 1994; Lin et al. 2004; Fray & Schmitt 2009), NH3 (Fray & Schmitt 2009; Lide 2009), and NH4SH (Stull 1947). We used a very simplified Ca–Al–Ti chemistry compared to, for instance, Lodders (2002), but the impact on our simulations was expected to be limited. We simulated the formation of Fe–Ni alloys by treating the Fe independently while keeping the Ni(g)/Fe(g) ratio constant. We considered the formation of H3 PO4 instead of P4 O6, based on the results of Wang et al. (2016).

We included non-equilibrium processes in the CH4–CO, CO–CO2, N2 –NH3–HCN chemical systems following the approach of Zahnle & Marley (2014). The chemical quenching level is determined by equating a reaction timescale to the mixing time H2∕Kzz, where H is the atmospheric scale height and Kzz the eddy mixing coefficient. Below the quenching level, the abundance profiles of the relevant species are governed by thermochemical equilibrium while, above this level, the mixing ratios of the quenched species are held constant. We used the chemical reaction timescales given by Eqs. (12)–(14) for CH4–CO, Eq. (44) for CO–CO2, Eq. (32) for NH3–N2, and Eq. (40) for HCN–NH3–N2 in Zahnle & Marley (2014). In the case of PH3, we assume that its conversion to H3 PO4 is inhibited, as observed in the giant planets of our solar system. We compared our chemical model in the case of K2-18b (taking identical temperature profile, elemental abundances, and eddy diffusion coefficient) against the model described in Venot et al. (2020)4. For CO, HCN, and NH3, we found VMR differences lower than 20% with our model in the HST sensitivity region, and lower than 2% in the case of CH4 and H2 O. The difference is larger for CO2, with our VMR being 25% higher than the value found in O. Venot’s model. We consider these differences to be satisfactory.

Table 6

Model grid.

3 Methodology

The star and planetary parameters used in our simulations are displayed in Table 1. A summary of our model grid is displayed in Table 6. We took a planetary radius at 105 Pa of Rp = 16 400 km, which is slightly different from the value used by Benneke et al. (2019; who define Rp as the radius at 103 Pa), because it is closer to the radius we find when fitting the observed data. To compare our model spectrum with the data, we applied an offset on Rp in the calculation of the transit depth, such that the χ2 is minimised.

The internal temperature resulting from the residual heat of formation of the planet was calculated following (Rogers & Seager 2010): (2)

where Lint is obtained from (Rogers & Seager 2010): (3)

where a1 = −12.46 ± 0.05, , , , tp is the age of the planet, and the astronomical constants are defined by the International Astronomical Union (IAU, Mamajek et al. 2015). According to a gyrochronological model from Guinan & Engle (2019), the age of K2-18 can be estimated at 2.4 ± 0.6 Gyr. Assuming that K2-18b ended its formation a few Myr after the formation of its star (so that t* ≈ tp), like what happened in our solar system, we obtained from Eq. (2) K. Rounding down this value, we chose 80 K as our nominal internal temperature.

The net fluxes are calculated from 40 to 30 000 cm−1, with a step of 20 cm−1. The atmospheric grid consists of 81 levels equally spaced in the log-space between 0.1 and 107 Pa. We imposeda correlation length of 0.5 pressure scale height to the solution temperature profile in order to avoid non-physical oscillations (see Baudino et al. 2015). We considered cloud radiative effects and scattering only from H2 O clouds. The NH3 and NH4SH clouds would form in atmospheres that are colder than what is expected for K2-18b, while NH4Cl clouds are too thin to have a significant impact, and other clouds are condensing too deeply into the atmosphere (see Sect. 4.2). Cloud vertical mass mixing ratios are calculated assuming equilibrium between the downward flux of falling particles and the upward flux of gas and cloud particles due to advection and turbulent mixing (see Charnay et al. 2018).

When we change the irradiance of the planet by a factor of k, we do so by multiplying ap by (see Eq. (1)). We proceed this way to make model comparison easier. Indeed, ap impacts only the irradiance, while the two other parameters, T*, eff and R*, affects respectively the stellar spectrum and both the stellar spectrum and the transmission spectrum.

We estimated Kzz using the values derived by Exo-REM (see Charnay et al. 2018). We found that typical Kzz values for K2-18b range from 106 to 109 cm2 s−1, with the highest values found in the convective layers. Given the uncertainties on the estimation of Kzz, we enlarged this range to 105 –1010 cm2 s−1. We also found that a quenching of our species often occurs just below the uppermost convective layer (see Sect. 4.3), hence, we set our nominal Kzz value at 106 cm2 s−1.

We used a fractional area covered by clouds fc = 0.15 for the calculation of the temperature profile, and fc = 1.0 for the calculation of the transmission spectra. To ensure numerical stability, the cloud mean particle radius, r, was fixed at a constant value (50 μm) that roughly corresponds to the maximum one predicted by the Exo-REM self-consistent cloud model at 300 (Z/H)⊙ and nominal irradiation (see Charnay et al. 2018, and Sect. 4.2), unless stated otherwise. We call this our nominal model. We explore the effects of changing some of these parameters in the following sections.

thumbnail Fig. 2

Best fit of our model to the dataset of Benneke et al. (2019) at nominal irradiation (1368 W m−2). Black: K2, HST and Spitzer data from Benneke et al. (2019). Grey: HST data from Tsiaras et al. (2019). Red: Exo-REM transmission spectrum at 175 (Z/H)⊙, with an offset of the 105-Pa level of −5 km. The χ2 of this spectrum against Benneke et al. (2019) data is indicated in parentheses.

4 Results

4.1 Metallicity and irradiance

Here, we assume that K2-18b has retained a relatively thick H2 –He atmosphere. This atmosphere can be enriched in heavy elements – compared to the initial proto-stellar nebula – during the formation of the planet via collisions with planetesimals (Fortney et al. 2013). There could also be a contribution from the erosion of the primordial core (Iaroslavitz & Podolak 2007).

In Fig. 2, we show our best fit for the dataset of Benneke et al. (2019) at nominal irradiation, along with the HST data reduced by Tsiaras et al. (2019). The contributions of the different opacity sources to this spectrum are represented in Fig. 3. The temperature profile and VMR of this spectrum are represented respectively in Figs. 4 and 5. Because the stellar spectrum we used is slightly more intense than a black body between ≈1 and 2.5 μm, the temperatures we obtain in the upper atmosphere are slightly higher (≈5 K at 1 kPa) than what would be obtained with a black body at the effective temperature of the star. The prevalence of CH4 absorptions over H2O absorptions is discussed in Sect. 5.1.

In Fig. 6, we display the χ2 of our nominal models for the datasets from Benneke et al. (2019) and Tsiaras et al. (2019), including their respective reduction of the same HST raw data, as well as in both cases the K2 and Spitzer data from Benneke et al. (2019). We performed simulations of K2-18b for 1, 3, 10, 30, 50, 75, 100, 125, 150, 175, 200, 300, 400, 500, and 1000 times (Z/H)⊙ and for irradiations between 0.5 and 1.5 time the nominal irradiation, with a step of 0.1. We consider a model as statistically accurate if it can reproduce the data within the 1σ confidence level (68%). Since there are 20 data points, we have 20 to 1 degrees of freedom, so the 1σ confidence level corresponds to a χ2 of 21.36. While the dataset from Tsiaras et al. (2019) allows, according to our interpolation, for metallicities ≥ 65 (Z/H)⊙ within K2-18b irradiation 3σ uncertainties, the dataset from Benneke et al. (2019) is much more restrictive, allowing only a metallicity between 100 and 200 (Z/H)⊙ at 0.9 times the nominal irradiation, or a metallicity between 150 and 200 (Z/H)⊙ at or above the nominal irradiation. These values are within the most common range of metallicity predicted by Fortney et al. (2013) for planets with radius in the 2–4 R⊕ range (between 100 and 400 + (Z/H)⊙). At nominal irradiation, our best fit against Benneke et al. (2019) and Tsiaras et al. (2019) data is located respectively at 175 (Z/H)⊙ (χ2 = 21.07) and 150 (Z/H)⊙ (χ2 = 17.24). In both cases, there is no H2O cloud formation, and we derive a Bond albedo of respectively 0.017 and 0.018.

The relationship between temperature, metallicity, and goodness of fit presented in Fig. 6 can be explained as follows. Theamplitudes of the features in the transmission spectra are correlated with the abundance of absorbers in the atmosphere as well as to its scale height, which can be written as (4)

where R is the gas constant, and T, μ and, g are respectively the temperature, the molar mass and the gravity in an atmospheric layer. As the metallicity decreases, two concurrent effects are occurring. The diminution of the VMR of absorbers in the atmosphere will of course decrease thespectrum amplitude. However, decreasing the VMR of heavy species will also leads to a decrease of μ, and thus to an increase of the amplitudes. This latter effect is dominant in our simulations for metallicities ≥ 10 (Z/H)⊙, as illustrated in Fig. 7. At lower metallicities, heavy elements contribute to less than 10% of μ, so the effect of less intense gas absorptions becomes dominant, flattening the transmission spectrum. If we decrease irradiation, the atmosphere gets colder and H2O clouds start to form and to thicken, flattening the spectrum at transparent wavelengths where part of the transmitted stellar light reaches the condensation level. As stellar irradiation decreases, the cloud forms lower and lower in the atmosphere, eventually ending below the region probed by transit spectroscopy where it can no longer directly affect the spectrum. The condensation also removes more gaseous H2 O, decreasing μ and H2O absorption. On the other hand, increasing the irradiation leads to an increase of the temperatures, and, hence, of the scale height and of the absorption amplitudes, as shown in Fig. 8.

thumbnail Fig. 3

Contributions of the absorbing species to our best-fit transmission spectrum at nominal irradiation within the K2 to Spitzer spectral range. The H2O cloud is not forming in this case, so there is no cloud contribution. The spectral contribution of individual species takes the CIA and Rayleigh scattering (represented as a dotted curve) into account. The spectral contribution of FeH, TiO, and VO are not represented here as they are insignificant.
thumbnail Fig. 4

Temperature profile of our best fit atmospheric model against the dataset from Bennekeet al. (2019) at nominal irradiation. Solid line: temperature profile. Dotted lines: condensation profiles of selected species. Red line: convective layers. Blue dot: H2 O ice-Ih–liquid–gas triple point.

4.2 Clouds

In Fig. 9, we show that H2O cloud can form at nominal irradiation and metallicities above 200 (Z/H)⊙. At 300 (Z/H)⊙, the cloud forms at 1 kPa. However, this cloud cannot form even at 1000 (Z/H)⊙ if the irradiation is more than 10% higher than the nominal irradiation. As a general tendency, the H2O cloud forms lower in the atmosphere as the irradiation decreases and the metallicity increases. This is because as the irradiation decreases, the temperature decreases, so the temperature profile crosses the H2 O condensation profile at higher pressures. Also, as the metallicity increases, the partial pressure of H2 O increases, shifting its condensation profile towards higher temperatures.

The effect of H2O clouds at nominal irradiation on the temperature profile is insignificant, but its effect on the transmission spectrum is quite important, as shown in Fig. 10. The cloud shields the transit spectrum from the star light passing below its altitude, with a particularly strong effect in the visible range. The cloud also has a slight effect on the Bond albedo, as shown in Fig. 11. Under 200 (Z/H)⊙, increasing the metallicity decreases the Bond albedo, because there is more and more gas absorption. At 300 (Z/H)⊙, adding clouds increases the Bond albedo from 0.017 to 0.018. This effect then increases with metallicity, because the cloud thickens and, thus, it is able to reflect more light. We note that above 300 (Z/H)⊙, the Bond albedo increases with metallicity even without clouds. This is due to the removal of gaseous H2 O from the atmosphere.

We were unable to make Exo-REM simulations converge at nominal irradiation when using its self-consistent cloud mode (Charnay et al. 2018). In this mode, the mean particle radius is determined from a comparison of the timescales of the microphysical processes governing the formation and growth of cloud particles. This is likely due to the temperature profile barely crossing the H2 O saturation profile, so that a slight variation of temperature in the solution profile can have a major effect on cloud formation. Nevertheless, we were able to determine that the H2 O cloud particles, according to Exo-REM, should be between 20 and 50 μm at 300 (Z/H)⊙, for a calculated Kzz of ≈ 9 × 107 cm2 s−1 at the layer of condensation. The result is essentially the same if we fix the sedimentation parameter (fsed) between 1 and 5, which is typical of the clouds in the solar system (see Charnay et al. 2018). If we fix Kzz at 1010 cm2 s−1 and fsed at 5, we obtain a maximum mean cloud particle radius of 600 μm. In Fig. 12, we show the effect of the mean cloud particle radius on the transmission spectrum. As the mean radius increases, the effect of the cloud is less and less visible. At 600 μm, it is almost indistinguishable from the spectrum without clouds.

We found that in our simulations, liquid H2O clouds can form at irradiations lower than 80% of the nominal irradiation, which is slightly lower than the 3σ lower uncertainty of K2-18b irradiation. However, according to Fig. 6, this case is not favoured by the data, and is above the 1σ confidence level against Benneke et al. (2019) data. It is not surprising that our results for the likelihood of liquid H2 O differ from Benneke et al. (2019), even though they also used a self-consistent model in their demonstration. Indeed, they assumed an albedo of 0.3, which is much higher than what we found. They probably also included a much lower amount of CH4, which significantly reduces the stellar heating.

We note from Fig. 4 that aside from H2 O, clouds of NH4Cl, KCl, ZnS and Na2S are condensing within our pressure grid. Clouds below KCl (included) form too low in the atmosphere to significantly impact the temperature profile above the convective layers, so they can safely be ignored. This is not the case for NH4Cl. According to our simulations, NH4Cl removes ≈ 30% of the totalamount of Cl in the atmosphere, with the remaining 70% being removed mainly via the condensation of KCl, RbCl, and CsCl. The latter two are not included in Exo-REM, the condensation of RbCl and CsCl removing less than 0.15% of the Cl, assuming a standard composition. While the NH4Cl cloud is, to our knowledge, rarely mentioned in the exoplanet literature, it offers an explanation to the lack of Cl-bearing species in the upper atmosphere of the giant planets of our solar system (Fouchet et al. 2004; Teanby et al. 2006). Nevertheless, we found that with our nominal model, NH4Cl clouds with r = 5 μm – roughly the minimum value found by the Exo-REM self-consistent cloud model on the layer with the most particles – have essentially no impact on the temperature profile and a marginal impact on the visual part of the transmission spectrum5, lowering the χ2 against Benneke et al. (2019) data by ≈0.2. Therefore, they can probably be neglected for K2-18b, but they might play a more important role in the transmission spectra of slightly hotter planets.

thumbnail Fig. 5

VMR of our best fit atmospheric model against the dataset from Benneke et al. (2019) at nominal irradiation (175 (Z/H)⊙, Kzz = 106 cm2 s−1). For clarity, only absorbing species and N2 are represented. FeH, TiO, and VO are not represented as their respective abundances are, respectively, < 10−12, < 10−35, and < 10−36 within our pressure grid. The bump of the PH3 VMR at 200 kPa is due to the chemical equilibrium between PH3 and P2 peaking in favour of P2 at this pressure.
thumbnail Fig. 6

Goodness of fit of our nominal models, varying only irradiation and metallicity. Left: data from Benneke et al. (2019), including HST, K2, and Spitzer data. Right: data from Tsiaras et al. (2019), including HST, K2, and Spitzer data. The white colour corresponds to the value of χ2 of 21.36, indicating a fit of the datasets at the 1σ confidence level. Accordingly, the blue colour indicates overfit and the red colour underfit of the datasets. The solid lines represents the 1 and 2σ confidencelevels of our fits. The dashed and dotted lines represents, respectively, the 1 and 3σ uncertainties on K2-18b irradiation.
thumbnail Fig. 7

Effect of metallicity on the transmission spectrum at nominal irradiation. The corresponding value of χ2 against Benneke et al. (2019) data is indicated in parentheses.
thumbnail Fig. 8

Effect of irradiation on the transmission spectrum at 10 (Z/H)⊙. The corresponding value of χ2 against Benneke et al. (2019) data is indicated in parentheses.
thumbnail Fig. 9

H2O condensation pressure level of our nominal models, varying only irradiation and metallicity. White squares: no cloud formation. Strips: simulation where the H2 O clouds are optically thick (normal optical depth τ > 1 at 1 μm). Dots: simulations where H2O condenses into its liquid phase. Dotted line: 3σ lower limit of K2-18b irradiation.
thumbnail Fig. 10

H2O cloud and metallicity effect on the transmission spectrum, at nominal irradiation. The χ2 against Benneke et al. (2019) data is indicated in parentheses.
thumbnail Fig. 11

H2O cloud and metallicity effect on the Bond albedo of K2-18b, at nominal irradiation. Cloud condensation occurs only above 200 (Z/H)⊙.

4.3 Eddy diffusion coefficient and internal temperature

In Fig. 13, we represented the temperature profiles and VMRs obtained with our upper and lower uncertainty boundaries on the internal temperature of K2-18b (assuming tp between 1 and 13 Gyr), as well as our selected nominal value (80 K). It appears that the internal temperature has only a minor effect on the temperature profile, except in the deepest layers of the atmosphere. This effect was also noted by other authors such as Morley et al. (2017). This is because the received flux (the stellar irradiance) is two to three orders of magnitude higher than the internal flux. Hence, above the convective layers (here at ≈ 20 kPa), the difference of temperature between our two extreme simulations is less than 40 K (8 K on average). However, below the convective layers (≈ 200 kPa), the atmosphere is opaque, thus, the stellar irradiance can no longer heat the atmosphere and the internal heating become prominent. Consequently, at the bottom of our pressure grid, the difference in temperature between our two extreme models reaches ≈ 500 K. This affects the species abundances, in the upper atmosphere as well if the quench level of a species is at a pressure level below the region of stellar heating.

In Fig. 14, we show the effect of the eddy diffusion coefficient on the VMR of the main C-, N-, and O-bearing species at 175 (Z/H)⊙. The VMR at thermochemical equilibrium are also displayed. The main effect of increasing Kzz is to shift down the quench levels of the species.

In Fig. 15, we display the VMR of CH4, CO, and H2O for Kzz = 106, 108, and 1010 cm2 s−1 as a function of the age of the planet and the corresponding internal temperature (see Eq. (2)). As theinternal temperature decreases (or as the planet gets older), the impact of the Kzz on these VMR gets smaller. If the Kzz is low (≾ 106 cm2 s−1), the chemical composition of the atmosphere essentially does not change. On the other hand, under vigorous mixing, the chemical composition of the upper atmosphere can be significantly altered over time. With Kzz ≿ 1010 cm2 s−1, if the planet is young (≾1 Gyr or Tp, int ≿ 100 K), CO could become the dominant O-bearing species and remain more abundant than CH4 for several Gyr. In parallel, Kzz could also decrease over time, as the heat to be dissipated by convection decreases. Moreover, if the quench level switches from a convective layer to a radiative one, Kzz could drop by several orders of magnitude, rapidly changing the chemistry of the upper atmosphere if this occurs when the planet is still young.

We compared our simulations within our range of Kzz, at Tp, int of 45, 80, and 115 K and at metallicities of 50, 75, 100, 125, 150, 175, 200, and 300 (Z/H)⊙. The best-fit parameters against Benneke et al. (2019) data for each tested internal temperature can be found in Table 7. We found no significant differences at constant metallicity in terms of goodness of fit for the tested Tp, int, the χ2 value varying by less than 0.7 across all Kzz values. In each case, our best fit is located at 175 (Z/H)⊙. The slight difference in goodness of fit at Tp, int = 45 K compared to the other tested internal temperatures is due to H2O cloud condensation occurring in the former case and not in the latter. χ2 values above 200 and below 150 (Z/H)⊙ are systematically above the 1σ confidence level. We note that even at 175 (Z/H)⊙, Tp, int = 115 K and Kzz = 1010 cm2 s−1, we retrieve a CH4 VMR above 1.5% in the upper atmosphere. Hence, internal heating from a relatively young (tp ≈ 1 Gyr) K2-18b and vigorous vertical mixing alone cannotexplain a CH4-depleted atmosphere.

thumbnail Fig. 12

H2O cloud particle radius effect on the transmission spectrum, at 300 (Z/H)⊙ and at nominal irradiation. The χ2 against Benneke et al. (2019) data is indicated in parentheses.
thumbnail Fig. 13

Effects of Tp, int on K2-18b, at 175 (Z/H)⊙, nominal irradiation, and Kzz = 106 cm2 s−1. Left: temperature profiles. The red area represents the convective layers. Right: VMR of the most abundant absorbers. Dotted: Tp, int = 45 K. Dashed: Tp, int = 80 K. Solid: Tp, int = 115 K.
thumbnail Fig. 14

Effects of Kzz on the VMR of K2-18b main C, N, O bearing species, at 175 (Z/H)⊙. Dashed: Kzz = 106 cm2 s−1. Dotted-dashed: Kzz = 108 cm2 s−1. Solid: Kzz = 1010 cm2 s−1. Dotted: thermochemical equilibrium.
thumbnail Fig. 15

Chemical evolution of the upper atmosphere of K2-18b for the most abundant C- and O-bearing species over time, depending on the value of Kzz, at 175 (Z/H)⊙. Dotted: Kzz = 106 cm2 s−1. Dashed: Kzz = 108 cm2 s−1. Solid: Kzz = 1010 cm2 s−1.
Table 7

Best fits against Benneke et al. (2019) dataset as a function of internal temperature.

4.4 C/O ratio

We simulated atmospheres with 30, 50, 75, 100, 125, 150, 175, 200, 300, 400, and 500 times the solar system metallicity for all elements except C and the noble gases, and with a carbon-to-hydrogen abundance ratio (C/H) of 0.3, 1, 2, 5, 10, 30, 50, 75, 100, 150, 200, 300, 400, and 500 times the solar ratio. The χ2 of these simulations against Benneke et al. (2019) data are displayed in Fig. 16. We call that our ’free C/H’ scenario. Our best fit against Tsiaras et al. (2019) data (χ2 = 17.07) is located at 125 (Z/H)⊙ and 50 (C/H)⊙ (C/O = 0.226), while our best fit against Benneke et al. (2019) data (χ2 = 18.86) is located at 125 (Z/H)⊙ and 30 (C/H)⊙ (C/O = 0.13). A discussion of these results is available in Sect. 5.2.2.

5 Comparisons with previous studies and discussions

5.1 Considering whether H2O or CH4 is the dominant absorber in the HST transit spectra

Our nominal model, presented in Sect. 4, shows that CH4 should be the main C-bearing species in the atmosphere, and, despite H2 O being more abundant than CH4, the latter should be the main contributor to the absorption features of the transmission spectrum. In particular, CH4 and H2 O have an overlapping band at 1.4 μm in the HST-WFC3 spectral range. This is outlined in Bézard et al. (2020), who showed that at thiswavelength and for sub-Neptunes with effective temperatures lower than ≈ 600 K and high atmospheric metallicities, CH4 should be a stronger absorber than H2O in transit spectroscopy due to the numerous weak lines of CH4. However, this is in disagreement with Benneke et al. (2019) and Tsiaras et al. (2019) who claim that the HST data provide evidence for the presence of H2O.

To investigate this discrepancy, Bézard et al. (2020) compared their best fit Exo-REM spectrum (in which CH4 absorption dominates over H2O) and an H2O-only spectrum (in which all other absorbers were removed) to the HST data reduced by Tsiaras et al. (2019). The χ2 associated with the two models are essentially the same (0.91 and 0.93). In addition, applying the retrieval algorithm TauREx 3 (Al-Refaie et al. 2019) to these data, Bézard et al. (2020) found solutions favouring a CH4-rich atmosphere whose absorption at 1.4 μm is dominated by CH4. Thus, the apparent disagreement with the analysis of Tsiaras et al. (2019) simply arises from the fact that these authors did not consider an atmosphere with significant amounts of CH4 in the three scenarios they investigated.

Regarding the Benneke et al. (2019) dataset, we used the same technique as Bézard et al. (2020): we took the H2 O contribution of our best-fit spectrum, displayed in Fig. 3, and compared it with the data. In the same way, we made another spectrum with only the contributions of the absorbers included in the retrieval analysis of Benneke et al. (2019) except CH4 (i.e. H2 O, CO, CO2, and NH3, HCN contribution is negligible and was not included). We will refer to this latter spectrum as “no CH4”. We stress that none of these spectra are directly the results of self-consistent simulations. They are displayed in Fig. 17. We found that our H2O-only spectrum does provide a superior fit compared with our best fit Exo-REM model. The fit is improved in the 1.15–1.20 μm spectral range, around 1.6 μm, and for the Spitzer data points. With our “no CH4” spectrum, the fit is also better than our best-fit model, but less so than our H2 O-only spectrum due to the absorptions of NH3 in the HST spectral range and the CO and CO2 absorptions in the Spitzer spectral range. We note that our VMRs of NH3, N2, CO, and CO2 of respectively 0.04%, 1.02%, 2.73%, and 0.96% (see Fig. 5), are all lower than the “2σ (97.5%)” upper limits given by Benneke et al. (2019) for these species (respectively 13.5%, 10.9%, 7.45% and 2.4%). The spectrum with only H2O absorption is strongly overfitting the Benneke et al. (2019) data, with a χ2 of 13.01. On the other hand, our best-fit model, with a χ2 of 21.07, is stillwithin the 1-σ confidence interval (χ2 < 21.36). The “no CH4” spectrum ischaracterised by a slight overfitting, providing a better χ2 (18.83) compared to our best fit, but the removal of CH4 absorption is artificial.

We also applied TauREx 3 to the Benneke et al. (2019) data. This algorithm uses the nested sampling code Multinest (Feroz et al. 2009)to explore the parameter space and find the best fit corresponding to a given spectrum. In this retrieval analysis, we used 500 live points, an evidence tolerance of 0.5, and the cross-sections provided by the TauREx website7 at a resolution power of 15 000 and between 0.3 and 15 μm. We simulated the atmosphere of K2-18b using an isothermal temperature profile with 100 atmospheric layers between 10−3 and 106 Pa. We took into account molecular absorptions from CH4, CO, CO2, H2 O, and NH3, CIA from H2 –H2 and H2–He, the contributionto μ of N2, Rayleigh scattering and spectrally gray clouds. The planet radius was allowed to vary by ± 10% of its value determined by Benneke et al. (2019), while the cloud top was allowed to be between 106 and 1 Pa. We chose to set bounds on the temperature and the species VMR, respectively, between 200 and 400 K, and between 10−10 and 0.3. The corresponding posterior distributions are shown in Fig. 18. We also tested models removing all molecular absorptions except H2O or CH4, removing clouds, or including all absorptions except that of H2O. For comparison, we also computed a “no active gas” model including only the cloud contribution, that is, a flat spectrum. Finally, we calculated the χ2 and corresponding σ confidence of each model using the same technique as in Sect. 4 and the “binned” spectrum output from TauREx 3. Our results are summarised in Table 8. A small discussion on the results is available in Appendix B.

We found that solutions with low amounts of CH4 are unambiguously favoured: the “all absorbers with clouds” model we used is favoured over the “no H2 O with clouds” model at 273:1 (3.79σ, see Benneke & Seager 2013). The difference in log-evidence between our CH4-only models and our “H2O-included” models (i.e. all models including at least the molecular absorption of H2 O) is, at worst of 4.16, indicating that the latter are favoured at ≥64:1 (3.35σ). Moreover, the 2σ (95.4%) upper limit of CH4 is 0.009% with our “all absorbers with clouds” model, even lower than what was inferred by Benneke et al. (2019; 0.248%). This confirms that an H2 O-dominated spectrum is a far better fit to the Benneke et al. (2019) dataset than a CH4-dominated spectrum, but not that the latter must be rejected. Indeed, looking at the χ2, we find that all of the spectra obtained from our TauREx 3 models that do not include H2 O are slightly below the 1σ confidence level, with χ2 values and CH4 VMR close to our Exo-REM best fit. In contrast, all the H2O-included TauREx 3 spectra strongly overfit the spectrum, with χ2 ≤ 10.50, which is much less than the number of data points, and likely much less than the number of effective degrees of freedom. We interpret this as an H2O-only spectrum leading to a significant overfit of the dataset.

The overfit that we identified could be explained by an over-estimation of the HST error bars. However, the Benneke et al. (2019) and Tsiaras et al. (2019) HST spectra that we analysed come from independent reduction methods (Tsiaras et al. 2018)and both show similar transit depth uncertainties. Moreover, the custom HST pipeline used by Tsiaras et al. (2019) is described in Tsiaras et al. (2018) to exhibit “nearly photon-noise limited” performance, meaning that the spectral uncertainties should be close to the theoretical minimum. It can also be seen, for example, in Fig. 7 that the data points near 1.2 μm from the two datasets are not compatible at the 1 − σ confidence level, which points to additional systematic errors. These are probably linked to uncertainties in some parameters used in the data reduction (orbital parameters, limb-darkening coefficients, or calibration, as suggested in Tsiaras et al. 2018). For all of these reasons, we regard this possibility as unlikely.

To summarise, of the two datasets available to us, only the one from Benneke et al. (2019) strongly favours a CH4-depleted atmosphere. This favouritism is likely due to free retrieval algorithms searching for the solution that best fits the data. But if such a solution is strongly overfitting, as is the case for an H2 O-only spectrum against the dataset from Benneke et al. (2019), any other solution that is considered to be more physically acceptable may also appear as statistically unlikely. However, if a strong depletion of CH4 in the atmosphere of K2-18b is not a straightforward scenario from a chemical standpoint and ends up overfitting the data, it is indisputablya statistically valid one and, thus, it must be considered. In the next section, we try to find a model that could explain sucha depletion.

thumbnail Fig. 16

Goodness of fit of our models for Tp, int = 80 K, Kzz = 106 cm2 s−1 and nominal irradiation, as a function of the metallicity and C/H, compared to the data from Benneke et al. (2019; left) and Tsiaras et al. (2019; right). The white colour corresponds to the value of χ2 = 21.36, indicating a fit at the 1σ confidence level. The solid black line represents the 1σ confidence level. The dotted lines indicates the C/O ratio. The solid red line indicates which species dominates on average the transmission spectrum in the 1.355–1.415 μm range. The cases where C/O > (C/O)⊙ was not explored are represented as black rectangles.
thumbnail Fig. 17

Best fit against Benneke et al. (2019) dataset with transit spectra in which the spectral contributions of selected species has been removed. Blue: H2 O contribution of the transit spectrum displayed in Fig. 3, with an offset of the 105 -Pa level of +98 km. Cyan: H2O, CO, CO2, and NH3 total contribution of the transit spectrum displayed in Fig. 3, with an offset of the 105 -Pa level of + 57 km. Black: K2, HST and Spitzer from Benneke et al. (2019). The dotted red line is our best-fit spectrum with all absorptions (see Sect. 4). The χ2 of these spectra against the data is indicated in parentheses.
thumbnail Fig. 18

TauREx 3 atmospheric retrieval posterior distributions against the dataset from Benneke et al. (2019), with planetary radius (R♃), temperature (K), decimal logarithm of the VMR of H20, CH4, CO, CO2, N2, and NH3, and decimallogarithm of the cloud top pressure (Pa) as free parameters.
Table 8

Comparison of the Bayesian log-evidence of TauREx 3 against Benneke et al. (2019) dataset for different models.

5.2 CH4-depleted scenarios for K2-18b

While our standard models, presented in Sect. 4, are statistically able to reproduce the observed spectrum of K2-18b, here we test some scenarios allowing for CH4 VMR to be as small as retrieved by all the other teams who analysed the data so far. Indeed, our nominal model gives a CH4 VMR of (for metallicities between 65 and 500 (Z/H)⊙), while Benneke et al. (2019), Madhusudhan et al. (2020) and Scheucher et al. (2020) give, respectively, a 2σ upper limit of 0.248%, a 99% upper limit of 3.47%, and an upper limit of 460 ppm. Benneke et al. (2019) propose three explanations to this depletion: (i) a high internal temperature, either from residual heat of formation or tidal heating, (ii) a low C/O ratio resulting from planetary formation process, and (iii) a catalytic destruction of CH4 by photolysis. The latter cannot be simulated using Exo-REM, so we will focus on the first two possibilities.

5.2.1 High internal temperature

Following the trend presented in Sect. 4.3, it is possible to reach a CH4 VMR that is compatible with what was found by Benneke et al. (2019) by considering Tp, int as high as 200 K. Using this internal temperature, we simulated atmospheres at 30, 50, 75, 100, 125, 150, 175, 200, 300, 400, and 500 (Z/H)⊙ and within our selected range of Kzz. The goodness of fit for our models compared to the datasets from Benneke et al. (2019) and Tsiaras et al. (2019) is displayed in Fig. 19. We note that H2 O never dominates the transmission spectrum in the 1.355–1.415 μm window. Our best fit against Tsiaras et al. (2019) data (χ2 = 16.78) is located at 150 (Z/H)⊙ and 107 cm2 s−1, while our best fit against Benneke et al. (2019) data (χ2 = 20.70) is located at 150 (Z/H)⊙ and 1010 cm2 s−1. There is no H2O cloud condensation in the latter scenario. The goodness of fit is again slightly better than for our nominal model. The VMR are displayed in Fig. 20. In this configuration, we obtain VMR for CH4, CO, CO2, NH3, and H2 O in the upper atmosphere at levels of 0.283%, 6.589%, 0.280%, 0.006%, and 1.657%, respectively, which are all belowor close to the 2σ upper limit found by Benneke et al. (2019).

Reaching such a high internal temperature, however, would require K2-18b to be a very young planet (≪ 1 Gyr old according to Eq. (2)) or a planet experiencing strong tidal heating, or both. According to the internal temperature model we used and the age of the system determined by Guinan & Engle (2019), it seems unlikely that the residual heat of formation could be sufficient to explain a high internal temperature. To assess the possibility of strong tidal heating on K2-18b, we can make a few considerations. The K2-18 system is composed of at least another planet: K2-18 c (Cloutier et al. 2019), orbiting closer to its star than K2-18b. Both planets may have a moderately high eccentricity that is rapidly evolving due to the secular interactions between the two objects (Gomes & Ferraz-Mello 2020). In that configuration, a lot of orbital energy can indeed be dissipated, but only in the innermost planet8, that is, K2-18 c. Thus, there should be no tidal heating on K2-18b.

thumbnail Fig. 19

Goodness of fit of our models at Tp, int = 200 K and nominal irradiation as a function of metallicity and Kzz, compared tothe data from Benneke et al. (2019; left) and Tsiaras et al. (2019; right). The white colour corresponds to the value of χ2 = 21.36, indicating a fit at the 1σ confidence level. The best fit against Benneke et al. (2019) data (χ2 = 20.64) is located at (Z/H)⊙ = 125 and Kzz = 1010 cm2 s−1. The best fit against Tsiaras et al. (2019) data (χ2 = 16.56) is located at (Z/H)⊙ = 200 and Kzz = 105 cm2 s−1.
thumbnail Fig. 20

VMR of selected species of our best-fit model against Benneke et al. (2019) data at Tp, int = 200 K and nominal irradiation.

5.2.2 Low C/O

To test the possibility of a low atmospheric C/O ratio on K2-18b, we use the results from the simulations presented in Sect. 4.4. The VMR of the most abundant absorbers in our best fit against Benneke et al. (2019) are represented in Fig. 21. The amount of CH4 in the upper atmosphere we obtain in this best-fit simulation is 1.23%, incompatible with the 2σ upper limit found by Benneke et al. (2019; 0.248%), but compatible with the 99% upper limit found by Madhusudhan et al. (2020; 3.47%). However, simulations at metallicities between 50 and 150 (Z/H)⊙ at 5 (C/H)⊙ and at 20 (Z/H)⊙ and 1 (C/H)⊙ are statistically able to reproduce the data while allowing the CH4 VMR to be close or below 0.248%, with a minimum at 0.101%. There is no simulation within the 1σ confidence level for C/H ≾1 (C/H)⊙ (CH4 VMR ≾ 0.05% at 50 (Z/H)⊙). Accordingly, C/O must be ≿0.01 (0.02 (C/O)⊙) in order to obtain a satisfactory fit of the data. We note that C/H must be ≾ 5 (C/H)⊙ (or C/O ≾ 0.1 (C/O)⊙) in order for H2 O to dominate the transmission spectrum in the 1.355–1.415 μm interval.

Heavily C-depleted atmospheres (C/O ≾ 0.01) do not seem to agree with the data. Lowering the abundance of C in the atmosphere naturally decreases the contribution of CH4 to the transmission spectrum, and improves the fit of Benneke et al. (2019) data, as discussed in Sect. 5.2. A side effect, however, is that because CH4 has a strong greenhouse effect, the temperatures decrease significantly, favouring the condensation of H2 O, dampening the absorption features of the latter, and worsening the fit around 1.4 μm. This side effect explains why our results are quite different from the results of other teams who used free retrieval algorithms that do not include any physical constraint on the model atmosphere. Moderately-C-depleted atmospheres (C/O ≈ 0.10), on the contrary, are clearly favoured, but in that case CH4 still dominates the absorption features of the transmission spectrum. Moreover, the precise phenomenon that could lead to such a depletion has, to our knowledge, yet to be discovered.

thumbnail Fig. 21

VMR of selected species at Tp, int = 80 K, Kzz = 106 cm2 s−1, nominal irradiation, 125 (Z/H)⊙, and 30 (C/H)⊙.
thumbnail Fig. 22

Contributions within some spectral ranges of a selection of absorbing species to the transmission spectrum of our nominal model (175 (Z/H)⊙) at nominal irradiation and a resolution of 0.5 cm−1 (resolving power of 20 000 at 1 μm). The spectral contribution of individual species takes the CIA and Rayleigh scattering (dotted curve) into account.

6 Observational perspectives

For K2-18b, acquiring higher resolution spectra would permit to discriminate between the three scenarios we identified. We used our 0.5 cm−1 resolution k-coefficient set to derive high-resolution spectra of K2-18b and show some spectral ranges in Fig. 22 that could prove interesting in that regard. We note that NH3 dominates the 1.45–1.56 μm spectral range for any of the following three scenarios.

Firstly, in the nominal scenario, CH4 should dominate over most of the infrared (IR) transmission spectrum, with particularly strong lines around 1.6 μm (end of the H band), or in the 3.0–4.0 μm spectral range (L band). CO would dominate the spectrum in the 4.6–5.0 μm spectral range (M band) but not in the 2.3–2.5 μm spectral range (end of the K band). Isolating H2O lines would be challenging in the IR, but might be possible around 1.16 μm and 1.37 μm (beginning and end of the J band). However, it should dominate around 823 and 924 nm.

Secondly, in the low C/O scenario, contrary to the nominal scenario, H2 O should overall dominate the IR spectrum over CH4. CH4 may still dominate the spectrum at the end of the H band (1.6–1.8 μm), and in the K and L bands, depending on the intensity of the C-depletion. However, in contrast to the nominal scenario, it should have a minor contribution in the J band. CO might still dominate the M band, again depending on the extent of the of C-depletion.

Thirdly, in the scenario involving high internal temperature or incomplete chemistry, H2 O should again dominate most of the IR spectrum. CO should dominate in the M band and possibly even in the K band, in contrast to the low C/O scenario. The CH4 spectral contribution should be similar to the low C/O scenario.

The James Webb Telescope (JWST) and the Atmospheric Remote-sensing Infrared Exoplanet Large-survey (ARIEL) mission, with their broader spectral range and precision higher than the data we studied here, will probably be very helpful to discriminate between these scenarios. In complement, we could acquire higher resolution spectra from ground-based instruments – such as the CRyogenic InfraRed Echelle Spectrograph+ (CRIRES+) or the Echelle SPectrograph for Rocky Exoplanets and Stable Spectroscopic Observations (ESPRESSO) mounted on the Very Large Telescope (VLT) – or by using the incoming generation of large telescopes (Extremely Large Telescope, Thirty Meter Telescope, etc.) in order to detect the absorption of individual lines and, thus, to unambiguously detect species.

As an example, we calculated the signal to noise ratio (S/N) we could obtain for a detection of CH4 in our nominal case using CRIRES+ at a resolution power of 100 000 in the 1.5–1.8 μm spectral range. This spectral range has the advantage of presenting a lot of deep CH4 lines (50–100 ppm) in a band for which CRIRES+ has a high sensitivity9. We calculated a synthetic transmission spectrum at a resolving power of 100 000 using a line-by-line radiative transfer program, including absorption from H2 –H2 and H2 O–H2O CIA and lines from H2O, CH4, NH3 and CO. Wecalculated the S/N, while neglecting the terrestrial absorptions and the photon noise of the star, using the following: (5)

where λ is the discretised wavelength, λmin and λmax are respectively the minimum and maximum wavelengths of the considered spectral range, and δtot are, respectively, the transit depths that consider only the absorption from CH4 and the absorption from all absorbers; the bar indicates a smoothing by a boxcar of width 16 resolution elements (also used by Brogi et al. 2016, for their CRIRES observation of HD 189733), N is the number of calculated points per resolution element within the spectral range, and S/N* is the expected CRIRES+ S/N on K2-18 for one hour of observation at R = 100 000. From the H limiting magnitude given in eso.org, we derive a S/N of 2500 for K2-18 at R = 50 000 per spectral dispersion element for a 1-hr integration. At R = 100 000, assuming no loss of the star signal, the S/N should thus be . Hence, we can obtain a total CH4 signal of 1500 ppm and, thus, from Eq. (5), we have an optimistic value for of ≈ 2.7 for one hour of observation with CRIRES+ at R = 100 000. We note that the transit duration of K2-18b is approximately three hours. A full transit would give S/N, which is sufficient for detecting CH4. However, this estimate is probably optimistic as we did not take into account the telluric absorption and the photon noise from the star.

Observing CO from the ground would, on the other hand, prove to be extremely difficult. Even at R = 100 000, the CO lines overlap with those of the star. The Doppler shift between the planet and the star during the transit is in the range ± 0.7 km s−1, which is too low to use a cross-correlation method as described by, for example, Snellen et al. (2010).

thumbnail Fig. 23

Atmospheric metallicity as a function of mass for a selection of planets. The Jupiter (J), Saturn (S), Uranus (U), and Neptune (N) C/H are taken respectively from Wong et al. (2004), Fletcher et al. (2009), Sromovsky et al. (2011) and Karkoschka & Tomasko (2011) (values compiled in Atreya et al. 2018). The data for GJ436b and WASP-127b are taken respectively from Morley et al. (2017) and Skaf et al. (2020), while the data for all the other planets are taken from Welbanks et al. (2019). The black area represents the 1σ C/H against mass trend for Jupiter, Saturn, Uranus, and Neptune resulting from a linear fit. The fit is log10 ((Z/H)∕(Z/H)⊙) = −1.10 ± 0.21log10(M[kg]) + 30.5 ± 5.5. The error bars of our results correspond to the simulations that were within the 1σ confidence interval against Tsiaras et al. (2019) data in the “free C/H” scenario. Our C/H results were slightly shifted to the right for clarity. The O/H ratios in the other works were estimated assuming [H2 ] + [He] + [H2O] + [CH4 ] = 1 and a solar He/H.
thumbnail Fig. 24

Summary of our results and results of other teams on K2-18b metallicity. The error bars of our results corresponds to the simulations that were within the 1σ confidence level against Tsiaras et al. (2019) data. The O/H and C/H ratio of the other works were estimated assuming [H2 ] + [He] + [H2O] + [CH4 ] = 1 and a solar He/H. The C/H upper limit of Benneke et al. (2019) and Madhusudhan et al. (2020) corresponds respectively to the 2σ and the 99% upper limit. The dotted line represents the metallicity expected for the mass of K2-18 b using the linear fit in Fig. 23, based on solar-system C/H ratios.
thumbnail Fig. 25

Comparison of K2-18b transmission spectrum best-fits for three different scenarios. Red: 175 (Z/H)⊙, 175 (C/H)⊙, Kzz = 106 cm2 s−1, Tp, int = 80 K. Green: 125 (Z/H)⊙, 30 (C/H)⊙, Kzz = 106 cm2 s−1, Tp, int = 80 K. Blue: 125 (Z/H)⊙, 125 (C/H)⊙, Kzz = 1010 cm2 s−1, Tp, int = 200 K. The χ2 against Benneke et al. (2019) data is indicated in parentheses.

7 Summary and conclusions

We analysed the transmission spectrum of K2-18b between 0.43 and 5.02 μm using a combination of K2, HST, and Spitzer data retrieved from Benneke et al. (2019) and Tsiaras et al. (2019). We studied the effect of irradiation, metallicity, clouds, internal temperature, eddy diffusion coefficient, and C/O ratio using our self-consistent model Exo-REM and assuming an atmosphere primarily composed of H2 and He. We also analysed the Benneke et al. (2019) dataset with the retrieval algorithm TauREx 3, and provided a S/N estimation for the detection of CH4 using CRIRES+.

We found that the data are compatible with a highly metal-enriched atmosphere, between 65 and 500 (Z/H)⊙ against Tsiaras et al. (2019) data – or between 100 and 200 (Z/H)⊙ against Benneke et al. (2019) data – when assuming a solar C/O ratio, and ≿ 40 (Z/H)⊙ when assuming a sub-solar C/O ratio. According to our results, the atmosphere of K2-18b appears quite similar to that of Neptune or Uranus and seems to follow the C/H-mass relationship of the giant planets of our solar system, as shown in Fig. 23. Most of the extrasolar planets seem to have an O/H ratio below that expected from the C/H-mass relationship for solar-system giant planets. Part of the explanation could reside in the fact that for most of these planets, the amount of CO and CO2 is unknown, leading to an underestimation of the O/H ratio (as pointed out by e.g. Wakeford et al. 2017; Welbanks et al. 2019). A summary of our results and results of other teams is displayed in Fig. 24, and the spectra of our best fit for each of the explored scenarios are displayed in Fig. 25.

We also show that thick to no H2O-ice clouds are allowed by the data. Liquid H2O clouds are possible on planets similar to K2-18b but receiving at most 80% of its irradiation. With or without clouds, we found that the Bond albedo of the planet should be around 0.02. The other studied parameters are not well constrained by the data. We note that in most of the cases we studied, CH4 absorption should dominate or be on par with that of H2O in the HST spectral window, as first outlined by Bézard et al. (2020).

Combining the work of Bézard et al. (2020) and our own retrievals using TauREx 3, we show that the discrepancy between our self-consistent results and the results from retrieval algorithms can be explained by either the discarding of CH4 absorptions, as in the case of the Tsiaras et al. (2019) dataset, or a strong overfitting of the data, as in the case of the Benneke et al. (2019) dataset.

Accordingly, in addition to our nominal scenario (i.e. Z/H = 65–500 (Z/H)⊙ and a solar C/O), scenarios with a CH4-depleted atmosphere could satisfactorily fit the observed spectrum. These could be obtained with a high internal temperature (≿ 200 K) or a low C/O ratio (≿0.01), with H2 O being the dominant absorber in HST/WFC3 spectral band if C/O ≾0.1. However, the high internal temperature scenario seems very unlikely, and the C-depleted scenario requires the existence of an unknown process. Therefore, it seems that a CH4-depleted atmosphere, a scenario favoured by all other teams who have analysed the same data so far, is more difficult to defend than a nominal scenario with a standard abundance of CH4. Moreover, a spectrum dominated by H2O absorptions seems even more unlikely given the relatively small range of self-consistent solutions allowing for this scenario and the fact that it is favoured by retrieval algorithms only because of overfitting. Another possibility is that our thermochemical model does not accurately describe the chemical conversion between CO and CH4 in the deep atmosphere, as suggested by Benneke et al. (2019) in their atmospheric analysis of GJ3470b (noting that GJ3470b is a much warmer planet with an effective temperature of ≈ 800 K). It should be possible to discriminate between these scenarios by acquiring data in different spectral ranges and/or at higher resolution using the incoming JWST and ARIEL space telescopes. In complement, ground-based instruments such as CRIRES+ could be used to detect individual lines of species. We estimate that with our nominal scenario, we could detect CH4 with CRIRES+ in the H band with an “optimistic” S/N of 4.7 after three hours of transit observation.

This study highlights the critical need for more precise spectroscopic data on sub-Neptunes. This will allow us to better characterise these intriguing bodies, which have no analogue in our Solar System.

Acknowledgements

D. B. acknowledges financial support from the ANR project “e-PYTHEAS” (ANR-16-CE31-0005-01). We acknowledge support from the Programme National de Planétologie of INSU/CNRS co-funded by CNES. We thank M. Rey for providing the TheoReTS CH4 line lists up to 13400 cm−1 over the whole Exo-REM temperature grid. We thank O. Venot for providing us with the results from Venot et al. (2020)’s thermochemical model for K2-18b. We thank J. Leconte for his valuable insight on tidal heating.

Appendix A Note on the effect of the H2O CIA

The effect ofH2O CIA on transmission spectra and temperature profiles is shown in Fig. A.1. In our K2-18b simulations at 175 (Z/H)⊙, we found that neglecting the H2O CIA leads to an underestimation of the atmospheric temperature by ≈ 40 K at 10 kPa. This can have an effect on the nature of H2O condensation as well. Taking a model with 0.8 time the nominal irradiance of K2-18b, if we include H2 O CIA, we find that H2O condenses into its solid phase, whereas if we do not, it condenses into its liquid phase at slightly higher pressures (≈ 8 kPa instead of 3 kPa). Moreover, because the cloud is forming lower without the CIA, it has less effect on the transmission spectrum. As a consequence, the amplitude of the absorption features of the transmission spectrum is increased. However, the effect on the temperature profile is smaller at low metallicities, with a temperature difference at 10 kPa reaching 20 K at 30 (Z/H)⊙, and 2 K at 3 (Z/H)⊙, at ≈300 kPa for both cases.

thumbnail Fig. A.1

Effect of the H2O CIA at 175 (Z/H)⊙ on the simulated transmission spectrum and temperature profile of K2-18b. Left: temperature profiles. Solid black: nominal model. Dotted-dashed black: nominal model at 0.8 time the nominal irradiance. Red: nominal model without H2 O CIA. Dotted-dashed red: nominal model at 0.8 times the nominal irradiance and without H2 O CIA. The phase diagram of H2O is represented as dotted blue lines, the dot corresponding to the H2O triple point. Right: transmission spectra at 0.8 time the nominal irradiation, with (black) and without (red) H2 O CIA.

CIA areinduced by molecules colliding with themselves or other species. On planets with low metallicity, H2 –H2 and H2–He collisions are the only significant contributors to CIA because species other than H2 and He are present in low quantities. However, this no longer holds true at higher metallicity. As the abundance of other species increases, their CIA increase which has an overall warming effect on the atmosphere.

Appendix B TauREx 3 retrievals

From our TauREx 3 retrievals, we found an Atmospheric Detectability Index (ADI, see Tsiaras et al. 2018) of ≥ 1.17, including only CH4 (2.12σ, or 3:1 relative odds), and ≥6.37 (3.99σ, or 584:1 relative odds), including at least H2O, indicating “weak” to “strong” detections of atmospheric absorptions. Similarly to Madhusudhan et al. (2020), we find only a “weak” detection of clouds by comparing similar models with and without clouds: the difference in log-evidence is at most 1.04 (2.05σ or 2.8:1 relative odds) with our CH4-only models and at least 0.29 (1.48σ or 1.3:1 relative odds) with our H2O-only models. Comparing our “all absorbers with clouds” model with our “no H2 O” and “no CH4” models, we also find that H2O is detected at 273:1 (3.79σ, that is, a “strong” detection), while CH4 is detected at 1.30:1 (1.44σ, meaning a “not significant” detection).

Regarding our retrieval with all absorbers and clouds (see Fig. 18), our retrieved VMR of H2 O (0.072–5.3%) is similar to what Benneke et al. (2019) found (0.033–8.9%). Our 2σ (95.4%) upper limits for CO, CO2, NH3, and N2 are respectively 4.39%, 0.789%, 0.002%, and 2.73%, all lower than those from Benneke et al. (2019; 7.45%, 2.4%, 13.5%, and 10.9%), especially for NH3. The cloud top pressure we retrieve (1.8–66.1 kPa) is slightly higher than theirs (0.77–13.9 kPa). From our Exo-REM simulations, this level corresponds more to NH4Cl condensation than to H2O condensation, but we showed that including NH4Cl clouds does not significantly enhance our fits (see Sect. 4.2).

References

  1. Allard, N. F., Kielkopf, J. F., Spiegelman, F., Tinetti, G., & Beaulieu, J. P. 2012, A&A, 543, A159 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Allard, N. F., Spiegelman, F., & Kielkopf, J. F. 2016, A&A, 589, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Al-Refaie, A. F., Changeat, Q., Waldmann, I. P., & Tinetti, G. 2019, ArXiv e-prints [arXiv:1912.07759] [Google Scholar]
  4. Atreya, S. K., Crida, A., Guillot, T., et al. 2018, in Saturn in the 21st Century, 20, 5 [CrossRef] [Google Scholar]
  5. Azzam, A. A. A., Tennyson, J., Yurchenko, S. N., & Naumenko, O. V. 2016, MNRAS, 460, 4063 [NASA ADS] [CrossRef] [Google Scholar]
  6. Baudino, J. L., Bézard, B., Boccaletti, A., et al. 2015, A&A, 582, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Baudino, J. L., Mollière, P., Venot, O., et al. 2017, ApJ, 850, 150 [NASA ADS] [CrossRef] [Google Scholar]
  8. Benneke, B., & Seager, S. 2013, ApJ, 778, 153 [NASA ADS] [CrossRef] [Google Scholar]
  9. Benneke, B., Werner, M., Petigura, E., et al. 2017, ApJ, 834, 187 [NASA ADS] [CrossRef] [Google Scholar]
  10. Benneke, B., Knutson, H., Lothringer, J., et al. 2019, Nat. Astron., 3, 813 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bernath, P. F. 2020, J. Quant. Spectr. Rad. Transf., 240, 106687 [CrossRef] [Google Scholar]
  12. Bézard, B., Charnay, C., & Blain, D. 2020, Nat. Astron., submitted, [arXiv:2011.10424] [Google Scholar]
  13. Borysow, A. 2002, A&A, 390, 779 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  14. Borysow, A., J⊘rgensen, U. G., & Fu, Y. 2001, J. Quant. Spectr. Rad. Transf., 68, 235 [NASA ADS] [CrossRef] [Google Scholar]
  15. Brogi, M., de Kok, R. J., Albrecht, S., et al. 2016, ApJ, 817, 106 [Google Scholar]
  16. Brown, L. R., & Peterson, D. B. 1994, J. Mol. Spectr., 168, 593 [NASA ADS] [CrossRef] [Google Scholar]
  17. Burch, D. E., Gryvnak, D. A., Patty, R. R., & Bartky, C. E. 1969, J. Opt. Soc. Am., 59, 267 [NASA ADS] [CrossRef] [Google Scholar]
  18. Burrows, A., Marley, M. S., & Sharp, C. M. 2000, ApJ, 531, 438 [NASA ADS] [CrossRef] [Google Scholar]
  19. Charnay, B., Bézard, B., Baudino, J. L., et al. 2018, ApJ, 854, 172 [Google Scholar]
  20. Chase, M. 1998, NIST-JANAF Thermochemical Tables, Tech. rep. [Google Scholar]
  21. Chen, J., & Kipping, D. 2016, ApJ, 834, 17 [NASA ADS] [CrossRef] [Google Scholar]
  22. Cloutier, R., Astudillo-Defru, N., Doyon, R., et al. 2019, A&A, 621, A49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Coles, P. A., Yurchenko, S. N., & Tennyson, J. 2019, MNRAS, 490, 4638 [CrossRef] [Google Scholar]
  24. Feroz, F., Hobson, M. P., & Bridges, M. 2009, MNRAS, 398, 1601 [NASA ADS] [CrossRef] [Google Scholar]
  25. Fletcher, L. N., Orton, G. S., Teanby, N. A., Irwin, P. G. J., & Bjoraker, G. L. 2009, Icarus, 199, 351 [NASA ADS] [CrossRef] [Google Scholar]
  26. Fortney, J. J., Mordasini, C., Nettelmann, N., et al. 2013, ApJ, 775, 80 [NASA ADS] [CrossRef] [Google Scholar]
  27. Fouchet, T., Orton, G., Irwin, P. G. J., Calcutt, S. B., & Nixon, C. A. 2004, Icarus, 170, 237 [CrossRef] [Google Scholar]
  28. Fray, N., & Schmitt, B. 2009, Planet. Space Sci., 57, 2053 [NASA ADS] [CrossRef] [Google Scholar]
  29. Fulton, B. J., Petigura, E. A., Howard, A. W., et al. 2017, ApJ, 154, 109 [Google Scholar]
  30. Furi, E., & Marty, B. 2015, Nat. Geosci., 8, 515 [NASA ADS] [CrossRef] [Google Scholar]
  31. Gomes, G., & Ferraz-Mello, S. 2020, MNRAS, 494, 5082 [CrossRef] [Google Scholar]
  32. Guinan, E. F., & Engle, S. G. 2019, Res. Notes AAS, 3, 189 [CrossRef] [Google Scholar]
  33. Harris, G. J., Tennyson, J., Kaminsky, B. M., Pavlenko, Y. V., & Jones, H. R. A. 2006, MNRAS, 367, 400 [NASA ADS] [CrossRef] [Google Scholar]
  34. Hartmann, J. M., Boulet, C., Brodbeck, C., et al. 2002, J. Quant. Spectr. Rad. Transf., 72, 117 [NASA ADS] [CrossRef] [Google Scholar]
  35. Howett, C. J. A., Carlson, R. W., Irwin, P. G. J., & Calcutt, S. B. 2007, J. Opt. Soc. Am. B, 24, 126 [NASA ADS] [CrossRef] [Google Scholar]
  36. Iaroslavitz, E.,& Podolak, M. 2007, Icarus, 187, 600 [NASA ADS] [CrossRef] [Google Scholar]
  37. Karkoschka, E., & Tomasko, M. G. 2011, Icarus, 211, 780 [NASA ADS] [CrossRef] [Google Scholar]
  38. Kissel, A., Sumpf, B., Kronfeldt, H. D., Tikhomirov, B. A., & Ponomarev, Y. N. 2002, J. Mol. Spectr., 216, 345 [NASA ADS] [CrossRef] [Google Scholar]
  39. Kramida, A., Ralchenko, Y., eader, J., et al. 2019, NIST Atomic Spectra Database [Google Scholar]
  40. Langlois, S., Birbeck, T. P., & Hanson, R. K. 1994, J. Mol. Spectr., 167, 272 [CrossRef] [Google Scholar]
  41. Lellouch, E., Bézard, B., Fouchet, T., et al. 2001, A&A, 370, 610 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Lerot, C., Walrand, J., Blanquet, G., Bouanich, J. P., & Lepère, M. 2003, J. Mol. Spectr., 219, 329 [CrossRef] [Google Scholar]
  43. Lide, D. 2009, CRC Handbook of Chemistry and Physics, 90th, ed. D. R. Lide [Google Scholar]
  44. Lin, J.-F., Militzer, B., Struzhkin, V. V., et al. 2004, J. Chem. Phys., 121, 8423 [CrossRef] [Google Scholar]
  45. Lodders, K. 1999, J. Phys. Chem. Ref. Data, 28, 1705 [NASA ADS] [CrossRef] [Google Scholar]
  46. Lodders, K. 2002, ApJ, 577, 974 [NASA ADS] [CrossRef] [Google Scholar]
  47. Lodders, K. 2010, Formation and Evolution of Exoplanets (NewYork: John Wiley & Sons), 157 [Google Scholar]
  48. Lodders, K. 2019, ArXiv e-prints [arXiv:1912.00844] [Google Scholar]
  49. Madhusudhan, N., Nixon, M. C., Welbanks, L., Piette, A. A. A., & Booth, R. A. 2020, ApJ, 891, L7 [CrossRef] [Google Scholar]
  50. Mamajek, E., Prsa, A., Torres, G., et al. 2015, ArXiv e-prints [arXiv:1510.07674] [Google Scholar]
  51. Margolis, J. S. 1996, J. Quant. Spectr. Rad. Transf., 55, 823 [CrossRef] [Google Scholar]
  52. McKemmish, L. K., Yurchenko, S. N., & Tennyson, J. 2016, MNRAS, 463, 771 [NASA ADS] [CrossRef] [Google Scholar]
  53. Mlawer, E. J., Payne, V. H., Moncet, J. L., et al. 2012, Phil. Trans.: Math. Phys. Eng. Sci., 370, 2520 [Google Scholar]
  54. Morley, C. V., Knutson, H., Line, M., et al. 2017, ApJ, 153, 86 [NASA ADS] [CrossRef] [Google Scholar]
  55. Nemtchinov, V., Sung, K., & Varanasi, P. 2004, J. Quant. Spectr. Rad. Transf., 83, 243 [CrossRef] [Google Scholar]
  56. Niemann, H. B., Atreya, S. K., Carignan, G. R., et al. 1998, J. Geophys. Res.: Planets, 103, 22831 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  57. Otegi, J. F., Bouchy, F., & Helled, R. 2020, A&A, 634, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Pine, A. S. 1992, J. Chem. Phys., 97, 773 [NASA ADS] [CrossRef] [Google Scholar]
  59. Rey, M., Nikitin, A. V., & Tyuterev, V. G. 2017, ApJ, 847, 105 [CrossRef] [Google Scholar]
  60. Richard, C., Gordon, I. E., Rothman, L. S., et al. 2012, J. Quant. Spectr. Rad. Transf., 113, 1276 [NASA ADS] [CrossRef] [Google Scholar]
  61. Rinsland, C. P., Malathy Devi, V., Smith, M. A. H., et al. 2003, J. Quant. Spectr. Rad. Transf., 82, 343 [CrossRef] [Google Scholar]
  62. Robie, R., & Hemingway, B. 1995, US Geol. Survey Bull., 2131, 461 [Google Scholar]
  63. Rogers, L. A., & Seager, S. 2010, ApJ, 712, 974 [NASA ADS] [CrossRef] [Google Scholar]
  64. Rothman, L. S., Gordon, I. E., Barber, R. J., et al. 2010, J. Quant. Spectr. Rad. Transf., 111, 2139 [Google Scholar]
  65. Scheucher, M., Wunderlich, F., Grenfell, J. L., et al. 2020, ApJ, 898, 44 [CrossRef] [Google Scholar]
  66. Schwenke, D. W. 1998, Faraday Discuss., 109, 321 [NASA ADS] [CrossRef] [Google Scholar]
  67. Seager, S., Kuchner, M., Hier-Majumder, C. A., & Militzer, B. 2007, ApJ, 669, 1279 [NASA ADS] [CrossRef] [Google Scholar]
  68. Skaf, N., Bieger, M. F., Edwards, B., et al. 2020, ApJ, 160, 109 [CrossRef] [Google Scholar]
  69. Snellen, I., De Kok, R., De Mooij, E., & Albrecht, S. 2010, Nature, 465, 1049 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  70. Sousa-Silva, C., Al-Refaie, A. F., Tennyson, J., & Yurchenko, S. N. 2014, MNRAS, 446, 2337 [NASA ADS] [CrossRef] [Google Scholar]
  71. Sromovsky, L. A., Fry, P. M., & Kim, J. H. 2011, Icarus, 215, 292 [NASA ADS] [CrossRef] [Google Scholar]
  72. Stull, D. 1947, Ind. Eng. Chem., 39, 540 [CrossRef] [Google Scholar]
  73. Teanby, N. A., Fletcher, L. N., Irwin, P. G. J., Fouchet, T., & Orton, G. S. 2006, Icarus, 185, 466 [NASA ADS] [CrossRef] [Google Scholar]
  74. Tsiaras, A., Waldmann, I. P., Zingales, T., et al. 2018, ApJ, 155, 156 [Google Scholar]
  75. Tsiaras, A., Waldmann, I., Tinetti, G., Tennyson, J., & Yurchenko, S. 2019, Nat. Astron., 3, 1086 [NASA ADS] [CrossRef] [Google Scholar]
  76. Valencia, D., Guillot, T., Parmentier, V., & Freedman, R. 2013, ApJ, 775, 10 [NASA ADS] [CrossRef] [Google Scholar]
  77. Venot, O., Cavalié, T., Bounaceur, R., et al. 2020, A&A, 634, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  78. Wagner, W., & Pruss, A. 1993, J. Phys. Chem. Ref. Data, 22, 783 [CrossRef] [Google Scholar]
  79. Wagner, W., Saul, A., & Pruss, A. 1994, J. Phys. Chem. Ref. Data, 23, 515 [NASA ADS] [CrossRef] [Google Scholar]
  80. Wakeford, H. R., Sing, D. K., Kataria, T., et al. 2017, Science, 356, 628 [NASA ADS] [CrossRef] [Google Scholar]
  81. Wang, D., Lunine, J. I., & Mousis, O. 2016, Icarus, 276, 21 [NASA ADS] [CrossRef] [Google Scholar]
  82. Welbanks, L., Madhusudhan, N., Allard, N. F., et al. 2019, ApJ, 887, L20 [CrossRef] [Google Scholar]
  83. Wilzewski, J. S., Gordon, I. E., Kochanov, R. V., Hill, C., & Rothman, L. S. 2016, J. Quant. Spectr. Rad. Transf., 168, 193 [CrossRef] [Google Scholar]
  84. Wong, M. H., Mahaffy, P. R., Atreya, S. K., Niemann, H. B., & Owen, T. C. 2004, Icarus, 171, 153 [NASA ADS] [CrossRef] [Google Scholar]
  85. Woodfield, B. F., Shapiro, J. L., Stevens, R., et al. 1999, J. Chem. Thermodyn., 31, 1573 [CrossRef] [Google Scholar]
  86. Yurchenko, S. N. 2015, J. Quant. Spectr. Rad. Transf., 152, 28 [CrossRef] [Google Scholar]
  87. Zahnle, K. J., & Marley, M. S. 2014, ApJ, 797, 41 [NASA ADS] [CrossRef] [Google Scholar]
  88. Zeng, L., Jacobsen, S. B., Sasselov, D. D., et al. 2019, Proc. Natl. Acad. Sci. USA, 116, 9723 [NASA ADS] [CrossRef] [Google Scholar]
2

Not all the confirmed planets have a measured radius and/or mass, thus the total count in the histograms is less than 4266.

4

O. Venot, priv. comm. The comparison was made using the “nominal model” described in Sect. 4.1

5

NH4Cl optical constants were not available, and we used the parameters from NH4 SH. Both molecules show similar strong absorptions due to NH and minor contributions from Cl− or HS− (NIST chemical WebBook and Howett et al. 2007).

6

(C/O) (Lodders 2019).

8

J. Leconte, priv. comm.

All Tables

Table 1

Parameters of K2-18b and its star.

In the text
Table 2

Collision-induced absorption references.

In the text
Table 3

Species cross-section parameters.

In the text
Table 4

Line broadening references.

In the text
Table 5

Species included in thermochemical calculations.

In the text
Table 6

Model grid.

In the text
Table 7

Best fits against Benneke et al. (2019) dataset as a function of internal temperature.

In the text
Table 8

Comparison of the Bayesian log-evidence of TauREx 3 against Benneke et al. (2019) dataset for different models.

In the text

All Figures

thumbnail Fig. 1

Mass-radius distribution of the 42662 confirmed exoplanets to date (17 June 2020). The scatter plot shows only the 454 planets for which the mass and the radius are known within ± 15%. The mass-radius curves for an Earth-like interior (orange) and a pure water ice planet (blue) are taken from Seager et al. (2007).
In the text
thumbnail Fig. 2

Best fit of our model to the dataset of Benneke et al. (2019) at nominal irradiation (1368 W m−2). Black: K2, HST and Spitzer data from Benneke et al. (2019). Grey: HST data from Tsiaras et al. (2019). Red: Exo-REM transmission spectrum at 175 (Z/H)⊙, with an offset of the 105-Pa level of −5 km. The χ2 of this spectrum against Benneke et al. (2019) data is indicated in parentheses.
In the text
thumbnail Fig. 3

Contributions of the absorbing species to our best-fit transmission spectrum at nominal irradiation within the K2 to Spitzer spectral range. The H2O cloud is not forming in this case, so there is no cloud contribution. The spectral contribution of individual species takes the CIA and Rayleigh scattering (represented as a dotted curve) into account. The spectral contribution of FeH, TiO, and VO are not represented here as they are insignificant.
In the text
thumbnail Fig. 4

Temperature profile of our best fit atmospheric model against the dataset from Bennekeet al. (2019) at nominal irradiation. Solid line: temperature profile. Dotted lines: condensation profiles of selected species. Red line: convective layers. Blue dot: H2 O ice-Ih–liquid–gas triple point.
In the text
thumbnail Fig. 5

VMR of our best fit atmospheric model against the dataset from Benneke et al. (2019) at nominal irradiation (175 (Z/H)⊙, Kzz = 106 cm2 s−1). For clarity, only absorbing species and N2 are represented. FeH, TiO, and VO are not represented as their respective abundances are, respectively, < 10−12, < 10−35, and < 10−36 within our pressure grid. The bump of the PH3 VMR at 200 kPa is due to the chemical equilibrium between PH3 and P2 peaking in favour of P2 at this pressure.
In the text
thumbnail Fig. 6

Goodness of fit of our nominal models, varying only irradiation and metallicity. Left: data from Benneke et al. (2019), including HST, K2, and Spitzer data. Right: data from Tsiaras et al. (2019), including HST, K2, and Spitzer data. The white colour corresponds to the value of χ2 of 21.36, indicating a fit of the datasets at the 1σ confidence level. Accordingly, the blue colour indicates overfit and the red colour underfit of the datasets. The solid lines represents the 1 and 2σ confidencelevels of our fits. The dashed and dotted lines represents, respectively, the 1 and 3σ uncertainties on K2-18b irradiation.
In the text
thumbnail Fig. 7

Effect of metallicity on the transmission spectrum at nominal irradiation. The corresponding value of χ2 against Benneke et al. (2019) data is indicated in parentheses.
In the text
thumbnail Fig. 8

Effect of irradiation on the transmission spectrum at 10 (Z/H)⊙. The corresponding value of χ2 against Benneke et al. (2019) data is indicated in parentheses.
In the text
thumbnail Fig. 9

H2O condensation pressure level of our nominal models, varying only irradiation and metallicity. White squares: no cloud formation. Strips: simulation where the H2 O clouds are optically thick (normal optical depth τ > 1 at 1 μm). Dots: simulations where H2O condenses into its liquid phase. Dotted line: 3σ lower limit of K2-18b irradiation.
In the text
thumbnail Fig. 10

H2O cloud and metallicity effect on the transmission spectrum, at nominal irradiation. The χ2 against Benneke et al. (2019) data is indicated in parentheses.
In the text
thumbnail Fig. 11

H2O cloud and metallicity effect on the Bond albedo of K2-18b, at nominal irradiation. Cloud condensation occurs only above 200 (Z/H)⊙.
In the text
thumbnail Fig. 12

H2O cloud particle radius effect on the transmission spectrum, at 300 (Z/H)⊙ and at nominal irradiation. The χ2 against Benneke et al. (2019) data is indicated in parentheses.
In the text
thumbnail Fig. 13

Effects of Tp, int on K2-18b, at 175 (Z/H)⊙, nominal irradiation, and Kzz = 106 cm2 s−1. Left: temperature profiles. The red area represents the convective layers. Right: VMR of the most abundant absorbers. Dotted: Tp, int = 45 K. Dashed: Tp, int = 80 K. Solid: Tp, int = 115 K.
In the text
thumbnail Fig. 14

Effects of Kzz on the VMR of K2-18b main C, N, O bearing species, at 175 (Z/H)⊙. Dashed: Kzz = 106 cm2 s−1. Dotted-dashed: Kzz = 108 cm2 s−1. Solid: Kzz = 1010 cm2 s−1. Dotted: thermochemical equilibrium.
In the text
thumbnail Fig. 15

Chemical evolution of the upper atmosphere of K2-18b for the most abundant C- and O-bearing species over time, depending on the value of Kzz, at 175 (Z/H)⊙. Dotted: Kzz = 106 cm2 s−1. Dashed: Kzz = 108 cm2 s−1. Solid: Kzz = 1010 cm2 s−1.
In the text
thumbnail Fig. 16

Goodness of fit of our models for Tp, int = 80 K, Kzz = 106 cm2 s−1 and nominal irradiation, as a function of the metallicity and C/H, compared to the data from Benneke et al. (2019; left) and Tsiaras et al. (2019; right). The white colour corresponds to the value of χ2 = 21.36, indicating a fit at the 1σ confidence level. The solid black line represents the 1σ confidence level. The dotted lines indicates the C/O ratio. The solid red line indicates which species dominates on average the transmission spectrum in the 1.355–1.415 μm range. The cases where C/O > (C/O)⊙ was not explored are represented as black rectangles.
In the text
thumbnail Fig. 17

Best fit against Benneke et al. (2019) dataset with transit spectra in which the spectral contributions of selected species has been removed. Blue: H2 O contribution of the transit spectrum displayed in Fig. 3, with an offset of the 105 -Pa level of +98 km. Cyan: H2O, CO, CO2, and NH3 total contribution of the transit spectrum displayed in Fig. 3, with an offset of the 105 -Pa level of + 57 km. Black: K2, HST and Spitzer from Benneke et al. (2019). The dotted red line is our best-fit spectrum with all absorptions (see Sect. 4). The χ2 of these spectra against the data is indicated in parentheses.
In the text
thumbnail Fig. 18

TauREx 3 atmospheric retrieval posterior distributions against the dataset from Benneke et al. (2019), with planetary radius (R♃), temperature (K), decimal logarithm of the VMR of H20, CH4, CO, CO2, N2, and NH3, and decimallogarithm of the cloud top pressure (Pa) as free parameters.
In the text
thumbnail Fig. 19

Goodness of fit of our models at Tp, int = 200 K and nominal irradiation as a function of metallicity and Kzz, compared tothe data from Benneke et al. (2019; left) and Tsiaras et al. (2019; right). The white colour corresponds to the value of χ2 = 21.36, indicating a fit at the 1σ confidence level. The best fit against Benneke et al. (2019) data (χ2 = 20.64) is located at (Z/H)⊙ = 125 and Kzz = 1010 cm2 s−1. The best fit against Tsiaras et al. (2019) data (χ2 = 16.56) is located at (Z/H)⊙ = 200 and Kzz = 105 cm2 s−1.
In the text
thumbnail Fig. 20

VMR of selected species of our best-fit model against Benneke et al. (2019) data at Tp, int = 200 K and nominal irradiation.
In the text
thumbnail Fig. 21

VMR of selected species at Tp, int = 80 K, Kzz = 106 cm2 s−1, nominal irradiation, 125 (Z/H)⊙, and 30 (C/H)⊙.
In the text
thumbnail Fig. 22

Contributions within some spectral ranges of a selection of absorbing species to the transmission spectrum of our nominal model (175 (Z/H)⊙) at nominal irradiation and a resolution of 0.5 cm−1 (resolving power of 20 000 at 1 μm). The spectral contribution of individual species takes the CIA and Rayleigh scattering (dotted curve) into account.
In the text
thumbnail Fig. 23

Atmospheric metallicity as a function of mass for a selection of planets. The Jupiter (J), Saturn (S), Uranus (U), and Neptune (N) C/H are taken respectively from Wong et al. (2004), Fletcher et al. (2009), Sromovsky et al. (2011) and Karkoschka & Tomasko (2011) (values compiled in Atreya et al. 2018). The data for GJ436b and WASP-127b are taken respectively from Morley et al. (2017) and Skaf et al. (2020), while the data for all the other planets are taken from Welbanks et al. (2019). The black area represents the 1σ C/H against mass trend for Jupiter, Saturn, Uranus, and Neptune resulting from a linear fit. The fit is log10 ((Z/H)∕(Z/H)⊙) = −1.10 ± 0.21log10(M[kg]) + 30.5 ± 5.5. The error bars of our results correspond to the simulations that were within the 1σ confidence interval against Tsiaras et al. (2019) data in the “free C/H” scenario. Our C/H results were slightly shifted to the right for clarity. The O/H ratios in the other works were estimated assuming [H2 ] + [He] + [H2O] + [CH4 ] = 1 and a solar He/H.
In the text
thumbnail Fig. 24

Summary of our results and results of other teams on K2-18b metallicity. The error bars of our results corresponds to the simulations that were within the 1σ confidence level against Tsiaras et al. (2019) data. The O/H and C/H ratio of the other works were estimated assuming [H2 ] + [He] + [H2O] + [CH4 ] = 1 and a solar He/H. The C/H upper limit of Benneke et al. (2019) and Madhusudhan et al. (2020) corresponds respectively to the 2σ and the 99% upper limit. The dotted line represents the metallicity expected for the mass of K2-18 b using the linear fit in Fig. 23, based on solar-system C/H ratios.
In the text
thumbnail Fig. 25

Comparison of K2-18b transmission spectrum best-fits for three different scenarios. Red: 175 (Z/H)⊙, 175 (C/H)⊙, Kzz = 106 cm2 s−1, Tp, int = 80 K. Green: 125 (Z/H)⊙, 30 (C/H)⊙, Kzz = 106 cm2 s−1, Tp, int = 80 K. Blue: 125 (Z/H)⊙, 125 (C/H)⊙, Kzz = 1010 cm2 s−1, Tp, int = 200 K. The χ2 against Benneke et al. (2019) data is indicated in parentheses.
In the text
thumbnail Fig. A.1

Effect of the H2O CIA at 175 (Z/H)⊙ on the simulated transmission spectrum and temperature profile of K2-18b. Left: temperature profiles. Solid black: nominal model. Dotted-dashed black: nominal model at 0.8 time the nominal irradiance. Red: nominal model without H2 O CIA. Dotted-dashed red: nominal model at 0.8 times the nominal irradiance and without H2 O CIA. The phase diagram of H2O is represented as dotted blue lines, the dot corresponding to the H2O triple point. Right: transmission spectra at 0.8 time the nominal irradiation, with (black) and without (red) H2 O CIA.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.