Tidal heating and stellar irradiation of Hot Jupiters

TL;DR

This paper investigates how stellar irradiation and tidal heating interact in Hot Jupiters, revealing a thermomechanical feedback mechanism that can explain their observed radius inflation over billions of years.

Contribution

It introduces an analytical model demonstrating the coupling of vibrational modes with thermal evolution, explaining planetary bloating through resonant tidal heating and stellar irradiation effects.

Findings

Tidal models produce interior radiative zones and enhanced mode dissipation.

Stellar irradiation traps heat at depth, increasing planetary entropy and radius.

Thermally driven winds help maintain spherical atmospheric structure.

Abstract

We study the interaction between stellar irradiation and tidal heating in gaseous planets with short orbital periods. The intentionally simplified atmospheric model we employ makes the problem analytically tractable and permits the derivation of useful scaling relations. We show that many tidal models provide thermal feedback, producing interior radiative zones and leading to enhanced g-mode dissipation with a wide spectrum of resonances. These resonances are dynamically tuned by the thermal feedback, and so represent a novel form of thermomechanical feedback, coupling vibrational modes to the very slow thermal evolution of the planet. We then show that stellar irradiation allows the heat produced by these modes to be trapped at depth with high efficiency, leading to entropy increase in the central convective region, as well as expansion of the planet's radius sufficient to match observed swelling. We find that thermally driven winds play an essential role in this process by making the thermal structure of the atmosphere spherically symmetric within a few scale heights of the photosphere. We characterise the relationship between the swelling factor, the orbital period and the host star and determine the timescale for swelling. We show that these g-modes suffice to produce bloating on the order of the radius of the planet over $\mathrm{Gyr}$ timescales when combined with significant insolation and we provide analytic relations for the relative magnitudes of tidal heating and insolation.

Figures

Equations

$$\psi$$

$$R_0$$

$$r$$

$$P$$

$${\nabla }_{\mathrm{a}}$$

$$\mathcal{R}$$

$$\mathcal{M}$$

$$\kappa_{0}\equiv${10}^{-2}\text{\,}\mathrm{c}\mathrm{m}^{2}\,\mathrm{g}^{-1}$.$$

$$T_{\mathrm{ph}}={\left(\frac{F}{\sigma }\right)}^{1/4},$$

$$F=F_{\mathrm{i}}+F_{\mathrm{e}},$$

$$T_{\mathrm{ph}}={\left(\frac{F_{\mathrm{e}}}{\sigma }\right)}^{1/4}.$$

$${\nabla }_{\mathrm{r}}=\frac{3\kappa PL}{64\pi GM\sigma T^4},$$

$${\nabla }_{\mathrm{a}}={\nabla }_{\mathrm{r}}=\frac{3\kappa P_{\mathrm{t}}{L}_{\mathrm{i}}}{64\pi GM\sigma T_{\mathrm{t}}^4},$$

$$L_{\mathrm{i}}=4\pi R^2{F}_{\mathrm{i}}.$$

$$\kappa ={\kappa }_0{T}^a{P}^b$$

$${\nabla }_{\mathrm{a}}=\frac{3{\kappa }_0{P}_{\mathrm{t}}^{1+b}{L}_{\mathrm{i}}}{64\pi GM\sigma T_{\mathrm{t}}^{4-a}}.$$

$${\nabla }_{\mathrm{r}}={\nabla }_{\mathrm{a}}{\left(\frac{P}{P_{\mathrm{t}}}\right)}^{1+b}{\left(\frac{T}{T_{\mathrm{t}}}\right)}^{a-4}=\frac{d\ln T}{d\ln P}.$$

$$1-{\left(\frac{T_{\mathrm{ph}}}{T_{\mathrm{t}}}\right)}^{4-a}=\frac{4-a}{1+b}{\nabla }_{\mathrm{a}}\left(1-{\left(\frac{P_{\mathrm{ph}}}{P_{\mathrm{t}}}\right)}^{1+b}\right).$$

$$\frac{T_{\mathrm{ph}}}{T_{\mathrm{t}}}={\left(1-\frac{4-a}{1+b}{\nabla }_{\mathrm{a}}\right)}^{\frac{1}{4-a}},$$

$$1-\frac{4-a}{1+b}{\nabla }_{\mathrm{a}}\ge 0.$$

$$F_{\mathrm{i}}=\frac{16\sigma g{\nabla }_{\mathrm{a}}{T}_{\mathrm{t}}^4}{3\kappa P_{\mathrm{t}}}=\left(\frac{32g{\nabla }_{\mathrm{a}}}{3\kappa P_{\mathrm{t}}}\right)F_{\mathrm{e}}.$$

$$P_{\mathrm{t}}=${10}^{13}\text{\,}\mathrm{e}\mathrm{r}\mathrm{g}\,\mathrm{c}\mathrm{m}^{-3}$\left(\frac{R}{R_{0}}\right)\left(\frac{T_{\mathrm{t}}}{$\mathrm{K}$}\right)^{5/2}\left(\frac{\psi}{$8\text{\times}{10}^{-6}$}\right)^{-5/2}.$$

$$\psi$$

$$=\frac{1}{4}\left(\mathcal{R}-2+\sqrt{{\mathcal{R}}^2+4\mathcal{R}-4}\right),$$

$$\Sigma_{\mathrm{t}}\approx\frac{P_{\mathrm{t}}}{g}=$4\text{\times}{10}^{9}\text{\,}\mathrm{g}\,\mathrm{c}\mathrm{m}^{-2}$\mathcal{R}^{3}\mathcal{M}^{-5/3}\left(\frac{T_{\mathrm{t}}}{$\mathrm{K}$}\right)^{5/2}\left(\frac{\psi}{$8\text{\times}{10}^{-4}$}\right)^{-5/2},$$

$$\displaystyle\Sigma_{\mathrm{t}}\approx$2\text{\times}{10}^{6}\text{\,}\mathrm{g}\,\mathrm{c}\mathrm{m}^{-2}$\mathcal{R}^{3}\mathcal{M}^{-5/3}\left(\frac{F_{\mathrm{e}}}{F_{\sun}}\right)^{5/8}\psi^{-5/2}.$$

$$\frac{F_{\mathrm{i}}}{F_{\mathrm{e}}}$$

$$4\pi R^2{F}_{\mathrm{e}}=\frac{\pi R^2{L}_{\star }}{4\pi a_{\mathrm{orbit}}^2},$$

$$F_{\mathrm{e}}=\frac{L}{16\pi a_{\mathrm{orbit}}^2}.$$

$$\frac{F_{\mathrm{e}}}{F_{\text{sun}}}=\frac{1}{4}\left(\frac{L_{\star }}{L_{\text{sun}}}\right){\left(\frac{a_{\mathrm{orbit}}}{R_{\text{sun}}}\right)}^{-2}.$$

$$\frac{F_{\mathrm{i}}}{F_{\mathrm{e}}}$$

