Exploring the deep atmospheres of HD 209458b and WASP-43b using a non-gray general circulation model

Aaron David Schneider; Ludmila Carone; Leen Decin; Uffe Gråe Jørgensen; Paul Mollière; Robin Baeyens; Sven Kiefer; Christiane Helling

doi:10.1051/0004-6361/202142728

Home

All issues

Volume 664 (August 2022)

A&A, 664 (2022) A56

Full HTML

Free Access

Issue		A&A Volume 664, August 2022


Article Number		A56
Number of page(s)		20
Section		Planets and planetary systems
DOI		https://doi.org/10.1051/0004-6361/202142728
Published online		05 August 2022

A&A 664, A56 (2022)

Exploring the deep atmospheres of HD 209458b and WASP-43b using a non-gray general circulation model

Aaron David Schneider¹^,2, Ludmila Carone³, Leen Decin², Uffe Gråe Jørgensen¹, Paul Mollière³, Robin Baeyens², Sven Kiefer²^,4^,5 and Christiane Helling⁴^,5

¹ Centre for ExoLife Sciences, Niels Bohr Institute, Øster Voldgade 5, 1350 Copenhagen, Denmark
e-mail: aaron.schneider@nbi.ku.dk
² Institute of Astronomy, KU Leuven, Celestijnenlaan 200D, 3001 Leuven, Belgium
³ Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany
⁴ Space Research Institute, Austrian Academy of Sciences, Schmiedlstrasse 6, 8042 Graz, Austria
⁵ TU Graz, Fakultät für Mathematik, Physik und Geodäsie, Petersgasse 16, 8010 Graz, Austria

Received: 23 November 2021
Accepted: 18 February 2022

Abstract

Simulations with a 3D general circulation model (GCM) suggest that one potential driver behind the observed radius inflation in hot Jupiters may be the downward advection of energy from the highly irradiated photosphere into the deeper layers. Here, we compare dynamical heat transport within the non-inflated hot Jupiter WASP-43b and the canonical inflated hot Jupiter HD 209458b, with similar effective temperatures. We investigate to what extent the radiatively driven heating and cooling in the photosphere (at pressures smaller than 1 bar) influence the deeper temperature profile (at pressures between 1 to 700 bar). Our simulations with the new non-gray 3D radiation-hydrodynamical model expeRT/MITgcm show that the deep temperature profile of WASP-43b is associated with a relatively cold adiabat. The deep layers of HD 209458b, however, do not converge and remain nearly unchanged regardless of whether a cold or a hot initial state is used. Furthermore, we show that different flow structures in the deep atmospheric layers arise. There, we find that WASP-43b exhibits a deep equatorial jet, driven by the relatively fast tidally locked rotation of this planet (0.81 days), as compared to HD 209458b (3.47 days). However, by comparing simulations with different rotation periods, we find that the resulting flow structures only marginally influence the temperature evolution in the deep atmosphere, which is almost completely dominated by radiative heating and cooling. Furthermore, we find that the evolution of deeper layers can influence the 3D temperature structure in the photosphere of WASP-43b. Thus, dayside emission spectra of WASP-43b may shed more light onto the dynamical processes occurring at greater depths.

Key words: radiation: dynamics / radiative transfer / scattering / planets and satellites: atmospheres–planets / satellites: gaseous planets

© ESO 2022

1 Introduction

The James Webb Space Telescope (JWST) is set to reveal data of unprecedented quality for several hot to ultra-hot Jupiters. Thanks to their inflated structure and hotter emission temperatures, these bodies are the easiest to observe, making them prime targets for further characterization. Furthermore, JWST observations may even yield insights into one of the long-standing mysteries in exoplanet science, which asks why many of the hot to ultra-hot Jupiters have large radii of 1.5 R_Jup (inflated), while others have compact atmospheres (not inflated).

A number of processes have been proposed to explain hot Jupiter inflation, including a reduction in the cooling efficiency or injection of additional energy into the interior. However, mechanisms that reduce the cooling of exoplanets fail to maintain inflated radii over timescales of billion years (see the reviews of Zhang 2020; Fortney et al. 2021). Recent studies have disfavored the reduced cooling hypothesis and point, rather, to the injection of stellar energy into the interior as a viable explanation. Examples include evidence of re-inflation in two close-in hot Jupiters as a result of their receiving more irradiation from their host stars during the course of their evolution (Hartman et al. 2016) as well as a correlation between the amount of radius inflation and incident stellar flux (Thorngren & Fortney 2018; Thorngren et al. 2019).

There are currently two favored processes for energy injection as an explanation for hot Jupiter inflation: (a) Ohmic dissipation as initially proposed by Batygin & Stevenson (2010), where partially ionized wind flows are coupled with the planetary magnetic field, releasing heat in the interior; and (b) dynamical downward energy transport from the irradiated atmosphere into deeper layers (Tremblin et al. 2017; Sainsbury-Martinez et al. 2019,2021). In this study, we focus on the latter mechanism.

Understanding the vertical transport of stellar energy from irradiated atmospheres into deeper layers is computationally challenging as it requires investigating the temperature and the dynamical evolution of atmospheric layers below the observable photosphere (pressures much larger than 1 bar) with large runtimes of »10 000 simulation days (Mayne et al. 2019; Mendonça 2020; Wang & Wordsworth 2020). Recent studies have tried to tackle this computational issue by invoking additional assumptions at greater depths or by using a simplified radiative transfer treatment (e.g., Tremblin et al. 2017; Mayne et al. 2017; Mendonça 2020; Carone et al. 2020; Sainsbury-Martinez et al. 2019, 2021).

Mendonça (2020) performed a detailed analysis of angular momentum and heat transport using a Fourier decomposition of the angular momentum transport. The author found that the mean circulation deposits heat from the upper atmosphere into the deeper layers, where he concludes that additional energy sources are needed to explain the radius inflation. Using much longer integration times and a parametric radiative transfer scheme without heating in the deep atmosphere, Sainsbury-Martinez et al. (2019) instead found that radius inflation can be explained by dynamical heat transport if the deep atmosphere is not heated by irradiation. Both of these studies, along with that of Showman et al. (2020), offer the conclusion that a self-consistent radiative transfer treatment is needed for the investigation of the deep layers to further quantify the coupling of the upper and lower atmospheres. Thus, there is a clear need for an efficient coupling mechanism between hydrodynamics and radiative transfer that is capable of monitoring deep processes over long simulation runtimes. We thus present expeRT/MITgcm¹, a fast and versatile radiative transfer scheme coupled with the 3D hydrodynamical code MITgcm.

The implementation of expeRT/MITgcm is described in Sect. 2. We apply the model to the exoplanets HD 209458b and WASP-43b in Sect. 3 and we explore the influence of layers deeper than 1 bar on the 3D thermal structure. HD 209458b is a well-observed and frequently modeled hot Jupiter, and WASP-43b is a JWST early release science target (ERS) (Venot et al. 2020; Bean et al. 2018). More crucially, despite having similar effective temperatures between 1400-1500 K, the first planet is inflated, whereas the latter is not (e.g., Helling et al. 2021, Table 1 and Fig. 1). We compare the behavior of the temperature convergence of WASP-43b and HD 209458b in the deep atmosphere in Sect. 4. We show the emission and transmission spectra for WASP-43b and HD 209458b in Sect. 5. We discuss in Sect. 6 the need for GCMs that are capable of tackling different complex processes in hot Jupiters. Finally, we provide a summary and conclusion in Sect. 7.

2 Methods

In this section, we describe the implementation of radiative transfer into expeRT/MITgcm. We first describe the general circulation model in Sect. 2.1 which is followed by a description of the radiative transfer extension introduced in this work (Sect. 2.2). The initialization of our simulations is discussed in Sect. 2.3.

2.1 Dynamical Core Set-Up

General circulation models (GCMs) are numerical frameworks that solve a (sub)set of the 3D equations of hydrodynamics on a rotating sphere. We integrated expeRT/MITgcm with the deep wind GCM framework introduced in Carone et al. (2020), which uses the dynamical core of the MITgcm (Adcroft et al. 2004). The MITgcm uses an Arakawa C type cubedsphere grid with resolution C32² to solve the 3D hydrostatic primitive equations (see e.g., Showman et al. 2009, Eqs. (1)–(5)) assuming an ideal gas. We used the C32 grid with 24 cores in MPI multiprocessing, resulting in four 16 × 16 tiles for each of the six cubedsphere faces. As in Showman et al. (2009), we applied a fourth-order Shapiro filter with τ_shap = 25 s to smooth horizontal grid-scale noise. The use of such smoothing schemes is common in GCMs, but may have a severe impact on atmospheric flows. We refer to Heng et al. (2011); Skinner & Cho (2021) for a discussion on the impact of smoothing schemes. Our vertical grid uses 41 logarithmically spaced grid cells between 1×10⁻⁵ bar and 100 bar. We extended the logarithmic domain by six linearly spaced grid cells between 100 bar and 700 bar in order to avoid a large spacing in the deep atmosphere.

We stabilized our model against nonphysical gravity waves by using a “soft” sponge layer that is similar to the one used in THOR (Mendonça et al. 2018b; Deitrick et al. 2020). The zonal horizontal velocity u [ms⁻¹] is damped by a Rayleigh friction term towards its longitudinal mean baru [m s⁻¹] via: $\frac{d u}{d t} = - {\tilde{k}}_{top} (u - bar u),$ ${{{\rm{d}}u} \over {{\rm{d}}t}} = - {\tilde k_{{\rm{top}}}}\left( {u - {\rm{bar}}u} \right),$ (1)

where t [s] is the time and ${\tilde{k}}_{top} [s^{- 1}]$ ${\tilde k_{{\rm{top}}}}\left[ {{{\rm{s}}^{ - 1}}} \right]$ is the strength of the opposed friction and is calculated as a function of pressure p [Pa] by ${\tilde{k}}_{top} (p) = k_{top} \cdot \max {[0, 1 - {(\frac{p}{p_{sponge}})}^{2}]}^{2} .$ ${\tilde k_{{\rm{top}}}}\left( p \right) = {k_{{\rm{top}}}} \cdot \max {\left[ {0,1 - {{\left( {{p \over {{p_{{\rm{sponge}}}}}}} \right)}^2}} \right]^2}.$ (2)

The fudge parameters p_sponge [Pa] and k_top [s⁻¹] are used to control the position and the strength of the opposed Rayleigh friction. Throughout this paper, we use k_top = 20 d⁻¹ and p_sponge = 1 × 10⁻⁴ bar, which results in a friction term that slowly increases in strength for p ≤ 1 × 10⁻⁴ bar. We want to note here that we calculate Eq. (1) by deprojecting the cubed sphere grid onto its geographic direction, averaging the deprojected u over 20 latitudinal bins and projecting the resulting baru on the cubedsphere grid.

To parametrize the deep magnetic drag and to stabilize the lower boundary, we followed Carone et al. (2020). We imposed additional Rayleigh friction to the winds in the deep layers (p ≥ 490 bar). Within these layers we damp the zonal velocity u [m s⁻¹] and meridional velocity υ [m s⁻¹] by $\begin{matrix} \frac{d u}{d t} = - \frac{u}{{\tilde{τ}}_{deep}}, & \frac{d υ}{d t} = - \frac{υ}{{\tilde{τ}}_{deep}}, \end{matrix}$ $\matrix{ {{{{\rm{d}}u} \over {{\rm{d}}t}} = - {u \over {{{\tilde \tau }_{{\rm{deep}}}}}},} & {{{{\rm{d}}\upsilon } \over {{\rm{d}}t}} = - {\upsilon \over {{{\tilde \tau }_{{\rm{deep}}}}}},} \cr }$ (3)

where we set the strength of the drag ${\tilde{τ}}_{deep}$ ${\tilde \tau _{{\rm{deep}}}}$ [s] as $\frac{1}{τ_{deep}} = \frac{1}{{\tilde{τ}}_{deep}} \max (0, \frac{p - 490 bar}{700 bar - 490 bar}) .$ ${1 \over {{\tau _{{\rm{deep}}}}}} = {1 \over {{{\tilde \tau }_{{\rm{deep}}}}}}\max \left( {0,{{p - 490\,{\rm{bar}}} \over {700\,{\rm{bar}} - {\rm{490}}\,{\rm{bar}}}}} \right).$ (4)

We used the same parameters as in Carone et al. (2020)³ and set τ_deep = 1 d. Carone et al. (2020) showed that these measures are important for maintaining numerical stability and to prevent nonphysical boundary effects forWASP-43b. These authors found that models with and without the deep drag yield small differences and only in the temperature and dynamics.

We did not use a convective adjustment scheme to model convection on a subgrid scale (e.g., Deitrick et al. 2020; Lee et al. 2021) because we are explicitly interested in the formation and stability of a deep adiabat by vertical heat-transport. Sainsbury-Martinez et al. (2019) argued that such a deep adiabat should form due to the absence of radiative heating and cooling at layers higher up in the atmosphere compared to layers where the atmosphere would become unstable to convection according to the Schwarzschild criterion.

The gas temperature T [K] of an ideal gas relates the density with the pressure and therefore ultimately determines the dynamics. MITgcm treats the temperature by means of the potential temperature, Θ [K], which is given by $Θ = T {(\frac{p_{0}}{p})}^{R / c_{p}},$ ${\rm{\Theta }}\,{\rm{ = }}\,T\,\,{\left( {{{{p_0}} \over p}} \right)^{{R \mathord{\left/ {\vphantom {R {{c_p}}}} \right. \kern-\nulldelimiterspace} {{c_p}}}}},$ (5)

where c_p [J kg⁻¹ K⁻¹ ] is the heat capacity at constant pressure, p₀ [Pa] is the pressure at the bottom of the computational domain, and R [J kg⁻¹ K⁻¹] is the specific gas constant. The potential temperature Θ is then forced using the thermodynamic energy equation (e.g., Showman et al. 2009): $\frac{dΘ}{d t} = \frac{Θ}{T} \frac{q}{c_{p}},$ ${{{\rm{d\Theta }}} \over {{\rm{d}}t}} = {{\rm{\Theta }} \over T}{q \over {{c_{\rm{p}}}}},$ (6)

where q [W kg⁻¹] is the heating rate. The heating rate q is given by $q = g \frac{\partial F^{net}}{\partial p},$ $q = g{{\partial {F^{{\rm{net}}}}} \over {\partial p}},$ (7)

where g [m s⁻²] is the gravity and F^net [W m⁻²] is the net radiative flux. We note that this approach towards radiative heating and cooling is not necessarily flux conserving, yielding the possibility of net cooling or heating of the planet.

2.2 Radiative Transfer

We used the Feautrier method (Feautrier 1964) in combination with lambda iteration in order to solve the radiative transfer equation. These schemes are frequently used in works on stellar atmospheres (e.g., Gustafsson et al. 1975, 2008), protoplanetary disks (Dullemond et al. 2002) and have been applied in some 1D models of exoplanetary atmospheres (Mollière et al. 2015,2017,2019,2020; Piette & Madhusudhan 2020). Our implementation of radiative transfer adopts the radiative transfer scheme used in the 1D planet atmosphere model petitRADTRANS (Mollière et al. 2019,2020).

2.2.1 Fluxes

The change of intensity I_v [W m⁻² sr⁻¹ Hz⁻¹] along a ray passing through a planetary atmosphere in the direction of the unit vector, n, may be described by the radiative transfer equation $n \cdot \nabla I_{v} = α_{v}^{tot} (S_{v} - I_{v}),$ ${\bf{n}} \cdot \nabla {I_v} = \alpha _v^{{\rm{tot}}}\left( {{S_v} - {I_v}} \right),$ (8)

