# Modeling the day-night temperature variations of ultra-hot Jupiters:
confronting non-grey general circulation models and observations

Xianyu Tan1,2,3, Thaddeus D. Komacek4, Natasha E. Batalha5, Drake Deming4,
Roxana Lupu6, Vivien Parmentier7, and Raymond T. Pierrehumbert3
1Tsung-Dao Lee Institute, Shanghai Jiao Tong University, 520 Shengrong Road,
Shanghai, People’s Republic of China
2School of Physics and Astronomy, Shanghai Jiao Tong University, 800 Dongchuan
Road, Shanghai, People’s Republic of China
3Atmospheric Oceanic and Planetary Physics, Department of Physics, University
of Oxford, OX1 3PU, UK
4Department of Astronomy, University of Maryland, College Park, MD 20742, USA
5NASA Ames Research Center, Moffett Field, CA 94035, USA
6Eureka Scientific Inc, Oakland, CA 94602, USA
7Université Côte d’Azur, Observatoire de la Côte d’Azur, CNRS, Laboratoire
Lagrange, Nice, France E-mail<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Ultra-hot Jupiters (UHJs) are natural laboratories to study extreme physics in
planetary atmospheres and their rich observational data sets are yet to be
confronted with models with varying complexities at a population level. In
this work, we update the general circulation model of Tan & Komacek (2019) to
include a non-grey radiative transfer scheme and apply it to simulate the
realistic thermal structures, phase-dependent spectra, and wavelength-
dependent phase curves of UHJs. We performed grids of models over a large
range of equilibrium temperatures and rotation periods for varying
assumptions, showing that the fractional day-night brightness temperature
differences remain almost constant or slightly increase with increasing
equilibrium temperature from the visible to mid-infrared wavelengths. This
differs from previous work primarily due to the increasing planetary rotation
rate with increasing equilibrium temperature for fixed host star type.
Radiative effects of varying atmospheric compositions become more significant
in dayside brightness temperature in longer wavelengths. Data-model
comparisons of dayside brightness temperatures and phase curve amplitudes as a
function of equilibrium temperature are in broad agreement. Observations show
a large scatter compared to models even with a range of different assumptions,
indicating significantly varying intrinsic properties in the hot Jupiter
population. Our cloud-free models generally struggle to match all observations
for individual targets with a single set of parameter choices, indicating the
need for extra processes for understanding the heat transport of UHJs.

###### keywords:

hydrodynamics — methods: numerical — planets and satellites: atmospheres —
planets and satellites: gaseous planets

††pubyear: 2022††pagerange: Modeling the day-night temperature variations of
ultra-hot Jupiters: confronting non-grey general circulation models and
observations–B

##  1 Introduction

Ultra-hot Jupiters are gas giant planets in extremely close-in orbits around
their host stars, with full-redistribution equilibrium temperatures in excess
of 2200 K. The enormous stellar irradiation on their dayside places their
atmospheres between the cooler hot Jupiters and those of late-type stars,
making them ideal laboratories to understand the transition of stellar to
planetary atmospheres and the atmospheric physics in extreme conditions. Owing
to the ease of observability, blooming space-based (Kreidberg et al., 2018;
Mansfield et al., 2018; Coulombe et al., 2023; Mikal-Evans et al., 2022) and
ground-based (Ehrenreich et al., 2020; Borsa et al., 2021; Pino et al., 2022)
observations over the past years have revealed the phase-dependent atmospheric
conditions of these objects.

A notable difference between the atmospheres of hot and ultra-hot Jupiters is
expected to be the thermal dissociation of molecules which could significantly
affect the thermal structure and observations of the atmospheres (Kitzmann et
al., 2018; Lothringer et al., 2018; Parmentier et al., 2018). Key radiatively
active molecules could be thermally dissociated and thus could impact the
radiative balance and resulting emission spectra of these objects. Thermal
dissociation of the major hydrogen molecules could significantly affect the
dynamics of these atmospheres, determining the day-to-night heat transport
efficiency and the detailed 3D circulation patterns (Bell & Cowan, 2018;
Komacek & Tan, 2018; Tan & Komacek, 2019; Roth et al., 2021).

Comparisons between observations and general circulation models (GCMs) are
vital to uncovering physical processes shaping the three-dimensional (3D)
atmospheres of hot Jupiters. Two modeling strategies are commonly used in
literature. The first approach is to model individual planets that have
detailed observations and models can be fine-tuned with different assumptions
to match many aspects of the observations including spectrum, the amplitudes
of the phase curves, phase offsets, and the detailed phase-curve shapes. This
approach has been conducted on many individual targets, revealing in-depth
characterization of their individual atmospheres (e.g., Showman et al., 2009;
Zellem et al., 2014; Stevenson et al., 2017; Kreidberg et al., 2018; Arcangeli
et al., 2019; Steinrueck et al., 2019; Lines et al., 2019; Mansfield et al.,
2020; Carone et al., 2020; May et al., 2021; Mikal-Evans et al., 2022; Lee et
al., 2022; Beltz et al., 2022b; Schneider et al., 2022; Zamyatina et al.,
2023; Lee et al., 2023).

Another approach, instead of focusing on explaining data for individual
targets, aims to model key trends as a function of system parameters and to
identify the common physical mechanisms that possibly sculpt some key
behaviors of the planetary population by systematic data-model comparisons
over the parameter space. Perna et al. (2012), Perez-Becker & Showman (2013),
and Komacek et al. (2017) investigated the trend of day-night temperature
variations under dynamical wave adjustment, radiative damping, and possible
wind dissipation mechanisms, showing that the fractional day-night temperature
variation is expected to increase with equilibrium temperature. Tan & Komacek
(2019) modeled the ultra-hot Jupiter population using the updated GCM with the
semi-grey radiative transfer and dynamical effects of hydrogen dissociation
and recombination that are proposed to be important in ultra-hot Jupiters
(Bell & Cowan, 2018), and showed that the fractional day-night temperature
variation may no longer increase with equilibrium temperature in the ultra-hot
regime. While these studies are illuminating in understanding the dynamical
mechanisms, the simplified radiative transfer calculations were not sufficient
to precisely model the atmospheric thermal structures and so insufficient to
facilitate fair comparisons to phase-curve observations which are
intrinsically wavelength dependent. Parmentier et al. (2016) and Parmentier et
al. (2021) carried out large grids of GCM simulations using the non-grey GCM
originally developed in Showman et al. (2009) and were able to pinpoint the
roles of clouds and dynamics affecting the wavelength-dependent phase-curve
observations of the canonical hot Jupiters population. This work is now being
pursued by Roth et al. (submitted), which is looking at a grid of several
thousand models. However, neither the dynamical core of Parmentier et al.
(2021) nor of Roth et al. (submitted) include the important effects of
hydrogen dissociation and recombination and, therefore did not extend fully to
the ultra-hot-Jupiter regime. In summary, there is now a gap in the GCM
modeling of the ultra-hot Jupiter population that needs to be filled with more
complete GCM simulations.

In this work, using GCMs with an updated non-grey radiative transfer scheme,
we improve upon our previous work of Tan & Komacek (2019) and perform
simulations over grids of varying equilibrium temperature to understand the
trends of the day-night temperature variation of ultra-hot Jupiters. By
comparing the modeling results to the rich observations of the ultra-hot
Jupiter population, we wish to identify whether our model physics is adequate
to explain the observed trends, or if not, what’s likely missing in our
models. We use the SPARC/MITgcm, originally developed in Showman et al.
(2009), and couple it with the hydrogen dissociation and recombination
treatments in Tan & Komacek (2019) and Komacek et al. (2022a). The realistic
radiative heating and cooling rates in our updated models guarantee more
precise modeling of the thermal structures, dynamics and chemistry, as well as
more self-consistent comparisons to wavelength-dependent phase-curve
observations.

Our work is also complementary to a variety of existing GCM parameter surveys
of hot Jupiter’s atmospheres. Some of these parameter studies do not perform a
systematic model-data comparison but focus on the theoretical understanding of
various physical and chemical processes over a large parameter space. Roman et
al. (2021) included temperature-dependent clouds with radiative feedback to
assess the effects on phase curves as a function of equilibrium temperature.
Helling et al. (2022) computed cloud microphysics based on a large grid of GCM
outputs, showing cloud inhomogeneity over a wide range of planetary
parameters. Baeyens et al. (2021) conducted a grid of GCM with a Newtonian
cooling scheme and outputs are used to calculate chemical mixing, unveiling
trends in the atmospheric chemistry of hot Jupiters.

The outline of this work is as follows. Section 2 describes the setup of the
GCM simulations, radiative transfer post-processing using PICASO, and system
parameters. In Section 3, we present results from the main suites of the model
grids as a function of planetary equilibrium temperature assuming solar
atmospheric metallicity, including the basic climatology, spectral, and phase-
curve amplitude calculations. We summarize our key findings in Section 4. In
the Appendix B, motivated by the systematic mismatch of our results to the
observed TESS day-night temperature variations, we explore additional effects
of varying atmospheric metallicity and carbon-to-oxygen ratio on the phase-
curve amplitudes, using TOI 1431b as a template.

##  2 Model Setup

###  2.1 The 3D general circulation model: SPARC/MITgcm

The global atmospheric circulation and thermal structure of hot Jupiters were
simulated using the SPARC/MITgcm (Showman et al., 2009). The model solves the
global primitive equations of dynamical meteorology on a cubed-sphere grid.
This equation set is suitable to model the outermost observable layers of
giant planets and brown dwarfs. Radiative heating and cooling rates in each
column of the GCM are calculated using the non-grey radiative transfer model
of Marley et al. (1999) that solves the two-stream radiative transfer
equations and employs the correlated-k method. The correlated-k opacity tables
have been updated based on the opacities used in Marley et al. (2021). We
assume thermochemical equilibrium in the atmosphere. In addition to molecular
and atomic opacities in k coefficients listed in Table 1, the radiative
transfer model also includes atomic line absorption from the neutral alkali
metals, continuum opacity sources including pressure-induced opacities, bound-
free and free-free opacity from ${\rm H}^{-}$ and ${\rm H_{2}}^{+}$ and free-
free opacity from ${\rm H_{2}}^{-}$ and ${\rm He}^{+}$, as well as electron
scattering and Rayleigh scattering from hydrogen and helium. See Marley et al.
(2021) for a detailed description. The k-coefficient tables used in the GCM
use 11 frequency bins as described in Kataria et al. (2015).

The hydrogen dissociation and recombination (${\rm H_{2}}$-${\rm H}$) is an
important heating and cooling process in atmospheres of ultra-hot Jupiters
(Bell & Cowan, 2018; Komacek & Tan, 2018). We utilize the same numerical
schemes as idealized GCMs in Tan & Komacek (2019) that implement
thermodynamically active tracers to track the evolution of atomic hydrogen
mixing ratio. Our current implementation is the same as that in Komacek et al.
(2022a) in which inert helium is included as the mass budget at Solar
abundance. An alternative approach that does not need to invoke tracers is to
cast the heating/cooling as a buffer effect which effectively increases the
heat capacity at constant pressure when there is dissociation or
recombination, as shown in a pseudo-2D framework of Roth et al. (2021). That
approach assumes the instantaneous thermal equilibrium of hydrogen. In our
tracer-based ${\rm H_{2}}$-${\rm H}$ scheme, we applied a “soft" adjustment in
which the conversion occurs over a finite timescale rather than
instantaneously reaching local ${\rm H_{2}}$-${\rm H}$ equilibrium. This soft
adjustment helps to maintain the numerical stability of our models by slightly
reducing the extreme horizontal temperature gradients locally occurring in the
model. When the ${\rm H_{2}}$-${\rm H}$ conversion timescale is very close to
the dynamical time step (mimicking instantaneous adjustment), models become
numerically unstable. The numerical challenge might stem from rapid dynamical
instabilities associated with the large horizontal temperature and ${\rm
H_{2}}$-H contrasts and they are difficult to be resolved within reasonable
time steps. After some tests, we choose a conservative ${\rm H_{2}}$-H
relaxation timescale of 40 s 111 Realistic chemical calculations (Tsai et al.,
2018, Shang-Min Tsai, personal communication) suggest that, for example, at
$10^{-3}$ bar, the expected chemical timescale of H2 thermal dissociation is
between $10^{2}$ to $10^{3}$ s at 2000 K and between 10 to $10^{2}$ s at 3000
K. Our choice of relaxation timescale = 40 s is reasonable at relatively low
pressures. At $10^{-1}$ bar, the expected timescale is well below 1 s for
temperatures higher than 1800 K. At relatively high pressures, the reaction
should be considered instantaneous in GCM step timescales. Dynamical
timescales relevant in large scales can be estimated as $R/c$ where $R$ is the
planetary radius and $c$ could be either horizontal wind speed or Kelvin wave
phase speed (both are on the order of several ${\rm km~{}s^{-1}}$), and are
typically on the order of $10^{4}$ to $10^{5}$ s. Our choice of chemical
relaxation timescale = 40 s is much smaller than the dynamical timescales.
(2-4 times the dynamical time step of the model) to guarantee stability and
this is applied for all models for consistency. Assuming a peak horizontal
wind speed of $\sim 10\;{\rm km\;s^{-1}}$, the minimum characteristic
“disequilibrium" length scale of hydrogen is $10^{5}-10^{6}$ m which is much
smaller than the planetary radius. Note that the atomic hydrogen mixing ratio
in our tracer scheme which affects the dynamics and heating/cooling rates
could be slightly different from those calculated by the full equilibrium
chemistry calculation used in the radiative transfer scheme, but this is a
small trade-off in order to maintain numerical stability.

