Thermochemistry and vertical mixing in the tropospheres of Uranus and Neptune: How convection inhibition can affect the derivation of deep oxygen abundances

TL;DR

This paper investigates how convection inhibition due to water condensation affects the thermal structure and oxygen abundance estimates in Uranus and Neptune's atmospheres, using updated models and observations.

Contribution

It introduces a revised atmospheric model that accounts for convection inhibition, leading to more accurate deep oxygen abundance estimates in ice giant planets.

Findings

Convection inhibition significantly alters the thermal profile.

Revised models suggest different oxygen enrichment levels.

Current limitations of the modeling approach are discussed.

Abstract

Thermochemical models have been used in the past to constrain the deep oxygen abundance in the gas and ice giant planets from tropospheric CO spectroscopic measurements. Knowing the oxygen abundance of these planets is a key to better understand their formation. These models have widely used dry and/or moist adiabats to extrapolate temperatures from the measured values in the upper troposphere down to the level where the thermochemical equilibrium between H$_2$O and CO is established. The mean molecular mass gradient produced by the condensation of H$_2$O stabilizes the atmosphere against convection and results in a vertical thermal profile and H$_2$O distribution that departs significantly from previous estimates. We revisit O/H estimates using an atmospheric structure that accounts for the inhibition of the convection by condensation. We use a thermochemical network and the latest observations of CO in Uranus and Neptune to calculate the internal oxygen enrichment required to satisfy both these new estimates of the thermal profile and the observations. We also present the current limitations of such modeling.

Tables

Table 1. Table 1: Summary of observational data.

Planet	Jupiter	Saturn	Uranus	Neptune
$y_{\mathrm{He}}$(1)	0.1359$\pm$0.0027	0.118$\pm$0.025	0.152$\pm$0.033	0.149${}_{-0.022}{}^{}{}_{}{}^{+0.017}$
$y_{\mathrm{CH}}$4(2)	(2.04$\pm$0.50)$\times {10}^{-3}$	(4.7$\pm$0.2)$\times {10}^{-3}$	0.01-0.05	0.01-0.05
$y_{\mathrm{CO}}$(3)	(1.0$\pm$0.2)$\times {10}^{-9}$	10-9	2.1$\times {10}^{-9}$	(0.20$\pm$0.05)$\times {10}^{-6}$
$T^{\left(4\right)}$	207.1	132.9	102.9	93.1
$F^{\left(5\right)}$	5.44$\pm$0.43	2.01$\pm$0.14	0.042$\pm$0.047	0.433$\pm$0.046
$T_{\mathrm{int}}^{\left(6\right)}$	99.0	77.2	29.3	52.6

Equations

$${\mathrm{H}}_2\mathrm{O}+{\mathrm{CH}}_4=\mathrm{CO}+3{\mathrm{H}}_2.$$

$$X_{\mathrm{He}}=y_{He}^{\mathrm{top}}/(y_{{\mathrm{H}}_2}^{\mathrm{top}}+2y_{{\mathrm{CH}}_4}^{\mathrm{top}})\times (X_{\mathrm{H}}-2X_{\mathrm{O}})/2,$$

$$X_{\mathrm{C}}=X_{\mathrm{He}}\times y_{{\mathrm{CH}}_4}^{\mathrm{top}}/y_{\mathrm{He}}^{\mathrm{top}},$$

$$y_{{\mathrm{H}}_2\mathrm{O}}^{\mathrm{bottom}}=X_{\mathrm{O}}/(X_{\mathrm{H}}/2+X_{\mathrm{He}}-X_{\mathrm{C}}),$$

$$T_{i+1}=T_i\times \exp \left({\nabla }_{\mathrm{dry}}\frac{\ln (p_{i+1})}{\ln (p_i)}\right)$$

$${\nabla }_{\mathrm{dry}}=\frac{R}{\mu c_{\mathrm{p}}},$$

$${\nabla }_{\mathrm{moist}}=\frac{p}{p-p_{\text{v}}}\frac{(1-q_{\text{v}})R_{\text{a}}+\frac{q_{\text{v}}{L}_{\text{v}}}{T}}{q_{\text{v}}{c}_{\text{p,v}}+q_{\text{a}}{c}_{\text{p,a}}+q_{\text{c}}{c}_{\text{p,c}}+q_{\text{v}}\frac{L_{\text{v}}}{T}\frac{p}{p-p_{\text{v}}}\frac{\mathrm{d}\ln p_{\text{s}}}{\mathrm{d}\ln T}},$$

$$\varpi q_{\mathrm{v}}\frac{\mathrm{d}\ln p_{\mathrm{s}}}{\mathrm{d}\ln T}>1,$$

$${\nabla }_{\text{r}}=\frac{3}{16}\frac{\kappa p{T}_{\text{int}}^4}{g{T}^4},$$

$$f_{\mu }=\frac{1}{1+({\mu }_{\mathrm{a}}/M_{{\mathrm{H}}_2\mathrm{O}}-1)y_{{\mathrm{H}}_2\mathrm{O}}^{\mathrm{mass}}}$$

$$y_{{\mathrm{H}}_2\mathrm{O}}^{\mathrm{mass}}=\frac{M_{{\mathrm{H}}_2\mathrm{O}}}{\mu }{y}_{{\mathrm{H}}_2\mathrm{O}}^{\mathrm{bottom}}$$

$$K_{zz}\simeq {\left(\frac{k_{B}F}{\rho m{c}_{\mathrm{p}}}\right)}^{1/3}\times H.$$

Full text

††journal: Icarus

Received 9 December 2016 ; Revised 3 March 2017 ; Accepted 10 March 2017

DOI: 10.1016/j.icarus.2017.03.015

1 Introduction

One of the great mysteries in the Solar System is how the gas and ice giant planets formed from the protoplanetary disk 4.5 billion years ago. This question is even more relevant regarding ice giants after the discovery of the commonality of Neptune-class planets among the exoplanets detected by Kepler (Batalha et al., 2013; Fressin et al., 2013). Two scenarios have been proposed regarding the formation of giant planets: disk gravitational instability (Boss, 1997) and core accretion (Pollack et al., 1996). These scenarios differ not only in the time required to form planets (a few hundred years vs. several million years, respectively), but also in the final composition of the planets’ interiors. While gravitational instability should result in $\sim$solar abundances of heavy elements (except if a significant external source of heavy elements is incorporated after the planet formation), core accretion formation should lead to an enrichment in heavy elements increasing with heliocentric distance. The level of enrichment then would depend on how the ices of the planetesimal that formed the cores of these planets condensed. Here again, several competing scenarios exist: condensation in amorphous ices (Bar-Nun et al., 1988; Owen et al., 1999) and clathration (Lunine and Stevenson, 1985; Gautier et al., 2001; Hersant et al., 2004; Gautier and Hersant, 2005).

While measuring the D/H ratio in the giant planets (Lellouch et al., 2001; Feuchtgruber et al., 2013) gives us an insight on the origin and condensation temperature of the protoplanetary ices in the outer Solar System (Hersant et al., 2001), measuring the heavy element abundances can help constraining the ice condensation processes. Enrichment in heavy elements has been observed by the Galileo probe in Jupiter’s troposphere, with an average enrichment factor of $\sim$4$\pm$2 in C, N, S, Ar, Kr and Xe, except for O which was found depleted (Niemann et al., 1998; Atreya et al., 1999; Mahaffy et al., 2000; Wong et al., 2004; Owen and Encrenaz, 2006). The O measurement may only be a lower limit because the Galileo probe most likely entered a dry zone of Jupiter and did not reach the levels where H2O is well-mixed (Wong et al., 2004), though alternative scenarios explaining an oxygen depletion exist (Lodders, 2004; Mousis et al., 2012). The global enrichment in heavy elements reported by Galileo favors the core accretion scenario for Jupiter. This formation model, if applied to the other giant planets, predicts enrichment factors that increase with increasing heliocentric distance: $\sim$7, 30 and 50 for Saturn, Uranus, and Neptune, respectively (Owen and Encrenaz, 2003). At Saturn, the C, N, P and S abundances have been measured (Fletcher et al., 2009; Hersant et al., 2008; Mandt et al., 2015) and are also found enriched compared to the solar value. At Uranus and Neptune, the information is more sparse, with only the C abundance has constraints (Baines et al., 1995; Sromovsky and Fry, 2008; Karkoschka and Tomasko, 2009).