$$\Sigma_{\mathrm{t}}\approx\frac{P_{\mathrm{t}}}{g}=$8\text{\times}{10}^{5}\text{\,}\mathrm{g}\,\mathrm{c}\mathrm{m}^{-2}$\mathcal{R}^{3}\mathcal{M}^{-5/3}\left(\frac{L_{\star}}{L_{\sun}}\right)^{5/8}\left(\frac{a_{\mathrm{orbit}}}{R_{\sun}}\right)^{-5/4}\psi^{-5/2}.$$

$${\mathbf{F}}_{\mathrm{v}}={\nu }_{\mathrm{v}}\frac{\partial \text{varv}}{\partial r},$$

Full text

††pubyear: 2017††pagerange: Tidal heating and stellar irradiation of Hot Jupiters–References

1 Introduction

In recent data from Kepler and ground-based followups, there is evidence for a large population of hot Jupiters which are substantially inflated relative to their degenerate radii (Hellier et al., 2012; Weiss et al., 2013; Hartman et al., 2012). The radii and periods of known members of this population as well as of the broader Jupiter-sized population are shown in Fig. 1 (Rein, 2012). There is an apparent split in the observed population around periods of $10\text{\,}\mathrm{d}$, such that planets with longer periods are generally not inflated while those with shorter periods are often substantially inflated. Importantly, planets at or above $2R_{\mathrm{J}}$ must be inflated relative to their degenerate radii, otherwise their implied masses would make them stars (Stevenson, 1991). Importantly, Jupiter has approximately the maximum radius for an unheated gas giant, so planets at $2R_{\mathrm{J}}$ must be bloated regardless of their mass. In order to achieve this level of expansion, the central convection zone must be heated considerably relative to what would be expected as a result of the residual heat of formation (Lopez & Fortney, 2016), and there is evidence of planets re-inflating after cooling down (Hartman et al., 2016). Complicating this is the thermodynamic requirement that heat flows only from hot to cold, not in the reverse fashion, and absent an internal heat source there is no means for blocked heat transport to heat the interior. This, combined with the expectation that temperature increases towards the core of the planet, means that any change in temperature at depth must be due to heat generated at or deeper than the point of interest.

A variety of mechanisms have been suggested to generate deep heating, with ohmic processes (Batygin, Stevenson & Bodenheimer, 2011; Spiegel & Burrows, 2013) and tidal dissipation (Socrates, 2013; Miller et al., 2009) among the more popular models. Given that heat cannot flow from the surface of the planet into its core, the stellar flux is often neglected. However, somewhat surprisingly, the observed radii correlate strongly with the incident stellar flux, so that this flux may play a role in the inflation process (Lopez & Fortney, 2016). Confounding this analysis is the fact that stellar flux is not independent of orbital period. So a theory of hot Jupiter inflation must separately handle the effects of orbital period and incident flux, particularly when dealing with tidal heating.

We investigate the effects of stellar flux on the structure of an internally heated hot Jupiter, making few assumptions about the nature or profile of the heating and considering the effects of wind redistribution. We show that the stellar flux acts to modulate the rate at which heat escapes from the planet. We then investigate the feedback that this heating produces on the thermal structure of the planet and show that a wide variety of realistic heating profiles gives rise to interior radiative zones. These zones migrate within the planet on thermal timescales, giving a broad and dynamically tuned spectrum of g-mode resonances which dissipate heat tidally in the planet. We then show that these modes suffice to produce the observed bloating. Finally we predict the relation between stellar flux, orbital period and planetary radius.

The new thermomechanical feedback mechanism we propose, shown schematically in Fig. 2, underscores the importance of considering planets as dynamical objects with complex behaviours coupling wildly different timescales. Vibrational effects with periods ranging from seconds to days can have a tremendous impact on thermal evolution over millions of years, and that thermal evolution in turn feeds back into the vibrational modes, creating a dynamically tuned spectrum which can ultimately determine the large-scale structure of the planet.

2 Isotropic Planetary Structure

We discuss the structure of a planet isotropically illuminated by flux $F_{\mathrm{e}}$ from its host star. Fig. 3 shows the orbital configuration of the planet-star system and Fig. 4 shows the thermodynamic structure of the planet with the relevant variables defined schematically.

For the deep interior of the planet we adopt the analytic brown dwarf structure of Stevenson (1991):

[TABLE]

where $\psi$ is the electron degeneracy parameter, $\rho$ is the density, $T$ is the temperature, $R_{0}$ is the degenerate radius, $r$ is the radial coordinate, $P$ is the pressure, $k_{\mathrm{B}}$ is the Boltzmann constant, $E_{\mathrm{F}}$ is the Fermi energy, $M$ is the mass of the planet and $\nabla_{\mathrm{a}}$ is the adiabatic temperature gradient. These relations effectively parametrise a $\gamma=5/3$ adiabatic atmosphere, accounting for electron degeneracy at high pressures. We assume solar composition in this paper, so that the mean molecular weight of electrons $\mu_{\mathrm{e}}\approx 1.15$. In addition, we take $R$ to be the radius of the planet. For convenience, we define the parameters

[TABLE]

where $M_{\mathrm{J}}=1.838\times 10^{30}\mathrm{g}$ is the mass of Jupiter. We expect the gas line opacity to dominate in hot atmospheres, so we use this as fiducial and define (Stevenson, 1991)

[TABLE]

We connect the top of the convection zone to the photosphere with a radiative zone at transition pressure $P_{\mathrm{t}}$. The photospheric temperature is given by

[TABLE]

where $F$ is the total flux leaving the planet’s atmosphere and $\sigma$ is the Stefan-Boltzmann constant. This flux may be divided into two components, as

[TABLE]

where $F_{\mathrm{i}}$ is the heat arriving from the planet’s interior and $F_{\mathrm{e}}$ is the heat arriving from the host star. When $F_{\mathrm{e}}$ is large relative to the flux which would escape through the planet’s natural cooling the photospheric temperature is determined entirely by $F_{\mathrm{e}}$, such that

[TABLE]

The escaping flux $F_{\mathrm{i}}$ from the planet’s core is just the flux which escapes from the convection zone. This may be generated by gravitational contraction, radioactive decay or by a decrease in the interior entropy; here we generally assume the last of these to be the dominant source of interior flux.

The radiative gradient is given by

[TABLE]

where $L$ is the luminosity driving the gradient. At the convective-radiative boundary we have

[TABLE]

where

[TABLE]

If we take $\kappa$ to be a power law in both $P$ and $T$ of the form

[TABLE]

then

[TABLE]

Above the transition we have

[TABLE]

Integrating from the photosphere to the transition yields

[TABLE]

The photosphere pressure is generally much lower than the transition pressure so

[TABLE]

which has a solution if and only if

[TABLE]