A Rayleigh drag with a coefficient of drag timescale of $10^{5}$ s at the
bottom and linearly decreasing to zero with decreasing pressure until a
pressure of 20 bars is applied in all models. This formulation is similar to
that in Liu & Showman (2013) and mimics interactions with the interior. In
practice, this “basal” drag scheme helps with model convergence and numerical
stability. In addition, another drag with a globally uniform drag timescale
$\tau_{\mathrm{drag}}$ is applied in a subset of the models to test the
sensitivity of circulation and phase curves. In the main text, by "models
without drag", we mean models without the uniform drag but still retaining the
basal drag.

The temperature of the lowest model layer center (at a pressure of about 85
bars) in the GCM is relaxed towards a prescribed value of 3344 K over a
timescale of $10^{2}$ s. This deep atmospheric profile is roughly consistent
with those presented in Thorngren et al. (2019) for inflated hot Jupiters.
However, hot Jupiter inflation varies with planetary $T_{\rm eq}$ and the
bottom temperature should in principle vary with $T_{\rm eq}$. Our prescribed
value is somewhat conservative for ultra-hot Jupiters because we want to
minimize the impact of the (uncertain) interior heat flux on the day-night
temperature difference. This value is also held fixed in our grids to allow us
to study the influence of only rotation period, host star type, and
irradiation on the dynamics and heat transport. Our bottom boundary condition
differs from previous setups of hot Jupiters (e.g., Showman et al., 2009) in
which a uniform upwelling heat flux was prescribed. Our bottom boundary
condition mimics a rapid entropy mixing by vigorous interior convection. We
argue that this is more appropriate here because convection is expected to
homogenize the entropy (towards the uniform interior entropy) rather than heat
flux, and the bottom of our domain is expected to lie below the radiative-
convective boundary (Thorngren et al., 2019). Using a similar hydrogen
dissociation modeling framework as here, Komacek et al. (2022b) recently
showed that different bottom boundary treatments could result in differences
in the isobaric temperature and circulation near the photosphere.

###  2.2 Radiative transfer post-processing: PICASO

Based on GCM outputs, we calculate the time-dependent planet spectra and
observables using PICASO which is an open-source radiative transfer code for
computing the reflected, thermal, and transmission spectrum of planets and
brown dwarfs (Batalha et al., 2019; Mukherjee et al., 2023).222Publicly
available on Github https://github.com/natashabatalha/picaso The methodology
of PICASO to calculate the radiative transfer for thermal emission is the same
as those in Marley et al. (1999). Emission, reflection spectra and phase curve
calculations by PICASO based on 3D GCM outputs of hot Jupiters have been
carried out in a few studies (Adams et al., 2022; Robbins-Blanch et al.,
2022).

PICASO has the same opacity source as the SPARC/MITgcm and this guarantees
energetically consistent treatment between the GCM outputs and post-processing
and can rigorously test effects of heat transport and ${\rm H_{2}}$-${\rm H}$
dissociation on the spectra and wavelength-dependent phase curves. In this
work, we use the mode of correlated-k methods in PICASO with pre-calculated
k-coefficients and equilibrium chemistry. These k tables have a much higher
frequency resolution of 622 bins from about 0.26 to 267 $\mu$m than the
version used in the GCM.

Parameter | Value
---|---
Stellar effective temperature | [5500, 6000, 6500] K
Stellar mass |  [0.94, 1.08, 1.34] $M_{\rm sun}$
Stellar radius |  [0.91, 1.10, 1.50] $R_{\rm sun}$
Planet radius | $10^{8}{\rm m}$
Planet gravity |  10 m s-2
Planet equilibrium temperature | [1800, 2000, 2200, 2400, 2600] K
Planet rotation period | [1.05, 0.77, 0.58, 0.44, 0.35] days
| [1.69, 1.23, 0.93, 0.71, 0.56] days
| [3.07, 2.24, 1.68, 1.30, 1.06] days
Temperature at 85 bars | 3344 K
H tracer relaxation timescale | 40 s
Drag timescale ($\tau_{\mathrm{drag}}$) |  $\infty$, $10^{4}$ s
Metallicity and C/O ratio |  $1\times$ solar
Species in k-coefficient tables |  ${\rm H_{2}O}$, ${\rm CO_{2}}$,${\rm CO}$, ${\rm CH_{4}}$, ${\rm NH_{3}}$, ${\rm OCS}$,
|  ${\rm N_{2}}$, ${\rm HCN}$, ${\rm C_{2}H_{2}}$, ${\rm C_{2}H_{6}}$, ${\rm
H_{2}S}$,
|  ${\rm C_{2}H_{4}}$, ${\rm FeH}$, ${\rm CrH}$, ${\rm LiCl}$, ${\rm SiO}$,
${\rm MgH}$,
|  [ ${\rm TiO}$, ${\rm VO}$, ${\rm Fe}$, ${\rm H_{2}}$, ${\rm H_{3}^{+}}$,
${\rm LiH}$, ${\rm LiF}$ ]
Horizontal resolution | C32, C48
Vertical resolution | 53 layers
Lower boundary | 100 bars
Upper boundary |  $5\times 10^{-6}$ bar
Shapiro filter order | 4
Shapiro filter timescale | 50 s
Dynamical time step | 10 - 20 s
Radiative time step | 50 s
Table 1: Stellar, planetary, and numerical parameters assumed for the suite of
GCM simulations presented Section 3. The last bracketed line of species in
k-coefficient tables indicates parameter sensitivity tests over a subset of
models in which these species are not included in the radiative transfer.

###  2.3 Parameters of our grid simulations

Following up our previous approach using idealized GCMs (Tan & Komacek, 2019)
and many other GCM studies, we sweep a wide parameter space in order to
understand the trend in the ultra-hot Jupiter population.

We specify a stellar effective temperature, mass, and radius, and then a set
of orbital periods can be determined for a set of zero-albedo planetary
equilibrium temperatures. We assume three stellar effective temperatures
$T_{\rm star}$ of 5500 K, 6000 K, and 6500 K which should be representative of
most stellar hosts of hot Jupiters. Using stellar mass-radius-temperature
relations found by Eker et al. (2018) based on solar neighbourhoods, we set
the stellar masses $M_{\rm star}$ as 0.94$M_{\rm sun}$, 1.08$M_{\rm sun}$ and
1.34$M_{\rm sun}$, and the stellar radius $R_{\rm star}$ as 0.91$R_{\rm sun}$,
1.1$R_{\rm sun}$ and 1.5$R_{\rm sun}$, respectively. With a given planetary
equilibrium temperature $T_{\rm eq}$ assuming zero albedos, the orbital period
(therefore the rotation period $P_{\rm rot}$ assuming tidal locking) is given
by

$P_{\rm rot}=\frac{\pi}{\sqrt{2}}\left(\frac{T_{\rm star}}{T_{\rm
eq}}\right)^{3}\sqrt{\frac{R^{3}_{\rm star}}{GM_{\rm star}}},$ (1)

where $G$ is the gravitational constant. Within the planetary $T_{\rm eq}$
that we consider in this work, tidal synchronization is a reasonable
assumption (Guillot et al., 1996). As will be discussed below, the change of
rotation rate of tidally locked planets corresponding to different stellar
properties and $T_{\rm eq}$ is the dominant effect on the circulation and day-
to-night heat transport, rather than the nature of the incident stellar
spectra.

We model atmospheres with equilibrium temperatures from 1800 to 2600 K,
encompassing the majority of the ultra-hot Jupiter population. The sets of
rotation period for planetary equilibrium temperature of [1800, 2000, 2200,
2400, 2600] K for three stellar temperatures are [1.05, 0.77, 0.58, 0.44,
0.35] days, [1.69, 1.23, 0.93, 0.71, 0.56] days and [3.07, 2.24, 1.68, 1.3,
1.06] days, respectively with increasing stellar effective temperature. The
horizontal resolution of most simulations is C32, equivalent to about
$2.81^{\circ}\times 2.81^{\circ}$ in longitude and latitude; some rapidly
rotating, drag-free cases of $T_{\rm eq}=2200,2400$ and 2600 K with $T_{\rm
star}=5500$, and $T_{\rm eq}=2600$ K with $T_{\rm star}=6000$, use a C48
resolution (about $1.88^{\circ}\times 1.88^{\circ}$ in longitude and latitude)
to better resolve smaller horizontal structures.

In the main set of simulations of Section 3, we include key optical absorbing
agents including TiO, VO and Fe, which have been known to drive a strong
thermal inversion on the dayside atmosphere and can be detected on ultra-hot
Jupiters with high resolution spectroscopy (Hoeijmakers et al., 2018, 2020;
Ehrenreich et al., 2020; Kesseli et al., 2022; Pelletier et al., 2023). Models
without drag are integrated for 1350 to 2250 days, and models with
$\tau_{\mathrm{drag}}=10^{4}$ s are integrated for 625 to 1000 days, depending
on the timescale for each case to reach an energy equilibrium. Outputs of the
last 200 days time-averaged for the results shown in Section 3. Key parameters
of our models are summarized in Table 1.

Figure 1: Temperature maps at 55 mbar for drag-free models with varying
planetary equilibrium temperature from 1800 K to 2600 K (from top to bottom)
and varying tidally synchronized rotation periods implied by stellar effective
temperatures of 5500 K, 6000 K and 6500 K (from left to right). Horizontal
winds are represented by streamlines and arrows which are normalized
separately for each panel. The rotation period of each model is labeled above
each panel. These results were obtained by time-averaging over the last
200-day drag-free model outputs. Temperature color maps of the same planetary
equilibrium temperature (e.g., in each row) share the same color bar on the
right.  Figure 2: Same as Figure 1 but showing temperature maps at 0.6 mbar.
Figure 3: Zonal-mean zonal winds for drag-free models with varying planetary
equilibrium temperature from 1800 K to 2600 K (from top to bottom) and varying
stellar effective temperature of 5500 K, 6000 K, and 6500 K (from left to
right). The rotation period of each model is labeled above each panel. These
results were obtained by time-averaging over the last 200-day drag-free model
outputs.

##  3 Results of the main suits of models

###  3.1 Basic circulation patterns and spectra

This subsection describes the temperature maps, zonal jets, dayside and
nightside temperature-pressure (TP) structures, effect of ${\rm H_{2}}$-${\rm
H}$ dissociation on heat transport, and the dayside and nightside thermal
spectra resulting from the circulation. These results are based on the main
suites of simulations assuming 1x solar compositions described in Table 1.

####  3.1.1 Temperature fields and jets

Figure 1 shows the horizontal temperature maps at 55 mbar which is near the
thermal photosphere level for a subset of the models without drag. Equilibrium
temperature increases from 1800 K to 2600 K from the top to the bottom, and
stellar temperature increases from 5500 K to 6500 K from left to right.
Streamlines with arrows represent horizontal velocity fields. All models show
qualitatively consistent features, including a prominent equatorial eastward
jet and stationary global-scale waves which are shifted eastward relatively to
the sub-stellar longitude of $0^{\circ}$ and extend from the equator to the
poles. These are the canonical circulation patterns that emerge from the
strongly irradiated, moderately rotating, tidally locked hot Jupiter dynamical
regimes (e.g., Showman & Guillot, 2002; Showman et al., 2009; Heng et al.,
2011; Rauscher & Menou, 2012; Mayne et al., 2014; Cho et al., 2015; Mendonça
et al., 2016; Hammond & Pierrehumbert, 2018; Tan & Komacek, 2019; Carone et
al., 2020; Lee et al., 2021; Hammond & Lewis, 2021, see also reviews by Heng &
Showman, 2015, Zhang, 2020, and Showman et al., 2020). At a given stellar
temperature, the extent to which the global thermal pattern is shifted
eastward of the sub-stellar longitude lessens with increasing equilibrium
temperatures. This is due to both the increasing radiative heating efficiency
and the stronger rotation, both of which act against horizontal heat transport
(Showman et al., 2013a; Tan & Showman, 2020). This trend is the same for all
three stellar temperatures.

At a lower pressure of 0.6 mbar as shown in Figure 2, the dayside atmospheres
are hotter due to the strong thermal inversion caused by optical absorbers,
and the atomic hydrogen mixing ratio is at higher values due to thermal
dissociation. This larger sub-stellar-to-limb atomic hydrogen contrast results
in a higher horizontal heat transport efficiency within the dayside from the
heating due to hydrogen recombination on the nightside and limb and cooling
due to dissociation near the substellar point (Tan & Komacek, 2019). Compared
to those at 55 mbar, the hot area is more significantly filled in the dayside.
Furthermore, for the hotter cases, the hot areas extend substantially to the
nightside through the circulation mostly at high latitudes. The extent to
which the pattern is shifted eastward of the sub-stellar longitude is weaker
than those shown at 55 mbar, similar to that shown for canonical hot Jupiters
(Showman et al., 2013b; Lee et al., 2021). Rotation is vital to determine
where the heat is preferentially transported. For example, on the left column
where the rotation periods are shorter for a given equilibrium temperature,
flows at high latitudes are strongly confined by the Coriolis force and form
Rossby gyres, whereas on the right column where rotation periods are longer,
flows are less confined by the Coriolis force, with a weaker Rossby gyre and
greater dayside to nightside (divergent) flow (Hammond & Lewis, 2021) which
efficiently warms the high latitude regions on the nightside.