where v [Hz] is the frequency, S_v [W m⁻² sr⁻¹ Hz⁻¹] is the source function and $α_{v}^{tot} [m^{- 1}]$ $\alpha _v^{{\rm{tot}}}\left[ {{{\rm{m}}^{ - 1}}} \right]$ is the inverse mean-free path of the light ray. We solve this equation for photons originating within the planetary atmosphere or scattered out of the incoming stellar ray. The extinction of the incoming stellar intensity is separately modeled using an exponential decay. This is not an approximation; such a separate treatment is allowed given to the linear nature of the equation of radiative transfer. The full solution is then obtained from adding the two intensities (planetary and scattered stellar photons along with extincted stellar intensity). The inverse mean free path $α_{v}^{tot}$ $\alpha _v^{{\rm{tot}}}$ is given by the sum of absorption and scattering inverse mean free paths and may be written as $α_{v}^{tot} = α_{v}^{abs} + α_{v}^{scat} .$ $\alpha _v^{{\rm{tot}}} = \alpha _v^{{\rm{abs}}} + \alpha _v^{{\rm{scat}}}.$ (9)

The mean intensity J_v [W m⁻² sr⁻¹ Hz⁻¹] is the zeroth order radiative momentum and is given by angle integration (e.g., Gaussian quadrature) of the intensity field. The zeroth and first order radiative moments in plane-parallel are given by $[J_{v}, H_{v}] = \frac{1}{2} \int_{- 1}^{1} [1, μ] I_{v} (z,μ) dμ$ $\left[ {{J_v},{H_v}} \right] = {1 \over 2}\int_{ - 1}^1 {\left[ {1,{\rm{\mu }}} \right]\,} {I_v}\left( {{\rm{z,\mu }}} \right){\rm{d\mu }}$ (10)

where μ = cos θ is the angle between the atmospheric normal and the ray and H_v [W m⁻² Hz⁻¹] is related to the flux (see below).

We use the source function in the coherent isotropic scattering approximation for the planetary and scattered stellar photon field, given by $S_{v} = ϵ_{v} B_{v} (T) + (1 - ϵ_{v}) (J_{v}^{pla} + J_{v}^{★}),$ ${S_v} = {_v}{B_v}\left( T \right) + \left( {1 - {_v}} \right)\left( {J_v^{{\rm{pla}}} + J_v^ \star } \right),$ (11)

where B_v(T) [W m⁻² sr⁻¹ Hz⁻¹] is the Planck function at temperature T [K], while $J_{v}^{pla}$ $J_v^{{\rm{pla}}}$ [W m⁻² sr⁻¹ Hz⁻¹] and $J_{v}^{⋆}$ $J_v^ \star$ [W m⁻² sr⁻¹ Hz⁻¹] are the mean intensity of the planetary atmosphere and the stellar attenuated radiation field respectively, and ∊_v is the photon destruction probability given by $ϵ_{v} = \frac{α_{v}^{abs}}{α_{v}^{tot}} .$ ${_v} = {{\alpha _v^{{\rm{abs}}}} \over {\alpha _v^{{\rm{tot}}}}}.$ (12)

The mean intensity of the stellar attenuated radiation field in a plane-parallel atmosphere $J_{v}^{⋆} [W m^{- 2} {sr}^{- 1} {Hz}^{- 1}]$ $J_v^ \star \left[ {{\rm{W}}\,{{\rm{m}}^{ - 2}}\,{\rm{s}}{{\rm{r}}^{ - 1}}{\rm{H}}{{\rm{z}}^{ - 1}}} \right]$ (Mollière et al. 2015) is given by $J_{v}^{★} = \frac{I_{v}^{★} (p = 0)}{4 π} \exp (- τ_{v} / μ_{*}),$ $J_v^ \star = {{I_v^ \star \left( {p = 0} \right)} \over {4\pi }}\exp \left( {{{ - {\tau _v}} \mathord{\left/ {\vphantom {{ - {\tau _v}} {{{\rm{\mu }}_ * }}}} \right. \kern-\nulldelimiterspace} {{{\rm{\mu }}_ * }}}} \right),$ (13)

where $J_{v}^{⋆} (p = 0) [W m^{- 2} {sr}^{- 1} {Hz}^{- 1}]$ $J_v^ \star \left( {p = 0} \right)\left[ {{\rm{W}}\,{{\rm{m}}^{ - 2}}\,{\rm{s}}{{\rm{r}}^{ - 1}}{\rm{H}}{{\rm{z}}^{ - 1}}} \right]$ is the stellar intensity field at the top of the planetary atmosphere. Here, expeRT/MITgcm uses the PHOENIX stellar model spectrum (Husser et al. 2013) that is part of petitRADTRANS to calculate the stellar attenuated radiation field. The angle between the atmospheric normal and the incoming stellar light μ_⋆ for a tidally locked exoplanet is given by $μ_{★} = \cos ϑ \cos φ,$ ${{\rm{\mu }}_ \star } = \cos \,\,\vartheta \,\cos \varphi ,$ (14)

where ϑ is the latitude and φ is the longitude. We note that we do not extend the plane-parallel assumption by including a geometric depth dependence in μ_⋆, which would improve the treatment of the poles (see Appendix A). The optical depth, τ_v, in a hydrostatic atmosphere that is parallel to the atmospheric normal is given by $τ_{v} = \int \frac{κ_{v}^{tot}}{g} d p,$ ${\tau _v} = \int {{{\kappa _v^{{\rm{tot}}}} \over g}{\rm{d}}p,}$ (15)

where $κ_{v}^{tot} [m^{2} {kg}^{- 1}]$ $\kappa _v^{{\rm{tot}}}\left[ {{{\rm{m}}^2}{\rm{k}}{{\rm{g}}^{ - 1}}} \right]$ is the cross-section per unit mass (opacity) and depends on the gas density, ρ [kg m⁻³], as: $κ_{v}^{tot} = \frac{α_{v}^{tot}}{ρ} .$ $\kappa _v^{{\rm{tot}}} = {{\alpha _v^{{\rm{tot}}}} \over \rho }.$ (16)

We can express the received bolometric stellar flux F^⋆ [W m⁻²] using Eqs. (10) and (13) as $F^{★} = - 4 π \int_{v} H_{v}^{★} d v =4 π μ_{★} \int_{v} J_{v}^{★} d v .$ ${F^ \star } = - 4\pi \int\limits_v {H_v^ \star {\rm{d}}v} {\rm{ = 4}}\pi {{\rm{\mu }}_ \star }\int\limits_v {J_v^ \star {\rm{d}}v.}$ (17)

We solved Eq. (8) using the well established Feautrier method (Feautrier 1964) with accelerated lambda iteration (ALI, Olson et al. 1986) in order to get the planetary flux. Once a converged solution (see Sect. 2.2.2) of the planetary intensity field has been found, we can calculate the emerging bolometric fluxes from the integration of the planetary radiation field F^pla [W m⁻²] using Eq. (10): $F^{pla} = \int_{v} F_{v}^{pla} d v = - 4 π \int_{v} H_{v}^{pla} d v .$ ${F^{{\rm{pla}}}} = \int\limits_v {F_v^{{\rm{pla}}}} {\rm{d}}v = - 4\pi \int\limits_v {H_v^{{\rm{pla}}}} {\rm{d}}v.$ (18)

The total bolometric flux F^net [W m⁻²] can now be calculated by $F^{net} = F^{pla} + F^{★} .$ ${F^{{\rm{net}}}} = {F^{{\rm{pla}}}} + {F^ \star }.$ (19)

Since Eq. (7) requires finite differencing of F^net, we followed the approach of Showman et al. (2009) and evaluated the fluxes on the vertically staggered cell interfaces. Lee et al. (2022) find that the use of a quadratic Bezier interpolation scheme for the interpolation of the temperature to the vertical interfaces greatly improves the accuracy and stability of the radiative transfer routine. Our tests confirm this finding and we therefore chose to also use this interpolation scheme.

We calculated Eq. (19) for every second grid column and linearly interpolate the fluxes in between (see e.g., Showman et al. 2009). The placement of the interpolated cells is shown in Fig. 1, where the bright cells indicate the interpolated cells. We note that the C32 geometry prohibits a strict interpolation of every second column. This is caused by the tiling of the horizontal domain into individual computational domains. It is therefore necessary to calculate the radiative transfer at the tile borders, as seen in Fig. 1. We test the interpolation scheme in Appendix B.3 and find that the error introduced by the interpolation is below 2%.

Fig. 1

Horizontal map displaying the number of iterations to converge the source function (Eq. (11)) during runtime (at t = 12000 d) for the HD 209458b simulation. The interpolated cells (ip) are staggered to each other. The C32 grid divides the horizontal domain into 24 individual computational domains, as seen by the line structure.

2.2.2 Convergence of the Source Function

There are two modes that expeRT/MITgcm can operate in: with scattering and without scattering. We set ∊_v = 1 if scattering is neglected, in which case Eq. (11) can be reduced to the LTE source function S_v = B_v(T). We iterate over the radiative transfer solution (so-called lambda iteration) if ∊_v < 1 (i.e., scattering is included), since S_v depends on I_v (see Eqs. (10) and (11)). We note that the Feautrier method has an exact solution in the isotropic scattering case. However, its solution would require a full matrix inversion (see e.g., Hubeny 2017), which is computationally expensive. Here we chose to solve the Feautrier equation by inverting a tridiagonal matrix instead, which requires iterating on the source function. For this, we use accelerated lambda iteration instead (ALI, Olson et al. 1986) to speed up the convergence of the source function. We choose to call the source function converged if the relative change of the source function is smaller than 2% of the reciprocal local lambda approximator (Olson et al. 1986; Rutten 2003). We performed numerical experiments, which confirmed that 2% yields a good tradeoff between accuracy and model performance, in agreement with Auer (1991); Rutten (2003). We use the knowledge of the source function of the previous radiative time-step as an initial guess for the current time-step. The number of iterations that are needed to converge the source function depends on the amplitude of structural changes in the temperature with time. During the first radiative time-step (see Table C.2) we find that 50–100 iterations are needed to converge the source function, while later time-steps reduce the amount of necessary iterations significantly to orders of one to three iterations per time-step (see Fig. 1).

Our numerical implementation of the ALI method in expeRT/MITgcm is an optimized version of the flux routine (with scattering) in petitRADTRANS (Mollière et al. 2020), which has been tested and benchmarked (Baudino et al. 2017). The major difference between the implementation in petitRADTRANS and expeRT/MITgcm is that we chose to not include Ng acceleration (Ng 1974) since the inclusion of Ng acceleration would only help within the first few time steps, whereas convergence of the source function is reached within a few iterations at later time steps (usually less than 2, see also Fig. 1).

The flexibility of our radiative transfer approach combined with its accuracy yields no drawback in terms of performance. When using the 24 core model setup, we report fast model runtimes of approximately 600 or 1000 simulated days per day of model runtime for a radiative time-step of 100 s or 300 s, respectively. We note that we take a different approach than e.g., Rauscher & Menou (2012); Lee et al. (2021), who developed a less complex radiative treatment to shorten the model run-time. Instead, our approach uses a full radiative transfer treatment, but with a different implementation than the one used in previous GCMs (e.g., Showman et al. 2009; Amundsen et al. 2016).

2.2.3 Opacities

Generally, expeRT/MITgcm operates on a precalculated grid of correlated-k opacities. The opacity in each correlated-k (e.g., Goody et al. 1989) wavelength bin is sorted in size and its distribution is sampled by a 16 point subgrid. Opacities are precalculated offline on a grid of 1000 logarithmically spaced temperature points between 100 K and 4000 K for every vertical layer. An overview of the used opacity sources can be found in Table 1. Opacities were obtained from the ExoMolOP database (Chubb et al. 2021) where available.

Opacities are precalculated and combined using a wrapper around the low-resolution mode of petitRADTRANS. Molecular abundances, the value for the specific gas constant R = 3590 J kg⁻¹ K⁻¹ and the specific heat capacity at constant pressure, c_p = 12 766 J kg⁻¹ K⁻¹, are inferred using the equilibrium chemistry package of petitRADTRANS. We used solar metallicity and solar C/O ratios, however, formation models for hot gas giants may hint towards sub-solar C/O ratios and super-solar metallicities (e.g., Mordasini et al. 2016; Schneider & Bitsch 2021a,b). We note that our approach to the chemical composition of the gas is computationally fast, but implies a fixed chemical composition during runtime. Additionally, this treatment neglects a physically correct treatment of cold traps, where certain species would condense out and not be available at altitudes above, which could prohibit a stratosphere caused by TiO and VO.

We used expeRT/MITgcm with three different wavelength grids (see Tables 2 and 3). The default resolution S0 uses N_bins = 6 wavelength bins from 0.26 µm to 300 µm. We used finer grids with 11 or 30 bins (S1 and S2 hereafter) to benchmark our model to the wavelength grids used in SPARC/MITgcm (Kataria et al. 2013). Those tests are shown in Appendix B.

Table 1

Opacity table.

Table 2

Low-resolution wavelength grid (S0).

Table 3

Wavelength resolutions.

2.3 Model Initialization

We followed the suggestion of Lee et al. (2021) and homogeneously initialize the upper layers (pressures below 1 bar) of our models with the analytic temperature profile of Parmentier et al. (2015). We calculated the substellar irradiation temperature T_irr of HD 209458b and WASP-43b using the stellar radius, R_⋆ [m], semi major axis, a_p [m], and stellar effective temperature, T_⋆ [K], as $T_{irr} = T_{★} \sqrt{\frac{R_{_{★}}}{2 a_{p}},}$ ${T_{{\rm{irr}}}} = {T_ \star }\sqrt {{{{R_{_ \star }}} \over {2{a_{\rm{p}}}}},}$ (20)

where we use an effective stellar temperature, T_⋆, and a stellar radius, R_⋆, of $T_{⋆}^{HD 209458} = 6092 K$ $T_ \star ^{{\rm{HD}}\,\,{\rm{209458}}} = 6092{\rm{K}}$ and $R_{⋆}^{HD 209458} = 1.203 R_{⊙}$ $R_ \star ^{{\rm{HD}}\,\,{\rm{209458}}} = 1.203{R_ \odot }$ (Boyajian et al. 2015) for HD 209458b and then $T_{⋆}^{WASP-43} = 4520 K$ $T_ \star ^{{\rm{WASP - 43}}} = 4520\,{\rm{K}}$ and $R_{⋆}^{WASP - 43} = 0.667 R_{⊙}$ $R_ \star ^{{\rm{WASP}} - {\rm{43}}} = 0.667\,{R_ \odot }$ for WASP-43b (Gillon et al. 2012). We assumed both planets to be tidally locked on a circular orbit. We set the rotation period of HD 209458b to $P_{rot}^{HD 209458b} = 3.47 d$ $P_{{\rm{rot}}}^{{\rm{HD}}\,{\rm{209458b}}} = 3.47\,{\rm{d}}$ and the orbital separation to $a_{p}^{HD 209458b} = 0.04747 AU$ $a_{\rm{p}}^{{\rm{HD}}\,\,{\rm{209458b}}} = 0.04747\,\,{\rm{AU}}$ (Southworth 2008,2010) and for WASP-43b to $P_{rot}^{WASP - 43 b} = 0.8135 d$ $P_{{\rm{rot}}}^{{\rm{WASP}} - 43{\rm{b}}} = 0.8135\,{\rm{d}}$ and $a_{p}^{WASP - 43 b} = 0.01526 AU$ $a_{\rm{p}}^{{\rm{WASP}} - 43{\rm{b}}} = 0.01526\,{\rm{AU}}$ (Gillon et al. 2012). Using Eq. (20), we estimated $T_{irr}^{HD 209458b} = 1479 K$ $T_{{\rm{irr}}}^{{\rm{HD}}\,\,{\rm{209458b}}} = 1479{\rm{K}}$ and $T_{irr}^{WASP - 43b} = 1441 K$ $T_{{\rm{irr}}}^{{\rm{WASP}} - {\rm{43b}}} = 1441\,{\rm{K}}$ for HD 209458b and WASP-43b, respectively. We used the fit from Thorngren et al. (2019) with these irradiation temperatures to estimate the interior temperature resulting from radiation of the planetary interior. This procedure yields interior temperatures of 575 K and 549 K, respectively. These irradiation and interior temperatures are then used to calculate the analytic temperature profile of Parmentier et al. (2015). We note that these interior temperatures were only used to initialize the temperature of our simulations, since we do not enforce radiative (-convective) equilibrium in the atmosphere of our GCM models during runtime.

