# Storms and convection on Uranus and Neptune: impact of methane abundance revealed by a 3D cloud-resolving model

Noé Clément, Jérémy Leconte, Aymeric Spiga, Sandrine Guerlet, Franck Selsis, Gwenaël Milcareck, Lucas Teinturier, Thibault Cavalié, Raphaël Moreno, Emmanuel Lellouch, et al.

► To cite this version:

Noé Clément, Jérémy Leconte, Aymeric Spiga, Sandrine Guerlet, Franck Selsis, et al. Storms and convection on Uranus and Neptune: impact of methane abundance revealed by a 3D cloud-resolving model. Astronomy & Astrophysics - A&A, 2024, 690 (October), pp.A227. 10.1051/0004-6361/202348936. insu-04690581v2

HAL Id: insu-04690581

https://insu.hal.science/insu-04690581v2

Submitted on 10 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Storms and convection on Uranus and Neptune: Impact of methane abundance revealed by a 3D cloud-resolving model

Noé Clément<fnref n="1" />,*, Jérémy Leconte<fnref n="1" />, Aymeric Spiga<fnref n="2" />, Sandrine Guerlet<fnref n="2" />,<fnref n="3" />, Franck Selsis<fnref n="1" />,
Gwenaël Milcareck<fnref n="2" />,<fnref n="4" />, Lucas Teinturier<fnref n="2" />,<fnref n="3" />, Thibault Cavalié<fnref n="1" />,<fnref n="3" />, Raphaël Moreno<fnref n="3" />,†,
Emmanuel Lellouch<fnref n="3" />,‡, and Óscar Carrión-González<fnref n="3" />,‡

Received 13 December 2023 / Accepted 26 July 2024

## ABSTRACT

**Context.** Uranus and Neptune have atmospheres dominated by molecular hydrogen and helium. In the upper troposphere (between 0.1 and 10 bar), methane is the third main molecule, and it condenses, yielding a vertical gradient in CH₄. As this condensable species is heavier than H₂ and He, the resulting change in mean molecular weight due to condensation serves as a factor countering convection, which is traditionally considered as governed by temperature only. This change in mean molecular weight makes both dry and moist convection more difficult to start. As observations also show latitudinal variations in methane abundance, one can expect different vertical gradients from one latitude to another.

**Aims.** In this paper, we investigate the impact of this vertical gradient of methane and the different shapes it can take, including on the atmospheric regimes and especially on the formation and inhibition of moist convective storms in the troposphere of ice giants.

**Methods.** We developed a 3D cloud-resolving model to simulate convective processes at the required scale. This model is non-hydrostatic and includes the effect of the mean molecular weight variations associated with condensation.

**Results.** Using our simulations, we conclude that typical velocities of dry convection in the deep atmosphere are rather low (on the order of 1 m/s) but sufficient to sustain upward methane transport and that moist convection at the methane condensation level is strongly inhibited. Previous studies derived an analytical criterion on the methane vapor amount above which moist convection should be inhibited in saturated environments. In ice giants, this criterion yields a critical methane abundance of 1.2% at 80 K (this corresponds approximately to the 1 bar level). We first validated this analytical criterion numerically. We then showed that this critical methane abundance governs the inhibition and formation of moist convective storms, and we conclude that the intensity and intermittency of these storms should depend on the methane abundance and saturation. In the regions where CH₄ exceeds this critical abundance in the deep atmosphere (at the equator and the middle latitudes on Uranus and at all latitudes on Neptune), a stable layer almost entirely saturated with methane develops at the condensation level. In this layer, moist convection is inhibited, ensuring stability. Only weak moist convective events can occur above this layer, where methane abundance becomes lower than the critical value. The inhibition of moist convection prevents strong drying and maintains high relative humidity, which favors the frequency of these events. In the regions where CH₄ remains below this critical abundance in the deep atmosphere (possibly at the poles on Uranus), there is no such layer. More powerful storms can form, but they are also a bit rarer.

**Conclusions.** In ice giants, dry convection is weak, and moist convection is strongly inhibited. However, when enough methane is transported upward, through dry convection and turbulent diffusion, sporadic moist convective storms can form. These storms should be more frequent on Neptune than on Uranus because of Neptune's internal heat flow and larger methane abundance. Our results can explain the observed sporadicity of clouds in ice giants and help guide future observations that can test the conclusions of this work.

**Key words.** planets and satellites: atmospheres – planets and satellites: gaseous planets

## 1. Introduction

Uranus and Neptune are the two most distant planets of our Solar System, and thus they receive little insolation (3.7 W m⁻<fnref n="2" /> for Uranus and 1.5 W m⁻<fnref n="2" /> for Neptune, compared to 1366 W m⁻<fnref n="2" /> for the Earth), in addition to having long radiative timescales (more than 100 terrestrial years at 1 bar). As a result, weak atmospheric activity might be expected, yet observations have highlighted intense meteorology, showing numerous discrete

cloud features (presumably composed of methane ice crystals) evolving on short timescales (Karkoschka & Tomasko 2011) and long-lasting, powerful storms (Molter et al. 2019). Hueso et al. (2020) list several observations as candidate moist-convection features in Uranus and Neptune, but they conclude that most cloud activity observed so far is probably not convective. The record of observations demonstrating moist convective activity in Uranus and Neptune is almost nonexistent because frequent observations at very high spatial resolution – which are lacking today – would be required. These clouds and storms need

<footnote n="1">¹ Laboratoire d'Astrophysique de Bordeaux, Univ. Bordeaux, CNRS, B18N, allée Geoffroy Saint-Hilaire, 33615 Pessac, France</footnote>

<footnote n="2">² Laboratoire de Météorologie Dynamique (IPSL), Sorbonne Université, Centre National de la Recherche Scientifique, École Polytechnique, École Normale Supérieure, Paris, France</footnote>

<footnote n="3">³ LESIA, Observatoire de Paris, Université PSL, CNRS, Sorbonne Université, Université de Paris, 5 place Jules Janssen, 92195 Meudon, France</footnote>

<footnote n="4">⁴ Laboratoire ATmosphère Milieux Observations Spatiales/Institut Pierre-Simon Laplace (LATMOS/IPSL), Sorbonne Universités, UPMC Univ Paris 06, Université Paris-Saclay, Université de Versailles Saint-Quentin-en-Yvelines, Centre National de la Recherche Scientifique, 78280 Guyancourt, France</footnote>

<footnote>* Corresponding author; noe.clement@u-bordeaux.fr</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Temperature profiles and variability of mean molecular weight in ice giants. Left panel: Uranus and Neptune temperature profiles from the Voyager 2 spacecraft, retrieved from radio occultations by Lindal et al. (1987, 1990); Lindal (1992). Saturation vapor pressure and triple point are indicated for 6% of CH₄ (the maximum observed value), which corresponds to a specific concentration of 0.30 kg/kg. Theoretically, methane only exists in a solid or gas form, and can only sublimate and condense. A lower methane concentration would shift the gas/liquid/solid transition curves toward the bottom of the plot. Right panel: Mean molecular weight variability in ice giants. We consider a mix of [85% H₂ + 15% He] for the background atmosphere composition corresponding to a 2.3 g mol⁻¹ mean molecular weight. We assume that CH₄ specific concentration follows saturation and its specific concentration in the deep atmosphere is set to 0.20 kg/kg for Uranus (equivalent to a 3.6% volume mixing ratio) and 0.30 kg/kg for Neptune (6.2% volume mixing ratio).</figcaption></figure>

to be modeled in order to understand the atmospheric dynamics of ice giants. A particularly crucial open question is related to the mechanisms of activation and inhibition of convection in those storms.

Among the properties they share, Uranus and Neptune both have atmospheres dominated by molecular hydrogen and helium, where all condensable species (CH₄, H₂S, NH₃, H₂O) are heavier than the background mixture of H₂ and He. Observations show a high abundance of methane in the troposphere, with significant latitudinal variations: 1–4% in Uranus at 2 bar (Sromovsky et al. 2014; Sromovsky et al. 2019), 2–6% in Neptune at 4 bar (Irwin et al. 2021; Tollefson et al. 2019). In ice giants, where the troposphere is located below the 0.1-bar level, methane is expected to condense around 1–2 bar (Fig. 1). The methane mixing ratio decreases along with pressure between 1–2 bar and 0.1 bar. As a consequence, the mean molecular weight in the atmosphere also decreases along with pressure in these layers (Fig. 1). In this study, we explore the 0.03–10 bar pressure range. Because the other minor species condense much deeper in the troposphere (e.g., H₂O around 100 bar) or are much less abundant (e.g., H₂S, above the 10-bar level, Moses et al. 2020), we consider only methane among the condensable species.

An abundant and heavy condensable species such as methane is thought to play an important role in the convective storms that occur in the troposphere. In the traditional view of convection, temperature alone sets the rules. Relative to the adiabatic gradient, warmer air rises, and colder air sinks. By contrast, in ice giants, both temperature and mean molecular weight control convection. Air at depth is warmer yet heavier; air higher up is colder yet lighter. This makes convection in ice giants more complex than in the traditional view. Guillot (2019) showed that