All models show a dominant equatorial superrotating jet with speeds up to
about 10 ${\rm km~{}s^{-1}}$ and most of the jets extend from nearly the model
top boundary down to almost 10 bars where the basal drag starts to dissipate
the horizontal winds. Figure 3 shows the time-averaged zonal-mean zonal winds
for the same set of models shown in Figure 1. In models with a fixed stellar
temperature, the jet speed increases with increasing planetary equilibrium
temperature; in models with a specific planetary equilibrium temperature, the
jet speed decreases with increasing stellar temperature (and hence increasing
rotation period). These trends make sense intuitively. The first trend emerges
likely because the eddy velocities are expected to increase with increasing
equilibrium temperature in the absence of strong dissipation (Perez-Becker &
Showman, 2013; Komacek & Showman, 2016; Koll & Komacek, 2018); the
superrotating jet is driven by horizontal convergence of the eddies (Showman &
Polvani, 2011), therefore it is suspected that the jet speed is also
positively correlated with the eddy magnitude which has been shown by
numerical simulations and scaling analysis (Zhang & Showman, 2017; Hammond et
al., 2020). On the second trend, given the same equilibrium temperature (hence
the same expected eddy magnitude), the magnitude of eddy convergence is
expected to increase when the eddy meridional lengthscale is smaller which is
expected for faster rotators (Showman & Polvani, 2011; Pierrehumbert &
Hammond, 2019).

The meridional width of the jets was expected to decrease with decreasing
rotation period due to the shrinking meridional extent of the eddies that is
proportional to the equatorial deformation radius
$L_{d}=\sqrt{\frac{c}{\beta}}$, where $c$ is the gravity wave phase speed,
$\beta=df/dy$ at the equator, $f=4\pi\sin\phi/P$ is the Coriolis parameter,
$P$ is the rotation period and $\phi$ is latitude.333The meridional extent of
the eastward acceleration on the zonal flow by the stationary Matsuno-Gill
waves could be modified by the presence of a background equatorial jet, see
Figure 7 of Hammond & Pierrehumbert (2018). Therefore, one might wonder what
wave speed is relevant – the one Doppler shifted by a background jet or in a
rest background (i.e., with no Doppler shift). Regardless of which wave speed
should be used, numerical results of the jet widths in the moderate-to-fast
rotating regime from Tan & Showman (2020) are well matched to the scaling of
rotation period. Therefore, we would expect that the jet width is proportional
to $\sqrt{P}$. In GCM simulations for both the non-synchronously rotating hot
and warm Jupiters (Showman et al., 2015) and synchronously rotating close-in
giant planets (Tan & Showman, 2020), the decreasing meridional width of the
jets with decreasing rotation periods are both investigated. In the rapid-to-
moderate rotating regime, the meridional jet width showed a clear linear
scaling to the equatorial deformation radius which is proportional to
$\sqrt{P}$ when everything else is held fixed (Tan & Showman, 2020). Indeed,
at the same $T_{\rm eq}$ but varying stellar temperatures and rotation periods
(horizontal comparisons in Figure 3), we see this qualitative trend of
increasing meridional jet with increasing rotation period.

However, this is not the case in the vertical comparisons of Figure 3 in which
the jet meridional widths change little with varying rotation periods under a
constant $T_{\rm star}$. The vertical structure of the eddies might be
somewhat sensitive to the heating and cooling rates which could affect their
interactions with the mean flow (Tsai et al., 2014). A key varying factor in
the vertical comparisons of Figure 3 may be the change of eddy structures
along with the change of effective temperature, and the wave speed $c$
associated with the eddies is no longer invariant. Holding the opacity
structure and rotation rate the same, Tan & Komacek (2019) showed that the
eastward equatorial jets near the photosphere become weaker and slightly
narrower as $T_{\rm eq}$ increases, and in addition, the upper layers
transition to westward zonal-mean equatorial flow; however, we don’t see this
behavior in our new sets of models with more realistic heating/cooling rates
and varying rotation. The sensitivity of the circulation and jet to heating
and cooling rates have been shown by Lee et al. (2021) and Steinrueck et al.
(2023). In addition to the radiative heating and cooling rates, heating rates
caused by ${\rm H_{2}}$-${\rm H}$ also change with equilibrium temperature,
further complicating the eddy structure. Overall, the jet structure
(determined by the wave-mean-flow interaction mechanism) requires a more
thorough understanding when more realistic factors are considered.

Figure 4: Left panel: temperature-pressure profiles from the sub-stellar
region (solid lines) and the anti-sub-stellar region (dotted lines) from
models with $T_{\rm eq}$ from 1800 K to 2600 K and without drag. Different
colors represent different planetary equilibrium temperatures. Profiles near
the sub-stellar area were averaged within $\pm 20^{\circ}$ in longitude and
latitude (where $0^{\circ}$ is the sub-stellar point); profiles near the anti-
sub-stellar point were averaged over longitudes >$160^{\circ}$ and
<$-160^{\circ}$ and latitudes within $\pm 20^{\circ}$. Right panel: profiles
averaged over the west limbs (solid lines) and east limbs (dashed lines).
Profiles of the west limbs were averaged over longitudes >$-100^{\circ}$ and
<$-80^{\circ}$ and latitudes within $\pm 75^{\circ}$, and profiles of the east
limbs were averaged over longitudes >$80^{\circ}$ and <$100^{\circ}$ and
latitudes within $\pm 75^{\circ}$. The “opening angles" in the transit
observations should depend on $T_{\rm eq}$ (Caldas et al., 2019; Wardenier et
al., 2022), and here we fix the sampling longitudinal regions of the limb
profiles merely for a demonstration of the heat transport effects. These
results were obtained by time-averaging over the last 200-day simulations. The
stellar effective temperature of these models is 6500 K. Figure 5: Radiative
heating and cooling profiles (solid lines) and ${\rm H_{2}}$-${\rm H}$ heating
and cooling profiles (dotted lines) as a function of pressure for the set of
models shown in Figure 4. Dashed lines are the ${\rm H_{2}}$-${\rm H}$ heating
minus the radiative heating rates, and correspond mostly to the dynamical
heating rates by the dry components. The left panel contains the dayside
profiles averaged over a domain that is $[\pm 60^{\circ},~{}\pm 60^{\circ}]$
of longitude and latitude around the sub-stellar point, and those in the right
panel are similarly averaged over a domain around the anti-sub-stellar point.
Figure 6:  Isobaric maps at about 6 mbar for atomic hydrogen mass mixing ratio
(left), radiative heating and cooling rates (middle) and ${\rm H_{2}}$-${\rm
H}$ heating and cooling rates (right) as a function of longitude and
latitudes. Streamlines with arrows on the left represent isobaric winds.
Results are from models with $T_{\rm eq}=1800$, 2200 and 2600 K from top to
bottom. The host star has an effective temperature of 6500 K.

####  3.1.2 Temperature and heating rate profiles

The left panel of Figure 4 shows the TP profiles averaged over a small area
around the sub-stellar point (solid lines) and around the anti-sub-stellar
point (dotted lines) for a set of simulations with $T_{\rm eq}=1800$ to 2600
K, a stellar temperature of 6500 K, and without drag. The dayside thermal
inversion above about 0.1 bar, i.e., temperature increases with increased
altitude, is primarily driven by absorption of TiO and VO first at relatively
high pressure and then by absorption of Fe at low pressures (Lothringer et
al., 2018; Parmentier et al., 2018; Hoeijmakers et al., 2018). Radiative
cooling is efficient on the nightside to cool the atmosphere and the TP
profiles show decreasing temperatures with decreased altitude on the nightside
photospheres. The overall nightside temperature variations among $T_{\rm
eq}=1800$ to 2200 K are small, and this has been discussed in Parmentier et
al. (2021) regarding radiative-dynamical timescale arguments. The nightside
temperatures with $T_{\rm eq}>2200$ K show stronger $T_{\rm eq}$ dependence,
presumably related to the ${\rm H_{2}}$-${\rm H}$ heating transport.

Dynamical heat transport induces strong thermal inversion and east-west
asymmetry near the terminators. The temperature profiles averaged in the east
and west limbs are shown in the right panel of Figure 4. Similar thermal
inversion in the upper atmosphere near the terminators has been seen in other
GCMs of ultra-hot Jupiters (Lee et al., 2022), and the effect here is perhaps
even stronger because of the additional ${\rm H_{2}}$-${\rm H}$ heat
transport. More interestingly, in the drag-free cases, the east-limb profiles
are systematically hotter than the west-limb profiles, at some pressure levels
showing east-west-limb temperature differences of more than 1000 K. The
overall east-west-limb temperature differences decrease with decreasing
$T_{\rm eq}$. In terms of pressure dependency, the east-west-limb differences
also decrease with decreasing pressures in the upper layers; this happens at
about $3.5\times 10^{-3}$ bar at $T_{\rm eq}=2600$ K and is at lower pressures
for cooler cases (about $2\times 10^{-4}$ bar at $T_{\rm eq}=1800$ K). The
pressure-dependent east-west-limb differences are partly related to the
relatively stronger day-to-night circulation at low pressures shown in Figures
1 and 2. The east-west limb differences should be readily detectable in the
transit lightcurve asymmetry by JWST as has been demonstrated by recent work
for WASP-39b (JWST Transiting Exoplanet Early Release Science Team, in prep.)
as well as high-resolution transmission spectroscopy (e.g., Savel et al.,
2023).

With a correlated-k radiative transfer scheme which calculates realistic
radiative heating and cooling rates and TP structure, we can investigate
pressure-dependent ${\rm H_{2}}$-${\rm H}$ heat transport, whereas we were not
confident to do so in previous work with an idealized, semi-grey scheme (Tan &
Komacek, 2019). The efficiency of hydrogen heat transport relative to the
nominal heat transport by the “dry" atmosphere can be shown by the ${\rm
H_{2}}$-${\rm H}$ heating rate profiles and compare to the radiative heating
profiles. Figure 5 shows the radiative heating profiles (solid lines) and
${\rm H_{2}}$-${\rm H}$ heating profiles (dotted lines) for the set of models
shown in Figure 4. The left panel contains the dayside profiles averaged over
a domain that is $[\pm 60^{\circ},~{}\pm 60^{\circ}]$ of longitude and
latitude around the sub-stellar point, and the right panel is similarly
averaged over a domain around the anti-sub-stellar point.

On the dayside, ${\rm H_{2}}$-${\rm H}$ cooling profiles are small at both
high and low pressures and show monotonic peaks at altitudes of about 10 mbar
to 0.4 mbar depending on the model $T_{\rm eq}$. The small ${\rm H_{2}}$-${\rm
H}$ cooling rates at high pressures are because of the globally small atomic
hydrogen fraction there. The small ${\rm H_{2}}$-${\rm H}$ rates at very low
pressure are because of the dissociation saturation, i.e., hydrogen is almost
fully dissociated over the whole dayside region where the heating profiles are
sampled, and there are little regions with large gradients of the atomic
hydrogen radio either in the horizontal or vertical directions. The ${\rm
H_{2}}$-${\rm H}$ heat transport depends on the gradient of the atomic
hydrogen mixing ratio rather than just the ratio itself, and so the ${\rm
H_{2}}$-${\rm H}$ heat transport is halted at low pressures. Only in the
intermedium pressure regions which are also where the thermal photospheres
roughly are, the transitions from molecular to atomic hydrogen lead to
significant ${\rm H_{2}}$-${\rm H}$ cooling that almost balances the radiative
heating. The radiative heating profiles on the dayside show corresponding
shapes to the hydrogen dissociation cooling profiles, but at low pressure,
they are still significantly nonzero and are balanced by the cooling induced
by the “dry" atmospheric heat transport. The “dry" heat transport cooling rate
profiles are represented as the dashed lines in Figure 5. Outside of regions
where the radiative heating rates peak and are balanced largely by ${\rm
H_{2}}$-${\rm H}$ cooling, the dry heating rates are the dominant balance to
radiative heating.

On the much cooler nightside, the atomic hydrogen mixing ratio is generally
smaller than on the dayside. For the hottest models, the ${\rm H_{2}}$-${\rm
H}$ heating rates near the thermal photosphere are comparable to that on the
dayside. All models show strong ${\rm H_{2}}$-${\rm H}$ heating at very low
pressures above the thermal photosphere (the right panel of Figure 5) but
their effects should be invisible in the (low-resolution) emission spectrum.
Likewise, the radiative cooling profiles have similar shapes to those of ${\rm
H_{2}}$-${\rm H}$ heating but with slightly larger amplitudes.