If $b>-1$ and $a>4$ or $b<-1$ and $a<4$ or $b>\frac{3-2a}{5}$ and $a<4$ this is satisfied. There are other conditions under which it is satisfied, but these are the most relevant common cases. The exponents are generally of order unity and the result is raised to a small power so $T_{\mathrm{ph}}\approx T_{\mathrm{t}}$. This agrees with other analyses, which have found that in Jupiter-like planets, $T_{\mathrm{ph}}<T_{\mathrm{t}}<1.5T_{\mathrm{ph}}$ (Ginzburg & Sari, 2015). The precise temperature ratio depends on the nature of the opacity function, so we simply take $T_{\mathrm{t}}=2^{1/4}T_{\mathrm{ph}}$ (the Eddington closed grey body) as representative. If this holds and the flux from the star is dominant then

[TABLE]

Eliminating $\mu_{e}$ between equations 2.1 and  2.4, we may write $P_{\mathrm{t}}$ in terms of $T_{\mathrm{t}}$ at $r\approx R$ and $R_{0}$, such that

[TABLE]

From equation (2.3) and the definition of $R_{0}$ (equation (2.6)) we find

[TABLE]

where we have taken the positive root because $\psi>0$ and $R>R_{0}$. Note that this is always of order unity, so to good approximation the majority of the variation in $P_{\mathrm{t}}$ comes from the $R$ and $T_{\mathrm{t}}$ dependence in equation (2.22).

If the convective-radiative transition occurs at a shallow point in the atmosphere the corresponding column density is just

[TABLE]

again with $\mathcal{M}=M/M_{\mathrm{J}}$. Eliminating $T_{\mathrm{t}}$ in favour of $F_{\mathrm{e}}$ and using $T_{\sun}=5777\mathrm{K}$ we find

[TABLE]

Inserting this result into equation (2.21) we obtain

[TABLE]

As one final manipulation, we wish to put our equations in terms of the stellar luminosity and orbital radius. The stellar luminosity is related to the external flux $F_{\mathrm{e}}$ by

[TABLE]

where $L$ is the stellar luminosity, $a_{\mathrm{orbit}}$ is the orbital radius of the planet. The factor of $\pi R^{2}$ on the right-hand side is just the cross-section of the planet as seen from the star, while the factor of $4\pi R^{2}$ on the left-hand side reflects the definition of $F_{\mathrm{e}}$ as an average over the surface of the planet. So

[TABLE]

Comparing with the Sun we find

[TABLE]

Thus

[TABLE]

Importantly, the exponent on the luminosity is greater than $-1$. This means that while the ratio of escaping to incident flux decreases with increasing stellar flux, the total escaping flux increases. This conclusion is dependent primarily on how strongly the ratio $T_{\mathrm{t}}/T_{\mathrm{ph}}$ varies with $F_{\mathrm{e}}$, which in turn depends on the form of the opacity function. In particular, it does not generally hold at extremely high temperatures where the gas line opacity ceases to dominate and Kramers-like rules take over. For brown dwarfs and hot Jupiters, however, this variation is small and should not pose a problem. It is also useful to compute the transition column density

[TABLE]

This is small enough that the shallow approximation is not bad.

3 Angular Temperature Distribution

The planets under consideration are generally highly insolated. This can lead to significant temperature differences between the day and night sides, particularly if the planet is tidally locked. In this section we show that winds suffice to make the thermal structure of the atmosphere spherically symmetric at depth even when there is a large temperature difference at the photosphere. This allows us to treat the structure of the planet as spherically symmetric where tidal effects are most prominent.

Consider a wind driven from one side of the planet to the other along isobars with characteristic velocity $\boldsymbol{\varv}$. Suppose further that the character of this wind changes in the vertical direction over distances of order the pressure scale height $h$ and that it changes in the horizontal direction over distances of order the planet’s radius. The specific force due to shear in the vertical direction is

[TABLE]

where $\nu_{\mathrm{v}}$ is the viscosity for a circumferential flow shearing in the vertical direction. The corresponding power dissipated is

[TABLE]

Likewise, the force due to shear in the horizontal direction is

[TABLE]

where $\nu_{\mathrm{h}}$ is the viscosity for a circumferential flow shearing in the other circumferential direction and $\xi$ is a coordinate along the flow. The corresponding power dissipated is

[TABLE]

The total power dissipated is then

[TABLE]

where we have approximated the velocity derivatives with the velocity magnitude and the relevant scale heights, the pressure scale height $h$ in the vertical direction and the radius $r$ in the horizontal. We have also simplified the viscosity from a rank-4 tensor to two scalars, so this relation ought only to be interpreted as an order of magnitude of the power.

To determine $v$, we now match this power to the work which the wind may extract as a heat engine. We are interested in cases where the temperature difference between the two sides is large so the efficiency of the heat engine is of order unity even if diffusive losses make it irreversible. We may neglect diffusive losses because we have taken the microscopic thermal diffusivity to be small on the relevant scales. So we may write the specific rate of work as

[TABLE]

where

[TABLE]

for a monatomic ideal gas, is the specific heat at constant pressure. This is simply the specific heat which is transported from one side of the planet to the other. In our case

[TABLE]

so $T$ refers to the average temperature while $\Delta T$ refers to the temperature difference. By definition, $\Delta T/T\leq 2$. In the most extreme case this gives

[TABLE]

where

[TABLE]

is the adiabatic sound speed and $\mu$ is the mean molecular weight. Equating the rate of work and power gives

[TABLE]

To proceed further we must examine the forms of $\nu_{\mathrm{v}}$ and $\nu_{\mathrm{h}}$. The nature of the viscosity differs between stably stratified and buoyantly unstable zones, so we must determine which of these are relevant and treat them separately.

We begin with radiative zones. In a stably stratified region the two viscosities differ because of Richardson stabilisation, an effect which limits the scale of turbulence in the vertical direction by means of a buoyant restoring force (Galperin, Sukoriansky & Anderson, 2007). A straightforward prescription for the viscosities in this context is

[TABLE]

where $\alpha$ is the microscopic thermal diffusivity (Mathis et al., 2004). Generally we expect $\alpha$ to be small compared to $\nu_{\mathrm{h}}$ because horizontal radiative transfer is inefficient., so we may neglect $\alpha$ and write

[TABLE]

By the Schwarzschild criterion $\nabla<\nabla_{\mathrm{a}}$ in a stably stratified zone. In general we expect radiative transport to be efficient far from the zone boundaries, so we take $\nabla\ll\nabla_{\mathrm{a}}$ in most of such a zone. Using this we write

[TABLE]

Now making use of

[TABLE]

we find

[TABLE]

Inserting equations (3.13) and  (3.18) into equation (3.12) gives

[TABLE]

This may be rearranged to

[TABLE]

Solving gives

[TABLE]

The positive branch is the one of interest, because we have implicitly taken $\varv>0$ in writing it as a magnitude. In the upper regions of the planet’s atmosphere $r\gg h$ so

[TABLE]