The deep layers (pressures above 10 bar) are homogeneously set to a hot adiabatic temperature profile of the form $T (p > 10 bar) = Θ_{ad} \cdot {(\frac{p}{1 bar})}^{η},$ $T\left( {p > 10{\rm{bar}}} \right) = {{\rm{\Theta }}_{{\rm{ad}}}} \cdot {\left( {{p \over {1\,{\rm{bar}}}}} \right)^\eta },$ (21)

where we follow Sainsbury-Martinez et al. (2019) to guess high values of $Θ_{ad}^{HD 209458b} = 1800 K$ $\Theta _{{\rm{ad}}}^{{\rm{HD}}\,\,{\rm{209458b}}} = 1800\,\,{\rm{K}}$ and $Θ_{ad}^{WASP - 43b} = 1400 K$ $\Theta _{{\rm{ad}}}^{{\rm{WASP}} - {\rm{43b}}} = 1400\,\,{\rm{K}}$ as the temperature of the adiabat at 1 bar. The slope of the adiabat η is given by $η = \frac{R}{c_{p}} \approx 3.56.$ $\eta = {R \over {{c_{\rm{p}}}}} \approx 3.56.$ (22)

The intermediate pressure (1 bar ≤ p ≤ 10 bar) levels are linearly interpolated between the adiabat and the analytic model of Parmentier et al. (2015). This procedure allows for a hot initial state in the deep layers, as proposed in Sainsbury-Martinez et al. (2019). We assume a value of g^{HD 209458b} = 8.98 m s⁻² (e.g., Lee et al. 2021) and g^WASP−43b = 46.9 m s⁻² (Gillon et al. 2012) for the surface gravity, respectively. Differences in the model setup of WASP-43b compared to HD 209458b (see Table C.1) only occur in the surface gravity, rotation period and radiative time-step (see Table C.2).

We ran the model for 12 000 days with a dynamical timestep of Δt = 25 s⁴. Since the temperature in the upper layers might change significantly during the first 500 days, we recalculate the radiative fluxes in every fourth dynamical time-step, which we do by defining a radiative time-step, Δt_rad, which we set to Δt_rad = 100 s (see Table C.2). In the case of WASP-43b, we continue the simulation with the same radiative time-step, whereas we use a longer radiative time-step for the rest of the HD 209458b simulation of 300 s. A shorter radiative time-step is needed for the simulations of WASP-43b, since WASP-43b undergoes large structural changes throughout the model runtime, whereas HD 209458b reaches a steady state in the observable atmosphere within the first few hundred days (see Sect. 4).

Fig. 2

Maps of HD 209458b displaying the color coded temperature at different pressure levels (as indicated above the maps). The substellar point is located at (0°,0°).

3 State of the Upper Atmosphere

The number of planets modelled with fully coupled GCMs has increased a great deal since Showman et al. (2009) first introduced a fully coupled GCM for hot Jupiters. We chose to focus on HD 209458b and WASP-43b in this work. HD 209458b has been modeled by many of the relevant fully coupled hot Jupiter GCMs (Showman et al. 2009; Amundsen et al. 2016; Lee et al. 2021), making it a suitable target for comparisons between our model and previous works. Our non-gray 3D GCM WASP-43b simulations will be compared to simulations, at a similar level of complexity, of the same planet by Kataria et al. (2015). Comparisons with Baeyens et al. (2021), where a simpler radiative transfer scheme is used, are shown in Appendix D.

3.1 HD 209458b

We show the time-averaged temperature for different pressure layers in Fig. 2. The temperature has been time-averaged over the last 100 d. The temperature maps demonstrate the typical features of a tidally locked hot Jupiter climate: the strongly irradiated dayside transports heat to the nightside via a super-rotating jet (e.g., Showman & Polvani 2011), while shifting the hottest point eastwards. The upper layers (p < 1 × 10⁻³ bar) are strongly influenced by stellar irradiation, whereas advection is negligible (Fig. 2, top left panel). In contrast, the deeper layers (p > 1 × 10⁻² bar) are more efficient in redistributing heat by advection (Fig. 2, bottom left panel).

The zonal mean wind velocities are shown in Fig. 3. The superrotating jet of HD 209458b as seen in Kataria et al. (2016); Amundsen et al. (2016); Lee et al. (2021) is reproduced with expeRT/MITgcm. The magnitude of the jet (≈5000 m s⁻¹) agrees best with Kataria et al. (2016); Lee et al. (2021) and is lower than the value reported by Amundsen et al. (2016) of 7000 ms⁻¹.

The temperature maps of the HD 209458b simulations (Fig. 2) qualitatively match those of Lee et al. (2021, Fig. 5) showing that the general flow pattern as well as the day to night temperature difference are similar. Small asymmetric features in the temperature map highlight that the model has not yet equillibriated, even though the model was run for 12 000 days. We note that our simulations are initialized with a hot deep adiabat that is maintained throughout the runtime (see Sect. 4) while other simulations of HD 209458b (e.g., Amundsen et al. 2016; Showman et al. 2009) did not use such a hot interior.

The simulations of Showman et al. (2009); Lee et al. (2021) exhibit a temperature inversion due to the presence of TiO and VO around the substellar point, which results in an enhanced absorption of incoming stellar light. Our simulation of HD 209458b reproduces this temperature inversion, since we also include contributions of TiO and VO to the gas opacities that we use in our models (see Table 1). Small uncertainties in temperature and opacities may result in slightly hotter or cooler gas temperatures, which could result in the presence or absence of TiO and VO, since the dayside temperature of HD 209458b is very close to the TiO and VO condensation curves. This further highlights the importance of accurate interpolation schemes, such as the Bezier interpolation for the temperature (see Sect. 2.2.1), as well as the fine grained opacity grid with 1000 logarithmically spaced temperature points.

Fig. 3

Zonally averaged zonal wind velocities for HD 209458b. Positive values correspond to prograde flow. The white line depicts the change of sign.

3.2 WASP-43b

Similarly to Fig. 2 for HD 209458b, the temperature maps for WASP-43b are shown in Fig. 4. The results are well in line with the SPARC/MITgcm results for WASP-43b (Kataria et al. 2015) with TiO and VO. The day-to-night temperature differences as well as the flow patterns are comparable.

Unlike the case of Carone et al. (2020), our simulations do not exhibit a retrograde flow. As noted in Carone et al. (2020), the retrograde flow in their simulation was embedded in strong superrotation. Thus, both tendencies – the one of the prograde flow and that of the retrograde flow – were engaged in a “tug of war”. The main difference between the simulations of Carone et al. (2020) and this work is, however, that this work does not require artificial temperature forcing in deep layers. Our simulations therefore capture the full feedback between dynamics and radiation that was missing in Carone et al. (2020). This feedback seems to strengthen superrotation, apparently allowing for full superrotation to win the “tug of war” over the tendency for retrograde flow. Carone et al. (2020) found that the tendency toward retrograde flow is associated with a deep jet that can transport zonal momentum upwards. The transport of zonal momentum has been analyzed in several studies (e.g., Showman et al. 2015; Mendonça 2020; Wang & Wordsworth 2020; Carone et al. 2020), but it is still not clear what conditions would lead a deep equatorial jet toward retrograde flow. We will investigate the transport of zonal momentum in more detail in a follow-up publication and instead focus on the temperature evolution at greater depths in this work.

Carone et al. (2020) find that the depth of the equatorial jet of WASP-43b is likely linked to the rotation period. Their argumentation bases on two simulations of WASP-43b. Their nominal WASP-43b simulation with a rotation period of 0.81 days shows signs of a deep equatorial jet, which is not found in a simulation of WASP-43b with a rotation period of 3.5 days (as for HD 209458b). Our non-gray model of WASP-43b confirms this finding of Carone et al. (2020). WASP-43b indeed exhibits a deeper equatorial jet compared to HD 209458b, which can be seen in the zonal mean of the wind flow (see Fig. 5). We note that the simulations of Kataria et al. (2015) do not exhibit such a deep jet, which is likely caused by the short runtime of 50 d compared to the 12 000 d used here. Future works are warranted to unveil the physical reasons behind the depth of the equatorial jet.

4 Temperature Convergence of the Deep Layers

The GCM experiments of Sainsbury-Martinez et al. (2019) indicate that hot Jupiters might pump the energy input from irradiation into the deep atmosphere (p > 1 bar) by means of dynamic heat transport, leading to hot temperatures in their deep atmosphere. The equation of entropy conservation in steady state can be written as (e.g., Tremblin et al. 2017; Sainsbury-Martinez et al. 2019) $\frac{dΘ}{d t} - u \cdot \nabla Θ=0,$ ${{{\rm{d\Theta }}} \over {{\rm{d}}t}} - {\bf{u}} \cdot \nabla {\rm{\Theta = 0,}}$ (23)

where u is the velocity vector, ∇Θ is the spatial gradient of the potential temperature, Θ, and $\frac{dΘ}{d t}$ ${{{\rm{d}}\Theta } \over {{\rm{d}}t}}$ encapsulates dissipative process that cool or heat the atmosphere⁵. If we only consider irradiation, $\frac{dΘ}{d t}$ ${{{\rm{d}}\Theta } \over {{\rm{d}}t}}$ becomes very small in the deep layers of the atmosphere because most of the stellar light is absorbed in the upper atmosphere. However, this requires that either u or ∇Θ become very small. To fulfill the latter, the potential temperature must be close to constant and therefore follow an adiabat (T ∝ p^η).

The potential temperature of such an adiabat would therefore be mainly determined by the temperature of the convergence pressure, which is defined as the region where heating and cooling by radiation becomes very small. Furthermore, the convergence pressure lies higher up in the atmosphere compared to the pressure at which the atmosphere would be convective.

The weakness of radiative cooling and heating in the deep layers of the atmosphere requires long model runtimes. We can utilize the fast runtime of expeRT/MITgcm of 600–1000 simulated days per day of model runtime to study the temperature convergence in the deeper layers of the atmosphere. We investigate the temperature evolution of HD 209458b and WASP-43b in Fig. 6, where we show the temperature profile at the substellar point as a function of time. For HD 209458b, we show two additional simulations with different initial temperature profiles. The first one is identical to the nominal model, but uses Θ_ad = 2000 K for the initial deep adiabat, instead of Θ_ad = 1800 K for the nominal model (see Sect. 2.3); whereas the second simulation was initialized with a globally isothermal temperature of T = 2000 K. We performed a similar test for WASP-43b, where we ran a model of WASP-43b that has been initialized, with Θ_ad = 900 K, instead of the nominal Θ_ad = 1400 K. We chose to stop the extra HD 209458b and WASP-43b simulations with an adiabatic interior at t = 3000 d, since the only purpose of these simulations is to gain information about the dependence of the temperature convergence behavior on initial conditions. However, it is important to note that the convergence time needed to cool down the deep atmosphere from a hot state to a colder state is shorter than the time needed to heat the deep atmosphere from a cold state to a hotter state (Sainsbury-Martinez et al. 2019).

We find that for WASP-43b, the temperature in the deep layers is steadily evolving during the full 12 000 d simulation. The rate of temperature change drops from approximately 1.5 K d⁻¹ to approximately 0.05 K d^–1 at the end of the simulation, while slowly converging. This process of cooling results in an almost isothermal temperature profile in the photosphere of WASP-43b, which is subsequently changing into a more pronounced temperature inversion once the deep cooling becomes less efficient (see Fig. 6, right panel). The decaying rate of cooling in the deep layers forecasts a final cool adiabat of Θ_ad ≈ 700 K.

The temperature evolution of the initially cool WASP-43b simulation overlaps with the temperature evolution of the nominal simulation. However, the later states of the nominal WASP-43b including the thermal inversion matches the initially cold simulation fairly well. We therefore find that the evolution of the temperature of the initially cold WASP-43b simulation confirms the independence of the final state of the deeper atmosphere of WASP-43b on initial conditions.

HD 209458b, on the other hand, appears to maintain its initial temperature profile in the deep layers (see Fig. 6, left panel). The temperature constantly decreases with a rate varying around 0.1 ± 0.1 K d^–1. Furthermore, we see that the deep temperature evolution is similar for the hot initial temperature profile, yielding very low values for the temperature decay in time. A perfect initial guess by chance is therefore impossible. We therefore conclude that we could not reach temperature convergence in the deep layers of HD 209458b within 12 000 days. This outcome, while “negative”, is of significance for other studies that have already published temperature maps of HD 209458b.

The temperature in the isothermally initialized HD 209458b simulation does not evolve towards a hotter state in the deep layers. It is not possible to predict whether the final state of the temperature in the deep layers of these simulations will indeed be adiabatic. We can thus not confirm Sainsbury-Martinez et al. (2019) in that the deep layers of this simulation heat up towards such a state. Additionally, the temperature exhibits a shaky pattern that is likely linked to the deeper parts of the atmosphere. This pattern arises at the convergence pressure and could be explained by upward propagating gravity waves; a process that has also been postulated for the atmospheres of brown dwarfs (e.g., Freytag et al. 2010; Showman & Kaspi 2013). We therefore postulate that this pattern is linked to the strength of the dynamic processes in the deeper parts of the atmosphere. It remains an open question if this pattern disappears when convergence is reached or if it could be forced to a stably stratified temperature profile with the use of a convective adjustment scheme (e.g., Deitrick et al. 2020; Lee et al. 2021), which we do not use within this work. We plan to examine these gravity waves and the effect of a convective adjustment scheme in detail in an upcoming work.

In the following, we investigate the effect of surface gravity and planetary rotation on the convergence behavior of the temperature in the deep atmosphere. Radiative timescales scale reciprocally with the surface gravity (see e.g., Showman et al. 2020), leading to longer radiative timescales in HD 209458b at similar pressure levels compared to WASP-43b, even though both planets receive a similar amount of energy from their host star. The pressure level where the radiative heating and cooling becomes irrelevant is therefore different for both planets, where irradiation is penetrating into deeper pressure layers in the atmosphere of WASP-43b compared to HD 209458b. This will then influence the pressure level and temperature at which a deep adiabat could couple to the radiatively dominated part of the atmosphere (according to the explanations of Tremblin et al. 2017; Sainsbury-Martinez et al. 2019). Hence, we think that a colder interior of WASP-43b in our simulations is a natural consequence of the assumption of the value of the surface gravity.