Figure 6 displays isobaric maps at about 6 mbar for atomic hydrogen mixing
ratio, radiative heating rates and H2-H heating rates of models with $T_{\rm
eq}=1800$, 2200 and 2600 K, for a closer look at the H mixing ratio and
heating rate structures. The shapes of the radiative heating rates, either as
shown in the mean vertical profiles in Figure 5 or in horizontal maps, display
a good correlation to the ${\rm H_{2}}$-${\rm H}$ heating rates but the
magnitude of the latter is always smaller. This is expected when the hydrogen
is near equilibrium. Neglecting time dependency and frictional dissipation of
heat, the thermodynamic equation can be written as (Equation [17] in Tan &
Komacek, 2019)

$(\mathbf{v}\cdot\nabla_{p}+\omega\frac{\partial}{\partial
p})(\bar{c}_{p}T+\mathcal{L}_{h}q)-\frac{\omega}{\bar{\rho}}=g\frac{\partial
F}{\partial p},$ (2)

where $\mathbf{v}$ is the isobaric wind vector, $p$ is pressure, $T$ is
temperature, $q$ is mass mixing ratio of atomic hydrogen H relative to total
mass, $\omega$ is the vertical velocity at pressure coordinate, $\nabla_{p}$
is the horizontal gradient at pressure coordinate, $\bar{c}_{p}$ is the mean
specific heat at constant pressure, $\mathcal{L}_{h}$ is the latent heat
associated with hydrogen dissociation, $\bar{\rho}$ is the mean gas density,
and $F$ is the radiative energy flux. If horizontal transport is the dominant
process, the above equation is simplified as

$\mathbf{v}\cdot(\bar{c}_{p}+\mathcal{L}_{h}\frac{\partial q}{\partial
T})\nabla_{p}T\sim g\frac{\partial F}{\partial p}.$ (3)

Heat transport is done by the “dry" and “latent heat" components which are
represented by the $\bar{c}_{p}$ and $\mathcal{L}_{h}\frac{\partial
q}{\partial T}$ (equivalently, the latter can be seen as a buffer in the heat
capacity, Roth et al., 2021). In nearly chemical equilibrium, i.e.,
$\frac{\partial q}{\partial T}\approx\frac{\partial q_{\rm eq}}{\partial T}$,
$\mathcal{L}_{h}\frac{\partial q_{\rm eq}}{\partial T}$ at the isobaric
surface is a function of temperature alone and is usually a factor of a few
larger than $\bar{c}_{p}$ at relatively high temperature before saturation.
Radiative heating is mostly balanced by the latent heating in these regions,
therefore showing good spatial correlations between radiative and latent
heating.

Figure 7: Thermal emission spectra (solid lines) of models with different
equilibrium temperatures and a stellar temperature of 6500 K on the dayside
(left) and on the nightside (right). The models shown here are the same as
those shown in Figure 4. The dotted lines are Blackbody spectra with total
energy corresponding to those of the model output spectra integrated over the
full wavelength range; they aid visualization of the spectral distributions.

####  3.1.3 Emission spectra

GCM outputs of the last 200 days were time-averaged and fed into the radiative
transfer code PICASO to post-process thermal emission spectra as a function of
phase as the planet rotates. Instantaneous time variability during the orbital
motion is therefore not considered in our post-processing.

Figure 7 shows the dayside (left) and nightside (right) thermal emission
spectra for the same suite of models with a stellar temperature of 6500 K
shown in Figure 4. Spectra on the dayside all show emission features that are
associated with the strong thermal inversions on the dayside, especially at
wavelength regions sculpted by the TiO, H2O, and CO molecules. The strength of
the features progressively decreases with increased temperature of the
atmosphere, primarily due to the thermal dissociation of the key molecules at
high altitudes and therefore the decreasing molecular scale heights probed by
the emission (e.g., Parmentier et al., 2018). Interestingly, the sharp and
narrow emission features at about 1.6 $\mu$m are caused by the CO molecules,
but this feature is only prominent at high temperatures. This is because CO
molecules are more resistant to higher temperatures compared with H2O, and
its’ spectral features at 1.6 $\mu$m appears more robustly when H2O becomes
significantly dissociated at high temperatures ($\gtrapprox 3000$ K). On the
nightside, the near-IR and IR spectral features change to absorption because
of the decreasing temperature with increasing altitude (see Figure 4). The
feature strength also slightly decreases with increased temperature but is not
as obvious as that on the dayside. The nightside spectra exhibit significant
excesses of fluxes in the visible wavelength compared to the Blackbody
spectra, and those are emitted mostly from the terminator regions when the
planetary nightsides face toward the observer according to the spatially
resolved flux map (not shown). The emitting agents at the visible wavelengths
are primarily TiO and VO, similar to the dayside, and those fluxes trace the
strong heated thermal layers near the terminators (right panel of Figure 4).
The nightside emission features of TiO and VO should be readily detectable in
the short wavelength range by JWST, and the near-limb heat transport could be
tested.

Thermal spectra of models with stellar temperatures of 5500 K and 6000 K have
qualitatively identical features except that they generally have higher
dayside emission and lower nightside emission. This is due to the poor day-
night heat transport when the rotation rate increases, which is more visible
in the phase curves discussed below.

####  3.1.4 Cases with a strong drag

Similar to previous GCM parametric surveys of Perez-Becker & Showman (2013);
Komacek & Showman (2016); Komacek et al. (2017), we performed the same set of
models as those partly shown in Figures 1 and 2 but with a strong drag of
$\tau_{\mathrm{drag}}=10^{4}$ s. Including the drag is motivated by the
presence of possible magnetohydrodynamic dissipation (e.g., Perna et al.,
2010; Rogers & Komacek, 2014; Beltz et al., 2022a) or mixing by small-scale
turbulence that cannot be resolved in global models (e.g., Li & Goodman,
2010). The simple drag scheme crudely mimics the dissipation of kinetic energy
of the flow over a global scale in both the zonal and meridional directions.
The strong drag typically reduces the wind speed, damps out the superrotating
jet and the global-scale Matsuno-Gill wave pattern, and increases the day-to-
night temperature variations. The horizontal pattern is instead mostly
characterized by the day-to-night flow, except for some rapidly rotating cases
in which the magnitude of the Coriolis force is comparable to the drag force
and global wave patterns still appear (Showman et al., 2013a). Therefore,
compared to the strong-drag cases in Tan & Komacek (2019) in which the primary
model sets assumed a fixed 2.43-day rotation period, the relatively rapid
rotating, strong-drag cases in this study exhibit more obvious wave patterns.
For more discussion of the circulation and temperature maps, readers are
referred to Tan & Komacek (2019) for results relevant to ultra-hot Jupiters
and we do not show them here.

On the dayside, as far as the vertical temperature structures, heating
profiles, and the resulting thermal spectra, those of the strong-drag models
are qualitatively similar to those of drag-free models (for temperature
structures, see the left panel of Figure 4 and the upper left panel of Figure
9). The primary quantitative difference is that the strong-drag models have
greater dayside fluxes and reduced nightside fluxes due to their poorer day-
to-night heat transport; the spectral features are almost identical otherwise.
This can be seen in the comparison in the left panel of Figure 8 for models
with or without drag, both of which have $T_{\rm eq}=2600$ K and a stellar
temperature of 6500 K.

On the nightside, the total spectral energy flux is significantly reduced in
the strong-drag models, as expected from the poorer heat transport. However,
in the case shown in the right panel of Figure 8, the nightside comparison
becomes much more subtle than the dayside comparison, showing much more muted
absorption features in the strong-drag case in near-IR and IR wavelengths; and
some regions seem to be blended with a weak emission feature within a broader
absorption feature, for example, around the water 1.4 $\mu$m window; as well
as comparable or even greater flux at short wavelengths $<0.7~{}\mu$m. The
nightside “nonstandard" near-IR and IR spectral properties of the strong-drag
model come from two distinct contributing areas, 1) near-limb regions where
there are strong thermal inversions and fluxes are generally emitted at low
pressures between 1 to 10 mbar, and 2) the broad area around the anti-
substellar point where there is no inversion and the photospheric pressures
are higher (pressure slightly less than 1 bar). Whereas in the drag-free
model, the nightside fluxes at near-IR and IR are primarily from the broad
anti-substellar areas (with much less flux contributed from the limb, low-
pressure regions) and therefore reflect the non-inverted temperature
structures of that region. In the visible wavelengths, both models probe the
low-pressure limb regions, and so they show similar flux. The above findings
were from constructing the 3D contribution functions by perturbing the 3D
temperature field and rerunning the radiative transfer post-processing every
time a local temperature is perturbed. Our results highlight that the
nightside spectrum of some UHJs may show exotic features depending on the
detailed circulation pattern, and JWST should be able to detect these possible
features (Mikal-Evans et al., 2023).

Figure 8: Modeled thermal emission spectra from the dayside (left panel) and
from the nightside (right panel). Red lines are from the model of a planetary
$T_{\rm eq}=2600$ K, a stellar temperature of 6500 K and without drag; blue
lines are from a similar model except for a strong drag with a timescale
$\tau_{\mathrm{drag}}=10^{4}$ s.

####  3.1.5 Cases without TiO, VO, and Fe

Figure 9: Temperature-pressure profiles on the dayside (upper left panel) and
nightside (upper right panel) for cases either including (solid lines) or
excluding TiO, VO, and Fe (dotted lines). Corresponding thermal spectra on the
dayside and nightside are shown in the bottom panels. All GCMs have a drag
timescale of $\tau_{\mathrm{drag}}=10^{4}$ s and a stellar temperature of 6500
K.

Optical absorbers such as TiO and VO have long been known to significantly
affect the atmospheric structure, chemistry, and dynamics of hot Jupiters
(e.g., Hubeny et al., 2003; Fortney et al., 2008; Showman et al., 2009). The
search for these species using emission and transmission spectroscopy so far
has not yielded a common existence of TiO and VO in the atmospheres of hot
Jupiters, potentially related to the cold trap of these species (Spiegel et
al., 2009; Parmentier et al., 2013). More recently, additional metal species
along with wind velocities have been found to exist in the atmospheres of some
ultra-hot Jupiters via high-resolution spectroscopy (e.g., Yan et al., 2019;
Gibson et al., 2020; Hoeijmakers et al., 2020; Pai Asnodkar et al., 2022) and
should contribute to significant heating in the upper atmospheres (e.g.,
Lothringer et al., 2018). Fe is one of the most important species among them,
which is also included in our main suite of GCMs. Similar to TiO and VO, Fe is
also condensible in atmospheres of hot Jupiters (Helling, 2019; Gao et al.,
2021) and might be subjected to the cold trap process. Notably, observations
of the ultra-hot Jupiter WASP76b indicate the partial rain-out of Fe
(Ehrenreich et al., 2020), but the transit asymmetry in Fe absorption could
instead be due to clouds or scale height variations (Savel et al., 2022;
Wardenier et al., 2021).

Here we explore the effects of excluding the strong optical absorbers on the
thermal structure, heat transport, and resulting thermal spectra of ultra-hot
Jupiters. Several species are excluded (see Table 1), most importantly TiO, VO
and Fe. For this test, we include a strong drag with
$\tau_{\mathrm{drag}}=10^{4}$ s to accelerate the convergence of the models
and study a limiting case with reduced day-night heat transport. Figure 9
summarizes our findings. The dayside and nightside temperature profiles are
shown in the top row. Near the infrared photosphere on the dayside, there is
no strong thermal inversion in cases without strong optical absorbers except
in the hottest case in which water cooling is weak. Correspondingly, the
dayside thermal emission spectra of these cases (dotted lines in the lower
left panel) display either a near blackbody in the hotter cases or absorption
features in the molecular bands in cooler cases, in contrast to those of the
canonical suite of models that show strong emission features (solid lines).
Because of the weak absorption in the upper atmosphere, more stellar flux can
be deposited in more depths, leading to hotter infrared photospheres in models
without TiO, VO and Fe due to the stronger infrared greenhouse (Guillot,
2010). For more in-depth discussion of the impact of metallic absorbers on the
temperature structure and emission spectra of ultra-hot Jupiters, see
Parmentier et al. (2018) and Lothringer et al. (2018).

The effect on the nightside is also profound. The nightside thermal
photospheres of models without TiO, VO and Fe are also hotter due to more heat
being transported from their hotter dayside thermal photosphere, which is
consistent with the finding in the pseudo-2D framework of Roth et al. (2021).
The nightside spectra shown in the lower right panel of Figure 9, as well as
those shown in Figures 7 and 8, demonstrate that even on the nightside, strong
TiO and VO spectral features can be readily detectable in models with TiO and
VO included. Interestingly, cases with $T_{\rm eq}=2400,2600$ K have obviously
hotter upper atmospheres than those including strong optical absorbers. The
optical absorbers of TiO, VO and Fe act as cooling agents on the nightside
which would in principle radiate away heat transport if they were included,
but in models without those absorbers the upper-layer temperature remains high
at low pressures due to the efficient heat transport by hydrogen recombination
and inefficient thermal cooling. The lack of thermal cooling by TiO and VO can
be seen in the lower right panel in which strong emission is absent in the
optical wavelengths in cases without TiO, VO and Fe. We also post-process the
GCM outputs without TiO, VO and Fe using k-coefficient tables that do include
those species, and excessive flux appears on the optical wavelength,
indicating that if there were those species, the upper atmosphere should be
effectively cooled down. High-resolution spectroscopy observations might be
able to measure the nightside-to-limb upper thermal structures and provide
constraints on the cooling role of these species.

