Evidence of radius inflation in radiative GCM models of WASP-76b due to the advection of potential temperature

the discrepancy between the radii of observed hot Jupiters and standard ‘radiativ e-conv ectiv e’ models


INTRODUCTION
Observations of hot Jupiters (Laughlin et al. 2011) and hot brown dwarfs (see Fig. 4 of Casewell et al. 2020) have revealed a significant discrepancy between standard 'radiative-convective' single column (1D) atmospheric models and the properties of observed objects: observed radii of highly irradiated objects tend to be significantly larger than 1D atmospheric models suggest (see, for example, Figure 1 of Komacek & Youdin 2017).This indicates that said 1D models are likely failing to capture some key physics or dynamics which drive the observed radius discrepancy (i.e.inflation).In 1D models this discrepancy is 'solved' via the inclusion of an intrinsic/internal temperature, which essentially acts to heat the deep atmosphere ★ E-mail: f.sainsbury-martinez@leeds.ac.uk (internal adiabat) to a more physical value, allowing for atmospheric retrievals of transit observations, without actually elucidating on exactly what dynamics drives this heating other than typically claiming that it is linked with thermal escape from the interior (Guillot & Showman 2002;Baraffe et al. 2003;Sudarsky et al. 2003;Chabrier et al. 2004;Thorngren et al. 2019).In a collective effort to understand this deep heating/radius inflation problem, a vast array of different physical mechanisms have been suggested as possible causes/solutions (see Baraffe et al. 2009;Fortney & Nettelmann 2010;Baraffe et al. 2014 for a more in-depth overview of many of the proposed mechanisms) including tidal heating and dissipation (Arras & Socrates 2010;Lee 2019), the ohmic dissipation of electrical/magnetic energy (Batygin & Stevenson 2010;Perna et al. 2010;Rauscher & Menou 2012;Helling et al. 2021;Knierim et al. 2022), the deep deposition of kinetic energy (Guillot & Showman 2002), enhanced opacities which inhibit interior cooling (Burrows et al. 2007), double-diffusive convection which hampers convective heat transport (Chabrier & Baraffe 2007), or the vertical advection of potential temperature (first proposed and studied in 2D by Tremblin et al. 2017and studied in 3D by Sainsbury-Martinez et al. 2019, 2021).Fortunately, observations can help us to narrow down which of the above mechanisms might be responsible for the observed radius inflation.Specifically, observational studies of hot Jupiters and hot brown dwarfs (e.g.Demory & Seager 2011;Laughlin et al. 2011;Lopez & Fortney 2016;Sestovic et al. 2018;Casewell et al. 2020) have revealed a clear trend: a general increase in the observed radius of highly irradiated gaseous planets with stellar irradiation, except in the case of very-highly irradiated objects in very short orbits (e.g.SDSS1411B - Casewell et al. 2018) where little to no inflation is observed.One such mechanism which can explain this trend without the inclusion of model-dependent fine tuning is the vertical advection of potential temperature (i.e.enthalpy).Briefly this mechanism can be understood as follows: for a tidally locked, gaseous exoplanet, the strong stellar irradiation leads to a very hot outer atmosphere paired with a very strong super-rotating equatorial jet.This driving can be understood via a 2D stationary circulation model, in which, due to mass and angular momentum conservation, significant vertical winds arise (as proposed/seen in Tremblin et al. 2017;Sainsbury-Martinez et al. 2019).Note that such a view is opposed by Showman & Polvani (2011), who assume/propose that only the irradiated layers of the atmosphere are meteorologically active, and that deeper layers are either quiescent or purely convective (not that the latter would have any negative implications for our mechanism, beyond the temperature of the adiabat).Not only do our results disagree with this view (see Figure 2), but other studies, such as Carone et al. (2020); Schneider et al. (2022a) have shown that significant wave activity and zonal/vertical winds can occur in these deep atmospheric layers.If this holds true, and we propose that it does, these aforementioned vertical winds carry high potential temperature fluid parcels from the hot (radiative) outer atmosphere deep into the interior (where radiative effects tend to zero -as shown in Figure 3), driving the formation of a non-convective (i.e.advective) adiabat at lower pressures than 1D models (without an artificially increased internal/intrinsic temperature) would predict.Because this adiabat forms at lower pressures , and because the radiative, advective and deep convective (i.e.interior) regions must smoothly connect, the internal adiabats temperature temperature is significantly increased when compared to a model which lacks advection and considers a radiative-convective boundary alone.In turn, this increase in the temperature of the internal adiabat, leads to an increase in the internal entropy, and hence an inflated radius.This is very similar to what occurs in a 1D model when the internal/intrinsic temperature is increased, although here it is occurring due to fundamental physics.An example of this can be seen in Thorngren et al. (2019), who find a clear link between the pressure of the RCB (radiative-convective boundary), i.e.where the outer atmosphere connects with the interior adiabat, and the intrinsic temperature, i.e. internal heat flux that their models impose.However the heating which drives the formation of this non-convective adiabat has nothing to do with heat transport from the interior.Rather it is heating associated with the irradiated outer atmosphere, which should, at steady-state, balance any outwards heat transport from the interior, stalling any internal cooling and leading to a net zero internal flux (i.e.no heating from the interior), a stable, inflated, radii, and a natural link between radius inflation and surface irradiation.
It is important to note that this mechanism is distinct from the kinetic energy transport and deposition mechanism proposed by Guillot & Showman (2002).In their mechanism, stellar irradiation is converted to kinetic energy in the outer atmosphere (by atmospheric pressure gradients), this energy is then somehow transported down towards the interior (possibly by, for example, Kelvin-Helmholtz instabilities, vertical advection, or waves), where it then dissipates, heating the deep atmosphere and warming the internal adiabat.Rather, the mechanism we (and Tremblin et al. 2017;Sainsbury-Martinez et al. 2019) propose skips these uncertain energy conversion steps, and instead we directly transport hot (high enthalpy) material from from the outer atmosphere to the deep atmosphere via already present flows and circulations.
Recently, Schneider et al. (2022b), called into question the validity of vertical potential temperature advection as a possible explanation for the radius inflation of the ultra-hot Jupiter WASP-76b (West et al. 2016;Seidel et al. 2019;Ehrenreich et al. 2020;Kesseli et al. 2022), arguing that their (hot-start) 3D atmospheric models, calculated with expeRT/MITgcm, and including a self-consistent, non-grey radiative transfer model (see Schneider et al. 2022a for a detailed discussion of this code), suggested that coupling between radiation and dynamics alone is not sufficient to explain the inflated radii of highly-irradiated, gaseous, exoplanets.
Here, we intend to investigate this claim in more detail, performing additional analysis of the nominal WASP-76b simulation discussed in Schneider et al. (2022b) along with additional, cooler-start (i.e.cooler initial deep adiabats) calculations that were run exclusively for this work.Specifically, we intend to investigate the vertical mass and enthalpy (i.e.potential temperature) transport in these models, confirming if vertical advection plays a significant role in the dynamics, before comparing the steady-state 3D models with internal-structure models based upon the work of Baraffe et al. (2010) in order to confirm how much, if any, of the inflated radius of WASP-76b potential temperature advection alone can explain.In section 2, we start with a brief overview of expeRT/MITgcm before introducing the models discussed as part of this work.This is followed, in section 3 with our analysis, focusing on the vertical transport of potential temperature and its implications for the steady-state deep atmosphere of our WASP-76b models.We finish, in section 4 by discussing the implications of our results, with a particular focus on the sustainability of potential temperature advection as an explanation for the inflated radii of highly irradiated, tidally locked, gaseous exoplanets.