We performed an additional simulation of a WASP-43b-like planet with a rotation period equal to the rotation period of HD 209458b (P_rot = 3.47 d), instead of the nominal short rotation period of P_rot = 0.8135 d. The goal of this additional simulation is to consider the dependence of the temperature evolution on rotation. The zonal wind velocity in the case of slow rotation is significantly different to our nominal model of WASP-43b, where we find that unlike the nominal simulation (Fig. 5), the slowly rotating WASP-43b simulation does not exhibit a deep equatorial jet (Fig. 7). This finding agrees well with the work of Carone et al. (2020), who predicted that WASP-43b exhibits a deep superrotating jet, which is linked to the fast rotation of WASP-43b.

Recently, Sainsbury-Martinez et al. (2021) proposed that the efficiency of downward heat transport of energy depends on the presence or absence of rotational dynamics over advective dynamics. If rotational dynamics dominate over advective dynamics, which is the case for a deep equatorial jet, we would expect a lower efficiency for the downward heat transport of energy. We therefore compared the temperature in the model of WASP-43b with a deep jet (nominal) to the model without a deep jet (P_rot = 3.5 d). We show the resulting substellar point temperature profiles in Fig. 8, where we find that both simulations exhibit a similar temperature evolution in the deep atmosphere. The independence of the deep temperature profile on the rotation period leads us to the conclusion that the fast temperature convergence of WASP-43b can be almost exclusively explained by the continuously strong radiative heating and cooling in the deep atmosphere caused by the high value of surface gravity.

The requirement of zero radiative heating and cooling (see Eq. (23)) for the theory of downward transport of energy (Sainsbury-Martinez et al. 2019) can only be monitored with a climate model that includes radiative feedback. Using our model, we find that radiative heating and cooling can still be important at greater depths, bringing on the question of how the requirements for the formation of a hot deep adiabat might be matched. However, we do note that our simulations demonstrate that our non-gray GCM succeeds in maintaining an adiabatic deep atmosphere, without the need for a convective adjustment scheme. We propose that a further detailed analysis is needed to settle the question of whether a downward transport of energy is possible.

Fig. 4

Maps of WASP-43b displaying the color coded temperature at different pressure levels (like Fig. 2).

Fig. 5

Zonally averaged zonal wind velocities (like Fig. 3) for WASP-43b.

Fig. 6

Evolution of the temperature profile at the substellar point for HD 209458b (left) and WASP-43b (right). Left panel: includes the evolution of two additional HD 209458b simulations with a colder and a hotter initial temperature profile. Right panel: includes one additional WASP-43b simulation with a cooler initial temperature profile. The additional simulations initialized with an adiabat ran for 3000 days, whereas the cool HD 209458b ran for 5000 days. The temperature time averaging and output interval is set to 100 days, therefore not resolving the fast semi-convergence of the temperature in the photosphere during the first 100 days. The horizontal orange dashed line illustrates the approximate altitude at which radiative heating or cooling becomes very small.

Fig. 7

Zonally averaged wind velocities (like Fig. 5) for a model of WASP-43b with a rotation period equal to the rotation period of HD 209458b.

Fig. 8

Evolution of the temperature profile at the substellar point (like Fig. 6). Three time steps (100 d, 5000 d and 12 000 d) of the nominal model of WASP-43b are compared to a model of WASP-43b with a rotation period equal to the rotation period of HD 209458b. The color scale is the same as in Fig. 6.

5 Synthetic Spectra

We post-processed our GCM results using petitRADTRANS to obtain emission and transmission spectra. We followed Baeyens et al. (2021) and calculated the transmission spectrum using the average temperature in latitude ±20° around the west and east terminators (φ = −90° and φ = 90°). Such an isothermal treatment of the terminator region is possible since temperature variations in the terminator region are very small for pressure layers below 1 × 10⁻² bar, which are most relevant to the transmission spectrum. The transmission spectrum is then calculated using the radius and gravity values of Table C.1 for a reference pressure of 0.01 bar. As in Baeyens et al. (2021), we determined the final transit spectrum by quadratically averaging the east and west terminators. For WASP-43b, we show the Kreidberg et al. (2014) HST/WFC3 data around the ≈ 1.4 µm H₂O feature, whereas for HD 209458b, we show the data of Crossfield et al. (2012); Sing et al. (2016); Evans et al. (2015). We calculate the mean of the observed values of the wavelength dependent radius as well as of the calculated transmission spectra for each observed spectral window to determine the offset in radius between observations and our model. We then corrected the observed radius by subtracting this offset. The offset correction is necessary due to differences in the normalization of the radius (see e.g., Lee et al. 2021).

For the emission spectrum, we chose to regrid our GCM output onto a low resolution rectangular grid with a 15° longitude-latitude resolution (288 horizontal grid cells). We then calculated the angle dependent planetary intensity at the top of each vertical column using a derivative of the flux routine of petitRADTRANS for 20 Gaussian µ quadrature points. We then integrated the intensity field over the solid angle to obtain the flux in the direction of the observer (see e.g., Seager 2010; Amundsen et al. 2016). The general procedure is identical to the one used in Carone et al. (2020) and only differs in its numerical implementation, where we use petitRADTRANS instead of petitCODE to obtain the intensity⁶. We tested the impact of the low-resolution spatial grid on the resulting emission spectrum by testing an even lower spatial resolution of 30° (72 grid cells) and found that only minimal differences occur, validating the 288 grid cells approach outlined above. We chose to plot the data of Stevenson et al. (2014) for WASP-43b and the data of Line et al. (2016) for HD 209458b. We note that we do not correct the observed emission spectrum for systematic errors, but instead plot the published values.

We investigate whether the difference in the magnitude of the temperature inversion caused by the ongoing convergence of the deeper atmosphere of our WASP-43b simulation would impact predictions for observables. The resulting transmission spectra for HD 209458b in Fig. 9 and WASP-43b in Fig. 10 show good agreement with the observed values. However, we find that the pressure region between around 1 bar contributes most to the spectrum for wavelengths between 0.4 and 2 µm, whereas the pressure range between 0.1 and 0.001 bar dominates the infrared spectrum. This would explain, why we only see a (very small) difference between the transmission spectra for WASP-43b at 500 d and 12 000 d in the infrared but not in the optical (see Fig. 9).

In examining the dayside emission spectra, we find that the shape of the water line feature of WASP-43b and HD 209458b does not seem to agree between the observed water feature and the model spectra. The reason for the difference in the shape of the water feature is found in the inclusion of TiO and VO opacities, which cause a temperature inversion (Showman et al. 2009). A similar trend can be found in Lee et al. (2021, Fig. 7, right panel), where the authors show the spectral differences of a semi-gray model without temperature inversion compared to more advanced models that include the temperature inversion. A model without TiO and VO in the upper atmosphere would therefore clearly fit the observed spectra much better. Further, the lack of clouds and the assumption of solar metallicity in our model affects the emitted flux and also effect the strength of superrotation (e.g., Kataria et al. 2015; Drummond et al. 2018a; Parmentier et al. 2021).

Moreover, we find that the 500 d and 12 000 d emission spectra of the WASP-43b model are clearly impacted by the state of the deep atmosphere (see Fig. 6). The flux differences in the emission spectrum thus highlight the effect of deep temperature convergence on the observable atmosphere. These differences are not seen for HD 209458b because the observable atmosphere seems to be unaffected by the state of the deep atmosphere (see Fig. 6).

Fig. 9

Synthetic transmission spectrum (top) and dayside emission spectrum (bottom) for HD 209458b calculated from simulation results after 12 000 days, compared to the observations. The thermal inversion, caused by the TiO and VO opacities, overcasts the observed water feature at ≈ 1.4 µm.

Fig. 10

Synthetic transmission spectra (top) and dayside emission spectra (bottom) for WASP-43b calculated from simulation results after 500 and 12 000 days, compared to the observations. The thermal inversion, caused by the TiO and VO opacities, overcasts the observed water feature at ≈ 1.4 µm.

6 Discussion

In this work, we employ expeRT/MITgcm, a new 3D GCM with non-gray radiative transfer specifically built to investigate the processes at work in greater depths in more detail. The first application of this model to WASP-43b and HD209458b shows interesting results in flow and temperature at such depths and may thus shed light onto the unsolved question of hot Jupiter inflation.

6.1 Versatile Fast 3D GCM with non-gray Radiative Transfer

There are currently several 3D exoplanet atmosphere models available, which vary in their complexity of radiative transfer (Showman et al. 2009; Rauscher & Menou 2012; Dobbs-Dixon & Agol 2013; Amundsen et al. 2016; Mendonça et al. 2018b; Lee et al. 2021) and even in the basic implementation of the Navier Stokes equation to treat fluid dynamics, ranging from hydrostatic primitive equations (e.g., Showman et al. 2009) to full (e.g., Mayne et al. 2014; Mendonça et al. 2018b; Deitrick et al. 2020). In addition, several 3D climate models include additional physical processes, which are thought to be important for the understanding of exoplanet atmospheres: clouds (Lee et al. 2016; Lines et al. 2019; Parmentier et al. 2021; Roman et al. 2021), photochemical hazes (Steinrueck et al. 2021), disequilibrium chemistry (Agúndez et al. 2014a; Drummond et al. 2018b; Mendonça et al. 2018a; Baeyens et al. 2021), and magnetic fields (Rogers & Showman 2014). Cloud modeling in 3D climate models ranges from detailed kinetic modeling (e.g., Lee et al. 2016) to simple but efficient parametrization (e.g., Roman et al. 2021; Christie et al. 2021).

In addition to cloud modeling, disequilibrium chemistry, and magnetic field interaction, there is a need for GCMs that are capable of investigating processes at greater depths such as deep wind flow and deep energy transport. In this work, we introduce a new non-gray GCM formalism that is numerically efficient and thus ideally suited to tackle the challenging processes at greater depths. The main reason for these advances is the use of the Feautrier method in combination with accelerated lambda iteration to solve the radiative transfer equation (Sects. 2.2.1 and 2.2.2). We find that these methods allow us to run non-gray GCM models for 1000 simulation days within just a single day of run time and using 24 cores. Applying the code to WASP-43b and HD 209458b shows that the upper atmosphere agrees well with the results of other non-gray GCMs (Sect. 3) for these planets (Kataria et al. 2015; Showman et al. 2009; Amundsen et al. 2016; Lee et al. 2021). The differences with Newtonian cooling as used in Carone et al. (2020) are shown in Appendix. D. The computational efficiency of expeRT/MITgcm allows extending the work of Sainsbury-Martinez et al. (2019) to investigate temperature convergence for very long simulation times at deeper layers without artificial external forcing mechanisms (Sect. 4).

6.2 Importance of the Temperatures in the Deep Atmosphere

The temperature in layers below 10 bar can influence the observable chemistry in several ways: for instance, iron and magnesium could condense in these layers and consequently would not be present in the gas phase of the upper atmosphere, where they would be observable via high-resolution spectroscopy (Sing et al. 2019, Fig. 13). Other possibilities include the quenching of methane via disequilibrium chemistry (Agúndez et al. 2014b; Baeyens et al. 2021) and the thickness of clouds also depends on the temperature at deeper layers (Fortney et al. 2020). Silicate cloud formation can remove or inject oxygen-bearing molecules from the gas phase (Helling et al. 2019,2021). Therefore, deep cloud condensation might affect local C/O ratios in cloud forming regions even more than thinner clouds. At the same time, the condensation of important cloud particles such as magnesium at greater depths may affect the composition of clouds higher up, as pointed out by Helling et al. (2021).

In this work, we focused on energy transport and thermal evolution at greater depths, because it is still an open question whether energy can be injected via 3D circulation from the irradiated photosphere into deeper layers in hot to ultra-hot Jupiters that could at least partly explain hot Jupiter inflation (Mayne et al. 2014,2017; Tremblin et al. 2017). For WASP-43b, we find that higher surface gravities compared to HD 209458b already naturally lead to a deeper convergence pressure. Thus, setting the surface gravity and assuming long-enough run times drives the simulation at greater depths naturally to a colder state. The surface gravity of WASP-43b, however, is higher because it is not inflated. In other words: the state of inflation sets the surface gravity, which, in turn, strongly affects the convergence behavior towards the final state. This dependence of the convergence behavior on the state of inflation is obviously not optimal, since it creates a “chicken-and-the-egg problem” when investigating energy transport at greater depths. However, we stress that our numerical efficient non-gray GCM has the potential to explore this problem in more detail, where we will also take additional dissipative processes at greater depths into account in the future.

6.3 Effect of Deep Dynamics on the Observable Atmosphere

Only a few GCM studies of hot Jupiters have been conducted thus far with the aim to explore dynamical effects at greater depths owing to the computational challenges that need to be resolved in the deeper atmosphere. Mayne et al. (2017) explored the possibility of deep wind flow when horizontal temperature differences at greater depths are imposed. Carone et al. (2020) showed that fast rotators like WASP-43b (orbital periods less than 1.5 d) naturally evolve equatorial jets that can extend very deep into the interior. In the Solar System, the fast-rotating gas giants are also shown to have deep jets (Kaspi et al. 2018; Iess et al. 2019). In this work, the deep equatorial jets in WASP-43b and shallow equatorial jets in HD 209458b are reproduced. Thus, we stress, as does Carone et al. (2020), that some extrasolar planets like the JWST ERS target WASP-43b (Venot et al. 2020; Bean et al. 2018) can be dynamically active at greater depth greater than 200 bar, which is the typical lower boundary assumed for hot Jupiter GCMs (e.g., Showman et al. 2009; Amundsen et al. 2016; Deitrick et al. 2020; Parmentier et al. 2021; Lee et al. 2021). More importantly, we find in this work that resolving dynamics at greater depths is important to allow the deep atmosphere of WASP-43b to self-consistently evolve towards cool adiabat (Θ_ad ≈ 700 K).

It is not only the energy transport at greater depths that is influenced by deep dynamics. Previous works indicate that the presence of a deep flow may impact the observable wind and temperature structure of hot Jupiters. Mayne et al. (2017) showed that the deep wind flow can greatly diminish superrotation in the observable atmosphere. Wang & Wordsworth (2020) found that the simulated wind flow structure in their mini-Neptune simulation can abruptly change from two jets to one jet, which could be observable via changes in the hotspot location. This change happens without forced conversion after 50 000 simulation days. Furthermore, Carone et al. (2020) find that in their simulations, superrotation can be completely disrupted at the dayside, leading to a local retrograde flow. In this work, we do not find this retrograde flow, but instead unperturbed superrotation. This can partly be attributed to the fact that our models do not rely on artificial temperature forcing in the deep atmosphere, but instead, they can capture the feedback between dynamics and irradiation that tends to strengthen superrotation. Carone et al. (2020) also noted that superrotation is never completely removed, and that there are apparently two tendencies at play in the upper atmosphere: one leading to superrotation and one leading to other flows.

Parmentier et al. (2021) noted that the strength of superrotation in the observable atmosphere at the dayside can be diminished by nightside cloud coverage. Beltz et al. (2022); Hindle et al. (2021) indicated that magnetic field interaction can completely suppress superrotation. Observationally, it was found that CoRoT-2b appears to have a retrograde flow at the dayside (Dang et al. 2018) and HAT-P-7b may even have a dynamically shifting hot spot (Armstrong et al. 2016). Thus, future works that includes large-scale flow at greater depths is warranted in order to investigate the interplay among more complex physics and its effect on the observable atmosphere.