####  3.1.6 Circulation influenced purely by different stellar spectra

Figure 10:  Temperature-pressure profiles on the dayside (solid lines) and
nightside (dotted lines) for models with $T_{\rm eq}=2400$ K and all other
parameters the same, including the same rotation period of 1.3 days, but only
the stellar effective temperature differs from 5500 K to 6500 K. These
profiles are from spatial averaging of the same manner as those shown in
Figure 4. All GCMs have a drag timescale of $\tau_{\mathrm{drag}}=10^{4}$ s.

The temperature structures of ultra-hot Jupiters could be sculpted by
different stellar spectra as modeled by static 1D models of Lothringer &
Barman (2019), but their role in circulation is unclear. Results in our main
suites of models are mixed with effect of varying rotation rates. To isolate
the rotational effect, we also performed a few additional models with $T_{\rm
eq}=2400$ K and all other parameters the same, including the same rotation
period of 1.3 days, but only the stellar effective temperature differs from
5500 K to 6500 K. The temperature profiles are shown in Figure 10. On the
dayside, the upper atmosphere is hotter and the infrared photosphere is cooler
when the stellar temperature is higher. The nightside temperature structures
show the same trend. This trend is in line with Lothringer & Barman (2019),
and is because hotter stars have a higher fraction of the stellar energy in
the shorter wavelengths that could be absorbed in the upper atmosphere by
metals. However, the magnitude of the differences, especially in the thermal
photosphere, is not significant in our work. This is not surprising partly
because the stellar temperature variations in this work are far less extreme
than those explored in Lothringer & Barman (2019) with stellar effective
temperatures from 5700 K to 10500 K. Dynamics might play a secondary role in
regulating this difference by heat transport. Because of the relatively small
effect on the thermal structure, the overall circulation is barely affected by
the different stellar temperatures when other parameters are held fixed.
Another difference is that models used in Lothringer & Barman (2019) reach
shorter wavelengths and capture more UV fluxes than our GCM and PICASO whose
shortest wavelength reaches about 0.26 $\mu$m. The UV absorption is important
for the hottest stars in the UHJ samples and drives stronger thermal inversion
at very low pressures, but this difference is expected to be minor here
because of the less extreme stellar temperatures we considered here and that
the low-resolution emission spectrum is generally less sensitive to structures
at very low pressures.

###  3.2 Comparisons to observations: dayside brightness temperatures and
phase-curve amplitudes

Figure 11: Dayside brightness temperature minus the planetary equilibrium
temperature as a function of planetary equilibrium temperature. Dots with
error bars are Spitzer data presented in Deming et al. (2023) with different
panels representing analysis at different Spitzer bandpasses, and HST WFC3
data collected in Mansfield et al. (2021). Lines are from different sets of
models: black, red and blue lines are from models with stellar temperatures of
5500 K, 6000 K and 6500 K, respectively; solid lines are drag-free and dashed
lines are with $\tau_{\mathrm{drag}}=10^{4}$ s; the dotted line is from models
with $\tau_{\mathrm{drag}}=10^{4}$ s, the stellar temperature of 6500 K and
without TiO, VO and Fe. Figure 12: Day-to-night brightness temperature
differences as a function of equilibrium temperature. Lines are from different
sets of models that are described in Figure 11, and dots with error bars are
observations described below. They are separately shown in different
instrument bandpasses as indicated above each panel. Data are updated from
those summarized in Tan & Komacek (2019) and Wong et al. (2021). At Spitzer
4.5 $\mu$m, most data were from analysis by Bell et al. (2021), WASP-76b is
from May et al. (2021), HD 149026b is from Zhang et al. (2018). At Spitzer 3.6
$\mu$m, data are from Wong et al. (2015, 2016), Maxted et al. (2013), Cowan et
al. (2012), Beatty et al. (2019) May et al. (2021) and Kreidberg et al.
(2018). HST WFC3 data are from Kreidberg et al. (2018); Arcangeli et al.
(2019); Mikal-Evans et al. (2022). TESS observations are obtained from Mancini
et al. (2022), Owens et al. (2021), Addison et al. (2021), Jansen & Kipping
(2020), Von Essen et al. (2020), Daylan et al. (2021), and Eftekhar & Adibi
(2022). For WASP-121b TESS phase curve, Bourrier et al. (2020) and Eftekhar
(2022) also show roughly consistent dayside and nightside brightness
temperatures as Daylan et al. (2021).  Figure 13: Normalized day-to-night
brightness temperature variations at bandpasses of different instruments. The
black solid lines are from a set of simulations with a fixed rotation period
of 2.24 days and different $T_{\rm eq}$ around a 6500 K star. The dotted lines
are models with the same parameters except that the rotation periods are
[3.07, 2.24, 1.68, 1.3, 1.06] days for $T_{\rm eq}=$ [1800, 2000, 2200, 2400,
2600] K. The thick grey lines are normalized day-to-night temperature
variations from a set of semi-grey GCMs with ${\rm H_{2}}$-H conversion and
$\tau_{\mathrm{drag}}=10^{7}$ s in Tan & Komacek (2019) (TK2019). TK2019
models assumed a fixed rotation period of 2.43 days, and these results were
derived from the bolometric fluxes and do not distinguish between different
instrument bandpasses (grey lines are identical in all panels). The TK2019
results are included for comparisons between GCMs with realistic non-grey and
idealized semi-grey radiative transfer schemes.  Figure 14: Phase offsets of
drag-free models with different equilibrium temperatures for Spitzer 4.5$\mu$m
and comparisons to those presented in the uniform analysis of Bell et al.
(2021). Different lines represent different stellar effective temperatures.

Using post-processed phase-resolved emission spectra by PICASO, we are able to
directly compare the model-predicted dayside emission and phase-curve
amplitudes with observed values over the population of ultra-hot Jupiters.

We first compare the dayside photometric brightness temperature over the
Spitzer  3.6, 4.5 and 8 $\mu$m bands as well as the HST WFC3 band in Figure
11. Dots with uncertainties in the three Spitzer bands are from the uniform
reanalysis of Deming et al. (2023), which is based on a similar analysis of
Garhart et al. (2020) but with more updated targets. Data in the HST WFC3 band
are summarized in Mansfield et al. (2021). The observational data sets extend
to lower effective temperatures but here we only show those above about
$T_{\rm eq}=1700$ K for comparisons with our model grids. Lines in the plots
represent model outputs, including cases with three stellar temperatures,
drag-free, and $\tau_{\mathrm{drag}}=10^{4}$ s, as well as one set that
excludes TiO, VO and Fe. The vertical axis is the brightness temperature at
that bandpass minus the $T_{\rm eq}$ of the planets.

Regarding just the model-predicted dayside brightness temperatures, there is a
general trend of rising brightness temperature with increasing $T_{\rm eq}$ at
all wavelengths and all model assumptions. Models with a strong drag of
$\tau_{\mathrm{drag}}=10^{4}$ s have overall higher brightness temperatures
than drag-free ones due to the inefficient day-to-night heat transport. At a
given drag setup, models with a cooler stellar temperature exhibit higher
dayside brightness temperature because these planets rotate faster which
inhibits day-to-night heat transport more efficiently (Tan & Showman, 2020).
In Section 3.1.6, we have discussed that results are quantitatively similar
among models with different stellar types but all other parameters are the
same. Here we demonstrate that with a reasonable stellar type range, the
effects of rotation dominate over the effects of the different stellar spectra
in the circulation, emission temperature, and day-night flux differences. The
dayside brightness data (dots in Figure 11) show a general agreement to models
with a wide range of parameters in terms of the trends with increasing
equilibrium temperature; however, the scatter of data is much wider than the
model predictions, indicating a much wider intrinsic parameter scattering
among the hot Jupiter population or that our model is missing some other key
processes such as the formation of clouds and their effect on the heat
transport (e.g., Roman et al., 2021; Parmentier et al., 2021).

The dayside brightness temperature at a certain wavelength of hot Jupiters is
affected by both the radiative and dynamical processes, and the question is
which process is more important at a certain wavelength. Figure 16 shows the
same model predictions as Figure 11 but without observational data for better
visualization. We compare the dayside brightness temperature variations caused
by different orbital and drag parameters (but with the same atmospheric
compositions) to those caused by different atmospheric compositions (but the
same orbital and drag parameters), and find that at relatively shorter
wavelengths, for example, the TESS and HST WFC3 bandpasses, day-night heat
transport is as important as the radiative processes. But towards longer
wavelengths from Spitzer 3.6 to 8 $\mu m$, the radiative process becomes more
and more important than the day-night dynamical heat transport.

Next, we turn our attention to the day-to-night variations of brightness
temperatures. Figure 12 shows the normalized day-to-night brightness
temperature variation for all models, together are observed values for ultra-
hot Jupiters at various bandpasses. Lines are modeled outputs and are labeled
the same as those in Figure 11. The peak and trough thermal flux over
individual bands are converted to brightness temperatures, which are then
normalized to get those values. The day-night brightness temperature
variations are less sensitive to whether there are strong visible absorbers,
in contrast to the dayside brightness temperature shown in Figure 11. Given a
fixed host star, the fractional day-night brightness temperature variations at
different wavelengths all exhibit either a nearly constant or slight increase
with increasing equilibrium temperature, which is in contrast to analytic
analysis of Komacek & Tan (2018) and GCMs with simplified radiative transfer
(Tan & Komacek, 2019), in which the fractional variations first increase but
then decrease with increasing equilibrium temperature. The previous work
assumed a fixed rotation period of slightly more than 2 days for models with
different equilibrium temperatures. This difference is due to the increased
rotation rate with increasing equilibrium temperature around a given star in
this work, and the increasing rotation rate tends to reduce the day-night heat
transport efficiency as discussed above.

We can notice a few obvious points from the data-model comparisons. First, the
range of our model predictions is still less than the observed values of the
day-night variations especially in the Spitzer 4.5 $\mu$m band and the TESS
band, although our models have covered a wide range of drag setup, rotation
periods and stellar properties. Second, in the TESS band, there seems to be a
systematic overestimate of the day-night temperature variation compared to the
observed values. This is a tougher issue than increasing the day-night
temperature contrasts because the latter can be potentially achieved in
several self-consistent ways such as nightside clouds and strong
magnetohydrodynamic drag. In contrast, some of our GCMs are nearly drag-free
(apart from numerical dissipation which is necessary to maintain numerical
stability) and include hydrogen dissociation and recombination which are
already prompted for the most efficient day-night heat transport. One remote
possibility is that the planets have enormous internal heat flux (even much
larger than those estimated by Thorngren et al., 2019) that is comparable to
the incoming stellar flux due to unknown processes, but it is implausible to
argue that the TESS observed population all have that amount of internal heat
flux. The scenario of large internal heat flux is also inconsistent with their
present-day radius. The systematic overestimation of the modeled day-night
temperature differences in the TESS band motivates us to look for alternative
solutions to enhance heat transport in the Appendix B. Third, for objects with
phase-curve observations over multiple wavelengths, no single model setup can
simultaneously explain all the observations. One notable example is WASP-103b,
for which the Spitzer 4.5 $\mu$m phase-curve amplitude is much larger than our
grid; the Spitzer 3.6 $\mu$m one lies in the $\tau_{\mathrm{drag}}=10^{4}$ s
branch of models; lastly the HST WFC3 one is in the drag-free branch of
models. Others include WASP-18b, WASP-14b, WASP-19b, WASP-121b, and WASP-33b.

We conducted an additional set of simulations with a fixed rotation period of
2.24 days and different $T_{\rm eq}$ around a 6500 K star. Results are shown
as the black solid lines in Figure 13. Dashed lines are included for
comparisons and they are from models with the same parameters except that the
rotation periods are [3.07, 2.24, 1.68, 1.3, 1.06] days for $T_{\rm eq}=$
[1800, 2000, 2200, 2400, 2600] K as described in Section 3.1. The comparisons
between the two sets of models immediately reveal the importance of rotation
rate in controlling the day-night heat transport. Interestingly, at Spitzer
3.6 and 4.5 $\mu$m bands, the fixed-rotation-rate models show a fractional
variation turning point at about 2200 K which is in good agreement with
previous work of Komacek & Tan (2018) and Tan & Komacek (2019) that used more
idealized models; however, the fractional variations at TESS and HST WFC3
bandpasses keep increasing with increasing temperature, although the trend
slows down after $T_{\rm eq}=2200$ K, unlike the IR bands and those predicted
in Komacek & Tan (2018) and Tan & Komacek (2019). Finally, results from a set
of semi-grey models with ${\rm H_{2}}$-H conversion,
$\tau_{\mathrm{drag}}=10^{7}$ s and a fixed rotation period of 2.43 days in
Tan & Komacek (2019) are included as thick grey lines in Figure 13 for
comparisons between realistic and simplified models. The qualitative trend of
the semi-grey curve agrees with the non-grey results in the Spitzer bands when
rotation periods are held fixed (grey lines vs the black solid lines) to some
extent. However, the non-grey results exhibit both qualitative and
quantitative differences from the semi-grey results through the influences of
the thermal structure by non-grey heating/cooling rates and also through non-
grey opacity (meaning that different wavelengths probe different atmospheric
layers). For example, the trends are qualitatively different in the WFC3 and
TESS bands with significant offsets between different modeling frameworks, and
the absolute values differ obviously in the Spitzer  bands. These again
highlight the wavelength-dependent nature of the thermal phase-curve
amplitudes for the hot Jupiter population (Parmentier et al., 2021).