The key measurement that would enable differentiating the condensation processes of the planetesimal ices and hence how other heavy elements were trapped is the deep water abundance. Indeed, the clathration scenario needs a larger amount of water at the time of the condensation of the planetesimal ices than the amorphous ice scenario (Owen, 2007), especially if the efficiency of the clathration process was lower than 100%. While Galileo probably failed to measure the Jovian deep water abundance below the water cloud (Niemann et al., 1998) because it entered a dry hot spot and probably did not reach the base of the water cloud, Juno (Matousek, 2007) should shed light on this long lasting question. However, there is no such mission planned in the near future to measure the deep water abundance in the other giant planets. A few mission concepts are being developed for Saturn (Mousis et al., 2012, 2014; Atkinson et al., 2016) and the ice giants (Arridge et al., 2012, 2014; Masters et al., 2014), but these challenging missions require probes that would have to survive high pressures to reach below the water cloud (up to $\sim$100 bars in the ice giants). In principle, cm-waves can probe down to several tens of bars and could probe below the water cloud (de Pater et al., 1991; Hofstadter and Butler, 2003; Klein and Hofstadter, 2006). However, the opacity at these levels can also be caused by NH3 and H2S and the degeneracy is difficult to waive (e.g., de Pater et al. 2005). Therefore, it is important to find other ways to constrain the deep water abundance. One interesting way is to take advantage of the chemical quenching of species like CO. The observed abundance of CO is in chemical disequilibrium and is inherited from deeper layer, where its abundance is in thermochemical equilibrium with water, through the combination of fast vertical mixing and slow chemical kinetics.

Following the first detection of CO in the atmosphere of Jupiter by Beer (1975), a simple model was proposed to constrain the deep water abundance in Jupiter by studying the tropospheric thermochemistry and vertical transport, and in particular the following thermochemical equilibrium reaction:

[TABLE]

This model, first developed by Fegley and Prinn (1988) is based on the approximation that the tropospheric mole fraction of CO is fixed at a so-called “quench” level, where the chemical timescale of the conversion of CO into H2O becomes longer than the timescale for its vertical transport by convection. This kind of model relies on the determination of the rate-limiting reaction of the conversion scheme. Therefore, assuming the kinetics of this rate-limiting reaction is known, the kinetics of the whole conversion scheme is constrained and the measured upper tropospheric mole fraction of CO can be linked to the deep water abundance. Prinn and Barshay (1977) initially proposed this reaction to be H2CO$+$H2$\rightarrow$CH3$+$OH. Later, Yung et al. (1988) proposed a two-step reaction scheme in which the rate-limiting reaction was H$+$H2CO$+$M$\rightarrow$CH3O$+$M. This “quench-level” approximation model was then used by Lodders and Fegley (1994) to constrain the atmospheric O/H ratio in all giant planets. However, Smith (1998) showed that the assumptions of these modelers on diffusion timescales were incorrect. His work was then applied by Bézard et al. (2002) to Jupiter, by Visscher and Fegley (2005) and Cavalié et al. (2009) to Saturn, and by Luszcz-Cook and de Pater (2013) to Neptune using either the Prinn and Barshay (1977) or Yung et al. (1988) limiting reactions, which renders any comparison between these results hazardous.

Visscher et al. (2010) first applied a comprehensive thermochemical and diffusive transport model to the troposphere of Jupiter. They have evaluated the Jovian deep water mole fraction to be (0.25-6.0)$\times 10^{-3}$, corresponding to an enrichment of 0.3-7.3 times the protosolar abundance (9.61$\times 10^{-4}$). Similar thermochemical models have been applied ever since to investigate the thermochemistry in exoplanets (Visscher and Moses, 2011; Moses et al., 2011, 2013; Venot et al., 2012, 2014, 2015) and in Solar System giants (Cavalié et al., 2014; Mousis et al., 2014; Wang et al., 2016). While all these models improve on the modeling of deep tropospheric chemistry compared to the “quench-level” approximation studies, they still have to rely on assumptions made on tropospheric vertical mixing and temperatures. Interestingly, recent theoretical work enables progress in the determination of these quantities:

Vertical mixing: vertical mixing in the troposphere of giant planets is caused by convection, and it is usually modeled by an eddy diffusion coefficient. It is estimated within an order of magnitude from the mixing length theory (Stone, 1976). Following rotating tank experiments, a recent paper by Wang et al. (2015) proposes a new formulation to estimate this coefficient with a narrower error bar. The authors even show that this vertical mixing is latitude-dependent, with a stronger magnitude at low latitudes.

• 2.

Tropospheric temperatures extrapolation: Until now, dry and/or moist adiabats have been used to extrapolate the thermal profiles of Solar System giant planets from observations around 1-bar to deeper levels (e.g., Luszcz-Cook and de Pater 2013). In a new paper, Leconte et al. (2017) propose a new criterion to compute the thermal gradient in the giant planet tropospheres. It takes into account not only dry and moist processes, but also the effect of the mean molecular weight gradient associated with the condensation of H2O and can produce a jump in temperature at the H2O condensation level that is caused by a thin radiative layer.

Both are composition-dependent and remain therefore quite uncertain. In addition, the possible radiative gradient at the H2O condensation level is poorly constrained and the CH4 tropospheric abundance is only known within a factor of two in Uranus and Neptune (Baines et al., 1995; Sromovsky and Fry, 2008; Karkoschka and Tomasko, 2009). We have therefore chosen to study the parameter space – which consists of the following four parameters: O/H, C/H, temperature, $K_{zz}$ – with a thermochemical and diffusion model of Uranus and Neptune to better evaluate the parameter space that is compliant with observations of CO in their upper tropospheres.

In this paper, we apply to the Solar System Ice Giants the thermochemical and diffusion model of Venot et al. (2012) along with the prescriptions of Leconte et al. (2017) to compute tropospheric thermal profiles. We present our models in Section 2, and review the various observational constraints on composition and temperature in Section 3. We detail the results of the 4D parameter space simulations and our nominal case in Section 4. We discuss our results in Section 5, and detail other sources of uncertainties that currently limit these kind of models in Section 6. We give our conclusions in Section 7.

2 Models

We have adapted the thermochemical and diffusion model of Venot et al. (2012), initially developed for the atmospheres of warm exoplanet atmospheres, to the giant planet tropospheres to constrain their deep oxygen abundance from CO observations. Obtaining the tropospheric composition from thermochemical and diffusion calculations requires some knowledge on the troposphere temperature, composition, and vertical mixing. In the following sections, we illustrate how we have run sequentially various steps to eventually link the observed tropospheric CO abundances to the deep oxygen abundances in Uranus and Neptune. These steps include:

temperature profile extrapolation to deep levels where thermochemistry prevails, assuming a given composition of the main compounds (H2, He, CH4 and H2O) 2. 2.

thermochemical equilibrium calculation using the temperature profile and a given elemental composition, to derive an initial composition state for the next round of computations 3. 3.

thermochemistry and diffusion calculations using the temperature profile and initial composition state of step 1 and step 2 and assuming a given vertical mixing $K_{zz}$, to obtain the steady state tropospheric composition 4. 4.

$K_{zz}$ assumption cross-check with theoretical estimates

To run the thermochemistry and diffusion model, we need at minimum to set as initial conditions the internal elemental abundances of H, He, C, N, and O, a thermal profile, and an eddy diffusion coefficient. The case of N will not be discussed any further as nitrogen chemistry has no significant impact on carbon and oxygen chemistries (Cavalié et al., 2014), the deep oxygen abundance is an assumption of the model, and other elemental abundances can be estimated a priori.

2.1 Step 0 - Estimating the internal composition by neglecting chemistry

As will be shown in what follows, the computation of the thermal profile needs a priori knowledge of the deep composition in terms of H2, He, CH4 and H2O. To estimate these molecular abundances as well as the elemental abundances of H, He and C, we start from the upper tropospheric abundances of H2, He, CH4, and an initial condition on the deep O abundance. Assuming that chemistry and the abundance of compounds other than H2, He, CH4, and H2O, can be neglected, and that the upper tropospheric abundance of H2O is negligible (because of condensation), we can derive a set of equations to link the upper tropospheric and deep tropospheric composition:

[TABLE]

[TABLE]

[TABLE]

where $X_{i}$ is the abundance of element $i$, $y_{i}^{\mathrm{top}}$ is the mole fraction of compound $i$ in the upper of the troposphere (but below the CH4 cloud), and $y_{i}^{\mathrm{bottom}}$ is the mole fraction of compound $i$ in the deep troposphere. By convention, $X_{\mathrm{H}}=10^{12}$.

In practice, we have measurements from which we can determine $y_{\mathrm{He}}^{\mathrm{top}}$, $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$, and $y_{\mathrm{H}_{2}}^{\mathrm{top}}$ (see Section 3), and we assume a deep oxygen abundance $X_{\mathrm{O}}$. We can thus derive an estimate of $X_{\mathrm{He}}$, $X_{\mathrm{C}}$, and $y_{\mathrm{H}_{2}\mathrm{O}}^{\mathrm{bottom}}$, which are needed for the preparation of the thermochemistry and diffusion computations, by using sequentially equations 2, 3, and 4. We control the validity of these assumptions on $X_{\mathrm{He}}$ and $X_{\mathrm{C}}$ a posteriori, by checking that the final upper tropospheric abundances of He and CH4 fit the observations.