Using equation (3.17) the rate at which heat is transported may be written as

[TABLE]

The region of interest is shallow so $gr^{2}\approx GM$ and

[TABLE]

The depth, as measured by column density $\Sigma_{\mathrm{i}}$, over which the winds make the flux distribution spherically symmetric is

[TABLE]

Evaluating the sound speed at the photosphere gives

[TABLE]

where $m_{\mathrm{p}}$ is the proton mass. For comparison, the photosphere is at a depth of

[TABLE]

Thus the temperature distribution becomes spherically symmetric deeper than the photosphere but shallower than the convective transition. So we need not worry about the viscosity in convection zones.

This remains valid as long as the planet rotates slowly relative to $v/R$, such that the characteristic scale of circumferential motion remains $R$ and is not reduced by Coriolis effects. At short periods, where this condition is most in danger, the anisotropy is very large, such that $v\approx hc_{\mathrm{s}}/r$, and the surface temperature should be quite high because of insolation, such that $c_{\mathrm{s}}\approx$3\text{\times}{10}^{5}\text{,}\mathrm{c}\mathrm{m},\mathrm{s}^{-1}. In this regime, the period of a Jupiter-radius planet must be at least $30\text{\,}\mathrm{d}$ with $h/r\approx${10}^{-2} or $3\text{\,}\mathrm{d}$ with $l/r\approx${10}^{-1}$$ for the Coriolis effect to be negligible. Even at the shortest known periods of just under a day, the correction term is not too great and does not alter the conclusion that the temperature distribution becomes spherical above the convection zone, so we continue to use this approximation with the knowledge that it becomes worse as the period diminishes.

4 Heated Thermal Structure

In this section we work on timescales long compared to the adjustment of radiative or convective zones to thermal perturbations but short compared to the characteristic thermal timescale of the planet. This is the instantaneous equilibrium approximation. This separation of scales exists because the thermal timescale of the planet is set by the thermal content of the core, whereas the radiative and convection regions of interest are shallow zones with much less mass and at much lower temperatures.

The equations governing the luminosity of the planet as a function of mass coordinate are

[TABLE]

where $\varepsilon(m)$ is the specific energy generation by tides, radioactive decay and ohmic processes and the mass coordinate $m$ corresponds to the spherical shell containing mass $m$. Note that mechanical expansion and contraction can generate energy, but in this coordinate system that generation provides no net contribution because it does not alter the specific entropy.

Thermodynamic consistency imposes the condition that heat travels from hot regions to cool ones. Assuming that $T$ increases towards the core of the planet, this means that $L(m)\geq 0$ everywhere. In a convective atmosphere, the thermal gradient is almost independent of the luminosity. This follows because the luminosity is determined by the superadiabaticity of the thermal gradient, rather than by the gradient. When convection is efficient, the convective zone is nearly isentropic so $L\propto(\nabla-\nabla_{\mathrm{a}})^{3/2}$ (Kippenhahn & Weigert, 1990) and the atmosphere achieves significant scaling of luminosity with only small changes to $\nabla$. As a result, the conditions on $L$ cannot generally be satisfied. This means that radiative zones are generically needed as interfaces between convective regions. More formally, we work in the limit of perfectly efficient convection, such that

[TABLE]

We also make the assumption that the convective turnover time for any region of interest is much shorter than the time-scale over which thermal quantities change, such that convection may be assumed to enforce an instantaneous adiabatic law.

Now suppose that we perturb a planet by injecting luminosity $\Delta L$ somewhere below the radiative-convective boundary. For $\Delta L\ll L_{\mathrm{i}}$, we may solve equation (4.1) by simply reducing the luminosity escaping from deeper regions of the planet. That is, $L_{\mathrm{i}}$ goes unchanged but the luminosity in regions deeper than the injection depth is reduced by $\Delta L$. In the limit of very efficient convection (or large opacity), this adjustment holds until $\Delta L\approx L_{\mathrm{i}}$. For $\Delta L>L_{\mathrm{i}}$ the adjustment still occurs, with the deep luminosity falling to the radiative luminosity at the adiabatic gradient, the minimum needed to maintain convection. The difference is that in this case there is an excess of luminosity reaching the convective-radiative transition and this must be accounted for. At the boundary we must have

[TABLE]

which must remain satisfied when we perturb $L_{\mathrm{i}}$ so

[TABLE]

If the transition temperature is similar to the photospheric temperature and if the radius does not change substantially owing to the perturbation $\Delta\ln L=\Delta\ln(L_{\mathrm{i}}+L_{\mathrm{e}})=4\Delta\ln T_{\mathrm{t}}$. Because $F_{\mathrm{e}}\gg F_{\mathrm{i}}$ and $L_{\mathrm{e}}$ is fixed

[TABLE]

So we may neglect the change in $T_{\mathrm{t}}$ and find

[TABLE]

We generally expect that, at fixed temperature, $\kappa$ rises as $P$ rises. As a result, $P_{\mathrm{t}}$ must fall to satisfy this relation, so either the entropy of the central adiabat must rise or the adiabatic law must be broken somewhere in the planet. The central entropy cannot rise unless either heat is being added at the core or the photosphere is hotter than the core, because heat cannot be forced to move up the temperature gradient. Neither of these are generally the case so the adiabatic law must be broken. As a result the planet must form an interior radiative zone.

To characterise these radiative zones, let $P_{1}$ be the transition pressure between the surface radiative zone and the new convection zone, $P_{2}$ the transition pressure between this zone and the interior radiative zone and $P_{3}$ the transition pressure between this zone and the central adiabat. The perturbed and unperturbed pressure-temperature structures are shown in Fig. 5. Let $T_{j}$, $m_{j}$, $\kappa_{j}$ and $L_{j}$ be the corresponding temperature, mass coordinate, opacity and luminosity at each transition. The new convection zone is adiabatic, so

[TABLE]

Assuming that $m\approx M$, the condition (4.5) for transition between radiative and convective zones gives

[TABLE]

Finally, recalling that the central adiabat is fixed and that the new adiabat contains a point at the fixed temperature $T_{\mathrm{t}}$, we find

[TABLE]

where the subscripts $\mathrm{i}$ and $\mathrm{f}$ refer to the initial unperturbed and final perturbed system respectively. Equation (4.11) thus expresses the entropy difference between the two adiabats. Note that with these subscript definitions,

[TABLE]

Again with $\kappa\propto T^{a}P^{b}$ we may write equations (4.6) – (4.12) as a system of linear equations in the logarithms of temperature and pressure. Solving this system yields

[TABLE]

where

[TABLE]

The transition temperature $T_{1}\equiv T_{\mathrm{t}}$ is known from the unperturbed state so with equations (4.13) and (4.14) we may determine the remaining temperatures. Likewise the unperturbed transition pressure $P_{\mathrm{t,i}}$ is known from the unperturbed state. The perturbed transition pressure $P_{1}$ is related to the unperturbed by equation (4.8) so, with equations (4.15) and (4.16) we may determine the remaining pressures.