It has been further hypothesized that shocks, mechanical dissipation of strong winds, and turbulent vertical mixing at greater depths may also inject energy into the interior and contribute to inflation (Li & Goodman 2010; Heng 2012; Perna et al. 2012; Fromang et al. 2016; Menou 2019). It is, however, unclear whether shear flow instabilities at greater depths (Menou & Rauscher 2009; Rauscher & Menou 2010; Liu & Showman 2013; Carone et al. 2020) may affect the vertical heat flow. Tackling these processes requires higher resolution non-hydrostatic 3D atmosphere simulations (Menou 2019; Fortney et al. 2021).

We note that the question of the importance of dynamical processes at greater depths in atmospheres is neither new nor solely constrained to hot Jupiters. There are clear parallels between the dynamics and the heat transport at greater depths in brown dwarfs (Sainsbury-Martinez et al. 2021). Considerable work has also been done in understanding such processes both for the Earth’s climate and for stellar atmospheres. Such examples are clumpy dust formation triggered by shock waves via large scale convective flows at greater depths in AGB stars (Fleischer et al. 1992; Woitke & Niccolini 2005; Höfner & Freytag 2019) and the need to resolve downdrafts reaching deep into the interior (Freytag et al. 2017, 2019) or “convective over-shooting” at the boundary of convective and stably stratified regions in stars (e.g., Hanasoge et al. 2015; Bressan et al. 1981). The issue of “convective overshooting” seems to resemble the boundary problem between a convective interior and a stably stratified irradiated atmosphere in hot Jupiters. Similar processes are also relevant in Earth GCMs: the importance of dry and moist convection in the Earth’s atmosphere and surface friction was identified early on (Manabe et al. 1965; Peixóto & Oort 1984). The importance of convection and surface friction was also identified in GCMs of tidally locked exo-Earths (Carone et al. 2016; Koll & Abbot 2016; Sergeev et al. 2020) and for brown dwarfs (Tremblin et al. 2019). Thus, it comes as no surprise that processes at greater depths matter in hot Jupiters as well. More importantly, deep dynamical processes impact the large-scale circulation higher up and can potentially yield important insights into the evolution of exoplanets.

7 Summary and conclusions

Here, we introduce and explain expeRT/MITgcm, which extends the deep wind framework exorad/MITgcm of Carone et al. (2020) with full gaseous radiative transfer, including scattering. We demonstrate that the Feautrier method with accelerated lambda iteration and correlated-k is a fast and accurate approach to model the non-gray radiative transfer while still achieving model runtimes of hundreds of simulation days in a single runtime day.

We simulated the two hot Jupiters HD 209458b and WASP-43b which are almost equally irradiated, but which have different orbital periods $(P_{rot}^{HD 209458b} = 3.47 d and P_{rot}^{WASP - 43 b} = 0.8135 d)$ $\left( {P_{{\rm{rot}}}^{{\rm{HD}}\,\,\,{\rm{209458b}}} = 3.47\,\,{\rm{d}}\,\,{\rm{and}}\,\,\,P_{{\rm{rot}}}^{{\rm{WASP}} - 43{\rm{b}}} = 0.8135\,{\rm{d}}} \right)$ and surface gravities $(g^{HD 209458b} = 8.98 m s^{- 2} and g^{WASP - 43 b} = 46.9 m s^{- 2})$ $\left( {{g^{{\rm{HD}}\,\,{\rm{209458b}}}} = 8.98\,\,{\rm{m}}\,{{\rm{s}}^{ - 2}}\,\,{\rm{and}}\,\,{g^{{\rm{WASP}} - 43{\rm{b}}}} = 46.9\,\,{\rm{m}}\,{{\rm{s}}^{ - 2}}} \right)$ .

We find that our results for the observable atmosphere of HD 209458b and WASP-43b agree well with the findings of previous works. Our WASP-43b models were able to confirm a deep equatorial jet as well as a dependency of the temperature in the observable atmosphere on the state of the deep atmosphere as predicted by Carone et al. (2020). However, we did not find the retrograde equatorial flow. Future work is needed to investigate the reason for these differences.

We investigated the deep atmosphere (pressures between 1 and 700 bar) of WASP-43b and HD 209458b and found a difference in the temperature evolution. WASP-43b cools down to a cold temperature in the deep atmosphere, whereas HD 209458b successfully maintains its initial temperature profile within the runtime. This difference in the temperature evolution can be linked to the different values of the surface gravity, where we find that due to the high surface gravity of our WASP-43b simulations, radiative heating and cooling are very important at greater depths. We looked at the dependency of the temperature evolution on the rotation period and found that the difference in temperature evolution is not significantly enhanced by the presence or absence of a deep equatorial jet.

Last but not least, we caution that longer convergence time scales of at least 10 000 simulation days are needed to properly resolve the observable atmosphere temperature structure of non-inflated exoplanets such as WASP-43b. The observable atmosphere of inflated hot Jupiters such as HD 209458b seems to be decoupled from the deep atmosphere and, therefore, it can be expected to be resolved after 1000 days of simulation time.

More work is warranted to identify whether and how energy transport from the irradiated atmosphere into the interior contributes to inflation. The computationally efficient non-gray GCM expeRT/MITgcm that we introduce here is an ideal tool to self-consistently tackle deep wind flow and energy transport.

Acknowledgements

A.D.S., L.D., U.G.J, S.K. and C.H. acknowledge funding from the European Union H2020-MSCA-ITN-2019 under Grant no. 860470 (CHAMELEON). L.C. acknowledges support from the DFG Priority Programme SP1833 Grant CA 1795/3. L.D. acknowledges support from the KU LEUVEN interdisciplinary IDN grant IDN/19/028 and the FWO research grant G086217N. U.G.J acknowledges funding from the Novo Nordisk Foundation Interdisciplinary Synergy Program grant no. NNF19OC0057374. P.M. acknowledges support from the European Research Council under the European Union’s Horizon 2020 research and innovation program under grant agreement No. 832428-Origins. R.B. acknowledges support from a PhD fellowship of the Research Foundation - Flanders (FWO). The computational resources and services used in this work were provided by the VSC (Flemish Supercomputer Center), funded by the Research Foundation Flanders (FWO) and the Flemish Government – department EWI. We also thank an anonymous referee for their very useful comments that helped a lot in the interpretation of the obtained results.

Appendix A Importance of Sphericity in Stellar Attenuation

The incident angle for the stellar light onto the planetary atmosphere is accounted by the factor μ_⋆ In Eq. 14, μ_⋆ is derived assuming plane-parallel geometry, which is a broadly used approximation (see e.g., Showman et al. 2009; Amundsen et al. 2016; Lee et al. 2021). Because this paper is meant to benchmark the radiative transfer implementation by comparing it with previous work, we used plane-parallel atmospheres and did not account for the sphericity of the model atmospheres.

While the incident angle in the substellar region is independent of depth and agrees well with the plane-parallel assumption, substantial differences are expected at the polar regions, where the effective solar path length deviates significantly from the solar path length in plane-parallel, ultimately altering the angle of incidence. The expected deviation in μ_⋆ between a spherical and a plane-parallel approximation increases with pressure (and, therefore, optical depth) and is potentially more important at greater depths. On the other hand, we note that the stellar component of the flux F^⋆ decreases with optical depth.

A correct treatment of sphericity would involve rewriting Eqs. 13 and 15, which would lead to a double integral over the solid angle. Since this is not computationally feasible, some authors (Li & Shibata 2006; Mendonça et al. 2018a) who account for sphericity in the environment of a GCM have expanded the plane-parallel assumption with the correction of the effective path-length. The resulting incident angle μ_⋆,eff (z) would then be dependent on the height z and the planetary radius R_p and is in its simplest approximation given by (Li & Shibata 2006): $μ_{★,eff} (z) = \sqrt{1 - {(\frac{R_{p}}{R_{p} + z})}^{2} (1 - μ_{★}^{2}) .}$ ${{\rm{\mu }}_{ \star {\rm{,eff}}}}\left( z \right) = \sqrt {1 - {{\left( {{{{R_p}} \over {{R_p} + z}}} \right)}^2}\left( {1 - {\rm{\mu }}_ \star ^2} \right).}$ (A.1)

In conclusion, the plane-parallel assumption is suitable if only the uppermost part of the atmosphere is probed or if the error of the optical depth within the polar regions causes only negligible effects. If the atmosphere is to be probed deeper, an additional treatment of the sphericity is vital (e.g., Li & Shibata 2006). We plan to incorporate this aspect in an upcoming work.

Appendix B Verification of the Radiative Transfer

It has been shown in previous works (Showman et al. 2009; Amundsen et al. 2016; Lee et al. 2021) as well as in this work that the inclusion of a full gaseous radiative transfer in GCMs is vital to obtain a full understanding of the atmospheric dynamics of hot Jupiters. However, a correct treatment of radiative transfer in the framework of a GCM is computationally very expensive and simplifications need to be well understood and justified. The correlated-k approximation and the assumption of equilibrium chemistry both lead to significant simplifications.

We tested the validity of radiative transfer in our model in a three-step process. We first show a comparison of the fluxes and heating rates obtained with petitRADTRANS and expeRT/MITgcm in Sect. B.1 and then our testing of the wavelength binning, done by comparing the resulting temperature profile in Sect. B.2. We finally show our testing of the influence of the horizontal interpolation of fluxes on the temperature profile in Sect. B.3.

B.1 Verification of Fluxes and Heating Rates

We follow the suggestion of Amundsen et al. (2014) and implement their dayside test for a cloud-free atmosphere in expeRT/MITgcm. We show the fluxes and heating rates in Fig. B.1, where we compare the results obtained with expeRT/MITgcm in different wavelength resolutions with scattering turned on and off to results obtained with petitRADTRANS, which has been benchmarked in Baudino et al. (2017).

We find that our results agree well with petitRADTRANS, which uses a much higher wavelength resolution of approximately 1000 times more correlated-k wavelength bins compared to the nominal resolution of our models (see Table 3). We note that the employed method of random sampling to combine opacities of individual species is prone to errors in the order of a few percent. On top of that, there is a general loss of accuracy when using low-resolution opacities (e.g., Leconte 2021), which is also in the order of a few percent.

The differences in the bolometric stellar flux $(F_{v}^{⋆})$ $\left( {F_v^ \star } \right)$ can be explained by the numerical accuracy of the integration, which depends on the amount of frequency bins used to integrate the spectral flux (see Eq. 17). Reducing this error would require calculations of the stellar bolometric flux for a high wavelength resolution and a normalization of the stellar spectral flux to the high-resolution bolometric flux. The current low-resolution implementation of petitRADTRANS does not include such a normalization, so any future implementation of expeRT/MITgcm might improve on that.

Overall, there seems to be a general good agreement between Fig. B.1 and the simulations shown in Amundsen et al. (2014). We note that Amundsen et al. (2014) used the Kurucz stellar spectrum⁷, while we used the PHOENIX stellar spectrum (Husser et al. 2013) that is part of petitRADTRANS. Differences might also occur due to the different opacities, where we used the same species as mentioned in Amundsen et al. (2014), but from updated opacity sources (see Table 1). Furthermore, we used the equilibrium chemistry of Mollière et al. (2019), while Amundsen et al. (2014) used the analytic abundances of Burrows & Sharp (1999).

Fig. B.1

Test 3 of Amundsen et al. (2014) for a mixed dayside hot Jupiter atmosphere, where we compare fluxes and heating rates for different wavelength resolutions (see Table 3) in expeRT/MITgcm (SO, S1, S2) to petitRADTRANS (prt). The bottom panels display the residuals of the above panels. From left to right, we show the incoming stellar flux, the emitting planetary flux, the gradient of the net flux, and the resulting change of potential temperature.

B.2 S0 vs S1

We perform an additional simulation of HD 209458b using the S1 resolution in order to benchmark the S0 resolution used in the model (see Table 3). The wavelength resolution S1 is introduced in Kataria et al. (2013). This additional simulation is initialized with the results of the S0 simulation, shown in Sect. 3.1 at 10000 days. We then run the model for 2000 days and compare the resulting temperature profile in Fig. B.2.

We find that differences between the two wavelength resolutions are on the order of a few percent in the observable part of the atmosphere, while the relative error decreases with pressure. The hotter temperature of the upper layers in the S0 simulation can partly be accounted for by the stronger incoming flux, as discussed above. We note that the smaller error in the deeper layers of the atmosphere demonstrates that the deep layers are not radiatively driven, and can only be driven by other mechanisms such as dynamics at greater depths.

Fig. B.2

Additional test demonstrating the performance of SO vs S1 by looking at the final temperature differences (at 12000 days) for several longitudes. The S1 simulation is initialized at 10000 days using the state of the nominal SO simulation.

B.3 Performance of the Interpolation Scheme

We performed an additional short term simulation (run for 100 days) without the horizontal interpolation of fluxes, where we solved the radiative transfer equation in every grid cell instead of every second grid cell (see Fig. 1). These models are identically to the cold HD 209458b model presented in Sect. 4, where we use a globally isothermal temperature of 2000 K to initialize the model. We show the resulting time averaged temperature after 100 days of the model with and without interpolation in Fig. B.3 and find that both agree within 2%. We can therefore conclude that the small wavelength resolution introduces a larger error compared to the error introduced by the horizontal interpolation of the fluxes.

Fig. B.3

Additional test demonstrating the performance of the flux interpolation scheme (see Fig. 1 and Sect. 2.2.1) by looking at an isothermally (T=2000 K) initialized HD 209458b simulation with and without the horizontal interpolation of fluxes.

Appendix C Parameters

We show the parameters that we employ in our models of HD 209458b and WASP-43b in Table C.1. We update the radiative fluxes of WASP-43b more often than we do for HD 209458b, since our models of WASP-43b keep evolving in temperature even after hundreds of days, whereas our model of HD 209458b shows no sign of temperature evolution after a few hundred days. The corresponding numerical flux update rates (radiative time-steps) are shown in Table C.2.

Table C.1

Parameters

Table C.2

Radiative time-steps

Appendix D Simplified Forcing vs Full Radiative Transfer

Our simulations of HD 209458b and WASP-43b with fully coupled radiative feedback are compared to similar simulations from a grid of 3D GCM models taken from Baeyens et al. (2021) with different parameters. The grid models are forced using the Newtonian cooling formalism and use the same model as in Carone et al. (2020). Specifically, to match the climate properties of HD 209458b a model with T_eff = 1400 K, g = 10 m s⁻² and host star type G5 is selected (Figs D.1 and D.2). For the WASP-43b-like case, we compare our simulation to a model with T_eff = 1400 K, g = 100 m s⁻² and host star type K5 (Figs D.3 and D.4). Additional information about these grid models can be found in Baeyens et al. (2021). Because the Newtonian-cooled grid models and the fully coupled models in this work make use of the same dynamical core (MITgcm) and similar radiative transfer module (petitRADTRANS), their comparison highlights the important role of radiative feedback in hot Jupiter climates. We note that the grid models do not use TiO and VO opacities, while the fully coupled models do include TiO and VO, leading to a temperature inversion in the fully coupled models, which is not to be expected in the Newtonian-cooled models.

In general, the Newtonian-cooled grid models show qualitative agreement with the simulations presented in this work, but quantitative differences are present. For the case of HD 209458b, the grid model (Fig. D.1) appears to be cooler than the fully coupled model (Fig. 2), in particular regarding the nightside of the planet, owing to the colder effective temperature of the grid model. The equatorial jet is not as fast in the grid model as in the fully coupled simulation. In the former, it reaches wind speeds up to 5.0 km s⁻¹, whereas in the latter, wind speeds up to 7.2 km s⁻¹ are attained (see Figs D.2 and 3). Finally, the jet persists up to the lowest pressures in the fully coupled simulation, whereas in the Newtonian-cooled model, the low-pressure regime is characterized by a thermally direct day-to-night circulation. The transition between a shallow and a deep jet seems to be for rotation periods between two and three days (e.g, compare first two columns of Fig. 3 from Baeyens et al. 2021), and because the grid model is not a perfect match to HD 209458b, it falls on the different side of this transition.