2.2 Step 1 - Extrapolating tropospheric temperatures

Formation models of the giant planets of the Solar System predict enrichment factors for heavy elements that increase together with heliocentric distance. The precise values for the O enrichment factor depend on the ice condensation scenario (amorphous ices vs. clathration). But even in the case where the O enrichment should be the lowest, i.e., in the scenario where ices condense and trap volatiles in amorphous form, the expected enrichment factors are 4, 7, 30, and 50 (respectively) for Jupiter, Saturn, Uranus, and Neptune (respectively), according to Owen and Encrenaz (2003). The high abundance of tropospheric CO in Neptune could even indicate that O is enriched by a factor $>$100 (Lodders and Fegley, 1994; Luszcz-Cook and de Pater, 2013). Actually, above an enrichment factor of $\sim$100-150 in Uranus and Neptune, we expect to see the abundance of O exceeding the abundance of He. Given the molecular mass of the main O carrier (i.e., water) compared to the molecular mass of the main atmospheric constituents (i.e., H2 and He), a significant mean molecular mass gradient $\nabla_{\mu}$ is expected around the layers where water condenses. Consequently, the Ledoux stability criterion (Ledoux, 1947; Sakashita and Hayashi, 1959) must be used instead of the Schwarzschild stability criterion in the computations of the temperature gradient $\nabla_{\textrm{T}}$, as double diffusive processes may appear. The temperature gradient can be significantly affected, as shown by Guillot (1995) in the case of the condensation of methane.

Hereafter, we will especially consider the cases of Uranus and Neptune, as these planets are the most “symptomatic” cases of mean molecular weight discontinuities around the water condensation level, because of the expected large oxygen enrichments.

2.2.1 Methodology

The thermochemistry and diffusion computations require a thermal profile for the troposphere down to the region where thermochemistry prevails over vertical diffusion, i.e. below the CO quench-level in our case.

Previous studies of the Solar System giant planets thermochemistry have used a variety of dry and/or moist adiabats, based on the application of the Schwarzschild and Härm (1958) stability criterion, to compute the tropospheric thermal structures (e.g., Lodders and Fegley 1994; Luszcz-Cook and de Pater 2013). In a new paper, Leconte et al. (2017) show that, in hydrogen-rich atmospheres, the mean molecular weight gradient around the cloud base (of a species heavier than H2) can be strong enough to stabilize the atmosphere against convection (i.e. inhibit moist convection) in this region. Although it shares some similarities with the processes described by Ledoux that work for non condensable species (Ledoux, 1947; Sakashita and Hayashi, 1959), this mechanism can suppress moist convection when the enrichment in the condensable species is higher than a critical threshold. Because the interior still needs to release its energy, a stable radiative layer with a steep temperature gradient develops. Guillot (1995) has already shown that such an effect can be produced by CH4 in Uranus and Neptune. We can expect this effect to happen with H2O too in Uranus and Neptune, in the transition zone between the H2O-rich region (deep troposphere) and the H2O-poor region (upper troposphere), and we have therefore implemented their stability criterion when extrapolating thermal profiles. We will refer to “3-layer profiles” for thermal profiles extrapolated using the prescription of Leconte et al. (2017).

In the computations we present in this paper, we start from the temperature measured at the 2-bar level (see Table 1). We will explore various cases for extrapolating deep tropospheric thermal profiles: (i) dry processes leading to a dry adiabat, (ii) latent heat effects leading to a moist adiabat, (iii) latent heat and mean molecular weight effects leading to what we will refer to as “3-layer profiles”. Finally, because magnitude of the temperature jump in the radiative layer of the 3-layer profiles is quite uncertain, we derive an extreme case (iv) from the Leconte et al. (2017) formulations in which the temperature jump of the radiative layer is limited by the deep water abundance.

2.2.2 Dry adiabat

Starting from the 2-bar level temperature, we extrapolate to deeper levels using the Schwarzschild and Härm (1958) criterion, which translates into the following:

[TABLE]

and

[TABLE]

where $\mu$ is the mean molecular weight, and $R$ is the ideal gas constant. We use the temperature-dependent expressions given by the NIST for the specific heat capacities of He, CH4 and H2O. For H2, we use the temperature- and pressure-dependent data from the Cryogenic Data Handbook111https://www.bnl.gov/magnets/staff/gupta/cryogenic-data-handbook/Section3.pdf. We assume equilibrium hydrogen, in agreement with observations of Baines et al. (1995) of the upper tropospheres of Uranus and Neptune. Because we start from the 2-bar level in both Uranus and Neptune, the lowest temperatures we consider are sufficiently high so that we do not have to consider the effects of CH4 condensation (Guillot, 1995). We then compute the next levels with very small steps in $\ln(p)$ ($\sim 10^{-4}$), with $p$ in Pa.

This kind of profile is obviously the least realistic given the expected high abundance of H2O in the interiors of Uranus and Neptune.

2.2.3 Moist adiabat

In the region where H2O condenses, latent heat is released and the gas is heated. Convection is thus reinforced. Adding this effect to the previous stability criterion results in replacing $\nabla_{\mathrm{dry}}$ in equation (5) by the expression given in Leconte et al. (2013):

[TABLE]

where subscripts i$$=a,v,c refer to non-condensable gas, condensable gas, and condensed material (solid or liquid), respectively. $p_{i}$ is the partial pressure, $M_{i}$ the molar mass, $q_{i}$ is the mass mixing ratio and is the ratio of the mass of component $i$ over the mass of the gas, $c_{\textrm{p},i}$ the mass heat capacity, and $R_{i}=R/M_{i}$. $p_{s}$ is the saturation vapor pressure (which is equal to $p_{\textrm{v}}$ in our case). For H2O, we use the Tetens formula. $L_{\mathrm{v}}$ is the latent heat of vaporization. We assume here that all the condensed material falls into deeper layers. As a consequence, the term $q_{\textrm{c}}c_{\textrm{p,c}}$ disappears and the process can no longer be considered as strictly adiabatic (hence the term pseudo-adiabat sometimes used). We note that below the H2O cloud, $\nabla_{\mathrm{moist}}=\nabla_{\mathrm{dry}}$.

2.2.4 3-layer profiles

Condensation of H2O can result in a significant gradient in molecular weight when high oxygen enrichments are considered. The increase of molecular weight with pressure introduces a stabilizing effect. Leconte et al. (2017) derived a new stability criterion that applies in this situation. When

[TABLE]

where $\varpi\equiv(\mu_{v}-\mu_{a})/\mu_{v}$ and $q_{\mathrm{v}}\equiv q_{\mathrm{H}_{2}\mathrm{O}}$, convection is inhibited and a radiative layer is formed. This layer is stable against double-diffusive processes in 1D, as demonstrated by Leconte et al. (2017). The thermal gradient from (7) must then be replaced by a radiative gradient

[TABLE]

where $\kappa$ is the Rosseland mean opacity, $g$ the gravity and $T_{\textrm{int}}$ the temperature associated with the internal heat flux (see Table 1). For the Rosseland mean opacities, we use the parametrization from Valencia et al. (2013). This parametrization requires a metallicity, which we compute by accounting for the elemental abundances of oxygen and carbon.

Thus, if $q_{\mathrm{H}_{2}\mathrm{O}}$ exceeds a given threshold (see equation 17 in Leconte et al. 2017), a positive feedback appears in the condensation region. As the mean molecular weight gradient, that manifests itself as an increase of $q_{\mathrm{H}_{2}\mathrm{O}}$ in equation (8), stabilizes the layer, the thermal gradient can increase and exceed the moist gradient owing to (9). As a consequence, more H2O is vaporized and the H2O mole fraction at the level $i+1$ increases, resulting in an increase of the mean molecular weight gradient. For high enough oxygen elemental abundances, this creates a very localized radiative region ($\Delta z\sim$1 km, whatever the O/H beyond the threshold) that separates the water-poor and water-rich regions. This interface is associated with a jump in temperature due to inefficient heat transport, as demonstrated by Leconte et al. (2017). Profiles for various oxygen enrichments are displayed in Figs. 1 and 2 for Uranus and Neptune.

2.2.5 Extreme case and sampling range

The magnitude of the temperature jump in the 3-layer profiles is controlled by the magnitude of the Rosseland opacities. Those are not well-constrained for heterogeneously enriched mixtures (Valencia et al., 2013). Moreover, any possible opacity of the condensate itself has been neglected with the formulation we have used. Because the jump in temperature at the transition zone between a water-rich and a water-poor region may be governed by these opacities, we may overestimate or underestimate this jump.

