CAMEMBERT: A Mini-Neptunes General Circulation Model Intercomparison, Protocol Version 1.0.A CUISINES Model Intercomparison Project

Super-Earths and mini-Neptunes represent a demarcation, albeit a nebulous one, between the giant planets with their thick atmospheres dominated by hydrogen and helium and the terrestrial planets with thinner secondary atmospheres (Lopez & Fortney 2014). The planets that have retained their hydrogen-dominated atmospheres are believed to have undergone runaway accretion during their formation in order to accumulate a thick atmosphere (Pollack et al. 1996; Lee et al. 2014) but have also been able to retain some or all of that atmosphere in the presence of irradiative evaporation (Owen & Jackson 2012).

Simulating these planets with general circulation models (GCMs) presents a unique set of challenges not necessarily seen in the Earth sciences community. The primitive equations that assume hydrostatic balance, a thin atmosphere, and a constant gravitational acceleration with height, may not yield accurate results for cases where the thickness of the modeled atmosphere becomes significant relative to the radius of the planet, limiting their applicability to planets with thick atmospheres (see Tokano 2013; Tort et al. 2015, and Mayne et al. 2019 for discussions related to Venus and Titan, Earth, and mini-Neptunes, respectively). It has also been argued that simulations of mini-Neptunes specifically may have extremely long convergence times, potentially of 50,000 Earth days of model time or more (Wang & Wordsworth 2020), which raises questions about the accuracy of simulations of only a few 1000 days, as are common in the exoplanetary modeling community (e.g., Charnay et al. 2015; Mayne et al. 2017, 2019). These issues motivate a better understanding of the tools we use to model these planets.

While intercomparison studies have been somewhat common in the Earth sciences community (see, e.g., Eyring et al. 2016; Haarsma et al. 2016; Pincus et al. 2016; Ullrich et al. 2017), it has not been until recently that intercomparisons have been done with a focus on exoplanetary targets. Although not an intercomparison of multiple GCMs, Heng et al. (2011) performed a comparison of the spectral and finite-difference dynamical cores in the Geophysical Fluid Dynamics Laboratory Princeton Flexible Modeling System (FMS) using the hot Jupiter HD 209458 b as a test case. The first true intercomparison of GCMs used in exoplanetary science, Polichtchouk et al. (2014), looked at highly idealized configurations—a steady state jet, a baroclinic wave, and diabatic forcing—for five GCMs with the intention of better understanding their respective dynamical cores. Increasing the complexity, Yang et al. (2019) compared GCMs in the context of terrestrial planets, specifically focusing on the cases of an Earth-like planet orbiting a G star and a tidally locked planet around an M star. More recently, the TRAPPIST-1 Habitable Atmosphere Intercomparison (THAI; Fauchez et al. 2020) compared GCM models of TRAPPIST-1e, investigating the dynamics (Sergeev et al. 2022; Turbet et al. 2022) and the synthetic observations (Fauchez et al. 2022) resulting from the simulations. In the final paper, they propose a "GCM uncertainty error bar" of ∼50% when interpreting transmission spectra with the uncertainty explained mostly by the cloud differences found between GCMs. The ability to provide this form of context to synthetic observations highlights the value of projects like THAI.

Based on the success of the THAI project, we propose here the Comparing Atmospheric Models of Extrasolar Mini-Neptunes Building and Envisioning Retrievals and Transits (CAMEMBERT) intercomparison of GCMs modeling mini-Neptunes under the umbrella of the Climates Using Interactive Suites of Intercomparisons Nested for Exoplanet Studies (CUISINES) framework for intercomparisons for exoplanets (L. E. Sohl et al. 2022, in preparation). The broad objectives of CUISINES are twofold. First, it provides a metaframework to quantify, and potentially mitigate, differences between exoplanet model outputs. Second, it aims to assess how these output differences affect the synthetic observations that are used to predict the detectability of atmospheric constituents and to interpret data from ground and space telescopes. With the increased focus on mini-Neptunes with the launch of TESS (e.g., Trifonov et al. 2019; Burt et al. 2021; Lacedelli et al. 2021), CHEOPS (e.g., Bonfanti et al. 2021; Leleu et al. 2021), JWST (e.g., Greene et al. 2017; Bean et al. 2021; Hu & Damiano 2021), and in anticipation of increased efforts to model these planets, we believe the timing is appropriate for an intercomparison of GCMs modeling mini-Neptunes to provide a foundational understanding of how our models behave and how our model choices may impact the interpretation of observations.