In equilibrium, the luminosities are related by

[TABLE]

With these we can compute the luminosity ratios. A consequence of equation (4.15) is that the new convective zone is maintained by heat generation in between $P_{1}$ and $P_{2}$, or equivalently between $m_{1}$ and $m_{2}$, because this is what allows for $L_{1}\neq L_{2}$.

The minimum luminosity required for convection may be calculated from equation (4.5) as

[TABLE]

Both sides of this equation are functions of pressure. We generally expect that $\nabla_{\mathrm{r}}$ rises quickly towards the interior of the planet as convection becomes more efficient so $L_{\mathrm{min}}$ is a small fraction of $L_{\mathrm{i}}$. This is actually guaranteed by equation (2.20) so we expect that $L_{\mathrm{min}}$ is suppressed relative to $L_{\mathrm{i}}$ by a power-law in $P$ and may calculate

[TABLE]

where $P_{\mathrm{inject}}$ is the pressure inside which minimal luminosity is injected.

5 Expansion

The expansion associated with changing the temperature profile of the planet is given by

[TABLE]

In the limit where $\Delta R/R$ is small, $\Delta P/P$ is small at fixed $m$, so

[TABLE]

Substituting this into equation (5.1) we find

[TABLE]

where the coordinate $r$ refers to the unheated system. The integration proceeds up to $R$ as an approximation, once more in the limit where $\Delta R/R$ is small. When this is the case and when the majority of the heating occurs near the surface at $r\approx R$ this may be approximated by

[TABLE]

Now we may approximate $\Delta\ln T$ as $\ln(T_{2}/T_{1})$. The pressure depth over which this approximation (rather than $\Delta\ln T\approx 0$) is valid is $\Delta\ln P\approx\ln(P_{1}/P_{\mathrm{t,i}})$. This corresponds to a physical depth of $h\ln(P_{1}/P_{\mathrm{t,i}})$, because $h$ is the characteristic scale of the thermal properties of the planet and hence sets the scale of the radiative zone which forms. So we may approximate equation (5.4) in terms of the heating parameters by

[TABLE]

With equation (4.8) we find

[TABLE]

For very deep zones, the relevant scale height is that near the base of the zone rather than the top, because the majority of the contribution to the integral comes from this region. This may be taken into account by noting that the scale height at the base of the radiative zone is given by

[TABLE]

Inserting equation (4.13) and equation (4.14) we have

[TABLE]

Making use of $T_{1}\approx T_{\mathrm{ph}}$, equations (2.11) and (2.28) give us

[TABLE]

where $\tau$ is the orbital period. The expansion is therefore

[TABLE]

A factor of a few from the luminosity term is therefore sufficient to substantially inflate the planet at short orbital periods.

6 G-Modes

The existence of internal radiative zones raises the possibility that g-modes may contribute to tidal heating. This is particularly interesting because, if g-mode dissipation is the dominant form of tidal heating $\epsilon$ is actually a function of the thermal structure of the planet. That is because g-modes predominantly resonate in radiative zones. What this amounts to is a form of feedback between the thermal and mechanical structures of the planet.

6.1 Dynamical Tide

In principle there are two sources of dynamical tides, namely gravitational and thermal. We expect that thermal tides do not couple to the g-modes considered here. There are two reasons for this. First, the thermal tide is significant only in the upper layers of the atmosphere where insolation is significant. In particular, the tide damps as $e^{-\kappa\Sigma}$ (Arras & Socrates, 2010). The internal radiative zone begins at a comparable column density to the unperturbed radiative-convective transition. Equation (2.32) gives $\kappa_{0}\Sigma_{\mathrm{t}}\approx 5\times 10^{2}$, so the damping is on the order of $\exp(-5\times 10^{2})$, which suffices to make this effect negligible. Secondly, the thermal tide relies on timescale for redistributing heat being large relative to the orbital time. We have shown that the temperature distribution becomes spherical very near the photosphere and well above the convection zone, even for a tidally locked planet. This means that it will not reach even the upper convection zone. As a result we restrict our analysis to gravitational tides.

Due to their frequencies being small relative to the acoustic frequency, g-modes are unlikely to substantially compress material in the planet. As a result we must treat them in the incompressible limit. This may be done by separating the perturbing tidal potential into a hydrostatic equilibrium tide and a dynamical tide (Zahn, 1975). The associated radial displacements $\xi^{\mathrm{eq}}$ and $\xi^{\mathrm{dyn}}$ obey the relations (Goodman & Dickson, 1998)

[TABLE]

and

[TABLE]

where $l$ is the latitudinal quantum number, $\omega$ is the frequency, $\Phi$ is the unperturbed planetary gravitational potential, $\delta\Phi$ is the perturbing tidal potential due to the star and $N$ is the Brunt-Väisälä frequency, with $N^{2}$ positive in the radiative zone and negative in the surrounding convective regions. To analyse this equation we first solve the homogeneous version, with the right hand side set to zero, and then compute the overlap between the resulting modes and the forcing term given by the right-hand side.

6.2 Mode Profile

The homogeneous part of equation (6.2) is

[TABLE]

Defining

[TABLE]

we find

[TABLE]

We now wish to perform a change of variables which will eliminate the first order derivative of $\xi$. To do this, we note that

[TABLE]

This may be written as

[TABLE]

Using this, we pick

[TABLE]

which gives

[TABLE]

With this substitution, equation (6.5) becomes

[TABLE]

Qualitatively we expect $N$ to peak near the centre of the radiative zone and fall to zero at the edges. To fit this, we pick a quadratic form in our new coordinate $y$, such that

[TABLE]

where $r_{0}$ is the radial coordinate of the centre of the radiative zone and $2\delta\!y$ is the width of the zone in $y$. Defining

[TABLE]

the differential equation equation (6.10) becomes

[TABLE]

We now define

[TABLE]

such that

[TABLE]

Note that $q$ is positive and large because

[TABLE]

This follows because $\rho$ has characteristic scale $h$ and because $h\ll r$ except near the core of the planet.

It is now worth noting that the physical width of the zone in $r$ is

[TABLE]

This holds because for a thin zone, $\rho$ does not change too much across it. In thick zones there would be deviations from this which we neglect. For convenience we now define

[TABLE]

With this, equation (6.16) becomes

[TABLE]

This may also be written as

[TABLE]

which is the same as the equation for a quantum harmonic oscillator with energy $\Omega^{2}-q$, mass $\hbar^{2}\beta/2$, and zero-point energy $\Omega/\sqrt{\beta}$. The eigenvalues are therefore quantised in the form

[TABLE]

The sign of $\Omega$ does not enter into equation (6.10) so we may take whichever branch of the solutions to this equation that we choose. Taking the positive we see that

[TABLE]

These correspond to periods and frequencies of

[TABLE]