For illustration, let us take the example of Neptune around the nominal case that will be presented later in Section 4. Even if we divide the Rosseland opacities by a factor of 10, the upper tropospheric mole fraction of CO is unnoticeably altered. On the other hand, if the opacities are increased by a factor of 5, $y_{\mathrm{CO}}^{\mathrm{top}}$ is already multiplied by a factor of 10. At Uranus, the effect is less prominent, because of the smaller temperature jump implied by the lower nominal O/H: an increase of the opacities by a factor of 10 increases $y_{\mathrm{CO}}^{\mathrm{top}}$ by 25%. It is thus important to see if we can set limits to this source of uncertainty.

Interestingly, the formulation of Leconte et al. (2017) implies that the temperature jump cannot exceed a limit. This means that even if the opacities are significantly underestimated, the temperature jump cannot grow infinitely. This limit is given by the following $f_{\mu}$ factor:

[TABLE]

and $f_{\mu}$ is a multiplicative factor that is applied to the temperature of a moist profile at the H2O condensation level. In the previous expression, $\mu_{\mathrm{a}}$ is the mean molecular weight of the dry air and

[TABLE]

is the mass fraction of H2O. Indeed, if we assume that a zone, in which dry convection would occur despite the presence of H2O, can exist, then the virtual potential temperature of equation 8 in Leconte et al. (2017) should be approximately constant. If we now consider that the mole fraction of H2O is negligible in the upper troposphere, the deep temperature cannot exceed the temperature given by the virtual potential temperature of the upper troposphere. This gives the limiting factor of equation 10. In practice, the temperature at the bottom of the atmosphere should be always lower than this limit, because of diffusive phenomena and moist convection that transport energy in the form of latent heat. We will explore the range of allowed thermal profiles by sampling the allowed profiles between the coldest one, which is given by the moist adiabat, and the warmest one, which is given by a moist adiabat in which the temperature at the H2O condensation layer is multiplied by $f_{\mu}$. We emphasize that the 3-layer profiles always fall inside this range of profiles. In practice in this study, we compute 20 profiles from the coldest to the warmest, for each set of O/H, $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$ and $K_{zz}$.

2.3 Step 2 - Calculating the thermochemical equilibrium

Once the internal elemental composition is estimated and the thermal profile derived, we can compute the atmospheric composition as a function of altitude assuming thermochemical equilibrium, i.e. in the absence of vertical mixing. This thermochemical equilibrium state is calculated by minimizing the Gibbs energy and these calculations are based on the algorithm of Gordon and McBride (1994). This equilibrium code has been developed by Agúndez et al. (2014) and adapted for our purposes. The results of this code are then used as initial conditions of the thermochemical and diffusion model. This step enables us to speed up the latter computations by feeding the model with a non-empty atmosphere.

2.4 Step 3 - Calculating the thermochemistry and diffusion steady state

We use the thermochemical and diffusion model of Venot et al. (2012) adapted to the giant planet tropospheres to compute the steady state composition of the troposphere. This model has the advantage to be based on a chemical scheme that has been validated intensively by the combustion industry (Bounaceur et al., 2007). A full description of the model and its chemical scheme can be found in Venot et al. (2012) and can be retrieved from the KIDA database (Wakelam et al., 2012).

We start from an initial condition assuming chemical equilibrium as detailed in the previous section, assume a value for $K_{zz}$, and integrate the continuity equation over $\sim 10^{8}$ s usually. This enables reaching a steady state in all cases. As CO has also an external source in Uranus and Neptune (Cavalié et al., 2014; Lellouch et al., 2005), we have made sure that any reasonable change in the magnitude of the CO external source had no significant effect on the tropospheric CO mole fraction that is usually implicitly assimilated to an internal source. In Neptune, the external source is probably a relatively recent comet (Lellouch et al., 2005; Luszcz-Cook and de Pater, 2013) that should therefore not have had enough time to contaminate the upper tropospheric CO profile. In Uranus, the nature of the external source is still uncertain (Cavalié et al., 2014). However, the fact that the vertical mixing time through the stratosphere is more than two orders of magnitude higher than the integration time we need in our computations makes the tropospheric CO more sensitive to the internal source than to the external source. As a consequence, we can constrain the deep oxygen abundance with our model by fitting the observed upper tropospheric CO.

2.5 Step 4 - Tropospheric vertical mixing estimation and validation

In 1D thermochemical models, convective mixing is usually approximated by a vertical eddy diffusion coefficient $K_{zz}$. A thermochemical reaction can be quenched when the vertical diffusion timescale becomes shorter than the chemical timescale. Therefore, constraining the level at which the quenching of CO happens requires an estimate of vertical transport timescales in the tropospheres of giant planets.

As shown in previous studies, it appears that the final CO tropospheric mole fraction is sensitive to the Kzz at the CO quench-level only. This quench-level remains typically around the same temperature level (i.e. $\sim$900 K). So, we need to use a value relevant for this temperature level, and will assume uniform $K_{zz}$ coefficients. The validity of the latter assumption is further discussed in Section 6.1.

Following the free-convection and mixing-length theories (MLT) of Stone (1976) and Gierasch and Conrath (1985), Visscher et al. (2010) give a generic form for $K_{zz}$:

[TABLE]

This formulation, in which $k_{B}$ is the Boltzmann constant, requires the knowledge of the planet’s internal heat flux $F$, the atmospheric mean mass density $\rho$, the atmospheric mean molecular mass $m$, the atmospheric mass specific heat at constant pressure $c_{\mathrm{p}}$, and the atmospheric scale height $H$. The value of $F$ has been measured by Voyager 2 (see Table 1). This formulation is supposed to provide estimates valid within an order of magnitude.

The aforementioned formulation based on MLT thus depends on the composition (via $\rho$, $m$, $c_{\mathrm{p}}$, and $H$) and on the thermal structure (via $H$ and $c_{\mathrm{p}}$). The thermal structure in our models also depends on composition (via mean molecular weight effects and, to a lesser extent, $c_{\mathrm{p}}$). Therefore, we cannot compute $K_{zz}$ a priori and we will investigate a broad range of value (typically 5-6 orders of magnitude) in our thermochemical computations. Only after obtaining the thermal and composition results are we able to establish the validity range for $K_{zz}$ using the MLT formulation.

More recently, Wang et al. (2015) have used laboratory studies of turbulent rotating convection to derive a new formulation of $K_{zz}$ that provides estimates with a relative uncertainty of $\sim$25% only. We note that their formulation for low latitudes is similar to equation 12, only corrected by a scaling factor. So, again, their formulation exhibits a dependence on temperature and composition.

3 Observational constraints

3.1 Tropospheric abundances observed in giant planets

The thermochemical and diffusion model is initiated by specifying the deep abundances of the following elements : H, He, C, N, and O. We thus have to adjust the elemental abundances to ensure fitting the abundances of the main species in the upper troposphere, where their abundances have been measured. We have reviewed the abundances of species relevant to our model in all the Solar System giant planets and listed the values in Table 1.

The tropospheric CO mole fractions come from recent interferometric and space-based observations for Uranus and Neptune. On Uranus, Teanby and Irwin (2013) have obtained an upper limit of 2.1$\times 10^{-9}$, which is an improvement over Cavalié et al. (2008b) by almost an order of magnitude. The observations of Luszcz-Cook and de Pater (2013) that we use for Neptune have quite significantly revised the previously accepted value of 1.0$\times 10^{-6}$ (Rosenqvist et al., 1992; Marten et al., 1993, 2005; Lellouch et al., 2005, 2010; Fletcher et al., 2010), and their result is confirmed (and the error bar is narrowed) by the recent IRAM-30m and Herschel/SPIRE observations (R. Moreno, priv. com.) that we use here as our nominal value: 0.2$\times 10^{-6}$.

The helium mole fraction we use for Uranus and Neptune come from the Voyager 2 measurement of Conrath et al. (1987) and from Infrared Space Observatory observations (Burgdorf et al., 2003). Both indicate that there is 15% of helium in both atmospheres.

Contrary to Jupiter and Saturn, Uranus and Neptune have cold enough tropopauses for CH4 to condense and sharply decrease from their tropospheres to their stratospheres (Lellouch et al., 2015). Karkoschka and Tomasko (2009) and Sromovsky et al. (2011, 2014) for Uranus, and Karkoschka and Tomasko (2011) for Neptune, have recently shown that the tropospheric CH4 abundance is not uniform with latitude. The widely accepted measurements of Baines et al. (1995) ($\sim$2% in both planets) seems to be representative of high latitudes, while the CH4 equatorial tropospheric mole fraction is around 4%. We will take this latter value as our nominal abundance for CH4.

3.2 Upper tropospheric temperatures and internal heat fluxes

As a starting point for the temperature extrapolation in the troposphere, we take the temperature at 2 bar in Uranus and Neptune, as measured by Orton et al. (2014), and Lindal (1992), respectively.

The computation of $K_{zz}$ requires the knowledge of the internal heat flux. We take the Voyager 2 measurements from Pearl et al. (1990) and Pearl and Conrath (1991) for Uranus and Neptune.