In this first Paper, we outline the protocol for the CAMEMBERT model intercomparison project, providing both the motivations for the test cases as well as sufficient details to reproduce them, in the hope that the results of this intercomparison can be used as a calibration for future GCMs. While tests of the protocol were run using the Unified Model (UM) and Exo-FMS to gauge the viability of the protocol, the results from these tests are not presented here, with the intention being to present the results from all participating GCMs in one or more follow-up papers. The outline of this Paper is as follows: An overview of the GCMs currently participating in the intercomparison and a discussion of the target planets are found in Section 2. The protocol and associated test cases are described in Section 3. The outputs and diagnostics are found in Section 4 and a final summary and discussion is found in Section 5.

With the goal of the intercomparison being to understand the differences between GCMs, the progression of simulations can be divided up into two stages: First, a simple investigation of the dynamics and the dynamical cores at the heart of each of the GCMs. GCM simulations of GJ 1214b, for example, show up to three zonal jets, with differing amplitudes depending on the study (Menou 2012; Kataria et al. 2014; Charnay et al. 2015; Drummond et al. 2018; Mayne et al. 2019; Wang & Wordsworth 2020). In understanding the origins of these differences, it is essential to be able to disentangle differences between dynamical cores and differences between radiative transfer schemes. Once this baseline understanding is established, we progress to a second stage where we compare models with radiative transfer.

The GCMs currently participating, based on an expressed interest at the BUFFET workshop ²² and during the protocol development process, are listed in Table 1. Although these are the participants at the time of the publication of this protocol, participation remains open and other teams are welcome to join.

2.1. Choice of Targets

In this section we outline the simulations associated with the intercomparison and their motivations. In general, we seek to maintain consistency in parameters throughout the protocol to facilitate comparison as complexity is increased. While the various cases were initially envisioned as a progression, we encourage participants to join the cases their GCMs are capable of completing regardless of whether or not their GCMs are capable of completing cases earlier in the series.

For all models, we adopt a lower boundary pressure of 3 ×10⁶ Pa. While previous studies have studied higher pressures—Charnay et al. (2015) and Mayne et al. (2019) used 8 × 10⁶ Pa and 2 × 10⁷ Pa, respectively—to simplify the test cases and limit differences between the primitive and less simplified equations of the dynamics for mini-Neptunes (Mayne et al. 2019), we limit ourselves to 3 × 10⁶ Pa. We adopt an upper boundary pressure p ≤ 10 Pa for codes that employ a pressure-based vertical grid. For GCMs that use a height-based vertical grid (e.g., the UM), we require a domain height sufficient to include pressures of 10 Pa throughout the simulations; however, this may require varying the domain height on a per-case basis. This will also result in those GCMs using height-based grids having potentially significantly lower pressures on the nightside upper boundary, which may necessitate modifications to ensure stability. Following the lead of the THAI project, we do not place specific requirements on the time step or grid spacing, instead encouraging participants to adopt parameters they would commonly use for exoplanet studies as requirements for stability may differ between GCMs.

For simplicity, we assume friction-free, impermeable boundaries to avoid complicating the tests with boundary-layer friction or mass exchange, and heat exchange is limited to a fixed internal heat flux with effective temperature T_int. We do not, however, exclude forms of dissipation that may be required for numerical stability (e.g., sponge layers, artificial viscosity).

As mini-Neptunes may have enhanced atmospheric metallicities relative to solar (Fortney et al. 2013), we adopt a value of 100× solar for K2-18b and use parameters and profiles consistent with this throughout the intercomparison. While there may be similar motivation to use this value for GJ 1214b, due to the number of previous simulations using solar metallicity and the differences between simulations already shown in the literature, specifically in the number and speed of the zonal jets, we adopt the solar value as it allows us to probe a part of the parameter space where we know differences between GCMs already occur.