the potential temperature increase required to both compensate for the mean molecular weight change and have a density profile that is neutrally stable to convection is 20 K, though this does not include possible latent heat effects. During flybys of the ice giants, the Voyager 2 spacecraft measured vertical temperature profiles (Fig. 1), showing near dry-adiabatic gradients in the troposphere for both Uranus and Neptune (Lindal et al. 1987, 1990; Lindal 1992). Since the thermal gradient is close to the dry-adiabatic gradient, it is relevant to address the impact of the mean molecular weight vertical gradient, as it can completely change the (in)stability of the atmosphere.

Global climate models (GCMs) have been used to model the atmosphere of Jupiter (Boissinot et al. 2024) and Saturn (Guerlet et al. 2014; Li & Ingersoll 2015; Spiga et al. 2020; Bardet et al. 2021), highlighting large-scale phenomena but also the need to include mesoscale processes that GCM cannot directly resolve. Indeed, GCMs for ice giants have a horizontal resolution of 1° in latitude at best (equivalent to 400 km), in addition to being vertically hydrostatic. Because of this assumption of hydrostaticity, a GCM cannot be used to study convective storms. Solving and studying convection requires models with a higher spatial resolution and the capacity to solve the vertical momentum equation without making the hydrostatic approximation.

Sugiyama et al. (2014) and Li & Chen (2019) have studied moist-convection on Jupiter with cloud-resolving models. On Jupiter, the condensable species (H₂O, NH₃, H₂S) also have a molecular weight higher than H₂ and He. Sugiyama et al. (2014) conclude that stable layers associated with condensation act as effective dynamic boundaries; intense cumulonimbus clouds develop with clear temporal intermittence, with a period that is roughly proportional to the deep abundance of H₂O gas; the active transport associated with these clouds leads to the establishment of mean vertical profiles of condensates and condensable gases that are markedly different from the usual three-layer structure. However, in none of their simulations do condensable species (H₂O, NH₃, H₂S) exceed their critical specific concentrations (see the next section) for moist convection inhibition. The simulations of an idealized Jovian atmosphere in radiative-convective equilibrium presented by Li & Chen (2019) show that the temperature gradient is super-adiabatic near the water condensation level because of the change of mean molecular weight.

Using non-hydrostatic simulations applied to Uranus and Neptune, Ge et al. (2024) show that CH₄ and H₂S condensation induces two stably stratified layers at about 1 bar and 8 bar when the abundance of these elements range from 30 times solar to 50 times solar. They find that, in these stable layers, the temperature profile is super-adiabatic and convection is inhibited because of the compositional gradient in sub-saturated weather layers. More generally, they find that weakly forced giant planets are less cloudy than previously expected and that moist convection is limited by the planetary heat flux.

To study how the mean molecular weight variability affects convection, the challenge is to build a model that can resolve convection and account for the effect of condensation on the mean molecular weight. In this paper, using a 3D cloud-resolving model, we investigate the impact of vertical mean molecular weight gradients induced in particular by methane condensation on the formation and inhibition regimes of convective storms in ice giants.

In Sect. 2 we review the criteria for convection inhibition in ice giants. In Sect. 3, we present the model, how we adapted it for this study, and which challenges it implies. In Sect. 4 we describe the simulated atmospheric structures, and we analyze their
temporal evolution in Sect. 5. In Sect. 6, we discuss the limitations of our model and several open questions about the understanding of ice giants, and we attempt to give an overall scenario for storm formation based on the results of our simulations.

## 2. Convection inhibition criteria in ice giants

Guillot (1995) and Leconte et al. (2017) calculated analytical criteria for convection inhibition caused by a vertical gradient of the mean molecular weight. To determine if convection can occur for a given thermal gradient, the density of a theoretical rising parcel should be compared to its surrounding environment. The rising parcel follows the adiabatic gradient, dry if the pressure level is not saturated, and moist if it is. In a dry and well-mixed environment, the parcel continues to rise (i.e., convection occurs) if the thermal gradient is higher in absolute value than the dry-adiabatic gradient. The parcel remains warmer than the environment, with lower density and positive buoyancy. In the presence of a gradient of methane (i.e., there is more methane at the bottom than at the top), convection can be inhibited even if the thermal gradient is super-adiabatic.

At dry levels, if there is a gradient of methane, a rising parcel coming from the bottom, contains more methane than the air at the top. The density of the parcel is increased relative to the density of the environment. If the methane gradient is strong enough, the density of the parcel becomes higher than the density of the environment, even if the parcel is warmer than the environment. The buoyancy of the rising parcel is negative, stopping convection. The criterion for dry convection inhibition is the Ledoux criterion:

$$ \nabla_T < \nabla_{ad} + \nabla_{\mu}, \quad (1) $$

where $\nabla_T = \frac{d\ln T}{d\ln P}$ is the thermal gradient of the atmosphere, $\nabla_{ad} = \frac{\partial \ln T}{\partial \ln P}|_{ad}$ the dry-adiabatic gradient and $\nabla_{\mu} = \frac{d\ln \mu}{d\ln P}$ the mean molecular weight gradient of the atmosphere, $(T, P)$ being the temperature versus pressure profile of the atmosphere. Here the mean molecular weight gradient and the thermal gradient come as additive factors: the greater the mean molecular weight gradient, the greater the thermal gradient required to trigger convection.

In the “moist” atmosphere, methane condenses, producing a vertical methane gradient. The resulting methane profile is close to the saturation vapor curve. A rising parcel follows the moist-adiabatic gradient: it is becoming cooler, with a decreasing methane abundance because of condensation. When the mean molecular weight is vertically constant (i.e., in an atmosphere where the condensable species would have the same weight as the background atmosphere), if the thermal gradient is steeper than the moist-adiabatic gradient, then the parcel continues to rise (i.e., convection occurs). (Warming of the rising parcel by latent heat release is included in the moist-adiabatic gradient, which is lower than the dry-adiabatic gradient).

When the mean molecular weight decreases with pressure, if the thermal gradient is steeper than the moist-adiabatic gradient, the abundance of methane also decreases with pressure faster in the environment than in the parcel.

A super moist-adiabatic thermal gradient leads to an inner competition. On the one hand, the steep temperature gradient and the release of latent heat favor convection, and on the other hand varying mean molecular weight prevents convection. The rising parcel is warmer but also heavier than the surrounding environment. The abundance of methane determines which factor

dominates, and convection is inhibited if the methane abundance exceeds a critical specific concentration $q_{cri}(T)$. Leconte et al. (2017) provide the criterion. Moist convection is inhibited if

$$ (\nabla_T - \nabla_{ad}^*)(q_v - q_{cri}(T)) > 0, \quad (2) $$

where $\nabla_T$ is the thermal gradient of the atmosphere, $\nabla_{ad}^*$ is the moist-adiabatic gradient, $q_v$ (kg/kg) is the specific concentration of vapor, and $q_{cri}(T)$ is the critical specific concentration:

$$ q_{cri}(T) = \frac{1}{1 - \frac{M_{gas}}{M_{cs}}} \frac{R}{\Delta H_{cs}} T = 0.078 \frac{T}{80 K} (\text{kg/kg}), \quad (3) $$

where $M_{gas} = 2.3 \text{ g mol}^{-1}$ is the mean molecular weight of the non-condensing atmosphere (here a mix of 85% of H$_2$ and 15% of He), $M_{cs} = 16.04 \text{ g mol}^{-1}$ the molecular weight of the condensable species (methane in this study), $\Delta H_{cs} = 10 000 \text{ J mol}^{-1}$ is the latent heat of sublimation (or vaporization depending on the pressure-temperature range) of the condensable species, methane in this study, $R = 8.314 \text{ J K}^{-1} \text{ mol}^{-1}$ is the perfect gas constant.

If $q_v$ exceeds $q_{cri}(T)$, moist convection is inhibited, even if the thermal gradient is stronger than the moist-adiabatic gradient. In this case, the vapor abundance does not come as an additional factor, as it was through the gradient $\nabla_{\mu}$ in the criterion for inhibition of dry convection, but as a multiplicative factor ($q_v - q_{cri}(T)$). Contrary to dry convection inhibition, no thermal gradient, however strong, can produce moist convection. This critical specific concentration is linearly dependent on temperature. As the condensation level is around 80 K in ice giants, we can keep in mind the value of 0.078 kg/kg. This corresponds to a volume mixing ratio of 1.2%. In Leconte et al. (2017) the value of 0.10 kg/kg at 80 K is proposed. This value was calculated with the latent heat of vaporization (instead of sublimation). As shown by Fig. 1, considering solid-gas equilibrium (i.e., using the latent heat of sublimation in the formula) is more adequate. In Fig. 2, we plot several possible methane vapor gradients constrained by observations and saturation curves. We can see that the criterion of moist convection inhibition should apply at some depth and latitude on both planets.

## 3. A 3D cloud-resolving model

The Generic Planetary Climate Model (Generic PCM) gathers different versions of a common structure to model the atmospheres of planets and moons of our Solar System (Spiga & Forget 2009; Lefèvre et al. 2017; Bardet et al. 2021; Boissinot et al. 2024), as well as exoplanets (Turbet et al. 2022). PCM is the new name of the planetary versions of the model known so far as the Laboratoire de Météorologie Dynamique (LMD) model. Our model is one of these versions.

To build it, two components are coupled:

- a dynamical core, which solves the Euler equations of motion, and