METHODS
The methodology and models used in this work are based on the work of Schneider et al. (2022a) and Schneider et al. (2022b).Here we give a brief overview of the GCM used to calculate the WASP-76b models considered here, before giving a more in depth description of said WASP-76b models setup.

expeRT/MITgcm
Briefly, expeRT/MITgcm (Carone et al. 2020;Schneider et al. 2022a) builds on the dynamical core of the MITgcm (Adcroft et al. 2004), pairing said core with the petitRADTRANS (Mollière et al. 2019) radiative transfer (RT) model in order to enable the long model integration times required to explore the steady-state dynamics of the deep atmospheres of hot Jupiters, whilst also maintaining the accuracy of a multi-wavelength radiation scheme.expeRT/MITgcm solves the primitive equations of meteorology (Vallis 2006;Showman et al. 2009), for an ideal gas, on an Arakawa C-type cubed sphere (designed to avoid numerical issues near the poles which occur due to singularities in the coordinate system; for more details of this grid, see, for example, Miller 1984) with a horizontal resolution C32 1 and a vertical grid that contains a combination of 41 linearly in log() (i.e.log-pressure) spaced layers between 1 × 10 −5 and 100 bar, paired with 6 linearly in  spaced layers between 100 and 700 bar.As in Showman et al. (2009); Carone et al. (2020) this model includes a horizontal fourth-order Shapiro filter (with  = 25 s) in order to smooth grid-scale noise.Additionally, expeRT/MITgcm includes a linear Rayleigh-drag (which is also known as a linear-basel drag scheme -see Carone et al. 2020, particularly Section 2.3 and Appendix A for a discussion of this dynamics preserving approach as well as comparisons with other drag-schemes) at the bottom of the atmosphere (between 490 and 700 bar) and a sponge layer at the top of the atmosphere (for  < 1 × 10 −4 bar).We discuss the implications of this Rayleigh-drag on the vertical advection of potential temperature, and hence radius inflation, in more detail in section 3. Note: we selected 700 bar as the maximum pressure of our simulation domain in-order to balance modelling a sufficient portion of the deep atmosphere with the increasing computational costs of modelling high-pressure regions (due to their increased dynamical timescales).Radiatively, the outer atmosphere is heated and cooled using a runtime (i.e.coupled), multi-wavelength, RT scheme based upon petitRADTRANS.Specifically, the radiative dynamics are updated every 100 seconds, quadruple the dynamical time step (Δ dy = 25 s), with the radiative transport calculated using a correlated-k approach that includes 5 wavelength bins each of which contain 16 Gaussian quadrature points (see Goody et al. 1989 for an introduction to the correlated-k approach to RT, and Appendix B of Schneider et al. 2022a for a discussion of the accuracy of the limited wavelength bin approach).Note that opacities for the RT scheme are based on a pre-calculated pressure-temperature grid, assume local chemical equilibrium, and include the following gas absorbers (with data taken from the ExoMol 2 database): H 2 O, CO 2 , CH 4 , NH 3 , CO, PH 3 , H 2 S, TiO, VO, HCN, Na, K and FeH.Additionally, the RT model includes Rayleigh-scattering for both H 2 and He, and collision-induced-absorption for H 2 − H 2 , He − He and H − (see Schneider et al. 2022b for more details).We do not include equilibrium condensation since, assuming that the latent heat release is low (Woitke & Helling 2003;Helling 2019;Helling et al. 2019), it should have little effect on the photosphere, especially for WASP-76b whose day-side can be assumed to be cloud free, and is simply too hot for condensation to occur.
Finally, the inclusion of an artificial Rayleigh-drag scheme in the deep atmosphere implies that an additional energy source term must be added to the deep atmosphere to account for the conversion of energy lost from drag to heat (Rauscher & Menou 2013;Carone et al. 2020;Schneider et al. 2022a), which is then locally returned to 1 C32 is comparable to a resolution of 128 × 64 in longitude and latitude 2 www.exomol.comand Tennyson et al. ( 2016); Chubb et al. (2021) the atmosphere.This takes the form: where  is the local temperature of the atmosphere,  is the horizontal (zonal plus meridional) wind speed,    is the Rayleigh-drag timescale at the bottom of the atmosphere, and   is the heat capacity at constant pressure.