4 Results

Hereafter, we present the results of our exploration of a 4D parameter space (O/H, C/H, temperature, $K_{zz}$) in thermochemical and diffusion computations in an attempt to constrain the deep oxygen abundance of Uranus and Neptune by fitting the observed upper tropospheric CO mole fraction. From this set of results, we also show how the O/H ratios derived from a thermochemical model that considers “3-layer” thermal profiles as nominal compare to previously published results in which dry or wet adiabats had been assumed. We want to keep the reader aware of the fact that the “3-layer” thermal profiles have been shown to be stable in 1D (Leconte et al., 2017), but may not be in a 3D treatment. The model is affected by sources of uncertainties other than those investigated in this 4D study and they are described in Section 6. So, in the following sections, we present the results of our 4D grid computations, and we provide the reader with the deep O/H obtained for Uranus and Neptune with our model, when assuming our nominal input parameters.

In the text that follows, the O/H and C/H ratios will be presented as a function of the O and C solar abundances, i.e. 8.73 dex for O and 8.39 dex, as reported by Lodders (2010).

4.1 $K_{zz}$ a posteriori validation

The estimation of the deep composition at step 0 and the calculation of the thermal profile at step 1 of our modeling enables us to compute estimates of $K_{zz}$ using the MLT or the Wang et al. (2015) formulation. For simplicity, we take the MLT formulation. The results presented below include an a posteriori consistency check between the estimated $K_{zz}$ and the value assumed in the thermochemical modeling (step 3). This is the step 4 of our analysis. For Uranus and Neptune, we find $K_{zz}\simeq 10^{8}$  cm${}^{2}\cdot$s*-1*.

4.2 Exploring the 4D parameter space

4.2.1 Methodology

The $y_{\mathrm{CO}}^{\mathrm{top}}$ we obtain from our computations depends on the following quantities: the deep O abundance, $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$ (as it controls $X_{C}$), the deep tropospheric temperature, and $K_{zz}$. To properly explore this 4D parameter space, we have run our models detailed in steps 0-3 of Section 2 in a loop format, by varying these quantities.

For the deep oxygen abundance we go from 10 to 300$\Sun$ for Uranus and from 10 to 600$\Sun$ for Neptune, and for the $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$ mole fraction we go from 1.5% to 4.5%. The temperature profiles are computed once these quantities are fixed.

The link between the upper tropospheric CO mole fraction and the deep O abundance depends on the level at which the quenching of CO occurs, i.e. the level at which the chemical lifetime of CO, which is mainly governed by temperature, equals its diffusion timescale given by $K_{zz}$. This level is usually located below the H2O cloud base, i.e. in a region where a dry adiabat prevails in any case. Because the behavior of the thermal profile remains poorly constrained around the H2O condensation level, and because the thermal profile should follow a simple dry adiabat below this level, we can explore the various possible deep thermal profiles by starting a dry adiabatic extrapolation from below the H2O cloud base. The starting temperature range from the temperature given by the moist adiabat, i.e. the coldest case, and the extreme case derived from Leconte et al. (2017). The 3-layer profiles fall, by definition, in this range. When computing the extreme thermal profiles, we apply the $f_{\mu}$ factor from equation 10 to the temperature at the level where the temperature jump occurs in the 3-layer profiles. We thus have the possibility to limit the range of possible temperature jumps, from a situation in which the opacities would be so low that no jump would be formed and the temperature profile would finally follow the moist adiabat to a situation in which the opacities would be high enough so that the temperature jump would reach its allowed limit. We sample this range with 20 profiles.

Finally, we vary the values of $K_{zz}$ over 5-6 orders of magnitude around the value that provides a fit to $y_{\mathrm{CO}}^{\mathrm{top}}$ in the “3-layer” case.

4.2.2 Effect of $K_{zz}$

This case is obvious as the deep O abundance required to fit the upper tropospheric CO decreases with increasing $K_{zz}$. The relationships for Uranus and Neptune are presented in Figs. 3 and 4. They show the range of values for ($K_{zz}$,O/H) that are compliant with observations. We remind the reader that the observations at Uranus have only enabled setting an upper limit on $y_{\mathrm{CO}}^{\mathrm{top}}$, and all values of O/H below the black curve shown in the figure are thus not ruled out.

4.2.3 Effect of $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$

Cavalié et al. (2014) investigated the effect of the tropospheric CH4 mole fraction on the deep O abundance in Uranus. They found that an increase of the CH4 abundance implied a decrease of the required O abundance to fit the tropospheric CO, essentially because increasing the C abundance implies an increase C atoms available to form CO. However, they used strictly similar thermal profiles for their two cases and neglected the influence of CH4 on the thermal profile.

In this work, we take into account the influence of the mean molecular weight $\mu$ and of the $c_{p}$ of CH4 and H2O when computing the thermal profile. Increasing the C abundance now also implies increasing $\mu_{a}$ in equation 8 and thus having a temperature jump located deeper in the troposphere. Moreover, the thermal gradient is lower because of a higher $\mu c_{p}$. Consequently, the deep tropospheric temperatures are colder when $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$ increases, if other quantities are kept constant. However, the increase of available C atoms to form CO remains the dominant effect when $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$ (and thus $X_{C}$) increases. So, if $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$ (and thus the deep C abundance) increases, a lower deep O abundance is needed to fit the tropospheric CO.

4.2.4 Effect of the temperature jump magnitude in the H2O condensation region

We have investigated 20 profiles for each set of O/H, C/H and $K_{zz}$ that range from the coldest possible profile, i.e. a wet adiabat, to the warmest one, i.e. a wet adiabat in which the multiplication factor $f_{\mu}$ is applied to the temperature in the layer where H2O condenses. The latter is the extreme profile of the “3-layer” case and the multiplication factor is given by equation 10. The nominal 3-layer profiles always fall within this range.

The results are quite difficult to present, in the sense that each individual O/H results in a different $f_{\mu}$ factor. A first indication can be obtained by comparing the results obtained with wet adiabats and nominal 3-layer profiles. Fig. 5 essentially shows that the 3-layer profiles require significantly lower O/H ratios to produce a given upper tropospheric CO mole fraction. To illustrate the full range of 3-layer profiles that are allowed and their impact on the results, we choose to normalize the results with respect to the maximal magnitude of multiplication factor for each O/H, noted $f_{\mu}^{\mathrm{max}}$. We thus compute the curves for f_{\mu}-1$$/$$f_{\mu}^{\mathrm{max}}-1 as a function of O/H. The values range from the wet adiabat, for which f_{\mu}$$=1 and f_{\mu}-1$$/$$f_{\mu}^{\mathrm{max}}-1=0, to the maximum amplitude temperature jump, for which f_{\mu}$$=$$f_{\mu}^{\mathrm{max}} and f_{\mu}-1$$/$$f_{\mu}^{\mathrm{max}}-1=1. Figs. 6 and 7 display the results for both planets. In this case, $K_{zz}$ is fixed to 108  cm${}^{2}\cdot$s*-1*. The downside is that the absolute magnitude of $f_{\mu}$ is not directly accessible to the reader. It essentially shows that a much higher O/H would be needed to fit the data if a pure wet adiabat was used in Uranus, and that even the warmest profile cannot produce enough CO in the upper troposphere if the deep O abundance is lower than $\sim$80 times solar. However, the tropospheric CO value is an upper limit in Uranus and this constraint is therefore useless. In Neptune, the fraction of the maximal multiplication factor sufficient to produce enough CO in the troposphere decreases quickly with increasing O/H, as this multiplication factor is already huge for O/H of several hundreds and a small fraction of it produces a big temperature jump in the zone where H2O condenses. The results also show that a minimum of O/H$\sim$190 times solar seems to be a lower limit in Neptune with our model.

4.3 Uranus nominal case

“3-layer” thermal profiles of Uranus present a less spectacular transition zone than in Neptune, but still modify significantly the upper limit on the O/H ratio when compared to cases in which we would use dry or wet adiabats. We find that reproducing the observed upper limit of 2.1$\times 10^{-9}$ of CO and the reported He and CH4 abundances (see Table 1) requires C/H and O/H ratios of 75, and 160 times the solar value. The resulting abundance profiles are displayed in Fig. 8. As in the Neptune case, using dry or wet adiabats would lead to much lower amounts of CO when starting from the same elementary abundances.

4.4 Neptune nominal case

As stated above, Neptune is expected to be the most “symptomatic” case demonstrating the impact of the mean molecular weight gradient in the transition zone from the water-rich to the water-poor region in its troposphere. It leads to a “3-layer” thermal profile with a significant departure from the purely dry or wet adiabats that had been used previously in the literature (e.g., Lodders and Fegley 1994; Luszcz-Cook and de Pater 2013). To fit the nominal 0.2$\times 10^{-6}$ of tropospheric CO (see Table 1), we had to set the O/H ratio to $\sim$540 times the solar value. The corresponding “3-layer” thermal profile is displayed in Fig. 10. The nominal C/H ratio is 40 times the solar value. This value is rather low, if compared to a value that would simply be deduced from the upper atmospheric abundances of CH4. The difference is caused by the high internal abundance of O and thus the high deep H2O abundance. Indeed, we find in this calculation that H2O is more abundant than H2 in Neptune’s interior.