The radial profiles of the solutions are the product of an exponential with an Hermite polynomial. For $l=2$, the dominant tidal mode, these are

[TABLE]

where $Y_{lm}$ are the spherical harmonics, $n\in\left\{0,1,2,...\right\}$ and $m\in\left\{-2,1,0,1,2\right\}$. The modes are normalised so that

[TABLE]

and we take $\rho$ and $r$ as constants throughout the radiative zone to compute this normalisation. This is consistent with approximations we make elsewhere.

6.3 Overlap Integral

In order to compute the tidal forcing $F_{n,m}(\omega)$, we must say something about the origin of the tidal potential. There are two potential sources, rotational asynchronisation and orbital eccentricity. In the former, the tidal forcing occurs at a frequency $\omega_{\mathrm{rotation}}-\Omega_{\mathrm{orbit}}$, while in the latter it occurs at a frequency of $\Omega_{\mathrm{orbit}}$. In both cases, working in the frame corotating with the planet’s orbit,

[TABLE]

and further involves a sum of $l=2$ spherical harmonics. Beyond this the two cases differ significantly because the eccentricity case has $\delta\Phi\propto e$ while the asynchronous case has no such factor. To capture both cases, we write

[TABLE]

where $\phi_{m^{\prime}}$ are phase factors, the factors $k_{m^{\prime}}$ capture the magnitudes of the various harmonics and sum in quadrature to unity and $\Pi$ is a dimensionless factor of order unity in the asynchronous case and of order $e$ in the eccentric case. From this form and equation (6.1) we may write the equilibrium tide as

[TABLE]

The driving term associated with this equilibrium tide is the right-hand side of equation (6.2), given by

[TABLE]

In computing the overlap of this with the eigenmodes of the homogeneous equation, we may treat factors of $r$ as constant, because the radiative zone ought to be thin on the scale of the planetary radius. As a result, the projection is

[TABLE]

where we have centred the integral on $r_{0}$, the radial coordinate corresponding to $y_{0}$, and defined

[TABLE]

We have also extended the integration bounds to infinity to make the computation easier because the exponential suppression in $x$ makes the precise bounds irrelevant. Note that changing variables from $r^{\prime}$ to $x$ is formally quite complicated, though in the approximation where $\rho$ changes little over the course of the zone it just produces a prefactor of $l_{r}$.

The integral may now be evaluated with the generating function of the Hermite polynomials,

[TABLE]

so

[TABLE]

As a result,

[TABLE]

for even $n$ and vanishes for odd $n$. There is complete degeneracy in both the dissipation and oscillation over $m^{\prime}$, so we may form a linear combination of spherical harmonics which precisely matches the forcing term. This amounts to summing the right hand side of equation (6.42) in quadrature over $m^{\prime}$ and taking the square root, which gives

[TABLE]

This expression gives the amplitude of the resonance. From this stage we take it as given that $l=2$ and drop the label on $\psi$.

6.4 Dissipation

The square of the displacement, which is proportional to the dissipation, has maxima at a distance of order $\pm a$ from the centre of the radiative zone so, even if the dampening were uniform, we would expect the dissipation to be greatest near the edges of the zone. In practice, convective turbulence increases the dissipation just outside the zone and this assertion is even stronger. To evaluate the strength of this effect we turn to various linear dissipation mechanisms. Both radiative and viscous damping are potentially relevant. For each of these we may calculate a quality factor $Q$, giving the number of undriven cycles required for an $e$-fold reduction in strength. These combine as

[TABLE]

We begin with radiative damping at finite opacity. The quality factor of mode $n$ is of order

[TABLE]

where $c$ is the speed of light and $\tau_{n}$, $\omega_{n}$ and $\lambda_{n}$ are the lifetime, frequency and wavelength corresponding to mode $n$ (Press, 1986). The wavelength is given by

[TABLE]

which just comes from the fact that mode $n$ has $n+1$ nodes over the zone width of $2l_{r}$. Thus

[TABLE]

Let

[TABLE]

an approximate zone mass. We find

[TABLE]

From equation (5.7) we know that the scale height $h$ is proportional to $T$, so

[TABLE]

We are interested in regions which are sufficiently shallow so that $g$ is nearly constant and so

[TABLE]

and

[TABLE]

so that

[TABLE]

The radiative zone is stably stratified so we expect turbulent damping to be limited to the evanescent part of the mode which leaks into the neighbouring convective zones. The Navier-Stokes equation with a simple viscosity term is

[TABLE]

We now neglect the non-linear term because we are interested in understanding the linear growth and decay of modes and expect the absolute velocities involved to be small. Without the nonlinear term

[TABLE]

The balance between gravity and the pressure gradient is what gives us g-modes, so we may write the balance as

[TABLE]

Now let

[TABLE]

so that

[TABLE]

The kinetic energy density

[TABLE]

where $\boldsymbol{\varv}^{*}$ is the complex conjugate of $\boldsymbol{\varv}$, and evolves as

[TABLE]

where $\nabla^{\dagger}$ is the adjunct gradient operator. When the spatial derivatives are greatest in the radial direction the Laplacian just produces a factor of $(2\pi/\lambda_{n})^{2}$ so

[TABLE]

Integrating this equation over the whole planet we find

[TABLE]

The viscosity $\nu$ is only significant in the convection zones on either side of the radiative zone, where turbulent viscosity dominates. The kinetic energy density $K$ is only significant inside the radiative zone and within a few wavelengths of the zone edges on either side. The damping integral is dominated by the region in which neither is small, so

[TABLE]

where $\nu$ is evaluated in the convecting regions and $2(l_{r}+\lambda_{n})$ is roughly the radial extent of the region where the kinetic energy density is significant. Similarly

[TABLE]

where $K$ on the right hand side is the average kinetic energy density in the zone. Then

[TABLE]

Because $m_{\mathrm{z}}$ is constant we find that

[TABLE]

and the damping timescale is

[TABLE]

with related quality factor

[TABLE]

The turbulent diffusivity is of order $v_{c}h$ when the convective turnover is on a timescale shorter than the forcing frequency $\omega$. It is $\omega$ rather than $\omega_{n}$ that matters here because the oscillation physically takes place at the driving frequency, not the mode period. The relevant turbulent frequency for motion over length-scale $l_{\mathrm{t}}$ is

[TABLE]

Taking $v_{c}$ to be given by a Kolmogorov spectrum, we find

[TABLE]

It follows that the relevant diffusive motions are on a scale

[TABLE]

and corresponding diffusivity is

[TABLE]

This is just the result of Goldreich & Keeley (1977) and yields a quality factor

[TABLE]

High-frequency driving leads to a high quality factor and $Q_{n}$ scales as $\omega_{n}^{2}$ because $\lambda_{n}\propto\omega_{n}$. This is a weaker scaling than the radiative $Q$, which goes as $\omega_{n}^{3}\propto(n+1)^{-3}$. Thus at large $n$ the convective mechanism dominates.