Models of WASP-76b
WASP-76b is a tidally locked ultra-hot Jupiter-like planet ( = 0.92 ± 0.03M J ) that orbits its host star at a distance of 0.033 AU, corresponding to an orbital period of 1.81 Earth days, and which appears to exhibit significant radius inflation, with an observed radius of 1.83 ± 0.06R J (West et al. 2016).The host star, WASP-76, is a hot yellow-white (F7V) main-sequence star with an effect temperature of  eff = 6250 ± 100 K and a radius of  * = 1.73 ± 0.04R ⊙ (Gaia Collaboration et al. 2018).Further, all our models assume a fixed specific heat capacity,   = 13784 J kg −1 K −1 and a fixed specific gas constant,  = 3707 J kg −1 K −1 , which corresponds to a adiabatic index  ≃ 1.36 (these values have been extracted from petitRAD-TRANS).However the Rayleigh-drag timescale does vary, with the majority of our models setting  drag = 1 day, and a low-drag model setting  drag = 1000 days.Finally, we include zero heat flux from the interior, meaning that any deep atmospheric heating is purely due to downwards enthalpy advection from the irradiated outer atmosphere.
Here we consider 6 models of WASP-76b, five of which only differ in the temperature profile used to initialise them, and one in which the strength of the deep Rayleigh-drag has been reduced (as previously mentioned).For the former models, the initialisation profile is a combination of an isotherm, based upon the stellar irradiation, in the outer atmosphere (i.e. for  < 1 bar), and an adiabat, with a reference temperature () taken at 1 bar, throughout the deep atmosphere (i.e. > 10 bar), with a linear interpolation between the two profiles between 1 and 10 bar.Here we consider reference temperatures of  = 4000 K (i.e. the nominal model which was first presented in Schneider et al. 2022b, but which has been further evolved as part of this work), 2500 K, 1800 K, 1400 K and 1000 K, which range from hotter than the adiabat of the final nominal model of Schneider et al. (2022b) to cooler -thus allowing us to explore models in which the deep atmosphere is both heating and cooling.These initial profiles can be seen in Figure 1, where we plot the initial profile of each variable initialisation model as a dashed line.On the other hand, the low-drag model (with  drag = 1000 days) is initialised from a snapshot of the nominal model (with  = 4000 K) taken after 40,000 days of simulation time.Note that, other than the nominal model, all the models featured here were performed as part of this work.