5 Discussion

This section will be divided in two parts. In the first one, we will discuss the results of the nominal model, as presented in Section 4. In the second part, we will discuss the microwave spectral counterpart of our model results to further constrain their validity.

5.1 Nominal results

We find nominally that Uranus has a C/H ratio of 75 times the solar value and a O/H ratio below 160 times the solar value. For Neptune, we get a C/H 40 times the solar value and a O/H 540 times the solar value.

5.1.1 A deep compositional difference between Uranus and Neptune?

One of the striking compositional differences in Uranus and Neptune atmospheres is their CO abundances, while their He and CH4 abundances and D/H ratios are relatively similar (Baines et al., 1995; Karkoschka and Tomasko, 2009, 2011; Sromovsky et al., 2011, 2014; Feuchtgruber et al., 2013). Indeed, they differ by at least two orders of magnitude in the troposphere (Teanby and Irwin, 2013; Luszcz-Cook and de Pater, 2013). However, previous evaluations of the O/H ratios in these bodies using the quench-level approximation and dry/wet adiabats never differed by more than a factor of two (Lodders and Fegley, 1994; Luszcz-Cook and de Pater, 2013; Cavalié et al., 2014). Accounting for the stability criterion of Leconte et al. (2017) in the water condensation layers has a significant effect on the planets thermal profile and hence on the derived O/H values. Our nominal values show that O/H$<$160 and O/H$=$540 times the solar value for Uranus and Neptune, respectively. So, there is a factor of $\geq$3.3 difference between these two values. The two-order of magnitude difference between the upper tropospheric CO mole fraction on Uranus and Neptune is still caused by a much smaller difference between their deep H2O abundances. This can be understood from the equilibrium equation 1. An increase on H2O implies a linear decrease of H2. As the power on H2 is 3 (and the CH4 abundance is relatively constant), CO has to increase non-linearly to equilibrate the reaction. This effect is shown for both planets on Fig. LABEL:CO_H2O, where we plot the relationship between the CO and H2O abundances at the quench level for our nominal cases for $K_{zz}$ and $y_{\mathrm{CH}_{4}}^{\mathrm{top}}$. Another interesting feature of our results is that a similar CH4 abundance in Uranus and Neptune does not reflect in a similar C/H ratio. This can be understood from equations 2 and 3. Indeed, the very high O/H ratio of Neptune implies a lower C/H ratio than in Uranus to reproduce both the CO and CH4 abundances in the upper troposphere. In Uranus, the deep H2O abundance is comparable to the He abundance, while the deep H2O abundance even exceeds the H2 abundance in Neptune.

We confirm that a large difference in terms of CO abundance (two orders of magnitude) does not imply such a large difference in terms of O/H ratios in the two planets in our nominal model. This rather small difference seen in their O/H ratios, compared to their CO abundances, is then mostly caused by the 10 K difference at 2 bar, where we initiate our thermal profiles.

5.1.2 Consistency with interior and formation models

While the tropospheric mole fractions of CO are very different in Uranus and Neptune, their D/H ratios are quite similar, according to Feuchtgruber et al. (2013). One of the possible explanations presented by these authors, based on the interior models of Helled et al. (2011) and Nettelmann et al. (2013), is that the cores of the planets may be more rocky than icy, with rock fraction of 68-86%, for their D/H to be still representative of the one seen in Oort cloud comets. The O/H in core ices implied in such cases ranges from 79 to 172 times the solar value for Uranus and from 68 to 148 times the solar value for Neptune (Feuchtgruber et al., 2013). These numbers, relevant only if the planets went through complete mixing at least once in their histories, are compatible with the O/H in the outer envelope of Uranus but seem nevertheless too low for Neptune’s core ices to result in a sufficiently high outer envelope O/H ratio compared to our nominal results.

A scenario has been proposed recently to reconcile the D/H measurements at Uranus and Neptune with the Oort cloud comet value by Ali-Dib et al. (2014). They propose that Uranus and Neptune formed on the CO snowline and that their ices were mostly composed of CO ices rather than H2O ices. Therefore, the present day atmospheric D/H would be representative of the dilution of the small fraction of D-enriched Oort-comet-like H2O ices in the H2O coming from the thermochemical conversion of the more abundant primordial CO. This scenario then predicts that oxygen should be enriched in a similar way as carbon in both planets. Our simulation results show that this is probably neither the case for Uranus nor for Neptune and seem to contradict the model of Ali-Dib et al. (2014). We find C/O ratios of $>$0.21 and $\sim$0.03, respectively, for the two planets. In Neptune, the C/O ratio could be brought back to one if CO tropospheric mole fraction was several orders of magnitude below its measured value, i.e., at the lower end of the measurement range of Luszcz-Cook and de Pater 2013). However, recent Herschel-SPIRE and IRAM-30m observations show that the CO tropospheric mole fraction is rather in the upper part of the measurement range of Luszcz-Cook and de Pater (2013), with a measured value of $y_{\mathrm{CO}}$ of 0.20$\pm$0.05 ppm (R. Moreno, priv. com.). In any case, the model of Ali-Dib et al. (2014) seems not to be compatible with the Nice model (Gomes et al., 2005; Morbidelli et al., 2005; Tsiganis et al., 2005) and it would anyway need to be improved by accounting for the protoplanetary disk temporal evolution as well as ice and core formation kinetics to better constrain the validity of their assumption that core formation on the CO snowline is possible.

Our results can further be used for direct comparison with formation model predictions. For instance, Hersant et al. (2004) predicted the O/H enrichment in Uranus and Neptune by exploring the range of possible C/H values when assuming that heavy elements were trapped by clathration. For instance, if C was carried by only CO during the formation of planetesimals, then 5.75 H2O molecules were required to trap CO in clathrates. Then, the resulting C/O ratio is 1/(1$+$5.75)$=$0.15. However, their prediction is less extreme and is as follows: for C/H ratios of 40, 60, and 80 times the solar value, the O/H ratio should be 87, 130, and 176 times the solar value (based on solar abundances of Anders and Grevesse 1989). In this case, half of the carbon comes from clathrated CO and half from pure CO2 ice, and the resulting C/O ratio is 1/(2$\times$1/2 + 5.75$\times$1/2 + 1/2)$=$1/4.375$=$0.23. For reference, the solar C/O ratio is 0.46, according to Lodders (2010). The lower limit we obtain on the C/O ratio of Uranus is at the limit of compatibility with the case in which carbon equally comes from clathrated CO and CO2 ice. On the other hand, the Neptune C/O ratio we derive would become compatible with the case in which all carbon comes from clathrated CO, as long as the efficiency of clathration was about 20%.

The rather large difference in terms of C/O ratio in Uranus and Neptune seems to imply different condensation environments, which may be difficult to explain. The question then is, whether the Uranus upper tropospheric CO value is representative for the CO quench level value? Actually, convection could be much less efficient than assumed in our computations. We have indeed based our computations of the tropospheric $K_{zz}$ of Uranus on the nominal value of $F$ reported in Pearl et al. (1990). However, the error bars on $F$ are such that even a zero value is compatible with the Voyager observations. A lower $F$ implies a lower $K_{zz}$, as $K_{zz}\propto F^{1/3}$ (see equation 12). We can then reversely use our model to constrain $K_{zz}$ to bring the O/H of Uranus to an O/H compatible with the clathration scenario of Hersant et al. (2004). If the O/H of Uranus is set to $\sim$450 times the solar value, so that C/O$\sim$0.15, then C/H is $\sim$50 times solar and K_{zz}$$\sim105 cm${}^{2}\cdot$s*-1*. Bringing the Uranus O/H further up to the Neptune nominal value implies $K_{zz}$ also $\lesssim$105 cm${}^{2}\cdot$s*-1*. It is thus possible that Uranus and Neptune are more similar regarding their oxygen reservoir than their tropospheric CO seem to tell us, provided that convection is strongly inhibited in the interior of Uranus. This seems to be supported by recent internal structure and evolution modeling. Nettelmann et al. (2016) attempt to fit both the low luminosity of Uranus and its gravitational data with an internal structure and evolution model. Achieving this requires a thermal boundary around 0.1 Mbar that breaks the adiabat. They also find that fitting $J_{4}$ for Uranus requires an O/H of less than 30, while C/H is about 80. However, they do not consider C/O$>$1 to be a good assumption for the core, implying that deep stratification is the cause for a lower O abundance at observable levels.