Lastly, we compare the model-predicted phase curve offsets with those from the
reanalysis of Spitzer  4.5 $\mu$m phase curves (Bell et al., 2021) in Figure
14. Only models with no drag are included; the phase curve offsets of models
with $\tau_{\mathrm{drag}}=10^{4}$ s are generally very small (but are still
positive) due to the suppression of jets and waves. The modeled offsets are in
general smaller than $\sim 22^{\circ}$ and larger than $\sim 15^{\circ}$
(positive values indicate the eastward displacements of the hot regions
relative to the sub-stellar point). This narrow variation is in large contrast
to the observed values that range from almost $-60^{\circ}$ to $30^{\circ}$.
In addition, many observed offsets are negative which is not predicted by many
GCMs of hot Jupiters (see a review of Showman et al., 2020). The factors
affecting the thermal phase-curve offsets include magnetohydrodynamics (e.g.,
Rogers, 2017) and clouds (e.g., Parmentier et al., 2021; Roman et al., 2021),
for which we do not include in our models. Given the incompleted physics in
our models and the potential sensitivity to data analyzing methods, we may say
that the drag-free and cloud-free models are likely insufficient to explain
the observed phase-curve offsets. Lastly, the phase offsets can be sensitive
to different data processing methods, for example, some data collected in Tan
& Komacek (2019) and presented in May et al. (2022) are different from those
reported in Bell et al. (2021) for the same objects. Even in Bell et al.
(2021), for many planets with negative phase offsets, the phase offsets could
be either positive or negative depending on the data reduction methods (see
their Figure 5). JWST will help to break the current conflicting phase-curve
observations of some targets (e.g, WASP-121b) with exquisite precision and
(simultaneously) large wavelength coverage for hot Jupiter atmospheres (Mikal-
Evans et al., 2023).

###  3.3 Variability

Searching for time variability of hot Jupiters’ atmospheres is a timely topic,
helping to constrain the dynamical mechanisms affecting the meteorology of
these objects (e.g., Armstrong et al., 2016; Jackson et al., 2019; Wong et
al., 2022). Komacek & Showman (2019) investigated the variability of hot
Jupiters using a suite of GCMs and found that about 2% variability can be
expected but the origin of the variability was mostly confined near the
equator. Tan & Komacek (2019) and Tan & Showman (2020) similarly found
variability in GCMs with short rotation periods, and the mechanism is similar
to those shown here. In Tan & Showman (2020), the zonal-mean meridional
potential vorticity is shown to have reversed signs in the vertical direction
in relatively rapid rotators of tidally locked hot Jupiters and the baroclinic
instability is the most likely mechanism here. Some GCM studies with ultra-
high spatial resolution and low numerical dissipation showed that significant
evolving storms may occur even under numerical setups for canonical hot
Jupiters (e.g., Skinner & Cho, 2022), in contrast to many other GCMs (Showman
et al., 2020).

Figure 15: Variability of the planet-to-star emission ratio relative to the
time-mean during the output time from the model with an equilibrium
temperature of 2000 K around the 5500 K star. Blue lines are from the dayside
spectra and red lines are from the nightside spectra. The flux ratios through
different bandpasses are shown in different panels. The left vertical axis of
each panel is variability measured by part per million, and the right vertical
axis represents the fractional variation relative to the time mean.

Our drag-free simulations show variability in their dayside and nightside
temperatures. Figure 15 shows the dayside and nightside emission as a function
of time from the case with $T_{\rm eq}=2000$ K around the 5500 K star. Long-
term variation modulations with timescales of tens of days dominate the
variability, and variability with much lower amplitudes also presents over the
orbital timescale. The amplitude of the hemispherically integral emission can
reach over 60 ppm at wavelengths near 3.6 $\mu$m and over 40 ppm at
wavelengths near 4.4 $\mu$m; but the variability amplitudes are relatively
small at visible and near-IR wavelengths where they probe higher-temperature
layers and thus the fractional flux perturbations caused by certain
temperature variations are smaller. A visual inspection of the temporal
evolution of the global temperature maps of these simulations suggests that
the long-term variation is driven mostly by off-equatorial eddies resulting
from atmospheric instabilities on the dayside and then subsequently migrating
longitudinally to the nightside and dissipated. Although the fractional
variability is similar at different bandpasses, their amplitudes at planet-to-
star flux ratio differ simply because of the different basic planet-to-star
ratios at different wavelengths. Therefore, the emission variability on the
dayside or nightside would be best tested in infrared wavelengths compared to
optical or near-IR wavelengths. JWST has yielded high-precision secondary
eclipse observations of (ultra) hot Jupiters’ atmospheres Coulombe et al.
(2023); the currently predicted variability amplitudes of over 60 ppm at about
3.6 $\mu$m wavelengths on the nightside and about 30 ppm on the dayside should
be within the limit of JWST single-eclipse precision. For optimal objects and
at infrared wavelengths, JWST should be able to put constraints on the
meteorological variability mechanisms proposed on hot Jupiters.

##  4 Summary and prospects

###  4.1 Summary

The day-night temperature variations inferred from phase-curve observations of
ultra-hot Jupiters provide rich information about the mechanisms controlling
the atmospheric circulation of these planets, including radiative heating and
cooling, dynamics of jets and waves, heat transport aided by hydrogen
dissociation and recombination, cloud formation and their radiative feedback,
and magnetohydrodynamics. In this study, we have upgraded our GCM for ultra-
hot Jupiters (Tan & Komacek, 2019; Komacek et al., 2022a) to couple the
dynamical core to a non-grey correlated-k radiative transfer scheme. This
allows us to capture the more realistic atmospheric thermal structure and
radiative heating and cooling rates. Using this model, we explore a wide range
of parameter space including different planetary equilibrium temperatures from
$T_{\rm eq}=$1800 to 2600 K, different tidally synchronized planetary rotation
periods associated with combinations of planetary $T_{\rm eq}$ and varying
stellar effective temperatures from 5500 K to 6500 K, and different frictional
drag setup in Section 3. To summarize our results in Section 3, we find that:

  * •

Assuming a fixed solar composition of the atmosphere, the atmospheric dynamics
of ultra-hot Jupiters depends sensitively on equilibrium temperature and host
star type, and thus rotation rate. For a fixed equilibrium temperature, hot
Jupiters orbiting earlier-type host stars have smaller day-night temperature
contrasts due to their longer rotation periods assuming spin-synchronization.
We also find that ultra-hot Jupiters with larger equilibrium temperatures
orbiting a fixed host star type have greater day-night temperature contrasts
and thus a faster speed of the equatorial jet. In agreement with our previous
double-grey models, we find that hotter ultra-hot Jupiters will have smaller
hot spot offsets. We also similarly find that strong drag further increases
day-night temperature differences and reduces hot spot offsets.

  * •

Our non-grey GCMs of ultra-hot Jupiters across a range of $T_{\rm eq}$
ubiquitously have strong dayside thermal inversions. In cases without strong
drag, these thermal inversions persist to both the east and west limbs, with
west limbs typically having higher-altitude thermal inversions than the east
limbs in these cloud-free simulations. The dayside inversions are driven by
shortwave absorption, which is largely balanced by cooling from H2
dissociation on the dayside – conversely, H2 recombination is the dominant
balance of radiative cooling on the nightside. By post-processing these models
to predict emergent spectra, we confirm with 3D non-grey GCMs the expectation
from 1D models that ultra-hot Jupiters will have dayside emission features
with emission strength decreases with increasing equilibrium temperature. Our
cloud-free, drag-free simulations also generally display nightside absorption
features, as nightside thermal inversions are limited to the upper atmosphere
in our hottest cases. The dayside spectral properties are not sensitive to our
model drag strength, but the nightside spectral shape could show significant
differences between the drag-free and strong-drag cases.

  * •

Our GCM simulations of ultra-hot Jupiters that include the thermodynamic
effect of hydrogen dissociation and recombination coupled to a non-grey
radiative transfer scheme predict a weaker temperature-dependence of dayside
temperatures and day-night temperature contrasts than previous non-grey
models. We find that given the orders of magnitude uncertainty in
$\tau_{\mathrm{drag}}$, our non-grey GCMs can broadly reproduce the observed
dayside temperatures and day-night temperature contrasts from Spitzer and HST.
In general, we require strong drag to explain observed ultra-hot Jupiter day-
night temperature contrasts at Spitzer Channel 1 and in the HST/WFC3 bandpass.
In tandem, the phase curve offsets of many ultra-hot Jupiters are near zero or
even negative (i.e., implying that peak flux occurs after secondary eclipse).
This may imply that MHD effects shape the observable properties of ultra-hot
Jupiters. Weather TiO and VO are included has a strong impact on the dayside
brightness temperatures especially for those at longer IR wavelengths. Besides
the dayside spectra, our modeled nightside spectra also displace emission
features of TiO and VO at visible wavelengths and they originate from
relatively low-pressure layers near the limb regions, where strong thermal
inversions occur and penetrate towards the nightside due to the day-to-night
heat transport.

In addition, for one particular test case at a given planetary $T_{\rm eq}$,
we explore the effects of varying atmospheric metallicity and carbon-to-oxygen
ratio in Section B, finding that:

  * •

At a given set of system parameters appropriate for TOI 1431b, we investigated
the effects of varying metallicity from 0.1 to 3 of the solar value, and C/O
ratio from 0.5 to 2.5 of the solar value, on the circulation and day-night
temperature variation in the TESS band. Reducing metallicity reduces the day-
night temperature variation and this trend is robust over different C/O
ratios. At a given metallicity, varying C/O ratio changes the day-night
temperature variation but the trend is less clear. The basic idea is the
efficiency of radiative cooling decreases with decreasing metallicity and in
some cases also the C/O ratio. Despite a large variation over this parameter
space, the lowest day-night temperature variation in the TESS band in our
models is still significantly greater than that observed. Further development
of physical processes in the GCM and/or additional observations should be
conducted to determine the source of this tension.

###  4.2 Prospects and future work

Our results serve as a baseline study for comparisons with those including
clouds and varying parameters such as atmospheric composition and gravity.
Clouds appear to be ubiquitous in atmospheres of close-in exoplanets. Their
effects on transmission, reflection, and thermal emission have been explored
in various contexts but their dynamic feedback in the general circulation of
close-in exoplanets is comparatively less understood. If clouds show a
significant non-grey effect on the thermal emission, phase curves at different
wavelengths are essential to test the hypothesis that clouds are present
(Parmentier et al., 2021). In parallel, mechanisms of cloud feedback in self-
luminous planetary atmospheres have been shown to be important to drive
variability and circulation of brown dwarfs and isolated giant planets (Tan &
Showman, 2021). Despite that energy of the circulation in hot Jupiter’s
atmospheres is injected on the largest scale, we aim to explore whether
mechanisms in self-luminous atmospheres are partly responsible for shaping the
cloud structure and circulation. In Komacek et al. (2022a), we started to
incorporate the cloud radiative feedback in a hot Jupiter GCM using a semi-
grey framework and showed that cloud patchiness could feed back on the
dynamics and thermal structure. It is important to test the results in a more
realistic (e.g., non-gray) modeling framework and over a wider planetary
parameter space in future work.

Another great challenge is the representation of MHD on hot close-in gas giant
planets. Self-consistent MHD models are required to fully understand the
dynamics of hot and ultra-hot Jupiters (Rogers & Komacek, 2014; Rogers, 2017).
However, there have been relatively few non-ideal and self-consistent MHD
models in the field to push the efforts further and confront observations,
especially those of westward hot spot offsets. A more commonly seen approach
is to parameterize the effects of (assumed) magnetic fields on the atmospheric
flow in hydrodynamic models (Perna et al., 2010; Rauscher & Menou, 2013; Beltz
et al., 2022a). In the future, adopting such an “active” kinematic magnetic
drag scheme in our model could be an immediate next step to test our GCMs
against observations. A long-term plan would be to move our current schemes of
radiative transfer and thermodynamics treatments of hydrogen dissociation and
recombination to self-consistent MHD models.

## Acknowledgements