RESULTS
A broader analysis of the nominal model, after 86,000 days of runtime, is presented in Schneider et al. (2022b).Instead, here, we focus our analysis on the vertical advection of potential temperature, including what drives this advection, what effect it has on the deep atmosphere, and how much, if any, of WASP-76b radius inflation can be attributed to it.
We start our analysis with the nominal model, which, after over 155,000 days of simulation time (which corresponds to over 10,000 Horizontal mean Temperature-Pressure profiles for our five WASP-76b atmospheric models with different initialisation temperatures.For each of the variable initial temperature models considered here, i.e. the nominal (4000 K), 1000 K, 1400 K, 1800 K, and 2500 K start models, we include a profile near initialisation (dotted) and a profile at the end of the models runtime (solid).Note that the nominal model has been run for significantly longer than the other models (Table 1), and hence is likely to represent the steady state that all aforementioned models are converging towards.
advective turnover timescales in the deep atmosphere -see Figure 4), is approaching steady-state at almost all simulated pressures.Here we find that the strong day/night temperature difference associated with the combination of both tidal-locking and a hot host-star has resulted in the formation of a rapid super-rotating jet (see Figure 2a, which plots the zonal mean zonal wind at 155,000 days) that extends significantly into the deep atmosphere: at the equator the region in which  zonal > 1000 m s −1 extends to pressures greater than 10 bar.Such deep jets were already predicted in Carone et al. (2020) and confirmed in Schneider et al. (2022a,b).Here, we emphasise that these deep jets facilitate the formation of an advective adiabat at the same depths as Sainsbury-Martinez et al. (2019) propose to explain the inflated radius of HD209458b.In turn, strong latitudinal and vertical flows also develop, as can be seen in the meridional mass streamfunction (i.e. the meridional circulation profile -Equation 16of Sainsbury-Martinez et al. 2019).
In Figure 2b we plot the meridional circulation profile for the nominal model at near steady-state, with clockwise circulations shown in red and anticlockwise circulations shown in blue.Here we find that, at the equator, the strong stellar irradiation on the day-side leads to a general upwelling between 10 −5 and ∼ 1 bardriven by the combination of a clockwise circulation in the northern hemisphere and an anti-clockwise circulation in the south, both of which also drive material away from the substellar point/equator in the outer atmosphere.However, as we move deeper into the atmosphere, where the radiative time-scale is longer and hence advective effects can start to play a more significant role, we find that the sense of the meridional circulations has changed, likely due flows associated with the super-rotating jet taking over the vertical driving, leading to a strong downflow at the equator balanced by mass-conserving upflows at mid-latitudes (i.e.around 45 • -i.e. at the edge of the super-rotating jet).A similar circulation pattern was found by Sainsbury-Martinez et al. ( 2021) for Kepler-13Ab, a hot brown dwarf with a very hot (A-class) host star, and was shown to be sufficient to drive significant deep heating.
We next explore if this is also the case for our WASP-76b models.Specifically, we start by investigating the vertical transport of enthalpy.We first recall briefly how this quantity impacts the averaged energy transport in the atmosphere.Assuming the density is near steady-state (a similar assumption to the anelastic approximation), the mass and energy conservation equations are given by (2) where , , , and  are the atmospheric density, pressure, internal, and total energy;  the velocity of the flow;  the gravitational potential and  rad the radiative flux (including the irradiation from the host star).We will assume that the flow is low Mach in the deep atmosphere and therefore neglect the contribution of the kinetic energy.Furthermore we rewrite the energy flux as a function of the enthalpy  +  =   , with   specific heat capacity at constant pressure and  the temperature.By averaging Equation 2 in 2D over the full sphere (Ω), we get assuming there is no mass flux out of the domain of interest in a planeparallel approximation.The gravitational potential does not depend on latitude/longitude, therefore, because of mass conservation, its contribution to the energy flux is zero.Only the contribution of the enthalpy and the radiative flux remain: If the temperature is uniform, e.g. a 1D model, the contribution of the enthalpy is zero similarly to the contribution of the gravitational potential.If not, e.g. a 3D GCM, cold downflows and hot upflows will tend to cool the deep atmosphere whereas hot downflows and cold upflows will tend to warm the deep atmosphere.This is how the circulation can transport energy from the irradiated hot top layers to the deep atmosphere, even in the absence of convective processes.
This split between upflows and downflows can be seen in Figure 3, where we plot the longitudinal variation of the latitudinal-mean vertical enthalpy (top) and the horizontal-mean vertical enthalpy (bottom) for three models, two of which are near-steady-state (left - In the top row figures, we plot the vertical enthalpy flux profiles at 6 different longitudes, ranging from just east of the anti-stellar point to just west of the substellar point, as well as the global mean vertical enthalpy flux profile.However, since the mean flux is significantly smaller than the local fluxes, we replot the mean profiles in the bottom row in order to better demonstrate the vertical variations in enthalpy transport, focusing on the advection into the deep atmosphere.Here we also include the horizontal mean stellar (incoming) and planetary (outgoing) fluxes in order to reinforce that the deep atmosphere is radiatively quiescent.Note: the 2500 K profiles are calculated near the start of the simulation when the cooling is strongest -similar results can be found near initialisation for the nominal and other hot-start models.2004), of a hot Jupiter with a mass of 0.9M J and an inflated radius of 1.98R J (purple).1000 K initialisation -centre -4000 K, i.e. nominal, initialisation) and one which was initialised with a hot deep adiabat that is still rapidly cooling at the time of the snapshot (right -2500 K initialisation after 200 days).Starting with the longitudinal variations of the latitudinal-mean vertical enthalpy (top -Figure 3), it is clear that the direction of enthalpy transport varies significantly across the planetary surface.This was to be expected as tidally-locked thermal and wind dynamics, particularly in the outer atmosphere, are highly spatially inhomogeneous.However, this is further complicated by the effect that the temperature of the initial deep adiabat has on the overall dynamics -when a model is initialised with a deep atmosphere that is hotter than its final steady-state, excess energy must leave the deep atmosphere and, since radiative time-scales in the deep atmosphere are long, this typically occurs via changes in the wind structure and hence vertical enthalpy transport.An example of this effect can be found when comparing the models shown in Figure 3: for the hot-start (2500 K) model near initialisation, Figure 3c, we find that vertical enthalpy transport is primarily outwards, other than over a limited longitude and latitude range associated with a mass-conserving downflow.Almost the exact opposite scenario is found for a cool (1000 K) initialisation model (throughout its runtime), Figure 3a, where we find that the enthalpy flow is directed downwards at most longitudes, albeit, once again, with a mass conserving counter flow.Finally the nominal model, Figure 3b, represents a mix of the two regimes, with dynamics that can be linked to a combination of its very hot initialisation, leading to significant initial cooling, and long-run-time, leaving the model close to steady-state (although still warming in the deeper regions of the atmosphere due to the very-long dynamical times required to heat high-pressure regions of a hot Jupiter -see the isothermal model of Knierim et al. (2022)).
This difference in regime is also reflected in the horizontal-mean vertical enthalpy profiles (bottom -Figure 3): both the 1000 K and nominal models reveal a net downwards enthalpy flux, extending from the outer atmosphere all the way to the bottom of the simulation domain.Furthermore, this peak in the downwards flux is married with the radiative flux (both outwards and inwards) tending towards zero, as required in the potential temperature advection mechanism (Tremblin et al. 2017).Note that the vertical extent of the enthalpy downflow is reduced in the nominal model when compared with the 1000 K model, which is due to the nominal model being closer to steady-state and hence heating being limited to the deepest regions of the simulation domain (see Figure 7 of Sainsbury-Martinez et al. 2019 for an example of this top-down evolution -similar top down heating can be found in the 1400 K model as it warms back up from the initial cooling that occurred during model initialisation).This effect (i.e. a switch from radiative to advective dynamics) can also be seen when comparing the vertical-advective and global-mean-radiative timescales: as we move deeper into the atmosphere, the dynamics switch from being radiatively dominated to advectively driven, at around the same pressure as the deep adiabat forms.However, it is important to note that this is a 1D view of an inherently 3D problem -between the tidally located nature of the planetary irradiation (i.e. the xied day-side and nigth-side), and the strong longitude and latitude dependence of the vertical winds, the exact pressure at which the atmosphere changes dynamical regimes is likely to be highly localised.Yet it is reassuring to confirm that, on a global scale, the regime transition occurs about where we would expect and as required for our mechanism to work.On the other hand, early outputs of the 2500 K model reveal, as expected, a strong enthalpy upflow throughout most of the deep atmosphere, although as the simulation evolves and the deep atmosphere finishes cooling, this slowly evolves towards the deep heating seen in the 1000 K and nominal models.Hints of this evolution towards deep heating can be seen around 1 bar where a weaker net downflow has started to develop.As shown in Table 1, the global steady-state vertical enthalpy flux is generally independent of the initialisation temperature.That is to say that, given enough time, almost all the models here should settle onto the same steady-state profile, with the initialisation temperature only affecting the time taken to reach that profile.The only exception to this rule is the model in which we have modified the deep Rayleigh-drag.
As shown in Table 1, the model with slower deep Rayleigh-drag exhibits a significantly stronger peak and importantly mean vertical enthalpy flux than the models with fast drag ( drag = 1 day), even when models are compared at the same point in time (∼ 50, 000 days).This difference in vertical heating rate, and hence the temperature of the steady-state deep atmosphere, can be understood through using the vertical advective timescale ( adv ∼    ): if we consider the scale height., to be on the order of the radius WASP-76b and the velocity to be the global mean vertical velocity (  = 734 m s −1 ), we find that  adv ∼ 2.06 days (see Figure 4).I.e., for most of the models considered, the advective and drag time-scales in the deep atmosphere are of the same order of magnitude, leading to a noticeable reduction in the vertical wind speed, and hence vertical enthalpy flux in these models when comparing them with a no/low drag model in which the advective time-scale is significantly shorter than the drag time-scale (see Table 1).A similar effect can be seen in the zonal-mean zonal-wind, with the Low Drag model exhibiting a jet that extends significantly deeper into the atmosphere than the nominal model it is based upon.We discuss the implications of this result on the expected level of advective radius inflation in section 4.
Finally we compare, in Figure 5, the near-steady-state temperaturepressure profile of the nominal model with both a 1D model of the outer and deep atmosphere calculated using ATMO (see Tremblin et al. 2015 for an overview of the ATMO model) and an internalstructure model (which extends down to over 10 7 bar), based upon the work of Baraffe et al. (2003); Chabrier et al. (2004), of a hot Jupiter with a mass of 0.9M J .This internal-structure model is rather unique, as it is very difficult to generate a model with such a large radius.In order to do so, a large amount of thermal energy (corresponding to a luminosity of 2 × 10 28 erg s −1 ) must be deposited deep enough into the planetary interior to modify the internal adiabat (i.e.inflate the radius).As a consequence, the radius of the planet becomes essentially constant with time from early ages and the evolution is stalled (see Figure 4 of Chabrier et al. (2004)).Note that the input physics and equation of state of these internalstructure models differs from that considered in expeRT/MITgcm (typically GCMs use simpler equations of state for computational efficiency reasons, and because they are focused upon relatively lowdensity dynamics).As such, the adiabatic index of our models and the internal-structure models also differ, complicating a direct comparison between the deep atmospheric temperature-pressure profiles in the two models.Instead, in order to divine which structure-model is the closet match to our steady-state GCM model, and hence calculate the level of radius inflation exhibited, we follow standard practice and perform the model comparison at a fixed pressure of 100 bar (i.e. at a reference-pressure which is sufficiently deep so that the atmosphere is optically thick and hence either convectively or advectively driven).The result of this comparison is the selection of a internal-structure model with a radius of  = 1.98RJ being chosen as the best 'match' to our steady-state atmospheric model.This radius is broadly compatible with the observed radius of WASP-76b,  = 1.83 ± 0.06R J , suggesting that potential temperature advection alone is enough to explain the radius inflation of WASP-76b.A conclusion that is further reinforced by the partially evolved T-P profiles found in our alternate initialisation temperature models (see the solid lines in Figure 1).Despite the shorter run time of the alternative start models, Figure 1 clearly shows that all of the models are converging towards the same, inflated, deep T-P profile found in the nominal model, albeit at differ-ent rates due to differences in the efficiency of deep cooling versus heating (see Sainsbury-Martinez et al. 2019), i.e. the slow heating of the 1000 K model in Figure 1.This suggests that our conclusion of advection alone being sufficient to explain the radius inflation of WASP-76b is fairly robust.