We can also compute the outer envelope heavy element mass fraction from the deep abundances of the species accounted for in the model and check whether they are in agreement with the internal structures inferred by Nettelmann et al. (2013). For Uranus and Neptune, we find Z_{1}$$<55% and $=$80%, respectively. The Neptune value is above the limit of $\sim$65% allowed by their model (see figure 3 in Nettelmann et al. 2013). Only if we consider the most extreme case of temperature jump in our simulations, in which the O/H of Neptune cannot be lower than 190 times solar, can we get within the limit of Nettelmann et al. (2013), with $Z_{1}$ equal to 58%. Our Uranus O/H value, being an upper limit higher than the range of values allowed by their model, is not constraining. Figure 3 of Nettelmann et al. (2013) indicates an upper limit for $Z_{1}$ in Uranus of $\sim$10%, which in turn would imply an O/H ratio much lower than our upper limit. Carbon which is mainly carried by CH4 already represents a mass fraction $\geq$10% in the outer envelope, and a lower O/H would result in a higher C/H and C mass fraction. Thus, it seems difficult to reconcile the model of Nettelmann et al. (2013) of Uranus with upper tropospheric measurements and their implication on deep tropospheric elemental composition. The case of Uranus thus remains puzzling, and more generally, the formation of Uranus and Neptune is difficult to explain (Helled and Lunine, 2014). It thus underlines the need for new data. In this sense, dedicated orbital missions with atmospheric descent probes are highly desirable (Arridge et al., 2012, 2014; Masters et al., 2014; Turrini et al., 2014).

Finally and as will be discussed in Section 6.2, our model may underestimate the production of CO. In this case, the model would require lower O abundances and bring the results closer to an agreement with the predictions from Nettelmann et al. (2013) for both planets.

5.2 Radio spectrum of Uranus and Neptune

To further constrain the validity of the results of this paper, we present radiative transfer simulations in the microwave range with the temperature and abundance profiles of our thermochemical computations. It is indeed worth checking whether the temperature jump implied for high deep O abundances and the sharp composition transition in the H2O condensation zone translate into spectra that are in agreement with observations in the mm-cm range (e.g. de Pater et al. 1991).

We thus take our results in terms of temperature profiles and corresponding abundance profiles. We use the radiative transfer model already described in Cavalié et al. (2008a, 2013) to compute synthetic spectra of Uranus and Neptune in the mm-cm range. We extend our tropospheric thermal profiles of Uranus and Neptune from 2 bar to 10 mbar with the data from Orton et al. (2014) and from the Herschel Science Centre “ESA5” model222http://archives.esac.esa.int/hsa/legacy/ADP/PlanetaryModels/ profiles for Uranus and Neptune, respectively. We add profiles of H2S and NH3 based on DeBoer and Steffes (1994) to account for their respective opacities, in addition to the microwave opacity of the H2O line at 22 GHz. The NH3 mole fraction between the NH4SH and NH3 clouds is set to 6$\times 10^{-6}$ and 1$\times 10^{-5}$ in Uranus and Neptune (Moreno 1998 and “ESA5” model, resp.). Below the NH4SH cloud, the NH3 mole fraction is set to 1$\times 10^{-4}$ and the deep H2S mole fraction is set such that all H2S is consumed to form the NH4SH cloud. We use Bellotti et al. (2016) for the H2O microwave opacity, Moreno (1998) for the NH3 microwave opacity, Liebe and Dillon (1969) for H2- and He-broadening and temperature dependence parameters for the 22 GHz line of H2O, Fletcher et al. (2007) for the H2-broadening and temperature dependence parameters of NH3 lines, DeBoer and Steffes (1994) for the H2S microwave opacity and broadening parameters, and de Pater et al. (1991) and Fray and Schmitt (2009) for H2S condensation laws. Other spectral line parameters come from the JPL Molecular Spectroscopy Database (Pickett et al., 1998). The collision-induced absorption caused by H2-H2, H2-He, and H2-CH4 pairs is computed following Borysow et al. (1985, 1988) and Borysow and Frommhold (1986).

The resulting spectra for Uranus and Neptune in the centimeter wavelength range, when considering the nominal cases of Section 4, are displayed in Fig. 11. Even though we make no particular attempt to improve the fit to the data (especially around 1 cm) by modifying the NH3 and/or H2S abundances, the obtained spectra are in good agreement with observations presented in Gulkis et al. (1978), Cunningham et al. (1981), Muhleman and Berge (1991), Griffin and Orton (1993), Greve et al. (1994), Klein and Hofstadter (2006), and Weiland et al. (2011), for Uranus, and Orton et al. (1986), Romani et al. (1989), de Pater and Richmond (1989), de Pater et al. (1991), Muhleman and Berge (1991), Griffin and Orton (1993) , Hofstadter (1993), Greve et al. (1994), and Weiland et al. (2011), for Neptune. Even the extreme cases, in which the O abundance is lower but compensated by the biggest temperature jump allowed by the model (see Section 4.2.4), the cm-wave spectrum of both planets remains in agreement with observations. This comes from the fact that the region where the jump in temperature and H2O abundance occur, i.e. the radiative layer, is only marginally probed by the H2O absorption at 20 cm wavelength, as shown in DeBoer and Steffes (1996).

6 Limitations and other questions

In this section, we will detail the limitations of our model and stress the data that are desirable for increasing the predictability of such models.

6.1 The radiative layer as an insulation layer?

One of the novel aspects of our work is the inclusion of the stability criterion of Leconte et al. (2017) when extrapolating measured tropospheric temperatures to deeper levels where the atmosphere is in thermochemical equilibrium. Accounting for the mean molecular weight gradient when extending thermal profiles from observed levels in the upper troposphere to deeper levels implies the formation of a thin radiative layer ($\Delta z\sim$1 km) in Uranus and Neptune at the levels where H2O condenses. If such a stable layer exists, it should act as an insulation layer for the both the transport of heat and the chemical mixing, making our assumption of a vertically uniform $K_{zz}$ profile questionable. There is no estimate for chemical mixing under such conditions, but within this layer $K_{zz}$ could, in principle, be as low as that of molecular diffusivity $D$, resulting in a gradient for CO, a species that has a non-zero flux at the boundaries of the model. According to our model, molecular diffusivity is in the range of 10*-2*-10*-3* cm${}^{2}\cdot$s*-1*, where the radiative layer in both planets resides. We estimate the variation in the CO mole fraction between the two limits of the radiative layer using \phi_{\mathrm{CO}}$$=$$-nD\frac{d{y_{\mathrm{CO}}}}{d{t}}, where $\phi_{\mathrm{CO}}$ is the flux of CO at the upper boundary of the model (i.e., the top of the troposphere) and $n$ is the total number density. As CO is relatively stable in the stratospheres of Uranus and Neptune, a good approximation for $\phi_{\mathrm{CO}}$ is given by the input fluxes caused by external sources, and we take 105 and 108 cm${}^{-2}\cdot$s*-1* respectively (Cavalié et al., 2014; Lellouch et al., 2005). We confirm such values running test cases in the stratosphere with the photochemical model of Dobrijevic et al. (2010). In this case only, the variation in the CO mole fraction between the two limits of the radiative layer would be of the order of magnitude of the observed $y_{\mathrm{CO}}^{\mathrm{top}}$. This means that the tropospheric CO could be partly caused by the external source because of the transport boundary caused by the radiative zone. Then, it would be impossible to constrain the deep oxygen abundance from the sole observation of $y_{\mathrm{CO}}^{\mathrm{top}}$ and thermochemical modeling, as presented in this paper. A full model accounting for thermochemistry in the troposphere, photochemistry in the stratosphere and external sources would be required to do so. This potentially is a significant limitation that should be kept in mind, although it only applies in the case where $K_{zz}$ is as low as $D$ in the radiative layer and the flux of CO from the stratosphere is as strong as the external source of CO. As long as $K_{zz}\geq 10^{-1}$, the gradient disappears and $y_{\mathrm{CO}}^{\mathrm{top}}$ can be taken as a probe to the deep oxygen abundance. In 3D, overshooting across the 1 km radiative layer is plausible and would help in this sense.

6.2 Chemical scheme

Moses (2014) compared her chemical scheme with the one we have adopted in this work (from Venot et al. 2012) for studying the composition of HD189733b, a hot Jupiter. In her chemical scheme, she used the reaction rates updated by Moses et al. (2011), based on ab initio calculations. She found that the main difference between the two schemes resided in the reaction rate used for CH3OH $+$ H $=$ CH3 $+$ H2O. Moses et al. (2011) then proceeded with a much slower reaction than Venot et al. (2012). The difference between these two reaction rates, when used in the Venot et al. (2012) model has been shown for Jupiter and Saturn by Wang et al. (2016). The scheme we use nominally in this study produces 10-30 times less CO than the same scheme with the Moses et al. (2011) rate for the CH3OH $+$ H $=$ CH3 $+$ H2O reaction (see also their figures 4 and 5).