We are often interested in the lowest $n$ because this mode is the least suppressed by overlap factors. The convective flux of interest is generally $F_{\mathrm{i}}$, which equation (2.31) shows is on the order of $10^{-5}F_{\mathrm{e}}$. The external flux is typically about $10^{-2}F_{\sun}$, so the relevant convective flux is on the order of $5\text{\times}{10}^{3}\text{\,}\mathrm{e}\mathrm{r}\mathrm{g}\,\mathrm{c}\mathrm{m}^{-2}\,\mathrm{s}^{-1}$. For a density of ${10}^{-1}\text{\,}\mathrm{g}\,\mathrm{c}\mathrm{m}^{-3}$ the convection speed is then

[TABLE]

So for a scale height of ${10}^{9}\text{\,}\mathrm{c}\mathrm{m}$, $v_{\mathrm{c}}/h\approx$3\text{\times}{10}^{-8}\text{,}\mathrm{H}\mathrm{z}$$. We show later that we are interested in frequencies on the order of ${10}^{-6}\text{\,}\mathrm{H}\mathrm{z}$. This means that the factor $[h\omega/v_{\mathrm{c}}(h)]^{2}$ accounting for the eddy time is of order $10^{5}$, so the convective quality factor is $Q_{1}\approx 300$. By comparison, the fiducial radiative quality factor with the same assumptions is $Q_{1}\approx 1$. Thus we expect radiative damping to dominate by a reasonable margin unless the fluxes involved are many orders of magnitude larger or the frequencies are several orders of magnitude smaller.

6.5 Boundaries

We have shown that, when g-modes dominate the dissipation, $\varepsilon$ is significant primarily near the edges of the radiative zone. It is also straightforward to show that the dissipation is symmetric because the squared mode profiles are even. So $\varepsilon$ is an even function, the integral of which is dominated by the regions just outside the zone boundaries. Suppose that the total luminosity produced by tides is $L_{\mathrm{t}}$ and that a fraction $f$ of this is produced inside the radiative zone. In steady-state

[TABLE]

and

[TABLE]

Using equation (4.21) we find

[TABLE]

If the heating is large relative to $L_{\mathrm{i}}$ and $f$ is small then we expect

[TABLE]

so

[TABLE]

Using equations (4.13) through (4.16) we find

[TABLE]

With equations (6.81) and (6.82) equation (5.11) yields

[TABLE]

Recall from equation (6.25) that the resonant frequencies are

[TABLE]

The right hand side has two characteristic regimes, one in which $l_{r}/h$ is large and one in which it is small or of order unity. In the former case, $\omega_{n}\propto h$, while in the latter $\omega_{n}\propto l_{r}$. Both $h$ and $l_{r}$ increase with $P_{3}/P_{2}$, so the resonant frequency goes up as the zone width increases. This means that the resonance shifts up in frequency as the luminosity increases. That is,

[TABLE]

So an increase in luminosity tends to tune the system towards resonance if it is being driven above resonance and pushes it away from resonance otherwise. The net result is that there is thermomechanical feedback which tends to bias systems towards resonance, particularly when their resonant frequency is below the driving frequency.

6.6 Power Production

Each mode may be treated as a separate damped and forced harmonic oscillator. Let $\xi_{n}$ be the amplitude for mode $n$. Then

[TABLE]

We may solve this differential equation in steady-state and fix the reference phase to find

[TABLE]

The power dissipated is

[TABLE]

With equation (6.43) this becomes

[TABLE]

where

[TABLE]

If the tides are driven by the rotational energy of the planet then $\omega$ is the planet’s rotation frequency. In many cases however the tides are driven by either orbital eccentricity, in which case the forcing frequency is just the orbital frequency (Arras & Socrates, 2010). In this more generally interesting case, the driving frequency is

[TABLE]

To compare, the highest resonant frequency occurs when $n=0$ and is

[TABLE]

equation (6.25). Because $\beta q\approx 1$,

[TABLE]

for $l=2$. Now $N_{0}$ is the peak Brunt-Väisälä frequency in the radiative zone where the temperature gradient is substantially subadiabatic, so

[TABLE]

This gives

[TABLE]

Now $l_{r}\approx h\ln(P_{1}/P_{\mathrm{t,i}})$ so

[TABLE]

This is quite close to the orbital frequency so we expect that resonances are not uncommon. Note also that it is somewhat greater than the orbital frequency, so there will generally be modes with frequencies lower than the orbital frequency which are pulled upward towards it by thermomechanical feedback.

The precise shape and spacing of the resonances depends on our ansatz for the Brunt-Väisälä frequency in the radiative zone and so it is not useful to make predictions which depend on our chosen form. If we instead average over resonances, we note that the forcing integral falls exponentially in $n$ while the resonances fall off as a power law in $n$ so that the $n=0$ resonance always dominates on average. Individual systems far from the $n=0$ resonance may exhibit significant dissipation by virtue of sitting directly on a higher resonance but we expect this to be rare. So we assume that $\omega$ and $\omega_{n}$ are of the same order and write the net power

[TABLE]

This may be converted to a flux as

[TABLE]

7 Equilibrium Radius

To first order, suppose that $l_{r}\approx h$. The tidal flux is then given by

[TABLE]

Inserting equation (6.54) yields

[TABLE]

If the radiative zone is shallow but dominates the mass above its base then $m\approx M$. In addition,

[TABLE]

so

[TABLE]

Noting that

[TABLE]

we find

[TABLE]

Recalling equation (2.32) and neglecting the logarithmic correction to $\Sigma_{\mathrm{t}}$ owing to the motion of the zone boundary we write

[TABLE]

With the fiducial values and $h\approx 0.2r$ this flux produces an expansion at a rate

[TABLE]

which is sufficient to produce expansion of order $R_{\mathrm{J}}$ over million-year timescales.

As discussed in section 5 the expansion eventually increases the escaping flux to match the generated flux, so the expansion does not continue forever. The relevant dimensionless parameter for this equilibrium is the ratio $F_{\mathrm{t}}$ to $F_{\mathrm{i,i}}$, which is the unperturbed $F_{\mathrm{i}}$. Recall from equation (2.31) that

[TABLE]

Inserting equation (2.30) yields

[TABLE]

Thus

[TABLE]

Inserting equation (5.10) we get

[TABLE]

If the expansion is small then the luminosity ratios may be replaced by flux ratios. If tidal heating is significant, the perturbed flux escaping from the interior of the planet $F_{\mathrm{i,f}}\approx F_{\mathrm{t}}$, so

[TABLE]

Using equations (6.82) and (6.83) we find

[TABLE]

Note that $R$ and $\mathcal{R}$ must be the equilibrium radii here, not pre-expansion radii. Note that the absolute magnitude of the opacity does not enter our final expression. This cancellation is due to the assumption of efficient convection.

The quantity

[TABLE]