DISCUSSION
In this work, we have performed additional analysis on extended and derivative versions of the WASP-76b models of Schneider et al. (2022b), focusing our analysis on the vertical advection of potential temperature, and its ability to heat the deep atmosphere with respect to 1D atmospheric models, leading to radius inflation with respect to these 1D models (as introduced by Tremblin et al. 2017 and explored, in a parametrised 3D model, by Sainsbury-Martinez et al. 2019, 2021).Importantly, thanks to the inclusion of a robust radiative transfer scheme (based upon petitRADTRANS) in expeRT/MITgcm, these models also allow us to complete the 'wish' of Sainsbury-Martinez et al. (2019): exploring the steady-state atmosphere of a hot Jupiter with a self-consistent radiative transfer scheme (in the outer atmosphere) so that a comparison between a atmospheric model and an internal-structure model can be made, thus quantifying, almost, the exact level of radius inflation that potential temperature advection alone can explain.
We started by exploring the zonal-mean zonal and meridional dynamics (Figure 2), with the aim of confirming the presence of a strong super-rotating jet that drives an equatorial downflow between the irradiated outer atmosphere and the advective deep atmosphere.This analysis was performed for six models, five of which have different initial deep adiabat temperatures ranging from significantly hotter to cooler than the expected steady-state deep atmosphere (see the dashed lines in Figure 1), and one which extends the nominal model of Schneider et al. (2022b), but with slower deep Rayleigh-drag, and which we include in order to explore the robustness of our results.For all five WASP-76b models with varying deep initialisation temperatures, we found that, once any deep atmospheric cooling had slowed/stopped, the strong super-rotating jet extends to  > 1 bar and drives a meridional circulation profile that includes a zonal-mean downflow that connects the radiative outer atmosphere with the advective deep atmosphere.This implies that high potential temperature fluid parcels from the outer atmosphere can indeed be transported vertically downwards, potentially heating the deep atmosphere.
Next, we explored if this was indeed the case, investigating how the mean vertical enthalpy advection ( H (, , ) =     (, , )) varies with both longitude and pressure (see Figure 3 and Table 1).This analysis revealed a number of trends which line up with the dynamics of the atmosphere.For example, for models that are initialised with an overly hot deep adiabat, and hence exhibit significant initial, deep cooling, the primary direction of enthalpy transport is from the deep to the outer atmosphere where it can be radiated away.However as such a model evolves, and the deep atmosphere cools towards (and maybe overshoots -an effect seen in the hot initialisation models of Sainsbury-Martinez et al. 2021) steady-state, we find that all of our models exhibit a net downwards flow of enthalpy, with the strength and pressure range of the downwards transport decreasing as the deep atmosphere very slowly equilibrates (a process that can take many hundreds to thousands of Earth years for  > 100 bar; Sainsbury-Martinez et al. 2019, 2021).Of course that is not to say that the vertical enthalpy transport lacks horizontal structure.As with the wind that drives it, differences in the vertical enthalpy transport are primarily linked with the differences in the day-side and night-side forcing, leading to a near global overturning circulation pattern that drives upwards vertical enthalpy transport on the day-side and downwards transport on the night-side, where divergent and wave driven circulations converge.We intend to explore the structure of the horizontal and vertical wind and enthalpy flux in more detail as part of a future study, including investigating how rotation impacts the dynamics (and hence may effect which hot Jupiters are inflated and which are not).Overall we find that, regardless of the initial conditions (i.e. with enough time), all of our fast drag models exhibit comparable peak and mean vertical enthalpy transport into the deep atmosphere.Furthermore this vertical enthalpy transport is also comparable, if not slightly stronger than that found in a reanalysis of the HD209458b models of Sainsbury-Martinez et al. (2019), reinforcing the idea that vertical potential temperature advection alone can explain the inflated radii of highly irradiated, gaseous, exoplanets.
We further confirm that this is the case via a comparison of our nominal WASP-76b models near-steady-state T-P profile (a T-P profile that all WASP-76b models appear to be converging towards -see Figure 1 -albeit at varying rates due to differences in the efficiency of cooling versus heating in the deep atmosphere) with an internal-structure model based upon the work of Baraffe et al. (2003); Chabrier et al. (2004).As shown in Figure 5, the closest match to the nominal model is an internal-structure model with a mass of 0.9M J and an inflated radius of 1.98R J , which is more than large enough to fully explain the observed radius of WASP-76b ( = 1.83 ± 0.06R J ).Note however that this comparison was performed by only considering the temperature at 100 bar (a fairly standard pressure at which atmospheric and internal-structure model comparisons are performed), a necessary approximation given the rather different adiabatic indexes found in our models and the internal-structure models considered here.Briefly, this difference occurs due to differences in the physics and specifically the equation of state considered in the models, with expeRT/MITgcm using a relatively simplified EOS (for both computational efficiency reasons as well as the GCMs focus upon modelling relatively low-density regions of the atmosphere) in comparison to that used in Baraffe et al. (2003); Chabrier et al. (2004).As such, an exact calculation of the level of radius inflation found in our model is beyond the current generation of GCMs, although work is in the pipeline to develop next-generation GCMs with updated dynamics and physics that will allow for even more robust comparisons with internal-structure models.However this does not mean that our calculation is without value, or that our results are far from the exact radius of our atmospheric model.For example, an internal-structure model with R = R J is a very poor fit to our atmospheric model with deep temperatures at 100 bar that are a order of magnitude cooler than than found with expeRT/MITgcm, reinforcing our inference that this model exhibits significant, advectively driven, radius inflation.However this is the not only effect that drives uncertainty in the exact level of radius inflation that advective heating can drive.For example, Mayne et al. (2019), showed that the dynamics of small-Neptunes and super-Earths varied significantly between models which solved the primitive equations of meteorology the the full Navier-Stokes equations.Other model choices can also affect the strength of the deep heating, such as the strength of any grid-scale smoothing (i.e. the inclusion of a Shapiro filter, which can affect the strength of the zonal jet and hence the vertical wind and advection -see Koll & Komacek 2018;Skinner & Cho 2021;Hammond & Abbot 2022), the atmospheric chemistry considered (e.g.equilibrium vs non-equilibrium chemistry), or the sources of opacity included (for example the inclusion of SiO, Fe and FeII opacity may affect atmospheric heating and the depth to which radiation penetrates, changing the T-P profile slightly.See, for example Lothringer et al. 2020).Here we investigate one of these possible sources of uncertainty: the inclusion, and thus strength, of deep Rayleigh-drag.This uncertainty can be seen by comparing the nominal model, with  drag = 1.0 days, with a model in which the deep Rayleigh-drag has been significantly slowed, such that  drag = 1000 days.i.e. a model in which the drag time-scale is significantly slower than the vertical advective timescale, which is of the order of 2 days for WASP-76b.Starting with the zonal-mean zonal-wind, our analysis indicates that the equatorial jet extends significantly deeper than in the nominal model.In turn, this drives stronger vertical mixing which results in a vertical enthalpy flux that is notably enhanced with respect to the nominal model.If we them compare the nominal model after 50,000 days with the low drag model after 40,000+10,000 days, we find that the deep T-P profile in the slow drag model is a little warmer, suggesting a slightly larger inflated radius.Comparing the vertical enthalpy flux at this time, confirms that the low drag model exhibits significantly enhanced deep heating.As such, and without a more complete understanding of how much, if any, Rayleigh-drag should be included in the deep atmosphere of hot Jupiter models, there will always remain an uncertainty on the exact level of radius inflation that vertical advection can drive.However, given that a) the Rayleigh-drag is confined to the highest pressure regions of the atmosphere (allowing for advective heat transport into the outer deep atmosphere, and then adiabatic mixing to carry heat deeper), and b) that the strength of the vertical advective transport is more than enough to explain the observed radius inflation, even in the nominal model with 'strong' drag, we are confident in our conclusion that the vertical advection of potential temperature alone is enough to explain the radius inflation of many hot Jupiters (and hot brown dwarfs), including WASP-76b.