Simulations are run for a fixed number of Earth days instead of specifying a convergence condition. The chosen simulation lengths do not ensure convergence in simulations that include the deep atmosphere; this may not be sufficient time for the deep atmosphere to have converged to a steady state (Mayne et al. 2017). However, as we are limited to regions with p ≤ 3 × 10⁶ Pa we should not be significantly impacted. While Wang & Wordsworth (2020) have found in their simulations of GJ 1214b that over integration times of 50,000 Earth days or more their model atmospheres transition from two off-equatorial jets to a single equatorial jet; including such long integration times in an intercomparison would likely prove computationally prohibitive and limit participation. We instead focus initially on comparing GCMs over shorter timescales, with hopes of extending the work in the future to look at these longer timescales. For cases 1 and 3, we run the simulations for 4000 Earth days. For case 2, however, we run each simulation for 10,000 Earth days as the dual-gray case is the case investigated by Wang & Wordsworth (2020). Although this does not approach the long integration times of Wang & Wordsworth (2020), an integration time of 10,000 Earth days may be sufficient to understand differences between GCMs as Menou (2012) observed the formation of the central equatorial jet in their 7800 Earth day simulation. It may be inevitable, however, that an understanding of a possible delayed formation of an equatorial jet may have to wait for a follow-up study with fewer participants investigating longer simulation times.

Although mini-Neptunes have the potential for clouds, and GJ 1214b in particular has been shown to have strong signs of clouds or hazes (Kreidberg et al. 2014), we do not include a cloudy benchmark as a part of the protocol as currently there are an insufficient number of GCMs capable of participating. We do hope that follow-up studies and intercomparisons will be able to include a cloud component as clouds will undoubtedly represent an important constituent of future mini-Neptune models.

3.1. Initial Conditions

To initialize our simulations, we use 1D pressure-temperature profiles (see Figure 1) with no initial winds. For K2-18b, we use a profile generated using Exo-REM from Charnay et al. (2021) and for GJ 1214b we use a profile generated using ATMO from Drummond et al. (2018). For all cases, we do not include any initial latitudinal or longitudinal variation. Each of these profiles, along with the profiles for the chemical abundances needed for case 3, are publicly available as an ASCII text file in the CAMEMBERT repository (see Section 4.2).

Zoom In Zoom Out Reset image size

3.2. Case 1: Temperature Forcing

The first case investigated as a part of the intercomparison is similar to the temperature-forced benchmark of Held & Suarez (1994). The motivation is to compare the robustness of the dynamical cores without also comparing differing radiative transfer schemes. In this case, we opt to use a Newtonian cooling prescription where the temperature T is forced to an equilibrium temperature profile T_eq on a given radiative timescale τ_rad

$\begin{eqnarray}&&\displaystyle \frac{{dT}}{{dt}}=\displaystyle \frac{{T}_{\mathrm{eq}}-T}{{\tau }_{\mathrm{rad}}},\end{eqnarray} \tag{ 1 }$

where t is time. We adopt a radiative timescale τ_rad for a hydrogen-dominated atmosphere (Zhang & Showman 2017)