If we proceed as Wang et al. (2016) and alter in the Venot et al. (2012) scheme the reaction rate for CH3OH $+$ H $=$ CH3 $+$ H2O by using the value of Moses et al. (2011), and still use “3-layer” thermal profiles, we find that the (O/H,C/H) ratios in Uranus and Neptune become $\sim$ ($<$55,85) and $\sim$ (280,65) times solar, respectively. These values are significantly lower than with our nominal chemical scheme, as expected. However, we emphasize that these values are just presented for the sake of comparison. While the ab initio calculations of Moses et al. (2011) are correct, including this rate in the scheme of Venot et al. (2012) as prescribed by Moses (2014) in the chemical scheme does not enable reproducing of the experimental data of different combustion studies (Cathonnet et al., 1982; Hidaka et al., 1989; Held and Dryer, 1998) (R. Bounaceur, priv. comm. June 2016). Using such a reaction rate is thus not validated in the frame of our chemical scheme. In any case, it certainly underlines the urgent need for a better understanding of the kinetics of CH3OH at high temperatures and pressures.

6.3 Tropospheric equation of state

In this work, we have implicitly used the ideal gas law. Significant water abundances at high pressures obviously render such a choice questionable. Nettelmann et al. (2008) proposed a new equation of state for H2-He-H2O mixtures, but for a higher pressure/temperature regime than investigated in this paper. More recently, Karpowicz and Steffes (2013) have proposed a tropospheric equation of state for Jupiter’s troposphere. The effect on pressure and temperature seems rather limited (2% decrease of T, compared to the ideal gas law, at their lower boundary of 200 bars), when compared to the temperature increases caused by the mean molecular weight gradient effect shown in this paper. Nevertheless, it remains to be seen what the effect of applying such equation of state to Uranus and Neptune would be, as they have presumably much more H2O in their deep tropospheres than Jupiter. Unfortunately, the authors underline that their “equation of state will likely have decreased accuracy under jovian conditions for pressures exceeding much beyond 100 bars, and may be invalid at pressures exceeding 2500 bars”. New laboratory experiments and theoretical simulations are thus needed to enable applying reliable equations of state for H2-H2O-He in the pressure/temperature range we are interested in (see Baraffe et al. 2014 for a review). However, even if we had a robust non-ideal equation of state relevant for our conditions, our chemical model would need to be updated. Indeed, though we explicitly assume an ideal gas law in many calculations of our models, we also implicitly assume it in the kinetic rates.

6.4 Latent heat release by other condensates

In this work, we do not account for the release of latent heat caused by the formation of various clouds besides the H2O clouds. Other expected clouds, from the deep troposphere to the upper troposphere, are supposedly composed of NH4SH, NH3 or H2S (depending on which one is the most abundant in the deep troposphere and “survives” the formation of the NH4SH cloud), and CH4. Our choice to start our thermal profiles below the CH4 cloud deck is meant so as not to have to account for the CH4 cloud formation, related latent heat release, and even possible mean molecular weight gradient effect (Guillot, 1995).

7 Conclusion

In this paper, we have modeled the thermochemistry in the tropospheres of Uranus and Neptune, in an attempt to constrain their deep oxygen abundance from upper tropospheric observations of CO. The derivation of the deep O/H ratio in Giant Planet tropospheres from thermochemical computations requires precise measurements of upper tropospheric CO and CH4, knowledge of oxygen species kinetics, as well as a good knowledge of tropospheric temperatures and vertical mixing.

We have shown that the transition between a water-rich and a water-poor region in giant planet tropospheres results in a significant gradient of mean molecular weight that impacts, in turn, the shape of the thermal profile. Accounting for this gradient in the condensation zone of H2O by applying the stability criterion of Leconte et al. (2017) implies a significant departure of temperature from the dry/wet adiabats that have been used so far in previous papers. Our results show that large oxygen enrichments can produce a thin radiative layer ($\Delta z\sim$1 km) where H2O condenses. The thermal gradient, stabilized by the downward increasing mean molecular weight, gets very strong and becomes radiative in this zone. This results in higher temperatures in deep levels compared to dry/wet adiabats.

While the mean molecular weight gradient effects are beginning to be applied to internal and formation models of giant planets (Leconte and Chabrier, 2012; Vazan et al., 2015), our model shows for the first time the importance of accounting for these effects already at upper tropospheric levels for oxygen-rich giant planets (Leconte et al., 2017). Our results show that these new profiles lower the required oxygen abundances to reproduce the CO observations. However, our estimates are affected by the current limitations of the model, like the lack of Rosseland opacities for the considered mixtures. Another source of uncertainty is related to the high pressures at the quench level which would require non ideal equation of states and modified kinetics. While a better equation of state can be expected soon and is rather simple to implement, its impact on the kinetics might be extremely difficult to derive. Finally, it is noticeable that if (but only if) vertical mixing is as low as molecular diffusivity in the $\sim$1 km thin radiative layer, then it is much more complicated to constrain the deep oxygen abundance in Uranus and Neptune because the CO external source then needs to be precisely characterized as well.

While using such temperature models should lead to lower O/H ratios than when using dry/wet adiabats, our nominal value for Neptune’s O/H ratio is close to the value derived by Luszcz-Cook and de Pater (2013) (they have used wet adiabats). This is caused by differences in oxygen chemistry. This certainly underlines the need for improved knowledge of the oxygen species thermochemistry. For instance, with kinetics that reduce significantly the destruction of CH3OH, as proposed by Moses et al. (2011) and applied by Wang et al. (2016), then the O/H required to fit the upper tropospheric observations of CO significantly decreases. However, using such different kinetics first require validation within a complete scheme over the temperature/pressure range relevant to such studies. We believe this is the main strength of the chemical scheme we are using, but we do not underestimate the work that still needs to be done to improve our knowledge of CH3OH kinetics.

In our nominal model, in which we fix the temperature, the CH4 abundance, and the tropospheric mixing to our best estimates, the O/H in Uranus and Neptune are $<$160 and 540 times the solar value, respectively. While Uranus’ and Neptune’s ices formation and heavy element enrichment could result from clathration, this scenario would imply different mixtures of CO/CO2 at the time of planetesimal condensation. The difference seen in their C/O ratios ($>$0.23 for Uranus and $\sim$0.03 for Neptune) is actually rather puzzling and may therefore be an indication that convection is strongly inhibited in the troposphere of Uranus. Indeed, this appealing explanation could enable reconciling formation models for the two planets with observations of their tropospheric composition. Among other things, the inconsistency between Uranus’ and Neptune’s O/H remains puzzling for now, and new data from dedicated orbiters with atmospheric descent probes are highly desirable to study these mysterious worlds (Arridge et al., 2012, 2014; Masters et al., 2014; Turrini et al., 2014). In the meantime, James Webb Space Telescope observations will enable measurements that will better constrain the deep abundances of several heavy elements, like Ge, As, and P (Norwood et al., 2016b, a).

Strong uncertainties remain in Uranus and Neptune on the upper tropospheric CH4 mole fraction, on the tropospheric mixing and on the tropospheric temperature profile. Because these are central to such thermochemical modeling, we have explored this whole parameter space with a chemical model that uses a scheme validated in the temperature/pressure range of interest here. When more precise measurements of temperature, mixing and CH4 abundance are obtained, our results can be used to find the corresponding O/H required to fit the tropospheric CO in both planets.

In the case of the gas giants, the lower oxygen enrichments that are expected should result in a much less significant mean molecular weight gradient in the condensation region of H2O and thus a less spectacular effect on the thermal profile and on the derivation of their deep O/H ratio. In this sense, it will be interesting to compare our thermochemical model results with the measurements that will be provided by Juno (Matousek, 2007; Helled and Lunine, 2014), which successfully performed its Jupiter orbit insertion on July 4, 2016. At Saturn, a descent probe has recently been proposed for an ESA M5 mission (Mousis et al., 2014, 2016). While this probe may not be able to reach the well-mixed region of H2O and thus measure Saturn’s O/H, it would provide ground truth as to the abundances of other heavy elements and disequilibrium species and shed light on the formation processes of Saturn.

Acknowledgements

T. Cavalié thanks the CNRS/INSU Programme National de Planétologie (PNP) for funding support. O. Venot acknowledges support from the KU Leuven IDO project IDO/10/2013, from the FWO Postdoctoral Fellowship program, and from the Centre National d’Études Spatiales (CNES). T. Cavalié and O.Venot thank R. Bounaceur (LRPG, Nancy, France) for useful discussions on the kinetics of CH3OH. T. Cavalié thanks I. de Pater (Univ. of California, Berkeley, USA) for providing her H2S condensation laws. We thank Dr. Dong Wang and an anonymous reviewer for they constructive comments that help us improve our paper.