- a physical package (by physical, we mean everything that is not related to dynamics), which calculates the tendencies of the relevant physical phenomena.

### 3.1. The dynamical core – the Weather Research and Forecasting Model

In this study, the dynamical core used was adapted from the Weather Research and Forecasting (WRF) Model. This model
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Possible methane vapor profiles in ice giants. In red we plot several possible methane vapor profiles constrained by observations from Sromovsky et al. (2014) and Irwin et al. (2021) (we plot the range of the latitudinal variations) and by saturation curves corresponding to Voyager 2 temperature profiles. Criterion curves in black indicate when moist convection is inhibited: when the specific concentration is higher than the critical specific concentration and the level is saturated, then the criterion applies. These possible methane vapor profiles are inspired by the tested profiles in Sromovsky et al. (2019) who talk about “descended profiles”.</figcaption></figure>

has been developed for a few decades for meteorological applications on Earth (large-eddy simulations) and used by a lot of meteorological research institutes. The version we used is the 4th one (Skamarock et al. 2021). To build our own model, we removed the physical package of the WRF model to keep only the dynamical core, which was coupled to our own physical package.

The dynamical core discretizes and solves the Euler equations of motion in a 3D rectangular grid. In addition to the classic fluid mechanics terms, these equations contain the physical tendencies of the relevant phenomena in the atmosphere, such as radiative transfer, for example, calculated by the physical package. As it is a cloud-resolving model designed for large-eddy simulations, with this chosen dynamical core, we can solve cloud formation. Clouds (identified here as saturated levels), being bigger than the resolution of the model, are spread over several grid cells. The WRF dynamical core solves dynamics and transport. It has several characteristics that are particularly interesting in our case:

- It uses the formalism by Laprise (1992) for the coordinates. The hydrostatic component of dry air pressure is used as the vertical coordinate.

- It is non-hydrostatic, so it can solve convection.

- It takes into account, in its equations, the variability of the mean molecular weight by including the variability of air density due to moisture (Leconte et al. 2024).

Some studies have already been done on Venus and Mars, using the same dynamical core with the adequate physical package from the PCM for each planet. Simulations on Mars (Spiga & Forget 2009; Spiga et al. 2017) demonstrated that localized convective snow storms can occur, and simulations on Venus (Lefèvre et al. 2017, 2018) reproduced the vertical position and thickness of the main convective cloud layer. In parallel with this study, the same association of the WRF dynamical core with the physical package we use here was made for K2-18b by Leconte et al. (2024) who give more details on the dynamical equations.

## 3.2. The physical package

The physical package is a “generic” physical package, which can be used for giant planets (Bardet et al. 2021; Guerlet et al. 2020),

paleoclimates of telluric planets and exoplanets (Turbet et al. 2022). Our physical package takes into account the following phenomena:

- Radiative transfer (absorption and emission by the gases with the correlated-k formalism; absorption, emission and scattering by aerosols layers; Collision-Induced Absorption (CIA); Rayleigh scattering). Considered species are H₂, He, CH₄, C₂H₂, C₂H₆.

- Methane thermodynamic cycle as a condensable species: condensation and sublimation with latent heat release, condensates precipitation.

- For Neptune, an internal heat flow (0.43 W m⁻², Guillot & Gautier 2015), which is a residual heating of the gravitational contraction of the planet, and is taken into account as a heat source at the bottom of the model.

For the radiative transfer, we used the parameters prescribed for ice giants by Milcareck et al. (2024), which integrate the aerosol scenario described by Irwin et al. (2022). Radiative-convective equilibrium 1D simulations produced temperature profiles close to the observed ones. Contrary to Milcareck et al. (2024) who used a fixed methane vertical profile, the methane abundance can vary in our simulations (due to condensation, sublimation, and transport). Radiative transfer calculations are updated during the simulation to take into account these methane variations.

The treatment of the methane cycle is generic and can be used for any condensable species, as it was done to study cloud formation on the Hot-Jupiter WASP-43b by Teinturier et al. (2024). This numerical thermodynamic scheme ensures that the methane abundance never exceeds the saturated value, which is derived from saturation pressure calculated with the Clausius-Clayperon formula. As soon as there is too much vapor at a given pressure-temperature level, the scheme condenses enough methane to bring the vapor amount below the saturation. The so-formed condensates then precipitate whenever their specific concentration exceeds a threshold that we kept arbitrarily small in this first study, so as not to have to take cloud radiative feedback into account. The condensates, which are ice crystals, precipitate and are transported instantaneously to deeper unsaturated layers where they sublimate. To do so, at each time step, the routine starts from the level where the condensates have formed
<table><caption>**Table 1.** Settings shared by all simulations.</caption><thead><tr><th>Parameters</th><th>Values</th></tr></thead><tbody><tr><td>Pressure range</td><td>0.03–10 bar</td></tr><tr><td>Vertical levels</td><td>200</td></tr><tr><td>Horizontal resolution (dx,dy)</td><td>2 km</td></tr><tr><td>Horizontal grid points</td><td>50×50</td></tr><tr><td>Insolation</td><td>0.93 W m<sup>-2</sup> for Uranus<br>0.38 W m<sup>-2</sup> for Neptune</td></tr><tr><td>Internal heat flow</td><td>None for Uranus<br>0.43 W m<sup>-2</sup> for Neptune</td></tr><tr><td>Heat capacity</td><td>10 200 J kg<sup>-1</sup> K<sup>-1</sup></td></tr><tr><td>Horizontal boundary conditions</td><td>cyclic</td></tr><tr><td>Duration</td><td>1 terrestrial year<br>(365×86 400 s)</td></tr><tr><td>Dynamical time step</td><td>5 s</td></tr></tbody></table>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Absorbed flux by Uranus' atmosphere. We tested two different methane abundances in the deep atmosphere (0.8% and 3.6%). We calculated the cumulative sum (from bottom to top) of the flux. The value at the top is the total absorbed flux, the so-called ASR.</figcaption></figure>

and carries them downward. At each unsaturated layer, the mass of methane that is needed to bring the layer back to saturation is computed accounting for the thermal effect of the sublimation. This mass is then sublimated into gas before repeating the process in the next layer until no more precipitations remain. We constrained the methane profile by a fixed value in the deep atmosphere ($q_{\text{deep}}(\text{kg/kg})$), and its value was chosen among the possibilities allowed by observations. This boundary condition works as an infinite source and sink of CH$_4$ that always maintains its abundance at the bottom of the model at this fixed value.

### 3.3. Simulation settings