CONCLUDING REMARKS
Overall, our analysis of the vertical mixing and vertical transport of potential temperature in an extended sample of the the WASP-76b models of Schneider et al. (2022b) has revealed that, contrary to their conclusions, the vertical advection of potential temperature alone is more than enough to explain the radius inflation of WASP-76b.This difference in conclusion arises for a number of reasons.The first is simply that the nominal model of Schneider et al. (2022b) was not run for long enough, and that their approach to avoid the computational expense of evolving a radiative GCM to steady-state in the deep atmosphere (i.e. the steroids model) made a number of assumptions about the deep dynamics which limit the applicability of such a extrapolative approach.Specifically, when extrapolating the evolution of their nominal models deep P-T profile, they focused on the evolution of the temperature at 650 bar, which, for the time frame they considered, revealed near exponential cooling.However, as shown in the isothermal-start model of Sainsbury-Martinez et al. (2019), advective heating of the deep atmosphere starts in the lower pressure regions of the deep atmosphere (i.e. at the bottom of the radiatively dominated outer atmosphere) and slowly pushes deeper with time, with the time to heat the atmosphere only increasing as the heating moves deeper and the local density increases.Evidence for this top down heating in the steroids models can be seen in Figure 2 of Schneider et al. (2022b), with slow heating occurring between ∼ 5 and ∼ 100 bar, leaving the region around 650 bar to appear steady and hence evolved.Here, by evolving the nominal model for an additional 69000 days of simulation time, we are approaching a true steady-state that is significantly hotter than the steroids model.Furthermore, when compared with an internal-structure model from Baraffe et al. (2003); Chabrier et al. (2004), this steady-state corresponds to a radius of 1.98R J , more than large enough to explain the observed, inflated, radius of WASP-76b ( = 1.83 ± 0.06R J ).The second reason for our difference in conclusion can be linked to the wide use of intrinsic/internal temperatures in the exoplanetary communities.Briefly, radius inflation is simply the difference between the observed radii of a hot Jupiter and a standard 'radiativeconvective' 1D model its outer atmosphere.This difference is believed to occur because 1D atmospheric models lack some fundamental physics that drive deep heating, with suggestion ranging from ohmic dissipation to vertical heat transport, and is 'fixed' (or accounted for) in 1D models by including an artificial, intrinsic/internal temperature meant to represent the heating of the deep atmosphere.Commonly this is linked with excess energy loss from the interior (hence the name internal temperature), however Tremblin et al. (2017) and Sainsbury-Martinez et al. (2019) proposed that this deep heating instead occurs due to vertical heat transport, with no need for any energy transport from the interior to the outer/deep atmosphere (i.e.zero net deep flux).In essence, this intrinsic temperature acts as a 'fudge' factor designed to allow for direct comparisons between observations (such as transmission spectra) and 1D models, and relying upon it outside of those scenarios can lead to either over or under (as was the case in Schneider et al. 2022b) estimation of the level of radius inflation.Instead, as done here, comparisons must be done with internalstructure models, even when the accuracy of those comparisons is limited by the different equations of state used (i.e. by the simplified EOS used in GCMs -although work is in progress to change this).
Of course, many questions remain about the exact level of radius inflation that vertical advection can drive, and if it can fully explain the differences seen in radius inflation for the broader hot Jupiter community, including those unusual objects that are very highly irradiated and yet show little to no sign of inflation (for example WASP-43b or WASP-18b).expeRT/MITgcm now makes a radiatively robust study of these objects possible for the first time, and we look forward to the results of future work with this, and other next-generation, models.However there is now no doubt that potential temperature advection provides a robust explanation for some if not all of the observed radius inflation of hot Jupiters and hot brown dwarfs, and as such changes to how future GCM studies are performed are recommended.Previously, it has been recommended that future GCM studies of hot Jupiters be initialised with a adiabat at the bottom of their simulation domain and then be allowed to evolve to steady-state (Sainsbury-Martinez et al. 2019).However this remains computationally expensive and can lead to mistakes when models are not allowed to evolve sufficiently.As such, given how well potential temperature advection alone can explain the inflated radii of hot Jupiters, we suggest that future studies should initialise their deep atmosphere with an adiabat based upon the best fitting internal-structure model that corresponds to the inflated radii, albeit modified to match the adiabatic index of the GCM.