The disappearance of cold nightsides with a reduced day-night temperature contrast and a thermally direct wind flow at high altitudes, when full radiative feedback is taken into account, can also be seen in Amundsen et al. (2016). A more efficient day-night heat redistribution, when radiative feedback is taken into account, was also noted by Showman et al. (2009). In the models of Amundsen et al. (2016), however, the wind speed of the zonal jet stream is slightly slower in the fully coupled GCM than in the Newtonian-cooled GCM. In contrast, we find that the full radiative feedback yields a zonal jet stream that is slightly faster than the Newtonian-cooled case. Hence, further detailed comparisons between the two forcing mechanisms are necessary to establish their impact on the jet stream speed.

Finally, for the case of WASP-43b, the same quantitative changes discussed above can be observed when a Newtonian-cooled and fully coupled model are compared. Additionally, a climate transition is observed, in which the Newtonian-cooled model yields equatorial retrograde flow (Fig. D.3), but the fully coupled model does not (Fig. 4). A possible explanation for the disappearance of retrograde flow in the latter case could be that retrograde and superrotating equatorial flow are competing tendencies (Carone et al. 2020). The increase in jet strength in the fully coupled model, as established in the HD 209458b comparison in the previous paragraph, would then become the dominant tendency in WASP-43b. Given that the fully coupled model should be the more realistic of the two forcing mechanisms, future investigations will be needed to investigate whether retrograde equatorial wind flow can also be elicited in models with radiative feedback.

Fig. D.1

Maps of a model of HD 209458b with Newtonian cooling displaying the color coded temperature at different pressure levels (like Fig. 2).

Fig. D.2

Zonally averaged zonal wind velocities (like Fig. 3) for the model of HD 209458b with Newtonian cooling.

Fig. D.3

Maps of WASP-43b with Newtonian cooling displaying the color coded temperature at different pressure levels (like Fig. 4).

Fig. D.4

Zonally averaged zonal wind velocities (like Fig. 5) for the model of WASP-43b with Newtonian cooling.

References

Adcroft, A., Campin, J.-M., Hill, C., & Marshall, J. 2004, Mon. Weather Rev., 132, 2845 [NASA ADS] [CrossRef] [Google Scholar]
Agúndez, M., Parmentier, V., Venot, O., Hersant, F., & Selsis, F. 2014a, A&A, 564, A73 [Google Scholar]
Agúndez, M., Venot, O., Selsis, F., & Iro, N. 2014b, ApJ, 781, 68 [Google Scholar]
Allard, N. F., Spiegelman, F., Leininger, T., & Molliere, P. 2019, VizieR Online Data Catalog: J/A+A/628/A120 [Google Scholar]
Amundsen, D. S., Baraffe, I., Tremblin, P., et al. 2014, A&A, 564, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Amundsen, D. S., Mayne, N. J., Baraffe, I., et al. 2016, A&A, 595, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Armstrong, D. J., de Mooij, E., Barstow, J., et al. 2016, Nat. Astron., 1, 0004 [Google Scholar]
Auer, L. 1991, ASI Ser. C, 341, 9 [NASA ADS] [Google Scholar]
Azzam, A. A. A., Tennyson, J., Yurchenko, S. N., & Naumenko, O. V. 2016, MNRAS, 460, 4063 [NASA ADS] [CrossRef] [Google Scholar]
Baeyens, R., Decin, L., Carone, L., et al. 2021, MNRAS, 505, 5603 [NASA ADS] [CrossRef] [Google Scholar]
Barber, R. J., Strange, J. K., Hill, C., et al. 2014, MNRAS, 437, 1828 [NASA ADS] [CrossRef] [Google Scholar]
Batygin, K., & Stevenson, D. J. 2010, ApJ, 714, L238 [NASA ADS] [CrossRef] [Google Scholar]
Baudino, J.-L., Mollière, P., Venot, O., et al. 2017, ApJ, 850, 150 [Google Scholar]
Bean, J. L., Stevenson, K. B., Batalha, N. M., et al. 2018, PASP, 130, 114402 [NASA ADS] [CrossRef] [Google Scholar]
Beltz, H., Rauscher, E., Roman, M. T., & Guilliat, A. 2022, AJ, 163, 35 [NASA ADS] [CrossRef] [Google Scholar]
Borysow, A. 2002, A&A, 390, 779 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Borysow, A., & Frommhold, L. 1989, ApJ, 341, 549 [NASA ADS] [CrossRef] [Google Scholar]
Borysow, J., Frommhold, L., & Birnbaum, G. 1988, ApJ, 326, 509 [NASA ADS] [CrossRef] [Google Scholar]
Borysow, A., Frommhold, L., & Moraldi, M. 1989, ApJ, 336, 495 [NASA ADS] [CrossRef] [Google Scholar]
Borysow, A., Jorgensen, U. G., & Fu, Y. 2001, J. Quant. Spectr. Rad. Transf., 68, 235 [NASA ADS] [CrossRef] [Google Scholar]
Boyajian, T., von Braun, K., Feiden, G. A., et al. 2015, MNRAS, 447, 846 [Google Scholar]
Bressan, A. G., Chiosi, C., & Bertelli, G. 1981, A&A, 102, 25 [Google Scholar]
Burrows, A., & Sharp, C. M. 1999, ApJ, 512, 843 [NASA ADS] [CrossRef] [Google Scholar]
Carone, L., Keppens, R., & Decin, L. 2016, MNRAS, 463, 3114 [NASA ADS] [CrossRef] [Google Scholar]
Carone, L., Baeyens, R., Mollière, P., et al. 2020, MNRAS, 496, 3582 [NASA ADS] [CrossRef] [Google Scholar]
Chan, Y. M., & Dalgarno, A. 1965, Proc. Phys. Soc., 85, 227 [NASA ADS] [CrossRef] [Google Scholar]
Christie, D. A., Mayne, N. J., Lines, S., et al. 2021, MNRAS, 506, 4500 [CrossRef] [Google Scholar]
Chubb, K. L., Rocchetto, M., Yurchenko, S. N., et al. 2021, A&A, 646, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Coles, P. A., Yurchenko, S. N., & Tennyson, J. 2019, MNRAS, 490, 4638 [CrossRef] [Google Scholar]
Crossfield, I. J. M., Knutson, H., Fortney, J., et al. 2012, ApJ, 752, 81 [Google Scholar]
Dalgarno, A., & Williams, D. A. 1962, ApJ, 136, 690 [NASA ADS] [CrossRef] [Google Scholar]
Dang, L., Cowan, N. B., Schwartz, J. C., et al. 2018, Nat. Astron., 2, 220 [NASA ADS] [CrossRef] [Google Scholar]
Deitrick, R., Mendonça, J. M., Schroffenegger, U., et al. 2020, ApJS, 248, 30 [NASA ADS] [CrossRef] [Google Scholar]
Dobbs-Dixon, I., & Agol, E. 2013, MNRAS, 435, 3159 [NASA ADS] [CrossRef] [Google Scholar]
Drummond, B., Mayne, N. J., Baraffe, I., et al. 2018a, A&A, 612, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Drummond, B., Mayne, N. J., Manners, J., et al. 2018b, ApJ, 855, L31 [NASA ADS] [CrossRef] [Google Scholar]
Dullemond, C. P., van Zadelhoff, G. J., & Natta, A. 2002, A&A, 389, 464 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Evans, T. M., Aigrain, S., Gibson, N., et al. 2015, MNRAS, 451, 680 [Google Scholar]
Feautrier, P. 1964, C. R. Acad. Sci., 258, 3189 [NASA ADS] [Google Scholar]
Fleischer, A. J., Gauger, A., & Sedlmayr, E. 1992, A&A, 266, 321 [NASA ADS] [Google Scholar]
Fortney, J. J., Visscher, C., Marley, M. S., et al. 2020, AJ, 160, 288 [NASA ADS] [CrossRef] [Google Scholar]
Fortney, J. J., Dawson, R. I., & Komacek, T. D. 2021, J. Geophys. Res. Planets, 126, e06629 [NASA ADS] [CrossRef] [Google Scholar]
Freytag, B., Allard, F., Ludwig, H. G., Homeier, D., & Steffen, M. 2010, A&A, 513, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Freytag, B., Liljegren, S., & Höfner, S. 2017, A&A, 600, A137 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Freytag, B., Höfner, S., & Liljegren, S. 2019, IAU Symp., 343, 9 [NASA ADS] [Google Scholar]
Fromang, S., Leconte, J., & Heng, K. 2016, A&A, 591, A144 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Gillon, M., Triaud, A. H. M. J., Fortney, J. J., et al. 2012, A&A, 542, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Goody, R., West, R., Chen, L., & Crisp, D. 1989, J. Quant. Spectr. Rad. Transf., 42, 539 [NASA ADS] [CrossRef] [Google Scholar]
Gray, D. F. 2008, The Observation and Analysis of Stellar Photospheres (Cambridge: Cambridge University Press) [Google Scholar]
Gustafsson, B., Bell, R. A., Eriksson, K., & Nordlund, A. 1975, A&A, 42, 407 [NASA ADS] [Google Scholar]
Gustafsson, B., Edvardsson, B., Eriksson, K., et al. 2008, A&A, 486, 951 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Hanasoge, S., Miesch, M. S., Roth, M., et al. 2015, Space Sci. Rev., 196, 79 [NASA ADS] [CrossRef] [Google Scholar]
Hartman, J. D., Bakos, G. Á., Bhatti, W., et al. 2016, AJ, 152, 182 [NASA ADS] [CrossRef] [Google Scholar]
Helling, C., Iro, N., Corrales, L., et al. 2019, A&A, 631, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Helling, C., Lewis, D., Samra, D., et al. 2021, A&A, 649, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Heng, K. 2012, ApJ, 761, L1 [NASA ADS] [CrossRef] [Google Scholar]
Heng, K., Menou, K., & Phillipps, P. J. 2011, MNRAS, 413, 2380 [NASA ADS] [CrossRef] [Google Scholar]
Hindle, A. W., Bushby, P. J., & Rogers, T. M. 2021, ApJ, 922, 176 [CrossRef] [Google Scholar]
Höfner, S., & Freytag, B. 2019, A&A, 623, A158 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Hubeny, I. 2017, MNRAS, 469, 841 [NASA ADS] [CrossRef] [Google Scholar]
Husser, T. O., Wendevon Berg, S., Dreizler, S., et al. 2013, A&A, 553, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Iess, L., Militzer, B., Kaspi, Y., et al. 2019, Science, 364, aat2965 [Google Scholar]
Kaspi, Y., Galanti, E., Hubbard, W. B., et al. 2018, Nature, 555, 223 [Google Scholar]
Kataria, T., Showman, A. P., Lewis, N. K., et al. 2013, ApJ, 767, 76 [NASA ADS] [CrossRef] [Google Scholar]
Kataria, T., Showman, A. P., Fortney, J. J., et al. 2015, ApJ, 801, 86 [NASA ADS] [CrossRef] [Google Scholar]
Kataria, T., Sing, D. K., Lewis, N. K., et al. 2016, ApJ, 821, 9 [NASA ADS] [CrossRef] [Google Scholar]
Koll, D. D. B., & Abbot, D. S. 2016, ApJ, 825, 99 [NASA ADS] [CrossRef] [Google Scholar]
Kreidberg, L., Bean, J. L., Désert, J.-M., et al. 2014, ApJ, 793, L27 [Google Scholar]
Leconte, J. 2021, A&A, 645, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Lee, E., Dobbs-Dixon, I., Helling, C., Bognar, K., & Woitke, P. 2016, A&A, 594, A48 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Lee, E. K. H., Lothringer, J. D., Casewell, S. L., et al. 2022, MNRAS, submitted [arXiv:2203.09854] [Google Scholar]
Lee, E. K. H., Parmentier, V., Hammond, M., et al. 2021, MNRAS, 506, 2695 [NASA ADS] [CrossRef] [Google Scholar]
Li, J., & Goodman, J. 2010, ApJ, 725, 1146 [NASA ADS] [CrossRef] [Google Scholar]
Li, J., & Shibata, K. 2006, J. Atmos. Sci., 63, 1365 [NASA ADS] [CrossRef] [Google Scholar]
Li, G., Gordon, I. E., Rothman, L. S., et al. 2015, ApJS, 216, 15 [NASA ADS] [CrossRef] [Google Scholar]
Line, M. R., Stevenson, K. B., Bean, J., et al. 2016, AJ, 152, 203 [NASA ADS] [CrossRef] [Google Scholar]
Lines, S., Mayne, N. J., Manners, J., et al. 2019, MNRAS, 488, 1332 [NASA ADS] [CrossRef] [Google Scholar]
Liu, B., & Showman, A. P. 2013, ApJ, 770, 42 [NASA ADS] [CrossRef] [Google Scholar]
Manabe, S., Smagorinsky, J., & Strickler, R. F. 1965, Mon. Weather Rev., 93, 769 [NASA ADS] [CrossRef] [Google Scholar]
Mayne, N. J., Baraffe, I., Acreman, D. M., et al. 2014, A&A, 561, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Mayne, N. J., Debras, F., Baraffe, I., et al. 2017, A&A, 604, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Mayne, N. J., Drummond, B., Debras, F., et al. 2019, ApJ, 871, 56 [Google Scholar]
McKemmish, L. K., Yurchenko, S. N., & Tennyson, J. 2016, MNRAS, 463, 771 [NASA ADS] [CrossRef] [Google Scholar]
McKemmish, L. K., Masseron, T., Hoeijmakers, H. J., et al. 2019, MNRAS, 488, 2836 [Google Scholar]
Mendonça, J. M. 2020, MNRAS, 491, 1456 [CrossRef] [Google Scholar]
Mendonça, J. M., Malik, M., Demory, B.-O., & Heng, K. 2018a, AJ, 155, 150 [Google Scholar]
Mendonça, J. M., Tsai, S.-M., Malik, M., Grimm, S. L., & Heng, K. 2018b, ApJ, 869, 107 [CrossRef] [Google Scholar]
Menou, K. 2019, MNRAS, 485, L98 [Google Scholar]
Menou, K., & Rauscher, E. 2009, ApJ, 700, 887 [NASA ADS] [CrossRef] [Google Scholar]
Mollière, P., van Boekel, R., Dullemond, C., Henning, T., & Mordasini, C. 2015, ApJ, 813, 47 [Google Scholar]
Mollière, P., van Boekel, R., Bouwman, J., et al. 2017, A&A, 600, A10 [Google Scholar]
Mollière, P., Wardenier, J. P., van Boekel, R., et al. 2019, A&A, 627, A67 [Google Scholar]
Mollière, P., Stolker, T., Lacour, S., et al. 2020, A&A, 640, A131 [Google Scholar]
Mordasini, C., van Boekel, R., Mollière, P., Henning, T., & Benneke, B. 2016, ApJ, 832, 41 [Google Scholar]
Ng, K. C. 1974, J. Chem. Phys., 61, 2680 [NASA ADS] [CrossRef] [Google Scholar]
Olson, G. L., Auer, L. H., & Buchler, J. R. 1986, J. Quant. Spectr. Rad. Transf., 35, 431 [Google Scholar]
Parmentier, V., Guillot, T., Fortney, J. J., & Marley, M. S. 2015, A&A, 574, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Parmentier, V., Showman, A. P., & Fortney, J. J. 2021, MNRAS, 501, 78 [Google Scholar]
Peixóto, J. P., & Oort, A. H. 1984, Rev. Mod. Phys., 56, 365 [CrossRef] [Google Scholar]
Perna, R., Heng, K., & Pont, F. 2012, ApJ, 751, 59 [NASA ADS] [CrossRef] [Google Scholar]
Piette, A. A. A., & Madhusudhan, N. 2020, ApJ, 904, 154 [NASA ADS] [CrossRef] [Google Scholar]
Piskunov, N. E., Kupka, F., Ryabchikova, T. A., Weiss, W. W., & Jeffery, C. S. 1995, A&AS, 112, 525 [Google Scholar]
Polyansky, O. L., Kyuberis, A. A., Zobov, N. F., et al. 2018, MNRAS, 480, 2597 [NASA ADS] [CrossRef] [Google Scholar]
Rauscher, E., & Menou, K. 2010, ApJ, 714, 1334 [NASA ADS] [CrossRef] [Google Scholar]
Rauscher, E., & Menou, K. 2012, ApJ, 750, 96 [NASA ADS] [CrossRef] [Google Scholar]
Richard, C., Gordon, I. E., Rothman, L. S., et al. 2012, J. Quant. Spectr. Rad. Transf., 113, 1276 [Google Scholar]
Rogers, T. M., & Showman, A. P. 2014, ApJ, 782, L4 [NASA ADS] [CrossRef] [Google Scholar]
Roman, M. T., Kempton, E. M. R., Rauscher, E., et al. 2021, ApJ, 908, 101 [CrossRef] [Google Scholar]
Rutten, R. J. 2003, Radiative Transfer in Stellar Atmospheres (Utrecht: Lecture Notes Utrecht University) [Google Scholar]
Sainsbury-Martinez, F., Casewell, S. L., Lothringer, J. D., Phillips, M. W., & Tremblin, P. 2021, A&A, 656, A128 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Sainsbury-Martinez, F., Wang, P., Fromang, S., et al. 2019, A&A, 632, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Schneider, A. D., & Bitsch, B. 2021a, A&A, 654, A71 [Google Scholar]
Schneider, A. D., & Bitsch, B. 2021b, A&A, 654, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Seager, S. 2010, Exoplanet Atmospheres: Physical Processes (Princeton: Princeton University Press) [Google Scholar]
Sergeev, D. E., Lambert, F. H., Mayne, N. J., et al. 2020, ApJ, 894, 84 [NASA ADS] [CrossRef] [Google Scholar]
Showman, A. P., & Kaspi, Y. 2013, ApJ, 776, 85 [NASA ADS] [CrossRef] [Google Scholar]
Showman, A. P., & Polvani, L. M. 2011, ApJ, 738, 71 [Google Scholar]
Showman, A. P., Fortney, J. J., Lian, Y., et al. 2009, ApJ, 699, 564 [NASA ADS] [CrossRef] [Google Scholar]
Showman, A. P., Lewis, N. K., & Fortney, J. J. 2015, ApJ, 801, 95 [NASA ADS] [CrossRef] [Google Scholar]
Showman, A. P., Tan, X., & Parmentier, V. 2020, Space Sci. Rev., 216, 139 [Google Scholar]
Sing, D. K., Fortney, J. J., Nikolov, N., et al. 2016, Nature, 529, 59 [Google Scholar]
Sing, D. K., Lavvas, P., Ballester, G. E., et al. 2019, AJ, 158, 91 [Google Scholar]
Skinner, J. W., & Cho, J. Y. K. 2021, MNRAS, 504, 5172 [NASA ADS] [CrossRef] [Google Scholar]
Sousa-Silva, C., Al-Refaie, A. F., Tennyson, J., & Yurchenko, S. N. 2015, MNRAS, 446, 2337 [Google Scholar]
Southworth, J. 2008, MNRAS, 386, 1644 [NASA ADS] [CrossRef] [Google Scholar]
Southworth, J. 2010, MNRAS, 408, 1689 [Google Scholar]
Steinrueck, M. E., Showman, A. P., Lavvas, P., et al. 2021, MNRAS, 504, 2783 [NASA ADS] [CrossRef] [Google Scholar]
Stevenson, K. B., Désert, J.-M., Line, M. R., et al. 2014, Science, 346, 838 [Google Scholar]
Thorngren, D. P., & Fortney, J. J. 2018, AJ, 155, 214 [Google Scholar]
Thorngren, D., Gao, P., & Fortney, J. J. 2019, ApJ, 884, L6 [Google Scholar]
Tremblin, P., Chabrier, G., Mayne, N. J., et al. 2017, ApJ, 841, 30 [NASA ADS] [CrossRef] [Google Scholar]
Tremblin, P., Padioleau, T., Phillips, M. W., et al. 2019, ApJ, 876, 144 [Google Scholar]
Venot, O., Parmentier, V., Blecic, J., et al. 2020, ApJ, 890, 176 [NASA ADS] [CrossRef] [Google Scholar]
Wang, H., & Wordsworth, R. 2020, ApJ, 891, 7 [NASA ADS] [CrossRef] [Google Scholar]
Wende, S., Reiners, A., Seifahrt, A., & Bernath, P. F. 2010, A&A, 523, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Woitke, P., & Niccolini, G. 2005, A&A, 433, 1101 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Yurchenko, S. N., Amundsen, D. S., Tennyson, J., & Waldmann, I. P. 2017, A&A, 605, A95 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Yurchenko, S. N., Mellor, T. M., Freedman, R. S., & Tennyson, J. 2020, MNRAS, 496, 5282 [NASA ADS] [CrossRef] [Google Scholar]
Zhang, X. 2020, Res. Astron. Astrophys., 20, 099 [NASA ADS] [CrossRef] [Google Scholar]