The parameters common to all simulations are summarized in Table 1. We chose to set the bottom pressure level at 10 bar to have enough pressure levels below methane condensation. We set the top of the model at 0.03 bar to encompass the tropopause, which is just below the 0.1-bar level, and the lower stratosphere. We defined 200 vertical pressure levels from 10 bar to 0.03 bar, with an almost regular distribution in log pressure. This pressure distribution allows us to have enough levels where condensation occurs. It corresponds to an average vertical resolution of 0.75 km. The choice of the dynamical time step is constrained by Courant-Friedrichs-Lewy condition. In our case, 5 s was an optimized choice. With this choice, there is no need to set a vertical-velocity damping in the first 190 levels. In the top 10 levels, which correspond to the range 0.03–0.05 bar, we set an implicit gravity-wave “damping” layer ($\text{damp\_opt} = 3$, $\text{dampcoef} = 0.2$; these parameters are detailed in the Modeling System User's Guide of WRF whose reference is given in Skamarock et al. 2021) that smoothes out high speeds to prevent the model from collapsing. When velocity damping happens, energy is lost. However, as discussed in Sect. 2.4 of Leconte et al. (2024), this energy sink is rather small (less than 2% of the global budget). In our model, the stratosphere starts at the 0.2 bar level (see Sect. 4). Between the 0.2 bar level (the bottom of the stratosphere) and the 0.05 bar level (the bottom of the damping layer), there are 30 vertical grid points. These 30 points are enough for gravity waves to propagate. Stratospheric dynamics are not studied here, as this artificial damping reduces and modifies them, and because the chosen pressure range and discretization of the model are not suited to study them. The last layer of the model at the 0.03 bar level is directly subjected to the solar flux, without attenuation by the non-simulated

layers above, in order to have a complete radiative balance. The WRF dynamical core requires a few specific settings. We have decided to turn off any subgrid parameterization of diffusion ($khdif = 0$, $kvdif = 0$) and hyperviscosity ($c_s = 0$) to ensure that any mixing that we see in the simulations is due to the resolved dynamics of the dynamical core. The only diffusive processes left in the dynamical core were divergence damping ($\text{smdiv} = 0.1$) and external model filters ($\text{emdiv} = 0.01$) that are necessary to stabilize the simulations. These are canonical values recommended by the WRF documentation.

To set the heat capacity, which is constant in our model, we chose a value that gives the best fit of a dry-adiabatic temperature profile to Voyager 2's temperature profiles with an associated molecular weight of non-condensing air (85% H$_2$ + 15% He) $M_{\text{gas}} = 2.3 \text{ g mol}^{-1}$. We take $c_p = 10 200 \text{ J kg}^{-1} \text{ K}^{-1}$. The chosen horizontal resolution is 2 km. It is a good compromise between high resolution and sufficient horizontal coverage, with a 100 km-width domain (50×2 km). Given the horizontal size, the effects of the Coriolis force were not studied (they appear at larger scales). Diurnal and seasonal effects are turned off. Their impact is addressed in Sect. 6. We ran the simulations with the diurnal-averaged insolation on the planet, which corresponds to 1/4 of the solar flux reaching the planet's orbit. While Uranus has a very low internal heat flow ($0.042_{-0.042}^{+0.047} \text{ W m}^{-2}$, Pearl et al. 1990; Pearl & Conrath 1991; Guillot & Gautier 2015), Neptune has a higher one ($0.433 \pm 0.046 \text{ W m}^{-2}$). We set Uranus' internal heat flow to zero and Neptune's one to $0.43 \text{ W m}^{-2}$. Calculations of radiative transfer in our model show that a bit less than $0.1 \text{ W m}^{-2}$ of the incoming solar flux ($0.93 \text{ W m}^{-2}$) penetrates the atmosphere below the condensation level (see Fig. 3). For a low methane abundance in the deep atmosphere (0.8%), there is even about $0.05 \text{ W m}^{-2}$ absorbed below the model bottom. These fluxes are enough to trigger convection in the deep atmosphere. In the model, this energy is assumed to be absorbed at the model bottom and is sufficient to power convection from there.

### 3.4. Initialization and convergence

Ice giants are weakly forced systems because of their long radiative timescales and the little insolation they receive. We therefore expect long convergence times – decades to centuries – for the 3D simulations to reach a steady state. A simulation of one terrestrial year with our 3D model requires two weeks of computation on a cluster. Simulating a few terrestrial years would
<table><caption>**Table 2.** Varying parameters related to methane between simulations.</caption><thead><tr><th>Planet</th><th>q<sub>deep</sub> – methane specific concentration<br>In the deep atmosphere</th><th>Equivalent methane volume mixing ratio<br>In the deep atmosphere</th><th>Moist convection<br>Inhibition criterion</th></tr></thead><tbody><tr><td>Uranus</td><td>0.05 kg/kg</td><td>0.8%</td><td>Unsatisfied</td></tr><tr><td>Uranus</td><td>0.20 kg/kg</td><td>3.6%</td><td>Satisfied</td></tr><tr><td>Neptune</td><td>0.05 kg/kg</td><td>0.8%</td><td>Unsatisfied</td></tr><tr><td>Neptune</td><td>0.30 kg/kg</td><td>6.2%</td><td>Satisfied</td></tr></tbody></table>

**Notes.** The methane specific concentration (kg/kg) is the mass of methane divided by the total mass. The methane volume mixing ratio (%) is the number of methane molecules divided by the number of the other (non-condensing) molecules (H₂ and He).

consequently extend the duration of computation, but it would still remain a short time compared to the radiative timescales. We decided to limit ourselves to this duration of one terrestrial year, being aware that our simulations may not reach a statistical steady state. This limitation must be kept in mind when interpreting the results, as discussed Sect. 6.

To mitigate this issue, the simulations are started as close as possible to the envisioned equilibrium state. Leconte et al. (2024) simulated temperate exo-Neptunes for which thermal equilibrium can be reached much faster, because of shorter radiative timescales. They show that a simple 1D model can well predict the behavior of the thermal profile of 3D simulations. To initialize the 3D simulations with 1D thermal profiles, we did preliminary work on 1D simulations at radiative-convective equilibrium using their approach. These 1D simulations use the same single-column physics package as the 3D model, for the radiative and microphysical considerations. Concerning the dynamics of the 1D simulations, we parameterized convective adjustment and turbulent vertical diffusion. The convective adjustments are performed by two schemes. The 1D dry convective adjustment scheme brings back any decreasing profile of virtual potential temperature (this concept is introduced in the next section) to a constant profile, in dry layers, and mixes the methane accordingly. The 1D moist convective adjustment scheme is triggered in regions where (i) the thermal gradient is steeper than the moist-adiabatic gradient, (ii) methane is at saturation, and (iii) the methane abundance is lower than the critical abundance discussed in Sect. 2. In these regions, the thermal gradient is brought back to the moist one and methane is condensed accordingly.

The 3D simulations that we present in the sections hereafter are initialized with 1D simulations that use these parametrizations. To test the sensitivity to the initial conditions, 3D simulations initialized with 1D simulations that only use a dry convective adjustment were also run and are discussed in Sect. 6.

In order to represent a diverse but limited set of conditions, we chose different configurations for initial methane profiles. In 1D simulations, the parameterizations (convective adjustment and turbulent vertical diffusion) and the constraints due to the saturation vapor curve, build the methane profile. It results in a constant profile at the value set by $q_{deep}$, from the bottom of the simulation to the condensation level, and above this level, the profile follows the saturation vapor curve. In this 1D initialization profile, the methane mixing ratio above the cold trap is constant at its value at the cold trap. Some observations show a decreasing methane gradient in Uranus' stratosphere and an excess of methane in Neptune's stratosphere (Lellouch et al. 2015). We did not take these variations into account, as their study would require other simulation settings.

The parameter $q_{deep}$ is the one that allows us to test different configurations that are inspired by both observations and analytical criteria. Observations show a latitudinal minimum around the 2-bar level of about 1% methane on Uranus and 2% methane on Neptune. The critical mixing ratio for moist convection inhibition being 1.2% (at 80 K), we chose to include a case study with only 0.8% methane ($q_{deep} = 0.05$ kg/kg) in order to be below the critical mixing ratio at all pressures. In the case of Neptune, observations show that this configuration might not exist; however, we kept it as an experiment. We also simulated atmospheres with more methane in the deep atmosphere: 3.6% for Uranus ($q_{deep} = 0.20$ kg/kg) and 6.2% for Neptune ($q_{deep} = 0.30$ kg/kg). Although there might be a (high) pressure level at which CH₄ is latitudinally homogeneous, our simulations did not attempt to capture this aspect – this is discussed later. Finally, we ran four different simulations in ice giants (two on Neptune and two on Uranus). Table 2 synthesizes the varying parameters related to methane between these simulations.

After initialization with a 1D profile, we ran 3D simulations for one terrestrial year. They need a few simulated months to reach a steady state. From the one terrestrial year simulation, we kept only the last 200 days, from day 150 to day 350, for our study.

### 3.5. 3D validation tests

Before using our 3D model in specific configurations, we ran more theoretical 3D simulations to check the application of our model to the ice giants. The first test simulates a dry atmosphere, where we remove the thermodynamic cycle of methane as a condensing species, only keeping it as a radiatively active species. The atmosphere is then only a radiative/dry-convective atmosphere. As expected, this test simulation exhibits a dry convective layer driven by radiative heating at depth overlain by a stratosphere. The second test simulates a saturated moist atmosphere where we treat methane as a condensing species, but do not account for its different mean molecular weight. Again, as expected, this test simulation exhibits the standard three-layer structure, with a moist convective layer between the dry troposphere and the stratosphere. These two tests are conclusive and confirm that the model can be configured to simulate the ice giants.

## 4. Structure and dynamics in simulated atmospheres

In this section, we describe the thermal structures arising in our simulations as well as dynamics and convective activity, then we highlight the role of convection inhibition in those properties.
<figure><img src="image_5.png" /><figcaption>Fig. 4. Temporally and horizontally averaged temperature profile, temperature difference between two simulations, and potential and virtual poten- tial temperature profile from simulations.</figcaption></figure>

4.1. Mean temperature profiles, dynamics and convective activity

Temperature profiles allow us to check if the simulations are capable of capturing the phenomena we want to study and if the model provides results in line with observations. In this section, we look at the temporally and horizontally averaged temperature profiles of the 3D simulations.

Mean temperature profiles are close to Voyager 2 profiles and similar to those obtained by Milcareck et al. (2024). Above condensation levels, temperature profiles are moist-adiabatic. In simulations where we expect inhibition of moist convection and the formation of a stable layer (Uranus 0.20 kg/kg CH$_{4}$ and Neptune 0.30 kg/kg CH$_{4}$), temperature profiles are super moist- adiabatic at condensation levels between the 1 and 2-bar levels. However, they remain sub dry-adiabatic (see Fig. 4). This struc- ture is intimately related to the initialization profiles due to the long radiative timescale. Another recent study of the tro- pospheres of ice giants by Ge et al. (2024) exhibits a super dry-adiabatic profile at these levels. We come back to these dif- ferences in Sect. 6. Profiles in Ge et al. (2024) have similar temperature ranges for Neptune, although they are a few Kelvin colder.

The potential temperature of a fluid parcel at a given pressure is the temperature the parcel would reach if it were brought adi- abatically to a standard reference pressure (the 10-bar pressure at the bottom of the model in our case). In this study, the virtual potential temperature is more useful. Virtual potential tempera- ture θ$_{v}$ is the theoretical potential temperature θ of dry air that would have the same density as moist air:

$$
\theta_v = \left(1 - \left(1 - \frac{1}{\epsilon_{\text{CH}_4}}\right) q_v\right) \theta, \qquad (4)
$$

where ϵ$_{CH}$$_{4}$ = M$_{CH}$$_{4}$ M$_{gas}$ = 16.04 2.3 = 6.97 is the ratio between the molec- ular weight of the condensable species, in this case methane, and the molecular weight of non-condensing air; q$_{v}$ (kg/kg) is the vapor specific concentration. As ϵ$_{CH}$$_{4}$ is greater than 1, θ$_{v}$ is lower than θ.

Variations of q$_{deep}$ from one simulation to another show no difference on the average temperature profiles above the 1-bar level (see Fig. 4). Below the 1-bar level, simulations with a high q$_{deep}$ are colder than simulations with a low q$_{deep}$. These average temperature profiles are highly dependent on the initial 1D pro- file. The absorption and emission of methane, and its impact on radiative-convective equilibrium, may explain these differences in temperature gradients. Ge et al. (2024) also show that methane enrichment in the deep atmosphere leads to colder temperature

profiles. Virtual potential temperature profiles also show differ- ences, in line with expectations. For simulations with q$_{deep}$ = 0.20 and 0.30 kg/kg, condensation level is reached deeper than for simulations q$_{deep}$ = 0.05 kg/kg. We expect dry convection to extend higher in simulations with q$_{deep}$ = 0.05 kg/kg. All simu- lations share the same global three-layer structure (see Fig. 5 for Uranus and Fig. 6 for Neptune):

– a dry layer between 10 and 1–2 bar, where the virtual poten- tial temperature is almost constant (Fig. 4). In this layer methane abundance is almost constant, and equals the value we have fixed in the deep atmosphere. Regions with con- stant virtual potential temperature correspond to well-mixed layers (CH$_{4}$ remains at the same concentration) close to the dry-adiabatic gradient. Dry convection permanently occurs in this layer, with typical speeds of about 1 m/s, as shown by the typical situation plot (second panel of Figs. 5 and 6). Downdrafts and updrafts can be observed. This is a mixing layer. The upper part of this layer may be non-convective as a thin methane gradient is formed to link the two con- straints that delimit this layer. The constraint at the bottom is the q$_{deep}$ quantity at 10 bar that is transported upward by convection, and the constraint at the top is the condensation level.

– a moist layer between 1–2 bar and 0.1 bar, where methane abundance is close to the saturation vapor curve and vir- tual potential temperature slightly increases. This layer is the moist troposphere where moist convection episodically occurs. Most of the time, no convection occurs in this layer, as shown by typical situation plots. Sometimes, convective storms occur, as shown by “storm event” plots (fourth panel of Figs. 5 and 6). Convective storms are characterized by positive speeds in the moist layer (that can be higher than 10 m/s), and a strong downdraft in the dry layer (the large blue cell on the fourth panel of Figs. 5 and 6). Gravity waves can be identified above the storm. They transport energy. When a storm occurs, the moist layer is saturated and rel- ative humidity reaches 100% in most of the vertical levels of the moist layer where the storm is located.

– an upper layer above the 0.1-bar level up to the top of the model at 0.03 bar, where virtual potential temperature increases more strongly. No convection is expected in this layer but gravity waves can propagate upward.

The simulated atmosphere alternates between a typical situation structure and a storm event structure. In the typical situation structure, there is an expected local minimum of relative humid- ity in the moist layer. The temperature profile of this layer is confined between the dry-adiabatic temperature gradient and the moist-adiabatic temperature gradient. The temperature gradient
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Uranus – synthesis of simulations. Simulations with 0.05 kg/kg methane (~0.8% methane in number) and 0.20 kg/kg methane (~3.6% methane in number) in the deep atmosphere are presented. The first plot shows methane specific concentration averaged over the domain for the “typical situation” profile and taken at a given column of the model for the “storm event” profile. For the speed color map of the “typical situation”, we chose an arbitrary snapshot that illustrates well what the typical situation looks like. This typical situation is sometimes perturbed: a convective storm occurs and the atmosphere has the storm event structure. For the color map of the storm event structure, we chose an arbitrary snapshot of a convective storm among the ones occurring during the whole simulation. On the fourth panel, the color bar for the vertical speed is linear in the range [-1,1] and becomes logarithmic in the ranges [-100,-1] and [1,100].</figcaption></figure>

is not steep enough for dry convection and the layer is not saturated for moist convection to occur. This layer remains as such until it reaches saturation and convective storms occur. The condensation level below and the cold trap above surround this minimum of relative humidity.

## 4.2. Convection inhibition in simulations

The analysis of the dynamics and convective activity in the three-layer structure highlights that convection can be inhibited. In the dry layer, methane abundance is constrained at both ends of the layer: by the fixed methane abundance $q_{\text{deep}}$ in the deep atmosphere and by the saturation at the condensation level. In most of the dry layer, methane abundance is constant and is equal to $q_{\text{deep}}$. At the top of the dry layer, methane abundance starts to decrease before reaching the condensation level. This thin part of the dry layer with a methane gradient is a “non-convective

boundary layer” where the criterion for dry convection inhibition applies. This can be seen in the bottom panels of Fig. 7, which are zooms of the first panel of Figs. 5 and 6. Dry convection stops slightly before the condensation level, marking the top of the dry layer. This gradient at the top of the dry layer is about 0.5 bar thick. It could be more or less thick as shown by the possible methane gradients in Fig. 2. This non-convective layer acts as an obstacle for upward methane transport, limiting the rise of methane from the dry layer to the moist layer. This rise is required for convective storms.

The critical specific concentration for moist convection inhibition ($q_{\text{cri}}(T)$) is plotted (dotted line) on the first panels of Figs. 5 and 6. For simulations run with $q_{\text{deep}} = 0.05$ kg/kg, on both Uranus and Neptune, methane abundance is always below the critical value $q_{\text{cri}}(T)$ at all pressures. Consequently, if we examine the fourth panel of Figs. 5 and 6 (i.e., storm event plots) for 0.05 kg/kg CH$_4$, we can see that moist convection occurs in
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Neptune – synthesis of simulations. Simulations with 0.05 kg/kg methane (~0.8% methane in number) and 0.30 kg/kg methane (~6.2% methane in number) in the deep atmosphere are presented. Same layout as in Fig. 5.</figcaption></figure>

the entire moist layer. Moist convection is never inhibited. For simulations run with $q_{\text{deep}}$ = 0.20 and 0.30 kg/kg, the methane abundance may exceed the critical value. A new layer, between the dry layer and the moist layer, appears, with a relative humidity very close to 100%. This is a non-convective layer (Fig. 7), we call it the “Stable Layer”.

In the simulations with $q_{\text{deep}}$ = 0.20 and 0.30 kg/kg, that is, when the criterion for moist convection inhibition is satisfied, the three-layer structure becomes a four-layer structure with this stable moist and non-convective layer appearing between the dry layer and the moist layer. This non-convective layer is a bit more than 0.5 bar thick in our simulation of Uranus with 0.20 kg/kg methane and almost 1 bar thick in our simulation of Neptune with 0.30 kg/kg methane. The more methane in the deep atmosphere, the thicker this layer is. It acts both as an obstacle and a “reservoir”. As moist convection is inhibited in this layer, methane is difficult to transport to higher levels (to levels where moist convection is not inhibited and where convective storms can occur). In addition, the top of this reservoir (the level denoted with an X marker in Fig. 7 must be completely saturated in methane to let convective storms occur in the moist layer (see the next section).

**5. Formation of convective storms and intermittency**

While the previous section focused on averaged structures, we describe here their evolution: how convective storms are formed, what is their frequency and released kinetic energy.

**5.1. Advection of methane: from the dry to the moist layer**

Convective storms are made of rising saturated air. Figs. 5 and 6 show the structure of storms. In a few columns of the domain where the storm occurs, the moist layer is filled with methane and reaches 100% of relative humidity. As the moist layer is confined between the dry-adiabatic gradient and the moist-adiabatic gradient, dry convection can never happen and moist convection happens when saturation is reached. Convective storms form when enough methane from the deep troposphere is brought to lower pressures, at condensation level where inhibition criteria are not satisfied. When these conditions are gathered, a perturbation can thus create a convective storm. The mechanical energy of the convective storm is dispersed by gravity waves propagating upward. The convective storm ends with condensed methane precipitating deeper.
<figure><img src="image_8.png" /><figcaption>Fig. 7. Average structures observed in Uranus and Neptune simulations. First line: average relative humidity in Uranus and Neptune simulations. The markers (O) and (X) indicate the minimum and maximum relative humidity levels in the moist layer, which are studied in Sect. 5. Second line: zoom between 0.5 bar and 4 bar. Methane specific concentration is plotted with a linear scale (instead of log-scale as before). A three-layer structure appears for all simulations with a “Dry Layer” (dry troposphere), a “Moist Layer” (moist troposphere), and an “Upper Layer” (bottom of the stratosphere). A fourth “Stable Layer” appears when the methane abundance in the deep atmosphere exceeds the critical specific concentration. The analytical criteria explain the convective and non-convective parts of these layers. The zoom (bottom panels of Fig. 7) highlights the methane gradient that exists in the upper part of the dry layer just below condensation level and that forms a “non-convective boundary layer”. On the panels of zoom plots where $q_{\text{deep}} = 0.05$ kg/kg, we do not see the critical specific concentration because the plotted range is lower than it.</figcaption></figure>

Dry convection in the dry layer transports methane to lower pressures. Dry convection faces an obstacle at the top of the dry layer that we call the non-convective boundary layer. To cross this obstacle, methane has to be transported by slow eddy diffusion and has to “climb” the methane gradient. This obstacle slows down the transport of methane to the moist layer and thus controls the frequency of convective storms. In the simulations with $q_{\text{deep}}$ lower than $q_{\text{cri}}(T)$ (Uranus 0.05 kg/kg CH$_4$ and Neptune 0.05 kg/kg CH$_4$), this gradient below condensation level is slowly becoming smaller with time until it is small enough for a perturbation coming from deeper to reach the moist layer and produce a convective storm. In the simulations with $q_{\text{deep}}$ exceeding $q_{\text{cri}}(T)$, the stable layer is filling up with methane, acting as a reservoir. The reservoir becomes full when saturation reaches its top, where methane concentration is lower than its critical value. And moist convection can occur.

To quantify methane transport, we estimate the equivalent mixing coefficient (the so-called eddy diffusivity or $K_{zz}$) defined as

$$K_{zz} = \frac{\langle \rho q_v w \rangle}{\langle \rho \partial_z q_v \rangle}, \qquad (5)$$

where $\rho$ is the density (kg/m$^3$), $q_v$ the CH$_4$ vapor specific concentration (kg/kg), $w$ the vertical speed (m/s), $\partial_z$ the partial derivative along the vertical axis, $\langle \cdot \rangle$ the average over horizontal and temporal dimensions.

The $K_{zz}$ profiles shown in Fig. 8 are consistent with the dynamics observed in the simulations. In the mixing dry layer, methane transport is efficient and $K_{zz}$ is of the order of $10^2$ to $10^4$ m$^2$/s. In the moist layer, intermittent convection is illustrated by the peak in the profiles around 0.5 bar and $K_{zz}$ is on the order of 1 m$^2$/s. At the interface between these two layers, there is a low $K_{zz}$ of about $10^{-1}$ m$^2$/s. $K_{zz}$ is particularly low in the stable layer of the simulations where the criterion for moist convection inhibition is satisfied (Uranus 0.20 kg/kg CH$_4$ and Neptune 0.30 kg/kg CH$_4$). Using another mesoscale model, Ge et al. (2024) also calculated $K_{zz}$ profiles that are close to ours: a high $K_{zz}$ in the dry troposphere, a low $K_{zz}$ in the upper layers where convection is inhibited, and a higher $K_{zz}$ in the layer where moist convection can occur.

In the dry layer, below the 2-bar level, maximum downward velocities are higher in absolute value than maximum upward velocities and correspond to the strong downdrafts forming when moist convection occurs (Fig. 9). In the moist layer, maximum
<table><caption>**Table 3.** Number of storms and average period between two storms in simulations.</caption><thead><tr><td>Simulations</td><td>Number of storms in 200 days</td><td>Average period between 2 storms</td></tr></thead><tbody><tr><td>Uranus – 0.05 kg/kg CH<sub>4</sub></td><td>2</td><td>100 days</td></tr><tr><td>Uranus – 0.20 kg/kg CH<sub>4</sub></td><td>5</td><td>40 days</td></tr><tr><td>Neptune – 0.05 kg/kg CH<sub>4</sub></td><td>4</td><td>50 days</td></tr><tr><td>Neptune – 0.30 kg/kg CH<sub>4</sub></td><td>36</td><td>6 days</td></tr></tbody></table>

<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Vertical profile of the equivalent vertical mixing coefficient $K_{zz}$ derived from the simulation.</figcaption></figure>

upward velocities are higher in absolute value than maximum downward velocities and correspond to the convective storms. In the simulations run above the critical specific concentration (Uranus 0.20 kg/kg and Neptune 0.30 kg/kg), we can see a strong minimum in the root-min-square (RMS) profiles around 1 bar. This minimum corresponds to the stable layer where the criterion for moist convection inhibition is satisfied. There is no such minimum in simulations below the critical specific concentration. We could see the RMS decrease slightly before reaching the 1-bar level, thus highlighting the non-convective boundary layer at the top of the dry layer, but the local minimum is not very marked.

## 5.2. Temporal evolution and intermittence

The existence of two different structures – a typical situation structure and a storm event structure – implies intermittency. The previous subsection has explained how the transition between these two structures happens. It works as a cycle.

To study this cycle, we identified on the average structures of Fig. 7 two important levels in the moist layer: a minimum (indicated by an O marker) and a maximum (indicated by an X marker) of relative humidity. This maximum is also the top of the reservoir. These levels allow us to identify when a storm occurs. The level of the relative humidity minimum is suddenly filled with methane. In Fig. 10, we plot the temporal evolution of the relative humidity at these two levels. While minima (O) remain at a constant value for a long time, at the same time maxima (X) are filled until they reach 100%, and then moisture is transferred upward, filling the upper levels, including the minimum level, which is more or less the epicenter of the convective storm. In simulations with 0.20 and 0.30 kg/kg CH$_4$, the levels just below the maximum (X) level and that form what we call the reservoir are also getting filled. Relative humidity at the maximum level

(the top of the reservoir) is always close to 100% in the simulations exceeding the critical specific concentration, that is, Uranus 0.20 kg/kg CH$_4$ and Neptune 0.30 kg/kg CH$_4$.

The frequency of storms varies a lot from one simulation to another. On the studied window of the last 200 days of the one terrestrial year duration of the simulations, we can count and estimate the period between two storms (Table 3). More storms occur in simulations where the critical specific concentration is exceeded: 2.5 times more storms in the simulation of Uranus with 0.20 kg/kg CH$_4$ than in the simulation of Uranus with 0.05 kg/kg CH$_4$, 9 times more storms in the simulation of Neptune with 0.30 kg/kg CH$_4$ than in the simulation of Neptune with 0.05 kg/kg CH$_4$. Storm occurrence is cyclic, with a chaotic aspect (i.e., the cycle is sometimes irregular), which is well illustrated in Neptune’s simulation with 0.05 kg/kg CH$_4$ (see the third panel of Fig. 10), where the time between two storms varies a lot. Storms are triggered when a sufficiently strong perturbation comes from deeper levels for the moisture of the reservoir to be transferred from the levels where moist convection is inhibited to the levels just above where moist convection is no longer inhibited.

Simulations show a much lower frequency of convective storms on Uranus than on Neptune. We attribute this to the absence of internal heat flow. This internal heat flow warms the deep atmosphere in simulations. It is a source of heat at the bottom of the model. There is more energy flux to evacuate from deep layers in Neptune than in Uranus. One might think that an interesting test would be to turn off the internal heat flow in the settings of Neptune’s simulations or to add one in Uranus’ simulations. We did this, but it changed so much the thermal balance so that the temperature gradients, and the conclusions, became unrealistic. Model and observations agree that Neptune is much more active in producing new cloud systems and cloud variations than Uranus (Hueso et al. 2020).

## 5.3. Intensity of convective storms

The evolution of relative humidity on plots from Fig. 10 allow us to identify each convective storm and monitor their frequency. However, these plots do not provide information about the intensity of storms. We computed the kinetic energy $E_c$ per unit of mass on the whole domain of our simulations as follows:

$$E_c(J/kg) = \frac{1}{2} \sum_i (u_i^2 + v_i^2 + w_i^2) \frac{dm_i}{m_{tot}}, \quad (6)$$

where $u_i$, $v_i$, $w_i$ are the three components of velocity of the i-th cell of the domain, $dm_i$ is the mass of this cell, $m_{tot} = \sum_i dm_i$ is the total mass of the domain. Fig. 11 illustrates the temporal evolution of $E_c$ in our simulations.

The average kinetic energy of simulations with more methane, that is, Uranus 0.20 kg/kg CH$_4$ and Neptune 0.30 kg/kg CH$_4$, is of the same order of magnitude as the average kinetic
<figure><img src="image_10.png" /><figcaption>**Fig. 9.** Vertical velocity profiles showing the RMS (dotted), maximum upward (dashed), and maximum downward (solid) velocities. Averages and maxima values are computed over temporal and horizontal dimensions.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Temporal evolution of the relative humidity in the moist layer at the minimum and maximum levels, respectively identified by the (O) and (X) markers in Fig. 7. For these levels, we kept only the highest value over the domain (latitudinally and longitudinally).</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>Fig. 11. Kinetic energy variations in the simulations calculated over the whole domain.</figcaption></figure>

energy of simulations with less methane – Uranus 0.05 kg/kg CH₄ and Neptune 0.05 kg/kg CH₄. But when looking at the convective storms that can be identified by the peaks of kinetic energy, we see that the intensity of these peaks is higher in simulations with 0.05 kg/kg CH₄. The non-convective layers of these configurations are thinner, but more importantly, the criterion for moist convection inhibition is never met. The release of energy is thus facilitated. There is a strong correlation between intermittency and intensity. The rarer the storms are the more intense they are. In Fig. 11, the kinetic energy of a storm corresponds to the difference between the peak associated with this storm and the “background” kinetic energy. This background kinetic energy is mainly due to dry convection in the troposphere and is of the same order of magnitude from one simulation to another.

**6. Open questions, limitations, and suggested scenario**

Here we present the limitations of our model and the associated issues. Knowing these limitations, we propose a scenario for storm formation in ice giants.

**Condensates micro-physics.** In our model, the micro-physics of condensation is limited to a very basic scheme. Growth, sedimentation, and sublimation of condensates would require a more detailed model. In this first study, we kept these processes simple and controlled by as few parameters as possible, which at least allows us to test extreme cases that bracket the actual behavior of condensates. The chosen microphysical parameters limit the potential retention of condensates (e.g., in moist updrafts) to 10⁻¹⁰ kg of methane ice per kg of air to avoid having to introduce a complex microphysical model. In reality, some condensates can be retained and slow down the updraft. To test whether this would affect our conclusions, we ran simulations where 10⁻³ kg/kg of methane ice can be retained before precipitating. As expected, the kinetic energy of moist convective events is reduced, by a factor of 2, when more condensates are allowed, which further favors convection inhibition. At the same time, the frequency of moist convective events is multiplied by a factor of 2, demonstrating again that frequency and intensity are strongly correlated.

<figure><img src="image_13.png" /><figcaption>Fig. 12. Evolution of temperature anomalies with respect to the initial thermal profile. The temperature anomaly is defined as the difference between the horizontally averaged temperature profile at each time and the initial thermal profile.</figcaption></figure>

**Convergence.** As discussed in Sect. 3, the radiative timescale is too long to enable running the simulation until complete thermal equilibration. As a result, even though we remove the spin-up phase of the simulations (i.e., the first 150 days) and analyze only the part of the simulation where the frequency of storms is rather stable, the thermal structure still evolves slightly over this period. Fig. 12 shows the evolution of the temperature anomaly with respect to the initial thermal profile for the four simulations. These anomalies remain relatively low and are mostly confined to the stable layer. This is to be expected because the thermal gradient in this region is determined by turbulent diffusion, which is one of the most difficult factors to take into
**Diurnal cycle and seasonal variations.** Our simulations used a constant averaged solar flux and zenith angle and did not include day/night alternation nor the daily and seasonal variations of solar zenith angle, which could have an impact. The insolation diurnal cycle should have little impact, as radiative time constants are very long: more than 100 terrestrial years at 1 bar. Concerning seasonal variations of solar flux, the period that we studied (200 days) was small compared to the Uranian and Neptunian years, so the solar flux reaching the planet was almost constant during that time. We could run simulations on the same planet with different constant averaged solar flux corresponding to different latitudes and seasons on the planet. This would change the thermal profile by several Kelvin and the level of methane condensation. It could increase or diminish storm frequency. Because of its obliquity, the case of Uranus is really interesting, and the study of seasonal variations will be the subject of future studies. Our model was run under average insolation conditions, and we tested different values of methane abundance in the deep atmosphere. One could build a pseudo-2D model (i.e., as a function of altitude and latitude) to study latitudinal variations in temperature and methane.

**Thermal gradient in stable layers.** In stable layers, our temperature profiles are super moist-adiabatic but remain sub dry-adiabatic. Using less computationally expensive 2D non-hydrostatic simulations that were run for a longer time, Ge et al. (2024) find a super dry-adiabatic temperature profile in the stable region and Leconte et al. (2024) find the same structure with a 3D model but for warmer atmospheres that have shorter radiative timescales. We thus believe that, given more integration time, our thermal profiles would probably converge toward a super dry-adiabatic one. This, however, should not affect our conclusions on the occurrence of storms as we showed that they are driven mainly by the methane cycle and the turbulence in the stable layer, which is only mildly affected by the thermal gradient in the stable layer.

**Global climate.** Simulations with a global climate model that accounts for the large-scale dynamics and seasonal changes would provide crucial information. But we still face major unknowns that prevent us from performing fully consistent large-scale simulations. For instance, it is not known whether the latitudinal variations in methane abundance observed at 1–2 bar on Uranus and 4 bar on Neptune are still valid at 10 bar. We do not know if there is a pressure at which methane abundance is homogeneous nor how these strong horizontal gradients affect the dynamics. Ideally, small-scale convective simulations should have evolving boundary conditions fed by large-scale simulations, while large-scale simulations should include a convective sub-grid scheme derived from small-scale simulations. There is quite some work to be done before reaching such a sophisticated coupling between models. At this early stage, our approach was to develop the small-scale simulations independently and in parallel with large-scale circulation models but to vary the deep methane concentration in order to capture the variety of conditions met at different latitudes. We do believe that by doing so, we were able to capture the mechanisms controlling storms at the scale of our simulations and to provide parametrizations for moist convection inhibition that can be used in GCMs.

**Stratospheric methane concentrations.** Another debated question is how to explain the abundance of methane in the stratosphere. Observations of Lellouch et al. (2015) show a mixing ratio of CH₄ higher than its value at the cold trap in Neptune's stratosphere, while methane abundance strongly decreases with pressure in Uranus' stratosphere. Though we did not seek to explain those particularities in this study, we do not see any methane transport by moist convection to the stratosphere. The overshoots that we simulated do not bring methane into the stratosphere either. The 3D GCMs are more suitable for such a study, as those variations could probably also be caused by large-scale dynamics. Mesoscale simulations that include the effect of the Coriolis force could also bring new pieces of information.

Aware of these limitations, we propose a description of the formation and structure of convective storms for given methane abundances in the deep atmosphere. Our simulations highlight a stable moist and non-convective layer in the regions where the critical specific concentration is exceeded. The thickness of this layer, which is located around the 1-bar level (Fig. 13), depends on the methane abundance in the deep atmosphere ($q_{deep}$). The frequency and intensity of storms depend on the presence or absence of this stable layer.

If the abundance at saturated levels is lower than the critical specific concentration, this stable layer does not exist. This situation is encountered in Uranus, at the poles according to observations. It is illustrated in the two plots done for methane-poor regions, on the edges of Fig. 13 (Uranus). A very thin line near the 1-bar level (in yellow on those plots) marks the maximum relative humidity and the separation between the dry layer and the moist layer. Convective storms occur when enough methane is brought up from deeper levels and are intense because moist convection is never inhibited. These intense storms reduce the relative humidity in the levels close to condensation by evaporating methane, and it then takes a long time to “reload” them with methane. The frequency of the storms is lower.

If this abundance is higher than the critical specific concentration at saturated levels, a stable layer appears (the larger $q_{deep}$ the thicker this stable layer). In Uranus, this situation corresponds to the three center plots of Fig. 13 (Uranus), and should occur at mid-latitudes and the equator. In Neptune, this situation is encountered in all plots of Fig. 13 (Neptune), with the variation of $q_{deep}$ inducing a variation in the thickness of the stable
<figure><img src="image_14.png" /><figcaption>**Fig. 13.** Relative humidity (slices of one horizontal dimension) during convective storms, from methane-rich latitudes to methane-poor latitudes, taken from our simulations in ice giants and inspired by latitudinal and vertical variations shown by observations. The simulations were run with different theoretical methane concentrations in the deep atmosphere. Uranus: Observations show latitudinal variations of 1–4% in methane at 2 bar (Sromovsky et al. 2019), from pole to equator. The relative humidity is plotted from simulations with three different methane specific concentrations: 0.05, 0.125, and 0.20 kg/kg. These correspond to volume mixing ratios of 0.8, 2.1, and 3.6% respectively. Neptune: Observations show latitudinal variations of 2–6% in methane at 4 bar (Irwin et al. 2021), from pole to equator. The relative humidity is plotted from our simulations with three different methane specific concentrations: 0.125, 0.20, and 0.30 kg/kg. These correspond to volume mixing ratios of 2.1, 3.6, and 6.2% respectively. The simulation with $q_{\text{deep}} = 0.05$ kg/kg is no longer plotted since this configuration might not exist on Neptune. We see powerful but rare storms at methane-poor latitudes and weak but frequent storms at methane-rich latitudes. The more methane there is, the more frequent storms are, but also the weaker they are. To make these plots, we ran three more short simulations with $q_{\text{deep}} = 0.125$ kg/kg CH$_4$ in Uranus, $q_{\text{deep}} = 0.125$ and 0.20 kg/kg CH$_4$ in Neptune, in addition to our simulations already presented and detailed in this article.</figcaption></figure>

layer. In this saturated layer of variable thickness, moist convection is inhibited, thus retaining methane and acting as a reservoir. This reservoir needs to be completely full for moist convective storms to occur above. Because of moist convection inhibition, this reservoir is always almost full and the frequency of convective storms is high (about a few days in Neptune when $q_{\text{deep}} = 0.30$ kg/kg). Also because of moist convection inhibition, the intensity of the storms is weak: the moist non-convective layers can never be driven upward by the storm above them. The

bottom of the storm is thus never deeper than this level where moist convection is no longer inhibited. The more methane there is, the more frequent storms are, but also the weaker they are.

We previously explained that convective storms should be more frequent in Neptune than in Uranus, due to the presence of an internal heat flow in Neptune. This might also be reinforced by the fact that the methane abundance in Neptune exceeds the critical abundance for the inhibition of moist convection at all latitudes, while in Uranus the methane abundance at the poles
might be lower than this critical abundance. As there is more methane in Neptune than in Uranus, convective storms should be even more frequent (but weaker).

The cyclic occurrence of storms in our simulations may suggest a link with the study by Li & Ingersoll (2015) on the frequency of Saturn's giant storms. The long timescale (about 60 years) and the resolution of convection that they used in their model are quite different from our study. Temperature variations in Saturn's atmosphere are much larger than in the ice giants (we report very few variations in our simulations), and are the main driver of the cyclic occurrence of storms, whereas, in our simulations, it is the evolution of the methane profile that is the driver. The frequencies that we calculated and presented in Table 3 should be treated with caution. These frequencies depend very much on the assumptions we make about methane microphysics. A better estimate of the frequencies may be obtained by a parametric study of our microphysical parameters. In the current state of our model, even if the frequency of storms is highly uncertain, we can explain the sporadicity of clouds observed on Uranus (Palotai et al. 2022; de Pater et al. 2015) and on Neptune (Karkoschka & Tomasko 2011).

As for the other condensable species in ice giants (H₂O, NH₃, H₂S), we expect the same behavior as for methane. However, some species are less abundant than methane (e.g., NH₃, H₂S, Moses et al. 2020), and should never exceed their associated critical mixing ratios for inhibition of moist convection. The behavior of these species could be similar to the simulations presented in this article with 0.05 kg/kg methane. In Uranus and Neptune, water should be very abundant in the deep atmosphere below the 100-bar level (Cavalié et al. 2017; Venot et al. 2020; Moses et al. 2020) and should exceed its associated critical mixing ratio. Even if it is difficult to extrapolate our simulated convective regimes in the range of 0.03–10 bar to such high pressures, we could expect behavior similar to the simulations presented in this article with 0.20 and 0.30 kg/kg. In Jupiter and Saturn, water abundance may exceed the critical mixing ratio associated with water. At certain latitudes, we would then expect behavior similar to the simulations presented in this article with 0.20 and 0.30 kg/kg.

## 7. Conclusions

Using a 3D cloud-resolving model, we investigated the impact of the change in mean molecular weight due to methane condensation on the formation and inhibition regimes of convective storms in ice giants. With methane being heavier than the H₂/He background, its condensation can indeed inhibit convection and moist convective storms.

Observations show both latitudinal variations – 1–4% in Uranus at 2 bar (Sromovsky et al. 2014; Sromovsky et al. 2019) and 2–6% in Neptune at 4 bar (Irwin et al. 2021) – as well as vertical variations caused by condensation. Vertical variations at non-saturated levels (i.e., dry levels) strongly stabilize the atmosphere, and a super dry-adiabatic gradient – which is convective in a mixed atmosphere – can remain stable. The literature (Guillot 1995; Leconte et al. 2017) highlights the existence of a critical methane abundance at saturated levels. Whether or not this critical abundance is exceeded can inhibit or activate moist convection. Depending on the form that the methane gradient takes, saturated levels may or may not be above this critical abundance, which is 1.2% at 80 K and corresponds approximately to the 1-bar level. Latitudinal variations observed at this level are almost all above this critical abundance, but these levels have to

be saturated for the criterion of moist convection inhibition to apply.

After having shown in our simulations that this critical methane abundance indeed governs convective storm inhibition and formation, we conclude the following:

- Typical velocities of dry convection in the deep atmosphere are rather low (on the order of 1 m/s) but are sufficient to sustain upward methane transport;

- Moist convection at the methane condensation level is strongly inhibited;

- Convective storms can form regardless of methane abundance in the deep atmosphere, but they can only form in saturated layers where the methane abundance does not exceed the critical abundance;

- The formation of convective storms on Uranus and Neptune should be intermittent and follow a loading cycle;

- The intermittency and intensity of storms depend on the methane abundance:

- Where CH₄ exceeds the critical abundance in the deep atmosphere (at the equator and the middle latitudes on Uranus and at all latitudes on Neptune), more frequent but weak storms form;

- Where CH₄ remains below the critical abundance in the deep atmosphere (possibly at the poles on Uranus), storms are rarer but more powerful;

- Storms should be more frequent on Neptune than on Uranus because of the internal heat flow of Neptune and because there is more methane in Neptune than in Uranus;

- Methane-rich latitudes at saturated (or near-saturated) levels should act as a barrier, allowing little energy to be released in one storm, while methane-poor latitudes should allow much more energy to be released in one storm.

These conclusions can explain the sporadicity of clouds observed in ice giants. Further observations with the James Webb Space Telescope, which can track moist convective events over a longer observational period or provide new constraints on methane abundance, could help bring new insight to the conclusions proposed in this article.

*Acknowledgements.* The authors acknowledge the support of the French Agence Nationale de la Recherche (ANR), under grant ANR-20-CE49-0009 (project SOUND).

## References

Bardet, D., Spiga, A., Guerlet, S., et al. 2021, Icarus, 354, 114042

Boissinot, A., Spiga, A., Guerlet, S., Cabanes, S., & Bardet, D. 2024, A&A, 687, A274

Cavalié, T., Venot, O., Selsis, F., et al. 2017, Icarus, 291, 1

de Pater, I., Sromovsky, L. A., Fry, P. M., et al. 2015, Icarus, 252, 121

Ge, H., Li, C., Zhang, X., & Moeckel, C. 2024, Planet. Sci. J., 5, 101

Guerlet, S., Spiga, A., Sylvestre, M., et al. 2014, Icarus, 238, 110

Guerlet, S., Spiga, A., Fouchet, T., & Delattre, H. 2020, Icarus, 351, 113935

Guillot, T. 1995, Science, 269, 1697

Guillot, T. 2019, Exp. Astron., 54, 1027

Guillot, T., & Gautier, D. 2015, Giant Planets (Elsevier), 529

Hueso, R., Guillot, T., & Sánchez-Lavega, A. 2020, Philos. Trans. Roy. Soc. Lond. Ser. A, 378, 20190476

Irwin, P. G. J., Dobinson, J., James, A., et al. 2021, Icarus, 357, 114277

Irwin, P. G. J., Teanby, N. A., Fletcher, L. N., et al. 2022, J. Geophys. Res. (Planets), 127, e2022JE007189

Karkoschka, E., & Tomasko, M. G. 2011, Icarus, 211, 780

Laprise, R. 1992, Monthly Weather Rev., 120, 197

Leconte, J., Selsis, F., Hersant, F., & Guillot, T. 2017, A&A, 598, A98

Leconte, J., Spiga, A., Clément, N., et al. 2024, A&A, 686, A131

Lefèvre, M., Spiga, A., & Lebonnois, S. 2017, J. Geophys. Res. (Planets), 122, 134

Lefèvre, M., Lebonnois, S., & Spiga, A. 2018, J. Geophys. Res. (Planets), 123, 2773

Lellouch, E., Moreno, R., Orton, G. S., et al. 2015, A&A, 579, A121
Li, C., & Ingersoll, A. 2015, Nat. Geosci., 8, 398

Li, C., & Chen, X. 2019, ApJS, 240, 37

Lindal, G. F. 1992, AJ, 103, 967

Lindal, G. F., Lyons, J. R., Sweetnam, D. N., et al. 1987, J. Geophys. Res., 92, 14987

Lindal, G. F., Lyons, J. R., Sweetnam, D. N., et al. 1990, Geophys. Res. Lett., 17, 1733

Milcareck, G., Guerlet, S., Montmessin, F., et al. 2024, A&A, 686, A303

Molter, E., de Pater, I., Luszcz-Cook, S., et al. 2019, Icarus, 321, 324

Moses, J. I., Cavalié, T., Fletcher, L. N., & Roman, M. T. 2020, Philos. Trans. Roy. Soc. Lond. Ser. A, 378, 20190477

Palotai, C., Brueshaber, S., Sankar, R., & Sayanagi, K. 2022, Remote Sensing, 15, 219

Pearl, J. C., & Conrath, B. J. 1991, J. Geophys. Res., 96, 18921

Pearl, J. C., Conrath, B. J., Hanel, R. A., Pirraglia, J. A., & Coustenis, A. 1990, Icarus, 84, 12

Skamarock, W. C., Klemp, J. B., Dudhia, J., et al. 2021. A Description of the Advanced Research WRF Model Version 4.3 (No. NCAR/TN-556+STR), https://doi.org/10.5065/1dfh-6p97

Spiga, A., & Forget, F. 2009, J. Geophys. Res. (Planets), 114, E02009

Spiga, A., Hinson, D. P., Madeleine, J.-B., et al. 2017, Nat. Geosci., 10, 652

Spiga, A., Guerlet, S., Millour, E., et al. 2020, Icarus, 335, 113377

Sromovsky, L. A., Karkoschka, E., Fry, P. M., et al. 2014, Icarus, 238, 137

Sromovsky, L. A., Karkoschka, E., Fry, P. M., de Pater, I., & Hammel, H. B. 2019, Icarus, 317, 266

Sugiyama, K., Nakajima, K., Odaka, M., Kuramoto, K., & Hayashi, Y. Y. 2014, Icarus, 229, 71

Teinturier, L., Charnay, B., Spiga, A., et al. 2024, A&A, 683, A231

Tollefson, J., de Pater, I., Luszcz-Cook, S., & DeBoer, D. 2019, AJ, 157, 251

Turbet, M., Fauchez, T. J., Sergeev, D. E., et al. 2022, Planet. Sci. J., 3, 211

Venot, O., Cavalié, T., Bounaceur, R., et al. 2020, A&A, 634, A78