Figure 2 .Figure 3 .
Figure2.The zonal-mean zonal wind (left) and meridional circulation streamfunction (right) for the nominal GCM model of WASP-76b.In the zonal wind profile, easterly winds are positive and westerly winds are negative, whilst in the meridional circulation profile, we plot the streamfunction using a log scale in order to clearly illustrate the full circulation profile, especially in the outer atmosphere.Here, clockwise circulations are shown in red and anti-clockwise in blue -these circulations combine to reveal an equatorial upwelling in the outer atmosphere driven by the strong day-side irradiation, and an equatorial downflow in the deep atmosphere, which is linked with the downward advection of potential temperature.

Figure 4 .Figure 5 .
Figure 4. Comparison of the peak vertical advective time-scale (  = /  , where  is the atmospheric scale height and   is the maximum downward velocity) and the global-mean radiative timescale for the near steady-state nominal model.Note how, despite both timescales increasing with pressure, the rapid increase in optical depth means that advection dominates over radiation in the deep atmosphere.

Table 1 .
Schneider et al. 2022bb)l mean vertical enthalpy flux for five WASP-76b models in which either the temperature of the initial deep adiabat or the strength of the deep drag have been changed, along with the nominal model presented inSchneider et al. (2022b).Note: the Low Drag model has been run for 10,000 additional days using a snapshot of the nominal model after 40,000 days (i.e. the 'evolved' model ofSchneider et al. 2022b) of simulation time as an initial condition.Further, the mean vertical enthalpy flux for the nominal model at an equivalent timestep to that of the Low Drag model remains essentially unchanged.