Atmospheric Circulations of Hot Jupiters as Planetary Heat Engines
DDraft Modified May 14, 2018
Preprint typeset using L A TEX style emulateapj v. 12/16/11
ATMOSPHERIC CIRCULATIONS OF HOT JUPITERS AS PLANETARY HEAT ENGINES
Daniel D.B. Koll , and Thaddeus D. Komacek , Department of Earth, Atmospheric, and Planetary Sciences, Massachusetts Institute of Technology, Cambridge, MA, 02139;[email protected] Lunar and Planetary Laboratory and Department of Planetary Sciences, University of Arizona, Tucson, AZ, 85721;[email protected] The authors contributed equally to this work, and are listed alphabetically.
Draft Modified May 14, 2018
ABSTRACTBecause of their intense incident stellar irradiation and likely tidally locked spin states, hot Jupitersare expected to have wind speeds that approach or exceed the speed of sound. In this work we developa theory to explain the magnitude of these winds. We model hot Jupiters as planetary heat enginesand show that hot Jupiters are always less efficient than an ideal Carnot engine. Next, we demonstratethat our predicted wind speeds match those from three-dimensional numerical simulations over a broadrange of parameters. Finally, we use our theory to evaluate how well different drag mechanisms canmatch the wind speeds observed with Doppler spectroscopy for HD 189733b and HD 209458b. Wefind that magnetic drag is potentially too weak to match the observations for HD 189733b, but iscompatible with the observations for HD 209458b. In contrast, shear instabilities and/or shocks arecompatible with both observations. Furthermore, the two mechanisms predict different wind speedtrends for hotter and colder planets than currently observed. As a result, we propose that a widerrange of Doppler observations could reveal multiple drag mechanisms at play across different hotJupiters.
Subject headings: hydrodynamics — methods: analytical — methods: numerical — planets and satel-lites: atmospheres — planets and satellites: individual (HD 189733b, HD 209458b) INTRODUCTION
Hot Jupiters provide a unique laboratory for testingour understanding of planetary atmospheres. Showman& Guillot (2002) were the first to consider the atmo-spheric circulations of these planets. Using numeri-cal simulations, Showman & Guillot predicted that hotJupiters should develop strongly superrotating equato-rial jets, with wind speeds up to several kilometers persecond. This prediction was confirmed by subsequent ob-servations which showed that the thermal emission peakon many hot Jupiters is shifted eastwards from the sub-stellar point, consistent with heat being advected down-wind by a superrotating jet (e.g., Knutson et al. 2007;Crossfield et al. 2010).More recent observations have started to directly con-strain the wind speeds of these jets. High-resolutiontransmission spectra have found Doppler shifts in molec-ular absorption lines for HD 209458b (Snellen et al. 2010)as well as HD 189733b (Wyttenbach et al. 2015; Louden& Wheatley 2015; Brogi et al. 2016). The significant( ∼ several km s − ) blueshifts detected for both planetsimply rapid dayside-to-nightside winds that are broadlyconsistent with the wind speeds predicted by a range ofnumerical simulations (Showman & Guillot 2002; Show-man et al. 2009; Heng et al. 2011a; Showman et al. 2013;Komacek et al. 2017).Although it is qualitatively understood why hotJupiters develop equatorial jets, there is still no generaltheory that explains the jets’ magnitude. Hot Jupitersare very likely tidally locked. This orbital spin state cre-ates a strong day-night forcing which excites standingwaves that flux angular momentum towards the equatorand drive equatorial superrotation (Showman & Polvani 2011). The strength of superrotation should thereforedepend on the ratio between horizontal wave propaga-tion and radiative cooling timescales (Koll & Abbot 2015;Komacek & Showman 2016; Zhang & Showman 2017).This basic expectation is complicated, however, by re-sults which show that the jet’s state depends on bothhorizontal standing waves and vertical eddies (Tsai et al.2014; Showman et al. 2015), and it is still unclear how thetwo mechanisms jointly determine the jet’s magnitude.In this paper we constrain the wind speeds of hotJupiters by modeling their atmospheric circulations asplanetary heat engines. The utility of this approach haspreviously been demonstrated for hurricanes on Earth(Emanuel 1986) and rocky exoplanets (Koll & Abbot2016). Atmospheric circulations can be considered heatengines because parcels of fluid tend to absorb heat at ahigh temperature (e.g., on the dayside of a hot Jupiter)and emit heat at a low temperature (on the nightside).The differential heating and cooling allows parcels to gen-erate work, and thus kinetic energy, which in steady statehas to be balanced by the dissipation of kinetic energyvia friction.In contrast to hurricanes and the atmospheres of rockyexoplanets, however, it is still poorly understood how hotJupiters dissipate kinetic energy (Goodman 2009). Po-tential mechanisms include magnetic drag in partiallyionized atmospheres (Perna et al. 2010; Menou 2012;Rauscher & Menou 2013; Rogers & Showman 2014),shocks in supersonic flows (Li & Goodman 2010; Heng2012; Perna et al. 2012; Dobbs-Dixon & Agol 2013; Fro-mang et al. 2016), and turbulence induced by fluid in-stabilities such as the Kelvin-Helmholtz instability (Li &Goodman 2010; Fromang et al. 2016). a r X i v : . [ a s t r o - ph . E P ] D ec Koll & KomacekOur goal is to evaluate these proposed mechanisms andto test which of them are able to match current obser-vations. To do so we first describe our numerical simu-lations (Section 2). Next, we develop the heat engineframework and test it with the numerical simulations(Section 3). Finally, we apply our framework to observa-tions (Section 4) and state our conclusions (Section 5).Our results show that current observations favor shearinstabilities and/or shocks as the dominant drag mech-anism for HD 189733b, and motivate extending similarobservations across a wider range of planets. NUMERICAL SIMULATIONS
We compare our theory with the GCM simulationsthat were previously described in Komacek et al. (2017).In summary, the simulations use the MITgcm (Adcroftet al. 2004) to solve the atmospheric fluid dynamics equa-tions coupled to double-gray radiative transfer with plan-etary parameters relevant for a typical hot Jupiter, HD209458b. The double-gray approximation divides thespectrum into an incoming collimated and a thermal dif-fuse part. The absorption coefficients were chosen tomatch more detailed radiative transfer calculations; theabsorption coefficient for incoming stellar radiation is setto a uniform value, κ SW = 4 × − m − kg − , the ther-mal absorption coefficient varies approximately with thesquare root of pressure, κ LW = 2 . × − m − kg − × ( p/ . , where the power-law exponent comesfrom fitting the analytic model of Parmentier & Guillot(2014) and Parmentier et al. (2015) to radiative transfermodels with realistic opacities. With these values thephotosphere (where the optical thickness equals unity)for stellar radiation lies at about 0 .
23 bar and the pho-tosphere for thermal radiation lies at 0 .
28 bar.The model’s resolution is C32 in the horizontal(roughly corresponding to a global resolution of 128 × τ drag . Simulations with τ drag ≤ s use a timescale thatis spatially uniform. Simulations with τ drag > s ad-ditionally include a “basal” drag term that allows themodel to equilibrate within reasonable integration times.The basal drag strength increases as a power-law withpressure, from no drag at 10 bar to a timescale of 10days at 200 bar (Komacek & Showman 2016). Second,to enforce numerical stability, the model includes a afourth-order Shapiro filter that damps wind and temper-ature variations at the model grid scale. The Shapirofilter acts as numerical drag at small spatial scales and,in simulations without any other sources of drag, even-tually helps to equilibrate the kinetic energy of the flow.The potential issue with relying on numerical drag is thatit relies on parameters which are generally chosen formodeling convenience, not because they are physically − − ( dK/dt ) num , rms / ( dK/dt ) Rayleigh , rms − − − − P [ b a r s ] (a) τ drag s10 s10 s10 s10 s ∞ − − − − ( dK/dt ) rms [J kg − s − ]10 − − − − P [ b a r s ] (b) RayleighNumerical
Fig. 1.—
Kinetic energy dissipation in many of our GCM simu-lations is dominated by numerical drag. Panel (a) shows the ratiobetween the global root-mean-square rate of kinetic energy dissi-pation by numerical drag, ( dK/dt ) num , rms , versus the global root-mean-square rate of kinetic energy dissipation by explicit Rayleighdrag, ( dK/dt ) Rayleigh , rms , as a function of pressure for modelswith an equilibrium temperature of T eq = 1500 K. The coloredlines show simulations with different Rayleigh drag timescales, withdarker lines representing longer drag timescales. The dashed ver-tical line shows the divide between dissipation dominated by nu-merical drag (to the right of the line) and Rayleigh drag (to theleft). Except for short Rayleigh drag timescales, τ drag ≤ s, nu-merical dissipation dominates. Note that the case with τ drag = ∞ still includes basal drag, so the ratio of numerical to Rayleigh dragdissipation is not infinite at depth. Panel (b) shows the absolutecontribution of Rayleigh drag and numerical effects to the kineticenergy dissipation. Only a subset of the simulations are shown forvisual convenience. The dissipation rate increases with decreasingpressure, largely due to the stronger wind speeds at lower pressures. motivated. This raises the question of which source ofdrag is dominant in our simulations. − − − − ( dM/dt ) rms / ( dM/dt ) Coriolis , rms − − − − P [ b a r s ] τ drag s, Rayleigh10 s ∞ Numerical
Fig. 2.—
Numerical effects are small relative to physical termsin the zonal angular momentum budget of our simulations. Thisplot shows the global root-mean-square of the change in zonal mo-mentum due to Rayleigh drag (solid lines) and numerics (dashedlines) relative to the change in zonal angular momentum due to theCoriolis force (i.e. rotation). Plots have the same color scheme asin Fig. 1, for visual convenience we only show a subset of all simu-lations. The acceleration from numerics is smaller than either theCoriolis force (if τ drag ≥ s) or Rayleigh drag (if τ drag ≤ s).As a result, numerics do not significantly affect the angular mo-mentum budget of our simulations. We find that numerical drag can play a key role in ourGCM simulations. Although the potential importanceof numerical drag has repeatedly been pointed out in thehot Jupiter literature (Goodman 2009; Li & Goodman2010; Thrastarson & Cho 2010; Heng et al. 2011a; Liu& Showman 2013; Mayne et al. 2014; Polichtchouk et al.2014; Cho et al. 2015), no work has previously quantifiedits effect relative to explicitly parametrized drag. Fig-ure 1 compares the rates at which our GCM is dissipat-ing kinetic energy via numerical drag from the Shapirofilter versus the dissipation rate due to Rayleigh dragas a function of pressure. Figure 1 (a) shows the rela-tive global root-mean square-dissipation due to numeri-cal drag versus Rayleigh drag, while Figure 1 (b) showsthe absolute global root-mean-square value of kinetic en-ergy dissipated by both drag mechanisms. We com-pute the root-mean-square change in kinetic energy as( ∂K/∂t ) rms = (cid:104) ( ∂K/∂t ) (cid:105) / , where the angle bracketsdenote an area average. We find that all simulations withmoderately long Rayleigh drag timescales, τ drag ≥ s,dissipate most kinetic energy through numerical drag.Moreover, even in the simulations with the strongestRayleigh drag (yellow curve in Fig. 1a,b) numerical dragdominates the dissipation of kinetic energy near the topand bottom of the model domain. Although the modelincludes a basal drag, we find that it contributes lesstowards kinetic energy dissipation than numerical dragnear the bottom of the domain. This is likely due tothe Shapiro filter acting as a sponge for waves that areexcited in the upper atmosphere. However, wind speedsat pressures greater than 10 bar are small so kinetic en-ergy dissipation near the domain bottom contributes rel-atively little to the overall dissipation (see Fig. 1b).Though numerical drag is a dominant factor in how our GCM dissipates kinetic energy, atmospheric circula-tions additionally depend on how the GCM resolves theangular momentum budget. We do not expect a priori that numerical effects will dominate the global angularmomentum budget, because the Shapiro filter is designedto not affect large-scale flow (Shapiro 1971). To checkthis insight, we explicitly compute the change in zonalangular momentum by numerics and Rayleigh drag as inPeixoto & Oort (1992): ∂M∂t = ∂u∂t a cos( φ ) . (1)In Eqn. (1) M is the zonal angular momentum per unitmass, ∂M/∂t is the rate of change of angular momentumwhich we compute in our simulations from the accelera-tion ∂u/∂t due to the Shapiro filter or Rayleigh drag, a is the planetary radius, and φ is latitude. Rayleigh dragalways acts as a sink of angular momentum whereas theShapiro filter can accelerate parts of the flow so we com-pare both terms via the root-mean-square change in mo-mentum, ( ∂M/∂t ) rms = (cid:104) ( ∂M/∂t ) (cid:105) / , where the anglebrackets as before denote an area average.We find that numerical effects play a relatively minorrole in the zonal angular momentum budget. Figure 2shows the change in angular momentum from numericsand Rayleigh drag relative to the change in angular mo-mentum from the Coriolis force, as a function of pres-sure. We compare both terms against the Coriolis forcebecause it is a small term in the zonal momentum budgetof hot Jupiters due to their slow rotation and winds thatpeak at the equator (Showman & Polvani 2011; Showmanet al. 2015). In relative terms, we find that the numer-ical change in angular momentum becomes larger thanRayleigh drag once τ drag > s (blue curves). However,in absolute terms, the momentum change from numericsremains one to two orders of magnitude smaller than theCoriolis term at most pressure levels. We conclude thatnumerical effects likely do not play a dominant role inthe angular momentum budget of our simulations.Given that many published simulations of hot Jupitersdo not include Rayleigh drag, our results indicate thatmany of these simulations rely on numerical drag to equi-librate kinetic energy. Further work is needed to ensurethat this kind of dissipation in hot Jupiter GCMs is phys-ically motivated and that its effects are robust with re-spect to changes in numerical parameters. At the sametime, the angular momentum budget in our simulationsis not dominated by numerics. We therefore expect thatGCMs are robust in simulating the qualitative featuresof hot Jupiter circulations (e.g., equatorial jets), but thatthe absolute kinetic energy and thus wind speeds in thesesimulations might be affected by numerical details. Ourresults agree with previous work, which has shown thatthe equilibrated flows in hot Jupiter GCMs largely con-serve angular momentum, are independent of initial con-ditions, and that the magnitude of winds is only weaklysensitive to changes in numerical parameters (e.g. Henget al. 2011; Liu & Showman 2013; Mayne et al. 2014). Inthe remainder of this paper we focus on existing GCMsto test our theoretical framework. To do so we developa theory in the next section that can account for bothexplicit and numerical drag. HOT JUPITERS AS HEAT ENGINES
Koll & Komacek − − − Temperature [K] P r e ss u r e [ ba r s ] DaysideNightsideAdiabat a bcd
Fig. 3.—
A diagram of the Ericsson cycle, overlaid on dayside-and nightside-averaged temperature profiles of a reference simula-tion ( T eq = 1500 K , τ drag = 10 s) and an adiabatic profile. TheEricsson cycle works as follows: a parcel of fluid starts at depth onthe nightside (a), moves towards the dayside (b), where it rises (c),moves back towards the nightside (d), and sinks (a). We assumethat rising and sinking motions (b-c, d-a) are isothermal and thatmotions between hemispheres (a-b, c-d) are isobaric. The isother-mal assumption is motivated by the GCM profiles, which showthat hot Jupiters are much closer to vertically isothermal than toadiabatic. In steady state, the rate W at which a heat engineperforms work is given by W = ηQ, (2)where η is the engine’s thermodynamic efficiency and Q is the rate at which the engine absorbs heat.First, the heating rate Q is equal to the average ab-sorbed stellar flux, Q = σT eq , (3)where T eq is the planetary equilibrium temperature.Second, we constrain the work output rate W . Weassume that work goes entirely towards generating anddissipating kinetic energy. If Rayleigh drag dominates,the rate at which kinetic energy is dissipated equals W Rayleigh = (cid:90) dpg × (cid:28) v τ drag (cid:29) , (4)where v is the velocity vector and the angle bracketsdenote an area average. If numerical drag dominates,kinetic energy is dissipated by the Shapiro filter whichdamps the highest wavenumber components of the flow.Because the highest wavenumber in the GCM is set bythe model’s grid spacing ∆ x we scale the Shapiro filter’sdamping timescale as τ ∼ ∆ x/U . This means the rate atwhich numerical drag dissipates kinetic energy is equalto W num ∼ U ∆ x/U × pg = U ∆ x × pg . (5)Third, we constrain the efficiency η . Previous workon hurricanes and the atmospheres of rocky planets con-strained this quantity by modeling atmospheric circu-lations as Carnot cycles (Emanuel 1986; Koll & Abbot 2016). Unfortunately it is difficult to argue that hotJupiters should also resemble Carnot cycles. In a Carnotcycle parcels of fluid expand and contract adiabaticallybetween heating and cooling. This model is physicallymotivated by the fact that hurricanes and rocky planetsundergo convection, so fluid parcels move rapidly andquasi-adiabatically. In contrast, the upper atmospheresof hot Jupiters are strongly irradiated by their host stars.The irradiation creates a stable stratification and sup-presses convection, which means the vertical tempera-ture structure is approximately in radiative equilibriumand lapse rates are small (Iro et al. 2005; Guillot 2010).As the temperature profiles from a reference simulationin Figure 3 illustrate, temperatures are indeed far fromadiabatic, which underlines that the Carnot cycle is apoor model for hot Jupiters.Here we constrain the efficiency η by modeling hotJupiters as Ericsson cycles (McCulloh 1876). The Er-icsson cycle is shown in Figure 3: a parcel of fluidstarts deep in the nightside atmosphere (Fig. 3, pointa). It moves at constant pressure towards the dayside(b), where the stellar heating causes it to rise (c). Theparcel then moves to the nightside (d), before cooling andsinking back to its starting position (a). Even though theassumption of isothermal vertical motions is an idealiza-tion, Figure 3 shows that the Ericsson cycle provides aphysically motivated model for hot Jupiters.The efficiency of the Ericsson cycle is given by η = (cid:72) δQ (cid:82) ca δQ = (cid:72) T ds (cid:82) ca T ds . (6)Here δQ is a change in a parcel’s heat content, and ds is achange in entropy. From the first law of thermodynamics, T ds = c p dT − dpρ = c p dT − RT d ln p, (7)where we have used the ideal gas law in the second step.We can then evaluate the numerator (cid:72) T ds as (cid:90) ba c p dT − (cid:90) cb RT d ln p + (cid:90) dc c p dT − (cid:90) ad RT d ln p, = c p ( T day − T night ) − RT day ln( p lo /p hi )+ c p ( T night − T day ) − RT night ln( p hi /p lo ) , = R ( T day − T night ) ln( p hi /p lo ) . (8)Similarly the denominator (cid:82) ca T ds in Eqn. (6) is (cid:90) ba c p dT − (cid:90) cb RT d ln p, = c p ( T day − T night ) + RT day ln( p hi /p lo ) . (9)The ratio of these two terms gives the efficiency, whichwe write as η = T day − T night T day × ln (cid:2) ( p hi /p lo ) R/c p (cid:3) T day − T night T day + ln (cid:2) ( p hi /p lo ) R/c p (cid:3) . (10)Importantly, the efficiency η is always lower thanthe efficiency of a Carnot cycle, η Carnot = ( T day − T night ) /T day , which is the maximum efficiency a heat en-gine can reach. The lower efficiency arises because heatis radiated to space as a parcel passes from the dayside Scaled heat input τ drag ησT eq [J/m ]10 U r m s f r o m G C M [ m / s ] Rayleigh drag ∝ ( τ drag ησT eq ) / Numerical drag ∝ ( ησT eq ) / (a) 10 U predicted by combined scaling [m/s]10 U r m s f r o m G C M [ m / s ] (b) τ drag ∞ sec10 sec10 sec10 sec10 sec U from Eq.11U from Eq.12
Fig. 4.—
Our heat engine scaling captures the strength of wind speeds across a wide range of hot Jupiter GCMs. (a) Hot Jupitersimulations fall into two regimes in which bulk wind speeds either scale following Rayleigh drag or numerical drag (black lines show the twodifferent slopes). (b) Our combined scaling predicts the GCM wind speeds in both regimes. The y-axis corresponds to the root-mean-squarewind speed averaged over pressures less than 1 bar in different hot Jupiter simulations, the x-axis is the wind speed predicted from Equation13. Each dot represents a different GCM simulation with varying T eq = 500 − τ drag = ∞ , we use τ drag = 5 × s in the left panel instead. to the nightside (c-d). If, instead, this heat could bestored and used later to heat up the parcel as it passesback from the nightside to the dayside (a-b), the Ericssoncycle’s efficiency would equal that of a Carnot cycle .As an example we consider the efficiency of WASP-18b, whose phase curve is consistent with zero heat re-distribution from dayside to nightside (Maxted et al.2013). We assume that a parcel of fluid moves twoscale heights in the vertical every time it traverses theplanet horizontally so ln (cid:2) ( p hi /p lo ) R/c p (cid:3) ∼ R/c p . Inthis case WASP-18b’s Carnot efficiency would be unity, η Carnot = 1, whereas its actual efficiency is smaller bya factor of three, η = 0 .
36. Hot Jupiters can thereforebe thought of as comparable to, but less efficient than,ideal Carnot engines. Their efficiency can be reducedeven further by molecular diffusion and irreversible phasechanges (Pauluis & Held 2002), so Equation 10 should beconsidered an upper limit.We are now able to test the extent to which hotJupiters resemble heat engines. A key prediction of ourtheory is that wind speeds are sensitive to whether windsare damped by Rayleigh drag or numerical drag. Basedon Equations 2-5, we expect that winds should scale asthe square root of the modified heat input for Rayleighdrag, U ∝ ( τ drag ησT eq ) / , whereas they should scale asthe one-third power of the heat input for numerical drag, If the heat lost during (c-d) could be captured and used toheat the parcel during (a-b), then Eqn. (9) becomes (cid:82) ca δT ds = (cid:82) cb δT ds = RT day ln( p hi /p lo ) and Eqn. (10) becomes η = ( T day − T night ) /T day . A parcel travels a vertical distance d z ∼ WaU , where a isthe planet radius and W the vertical wind speed. Using char-acteristic values from a simulation with T eq = 1500 K and nodrag, W ∼ − , U ∼ ms − , and a = a HD209458b , we find d vert ∼ H , where H is the scale height. In agreement with thisestimate, we mapped streamfunctions in our simulations and foundthat the vertical extent of both zonal and meridional flows is nor-mally confined to ∼ − U ∝ (∆ xησT eq ) / . To compare both scalings in a singleplot and because ∆ x depends on numerical parameterswe first use the quantity τ drag ησT eq .Figure 4(a) shows that our simulations indeed exhibita dichotomy between Rayleigh and numerical drag. Thex-axis shows the scaled heat input τ drag ησT eq while they-axis shows the root-mean-square wind speed, U rms =( p − (cid:82) (cid:104) u + v (cid:105) dp ) / , where u and v are the zonal andmeridional wind speeds and where we average horizon-tally and over the meteorologically active region above p = 1 bar (see Fig. 1). To evaluate η we use the day-side and nightside brightness temperatures that wouldbe seen by an observer and assume that a parcel crossestwo scale heights, ln (cid:2) ( p hi /p lo ) R/c p (cid:3) ∼ R/c p .We find that wind speeds in most strongly damped sim-ulations with τ drag ≤ s increase according to Rayleighdrag (Fig. 4a). In contrast, winds in simulations with τ drag ≥ s increase more slowly and approximately fol-low the one-third slope predicted for numerical drag. Anotable exception to the Rayleigh scaling is given by thehottest simulations with τ drag = 10 s (yellow dots), inwhich winds increase with a one-thirds slope instead.This is due to the relative increase of numerical dissi-pation in strongly damped simulations. At τ drag = 10 swinds are so weak that Rayleigh drag, which is propor-tional to wind speed, becomes small relative to numer-ical drag in parts of the model domain. Similarly, ournumerical scaling performs worst for simulations with τ drag = 10 s (purple dots), in which wind speeds flat-ten out at high T eq even though the heat input keepsincreasing. Given that our theory performs well in thestrongly damped limit, deviations from it are likely dueto inaccuracies in our numerical scaling, which we discussbelow.We now constrain the wind speeds inside a hot Jupiteratmosphere. If the atmospheric circulation is primarilybalancing Rayleigh drag then wind speeds should scale Koll & Komacekas U Rayleigh = k (cid:18) τ drag ησT eq gp (cid:19) / , (11)whereas if the circulation is balancing numerical dragthen winds should scale as U num = k (cid:18) ∆ xησT eq gp (cid:19) / . (12)Here k and k are fitting constants of order unity thataccount for various approximations, in particular our as-sumption that temperature profiles are isothermal. Weuse k = 0 . k = 1 . T eq = 3000 K with τ drag = 10 s and τ drag = ∞ , re-spectively. We combine Eqns. 11 and 12 by demandingthat a GCM’s work output equals whichever is stronger,Rayleigh or numerical drag, so U = min( U Rayleigh , U num ) . (13)To evaluate Eqn. (12) we use the model’s grid spacingat the equator ∆ x ∼ πa/ a is the planetaryradius.We find that our theory matches the GCM simulationswell. Figure 4(b) compares our predicted winds with thesimulated root-mean-square wind speeds U rms , definedabove. As in Figure 4(a), we find that our scaling worksbest in the strongly damped limit, particularly for thesimulations with τ drag = 10 − s which our scalingmatches to better than 33%. These are also the simula-tions in which numerical drag is not dominant yet, andfor which we scale winds using Eqn. (11).Our scaling additionally matches the weakly dampedsimulations that are dominated by numerical drag( τ drag > s), even though the fit is less good thanin the strongly damped regime. This is likely due tothe approximations we made in deriving Eqn. (12). Totest this point we performed additional simulations inwhich we varied the model resolution and timestep. Wefound that Eqn. (12) over-predicts the sensitivity of windspeeds to numerical resolution (see Appendix). Furtherwork is needed to understand exactly how hot Jupitersimulations equilibrate through numerical drag.Nevertheless, given that our scaling captures the ba-sic dependence of wind speeds on a planet’s heat input(Fig. 4a) and additionally matches the GCM to betterthan a factor of two even when the models are domi-nated by numerical drag (Fig. 4b), we argue that themain shortcoming in Figure 4 is due to our imperfect de-scription of numerical drag, not due to the heat engineframework. We therefore sidestep the intricacies of nu-merical simulations and in the last section apply the heatengine framework directly to data. EVALUATING DRAG MECHANISMS WITHOBSERVATIONS
In this section we use the heat engine framework topredict how strong winds would have to be to balancethe two main proposed drag mechanisms on hot Jupiters,namely magnetic drag and shear instabilities. We thenevaluate our predictions by comparing them to observedwind speeds obtained from Doppler spectroscopy.For magnetic drag we combine Eqn. (11) with a kine-matic scaling for the effective Lorentz drag timescale U [ m / s ] Magnetic drag B eq [K]10 U [ m / s ] speed of sound H D b H D b Shear instability L = HL = 2 πa Fig. 5.—
Top: solid lines show the predicted wind speeds fromEqn. (11), assuming dissipation is caused by magnetic drag. Col-ored envelopes indicate that our theoretical scalings are subject touncertainty. The uppermost line for each magnetic field strengthshows the wind speed predicted for dissipation occuring at 1 bar,the lower line shows the wind speed predicted for dissipation oc-curing at 10 − bar, and the colored envelope shows intermediatepressures. Dots show wind speeds constrained via Doppler spec-troscopy for HD 189733b and HD 209458b (Snellen et al. 2010;Louden & Wheatley 2015). Bottom: solid lines show the pre-dicted wind speeds from Eqn. (15), assuming dissipation is causedby shear instabilities. Colored envelopes here indicate our esti-mated uncertainty for our heat engine scaling (see text). Windsfaster than the speed of sound (dashed black line a ) can also de-velop shocks. Magnetic drag can match both observations, butdoing so requires a large dipole field ( (cid:38) a We assume solar composition and that atmospheric temperatureis equal to the equilibrium temperature. (Perna et al. 2010; Menou 2012; Rauscher & Menou2013). To be consistent with Section 3, we use k = 0 . τ mag = 4 πH e ρB , (14)where B is the dipole field strength, H e the atmosphericelectrical resistivity, and ρ the gas density. The electricalresistivity is inversely related to the ionization fraction x e , H e ∝ √ T /x e , where x e is calculated from the Sahaequation (Perna et al. 2010). For hot Jupiters the ion-ized gas is largely potassium, for which we assume a so-lar abundance . We expect that most dissipation occurs For a planet with the equilibrium temperature of HD 209458b, somewhere between the upper levels probed by Dopplerobservations ( ∼ − bar) and the photosphere, so wecalculate winds over the range 10 − ≤ p ≤ L and damp the flow over a timescale L/U , sowind speeds scale as U shear = k (cid:18) LησT eq gp (cid:19) / . (15)For consistency we use k = 1 .
1, as in Section 3. We notethat Doppler observations probe the upper atmosphereonly whereas our theory constrains large-scale dissipationand thus should be representative of the bulk flow. Ob-servable wind speeds could potentially deviate from thebulk flow in atmospheres with large vertical shear. Nev-ertheless, we expect that the comparison between ourtheory and observations is warranted, given that a widerange of hot Jupiter GCMs produce equatorial jets thatare strongly vertically coherent (Showman et al. 2009;Heng et al. 2011a; Liu & Showman 2013; Mayne et al.2014; Polichtchouk et al. 2014; Cho et al. 2015).Figure 5 compares the observed wind speeds of1 . +0 . − . km s − for HD 189733b (Louden & Wheatley2015) and 2 ± − for HD 209458b (Snellen et al.2010) with our theoretical predictions for the two dragmechanisms . To indicate that our scalings aren’t ex-act, the colored envelopes in Figure 5 reflect the domi-nant sources of uncertainty in our scalings. For magneticdrag the uncertainty is dominated by the pressure atwhich dissipation is assumed to occur, for shear instabil-ities we use the remaining mismatch between theory andGCM simulations in Section 3. Because the magneticdrag timescale is relatively sensitive to both temperatureand pressure we additionally explored the impact of dif-ferent pressure-temperature profiles, and find that mostfeatures in Figure 5 are robust (see Appendix).First, we find that the observations for HD 189733bcan only be matched with a very strong dipole field of ∼ (cid:38) ∼
50G (R. Yadav, per-sonal communication). We conclude that magnetic drag T eq = 1450 K, the ionization fraction is x e = 4 . × − ,which is much smaller than the solar abundance of neutral potas-sium and thus consistent with the approximations made in Pernaet al. (2010). The corresponding magnetic resistivity is H e =2 . × cm s − . Note that these are 1 σ error bars and the detection itself wasonly significant at 2 σ . We assume p = 1 bar, η = 0 .
2, and g = 23m s − , with thelast two values motivated by the phase curve amplitude and mass-radius measurements of HD 189733b. We conservatively use 100% uncertainty (a factor of two) forwinds predicted with Eqn. (15). is a plausible drag mechanism for HD 209458b. In ad-dition, given the potentially large uncertainties in boththe Lorentz drag timescale (Equation 14) and dynamoscaling laws, magnetic drag cannot be ruled out for HD189733b, even though the required field strengths wouldbe larger than currently expected. Further theoreticalwork could help reduce these uncertainties. Our resultthat Lorentz forces are potentially unimportant for HD189733b but may be important for HD 209458b thereforeagrees with previous estimates that magnetic drag couldbecome significant at T eq (cid:38) T eq , in agreement with the observations. We also findthat the vertical scale height H , which has been pro-posed as the characteristic scale of Kelvin-Helmholtz in-stabilities in hot Jupiters (Goodman 2009; Li & Good-man 2010), would yield wind speeds that are an order ofmagnitude too slow to match the observed wind speeds.Instead, a damping length 2 πa , where a is the planet ra-dius, is needed to match the observed wind speeds. Sucha damping length could be either due to a horizontalKelvin-Helmholtz instability or due to the steepening ofday-night standing waves into shocks. We note that theshock-resolving simulations in Fromang et al. (2016) alsofound a dominant scale for horizontal shear instabilitiesof L ∼ πa/
5, and are thus consistent with our results.The upper end of our wind speed estimate is additionallyconsistent with the bulk flow becoming supersonic, andthus prone to dissipation via shocks (Fig. 5). CONCLUSION
We describe the large-scale atmospheric dynamics ofhot Jupiters by modeling them as planetary heat engines.Hot Jupiters are comparable to, but less efficient than,ideal Carnot engines because parcels lose heat to spaceas they move between dayside and nightside. Our theorysuccessfully captures the intensity of winds in a largenumber of hot Jupiter simulations (Fig. 4). Remainingdifferences between theory and simulations are likely dueto our imperfect understanding of numerical dissipationin the simulations, instead of a fundamental shortcomingin our theory.Applying our theory to observations, we find that ei-ther the magnetic dipole field of HD 189733b could bestronger than current estimates suggest, or that its atmo-sphere is dissipating kinetic energy via shear instabilitiesand/or shocks. For HD 209458b our results indicate thatboth drag mechanisms can plausibly match the observa-tions.Looking towards future observations, we expect thatmagnetic drag should become dominant on hotter exo-planets with T eq > APPENDIX
SENSITIVITY TO NUMERICAL PARAMETERS
Our scalings suggest that, for simulations that are dominated by numerical drag, large-scale wind speeds shouldbe sensitive to horizontal resolution (Eqn. 12). To explore this possibility we performed additional simulations inwhich we did not include any Rayleigh drag (including no basal drag) and kept the equilibrium temperature fixed to1500 K while varying different numerical parameters in the model. The two parameters we considered are the model’shorizontal resolution and its timestep dt . Table 1 summarizes the numerical parameter variations for this suite ofsimulations. The Shapiro filter timescale τ num was always kept equal to the timestep.Fig. 6 shows that wind speeds are largely independent of the GCM timestep. We only find a (cid:46)
3% variation in theRMS wind speed while changing dt (and thus also τ num ) over an order of magnitude. Given that Equation 12 predictswind speeds should be independent of dt , this implies a general agreement between our theory and our GCM results.In addition, Figure 6 shows that large-scale wind speeds are less sensitive to horizontal resolution than our scalingwould suggest. Following Equation 12, wind speeds should scale with resolution as U ∝ N − / x , where N x is thenumber of horizontal grid points. Our GCMs do not follow such a scaling and instead we find that the wind speedis independent of resolution to (cid:46)
10% over a factor of 4 change in horizontal resolution, going from C16 to C64.One potential explanation is that our weakly damped simulations develop a direct turbulent cascade of energy tosmaller scales, so that the large-scale kinetic energy of the flow becomes insensitive to the dissipation scale. Anotherexplanation is that hot Jupiter GCM simulations are prone to developing shocks (see Rauscher & Menou 2010; Pernaet al. 2012; Dobbs-Dixon & Agol 2013; Fromang et al. 2016), in which case the large-scale kinetic energy might be lesssensitive to how well the shock is being resolved than Eqn. (12) suggests.Our result is consistent with the suggestion of Heng et al. (2011a) that changes in numerics can change wind speeds inGCMs at the (cid:46)
10% level, but shows that our scaling does not adequately capture the dependence of large-scale GCMwind speeds on numerical resolution. As a result, a better description of numerical drag than our scaling is neededto capture how hot Jupiter GCMs converge with numerical drag. Nevertheless, although our scaling over-predicts thesensitivity to numerical parameters, it does correctly predict the sensitivity to physical parameters, such as equilibriumtemperature (see Fig. 4, left panel).
Physical Parameter Parameter Value(s) Unit
Equilibrium temperature T eq , 2000, 2500, 3000 KVisible absorption coefficient κ SW × − m − kg − Thermal absorption coefficient κ LW . × − × ( p/ . m − kg − Drag timescale τ drag , , , , , ∞ sGravity g − Rotation rate Ω 2 . × − s − Planet Radius a . × mHeat capacity C p . × J kg − K − Specific gas constant R − K − Numerical Parameter Parameter Value(s) Unit
Horizontal resolution ( N x ) C16 (64), C32 (128) , C64 (256) n/aVertical resolution N z
40 n/aTimestep dt sShapiro filter timescale dt num sShapiro filter length scale l num = 2 πa/N x πa/ π a/128 , 2 πa/
256 mShapiro filter order n TABLE 1Range of physical and numerical parameters used in our suite of simulations. Numerical parameters in bold show fiducialvalues used for our main suite of simulations with varying physical parameters, and physical parameters in bold highlightfiducial values used for our secondary suite of simulations with varying numerical parameters. Numbers in parenthesesfor horizontal resolution show the approximate number of horizontal grid points.SENSITIVITY OF MAGNETIC DRAG TIMESCALE TO TEMPERATURE-PRESSURE PROFILE
Because the magnetic drag timescale is highly sensitive to temperature (Perna et al. 2010; Menou 2012; Rauscher& Menou 2013), we explored the impact of the assumed temperature-pressure profile on our results in Section 4.In Section 4 we assume an isothermal atmosphere, here we constrain the vertical temperature structure using theanalytical solutions from Guillot (2010) as follows: we use Eqn. 29 from Guillot (2010) with parameters similar to
100 150 200 250Number of horizontal grid points N x U r m s f r o m G C M [ m s − ] . . . . . . dt = τ num [s] Fig. 6.—
Our scaling for how wind speeds depends on numerical parameters (Eqn. 12) matches the independence of U rms on timestepwell, but does not match the dependence of U rms on grid size. Shown are GCM results for U rms as a function of horizontal resolution (blackdots) and timestep (magenta dots) from simulations with T eq = 1500 K and no Rayleigh drag. In this set of simulations the Shapiro filtertimescale τ num is kept equal to the timestep. Dashed lines show our predicted dependence of U rms on timestep (magenta) and resolution(black), using a value of k such that the theory matches the intermediate GCM point. Eqn. 12 correctly predicts that the wind speed isindependent of timestep (accurate to the 3% level in our GCMs), but predicts that the wind speeds should decrease steeply with increasingresolution, which is not found in our GCM simulations. those used in that paper ( κ LW = 10 − cm g − , γ = 0 . T int = 100K, f = 0 . − bar, and compute wind speeds followingEqn. 11.Figure 7 shows that our conclusions from Section 4 are robust. The most significant difference in Figure 7 comparedto Figure 5 occurs above T eq (cid:38) T eq < H e has an exponential sensitivity to temperature, the absolute value of temperature varies less than a factorof two between 1bar and 10 − bar. This compares to a three order of magnitude change in pressure, which appears inboth density ( ρ ∝ p ) and resistivity ( H e ∝ x − e ∝ p / ) in Eqn. 14. eq [K]10 U [ m / s ] Magnetic drag B Fig. 7.—
Same as the top panel in Figure 5, but instead of an isothermal atmosphere we assume that temperature increases with pressurefollowing the analytic solutions in Guillot (2010). Solid lines are evaluated at 1bar, dashed lines are evaluated at 10 − bar. Compared withFigure 5, our main conclusions are robust to changes in thermal structure.REFERENCESAdcroft, A., Hill, C., Campin, J., Marshall, J., & Heimbach, P.2004, Monthly Weather Review, 132, 2845 Brogi, M., de Kok, R., Albrecht, S., Snellen, I., Birkby, J., &Schwarz, H. 2016, The Astrophysical Journal, 817, 1060 Koll & Komacek