$\begin{eqnarray}{\tau }_{\mathrm{rad},{{\rm{H}}}_{2}}(p)=\left\{\begin{array}{ll}{10}^{4}\,{\rm{s}} & p\leqslant {10}^{2}\,\mathrm{Pa}\\ {10}^{5/2}{p}^{3/4}\,{\rm{s}} & {10}^{2}\,\mathrm{Pa}\lt p\lt {10}^{6}\,\mathrm{Pa}\\ {10}^{7}\,{\rm{s}} & p\geqslant {10}^{6}\,\mathrm{Pa}\end{array}\right.,\end{eqnarray} \tag{ 2 }$

$\begin{eqnarray}&&{\tau }_{\mathrm{rad}}(p)={\tau }_{\mathrm{rad},{{\rm{H}}}_{2}}(p)\left(\displaystyle \frac{{c}_{p}}{7R/2}\right)\displaystyle \frac{2}{\mu },\end{eqnarray} \tag{ 3 }$

with the specific values of R, c_p , and μ found in Table 2.

Table 2. Input Parameters

	GJ 1214b	K2-18b
Common Planetary Parameters
Mass (kg)	4.88 × 1025	5.15 × 1025
Radius (m)	1.75 × 107	1.66 × 107
Orbital radius (au)	0.014 85	0.159
Orbital and rotation period (days)	1.58	32.94
Gravity (m s−2)	10.7	12.4
Common Stellar Parameters
Mass (kg)	2.98 × 1029	7.14 × 1029
Radius (m)	1.50 × 108	3.09 × 108
Teff (K)	3250	3457
Metallicity [Fe/H]	0.29	0.123
log(g)	5.026	4.8
Case 1
${\rm{\Delta }}{T}_{\mathrm{eq},\max }$ (K)	600	50
Specific gas constant R (J kg−1 K−1)	3.513 × 103	1.732 × 103
Specific heat capacity cP (J kg−1 K−1)	1.200 × 104	6.682 × 103
Mean molar mass μ (g mol−1)	2.367	4.801
Case 2
Shortwave absorption κsw (m2 kg−1)	1 × 10−4	2 × 10−5
Longwave absorption κlw (m2 kg−1)	3 × 10−3	1.4 × 10−2
Instellation (W m−2)	2.17 × 104	1.37 × 103
Tint (K)	100	90
Case 3
Stellar spectrum	3000 K BT-Settl	3500 K BT-Settl
	with [Fe/H] = 0.3	with [Fe/H] = 0.0
	and log(g) = 5	and log(g) = 5

Note. Planetary parameters for GJ 1214b are from Cloutier et al. (2021) and parameters for K2-18b are from Benneke et al. (2019) and Cloutier et al. (2019). Stellar parameters are for GJ 1214 and K2-18 are taken from Cloutier et al. (2021) and Benneke et al. (2017), respectively. To avoid ambiguities in the mass and radius values used in the intercomparison, the input values are quoted in units of kilograms and meters, respectively, with the conversion being done using the appropriate conversion factors in Prša et al. (2016).

Download table as:  ASCIITypeset image

The equilibrium 3D temperature profile T_eq is generated from the initial temperature profile T₀(p) and temperature difference ΔT_eq(p) intended to mimic qualitatively the day–night temperature contrast expected from a tidally locked planet

$\begin{eqnarray}{T}_{\mathrm{eq}}=\left\{\begin{array}{ll}{T}_{0}(p)+{\rm{\Delta }}{T}_{\mathrm{eq}}(p)\left(| \cos \phi \cos \lambda | -\displaystyle \frac{1}{2}\right) & \mathrm{dayside}\\ {T}_{0}(p)-\displaystyle \frac{1}{2}{\rm{\Delta }}{T}_{\mathrm{eq}}(p) & \mathrm{nightside}\end{array}\right.,\end{eqnarray} \tag{ 4 }$

where the longitude is λ and the latitude is ϕ. The temperature difference ΔT_eq is taken to be ${\rm{\Delta }}{T}_{\mathrm{eq},\max }$ at p ≤ 10 Pa and to decrease linearly with $\mathrm{log}p$ until p = 10⁶ Pa where is becomes zero. The range of equilibrium temperatures at a given pressure is shown in Figure 1. For GJ 1214b, the contrast parameter ${\rm{\Delta }}{T}_{\mathrm{eq},\max }$ is chosen to maintain consistency with the temperature forcing tests of Mayne et al. (2019). As published temperature forcing tests do not exist for K2-18b, we instead look to Charnay et al. (2021), which shows a more modest ∼50 K temperature contrast at the top of the atmosphere for their 100× solar metallicity atmosphere, motivating the chosen value of ${\rm{\Delta }}{T}_{\mathrm{eq},\max }$.

We note that the temperature-forced case of K2-18b has presented significant difficulty in terms of numerical stability for the UM and Exo-FMS in tests of the protocol. Rather than remove it, we retain it as a part of the protocol with a note of caution to participants.

3.3. Case 2: Gray Radiative Transfer

3.4. Case 3: Nongray Radiative Transfer with Fixed Abundance Profiles

To model atmospheric chemistry, we limit ourselves to H₂/He-dominated atmospheres with H₂O, CH₄, NH₃, CO, and CO₂ as well as H₂ and He collisionally induced absorption and Rayleigh scattering as opacity sources. Volume mixing ratios for each species as a function of gas pressure, taken from the same simulations that generated the initial conditions, are shown in Figure 2 and are provided as a part of the initial conditions archive (see Section 3.1). All participating GCMs are to use these abundance profiles, as this allows GCMs without coupled chemistry solvers to participate.

Zoom In Zoom Out Reset image size

For the stellar spectra, we use the model stellar spectra from PHOENIX BT-Settl (Allard et al. 2012), which closest matches the target star (see Table 2). These spectra are made available in the CUISINES repository along with the other required inputs. The specific line lists, calculation method, and spectral resolution are left to the individual groups.

In this Section we outline the procedures for formatting and archiving data. The goal is to standardize the output and storage of data as much as possible to facilitate not only the initial analysis but also future analysis by third parties.

4.1. Outputting and Formatting of Data

4.2. The CAMEMBERT Repositories

4.3. Simulated Observables

4.4. Environmental Impact

In this Paper, we have presented the protocol for the CAMEMBERT project, which seeks to compare GCMs used by the exoplanetary science community, with models of mini-Neptunes being the primary focus. Two benchmarks were chosen—the warm GJ 1214b and the relatively cooler K2-18b—based on the volume of prior modeling work and observational potential, and a series of simulations of increasing complexity are described to calibrate and compare the participating GCMs, with all of the requisite parameters provided here and in the CUISINES repository. Membership in CAMEMBERT remains open, and other groups interested in participating are invited to contact the authors. Collaboration meetings will be held in 2022 and beyond, and the aim is to present the results of the intercomparison in one or more follow-up papers. It is hoped that the results from this intercomparison will provide a strong foundation for follow-up studies exploring elements not included in this initial protocol such as chemistry, clouds, and convergence timescales with any model differences due to these elements being more easily isolated and interpreted. As CUISINES brings together researchers employing a diverse range of tools and approaches, we anticipate that it will act to create new collaborations and stimulate progress in understanding exoplanets. For CAMEMBERT specifically, connections with Modeling Atmospheric Lines By the Exoplanet Community (MALBEC), an intercomparison of radiative transfer codes, and Photochemical model Intercomparison for Exoplanet science (PIE), an intercomparison of 1D photochemistry codes, will begin in 2022 with the goal of comparing results and better informing our own models with the insights from the other projects.

CAMEMBERT belongs to the CUISINES metaframework, a Nexus for Exoplanet System Science (NExSS) science working group. The testing of the protocol was done using Met Office Software. The contributions of D.A.C., N.J.M., D.S., and M.Z. made use of the ISCA High Performance Computing Service at the University of Exeter and the DiRAC Data Intensive service at Leicester, operated by the University of Leicester IT Services, which forms part of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The equipment for DIRAC was funded by BEIS capital funding via STFC capital grants ST/K000373/1 and ST/R002363/1 and STFC DiRAC Operations grant ST/R001014/1. DiRAC is part of the National e-Infrastructure. The work of D.A.C., N.J.M., D.S., and M.Z. was also partly funded by the Leverhulme Trust through a research project grant (RPG-2020-82), a Science and Technology Facilities Council Consolidated Grant (ST/R000395/1), and a a UKRI Future Leaders Fellowship (grant No. MR/T040866/1). T.J.F. acknowledges support from the GSFC Sellers Exoplanet Environments Collaboration (SEEC), which is funded in part by the NASA Planetary Science Divisions Internal Scientist Funding Model. E.K.H.L and P.A.N. are supported by the SNSF Ambizione Fellowship grant (#193448). H.I. is supported by funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement No. 740963). The work of B.C. was granted access to the HPC resources of MesoPSL financed by the Region Ile de France and the project Equip@Meso (reference ANR-10-EQPX-29-01) of the program Investissements d'Avenir supervised by the Agence Nationale pour la Recherche. The authors also thank Aymeric Spiga for proposing the CAMEMBERT acronym.