¹

expeRT/MITgcm is an extension of the exorad/MITgcm scheme established in Carone et al. (2020).

²

C32 is comparable to a resolution of 128 × 64 in longitude and latitude.

³

There is a typo in the drag description in Carone et al. (2020): τ_fric and τ_bottom,fric should be $\frac{1}{τ_{fric}}$ ${1 \over {{\tau _{{\rm{fric}}}}}}$ and $\frac{1}{τ_{bottom,fric}}$ ${1 \over {{\tau _{{\rm{bottom,fric}}}}}}$ respectively.

⁴

Day refers to 24 h or 1 earth day.

⁵

$\frac{dΘ}{d t}$ ${{{\rm{d}}\Theta } \over {{\rm{d}}t}}$ is given by Eq. (6) in this work, but other processes such as Ohmic dissipation could also play a role.

⁶

The python package used to obtain and integrate the emission spectrum has been made publicly available and can be found at https://prt-phasecurve.readthedocs.io/

⁷

see http://kurucz.harvard.edu/stars/hd289458/

All Tables

Table 1

Opacity table.

In the text

Table 2

Low-resolution wavelength grid (S0).

In the text

Table 3

Wavelength resolutions.

Parameters

Radiative time-steps

All Figures

	Fig. 1 Horizontal map displaying the number of iterations to converge the source function (Eq. (11)) during runtime (at t = 12000 d) for the HD 209458b simulation. The interpolated cells (ip) are staggered to each other. The C32 grid divides the horizontal domain into 24 individual computational domains, as seen by the line structure.
In the text

	Fig. 2 Maps of HD 209458b displaying the color coded temperature at different pressure levels (as indicated above the maps). The substellar point is located at (0°,0°).
In the text

	Fig. 3 Zonally averaged zonal wind velocities for HD 209458b. Positive values correspond to prograde flow. The white line depicts the change of sign.
In the text

	Fig. 4 Maps of WASP-43b displaying the color coded temperature at different pressure levels (like Fig. 2).
In the text

	Fig. 5 Zonally averaged zonal wind velocities (like Fig. 3) for WASP-43b.
In the text

Fig. 6

Evolution of the temperature profile at the substellar point for HD 209458b (left) and WASP-43b (right). Left panel: includes the evolution of two additional HD 209458b simulations with a colder and a hotter initial temperature profile. Right panel: includes one additional WASP-43b simulation with a cooler initial temperature profile. The additional simulations initialized with an adiabat ran for 3000 days, whereas the cool HD 209458b ran for 5000 days. The temperature time averaging and output interval is set to 100 days, therefore not resolving the fast semi-convergence of the temperature in the photosphere during the first 100 days. The horizontal orange dashed line illustrates the approximate altitude at which radiative heating or cooling becomes very small.

In the text

	Fig. 7 Zonally averaged wind velocities (like Fig. 5) for a model of WASP-43b with a rotation period equal to the rotation period of HD 209458b.
In the text

	Fig. 8 Evolution of the temperature profile at the substellar point (like Fig. 6). Three time steps (100 d, 5000 d and 12 000 d) of the nominal model of WASP-43b are compared to a model of WASP-43b with a rotation period equal to the rotation period of HD 209458b. The color scale is the same as in Fig. 6.
In the text

	Fig. 9 Synthetic transmission spectrum (top) and dayside emission spectrum (bottom) for HD 209458b calculated from simulation results after 12 000 days, compared to the observations. The thermal inversion, caused by the TiO and VO opacities, overcasts the observed water feature at ≈ 1.4 µm.
In the text

	Fig. 10 Synthetic transmission spectra (top) and dayside emission spectra (bottom) for WASP-43b calculated from simulation results after 500 and 12 000 days, compared to the observations. The thermal inversion, caused by the TiO and VO opacities, overcasts the observed water feature at ≈ 1.4 µm.
In the text

Fig. B.1

Test 3 of Amundsen et al. (2014) for a mixed dayside hot Jupiter atmosphere, where we compare fluxes and heating rates for different wavelength resolutions (see Table 3) in expeRT/MITgcm (SO, S1, S2) to petitRADTRANS (prt). The bottom panels display the residuals of the above panels. From left to right, we show the incoming stellar flux, the emitting planetary flux, the gradient of the net flux, and the resulting change of potential temperature.

In the text

	Fig. B.2 Additional test demonstrating the performance of SO vs S1 by looking at the final temperature differences (at 12000 days) for several longitudes. The S1 simulation is initialized at 10000 days using the state of the nominal SO simulation.
In the text

	Fig. B.3 Additional test demonstrating the performance of the flux interpolation scheme (see Fig. 1 and Sect. 2.2.1) by looking at an isothermally (T=2000 K) initialized HD 209458b simulation with and without the horizontal interpolation of fluxes.
In the text

	Fig. D.1 Maps of a model of HD 209458b with Newtonian cooling displaying the color coded temperature at different pressure levels (like Fig. 2).
In the text

	Fig. D.2 Zonally averaged zonal wind velocities (like Fig. 3) for the model of HD 209458b with Newtonian cooling.
In the text

	Fig. D.3 Maps of WASP-43b with Newtonian cooling displaying the color coded temperature at different pressure levels (like Fig. 4).
In the text

	Fig. D.4 Zonally averaged zonal wind velocities (like Fig. 5) for the model of WASP-43b with Newtonian cooling.
In the text

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

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

Initial download of the metrics may take a while.

[1] Adcroft, A., Campin, J.-M., Hill, C., & Marshall, J. 2004, Mon. Weather Rev., 132, 2845 [NASA ADS] [CrossRef] [Google Scholar]

[2] Agúndez, M., Parmentier, V., Venot, O., Hersant, F., & Selsis, F. 2014a, A&A, 564, A73 [Google Scholar]

[3] Agúndez, M., Venot, O., Selsis, F., & Iro, N. 2014b, ApJ, 781, 68 [Google Scholar]

[4] Allard, N. F., Spiegelman, F., Leininger, T., & Molliere, P. 2019, VizieR Online Data Catalog: J/A+A/628/A120 [Google Scholar]

[5] Amundsen, D. S., Baraffe, I., Tremblin, P., et al. 2014, A&A, 564, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[6] Amundsen, D. S., Mayne, N. J., Baraffe, I., et al. 2016, A&A, 595, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[7] Armstrong, D. J., de Mooij, E., Barstow, J., et al. 2016, Nat. Astron., 1, 0004 [Google Scholar]

[8] Auer, L. 1991, ASI Ser. C, 341, 9 [NASA ADS] [Google Scholar]

[9] Azzam, A. A. A., Tennyson, J., Yurchenko, S. N., & Naumenko, O. V. 2016, MNRAS, 460, 4063 [NASA ADS] [CrossRef] [Google Scholar]

[10] Baeyens, R., Decin, L., Carone, L., et al. 2021, MNRAS, 505, 5603 [NASA ADS] [CrossRef] [Google Scholar]

[11] Barber, R. J., Strange, J. K., Hill, C., et al. 2014, MNRAS, 437, 1828 [NASA ADS] [CrossRef] [Google Scholar]

[12] Batygin, K., & Stevenson, D. J. 2010, ApJ, 714, L238 [NASA ADS] [CrossRef] [Google Scholar]

[13] Baudino, J.-L., Mollière, P., Venot, O., et al. 2017, ApJ, 850, 150 [Google Scholar]

[14] Bean, J. L., Stevenson, K. B., Batalha, N. M., et al. 2018, PASP, 130, 114402 [NASA ADS] [CrossRef] [Google Scholar]

[15] Beltz, H., Rauscher, E., Roman, M. T., & Guilliat, A. 2022, AJ, 163, 35 [NASA ADS] [CrossRef] [Google Scholar]

[16] Borysow, A. 2002, A&A, 390, 779 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[17] Borysow, A., & Frommhold, L. 1989, ApJ, 341, 549 [NASA ADS] [CrossRef] [Google Scholar]

[18] Borysow, J., Frommhold, L., & Birnbaum, G. 1988, ApJ, 326, 509 [NASA ADS] [CrossRef] [Google Scholar]

[19] Borysow, A., Frommhold, L., & Moraldi, M. 1989, ApJ, 336, 495 [NASA ADS] [CrossRef] [Google Scholar]

[20] Borysow, A., Jorgensen, U. G., & Fu, Y. 2001, J. Quant. Spectr. Rad. Transf., 68, 235 [NASA ADS] [CrossRef] [Google Scholar]

[21] Boyajian, T., von Braun, K., Feiden, G. A., et al. 2015, MNRAS, 447, 846 [Google Scholar]

[22] Bressan, A. G., Chiosi, C., & Bertelli, G. 1981, A&A, 102, 25 [Google Scholar]

[23] Burrows, A., & Sharp, C. M. 1999, ApJ, 512, 843 [NASA ADS] [CrossRef] [Google Scholar]

[24] Carone, L., Keppens, R., & Decin, L. 2016, MNRAS, 463, 3114 [NASA ADS] [CrossRef] [Google Scholar]

[25] Carone, L., Baeyens, R., Mollière, P., et al. 2020, MNRAS, 496, 3582 [NASA ADS] [CrossRef] [Google Scholar]

[26] Chan, Y. M., & Dalgarno, A. 1965, Proc. Phys. Soc., 85, 227 [NASA ADS] [CrossRef] [Google Scholar]

[27] Christie, D. A., Mayne, N. J., Lines, S., et al. 2021, MNRAS, 506, 4500 [CrossRef] [Google Scholar]

[28] Chubb, K. L., Rocchetto, M., Yurchenko, S. N., et al. 2021, A&A, 646, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[29] Coles, P. A., Yurchenko, S. N., & Tennyson, J. 2019, MNRAS, 490, 4638 [CrossRef] [Google Scholar]

[30] Crossfield, I. J. M., Knutson, H., Fortney, J., et al. 2012, ApJ, 752, 81 [Google Scholar]