X.T. and R.T.P. acknowledge support from the European community through the
ERC advanced grant EXOCONDENSE (PI: R.T. Pierrehumbert). R. L. acknowledges
support from the NASA ROSES XRP program with the grant 80NSSC22K0953, and from
STScI grants JWST-AR-01977.007-A and JWST-AR-02232.008-A. The authors would
like to acknowledge the use of the University of Oxford Advanced Research
Computing (ARC) facility in carrying out this work
(http://dx.doi.org/10.5281/zenodo.22558). Part of this work was completed with
resources provided by the University of Chicago Research Computing Center.

## Data availability

The data underlying this article, including GCM outputs and post-processed
phase-dependent spectra, is publicly available in Zenodo
https://doi.org/10.5281/zenodo.10121933. The 11-bin and 622-bin k-coefficient
tables used in this work can be found in Lupu et al. (2021, 2022).

## References

  * Adams et al. (2022) Adams D. J., Kataria T., Batalha N. E., Gao P., Knutson H. A., 2022, ApJ, 926, 157
  * Addison et al. (2021) Addison B. C., et al., 2021, The Astronomical Journal, 162, 292
  * Arcangeli et al. (2019) Arcangeli J., et al., 2019, A&A, 625, A136
  * Armstrong et al. (2016) Armstrong D. J., de Mooij E., Barstow J., Osborn H. P., Blake J., Saniee N. F., 2016, Nature Astronomy, 1, 1
  * Baeyens et al. (2021) Baeyens R., Decin L., Carone L., Venot O., Agúndez M., Mollière P., 2021, Monthly Notices of the Royal Astronomical Society, 505, 5603
  * Batalha et al. (2019) Batalha N. E., Marley M. S., Lewis N. K., Fortney J. J., 2019, ApJ, 878, 70
  * Beatty et al. (2019) Beatty T. G., Marley M. S., Gaudi B. S., Colón K. D., Fortney J. J., Showman A. P., 2019, The Astronomical Journal, 158, 166
  * Bell & Cowan (2018) Bell T. J., Cowan N. B., 2018, ApJ, 857, L20
  * Bell et al. (2021) Bell T. J., et al., 2021, Monthly Notices of the Royal Astronomical Society, 504, 3316
  * Beltz et al. (2022a) Beltz H., Rauscher E., Roman M. T., Guilliat A., 2022a, AJ, 163, 35
  * Beltz et al. (2022b) Beltz H., Rauscher E., Kempton E. M. R., Malsky I., Ochs G., Arora M., Savel A., 2022b, AJ, 164, 140
  * Borsa et al. (2021) Borsa F., Allart R., Casasayas-Barris N., et al. 2021, Astronomy & Astrophysics, 645, A24
  * Bourrier et al. (2020) Bourrier V., et al., 2020, Astronomy & Astrophysics, 637, A36
  * Caldas et al. (2019) Caldas A., Leconte J., Selsis F., Waldmann I., Bordé P., Rocchetto M., Charnay B., 2019, Astronomy & Astrophysics, 623, A161
  * Carone et al. (2020) Carone L., et al., 2020, Monthly Notices of the Royal Astronomical Society, 496, 3582
  * Cho et al. (2015) Cho J.-K., Polichtchouk I., Thrastarson H. T., 2015, Monthly Notices of the Royal Astronomical Society, 454, 3423
  * Coulombe et al. (2023) Coulombe L.-P., et al., 2023, arXiv e-prints, p. arXiv:2301.08192
  * Cowan et al. (2012) Cowan N. B., Machalek P., Croll B., Shekhtman L. M., Burrows A., Deming D., Greene T., Hora J. L., 2012, ApJ, 747, 82
  * Daylan et al. (2021) Daylan T., et al., 2021, The Astronomical Journal, 161, 131
  * Deming et al. (2023) Deming D., Line M. R., Knutson H. A., Crossfield I. J. M., Kempton E. M. R., Komacek T. D., Wallack N. L., Fu G., 2023, AJ, 165, 104
  * Dobbs-Dixon & Lin (2008) Dobbs-Dixon I., Lin D. N. C., 2008, ApJ, 673, 513
  * Eftekhar (2022) Eftekhar M., 2022, Rev. Mex. Astron. Astrofis., 58, 99
  * Eftekhar & Adibi (2022) Eftekhar M., Adibi P., 2022, The Planetary Science Journal, 3, 255
  * Ehrenreich et al. (2020) Ehrenreich D., et al., 2020, Nature, 580, 597
  * Eker et al. (2018) Eker Z., et al., 2018, Monthly Notices of the Royal Astronomical Society, 479, 5491
  * Fortney et al. (2008) Fortney J. J., Lodders K., Marley M. S., Freedman R. S., 2008, The Astrophysical Journal, 678, 1419
  * Gao et al. (2021) Gao P., Wakeford H. R., Moran S. E., Parmentier V., 2021, Journal of Geophysical Research (Planets), 126, e06655
  * Garhart et al. (2020) Garhart E., et al., 2020, The Astronomical Journal, 159, 137
  * Gibson et al. (2020) Gibson N. P., et al., 2020, Monthly Notices of the Royal Astronomical Society, 493, 2215
  * Guillot (2010) Guillot T., 2010, A&A, 520, A27
  * Guillot et al. (1996) Guillot T., Burrows A., Hubbard W. B., Lunine J. I., Saumon D., 1996, ApJ, 459, L35
  * Hammond & Lewis (2021) Hammond M., Lewis N. T., 2021, Proceedings of the National Academy of Sciences, 118, e2022705118
  * Hammond & Pierrehumbert (2018) Hammond M., Pierrehumbert R., 2018, The Astrophysical Journal, 869, 65
  * Hammond et al. (2020) Hammond M., Tsai S.-M., Pierrehumbert R. T., 2020, The Astrophysical Journal, 901, 78
  * Helling (2019) Helling C., 2019, Annual Review of Earth and Planetary Sciences, 47, 583
  * Helling et al. (2022) Helling C., et al., 2022, arXiv e-prints, p. arXiv:2208.05562
  * Heng & Showman (2015) Heng K., Showman A. P., 2015, Annual Review of Earth and Planetary Sciences, 43, 509
  * Heng et al. (2011) Heng K., Menou K., Phillipps P. J., 2011, MNRAS, 413, 2380
  * Hoeijmakers et al. (2018) Hoeijmakers H., Ehrenreich D., Heng K., et al. 2018, Nature, 560, 453
  * Hoeijmakers et al. (2020) Hoeijmakers H., Seidel J., Pino L., et al. 2020, Astronomy & Astrophysics, 641, A123
  * Hubeny et al. (2003) Hubeny I., Burrows A., Sudarsky D., 2003, The Astrophysical Journal, 594, 1011
  * Jackson et al. (2019) Jackson B., Adams E., Sandidge W., Kreyche S., Briggs J., 2019, The Astronomical Journal, 157, 239
  * Jansen & Kipping (2020) Jansen T., Kipping D., 2020, Monthly Notices of the Royal Astronomical Society, 494, 4077
  * Kataria et al. (2014) Kataria T., Showman A. P., Fortney J. J., Marley M. S., Freedman R. S., 2014, ApJ, 785, 92
  * Kataria et al. (2015) Kataria T., Showman A. P., Fortney J. J., Stevenson K. B., Line M. R., Kreidberg L., Bean J. L., Désert J.-M., 2015, ApJ, 801, 86
  * Kesseli et al. (2022) Kesseli A., Snellen I., Casasayas-Barris N., Molliére P., Sanchez-Lopez A., 2022, The Astronomical Journal, 163, 107
  * Kitzmann et al. (2018) Kitzmann D., et al., 2018, The Astrophysical Journal, 863, 183
  * Koll & Komacek (2018) Koll D. D. B., Komacek T. D., 2018, ApJ, 853, 133
  * Komacek & Showman (2016) Komacek T. D., Showman A. P., 2016, The Astrophysical Journal, 821, 16
  * Komacek & Showman (2019) Komacek T. D., Showman A. P., 2019, The Astrophysical Journal, 888, 2
  * Komacek & Tan (2018) Komacek T. D., Tan X., 2018, Research Notes of the American Astronomical Society, 2, 36
  * Komacek et al. (2017) Komacek T. D., Showman A. P., Tan X., 2017, ApJ, 835, 198
  * Komacek et al. (2022a) Komacek T. D., Tan X., Gao P., Lee E. K. H., 2022a, ApJ, 934, 79
  * Komacek et al. (2022b) Komacek T. D., Gao P., Thorngren D. P., May E. M., Tan X., 2022b, ApJ, 941, L40
  * Kreidberg et al. (2018) Kreidberg L., et al., 2018, The Astronomical Journal, 156, 17
  * Lee et al. (2021) Lee E. K. H., Parmentier V., Hammond M., Grimm S. L., Kitzmann D., Tan X., Tsai S.-M., Pierrehumbert R. T., 2021, MNRAS, 506, 2695
  * Lee et al. (2022) Lee E. K. H., Prinoth B., Kitzmann D., Tsai S.-M., Hoeijmakers J., Borsato N. W., Heng K., 2022, MNRAS, 517, 240
  * Lee et al. (2023) Lee E. K. H., Tsai S.-M., Hammond M., Tan X., 2023, A&A, 672, A110
  * Li & Goodman (2010) Li J., Goodman J., 2010, The Astrophysical Journal, 725, 1146
  * Lines et al. (2019) Lines S., Mayne N., Manners J., Boutle I., Drummond B., Mikal-Evans T., Kohary K., Sing D., 2019, Monthly Notices of the Royal Astronomical Society, 488, 1332
  * Liu & Showman (2013) Liu B., Showman A. P., 2013, ApJ, 770, 42
  * Lothringer & Barman (2019) Lothringer J. D., Barman T., 2019, ApJ, 876, 69
  * Lothringer et al. (2018) Lothringer J. D., Barman T., Koskinen T., 2018, The Astrophysical Journal, 866, 27
  * Lupu et al. (2021) Lupu R., Freedman R., Gharib-Nezhad E., Visscher C., Molliere P., 2021, Correlated k coefficients for H2-He atmospheres; 11 spectral windows and 1460 pressure-temperature points, doi:10.5281/zenodo.5590997, https://doi.org/10.5281/zenodo.5590997
  * Lupu et al. (2022) Lupu R., Freedman R., Gharib-Nezhad E., Visscher C., Molliere P., 2022, Correlated k coefficients for H2-He atmospheres; 622 spectral windows and 1460 pressure-temperature points, doi:10.5281/zenodo.6829310, https://doi.org/10.5281/zenodo.6829310
  * Mancini et al. (2022) Mancini L., et al., 2022, Monthly Notices of the Royal Astronomical Society, 509, 1447
  * Mansfield et al. (2018) Mansfield M., Bean J., Line M., et al. 2018, The Astronomical Journal, 156, 10
  * Mansfield et al. (2020) Mansfield M., et al., 2020, ApJ, 888, L15
  * Mansfield et al. (2021) Mansfield M., et al., 2021, Nature Astronomy, 5, 1224
  * Marley et al. (1999) Marley M. S., Gelino C., Stephens D., Lunine J. I., Freedman R., 1999, ApJ, 513, 879
  * Marley et al. (2021) Marley M. S., et al., 2021, The Astrophysical Journal, 920, 85
  * Maxted et al. (2013) Maxted P. F. L., et al., 2013, MNRAS, 428, 2645
  * May et al. (2021) May E. M., et al., 2021, The Astronomical Journal, 162, 158
  * May et al. (2022) May E. M., et al., 2022, AJ, 163, 256
  * Mayne et al. (2014) Mayne N. J., et al., 2014, A&A, 561, A1
  * Mendonça et al. (2016) Mendonça J. M., Grimm S. L., Grosheintz L., Heng K., 2016, ApJ, 829, 115
  * Mikal-Evans et al. (2022) Mikal-Evans T., et al., 2022, Nature Astronomy, 6, 471
  * Mikal-Evans et al. (2023) Mikal-Evans T., et al., 2023, ApJ, 943, L17
  * Mukherjee et al. (2023) Mukherjee S., Batalha N. E., Fortney J. J., Marley M. S., 2023, ApJ, 942, 71
  * Owens et al. (2021) Owens N., de Mooij E. J. W., Watson C. A., Hooton M. J., 2021, MNRAS, 503, L38
  * Pai Asnodkar et al. (2022) Pai Asnodkar A., Wang J., Eastman J. D., Cauley P. W., Gaudi B. S., Ilyin I., Strassmeier K., 2022, AJ, 163, 155
  * Parmentier et al. (2013) Parmentier V., Showman A. P., Lian Y., 2013, Astronomy & Astrophysics, 558, A91
  * Parmentier et al. (2016) Parmentier V., Fortney J. J., Showman A. P., Morley C., Marley M. S., 2016, ApJ, 828, 22
  * Parmentier et al. (2018) Parmentier V., et al., 2018, Astronomy & Astrophysics, 617, A110
  * Parmentier et al. (2021) Parmentier V., Showman A. P., Fortney J. J., 2021, Monthly Notices of the Royal Astronomical Society, 501, 78
  * Pelletier et al. (2023) Pelletier S., et al., 2023, arXiv e-prints, p. arXiv:2306.08739
  * Perez-Becker & Showman (2013) Perez-Becker D., Showman A. P., 2013, ApJ, 776, 134
  * Perna et al. (2010) Perna R., Menou K., Rauscher E., 2010, The Astrophysical Journal, 719, 1421
  * Perna et al. (2012) Perna R., Heng K., Pont F., 2012, ApJ, 751, 59
  * Pierrehumbert & Hammond (2019) Pierrehumbert R. T., Hammond M., 2019, Annual Review of Fluid Mechanics, 51, 275
  * Pino et al. (2022) Pino L., et al., 2022, A&A, 668, A176
  * Rauscher & Menou (2012) Rauscher E., Menou K., 2012, The Astrophysical Journal, 750, 96
  * Rauscher & Menou (2013) Rauscher E., Menou K., 2013, The Astrophysical Journal, 764, 103
  * Robbins-Blanch et al. (2022) Robbins-Blanch N., Kataria T., Batalha N. E., Adams D. J., 2022, ApJ, 930, 93
  * Rogers (2017) Rogers T., 2017, Nature Astronomy, 1, 0131
  * Rogers & Komacek (2014) Rogers T. M., Komacek T. D., 2014, The Astrophysical Journal, 794, 132
  * Roman et al. (2021) Roman M. T., Kempton E. M.-R., Rauscher E., Harada C. K., Bean J. L., Stevenson K. B., 2021, The Astrophysical Journal, 908, 101
  * Roth et al. (2021) Roth A., Drummond B., Hébrard E., Tremblin P., Goyal J., Mayne N., 2021, MNRAS, 505, 4515
  * Savel et al. (2022) Savel A., Kempton E., Malik M., et al. 2022, The Astrophysical Journal, 926, 85
  * Savel et al. (2023) Savel A. B., Kempton E. M. R., Rauscher E., Komacek T. D., Bean J. L., Malik M., Malsky I., 2023, ApJ, 944, 99
  * Schneider et al. (2022) Schneider A. D., Carone L., Decin L., Jørgensen U. G., Mollière P., Baeyens R., Kiefer S., Helling C., 2022, Astronomy & Astrophysics, 664, A56
  * Sheppard et al. (2017) Sheppard K. B., Mandell A. M., Tamburo P., Gandhi S., Pinhas A., Madhusudhan N., Deming D., 2017, The Astrophysical Journal Letters, 850, L32
  * Showman & Guillot (2002) Showman A. P., Guillot T., 2002, A&A, 385, 166
  * Showman & Polvani (2011) Showman A. P., Polvani L. M., 2011, The Astrophysical Journal, 738, 71
  * Showman et al. (2009) Showman A. P., Fortney J. J., Lian Y., Marley M. S., Freedman R. S., Knutson H. A., Charbonneau D., 2009, ApJ, 699, 564
  * Showman et al. (2013a) Showman A. P., Wordsworth R. D., Merlis T. M., Kaspi Y., 2013a, Comparative Climatology of Terrestrial Planets, 1, 277
  * Showman et al. (2013b) Showman A. P., Fortney J. J., Lewis N. K., Shabram M., 2013b, ApJ, 762, 24
  * Showman et al. (2015) Showman A. P., Lewis N. K., Fortney J. J., 2015, The Astrophysical Journal, 801, 95
  * Showman et al. (2020) Showman A. P., Tan X., Parmentier V., 2020, Space Sci. Rev., 216, 139
  * Skinner & Cho (2022) Skinner J., Cho J. Y., 2022, Monthly Notices of the Royal Astronomical Society, 511, 3584
  * Spiegel et al. (2009) Spiegel D. S., Silverio K., Burrows A., 2009, The Astrophysical Journal, 699, 1487
  * Steinrueck et al. (2019) Steinrueck M. E., Parmentier V., Showman A. P., Lothringer J. D., Lupu R. E., 2019, ApJ, 880, 14
  * Steinrueck et al. (2023) Steinrueck M. E., Koskinen T., Lavvas P., Parmentier V., Zieba S., Tan X., Zhang X., Kreidberg L., 2023, arXiv preprint arXiv:2305.09654
  * Stevenson et al. (2017) Stevenson K. B., et al., 2017, The Astronomical Journal, 153, 68
  * Tan & Komacek (2019) Tan X., Komacek T. D., 2019, ApJ, 886, 26
  * Tan & Showman (2020) Tan X., Showman A. P., 2020, ApJ, 902, 27
  * Tan & Showman (2021) Tan X., Showman A. P., 2021, MNRAS, 502, 2198
  * Thorngren et al. (2019) Thorngren D., Gao P., Fortney J., 2019, The Astrophysical Journal Letters, 884, L6
  * Tsai et al. (2014) Tsai S.-M., Dobbs-Dixon I., Gu P.-G., 2014, ApJ, 793, 141
  * Tsai et al. (2018) Tsai S.-M., Kitzmann D., Lyons J. R., Mendonça J., Grimm S. L., Heng K., 2018, The Astrophysical Journal, 862, 31
  * Von Essen et al. (2020) Von Essen C., Mallonn M., Borre C., Antoci V., Stassun K., Khalafinejad S., Tautvaišienė G., 2020, Astronomy & Astrophysics, 639, A34
  * Wardenier et al. (2021) Wardenier J., Parmentier V., Lee E., Line M., Gharib-Nezhad E., 2021, Monthly Notices of the Royal Astronomical Society, 506, 1258
  * Wardenier et al. (2022) Wardenier J. P., Parmentier V., Lee E. K., 2022, Monthly Notices of the Royal Astronomical Society, 510, 620
  * Wong et al. (2015) Wong I., et al., 2015, ApJ, 811, 122
  * Wong et al. (2016) Wong I., et al., 2016, ApJ, 823, 122
  * Wong et al. (2020) Wong I., et al., 2020, The Astronomical Journal, 160, 155
  * Wong et al. (2021) Wong I., et al., 2021, AJ, 162, 256
  * Wong et al. (2022) Wong I., Shporer A., Vissapragada S., Greklek-McKeon M., Knutson H. A., Winn J. N., Benneke B., 2022, The Astronomical Journal, 163, 175
  * Yan et al. (2019) Yan F., et al., 2019, Astronomy & Astrophysics, 632, A69
  * Zamyatina et al. (2023) Zamyatina M., et al., 2023, Monthly Notices of the Royal Astronomical Society, 519, 3129
  * Zellem et al. (2014) Zellem R. T., et al., 2014, The Astrophysical Journal, 790, 53
  * Zhang (2020) Zhang X., 2020, Research in Astronomy and Astrophysics, 20, 099
  * Zhang & Showman (2017) Zhang X., Showman A. P., 2017, ApJ, 836, 73
  * Zhang et al. (2018) Zhang M., et al., 2018, AJ, 155, 83

##  Appendix A Additional figure

Figure 16: Dayside brightness temperature from models as shown in Figure 11
but with additional outputs at the TESS band and without observational data
for better visualization. Solid lines are models with no drag and three
stellar temperatures; dashed lines are models with
$\tau_{\mathrm{drag}}=10^{4}$ s and three stellar temperatures; finally, the
dotted line is from models with $\tau_{\mathrm{drag}}=10^{4}$ s, a stellar
temperature of 6500 K but excludes TiO, VO and Fe.

##  Appendix B Varying metallicity and C/O ratio: a case study for TOI 1431b

Figure 17: Simulated phase curves in the TESS band from models appropriate for
TOI 1431b. In the left panel, we show models with a fixed 1 solar C/O ratio
but varying metallicity; in the right panel, we show those with a fixed 1
solar metallicity but varying C/O ratio. Figure 18: The fractional day-night
brightness temperature variations in the TESS band as a function of varying
metallicity and C/O ratios for models of TOI 1431b. Figure 19: Dayside and
nightside TP profiles along with equilibrium water vapor volume mixing ratio
in the temperature-pressure space for two models of TOI 1431b: on the left is
the 0.3 solar metallicity case and on the right is the 3 solar metallicity
case. Both cases have a 1 solar C/O ratio. Figure 20: Dayside (left column)
and nightside (right column) emission spectra from models of TOI 1431b with
varying metallicity and C/O ratio. Top row: models with a fixed 1 solar C/O
ratio but with varying metallicity. Bottom row: models with a fixed 1 solar
metallicity but with varying C/O ratio.

The systematic overestimation of the modeled day-night temperature differences
in the TESS band motivates us to look for alternative solutions to enhance
heat transport. The majority of TESS targets that have day-night temperature
contrasts much lower than model predictions do not have phase curve
observations in other bandpasses, and so in this section, we discuss TESS
phase curve modeling only. We vary the metallicity and C/O ratio to see if it
would cause poor nightside radiative cooling that could lead to better day-
night heat transport. Rather than modeling grids over planetary equilibrium
and stellar temperatures which require significant computational resources, we
focus on one planetary condition while changing the metallicity and C/O ratio.
We use TOI 1431b as a template. Our previous work using a semi-grey GCM with
radiatively active cloud tracers (Komacek et al., 2022a) has also modeled this
planet.

We carry out a grid of models with metallicity of 0.1 0.3, 1, and 3 solar, and
C/O ratio of 0.5, 1, 2, and 2.5 solar. Note that in this set of GCMs we did
not include opacities of Fe and a few other minor species in the k tables as
was done in Section 3. We have tried modeling this grid using the same k
tables as in Section 3 but many of these models are computationally unstable.
We then have to switch to a version of k tables without the Fe opacity that
guarantee the stability of the GCMs. The main discrepancy is that there is no
strong thermal inversion at very low pressure anymore which was supposed to be
mainly driven by Fe absorption, and the stellar flux reaching the deeper IR
photosphere is thus slightly more. Nevertheless, this is not expected to
affect our qualitative conclusions.

We first show the resulting phase curves in the TESS bandpass from subsets of
this grid in Figure 17. The normalized phase curves with a C/O ratio of 1
solar but varying metallicity are shown on the left panel. Higher-metallicity
cases show higher day-night flux contrast and smaller phase-curve offsets and
this trend is quite robust. In the right panel, with 1 solar metallicity,
varying the C/O ratio does not have a clear trend and the over effect is less
than varying the metallicity while holding the C/O ratio fixed as shown on the
left. The trend that high metallicity generally results in larger day-to-night
flux contrast and more aligned peak with the secondary eclipse is in agreement
with previous studies of Dobbs-Dixon & Lin (2008); Kataria et al. (2014);
Kataria et al. (2015).

The normalized day-night temperature variations in the TESS bandpass of this
grid are summarized in Figure 18, in which a grey shaded area represents the
observed range of TOI 1431b. As seen in the phase curves, there is a general
trend of increasing day-night temperature variations with increasing
metallicity, while the varying C/O ratio adds on certain complexity. While
low-metallicity GCMs significantly reduce the day-night temperature variation
to about 0.23, this value is still significantly larger than the upper limit
reported for TOI 1431b. A qualitatively similar trend of the day-night
temperature variations also holds in other wavelengths which we do not further
discuss here.

The trend could be understood via the competition between the day-night heat
transport and radiative cooling/heating processes. The former is on the zeroth
order held invariant for a given setup because the dominant heat transport
near the photosphere is the hydrogen dissociation and recombination which is a
priori independent of metallicity. Radiative heating/cooling efficiency is our
primary focus. Decreasing the metallicity generally results in a deeper
photosphere simply because the opacity is positively correlated to the
metallicity; therefore the longer radiative timescale at higher pressures
(Showman & Guillot, 2002) naturally leads to smaller day-night temperature
variations (Perez-Becker & Showman, 2013; Komacek & Showman, 2016). Examples
are shown in Figure 19 in which 0.3 solar and 3 solar metallicities are
displayed. The nightside structure at relatively low pressures is also
significantly heated in the low-metallicity case (see the left panel of Figure
19). Without sufficient cooling agents in low-metallicity cases, for example,
water content as shown in the color contours of Figure 19 which is one of the
most important coolants in giant planet atmospheres, the upper atmosphere is
driven toward a strong thermal inversion. The hotter it gets, the less water
to cool due to thermal dissociation, and this is positive feedback to drive a
strong thermal inversion on the nightside of ultra-hot Jupiters. The other
factor enhancing the nightside thermal inversion is the increasing heat
transport efficiency by hydrogen dissociation and recombination because of
increasing atomic hydrogen mixing ratio with decreasing pressure. Whereas in
the high-metallicity case on the right panel, water vapor is rich to
sufficiently cool on the nightside and we see a decreasing temperature with
decreasing pressure of the nightside profile that is typical in other hot
Jupiter GCMs.

Meanwhile, we are also cautious about the systematically low day-night
temperature variations observed in the TESS band. Some studies on the same
objects claimed only the upper limits of the nightside flux which are smaller
than those quoted in Figure 12 (e.g., Wong et al., 2020; Wong et al., 2022),
indicating that the measurements could also be sensitive to data reduction
methods. The significant mismatch of the GCM-derived day-night variations to
the observed TESS values requires a better understanding of the data reduction
processes in addition to searching for more theoretical possibilities.
Nonetheless, this issue motivated the exploration of the strong dependence of
hot Jupiter circulation on atmospheric compositions.

Lastly, we show the prediction of the dayside and nightside emission spectra
obtained from a subset of the models with varying metallicity and the C/O
ratio in Figure 20. With a 1-solar C/O ratio, varying metallicity in the top
panels merely changes the size of the emission features on the dayside and the
absorption features on the nightside. With a 1-solar metallicity, cases with
0.5 and 1 solar C/O ratios in the bottom panels show emission on the dayside
and absorption on the nightside. However, cases with 2 and 2.5 solar C/O
ratios are qualitatively different, lacking features in the water band, and
are significantly closer to the blackbody spectrum. This is presumably because
of the lack of sufficient water vapor in the atmosphere, and other opacity
sources such as that of H- contribute more to the near-IR wavelength. The
super solar C/O ratio to explain the much more muted features has been
proposed for the HST dayside spectrum and Spitzer photometry of WASP-18b
(Sheppard et al., 2017) based on 1D retrievals for only the dayside (although
recently been ruled out by JWST observations Coulombe et al., 2023). The
spectral behaviors for super solar C/O ratio appear also on the nightside
which is partly driven by the atmospheric heat transport and its effect on the
nightside thermal inversion, and this should provide a unique test case for
testing atmospheric heat transport with varying C/O ratio and metallicity.