is of key importance to the nature of the solution. If $f<1$ then the solution is stable, meaning that a system initially perturbed away from this equilibrium solution returns to it over time. What this means physically is that an increase in the generated flux increases the temperature at the base of the radiative zone enough that the flux which escapes increases by more, leading to a negative feedback loop. If $f>1$ then the solution is unstable, meaning that an increase in the generated flux increases the temperature in the radiative zone by less than what is required to allow that additional flux to escape, leading to a positive feedback loop. If a system has $f>1$ then either the initial perturbation has the radiative zone deep enough that the runaway process keeps increasing the flux until some of our assumptions break down or the initial perturbation has the radiative zone shallow enough that the runaway process prevents it from migrating inward causing it to stay where it initially forms. If the prefactor on the right side of equation (7.16) exceeds unity then the unstable branch is always the relevant one because the initial perturbation must yield a ratio of at least unity in order to cause a radiative zone to form. If the prefactor is less than one then any radiative zone formed simply stays where the initial perturbation produces it.

Despite the uncertainties in the precise numbers involved, our results are fairly robust because the right-hand side of equation (7.16) scales as $\tau_{\mathrm{orbit}}^{-10}$, and this scaling only becomes stronger once $f$ is taken into account. For orbital periods shorter than of about $10\text{\,}\mathrm{d}$, with some uncertainty, the flux generated may exceed the flux escaping from the centre of the planet. When the tidal flux dominates the expansion is given by equation (6.88) and may be approximated by

[TABLE]

Neglecting the logarithmic dependence, this may be combined with equation (7.16) to give

[TABLE]

Even though the dependence on the flux ratio is small, the ratio itself can be quite large, particularly at smaller periods. Many cases, such as $a=0,b=2$ or $a=1,b=1$, have $f>1$ and so orbital periods of order 30\text{\,}\mathrm{d}$\Pi^{1/5}$ suffice to cause expansion of order $R_{0}$. It is difficult to say more because many of our approximations break down at this point. If $f<1$, as can be achieved for example with $a=4,b=2$, $\Delta R/R\propto\tau_{\mathrm{orbit}}^{-64/3}$, orbital periods of order 20\text{,}\mathrm{d}$\Pi^{3/16}$ suffice to produce unit expansion. This is consistent with most of the known cases of highly inflated Jupiter-mass planets, assuming $\Pi\approx e\approx 0.1$. For low-mass planets, the lower surface gravity makes larger expansion more feasible. This has been recently observed (Bakos et al., 2016). Stronger claims are difficult to make analytically because the dependence of the flux ratio on the specifics of the opacity are quite severe and the detailed compositions of the atmospheres of exoplanets at intermediate depths remain largely unknown. Precision studies of this thermomechanical feedback will likely require numerical tools in all but the simplest cases.

At sufficiently low masses, large flux ratios become impossible to attain given the factor of $\mathcal{M}^{-3}$ in equation (7.16). At this point further expansion is impossible. Likewise at large enough radii the central adiabat disappears so that much of this analysis becomes invalid. We do not expect this to be a limiting factor, however. At short periods Roche lobe overflow becomes a substantial barrier. Substantial changes in the opacity may also occur, particularly if the Kramer regime becomes relevant, and this may invalidate much of the analysis too. In addition at large radii the neglected factors of $R$ in converting from luminosities to fluxes become relevant and these act to limit the expansion.

8 Energetic Timescales

If the tides are eccentricity-driven, it is important to consider the timescale over which the orbit circularises. It suffices to the level of accuracy of interest to note that the energy which may be extracted from an orbit of eccentricity $e<1$ is of order $e^{2}\varv_{\mathrm{orbit}}^{2}M$. The circularisation timescale is therefore

[TABLE]

where we have taken $\Pi\approx e$ and used our fiducial values for all parameters other than $h/r$, which we have taken to be $0.2$. From this it is clear that most systems of interest can be eccentricity-driven for billion-year timescales, even if they require shorter periods than the fiducial.

If the tides are driven by the planet’s rotation they generally have many orders of magnitude less energy to draw from, and so are not sustained on the timescales of interest. They may still produce bloating, but not for long enough to be easily observable.

9 Comparison

For comparison with other mechanisms it is useful to compute the Love number associated with these g-modes. The Love number is defined in terms of the power and perturbing tidal potential as

[TABLE]

(Ogilvie, 2014). Summing in quadrature over all $m$ and using $l=2$ and equation (7.8) we find

[TABLE]

For comparison, inertial waves result in a frequency-averaged value of

[TABLE]

(Ogilvie, 2013), where $R_{\mathrm{c}}$ is the core radius. Likewise viscoelastic dissipation in a solid core potentially yields $\Im(k)\approx 10^{-6}$, scaling with the elastic properties of the core (Guenel et al., 2014; Remus et al., 2012) and once more averaging over frequency. These mechanisms are therefore of comparable order, depending on precisely which fiducial value is chosen. The reason the g-mode mechanism inflates planets more readily despite having comparable or somewhat less power dissipation is that it heats primarily near the surface where the radius is more easily perturbed.

10 Conclusions

We have characterised the response of heavily insolated Jupiter-like planets to tidal heating for a wide range of tidal heating models. A necessary condition for significant bloating of these planets is deep heating. We find that tidal heating, either directly through tide–core interactions or indirectly through resonance-sensitive migration of radiative zones induced by g-modes, is of the right order of magnitude to induce the observed bloating if it is sufficiently deep. We have further shown that nearly every tidal heating model results in deep heating so long as the atmosphere is sufficiently irradiated. This explains the observation of substantially bloated hot Jupiters with a physically reasonable orbital period cutoff for such effects. The migration of interior radiative zones provides a natural explanation for the matching of tidal frequency with orbital frequency despite the observed wide range of orbital frequencies of hot Jupiters.

This entire analysis hinges on there being a luminosity perturbation to start. This luminosity then produces a self-sustaining interior radiative zone which dissipates substantially more heat. The initial perturbation may come from non-linear instabilities and so may provide an indirect probe of these effects. It may also come from planetary migration. If a planet migrates inward and if opacity falls as a result, the incident flux can temporarily force the creation of a radiative zone while the convection zone adjusts to the reduced flux it must carry. The g-mode hysteresis described in this paper may then prevent the zone from disappearing, even if its location shifts to better match resonance. Inflated planets may therefore carry a record of their migration histories.

Finally, the thermomechanical feedback mechanism we propose highlights the importance of considering dynamical effects across many timescales. Feedback is possible both from short timescales to long, as in tidal heating, and from long timescales to short, as in the dynamical tuning of g-modes. By their very nature couplings across so many scales are difficult to track down and so there may be many more which have yet to be discovered.

Acknowledgements

We thank Sterl Phinney for productive conversations on planetary atmospheres and heated planets. ASJ acknowledges support from the Goldwater scholarship and the Marshall scholarship. CAT thanks Churchill College for his fellowship.