[31] Dalgarno, A., & Williams, D. A. 1962, ApJ, 136, 690 [NASA ADS] [CrossRef] [Google Scholar]

[32] Dang, L., Cowan, N. B., Schwartz, J. C., et al. 2018, Nat. Astron., 2, 220 [NASA ADS] [CrossRef] [Google Scholar]

[33] Deitrick, R., Mendonça, J. M., Schroffenegger, U., et al. 2020, ApJS, 248, 30 [NASA ADS] [CrossRef] [Google Scholar]

[34] Dobbs-Dixon, I., & Agol, E. 2013, MNRAS, 435, 3159 [NASA ADS] [CrossRef] [Google Scholar]

[35] Drummond, B., Mayne, N. J., Baraffe, I., et al. 2018a, A&A, 612, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[36] Drummond, B., Mayne, N. J., Manners, J., et al. 2018b, ApJ, 855, L31 [NASA ADS] [CrossRef] [Google Scholar]

[37] Dullemond, C. P., van Zadelhoff, G. J., & Natta, A. 2002, A&A, 389, 464 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[38] Evans, T. M., Aigrain, S., Gibson, N., et al. 2015, MNRAS, 451, 680 [Google Scholar]

[39] Feautrier, P. 1964, C. R. Acad. Sci., 258, 3189 [NASA ADS] [Google Scholar]

[40] Fleischer, A. J., Gauger, A., & Sedlmayr, E. 1992, A&A, 266, 321 [NASA ADS] [Google Scholar]

[41] Fortney, J. J., Visscher, C., Marley, M. S., et al. 2020, AJ, 160, 288 [NASA ADS] [CrossRef] [Google Scholar]

[42] Fortney, J. J., Dawson, R. I., & Komacek, T. D. 2021, J. Geophys. Res. Planets, 126, e06629 [NASA ADS] [CrossRef] [Google Scholar]

[43] Freytag, B., Allard, F., Ludwig, H. G., Homeier, D., & Steffen, M. 2010, A&A, 513, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[44] Freytag, B., Liljegren, S., & Höfner, S. 2017, A&A, 600, A137 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[45] Freytag, B., Höfner, S., & Liljegren, S. 2019, IAU Symp., 343, 9 [NASA ADS] [Google Scholar]

[46] Fromang, S., Leconte, J., & Heng, K. 2016, A&A, 591, A144 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[47] Gillon, M., Triaud, A. H. M. J., Fortney, J. J., et al. 2012, A&A, 542, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[48] Goody, R., West, R., Chen, L., & Crisp, D. 1989, J. Quant. Spectr. Rad. Transf., 42, 539 [NASA ADS] [CrossRef] [Google Scholar]

[49] Gray, D. F. 2008, The Observation and Analysis of Stellar Photospheres (Cambridge: Cambridge University Press) [Google Scholar]

[50] Gustafsson, B., Bell, R. A., Eriksson, K., & Nordlund, A. 1975, A&A, 42, 407 [NASA ADS] [Google Scholar]

[51] Gustafsson, B., Edvardsson, B., Eriksson, K., et al. 2008, A&A, 486, 951 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[52] Hanasoge, S., Miesch, M. S., Roth, M., et al. 2015, Space Sci. Rev., 196, 79 [NASA ADS] [CrossRef] [Google Scholar]

[53] Hartman, J. D., Bakos, G. Á., Bhatti, W., et al. 2016, AJ, 152, 182 [NASA ADS] [CrossRef] [Google Scholar]

[54] Helling, C., Iro, N., Corrales, L., et al. 2019, A&A, 631, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[55] Helling, C., Lewis, D., Samra, D., et al. 2021, A&A, 649, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[56] Heng, K. 2012, ApJ, 761, L1 [NASA ADS] [CrossRef] [Google Scholar]

[57] Heng, K., Menou, K., & Phillipps, P. J. 2011, MNRAS, 413, 2380 [NASA ADS] [CrossRef] [Google Scholar]

[58] Hindle, A. W., Bushby, P. J., & Rogers, T. M. 2021, ApJ, 922, 176 [CrossRef] [Google Scholar]

[59] Höfner, S., & Freytag, B. 2019, A&A, 623, A158 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[60] Hubeny, I. 2017, MNRAS, 469, 841 [NASA ADS] [CrossRef] [Google Scholar]

[61] Husser, T. O., Wendevon Berg, S., Dreizler, S., et al. 2013, A&A, 553, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[62] Iess, L., Militzer, B., Kaspi, Y., et al. 2019, Science, 364, aat2965 [Google Scholar]

[63] Kaspi, Y., Galanti, E., Hubbard, W. B., et al. 2018, Nature, 555, 223 [Google Scholar]

[64] Kataria, T., Showman, A. P., Lewis, N. K., et al. 2013, ApJ, 767, 76 [NASA ADS] [CrossRef] [Google Scholar]

[65] Kataria, T., Showman, A. P., Fortney, J. J., et al. 2015, ApJ, 801, 86 [NASA ADS] [CrossRef] [Google Scholar]

[66] Kataria, T., Sing, D. K., Lewis, N. K., et al. 2016, ApJ, 821, 9 [NASA ADS] [CrossRef] [Google Scholar]

[67] Koll, D. D. B., & Abbot, D. S. 2016, ApJ, 825, 99 [NASA ADS] [CrossRef] [Google Scholar]

[68] Kreidberg, L., Bean, J. L., Désert, J.-M., et al. 2014, ApJ, 793, L27 [Google Scholar]

[69] Leconte, J. 2021, A&A, 645, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[70] Lee, E., Dobbs-Dixon, I., Helling, C., Bognar, K., & Woitke, P. 2016, A&A, 594, A48 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[71] Lee, E. K. H., Lothringer, J. D., Casewell, S. L., et al. 2022, MNRAS, submitted [arXiv:2203.09854] [Google Scholar]

[72] Lee, E. K. H., Parmentier, V., Hammond, M., et al. 2021, MNRAS, 506, 2695 [NASA ADS] [CrossRef] [Google Scholar]

[73] Li, J., & Goodman, J. 2010, ApJ, 725, 1146 [NASA ADS] [CrossRef] [Google Scholar]

[74] Li, J., & Shibata, K. 2006, J. Atmos. Sci., 63, 1365 [NASA ADS] [CrossRef] [Google Scholar]

[75] Li, G., Gordon, I. E., Rothman, L. S., et al. 2015, ApJS, 216, 15 [NASA ADS] [CrossRef] [Google Scholar]

[76] Line, M. R., Stevenson, K. B., Bean, J., et al. 2016, AJ, 152, 203 [NASA ADS] [CrossRef] [Google Scholar]

[77] Lines, S., Mayne, N. J., Manners, J., et al. 2019, MNRAS, 488, 1332 [NASA ADS] [CrossRef] [Google Scholar]

[78] Liu, B., & Showman, A. P. 2013, ApJ, 770, 42 [NASA ADS] [CrossRef] [Google Scholar]

[79] Manabe, S., Smagorinsky, J., & Strickler, R. F. 1965, Mon. Weather Rev., 93, 769 [NASA ADS] [CrossRef] [Google Scholar]

[80] Mayne, N. J., Baraffe, I., Acreman, D. M., et al. 2014, A&A, 561, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[81] Mayne, N. J., Debras, F., Baraffe, I., et al. 2017, A&A, 604, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[82] Mayne, N. J., Drummond, B., Debras, F., et al. 2019, ApJ, 871, 56 [Google Scholar]

[83] McKemmish, L. K., Yurchenko, S. N., & Tennyson, J. 2016, MNRAS, 463, 771 [NASA ADS] [CrossRef] [Google Scholar]

[84] McKemmish, L. K., Masseron, T., Hoeijmakers, H. J., et al. 2019, MNRAS, 488, 2836 [Google Scholar]

[85] Mendonça, J. M. 2020, MNRAS, 491, 1456 [CrossRef] [Google Scholar]

[86] Mendonça, J. M., Malik, M., Demory, B.-O., & Heng, K. 2018a, AJ, 155, 150 [Google Scholar]

[87] Mendonça, J. M., Tsai, S.-M., Malik, M., Grimm, S. L., & Heng, K. 2018b, ApJ, 869, 107 [CrossRef] [Google Scholar]

[88] Menou, K. 2019, MNRAS, 485, L98 [Google Scholar]

[89] Menou, K., & Rauscher, E. 2009, ApJ, 700, 887 [NASA ADS] [CrossRef] [Google Scholar]

[90] Mollière, P., van Boekel, R., Dullemond, C., Henning, T., & Mordasini, C. 2015, ApJ, 813, 47 [Google Scholar]

[91] Mollière, P., van Boekel, R., Bouwman, J., et al. 2017, A&A, 600, A10 [Google Scholar]

[92] Mollière, P., Wardenier, J. P., van Boekel, R., et al. 2019, A&A, 627, A67 [Google Scholar]

[93] Mollière, P., Stolker, T., Lacour, S., et al. 2020, A&A, 640, A131 [Google Scholar]

[94] Mordasini, C., van Boekel, R., Mollière, P., Henning, T., & Benneke, B. 2016, ApJ, 832, 41 [Google Scholar]

[95] Ng, K. C. 1974, J. Chem. Phys., 61, 2680 [NASA ADS] [CrossRef] [Google Scholar]

[96] Olson, G. L., Auer, L. H., & Buchler, J. R. 1986, J. Quant. Spectr. Rad. Transf., 35, 431 [Google Scholar]

[97] Parmentier, V., Guillot, T., Fortney, J. J., & Marley, M. S. 2015, A&A, 574, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[98] Parmentier, V., Showman, A. P., & Fortney, J. J. 2021, MNRAS, 501, 78 [Google Scholar]

[99] Peixóto, J. P., & Oort, A. H. 1984, Rev. Mod. Phys., 56, 365 [CrossRef] [Google Scholar]

[100] Perna, R., Heng, K., & Pont, F. 2012, ApJ, 751, 59 [NASA ADS] [CrossRef] [Google Scholar]

[101] Piette, A. A. A., & Madhusudhan, N. 2020, ApJ, 904, 154 [NASA ADS] [CrossRef] [Google Scholar]

[102] Piskunov, N. E., Kupka, F., Ryabchikova, T. A., Weiss, W. W., & Jeffery, C. S. 1995, A&AS, 112, 525 [Google Scholar]

[103] Polyansky, O. L., Kyuberis, A. A., Zobov, N. F., et al. 2018, MNRAS, 480, 2597 [NASA ADS] [CrossRef] [Google Scholar]

[104] Rauscher, E., & Menou, K. 2010, ApJ, 714, 1334 [NASA ADS] [CrossRef] [Google Scholar]

[105] Rauscher, E., & Menou, K. 2012, ApJ, 750, 96 [NASA ADS] [CrossRef] [Google Scholar]

[106] Richard, C., Gordon, I. E., Rothman, L. S., et al. 2012, J. Quant. Spectr. Rad. Transf., 113, 1276 [Google Scholar]

[107] Rogers, T. M., & Showman, A. P. 2014, ApJ, 782, L4 [NASA ADS] [CrossRef] [Google Scholar]

[108] Roman, M. T., Kempton, E. M. R., Rauscher, E., et al. 2021, ApJ, 908, 101 [CrossRef] [Google Scholar]

[109] Rutten, R. J. 2003, Radiative Transfer in Stellar Atmospheres (Utrecht: Lecture Notes Utrecht University) [Google Scholar]

[110] Sainsbury-Martinez, F., Casewell, S. L., Lothringer, J. D., Phillips, M. W., & Tremblin, P. 2021, A&A, 656, A128 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[111] Sainsbury-Martinez, F., Wang, P., Fromang, S., et al. 2019, A&A, 632, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[112] Schneider, A. D., & Bitsch, B. 2021a, A&A, 654, A71 [Google Scholar]

[113] Schneider, A. D., & Bitsch, B. 2021b, A&A, 654, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[114] Seager, S. 2010, Exoplanet Atmospheres: Physical Processes (Princeton: Princeton University Press) [Google Scholar]

[115] Sergeev, D. E., Lambert, F. H., Mayne, N. J., et al. 2020, ApJ, 894, 84 [NASA ADS] [CrossRef] [Google Scholar]

[116] Showman, A. P., & Kaspi, Y. 2013, ApJ, 776, 85 [NASA ADS] [CrossRef] [Google Scholar]

[117] Showman, A. P., & Polvani, L. M. 2011, ApJ, 738, 71 [Google Scholar]

[118] Showman, A. P., Fortney, J. J., Lian, Y., et al. 2009, ApJ, 699, 564 [NASA ADS] [CrossRef] [Google Scholar]

[119] Showman, A. P., Lewis, N. K., & Fortney, J. J. 2015, ApJ, 801, 95 [NASA ADS] [CrossRef] [Google Scholar]

[120] Showman, A. P., Tan, X., & Parmentier, V. 2020, Space Sci. Rev., 216, 139 [Google Scholar]

[121] Sing, D. K., Fortney, J. J., Nikolov, N., et al. 2016, Nature, 529, 59 [Google Scholar]

[122] Sing, D. K., Lavvas, P., Ballester, G. E., et al. 2019, AJ, 158, 91 [Google Scholar]

[123] Skinner, J. W., & Cho, J. Y. K. 2021, MNRAS, 504, 5172 [NASA ADS] [CrossRef] [Google Scholar]

[124] Sousa-Silva, C., Al-Refaie, A. F., Tennyson, J., & Yurchenko, S. N. 2015, MNRAS, 446, 2337 [Google Scholar]

[125] Southworth, J. 2008, MNRAS, 386, 1644 [NASA ADS] [CrossRef] [Google Scholar]

[126] Southworth, J. 2010, MNRAS, 408, 1689 [Google Scholar]

[127] Steinrueck, M. E., Showman, A. P., Lavvas, P., et al. 2021, MNRAS, 504, 2783 [NASA ADS] [CrossRef] [Google Scholar]

[128] Stevenson, K. B., Désert, J.-M., Line, M. R., et al. 2014, Science, 346, 838 [Google Scholar]

[129] Thorngren, D. P., & Fortney, J. J. 2018, AJ, 155, 214 [Google Scholar]

[130] Thorngren, D., Gao, P., & Fortney, J. J. 2019, ApJ, 884, L6 [Google Scholar]

[131] Tremblin, P., Chabrier, G., Mayne, N. J., et al. 2017, ApJ, 841, 30 [NASA ADS] [CrossRef] [Google Scholar]

[132] Tremblin, P., Padioleau, T., Phillips, M. W., et al. 2019, ApJ, 876, 144 [Google Scholar]

[133] Venot, O., Parmentier, V., Blecic, J., et al. 2020, ApJ, 890, 176 [NASA ADS] [CrossRef] [Google Scholar]

[134] Wang, H., & Wordsworth, R. 2020, ApJ, 891, 7 [NASA ADS] [CrossRef] [Google Scholar]

[135] Wende, S., Reiners, A., Seifahrt, A., & Bernath, P. F. 2010, A&A, 523, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[136] Woitke, P., & Niccolini, G. 2005, A&A, 433, 1101 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[137] Yurchenko, S. N., Amundsen, D. S., Tennyson, J., & Waldmann, I. P. 2017, A&A, 605, A95 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[138] Yurchenko, S. N., Mellor, T. M., Freedman, R. S., & Tennyson, J. 2020, MNRAS, 496, 5282 [NASA ADS] [CrossRef] [Google Scholar]

[139] Zhang, X. 2020, Res. Astron. Astrophys., 20, 099 [NASA ADS] [CrossRef] [Google Scholar]