Chaotic Tides in Migrating Gas Giants: Forming Hot and Transient Warm Jupiters via High-Eccentricity Migration

TL;DR

This paper investigates how chaotic dynamical tides influence high-eccentricity migration, leading to the formation of hot and warm Jupiters with specific orbital characteristics, and enhances understanding of their observed properties.

Contribution

It introduces the role of chaotic tides in high-eccentricity migration, showing their impact on orbit evolution, planet survival, and the formation of hot and warm Jupiters.

Findings

Chaotic tides can rapidly shrink planetary orbits during high-eccentricity phases.

This mechanism produces very eccentric warm Jupiters on short timescales.

Chaotic tides increase the formation fraction of hot Jupiters by preventing tidal disruption.

Abstract

High-eccentricity migration is an important channel for the formation of hot Jupiters (HJs). In particular, Lidov-Kozai (LK) oscillations of orbital eccentricity/inclination induced by a distant planetary or stellar companion, combined with tidal friction, have been shown to produce HJs on Gyr timescales, provided that efficient tidal dissipation operates in the planet. We re-examine this scenario with the inclusion of dynamical tides. When the planet's orbit is in a high-eccentricity phase, the tidal force from the star excites oscillatory f-modes and r-modes in the planet. For sufficiently large eccentricity and small pericentre distance, the mode can grow chaotically over multiple pericentre passages and eventually dissipate non-linearly, drawing energy from the orbit and rapidly shrinking the semi-major axis. We study the effect of such chaotic tides on the planet's orbital evolution. We find that this pathway produces very eccentric ($e\gtrsim 0.9$) warm Jupiters (WJs) on short timescales (a few to 100 Myrs). These WJs efficiently circularize to become HJs due to their persistently small pericentre distances. Chaotic tides can also save some planets from tidal disruption by truncating the LK eccentricity oscillations, significantly increasing the HJ formation fraction for a range of planet masses and radii. Using a population synthesis calculation, we determine the characteristics of WJs and HJs produced in this scenario, including the final period distribution, orbital inclinations and stellar obliquities. Chaotic tides endow LK migration with several favorable features to explain observations of HJs. We expect that chaotic tides are also important in other flavours of high-$e$ migration.

Figures

Tables

Table 1. Table 1: Properties of l = m = 2 𝑙 𝑚 2 l=m=2 f-modes and inertial modes for a γ = 2 𝛾 2 \gamma=2 polytrope and a γ = 5 / 3 𝛾 5 3 \gamma=5/3 polytrope. The mode frequencies ω α subscript 𝜔 𝛼 \omega_{\alpha} (in the rotating frame), ϵ α subscript italic-ϵ 𝛼 \epsilon_{\alpha} (equation 5 ), and σ α subscript 𝜎 𝛼 \sigma_{\alpha} (in the inertial frame), and the planetary rotation rate Ω s subscript Ω 𝑠 \Omega_{s} are in units of ( G ​ M p / R p 3 ) 1 / 2 superscript 𝐺 subscript 𝑀 𝑝 superscript subscript 𝑅 𝑝 3 1 2 (GM_{p}/R_{p}^{3})^{1/2} , i.e. ω α = ω ¯ α ​ ( G ​ M p / R p 3 ) 1 / 2 subscript 𝜔 𝛼 subscript ¯ 𝜔 𝛼 superscript 𝐺 subscript 𝑀 𝑝 superscript subscript 𝑅 𝑝 3 1 2 \omega_{\alpha}=\bar{\omega}_{\alpha}(GM_{p}/R_{p}^{3})^{1/2} . The tidal overlap integral Q α subscript 𝑄 𝛼 Q_{\alpha} is defined by equation ( 4 ).

Model Type	Mode	${\overline{\omega }}_{\alpha }$	${\overline{ϵ}}_{\alpha }$	${\overline{\sigma }}_{\alpha }$	$Q_{\alpha }$
$\gamma =2$ polytrope	f-mode	$1.22-{\overline{\Omega }}_s$	$1.22$	$1.22+{\overline{\Omega }}_s$	0.56
	i-mode 1 ($j=3$)	0.56 ${\overline{\Omega }}_s$	0.28 ${\overline{\Omega }}_s$	2.56 ${\overline{\Omega }}_s$	$0.015{\overline{\Omega }}_s^2$
	i-mode 2 ($j=3$)	-1.1 ${\overline{\Omega }}_s$	-0.55${\overline{\Omega }}_s$	0.9 ${\overline{\Omega }}_s$	0.01 ${\overline{\Omega }}_s^2$
$\gamma =5/3$ polytrope	f-mode	$1.46-{\overline{\Omega }}_s$	$1.46$	$1.46+{\overline{\Omega }}_s$	0.49

Table 2. Table 2: Outcomes of LK migration with chaotic tides for three planet models. The fraction of systems that undergo chaotic evolution is denoted by ℱ mig subscript ℱ mig \mathcal{F}_{\rm mig} ; these planets are either tidally disrupted ( ℱ dis subscript ℱ dis \mathcal{F}_{\rm dis} ) or exit chaotic tides as eccentric WJs ( ℱ WJ = ℱ mig − ℱ dis subscript ℱ WJ subscript ℱ mig subscript ℱ dis \mathcal{F}_{\rm WJ}=\mathcal{F}_{\rm mig}-\mathcal{F}_{\rm dis} ). The fraction of planets that circularize to e = 0.1 𝑒 0.1 e=0.1 with 1 Gyr is ℱ HJ subscript ℱ HJ \mathcal{F}_{\rm HJ} . We compare our results with semi-analytical predictions using the method described in Section  6 .

	1 $R_J$, 1$M_J$			1.6 $R_J$, 1$M_J$			1.6 $R_J$, 0.3$M_J$
	Rate	Prediction		Rate	Prediction		Rate	Prediction
${ℱ}_{\mathrm{mig}}$	14.0 %	13.7 %		15.4 %	14.2 %		16.7 %	15.2%
${ℱ}_{\mathrm{dis}}$	10.1%	$\le$11.9%		11.5 %	$\le$11.9 %		13.2 %	$\le$13.1 %
${ℱ}_{\mathrm{WJ}}$	3.9%	$\ge$1.8 %		3.9%	$\ge$2.3 %		3.6 %	$\ge$ 2.1 %
${ℱ}_{\mathrm{HJ}}$	3.9%	—		3.8%	—		2.2 %	—

Equations

$$U(\mathbf{r},t)=-G{M}_{\ast }\sum _{lm}\frac{W_{lm}{r}^l}{D^{l+1}}{e}^{-im\Phi (t)}{Y}_{lm}(\theta ,{\varphi }_i),$$

$$\left[\begin{array}[]{c}\mn@boldsymbol{\xi}\\ {\partial\mn@boldsymbol{\xi}/\partial t}\end{array}\right]=\sum_{\alpha}b_{\alpha}(t)\left[\begin{array}[]{c}\mn@boldsymbol{\xi}_{\alpha}(\bm{r})\\ -i\omega_{\alpha}\mn@boldsymbol{\xi}_{\alpha}(\bm{r})\end{array}\right].$$

$${\dot{b}}_{\alpha }+i{\omega }_{\alpha }{b}_{\alpha }=\frac{i{M}_{\ast }{W}_{lm}{Q}_{\alpha }}{2{ϵ}_{\alpha }{D}^{l+1}}{\text{e}}^{-im\Phi (t)+im{\Omega }_{s}t},$$

$$Q_{\alpha }=\int d^{3}x\rho \text{mn@boldsymbol}{\xi }_{\alpha }^{\star }\cdot \nabla (r^l{Y}_{lm})$$

$${ϵ}_{\alpha }={\omega }_{\alpha }+\int d^{3}x\rho \text{mn@boldsymbol}{\xi }_{\alpha }^{\star }\cdot (\text{i}{\mathbf{\Omega }}_{\mathbf{s}}\times {\mathbf{\xi }}_{\mathbf{\alpha }}).$$

$$\int d^{3}x\rho (\mathbf{r})|\text{mn@boldsymbol}{\xi }_{\alpha }(\mathbf{r}){|}^2=1.$$

$$b_{\alpha }(t)={\text{e}}^{-i{\omega }_{\alpha }t}{\int }_{t_0}^t\frac{i{M}_{\ast }{W}_{lm}{Q}_{\alpha }}{2{ϵ}_{\alpha }D(t^{\prime }{)}^{l+1}}{\text{e}}^{i{\sigma }_{\alpha }{t}^{\prime }-im\Phi (t^{\prime })}d{t}^{\prime }+b_{\alpha }(t_0),$$

$${\sigma }_{\alpha }={\omega }_{\alpha }+m{\Omega }_s$$

$$t_k=t_{k-1}+\frac{1}{2}\left(P_{k-1}+P_k\right),$$

$$\Delta b_{\alpha ,k}={\int }_{-P_{k-1}/2}^{P_k/2}\frac{i{M}_{\ast }{W}_{lm}{Q}_{\alpha }}{2{ϵ}_{\alpha }D(t^{\prime }{)}^{l+1}}{\text{e}}^{i{\sigma }_{\alpha }{t}^{\prime }-im\Phi (t^{\prime })}d{t}^{\prime }.$$

$$b_{\alpha }=b_{\alpha ,0}+{\text{e}}^{-\text{i}{\omega }_{\alpha }{t}_k}\sum _{j=1}^k{\text{e}}^{\text{i}{\sigma }_{\alpha }(t_{j-1}+P_{j-1}/2)}\Delta b_{\alpha ,j}.$$

$$b_{\alpha ,k}=b_{\alpha ,k-1}{\text{e}}^{-\text{i}{\omega }_{\alpha }(P_{k-1}+P_k)/2}+\Delta b_{\alpha ,k}{\text{e}}^{\text{i}(m{\Omega }_s{t}_k-{\sigma }_{\alpha }{P}_k/2)}.$$

$$c_{\alpha ,k}\equiv b_{\alpha ,k}{\text{e}}^{-\text{i}(m{\Omega }_s{t}_k+{\sigma }_{\alpha }{P}_k/2)}.$$

$$c_{\alpha ,k}=(c_{\alpha ,k-1}+\Delta c_{\alpha ,k}){\text{e}}^{-\text{i}{\sigma }_{\alpha }{P}_k},$$

$$E_{\alpha ,k}=2{\sigma }_{\alpha }{ϵ}_{\alpha }|c_{\alpha ,k}{|}^2,$$

$${\tilde{E}}_{\alpha ,k}\equiv |{\tilde{c}}_{\alpha ,k}{|}^2\equiv \frac{E_{\alpha ,k}}{|E_{B,0}|}.$$

$$\Delta {\tilde{E}}_k=\sum _{\alpha }\Delta {\tilde{E}}_{\alpha ,k}=\sum _{\alpha }(|{\tilde{c}}_{\alpha ,k-1}+\Delta {\tilde{c}}_{\alpha }{|}^2-|{\tilde{c}}_{\alpha ,k-1}{|}^2).$$

$${\tilde{E}}_{B,k}={\tilde{E}}_{B,k-1}-\Delta {\tilde{E}}_k={\tilde{E}}_{B,k-1}-\sum _{\alpha }{\tilde{E}}_{\alpha ,k},$$

$$\frac{P_k}{P_0}={\left(\frac{{\tilde{E}}_{B,0}}{{\tilde{E}}_{B,k}}\right)}^{3/2},$$

$$|\Delta {\tilde{c}}_{\alpha }{|}^2=\frac{\Delta E_{\alpha }}{|E_{\mathrm{B},0}|}.$$

$$\Delta E_{\alpha }$$

$$=\frac{G{M}_p^2}{R_p}{\left(\frac{M_{\ast }}{M_{\mathrm{tot}}}\right)}^2{\eta }^{-6}T(\eta ,{\sigma }_{\alpha }/{\Omega }_p,e),$$

$$\eta \equiv \frac{r_{\mathrm{p}}}{r_{\mathrm{tide}}},$$

$${\Omega }_p\equiv {\left(\frac{G{M}_{\mathrm{tot}}}{r_{\mathrm{p}}^3}\right)}^{1/2}={\left(\frac{G{M}_p}{R_p^3}\right)}^{1/2}{\eta }^{-3/2}.$$

$$r_{\mathrm{p}}\ge r_{\mathrm{p},\mathrm{dis}}=2.7r_{\mathrm{tide}},$$

$$T=2{\pi }^2\frac{{\sigma }_{\alpha }}{{ϵ}_{\alpha }}(Q_{\alpha }{K}_{lm}{)}^2,$$

$$K_{lm}=\frac{W_{lm}}{2\pi }{\int }_{-P/2}^{P/2}dt{\left(\frac{r_{\mathrm{p}}}{D(t)}\right)}^{l+1}{\text{e}}^{\text{i}{\sigma }_{\alpha }t-\text{i}m\Phi (t)}.$$

$$K_{22}\simeq \frac{2z^{3/2}\exp (-2z/3)}{\sqrt{15}}\left(1-\frac{\sqrt{\text{uppi}}}{4\sqrt{z}}\right){\eta }^{3/2}.$$

$$|\Delta {\hat{P}}_{\alpha }|={\sigma }_{\alpha }\Delta P$$

$$\simeq \frac{6\pi {\sigma }_{\alpha }/{\Omega }_p}{(1-e{)}^{5/2}}{\left(\frac{M_p}{M_{\ast }}\right)}^{2/3}{\eta }^{-5}T(\eta ,{\sigma }_{\alpha }/{\Omega }_p,e).$$

$${\Omega }_{\mathrm{s}}={\Omega }_{\mathrm{ps}}\equiv \frac{f_2(e)}{(1-e^2{)}^{3/2}{f}_5(e)}n,$$

Full text

††pubyear: 2018††pagerange: Chaotic Tides in Migrating Gas Giants: Forming Hot and Transient Warm Jupiters via High-Eccentricity Migration–References

1 Introduction

Despite over two decades of observations, the origins of hot Jupiters (HJs, giant planets with orbital periods $\lesssim 10$ days) remain puzzling (see Dawson & Johnson, 2018, for a review). The materials and conditions necessary to form such massive planets are not thought to exist so close to a protostar. Some have suggested that HJs could form in-situ when a protoplanet migrates inward and accumulates a gas envelope (Boley et al., 2016; Batygin et al., 2016). However, most theories propose that a fully formed gas giant migrated toward its host star either via interaction with the protoplanetary disk (e.g. Lin et al., 1996; Kley & Nelson, 2012) or through tidal decay and circularization of a high-eccentricity orbit, a process termed “high-eccentricity migration." Equally puzzling is the origin of warm Jupiters (WJs), giant planets with periods between about 10 days and 200 days. Although these close-in giant planets (HJs and WJs) represent a relatively small population of exoplanetary systems compared to the more abundant super-earths, their dynamical history can potentially shape the architecture of planetary systems.

High-eccentricity migration is an appealing avenue for HJ formation. In this scenario, a gas giant is excited into a highly eccentric orbit via interactions with other planets or with a distant stellar companion. Strong planet-planet scattering (Rasio & Ford, 1996; Chatterjee et al., 2008; Jurić & Tremaine, 2008), various forms of secular interactions (Wu & Lithwick, 2011; Hamers & Portegies Zwart, 2016; Petrovich, 2015b), or a combination of both (Nagasawa et al., 2008; Beaugé & Nesvorný, 2012), can produce very eccentric gas giants. An highly inclined stellar or planetary companion can excite “Lidov-Kozai” (LK) oscillations (Lidov, 1962; Kozai, 1962; Naoz, 2016, for a review), pushing the eccentricity of the giant planet to near unity (Wu & Murray, 2003; Fabrycky & Tremaine, 2007; Naoz et al., 2012; Correia et al., 2012; Petrovich, 2015a; Anderson et al., 2016). In all cases, the orbit circularizes and decays to a period of few days due to the tidal dissipation within the planet.

High-eccentricity migration can account for several puzzling characteristics of observed HJs. For instance, tidal orbital decay produces the observed pile-up of HJs at an orbital period of $\sim 3$ days (e.g. Santerne et al., 2016), corresponding to semi-major axes of a few times the Roche radii of the planets. Lidov-Kozai migration driven by an inclined stellar/planetary companion can naturally generate significant “spin-orbit” misalignments between the rotation axis of the host star and the orbital angular momentum axis of the planet, as observed in many HJ systems (Hébrard et al., 2008; Narita et al., 2009; Winn et al., 2009; Triaud et al., 2010; Albrecht et al., 2012; Winn & Fabrycky, 2015); the chaotic evolution of the stellar spin axis driven by the changing planetary orbit plays a dominant role in producing the observed stellar oblqiuities (Storch et al., 2014; Anderson et al., 2016; Storch et al., 2017). Other flavors of high-$e$ migration can also produce appreciable spin-orbit misalignments (e.g. Lithwick & Wu, 2014; Petrovich, 2015b; Teyssandier et al., 2018). In addition, direct RV and AO-imaging searches have shown that a large fraction of HJs have external massive planet companions at 5-20 au or distant stellar companions (50-2000 au) (Knutson et al., 2014; Ngo et al., 2015; Wang et al., 2015; Bryan et al., 2016), suggesting that dynamical interactions play a role in HJ formation. Huang et al. (2016) showed that most HJs do not have any detectable neighbors (while half of WJs are closely flanked by small neighbors), again suggesting that a significant fraction of HJs may have formed through high-$e$ migration.

A significant uncertainty in any high-$e$ migration scenario is tidal dissipation in the planet. Regardless of how the planet attains its high eccentricity, efficient tidal dissipation is necessary in order to circularize the planet’s orbit and to bring it from a semi-major axis of several au’s to $\lesssim 0.05$ au. So far, almost all studies of high-$e$ migration have relied on the parameterized weak-friction theory of equilibrium tides (Alexander, 1973; Hut, 1981). In order achieve tidal circularization within $\sim 5$ Gyrs in high-$e$ migration scenarios, the giant planet generally should be more dissipative than Jupiter by more than an order of magnitude (e.g. Socrates et al., 2012; Petrovich, 2015a; Anderson et al., 2016).

In general, the response of a fluid body (such as a giant planet) to tidal perturbation involves the excitation of internal waves and modes at different frequencies (see, e.g., Ogilvie, 2014). For a body in a sufficiently eccentric orbit, oscillatory modes excited at pericentre by the star’s tidal potential can grow chaotically in amplitude, rapidly draining energy from the orbit over multiple close passages (Mardling, 1995a, b). This effect is important for high-eccentricity migration because the same conditions that allow a planet to migrate quickly – small pericentre distance and large eccentricity – are those that result in chaotic tidal behaviour. Ivanov & Papaloizou (2004) first studied the possibility that chaotic dynamical tides could speed up the circularization of eccentric gas giants. They developed an iterative algebraic map to simply follow mode evolution over many orbits from one pericentre passage to the next. More recently, Vick & Lai (2018) used a similar map (generalized to include linear mode dampling) to quantify the various dynamical behaviours of the “mode + eccentric orbit” system (see also Section 2 below); they suggested that chaotic tides could quickly produce eccentric warm Jupiters (WJs), with a semi-major axis between 0.1 and 1 au, that efficiently circularize to HJs. Wu (2018) also studied the orbital evolution of a planet experiencing chaotic tides, and uncovered some of the key features of LK high-$e$ migration with chaotic tides – we will examine these features in more details in the later sections of this paper.

In this paper, we present a comprehensive study of high-eccentricity migration of giant planets, incorporating a detailed model of chaotic tidal evolution. We focus on migration via the LK effect with a stellar companion. Although this flavour of high-$e$ migration may only account for a fraction of the HJ population, it is most suitable for systematic study since the initial conditions (such as the orbital properties of the stellar perturber) are reasonably well-justified, and it serves as a benchmark for other high-$e$ migration scenarios. Indeed, several aspects of our results presented in this paper can be applied to other flavours of high-$e$ migration (see Teyssandier et al., 2018, which studies dynamical tides in the “secular chaos” scenario). Comparing to our previous study of LK migration with static tides (Anderson et al., 2016; Muñoz et al., 2016), we show that including chaotic tides leads to a number of “favorable” features for the formation of close-in giant planets: e.g., it naturally produces eccentric WJs on short timescales ($\lesssim 10-100$ Myrs) and speeds up the formation of HJs; it generates a wider HJ period distribution; and it prevents some gas giants from undergoing tidal disruption, thereby increase the HJ formation fraction. Overall, chaotic tides make high-$e$ migration a more promising mechanism for producing HJs

The structure of the paper is as follows. In Section 2 we revisit the derivation of the iterative map for the evolution of planetary oscillation modes (dynamical tides), starting from hydrodynamical equations and including the effect of planetray rotation; we examine the different dynamical behaviours of the “modes

• eccentric orbit” system and discuss the conditions for various modes in the planet to become chaotic. In Section 3 we describe our model for coupling LK oscillations with chaotic tides. Section 4 discusses the key features and signposts of LK migration with chaotic tides and includes a number of analytical results that characterize chaotic tidal migration. In Section 5 we present our population synthesis study for planets undergoing LK oscillations with chaotic tides. In Section 6 we give analytical calculations for the HJ formation and tidal disruption fractions from this migration mechanism and compare them with our numerical results. We conclude in Section 7 with a summary of results and discussion of their implications.

2 Chaotic Dynamical Tides in Giant Planets in Eccentric Orbits

Consider a planet (mass $M_{p}$ and radius $R_{p}$) in a highly eccentric orbit around a star (mass $M_{*}$). At each pericentre passage, the time-varying gravitational potential of the star excites oscillatory modes (e.g. f-modes and r-modes) in the planet. Over multiple pericentre passages, the orbit and the modes exchange energy and angular momentum. In general, this exchange can occur in either direction depending on the phase of the modes (Kochanek, 1992; Mardling, 1995a). Provided the oscillation amplitudes are sufficiently small, one can determine the dynamical behaviour of the system using linear hydrodynamics (e.g. Mardling, 1995a, b; Lai, 1996; Kumar & Goodman, 1996; Lai, 1997). However, following the evolution of the system over many orbits can be challenging because of the short timescale associated with the pericentre passage. Fortunately, when the orbit is eccentric enough, such that the mode and the orbital evolution mainly occur at pericentre, the full hydrodynamic solution of the “eccentric orbit + oscillation modes" system can be expressed as an iterative map (Ivanov & Papaloizou, 2004; Vick & Lai, 2018).

2.1 Hydrodynamics of Tidally Forced Oscillations & Iterative Map

In Vick & Lai (2018), we presented a simple iterative map to describe the long-term evolution of an “eccentric orbit + oscillation modes" system. Our map neglected the rotation of the planet and assumed that the change in the mode amplitude during each pericentre passage is constant throughout the orbital evolution. Giant planets can have substantial rotations, which affect their mode properties and introduce a few subtleties to the mapping. Also, during the long-term evolution of the planet’s orbit driven by an external companion (the Lidov-Kozai effect), the pericentre distance can change appreciably, which affects mode excitation. Here we present a general derivation of the iterative map from linear hydrodynamics, allowing for the planet to rotate and for the pericentre distance (and thus the change in mode amplitude at pericentre) to vary between passages.

The stellar gravitational potential that excites oscillations in the planet is given by

[TABLE]

where ${\bf r}=(r,\theta,\phi_{i}=\phi+\Omega_{s}t)$ is the position vector (in spherical coordinates) relative to the centre of mass of the planet (the azimuthal angle $\phi$ is measured in the rotating frame of the planet, and $\Omega_{s}$ is the rotation rate of the planet111We assume the spin axis of the planet is aligned with the orbital angular momentum axis throughout this paper.). In equation (1), $D(t)$ is the time-varying separation between the star and planet, $\Phi(t)$ is the orbital true anomaly, and $W_{lm}$ is a constant defined in Press & Teukolsky (1977). The dominant (quadrupole) terms have $l=|m|=2$ with $W_{2\pm 2}=\sqrt{3\uppi/10}$.

The response of the planet to tidal forcing is described by the Lagrangian displacement, $\bm{\xi}(\bm{r},t)$. A free oscillatory mode of frequency $\omega_{\alpha}$ (in the rotating frame) has the form $\bm{\xi}_{\alpha}(\bm{r},t)=\bm{\xi}_{\alpha}(\bm{r})\,\text{e}^{-\text{i}\omega_{\alpha}t}\propto\text{e}^{\text{i}m\phi-i\omega_{\alpha}t}$, where $\alpha$ is the mode index. We expand $\bm{\xi}(\bm{r},t)$ in terms of the eigenmodes in the phase space (Schenk et al., 2002):

[TABLE]

The linear fluid dynamics equations then reduce to a set of first-order differential equations (Lai & Wu, 2006),

[TABLE]

where

[TABLE]

is the dimensionless tidal overlap integral (in units where $G=M_{p}=R_{p}=1$), and

[TABLE]

Note that in the slow rotation limit, $\epsilon_{\alpha}$ is simply the mode frequency for a non-rotating planet, $\omega_{\alpha}(0)$. When $\epsilon_{\alpha}\gg\Omega_{s}$, the first-order correction to the mode frequency due to rotation is $\omega_{\alpha}-\omega_{\alpha}(0)\approx-mC_{\alpha}\Omega_{s}=-\int d^{3}x\,\rho\mn@boldsymbol{\xi}_{\alpha}^{\star}\cdot(\text{i}\bm{\Omega_{s}}\times\bm{\xi_{\alpha}})$ (e.g. Unno et al., 1989). In equations (3)-(5), the eigenmode is normalized according to

[TABLE]

The general solution to equation (3) is

[TABLE]

where

[TABLE]

is the mode frequency in the inertial frame.

For a highly eccentric orbit, we can assume that any energy transfer between a mode and the orbit occurs at pericentre. We can then manipulate equation (7) into a discrete form by defining $t_{k}$ as the time at apocentre after the $k$-th pericentre passage, i.e.

[TABLE]

where $P_{k}$ is the orbital period after the $k$-th pericentre passage. We also define

[TABLE]

Physically, this is the real change in mode amplitude during the $k$-th pericentre passage. Equation (7) becomes

[TABLE]

We assume the initial condition $b_{\alpha,0}=0.$ Equation (11) can be rewritten in an iterative form:

[TABLE]

We now define the mode amplitude in the inertial frame and shift the index $k$ to count pericentre passages using

[TABLE]

Physically, $c_{\alpha,k}$ is the mode amplitude just before the $(k+1)$-th pericentre passage. Equation (12) then becomes

[TABLE]

where $\Delta c_{\alpha,k}=\Delta b_{\alpha,k}$. Equation (14) has a straight-forward physical interpretation: when the planet passes through pericentre (the $k$-th passage), the mode amplitude changes by $\Delta c_{\alpha,k}$, and the orbital period changes to $P_{k}$; as the planet completes its orbit, the phase of the mode evolves due to the passage of time, and the complex mode amplitude gains a factor $\text{e}^{-\text{i}\sigma_{\alpha}P_{k}}$. Note that in general $\Delta c_{\alpha,k}$ depends on the “current" parameters of the system. If the pericentre distance and the shape of the orbit near pericentre remain unchanged over many passages, $\Delta c_{\alpha,k}$ is nearly constant from one passage to the next. In this case, equation (14) reduces to the result from Vick & Lai (2018) but without mode damping and with the phase evolution determined by the mode frequency in the inertial frame (see also Ivanov & Papaloizou, 2004).

The mode amplitude is directly related to the mode energy via

[TABLE]

where the factor of 2 arises because a mode with given $m\;(>0)$, and $\sigma_{\alpha}$ is physically identical to that with $(-m)$ and $(-\sigma_{\alpha})$ — we treat them as the same mode. We define the dimensionless mode energy, ${\tilde{E}}_{\alpha,k}$, in units of the initial binary (planetary) orbital energy, $|E_{B,0}|=GM_{*}M_{p}/(2a_{0})$ (where $a_{0}$ is the semi-major axis), and define the re-normalized mode amplitude $\tilde{c}_{\alpha,k}$ such that

[TABLE]

The total energy transfer to modes in the $k$-th passage is therefore

[TABLE]

From energy conservation, the orbital energy (${\tilde{E}}_{B,k}$) immediately after the $k$-th passage is given by

[TABLE]

and the corresponding orbital period is

[TABLE]

where ${\tilde{E}}_{B,0}=-1$, and $P_{0}$ is the initial period.

Equations (14) and (17)-(19) complete the map. The required inputs are the initial mode energy, ${\tilde{E}}_{\alpha,0}=|\tilde{c}_{\alpha,0}|^{2}$, and $\Delta\tilde{c}_{\alpha,k}$, the change in mode amplitude during each pericentre passage. This map accurately determines the evolution of the orbit and mode energies provided that the orbit is highly eccentric and the modes remain linear (i.e., the sum of the mode energies is much less than the binding energy of the planet).

The dimensionless change in mode amplitude during a pericentre passage (suppressing the subscript ‘k’), $\Delta\tilde{c}_{\alpha}$, depends on the orbital parameters (mainly the pericentre distance $r_{\rm p}$) and the mode properties. It is related to the energy transfer $\Delta E_{\alpha}$ by

[TABLE]

Note that in general the energy transfer to mode $\alpha$ in a pericentre passage depends on the amplitude and phase of “pre-existing" oscillations [see equation (17)]. In equation (20), $\Delta E_{\alpha}$ refers to the energy transfer when the mode has zero amplitude prior to the pericentre passage — in Vick & Lai (2018), we term this “the first passage." The energy transfer $\Delta E_{\alpha}$ for a parabolic encounter ($e\to 1$) was first derived by Press & Teukolsky (1977). This can be generalized to eccentric orbits and modified to include the effect of rotation (Lai, 1997; Fuller & Lai, 2012). For a pericentre distance $r_{\rm p}$, $\Delta E_{\alpha}$ can be written as (keeping only the $l=2$ terms)

[TABLE]

where $\eta$ is the ratio of $r_{\rm p}$ and the tidal radius $r_{\rm tide}$, i.e.

[TABLE]

with $M_{\rm tot}=M_{*}+M_{p}\simeq M_{*}$, and $\Omega_{\rm p}$ is the pericentre orbital frequency

[TABLE]

Note that $\sigma_{\alpha}/\Omega_{p}=\bar{\sigma}_{\alpha}\eta^{3/2}$, with $\bar{\sigma}_{\alpha}\equiv\sigma_{\alpha}(R_{p}^{3}/GM_{p})^{1/2}$. To avoid tidal disruption of the planet, we require (Guillochon et al., 2011)

[TABLE]

The dimensionless function $T$ is given by

[TABLE]

where $Q_{\alpha}$ is an overlap integral (see equation 4), and

[TABLE]

In equation (26), all quantities are in units where $G=M_{p}=R_{p}=1$ so that $K_{lm}$ is dimensionless.

To follow the evolution of the system over many orbits, it is useful to be able to compute $K_{lm}$ (and thus $\Delta E_{\alpha}$) efficiently. For the dominant $l=m=2$ prograde mode ($\sigma_{\alpha}>0$), Lai (1997) derived an approximation for $K_{22}$ that is accurate to within 2% for $(1-e)\ll 1$ and $z\equiv\sqrt{2}\sigma_{\alpha}/\Omega_{p}\gtrsim$ a few:

[TABLE]

As shown in Vick & Lai (2018), when the mode has no initial energy, i.e. $\tilde{E}_{\alpha,0}=0$, the important quantity that determines the dynamical behaviour of the system is $\Delta\hat{P}_{\alpha}=\sigma_{\alpha}\Delta P$, where $\Delta P=P_{1}-P_{0}\simeq-(3/2)P_{0}(\Delta E_{\alpha}/|E_{B,0}|)$ is the change in the initial orbital period ($P_{0}$) due to the energy transfer $\Delta E_{\alpha}$. From equation (21), we have

[TABLE]

Physically, $|\Delta\hat{P}_{\alpha}|$ is the phase change in the mode due to energy transfer at pericentre when there is no “pre-existing" mode energy.

2.2 Planetary Oscillation Modes

The gravitational potential of a star can excite many types of oscillation modes in the eccentric orbiting planet. Some of the modes will be more strongly excited than others and therefore more influential in long-term orbital evolution. The dominant modes are quadrupole, with $l=|m|=2$. The $m=0$ modes are generally much less important because $K_{l0}$ (equation 26) is small. In this paper, we adopt the convention $m>0$, so a mode with $\sigma_{\alpha}>0$ ($\sigma_{\alpha}<0$) is prograde (retrograde) with respect to the rotation of the planet in the inertial frame. For the rest of this paper, we only consider modes with $l=m=2$.

We adopt a simple giant planet model composed of a neutrally stratified fluid characterized by a $\gamma=2$ polytrope equation of state ($P\propto\rho^{2}$). Of the acoustic modes (with the restoring force being pressure), the f-mode (fundamental mode) has a frequency closest to the frequency of tidal forcing near pericentre, and the largest tidal overlap integral $Q_{\alpha}$. Other acoustic modes (p-modes) have higher frequencies and smaller tidal overlaps, and are much less important than the f-mode for energy transfer. For a rotating planet, the stellar potential can also excite i-modes [i.e. inertial modes, with the restoring force being the Coriolis force; see Xu & Lai (2017) and references therein]. In a gas giant that has stratification, g-modes (driven by buoyancy) may also be excited (Fuller, 2014). The presence of a large solid core also supports elastic modes which can “mix" with the f-mode (e.g. Fuller et al., 2014). We do not consider such planet models in this paper because there are considerable uncertainties in the stratification and core size of giant planets (e.g. Guillot, 2005).

Table 1 lists the key properties of the $l=m=2$ f-modes and the first two inertial modes for the $\gamma=2$ polytropic planet model. For comparison, the f-mode properties of the $\gamma=5/3$ polytrope are also given. The f-modes are calculated using the slow-rotation approximation, which gives $\omega_{\alpha}=\omega_{\alpha}(0)-mC_{\alpha}\Omega_{s}$, and we find $C_{\alpha}\simeq 0.5$ for both the $\gamma=2$ and $\gamma=5/3$ models. The i-mode result is from Xu & Lai (2017) (see their Table IV), based on calculations using a non-perturbative spectral code. Note that for $m=2$, the mode frequencies in the rotating frame ($\omega_{\alpha}$) and in the inertial frame ($\sigma_{\alpha}$) are related by $\sigma_{\alpha}=\omega_{\alpha}+2\Omega_{s}$.

To determine which of these modes produces the strongest dynamical tides, we calculate $|\Delta{\hat{P}}_{\alpha}|$ with equation (28) as a function of $r_{\rm p}$ for each of the modes (see Fig. 1). For definiteness, we assume that the planet rotates at the pseudosynchronous rate $\Omega_{\rm ps}$ given by the weak friction theory of equilibrium tides (e.g. Hut, 1981), i.e.

[TABLE]

where $n=(GM_{\rm tot}/a^{3})^{1/2}$ is the mean motion of the planet and

[TABLE]

For $(1-e)\ll 1$, $\bar{\Omega}_{\rm ps}=\Omega_{\rm ps}(GM_{p}/R_{p}^{3})^{-1/2}\simeq 1.17/\eta^{3/2}$, and the mode frequency in the inertial frame, $\sigma_{\alpha}$, is nearly independent of the orbital eccentricity. In this limit, $T(\eta,\sigma_{\alpha}/\Omega_{p},e)$ also has negligible eccentricity dependence, and $|\Delta{\hat{P}}_{\alpha}|\propto(1-e)^{-5/2}$ [see equation (28)]. The two panels of Fig. 1 shows that $|\Delta{\hat{P}}_{\alpha}|$ indeed increases with $e$ for all modes. We are most interested in systems where $r_{\rm p}$ is sufficiently small such that $|\Delta{\hat{P}}_{\alpha}|\gtrsim 1$ (the criterion for chaotic mode growth; see Section 2.3). For such systems, it is clear that the $l=m=2$ f-mode is most strongly excited at pericentre passages. For the remainder of the paper, we focus on the prograde $l=2$ f-mode and neglect the contributions to energy transfer from other modes.

2.3 Conditions for Chaotic Tides

When the planetary orbit is sufficiently eccentric, the change in mode amplitude in a pericentre passage, $\Delta c_{\alpha}$ (or equivalently the energy transfer in the “first" passage, $\Delta E_{\alpha}$, see equation 21), is constant over many pericentre passages.222This requires that the pericentre distance is constant, which in turn requires that the fractional change in the orbital angular momentum, $\Delta L/L$, remains small throughout orbital evolution. See equation (13) of Vick & Lai (2018). In this case, if the mode has negligible dissipation, the dynamical behaviour of the iterative map depends on ${\hat{P}}_{\alpha,0}=\sigma_{\alpha}P_{0}$, the initial mode energy $E_{\alpha,0}$, and $|\Delta{\hat{P}}_{\alpha}|$ (see equation 28), assuming that the initial mode energy $E_{\alpha,0}$ is negligible (see Section 2.4). In particular, a system with $|\Delta{\hat{P}}_{\alpha}|$ greater than a critical value, $\Delta\hat{P}_{\rm crit}$, behaves very differently from systems with $|\Delta{\hat{P}}_{\alpha}|\lesssim\Delta\hat{P}_{\rm crit}$. The exact value of $\Delta\hat{P}_{\rm crit}$ depends on the dimensionless orbital period, ${\hat{P}}_{\alpha,0}$. On average, $\Delta\hat{P}_{\rm crit}\sim 1$. The map exhibits three classes of behaviour (Vick & Lai, 2018):

(i) When $|\Delta{\hat{P}}_{\alpha}|\lesssim\Delta{\hat{P}}_{\rm crit}\sim 1$, the mode energy undergoes low-amplitude oscillations. A small amount of energy (of order $\Delta E_{\alpha}$) is transferred back and forth between the mode and the orbit over multiple pericentre passages. An example of this behaviour is shown in the top panel of Fig. 2.

(ii) When $|\Delta{\hat{P}}_{\alpha}|\lesssim\Delta{\hat{P}}_{\rm crit}$ and $\hat{P}_{\alpha,0}/(2\pi)$ is close to an integer (i.e., when $\sigma_{\alpha}$ is an integer multiple of the orbital frequency, $2\pi/P_{0}$), the mode energy exhibits larger-amplitude oscillations with a mean mode energy $\gg\Delta E_{\alpha}$ (see the middle panel of Fig. 2.)

(iii) When $|\Delta{\hat{P}}_{\alpha}|\gtrsim\Delta{\hat{P}}_{\rm crit}$, the mode amplitude evolves chaotically and can grow to very large values of order the initial orbital energy $|E_{B,0}|$; meanwhile, the orbit experiences significant decay in the semi-major axis. The bottom panel of Fig. 2 provides an example. Note that in this example the binding energy of the planet is of order the initial orbital energy ($GM_{p}^{2}/R_{p}=2.6\;|E_{\rm B,0}|$ for an initial semi-major axis $a_{0}=1$ au). For such chaotic systems, the mode energy will eventually grow to values comparable to the binding energy of the planet, at which point the linear treatment is no longer appropriate and the mode will dissipate energy due to non-linear effects. We will discuss the implications of non-linear dissipation in later sections.

These three different behaviours can be characterized by the maximum mode energy that is attained over many pericentre passages. The examples in Fig. 2 demonstrate that the mode energy remains a small fraction of the initial orbital energy for non-chaotic evolution. However, when the system behaves chaotically, the maximum mode energy, $\rm{max}(E_{\alpha,k})$, can be significant compared to $|E_{\rm B,0}|$.

In Fig. 3 we explore the boundary for chaotic behaviour in the $r_{\rm p}-e$ plane by recording $\text{max}(E_{\alpha,k})$ over $10^{3}$ orbits for different planet models with $M_{p}=1M_{J}$ and $R_{p}=1.6~{}R_{J}$ (the stellar mass is fixed at $M_{*}=1M_{\odot}$). The left and right panels of Fig. 3 assume that $\Omega_{s}$ is the pseudosynchronous rate from equation (29) while the middle panel uses a constant value of $\Omega_{s}=0.2(GM_{p}/R_{p}^{3})^{1/2}$. Small changes in $\Omega_{s}$ do not have a large effect on the chaotic boundary. The right panel of Fig. 3 uses a $\gamma=5/3$ polytropic planet model rather than our standard choice of $\gamma=2$. A comparison between the right panel and others suggests that changes in the structure of the planet can significantly shift the boundary for chaotic behaviour.

Figure 3 demonstrates that the condition $|\Delta\hat{P}_{\alpha}|\gtrsim 1$ generally captures the boundary for chaotic behaviour. This boundary is somewhat fuzzy because the precise value of $\Delta{\hat{P}}_{\rm crit}$ spans a wide range, from $<0.1$ to $>1$, depending on the value of ${\hat{P}}_{\alpha,0}$ (see Fig. 1 of Vick & Lai, 2018). In particular, $\Delta{\hat{P}}_{\rm crit}\ll 1$ when ${\hat{P}}_{\alpha,0}/(2\pi)$ is close to an integer (resonance). Figure 3 confirms that $\Delta{\hat{P}}_{\rm crit}\sim 1$ on average. When $|\Delta\hat{P}_{\alpha}|\gtrsim 1,$ the change in the mode phase is nearly a random number mod $2\uppi$, and each pericentre “kick" to the mode occurs at a random phase. Under this condition, the mode energy grows in a diffusive manner (Mardling, 1995a; Ivanov & Papaloizou, 2004; Wu, 2018). Note that this “diffusion" is approximate; as explained in Vick & Lai (2018) (see their Fig. 3), there exists an upper limit to the mode energy even in the absence of dissipation. For this reason, we prefer to call the mode growth and dissipation in the regime $|\Delta\hat{P}_{\alpha}|\gtrsim 1$ “chaotic tides."

From equation (28), we see that $|\Delta{\hat{P}}_{\alpha}|$ is larger for more eccentric orbits, and smaller $\eta$ (i.e. smaller $r_{\rm p}$ and larger $R_{p}$). A young gas giant that is pumped into a highly eccentric orbit can satisfy the condition for chaotic tidal evolution (see Fig. 3) and rapidly transfer orbital energy to the f-mode. This behaviour is most likely to occur for $e\gtrsim 0.95$. The critical condition for chaotic tides, $|\Delta{\hat{P}}_{\alpha}|\simeq\Delta{\hat{P}}_{\rm crit}\sim 1$, can be written (using equation 28) as

[TABLE]

Note that $T$ depends on $\eta$ and $\sigma_{\alpha}/\Omega_{\rm p}=\bar{\sigma}_{\alpha}\eta^{3/2},$ and is nearly independent of $e$ for large eccentricities. Equation (32) gives the critical eccentricity, for a given $\sigma_{\alpha}$, above which chaotic tides operate as a function of the dimensionless pericentre distance $\eta$ (equation 22).333The mode frequency $\sigma_{\alpha}$ depends on the planet rotation rate $\Omega_{s}$, which may depend on $e$ if the equilibrium rotation rate depends on $e$. For pseudosynchronous rotation, this dependence is very weak when $(1-e)\ll 1$ (for which $\Omega_{s}\simeq 1.17\Omega_{p}$); see equation (29). If we make further simplification for the function $T$ by replacing $K_{22}$ in equation (25) with the approximate fitting formula $K_{2,2}\approx 1.79\times 10^{4}z^{-6}\eta^{3/2}$ (Wu, 2018), equation (32) can be inverted to yield

[TABLE]

or equivalently,

[TABLE]

Equation (34) gives the critical pericentre distance for chaotic tides to operate as a function of eccentricity.

2.4 Conditions for Chaotic Tides When $E_{\alpha,0}>0$

In Section 2.3, we have discussed the condition for the onset of chaotic mode growth for an “eccentric orbit + oscillatory modes" system where there is no initial energy in the mode ($E_{\alpha,0}=0$). When dynamical tides are the only influence on the planet’s orbital evolution, and when the total energy of the system is conserved (i.e. the modes do not dissipate), ${\hat{P}}_{\alpha,0}$ and $|\Delta{\hat{P}}_{\alpha}|$ completely determine the behaviour of the system. For instance, a system with $\Delta{\hat{P}}_{\alpha}\lesssim\Delta{\hat{P}}_{\rm crit}\sim 1$ will never undergo chaotic mode growth, and one with $\Delta{\hat{P}}_{\alpha}\gtrsim\Delta{\hat{P}}_{\rm crit}$ will always be chaotic. However, when non-tidal effects alter the planet’s orbit (e.g. LK eccentricity oscillations driven by an external companion), a system can transition in and out of the chaotic regime as $r_{\rm p}$ is driven below and above $r_{\rm p,crit}$ [see equations (32)-(34) and Fig. 3]. Therefore, it is useful to know whether a system is “currently chaotic" when there is “pre-existing" energy in the mode ($E_{\alpha,0}>0$). Here, we generalize the criterion for the onset of chaotic tides to account for the “pre-existing" mode energy (see Mardling, 1995a; Wu, 2018).

As discussed in Section 2.3, the mode amplitude evolves chaotically when the mode phase at pericentre is nearly random from one passage to the next. This occurs when the phase-shift from pericentre energy transfer is of order unity. For a mode with no “pre-existing" energy, the pericentre phase shift is solely due to $\Delta E_{\alpha}$. When the mode already has energy $E_{\alpha,0}$ prior to the first passage, the pericentre phase shift can reach a maximum value of

[TABLE]

where $\text{max}(\Delta E_{\alpha,1})$ is the maximum possible energy transfer in the first pericentre passage given the initial mode energy $E_{\alpha,0}$, and is given by

[TABLE]

The general condition for chaotic behaviour is that the maximum pericentre phase shift is sufficiently large that that the mode phase is nearly random at pericentre, i.e.

[TABLE]

Note that $|\Delta{\hat{P}}_{\alpha}(E_{\alpha,0})|$ can be many times larger than $|\Delta{\hat{P}}_{\alpha}|=|\Delta{\hat{P}}_{\alpha}(0)|$ if $E_{\alpha,0}\gg\Delta E_{\alpha}$; in this case equation (37) reduces to equation (18) of Wu (2018).444Wu (2018) adopted $\Delta\hat{P}_{\rm crit}=2/3$. In reality, $\Delta\hat{P}_{\rm crit}$ depends strongly on $\hat{P}_{0}$; see Fig. 4.

Figure 4 shows how the dynamical behaviour of the system depends on ${\hat{P}}_{\alpha,0}$ and $|\Delta{\hat{P}}_{\alpha}(E_{\alpha,0})|$ for $E_{\alpha,0}=0,\Delta E_{\alpha}$, and $0.001GM_{p}^{2}/R_{p}$. The left panel is similar to Fig. 1 from Vick & Lai (2018). The middle panel corresponds to a system where the planet has an initial energy $E_{\alpha,0}=\Delta E_{\alpha}$. This is likely the case the time LK oscillations drive the planet’s orbital eccentricity to the regime for chaotic mode growth. In the right panel, $E_{\alpha,0}\gg\Delta E_{\alpha}$, which is relevant for determining whether a system undergoing chaotic mode growth will continue to behave chaotically (see Section 3.1). As in Fig. 3, the orange regions indicate chaotic tides. While the exact value of $\Delta{\hat{P}}_{\rm crit}$ depends on whether the mode frequency is near a resonance with the orbital frequency, in general a system with $|\Delta{\hat{P}}_{\alpha}(E_{\alpha,0})|\gtrsim 1$ is more likely than not to experience chaotic tides. This criterion works best for $\tilde{E}_{\alpha,0}\equiv E_{\alpha,0}/|E_{B,0}|\ll 1$; for larger $\tilde{E}_{\alpha,0}$, $\Delta\hat{P}_{\rm crit}$ gradually moves to larger values, as suggested by the right panel of Fig. 4.

We can also test the generalized criterion for chaotic tides over the range of orbital parameters relevant to migrating gas giants. Figure 5 shows the similar results as Fig. 3, but for multiple values of $E_{\alpha,0}$. We see that the condition $|\Delta{\hat{P}}_{\alpha}(E_{\alpha,0})|=1$ (the solid black line) matches well with the boundary for chaotic tidal behaviour. The critical pericentre distance for chaotic tides (for a given eccentricity $e$), $r_{\rm p,crit}$, increases with $E_{\alpha,0}$ while $e_{\rm crit}$ (for given $r_{\rm p}$) decreases.

In summary, chaotic tides are easier to achieve when planetary modes are already oscillating. The larger the mode energy, the more relaxed the conditions for chaotic tidal behaviour.

2.5 Conditions for Continued Chaotic Behaviour

As mentioned in Section 2.3, we expect the f-mode of a planet experiencing chaotic tides to dissipate a significant amount of energy when the mode amplitudes become non-linear. As a mode drains energy from the orbit, the orbit becomes more tightly bound and a larger amount of energy transfer is required to significantly alter the mode phase at pericentre. Eventually, a system undergoing chaotic tides can dissipate enough energy that the planetary mode enters a quasi-steady state. For dissipative systems, it is useful to know when the orbit is too tightly bound for the mode to continue chaotic behaviour.

To determine whether the system is chaotic after the $k$-th pericentre passage, we ask whether the current maximum pericentre phase shift, $|\Delta{\hat{P}}_{\alpha,k}(E_{\alpha,k-1})|$ is sufficiently large ($\gtrsim 1$). A natural way to think of this condition is to rescale $\Delta{\hat{P}}_{\alpha}({\tilde{E}}_{\alpha,0})$ to the current orbital energy and period, $E_{B,k}$ and $P_{k}$. Similar to equation (35), we have

[TABLE]

where

[TABLE]

When evaluating equations (38) and (39), one must use the “current" parameters of the system, e.g. the orbital period $P_{k-1}$ and mode frequency $\sigma_{\alpha,k-1}$ just before the $k$-th pericentre passage. Note that the mode frequency may change as the system evolves because it depends on the planet’s spin.

Equations (38) and (39) suggest that a system can escape the influence of chaotic tides in several ways. If a non-tidal effect is dominating the orbital evolution and $r_{\rm p}$ increases or $e$ decreases, $\Delta E_{\alpha}$ and $\text{max}(\Delta E_{\alpha,k})$ can become too small for chaotic tides to continue operating; if the orbital evolution is dominated by chaotic tides, $r_{\rm p}$ and therefore $\Delta E_{\alpha}$ are nearly constant, but the orbit can become too tightly bound (i.e. $|E_{B,k-1}|$ is too large and $P_{k-1}$ is too small), for the system to satisfy $|\Delta{\hat{P}}_{\alpha,k}(E_{\alpha,k-1})|\gtrsim 1$. Chaotic behaviour can also be suppressed if the mode suddenly dissipates a large amount of energy such that $E_{\alpha,k}\ll E_{\alpha,k-1}$.

3 Lidov-Kozai migration with chaotic tides

We now discuss our implementation of dynamical tides in giant planet migration driven by the Lidov-Kozai (LK) mechanism. In this scenario, an external stellar or planetary companion induces quasi-periodic oscillations in the eccentricity of the giant planet’s orbit; as the eccentricity attains a large value $(\sim 1)$, tidal dissipation in the planet at pericentre reduces the orbital energy, leading to circularization and decay of the planet’s orbit. As noted in Section 1, previous works on LK migration adopted the (parametrized) weak friction theory of static tides (ST). The grey curves of Fig. 6 give an example of this “standard" LK migration scenario.

The theory of LK migration with static tides involves two important timescales. The first is the timescale for quadrupole LK eccentricity oscillations, $t_{\rm LK}$, given by

[TABLE]

where $a_{0}$ is the initial semi-major axis of the planet’s orbit, $M_{b}$ is the mass of the (stellar) companion, and

[TABLE]

with $a_{\rm b}$ and $e_{\rm b}$ the semi-major axis and eccentricity of the companion’s orbit. The time that the planet spends near the maximum eccentricity, $e_{\rm max}$, is of order (e.g. Anderson et al., 2016)

[TABLE]

The second important timescale is that of orbital decay due to static tides given by (Alexander, 1973; Hut, 1981).

[TABLE]

where $k_{2p}$ is the tidal Love number of the planet, $\Delta t_{L}$ is the lag time, $f_{1}(e)=1+31e^{2}/2+255e^{4}/8+185e^{6}/16+25e^{8}/64$, and $f_{2}(e)$ and $f_{5}(e)$ are given in equations (30) and (31). (We have assumed that the planet has a pseudosynchronous spin rate.) Tidal dissipation is most efficient near the maximum eccentricity. Since the planet only spends a fraction ($\sim\sqrt{1-e_{\rm max}^{2}}$) of the time near $e_{\rm max}$, the effective orbital decay rate during LK migration is

[TABLE]

where $r_{\rm p,min}=a_{0}(1-e_{\rm max})$ is the minimum pericentre distance (in the second of the above equalities, we have used $e_{\rm max}=0.96$). Successful migration within a few Gyrs requires $r_{\rm p,min}\lesssim 0.025$ au, corresponding to a final (circularized) planet semi-major axis $a_{\rm F}\lesssim 0.05$ au. A planet would need to be more dissipative than Jupiter (larger $\Delta t_{\rm L}$)555For Jupiter, $\Delta t_{L}\simeq 0.1$ s, corresponding to $k_{2p}/Q_{p}\simeq 10^{-5}$ at a tidal forcing period of $0.5$ hrs. to become a HJ with a larger $a_{\rm F}$.

As we shall see, chaotic tides can change the standard LK migration scenario in several important aspects (see Fig. 6). Although the planet spends only a small fraction of time in the high-eccentricity phase of the LK oscillations, this time covers enough orbits to allow chaotic mode growth. Consider the example depicted in Fig. 6: the critical eccentricity (for $a_{0}=1.5$ au) for chaotic mode growth to operate is $e_{\rm crit}\simeq 0.985$ (see Fig. 3 and equation 32). So the time the planet spends around $e\gtrsim e_{\rm crit}$ is $\sim\sqrt{1-e_{\rm crit}^{2}}t_{\rm LK}\approx 0.17t_{\rm LK}\approx 0.23$ Myr, which is more than $1.2\times 10^{5}$ orbital periods. The planet’s oscillation mode can grow to a large amplitude within this time frame, stealing orbital energy in the process and allowing the orbit to decay within one or a few LK cycles.

3.1 Description of the Model

We now describe our method for coupling the evolution of the oscillation mode (the $l=m=2$ f-mode) of a gas giant with the evolution of its orbit driven by an external binary companion (LK oscillations). The planet starts out with eccentricity $e_{0}\simeq 0$ and zero mode amplitude ($E_{\alpha,0}=0$). The other planet and stellar properties as well as the initial values of $e_{\rm b}$, $a_{0}$, $a_{\rm b}$, $\theta_{\rm lb,0}$ (the mutual inclination of the orbits) and $\Omega_{0}$ (the longitude of the ascending node of the planet’s orbit) are taken as input parameters. Apart from dynamical tides, we evolve the orbital angular momentum vector $\bm{L}$ and eccentricity vector $\bm{e}$ of the planet, and the spin angular momentum vector, $\bm{S_{*}}$ of the host star in the same way as in Anderson et al. (2016) (hereafter ASL16), i.e.

[TABLE]

Here the LK terms are contributions (to the octupole order) from the binary companion that give rise to LK oscillations [equations A1-A2 of ASL16; from Liu et al. (2015)]; the SL terms arise from the spin-orbit coupling between the host star spin, $\bm{S_{*}}$, and the planet’s orbital angular momentum (equations 60-61 of ASL16); and the short-range force (SRF) terms account for periastron precession of the planet’s orbit due to general relativity (GR) and tidal/rotational distortions of the planet (sections A3 and A4 of ASL16). The stellar spin evolution includes the spin-down torque due to magnetic braking (MS), i.e.

[TABLE]

where $\alpha_{\rm MB}$ is taken to be $1.5\times 10^{-14}$ yr to model a solar type star, $I_{*}$ is the moment of inertia of the star, and $\Omega_{*}$ is its rotation rate. We assume an initial spin period of 2.3 days. We also evolve the angular momentum and eccentricity vectors of the outer orbit ($\bm{L}_{b}$ and $\bm{e}_{b}$) according to the octupole LK effect (equations 19 and 20 of Liu et al., 2015); these changes in the outer orbit are small because $M_{p}\ll M_{b}$. Throughout evolution, we take the planet’s rotation rate to be the pseudosynchronous rotation rate given by equation (29). Equations (45)-(46) include contributions from the dissipation of static tides (ST); these are parameterized by the tidal lag time $\Delta t_{\rm L}$, and are distinct from the effects of dynamical tides.

In the presence of dynamical tides, we must supplement equations (45)-(47) with the evolution of the mode amplitude. This is based on the iterative map discussed in Section 2.1. Here, we provide a summary of the implementation of this procedure. At a certain time in the orbital evolution, just before the $k$-th pericentre passage, the orbital parameters are $a_{k-1}$ and $e_{k-1}$, the planet spin rate is $\Omega_{s,k-1}$, and the mode amplitude is $c_{\alpha,k-1}$. To advance to the next orbit, we carry out the following steps:

(i) Calculate $\Delta E_{\alpha}$ from equations (21), (25), and (27) using mode parameters from Table 1 and the current orbital/planet parameters (i.e. $a_{k-1}$, $e_{k-1}$, $\Omega_{s,k-1}$, etc.) and obtain $\Delta c_{\alpha}\propto\sqrt{\Delta E_{\alpha}}$. In practice, we normalize mode energy by $|E_{B,0}|$ (the initial orbital energy), and thus $\Delta\tilde{c}_{\alpha}=\sqrt{\Delta{\tilde{E}}_{\alpha}}=\sqrt{\Delta E_{\alpha}/|E_{B,0}|}$.

(ii) Calculate the energy transfer in the $k$-th passage using $\Delta E_{\alpha,k}=|E_{B,0}|(|\tilde{c}_{\alpha,k-1}+\Delta\tilde{c}_{\alpha}|^{2}-|\tilde{c}_{\alpha,k-1}|)$ (see equation 17) and obtain the new orbital energy $E_{B,k}=E_{B,k-1}-\Delta E_{\alpha,k}$.

(iii) Obtain the new orbital semi-major axis and eccentricity (after the $k$-th pericentre passage) according to

[TABLE]

where we have assumed that the orbital angular momentum is conserved during the passage.

(iv) Evolve the (complex) mode amplitude over one orbit to obtain its value just before the $(k+1)$-th passage using $c_{\alpha,k}=(c_{\alpha,k-1}+\Delta c_{\alpha})\text{e}^{-\text{i}\sigma_{\alpha,k}P_{k}}$ (see equation 14), where $P_{k}$ is the orbital period corresponding to $a_{k}$ and $\sigma_{\alpha,k}$ is the mode frequency after the $k$-th passage (assuming a pseudosynchronous rotation rate).

Left unchecked, the mode energy in this model can grow to un-physically large values (as in the bottom panel of Fig. 2). We therefore assume that, when the mode energy reaches a significant fraction of the planet’s binding energy, non-linear effects dissipate nearly all of the mode energy within an orbital period. We parametrize the maximum energy that the mode can reach before non-linear effects become important, $E_{\rm max}$, and the residual energy in the mode after an episode of nonlinear mode dissipation and tidal heating, $E_{\rm resid}$. Thus, when $|\tilde{c}_{\alpha,k}|^{2}>E_{\rm max}/|E_{B,0}|$, the mode amplitude is immediately changed to $|\tilde{c}_{\alpha,k}|^{2}=E_{\rm resid}/|E_{B,0}|$. This method for handling dissipation was used in Wu (2018). The parameter $E_{\rm max}$ does not change the qualitative features of evolution due to chaotic tides (see Fig. 7), but can alter the details. The quantity $E_{\rm resid}$ is more important because chaotic behaviour is easier to excite when the f-mode already has some energy, as discussed in Section 2.4, but still does not change the overall behaviour of the model (see Fig. 8). Unless stated otherwise, we use $E_{\rm max}=0.1(GM_{p}^{2}/R_{p})$ and $E_{\rm resid}=0.01E_{\rm max}=10^{-3}(GM_{p}^{2}/R_{p})$ in the results presented in this paper.

Since dynamical tides must be implemented on an orbit by orbit basis, it is not practical to evolve the whole system over many Myrs. Because dynamical tides are effective only during the high-e phase of the LK cycle, when chaotic mode growth occurs, we only need to evolve the modes when the eccentricity is sufficiently high and the pericentre distance is sufficiently small. As discussed in Section 2.5, chaotic tides operate when $|\Delta\hat{P}_{\alpha,k}(E_{\alpha,k-1})|\gtrsim\Delta\hat{P}_{\rm crit}$ (see equation 38). In our implementation, we turn on dynamical tides only when $|\Delta\hat{P}_{\alpha,k}(E_{\alpha,k-1})|\gtrsim 0.01$. As we see from Fig. 4, this is a conservative choice, allowing us to capture the onset of chaotic behaviour for all reasonable parameters regardless of resonances and the value of “pre-existing" mode energy.

After the orbit has decayed through chaotic tides, the system will eventually satisfy $|\Delta\hat{P}_{\alpha,k}(E_{\alpha,k-1})|<\Delta\hat{P}_{\rm crit}\sim 1$ (see the right panel of Fig. 4) while the orbital eccentricity is still large. This typically occurs just after the mode energy is dissipated non-linearly from $E_{\rm max}$ to $E_{\rm resid}$. Thus, when the system satisfies

[TABLE]

(where we have assumed $E_{\alpha,k-1}\sim E_{\rm resid}\gg\Delta E_{\alpha}$) the f-mode will no longer behave chaotically and influence the orbital evolution. Instead, static tides drive the gradual decay and circularization of the planet’s orbit. To save computation time, when equation (51) is satisfied while the orbital eccentricity is large, we take note. If equation (51) is still satisfied in 30,000 orbits, we stop evolving the dynamical tides; if equation (51) is not satisfied after 30,000 orbits, we check again in another 30,000 orbits. The multiple checks are to account for the fact that $\Delta{\hat{P}}_{\rm crit}$ varies significantly depending on the orbital period and can be less than 1.

3.2 Examples of the Model

LK migration with chaotic tides proceeds in three steps (see the blue curves in Fig. 6). First, the external companion drives the planet’s orbital eccentricity beyond $e_{\rm crit}$ (see equation 32), where chaotic tides set in and begin to take over orbital evolution, freezing $r_{\rm p}$ due to angular momentum conservation. The blue lines in Figs. 7 and 8 show a zoom-in view. We see that $r_{\rm p}$ decreases due to the LK effect until settling to an approximately constant value around 2.1 Myr while $e$ and $a$ start to wander chaotically. In the second stage, chaotic tides dominate orbital evolution and quickly shrink the orbit until it is too bound for such tides to continue operating. In Figs. 7 and 8, it is clear that the planet undergoes multiple episodes of non-linear dissipation (where the mode energy rapidly dissipates to a small residual value of $E_{\rm resid}$ after climbing to maximum energy $E_{\rm max}$) while subject to chaotic tides. The result is a very eccentric WJ. At this new semi-major axis, the planet is decoupled from the influence of the stellar companion, and no longer experiences significant eccentricity oscillations. Finally, weak tidal friction circularizes and shrinks the orbit to form a HJ on a Gyr timescale (depending on the dissipation rate of static tides in the planet).

In some cases, a system will undergo multiple high-eccentricity phases before chaotic tides reach the threshold $E_{\rm max}$ and the f-mode dissipates energy. An example is shown in the blue lines of Fig. 9, where the system has the same parameters and initial conditions as in Fig. 6 but with a slightly larger companion semi-major axis ($a_{\rm b}=225$ au). We see that the planet’s eccentricity undergoes five LK oscillations before the f-mode energy climbs to $E_{\rm max}$. Before an episode of tidal dissipation, any energy transferred from the orbit to the mode can easily pass back to the orbit. For example, in the third high-eccentricity phase of Fig. 9 (shown in an inset of the bottom right panel) the f-mode loses a small amount of energy and the semi-major axis increases slightly. Once the mode energy reaches $E_{\rm max}$ and is dissipated, the orbit shrinks irreversibly and continues to rapidly decay until the orbit is too bound for chaotic tides to operate.

We can gain some insight into when a system requires multiple high-eccentricity phases to reach $E_{\rm max}$ by studying the rate of chaotic orbital decay predicted by the iterative map from Section 2.1. While $|\Delta\hat{P}_{\alpha}|>1$ (for $E_{\alpha,0}=0$), the mode amplitude grows in an approximately diffusive manner [Fig. 5 of Vick & Lai (2018); see also Mardling (1995a); Ivanov & Papaloizou (2004); Wu (2018)], with the averaged mode energy

[TABLE]

For a sufficiently eccentric orbit, $\Delta E_{\alpha}$ is roughly constant over many pericentre passages. We can define a timescale for orbital decay from the energy transfer and dissipation,

[TABLE]

If we use the largest value of $e$ attained by the system in Fig. 9, we find that $t_{\rm decay}=0.31$ Myr and the timescale for the system to remain at high eccentricity is $\Delta t(e=0.983)=0.13$ Myr (see equation 42). Because $t_{\rm decay}<\Delta t(e=0.983)$, the system must undergo multiple LK cycles for the mode energy to climb to $E_{\rm max}$. In contrast, the system in Fig. 6 reaches $e=0.989$, for which $t_{\rm decay}=3.9\times 10^{-2}$ Myr and $\Delta t(e=0.989)=0.15$ Myr. For this system, where the timescale for chaotic tides to alter the orbit is very short, the planet’s orbit decays within one high-eccentricity phase.

4 Key Features of Chaotic Tidal Migration: Analytical Understanding

LK migration with chaotic tides has a few characteristic features. First, it can prevent some gas giants from undergoing tidal disruption. Second, it naturally produces very eccentric WJs. Lastly, this population of WJs circularizes relatively quickly, due to weak tidal friction because the planets are “detached" from the companions’ perturbations. In the following, we discuss the necessary conditions for LK chaotic tidal migration, hallmarks of this process, and predictions for the stellar obliquities of HJs formed via this mechanism.

4.1 Condition for Chaotic Tidal Migration

As discussed in Section 2.3, a planet must have $e\gtrsim e_{\rm crit}$ to initiate chaotic tides, where $e_{\rm crit}$ as a function of $a_{0}$ (the initial semi-major axis of the planet) is implicitly given by equation (32). To reach such a high eccentricity, a system must have sufficiently large initial mutual inclination, $\theta_{\rm lb,0}$. In the idealized case of quadrupole LK oscillations, the maximum eccentricity is $e_{\rm max}=\left[1-(5/3)\cos^{2}\theta_{\rm lb,0}\right]^{1/2}.$ Accordingly, $e_{\rm max}\sim 1$ can be achieved only for $\theta_{\rm lb,0}\approx 90^{\circ}$. Including octupole terms complicates the relationship between the initial mutual inclination and the maximum attainable eccentricity and broadens the inclination window for extreme eccentricity excitation. A realistic treatment of the LK effect must also account for the short-range forces (SRFs) that cause the planet’s pericentre to precess. These limit the maximum attainable eccentricity in LK cycles to a value $e_{\rm lim}$, where $e_{\rm lim}$ is the maximum eccentricity achieved when $\theta_{\rm lb,0}=90^{\circ}$ in the quadrupole limit. In the presence of the octupole potential, the eccentricity still cannot exceed $e_{\rm lim}$, but this eccentricity can be realized for a wider range of initial inclinations (Liu et al., 2015). The relevant SRF’s for this problem are GR and the effect of static tides raised on the planet. When $(1-e_{\rm lim})\ll 1$, the limiting eccentricity is given by

[TABLE]

where $\epsilon_{\rm GR}$ and $\epsilon_{\rm Tide}$ measure the strengths of precession due to GR and the planet’s tidal bulge relative to the quadrupole LK effect, respectively:

[TABLE]

For our “standard" system (shown in Fig. 6) with $a_{0}=1.5$ au, $a_{\rm b}=200$ au, and $e_{\rm b}=0$, we find $e_{\rm lim}\approx 0.998$ and $e_{\rm crit}\approx 0.985$, so the planet can reach large enough eccentricities to experience chaotic tides, provided the initial mutual inclination angle, $\theta_{\rm lb,0}$, is sufficiently large.

In general, the necessary condition for chaotic tidal migration is

[TABLE]

or equivalently $r_{\rm p,lim}\equiv a_{0}(1-e_{\rm lim})\leq r_{\rm p,crit}$, where $e_{\rm crit}$ and $r_{\rm p,crit}$ are given by equation (32) (or more approximately by equation 34). Because $e_{\rm lim}$, given by equation (54), depends on $a_{\rm b,eff}$, the condition on $e_{\rm lim}$ translates to an upper limit on $a_{\rm b,eff}$ for chaotic migration to occur:

[TABLE]

The upper limit in equation (58) is shown as a blue line in Fig. 10.

We can also derive the maximum $r_{\rm p}$ that allows for standard LK migration with static tides within a stellar lifetime. Using equation (44) and requiring that the planet migrate within $t_{\rm mig}$ (i.e. $t_{\rm ST,LK}\lesssim t_{\rm mig}$), we find

[TABLE]

[$r_{\rm p,ST}$ is called $r_{\rm p,mig}$ in Muñoz et al. (2016)]. Substituting $r_{\rm p,ST}$ for $r_{\rm p,crit}$ in equation (58) yields an approximate maximum $a_{\rm b,eff}$ for standard LK migration (without chaotic tides). The result is shown as the dashed line in Fig. 10. For a $1M_{J}$, $1.6R_{J}$ gas giant, the conditions for chaotic tidal migration are typically more generous than for standard LK migration. This suggests that LK migration with chaotic tides yields more HJs than LK migration with static tides.

4.2 Evading Tidal Disruption

The planet can be tidally disrupted if its pericentre distance, $r_{\rm p}$, is less than the tidal disruption limit from equation (24),

[TABLE]

Therefore, when $r_{\rm p,lim}<r_{\rm p,dis}$, the companion can induce tidal disruption of the planet if $\theta_{\rm lb,0}$ is sufficiently large. Substituting $r_{\rm p,dis}$ for $r_{\rm p,crit}$ in equation (58) yields the maximum $a_{\rm b,eff}$ for tidal disruption to be possible. This is shown as the red line in Fig. 10.

In the standard LK migration scenario (with static tides), a large fraction of migrated giant planets are tidally disrupted because $r_{\rm p,dis}$ can be quite close to $r_{\rm p,ST}$; moreover, Saturn-mass planets can rarely migrate successfully via LK oscillations and static tides because $r_{\rm p,dis}\gtrsim r_{\rm p,ST}$ (ASL16, Muñoz et al., 2016). In contrast, with chaotic tides, there is always a region of parameter space where chaotic tidal migration is possible without disruption ($r_{\rm p,crit}>r_{\rm p,dis}$; see Fig. 3) for any reasonable values of $M_{p}$ and $R_{p}$ for giant planets.

In some cases, chaotic tides can shepherd to safety the planets that are otherwise destined for tidal disruption by acting to rapidly decrease $e$ and increase $r_{\rm p}$. This possibility was suggested in Wu (2018). Figure 11 shows one such example. Before the planet is pushed to extreme eccentricity and tidal disruption by the octupole potential, it is frozen into a high-eccentricity orbit with a small semi-major axis and becomes decoupled from the companion. Without further interference from the LK effect, the eccentric WJ circularizes in orbit over a Gyr timescale to eventually become a HJ.

We can understand the condition for chaotic tides to save a planet from disruption using a simple time-scale argument. Tidal disruption occurs when $r_{\rm p}\leq r_{\rm p,dis}$. At high eccentricity, the timescale for a planet to remain below a given $r_{\rm p}$ is $\Delta t(<r_{\rm p})\sim t_{\rm LK}\sqrt{1-e^{2}}\sim\sqrt{2}t_{\rm LK}(r_{\rm p}/a_{0})^{1/2}$ (see equation 42), while the timescale for chaotic tides to decay the orbit is $t_{\rm decay}$, given by equation (53). The planet reaches a minimum $r_{\rm p}$ when these timescales are roughly equal. For planets that are “just saved" from disruption, this minimum is nearly $r_{\rm p,dis}$. A system that can be saved from disruption must satisfy $t_{\rm decay}\lesssim\Delta t(<r_{\rm p,dis})$, i.e.

[TABLE]

where $\Delta E_{\alpha}$ is evaluated at $r_{\rm p}=r_{\rm p,dis}$ (see equation 21). Recall that $\Delta E_{\alpha}$ depends on $T(\eta,\sigma_{\alpha}/\Omega_{p},e)$. When $(1-e)\ll 1$, $T$ and the pseudosynchronous spin rate, $\Omega_{ps}$, have negligible dependence on $e$. We can then evaluate $T(\eta)$ by using equation (25) with the approximate $K_{22}$ from equation (27) and taking $\bar{\sigma}_{\alpha}\simeq\bar{\omega}_{\alpha}+1.17/\eta^{3/2}$ for the mode frequency in the inertial frame. At $r_{\rm p,dis}$, $T(\eta=2.7)=2.5\times 10^{-2}$.

By rearranging equation (61), we can find the minimum $a_{\rm b,eff}$ needed for chaotic tidal migration:

[TABLE]

with $r_{\rm p,dis}=2.7r_{\rm tide}$. The limit on $a_{\rm b,eff}$ is shown as a black line in Fig. 10.

4.3 Eccentric Warm Jupiter Formation

A giant planet that undergoes chaotic orbital decay becomes an eccentric WJ after reaching a semi-major axis where both chaotic tides and the LK effect have little influence on the continued evolution of the planet’s orbit. For example, in Figs. 6 - 8, all of the integrations produce a planet with a semi-major axis between $0.35$ and $0.65$ au and an eccentricity between 0.93 and 0.97. We can understand why this occurs by re-examining the criterion for chaotic tides to stop affecting orbital evolution (see equation 51). Typically, chaotic tides end just after a dissipation episode, where $E_{\alpha,k-1}\sim E_{\rm resid}$. Setting $|\Delta\hat{P}_{\alpha,k}(E_{\rm resid})|\sim 1$ with $P_{k}\propto a_{k}^{3/2}$ and $|E_{B,k}|\propto a_{k}^{-1}$, we find

[TABLE]

where $\Delta E_{\alpha}$ and $\sigma_{\alpha}$ are evaluated at $r_{\rm p}=r_{\rm p,WJ}$, the pericentre distance of the planet when chaotic tides stop operating. The value of $r_{\rm p,WJ}$ varies erratically for small changes of initial conditions, but must fall between $r_{\rm p,dis}$ ($\eta=2.7$) and $r_{\rm p,crit}$ (very roughly $\eta\sim 4$) for a system with $a_{0}<5$ au. Equation (63) gives the relationship between $a$ and $r_{\rm p}$ for WJs formed by chaotic tides. Parametrizing the residual mode energy by $E_{\rm resid}\equiv fGM_{p}^{2}/R_{p}$, and using equation (21), we can write $a_{\rm WJ}$ in terms of the planet properties and $\eta$:

[TABLE]

Equation (64) is plotted in Fig. 12 for two different planet models. For values of $\eta$ between $2.7$ and $3.6$, equation (64) yields semi-major axes in the range of $0.3-1.0$ au (for $R_{p}=1.6R_{J}$ and $M_{*}/M_{p}=1000$). We can also find the eccentricity of the WJs via

[TABLE]

For $\eta$ between $2.7$ and $3.6$, $(1-e_{\rm WJ})$ is between 0.027 and 0.060 (again for $R_{p}=1.6R_{J}$ and $M_{*}/M_{p}=1000$).

Figures 6, 9, and 11 show that after the WJ forms through chaotic tides, the eccentricity freezes at a high value and the planet’s orbit decouples from the companion. In general, LK eccentricity oscillations freeze when $\dot{\omega}_{\rm SRF}t_{\rm LK}\sqrt{1-e^{2}}\gtrsim 1,$ where $\dot{\omega}_{\rm SRF}$ is the rate of precession due to SRFs. This is equivalent to $\epsilon_{\rm GR}/\sqrt{1-e^{2}}\gtrsim 1$ or $\epsilon_{\rm Tide}/(1-e^{2})^{9/2}\gg 1$. We find that tidal effects play a more important role than GR in the “freezing" of eccentricity oscillations of WJs.

4.4 Fast Long-Term Orbital Decay and Hot Jupiter Formation

When the planet’s orbit is no longer strongly influenced by chaotic f-mode evolution and the LK effect, the orbit decays and eventually circularizes due to tidal friction from static tides raised on the planet. The final result is a HJ. This is the same mechanism that has been used to explain HJ formation in previous studies of migration via the LK effect (ASL16; Petrovich, 2015a). Here, the process occurs more quickly because the orbit is frozen into a small pericentre value rather than oscillating between high and low eccentricities. Figures 6 and 9 compare the long-term orbital evolution predicted by our model (in blue) with the standard calculation that does not include chaotic tides (in grey). In both cases, the orbit circularizes much more quickly when chaotic tides are included.

4.5 Spin-Orbit Misalignment and Final Planet-Binary Inclination

Recent studies (ASL16 Storch et al., 2014, 2017) have shown that a giant planet undergoing LK migration can induce complex dynamics in the spin axis of the oblate host star before the spin-orbit misalignment angle (stellar obliquity), $\theta_{\rm sl}$, becomes frozen. This freezing occurs when the precession rate, $\Omega_{\rm p*}$, of the stellar spin vector ($\bm{S_{*}}$) around $\hat{\bm{L}}$ (the unit orbital angular momentum vector of the planet) becomes much faster than the precession rate ($\Omega_{L}$) of $\hat{\bm{L}}$ around the binary angular momentum axis $\bm{\hat{L}_{b}}$. The ratio of precession rates is

[TABLE]

where

[TABLE]

and

[TABLE]

with $k_{*}$ as the stellar moment of inertia constant, $k_{q*}$ as the stellar rotational distortion coefficient, and $\Omega_{*}$ as the stellar spin rate (see ASL16). We take $k_{*}$ = 0.1 and $k_{q*}=0.05$. We call $\mathcal{A}$ the adiabaticity parameter. The spin-orbit angle freezes when $\mathcal{A}\gg 1$.

Figure 13 shows an example of the evolution of the spin-orbit misalignment angle $\theta_{\rm sl}$ during chaotic tidal migration. The stellar spin rate evolves according to the Skumanich law given in equation (48). When chaotic tides stop operating, around 4.5 Myr, $\mathcal{A}\sim 1$, and $\theta_{\rm sl}$ still oscillates with an amplitude of $\sim 16^{\circ}$. Because $\mathcal{A}$ depends on the stellar spin rate and $\theta_{\rm lb}$, the condition for the adiabatic “spin-orbit freeze" ($\mathcal{A}\gg 1$) is generally different from the condition for “LK oscillations freeze" [$\epsilon_{\rm Tide}/(1-e^{2})^{9/2}\gg 1$; see section 4.3]. Many of the planets that migrate via LK chaotic tides will not have a fixed $\theta_{\rm sl}$ when they become eccentric WJs.

The mutual inclination $\theta_{\rm lb}$ between the orbits of the planet and the stellar companion is also of interest. We see from Fig. 13 that $\theta_{\rm lb}$ can oscillate even at the end of the evolution (when the planet has become a HJ). This variation of $\theta_{\rm lb}$ arises from the coupling between $\bm{L}$ and $\bm{S_{*}}$ (with $|\bm{S_{*}}|$ non-negligible compared to $|\bm{L}|$). We define

[TABLE]

We find that the angle between $\bm{J}$ and $\bm{L_{b}}$, $\theta_{\rm jb}$, is constant when the planet’s orbit is decoupled from the influence of the stellar companion. To understand this, consider the time evolution of the planet’s orbital angular momentum axis $\hat{\bm{L}}$ and the stellar spin axis $\hat{\bm{S}}_{*}$. Both are rapidly precessing around $\bm{J}$. An external torque acting on $\bm{S}_{*}$ (from magnetic braking) or $\bm{L}$ (from the binary companion) also acts on $\bm{J}$. But after averaging over the fast timescale $2\pi/|\Omega_{p*}|$, the net effect is that $\bm{J}$ precesses around $\bm{L}_{B}$. As a result, $\theta_{\rm jb}$ is constant and $\theta_{\rm lb}$ oscillates around $\theta_{\rm jb}$ with amplitude $\theta_{\rm lj}$, where

[TABLE]

The top panel of Fig. 13 illustrates this behaviour. As the star gradually spins down due to magnetic braking, $\theta_{\rm lb}$ asymptotes to $\theta_{\rm jb}$.

5 Population Synthesis

In this section we carry out a population synthesis study of giant planets undergoing LK migration with chaotic tides. Our goal is to determine the production efficiencies of eccentric WJs and HJs, as well as their general properties as predicted by this scenario.

5.1 Setup and Method

All of our population synthesis calculations use the same stellar and companion masses $M_{*}=M_{b}=M_{\odot}$. The initial rotational period of the host star is $2.3$ days, and the star spins down according to the Skumanich law (equation 48). The planet model is a $\gamma=2$ polytrope. We consider two planet radii, $R_{p}=R_{J}$ and $R_{p}=1.6R_{J}$, and two planet masses $M_{p}=1M_{J}$ and $M_{p}=0.3M_{J}$.

While some of our earlier examples have assumed that $e_{\rm b}=0$ (so that the octupole terms vanish), in the population synthesis we consider a broad range of $e_{b}$ and incorporate the octupole LK effect. The importance of the octupole effect (relative to the quadrupole effect) is encoded in the dimensionless parameter

[TABLE]

We have explored a variety of initial conditions, uniformly spanning $a_{0}=[1,5]~{}\text{au}$, $a_{\rm b}=[10^{2},10^{3}]$ (au) (uniformly sampled in $\log{a_{\rm b}}$), $\cos(\theta_{\rm lb,0})=(-0.77,0.77)$, $e_{\rm b}=[0,0.8]$ and $\Omega_{0}=[0,2\uppi]$ (where $\Omega_{0}$ is the initial longitude of the ascending node of the planet’s orbit). Our choice of initial mutual inclinations, $\theta_{\rm lb,0}$, is limited to the range where the quadrupole LK effect can operate, i.e. $\cos^{2}\theta_{\rm lb,0}<3/5$. Systems that do not obey the stability criterion (Mardling & Aarseth, 2001),

[TABLE]

are discarded. Note that the range of initial conditions is identical to that of ASL16 for straightforward comparison between the HJ formation fraction from LK migration with and without chaotic tides.

The details of how the evolution of the planetary f-mode is calculated in combination with the LK effect are discussed in Section 3.1. In all calculations, we use $E_{\rm max}=0.1GM_{p}^{2}/R_{p}$ and $E_{\rm resid}=0.001GM_{p}^{2}/R_{p}$ to parametrize the maximum f-mode energy and the residual energy after an episode of non-linear tidal dissipation. The planet is assumed to be rotating at the pseudo-synchronous rate (equation 29). The effect of dynamical tides on the planet’s spin is not accounted for in this investigation, but could be incorporated into future studies that follow the evolution of multiple oscillatory modes in the planet.

Our study is structured to efficiently determine the likelihood of three possible outcomes: “No Chaotic Tidal Migration", “Tidal Disruption", and ‘Chaotic Tidal Migration". Each calculation is stopped when the integration time reaches 1 Gyr or when one of the following conditions is met:

(i) If $|\Delta\hat{P}_{\alpha,k}(E_{\alpha,k-1})|$ (equation 38) has not reached 1 after min$(500t_{\rm LK},5t_{\rm LK}/\epsilon_{\rm oct})$, the calculation is terminated and the system is labelled as “No Chaotic Tidal Migration". Note that a small percentage of such systems may migrate within 1 Gyr via standard LK migration without the assistance of chaotic tides. These are not identified as a separate outcome from planets that will not migrate at all.

(ii) If the system has evolved for more than $10^{7}$ planetary orbits in the chaotic phase (i.e. where the orbital evolution is coupled with the f-mode evolution, as discussed in Section 3.1) without reaching $E_{\alpha,k}=E_{\rm max}$, the integration is stopped and the system is classified as “No Chaotic Tidal Migration".

(iii) If, at any time, the pericentre distance is within the tidal disruption radius, i.e. $r_{\rm p}<r_{\rm p,dis}$, the outcome is classified as “Tidal Disruption", and the integration is stopped.

(iv) If the planet has experienced chaotic orbital evolution and its orbit has circularized to $e=0.1$ within 1 Gyr, the integration is stopped and the system is classified as having undergone “Chaotic Tidal Migration".

For each system, we begin by integrating the full equations of motion and incorporating the evolution of the planet f-mode (as described in Section 3.1). However, at different points in the giant planet’s evolution, it is possible to “turn off" various effects without losing accuracy. As discussed in Section 3.1, the effect of dynamical tides is only accounted for when it can significantly affect orbital evolution. Additionally, when the semi-major axis is small, LK oscillations are suppressed. As the planet’s orbit shrinks and circularizes on a long time-scale due to static tides [see equation (43) with $\Delta t_{L}=1$ s], the eccentricity of the planet’s orbit precesses on a much shorter time-scale than that of tidal decay. Following this precession is computationally expensive and does not affect the final properties of the migrating planet. When the LK oscillations are “frozen" [$\epsilon_{\rm Tide}/(1-e^{2})^{9/2}>30$] and the spin-orbit angle is safely adiabatic ($\mathcal{A}>50$), we continue the integration without LK and SRF terms.

5.2 Migration and WJ/HJ Formation Fractions

The results of our population synthesis are summarized in Table 2. For each combination of planetary mass and radius, we ran $10^{4}$ calculations to determine the fraction of systems that undergo chaotic tidal migration ($\mathcal{F}_{\rm mig}$). These migrated planets are either tidally disrupted (fraction $\mathcal{F}_{\rm dis}$) or survive as eccentric WJs (fraction $\mathcal{F}_{\rm WJ}=\mathcal{F}_{\rm mig}-\mathcal{F}_{\rm dis}$) that undergo further orbital decay and circularization due to static tides, eventually becoming HJs. The fraction of systems that evolve into HJs within 1 Gyr is denoted by $\mathcal{F}_{\rm HJ}$. Note that our population synthesis did not sample initial mutual inclinations with $|\cos(\theta_{\rm lb,0})|>0.77$, as such planets do not experience large excursions in eccentricity. The migration, disruption, WJ, and HJ fractions are calculated assuming a uniform distribution in $\cos(\theta_{\rm lb,0})$, e.g. $\mathcal{F}_{\rm mig}=0.77N_{\rm mig}/N_{\rm run}$ with $N_{\rm mig}$ the number of systems that displayed chaotic tidal migration and $N_{\rm run}$ the total number of runs. For $1M_{J}$, $1R_{p}$ planets, the HJ formation percentage from our population synthesis is larger than the value obtained for standard LK migration by a factor of $\sim 1.6$ (3.9 % vs 2.4 %); see(ASL16; Muñoz et al., 2016).

5.3 Parameter Space for WJ Formation

We can learn about the parameter space that leads to eccentric WJ formation by examining how the outcome of orbital evolution depends on the initial conditions and the planet properties (see Fig. 14; see also Fig. 18 from ASL16 for comparison with the standard LK migration results). It is clear that smaller values of $\theta_{\rm lb,0}$ do not yield migration. This is unsurprising as the eccentricities required for chaotic tides are very large, which necessitates large initial mutual inclinations. Another pronounced feature is that the systems that produce eccentric WJs are clustered in $a_{\rm b}/a_{0}$. The reason for this is discussed in Section 4 (see Fig. 10). In essence, only a narrow window in $a_{\rm b}/a_{0}$ can produce systems that will reach large enough eccentricities to undergo chaotic tides (see equation 58), but not so large that the planets are tidally disrupted (see equation 62). Figure 15 demonstrates that the systems that produce surviving WJs indeed satisfy equations (58) and (62). The range of $\epsilon_{\rm oct}$ that can produce eccentric WJs is capped by the limit on $a_{\rm b}/a_{0}$ (see equation 71). Lastly, Fig. 14 shows that changing the properties of the gas giant has little effect on the parameter space that produces eccentric WJs. Unlike standard LK migration with static tides, where systems with $M_{p}=0.3M_{J}$ produce hot Saturns at a rate of $0.5\%$ (assuming $\Delta t_{\rm L}=1$ s) after evolving for 5 Gyr (see ASL16, Muñoz et al., 2016), chaotic tides allow low-mass planets to survive high-eccentricity migration without suffering tidal disruption, and produce hot Saturns at a rate of about $2.2\%$ after only 1 Gyr.

5.4 Properties of WJs Formed by Chaotic Tidal Evolution

We now examine the orbital properties of planets produced by chaotic tides. The scatter plots in Fig. 16 show an overview of the WJ eccentricities and semi-major axes ($e_{\rm WJ}$ and $a_{\rm WJ}$) for the $M_{p}=1M_{J},R_{p}=1.6R_{J}$ planet model after the planetary f-mode has settled into a quiescent state (but before weak tidal friction circularizes the orbit over a long time-scale). Most planets have $a_{\rm WJ}$ between $0.2$ and $1$ au and $e_{\rm WJ}$ between 0.91 and 0.97. Systems that have undergone multiple LK cycles tend to have larger $a_{\rm WJ}$ and $e_{\rm WJ}$. This is expected because larger values of $a_{\rm WJ}$ and $e_{\rm WJ}$ correspond to systems that do not reach a small enough value of $r_{\rm p}$ for chaotic tides to significantly alter the orbit within the time-frame of a single high-eccentricity phase. In general, $a_{\rm WJ}$ and $e_{\rm WJ}$ increase with $\eta=r_{\rm p}/r_{\rm tide}$ [see equations (64) and (65) and Fig. 12].

In the top left panel of Fig. 16, the values of $a_{\rm WJ}$ seem to fall on a set of $a_{\rm WJ}$ - $a_{0}$ curves. This feature arises from our treatment of non-linear dissipation, where the f-mode rapidly dissipates energy after reaching $E_{\rm max}$. Each curve corresponds to some number $N$ of non-linear dissipation episodes, after which the orbital energy becomes $(E_{B,0}-NE_{\rm max})$. The systems with the most distant stellar companions can produce WJs with relatively large semi-major axes ($\gtrsim$1 au) after chaotic tidal evolution. This occurs because the minimum $\eta$ required for chaotic tides to compete with the LK effect is larger for systems with larger $a_{\rm b,eff}$. In Fig. 16, we can also see that most systems that exhibit chaotic behaviour have initial mutual inclinations around $90^{\circ}$, as expected. Lastly, as seen earlier in Fig. 14, large values of $\epsilon_{\rm oct}$ do not produce surviving WJs.

Figure 17 shows the histogram of the parameters of WJs produced by chaotic tidal migration. The distributions of planet properties are generally more sharply peaked for denser planets. The ranges of $a_{\rm WJ}$ and $e_{\rm p,WJ}$ are in agreement with equations (64) and (65), which predict that $a_{\rm WJ}$ scales linearly with $R_{p}$ while $e_{\rm WJ}$ has no dependence on the planet radius and a weak dependence on the planet mass (for a given minimum $\eta$ achieved during the orbital evolution). Combining equations (64) and (65), we can see that the distribution of $r_{\rm p,WJ}=a_{\rm WJ}(1-e_{\rm WJ})$ should scale directly with $r_{\rm tide}$. This would imply that the peak at $0.015$ au in the $r_{\rm p,WJ}$ distribution for the $M_{p}=1M_{J},R_{p}=1\;R_{J}$ planet model should appear near $0.024$ au and $0.036$ au in the $M_{p}=1M_{J},R_{p}=1.6\;R_{J}$ and $M_{p}=0.3M_{J},R_{p}=1.6\;R_{J}$ models, respectively, as is the case in Fig. 17. The timescale for chaotic evolution to shrink the orbit, $t_{\rm WJ}$, peaks near 10 Myr for all planet models.

The range of $r_{\rm p,WJ}$ is narrow for all three planet models. The lower edge of this distribution is determined by the tidal disruption radius. A planet cannot survive if $\eta<2.7$. The upper value of $r_{\rm p,WJ}$ is roughly set by the pericentre distance where the planet with the largest semi-major axis (in our simulation, $a_{0}=5$ au) crosses the chaos boundary, generally near $\eta\sim 4$. Figure 18 shows the relationship between $a_{\rm WJ}$, $e_{\rm WJ}$ and $r_{\rm p,WJ}$ for planets that have undergone chaotic evolution. The solid lines show equation (65) for each planet model. Recall that equation (65) yields the eccentricity where $|\Delta\hat{P}_{\alpha}(E_{\rm resid})|=1$ for a given $r_{\rm p}$. Therefore, the restrictions on $r_{\rm p,WJ}$ determine the range in the distributions of $a_{\rm WJ}$ and $e_{\rm WJ}$. This simple calculation explains the properties of planets that survive chaotic tidal evolution very well. The spread in the results from the population synthesis arises from the fuzziness in the chaotic tides boundary, i.e. variations in $\Delta\hat{P}_{\rm crit}$. These variations are larger for smaller $E_{\rm resid}$ (see Fig. 4). On average, systems with smaller $M_{p}$ and larger $R_{p}$ have smaller $E_{\rm resid}$. Accordingly, the spread in the results is largest for the $M_{p}=0.3M_{J},R_{p}=1.6\;R_{J}$ model and smallest for the $M_{p}=1M_{J},R_{p}=1\;R_{J}$ model.

5.5 Hot Jupiter Properties

The eccentric WJs formed by chaotic tidal migration continue to experience orbital decay and circularization due to static tides. As long as the static tide is sufficiently dissipative, the planet will circularize to

[TABLE]

For our population synthesis, most WJs circularize within 1 Gyr (assuming $\Delta t_{\rm L}=1$ s) to become HJs (see Table 2). Those that do not are almost exclusively from the low-density planet model with the largest $r_{\rm tide}$. The HJ period distributions for our population synthesis calculations are shown in Fig. 19. The peak lies between 3.5 and 4 days for our standard model with $M_{p}=1M_{J},R_{p}=1.6R_{J}$; this shifts to smaller periods for denser planets (with $M_{p}=1M_{J}$, $R_{p}=1R_{J}$). As discussed in Section 5.4, the range of $r_{\rm p,WJ}$ is set by the pericentre distance requirements for tidal disruption and for chaotic tidal behaviour. Both of these conditions vary directly with $r_{\rm tide}$, so the peak of the HJ period distribution scales as $r_{\rm tide}^{3/2}$.

The final spin-orbit misalignments ($\theta_{\rm sl}$) are also shown in Fig. 19. For all three planet models, the distribution in $\theta_{\rm sl}$ is bimodal, peaking near $30^{\circ}$ and $140^{\circ}$. These distributions are qualitatively similar to those obtained in ASL16 for standard LK migration with static tides. In general, the percentage of retrograde configurations ($\theta_{\rm sl}>90^{\circ}$) is larger for chaotic tidal migration than for standard LK migration. This may be due to the fact that, because chaotic tides dramatically speed up orbital decay, the star is rotating more rapidly when $\theta_{\rm sl}$ becomes “frozen," and the feedback torque from the star on the orbit is larger.

We have also provided the final values of $\theta_{\rm jb}$, the angle between $\bm{J}=\bm{S_{*}}+\bm{L}$ and $\bm{L}_{\rm B}$. In all of our calculations, the ratio $|\bm{S_{*}}|/|\bm{L}|$ is still of order unity near $t_{\rm HJ}$ (the time when the planet’s orbit has circularized to $e=0.1$), so the mutual inclination of the inner and outer orbits, $\theta_{\rm lb}$, is not fixed at $t_{\rm HJ}$. As the star continues to spin down, $\theta_{\rm lb}$ will approach $\theta_{\rm jb}$, as discussed in Section 4.5. The distributions of $\theta_{\rm jb}$ have two strong peaks around $65^{\circ}$ and $115^{\circ}$.

Lastly, Fig. 19 shows the distribution of “arrival times" for HJs, $t_{\rm HJ}$, when $e=0.1$. For the $M_{p}=1M_{J},R_{p}=1R_{J}$ planet model, all orbits circularize to $e=0.1$ within 1 Gyr. For the $M_{p}=0.3M_{J},R_{p}=1.6R_{J}$ model, only 65% of the WJs circularize within a Gyr timeframe, yet 97% can become HJs within 10 Gyr. In general, planets with larger $r_{\rm tide}$ finish chaotic evolution at larger $r_{\rm p,WJ}$ and take longer to circularize via static tides. Most planets that undergo chaotic tidal evolution can be expected to become HJs within the lifetimes of their host stars. As a result, the population of eccentric WJs formed by chaotic tides is transient (assuming $\Delta t_{\rm L}=1$ s).

6 Analytical Calculation of WJ Formation Fraction

Muñoz et al. (2016) developed an analytical method to calculate the HJ formation fraction in the standard LK migration (with static tides) scenario. In this section, we adapt this method to calculate the fractions of systems that experience chaotic tidal migration or tidal disruption. A system will become chaotic or suffer disruption when the planet’s pericentre distance is smaller than a critical value — for chaotic behaviour, $r_{\rm p,crit}$ from equation (34), for disruption, $r_{\rm p,dis}$ from equation (60). We define $\theta_{\rm lb,crit}$ as the minimum initial mutual inclination necessary for a system to reach $r_{\rm p,crit}$, or equivalently $e_{\rm crit}$ (see equation 32). Assuming uniform distribution in $\cos{\theta_{\rm lb,0}}$, the fraction of systems (for a given $a_{0}$, $a_{\rm b}$, $e_{\rm b}$) that experience chaotic tidal behaviour is then

[TABLE]

where $\cos{\theta_{\rm lb,crit}}=0$ when $e_{\rm lim}<e_{\rm crit}$. We refer to Muñoz et al. (2016) for the calculation of $\cos{\theta_{\rm lb,crit}}$ given $e_{\rm crit}$. By performing a Monte Carlo integration over $a_{0}=[1,5]$ (au), $a_{\rm b}=[10^{2},10^{3}]$ (au) (uniformly sampled in $\log{a_{\rm b}}$), and $e_{\rm b}=[0,0.8]$, we can predict the fraction of systems in our population synthesis where the planet undergoes chaotic tidal migration ($\mathcal{F}_{\rm mig}$) and fraction of systems where the planet is tidally disrupted ($\mathcal{F}_{\rm dis}$). The fraction of systems that become chaotic but are not disrupted is given by

[TABLE]

As discussed in Section 4.2, chaotic tides can compete with the LK effect and prevent planet tidal disruption. Our analytic calculations therefore provide an upper limit to the tidal disruption fraction ($\mathcal{F}_{\rm dis}$) and a lower bound on the fraction of WJs formed ($\mathcal{F}_{\rm WJ}$) in this scenario. Table 2 shows that our analytical results are in good agreement with the populaion synthesis results. Our analytical $\mathcal{F}_{\rm mig}$ is slightly smaller than the value obtained from population synthesis, likely because the boundary for chaotic tides is fuzzy and we have chosen a conservatively low value for $r_{\rm p,crit}$.

In addition to calculating $\mathcal{F}_{\rm mig}$ and $\mathcal{F}_{\rm dis}$, we can estimate the total migration fraction $\mathcal{F}^{\;\rm st+ct}_{\rm mig}$ where the planet either migrates via chaotic tides or migrates within 1 Gyr due to weak tidal friction (the standard LK migration scenario). For a given $a_{0},a_{\rm b}$, and $e_{\rm b}$, this is

[TABLE]

where $\cos{\theta_{\rm lb,st}}$ is the maximum value of $\cos{\theta_{\rm lb,0}}$ for which a system can reach $r_{\rm p,ST}$, given in equation (59). We can integrate $f^{\;\rm st+ct}_{\rm mig}$ over orbital parameters $a_{0}$, $a_{\rm b}$, and $e_{\rm b}$ (as before) to find $\mathcal{F}^{\;\rm st+ct}_{\rm mig}$. Because $r_{\rm p,ST}$ is generally less than $r_{\rm p,crit}$, the percentage of systems that are not chaotic but become HJs due to weak tidal friction is much less than $1\%$.

In our population synthesis calculations (Section 5) we have only considered planets with $M_{p}=M_{J}$ and $0.3M_{J}$ and radius $R_{p}=R_{J}$ or $1.6\;R_{J}$. With analytic calculations, we can easily predict the fractions of systems that survive chaotic tides or suffer disruption as a function of $M_{p}$ and $R_{p}$. The results are shown in Fig. 20. The fraction of systems that undergo chaotic tidal migration is nearly independent of $M_{p}$, $R_{p}$ ($\mathcal{F}_{\rm mig}\sim 13-15$%). The WJ formation fraction from LK chaotic tides is also insensitive to $M_{p}$ and $R_{p}$ because $r_{\rm p,crit}$ and $r_{\rm p,dis}$ have roughly the same scaling with $M_{p}$ and $R_{p}$. This result differs from the standard picture of LK migration, where Saturn-mass planets that migrate are destined for tidal disruption (see ASL16).

7 Summary and Discussion

7.1 Summary of Key Results

In this paper, we have systematically studied the role of dynamical tides in high-eccentricity gas giant migration via the LK effect. Our investigation has revealed (i) the conditions under which dynamical tides in the planet have a significant effect on migration, (ii) the characteristics of the planet population that results from LK migration with dynamical tides, and (iii) the expected HJ fraction from this formation path. The overall summary of the paper is already given in the abstract. Here we provide a guide to the key results of each section.

In Section 2, we explored the orbital parameters necessary for a gas giant to undergo chaotic tidal migration. When a planet’s orbit is very eccentric and its pericentre distance is sufficiently small, the tidal force from the host star can excite planetary oscillations at pericentre in such a way that these oscillations chaotically climb in amplitude over repeated close passages (an example is shown in the bottom panel of Fig. 2). We found that, for the $\gamma=2$ polytrope planetary model, the fundamental mode (f-mode) is more strongly excited than inertial modes, and is the most likely to experience chaotic growth over multiple orbits. The condition for these “chaotic tides” to operate is given by equation (32) and plotted in Fig. 3. In general, chaotic tides can occur in planets on highly eccentric orbits with small pericentre distances. We also explored (in Sections 2.4 and 2.5) how the conditions for chaotic tides become more generous when the f-mode already has some non-zero energy.

In Section 3, we developed a model to couple the LK migration with orbital changes due to dynamical tides. The model is fully described in Section 3.1, and can be readily adapted to other high-eccentricity migration mechanisms. Although the model is largely ab initio, there are two free parameters that characterize the (uncertain) nonlinear mode dissipation: $E_{\rm max}$ (the maximum energy the mode can reach before nonlinear dissipation sets in) and $E_{\rm resid}$ (the residual energy in the mode after an episode of nonlinear dissipation). We found that these parameters (with reasonable values) do not change the general features of evolution due to chaotic tides, but can alter the details (see Figs. 7-8).

In Sections 3.2 and 4, we presented example calculations of the LK chaotic tidal migration, and explained some key features of this migration mechanism:

• •

Chaotic tidal migration leads to rapid formation of WJs (see Fig. 6 for an example) on timescales a few to 100 Myrs (see Fig. 17). These WJs have eccentricities in the range $0.85\lesssim e_{\rm WJ}\lesssim 0.98$ and semi-major axes between 0.1 and 2 au (with lower values more likely), depending on the giant planet mass and radius (see Fig. 16).

• •

These eccentric WJs efficiently circularize to HJs via dissipation of static tides. This process is faster than standard LK migration because a planet that experiences chaotic tides quickly becomes decoupled from the influence of the perturber and the eccentricity oscillations are frozen (see Figs. 6, 9 and 11).

• •

Some giant planets that are otherwise fated for tidal disruption can be saved by chaotic tides (see Fig. 11). This occurs when the orbital decay timescale associated with chaotic tides is shorter than the timescale for the perturber to change the planet’s eccentricity (see Section 4.2).

In Section 5, we conducted a population synthesis study to determine the formation fractions and properties of eccentric WJs and HJs for a few combinations of planet mass and radius. Our study sampled the same range of initial orbital parameters as in Anderson et al. (2016) for the standard LK migration scenario. The calculated fractions for chaotic tidal migration ($\mathcal{F}_{\rm mig}$), tidal disruption ($\mathcal{F}_{\rm dis}$), WJ formation ($\mathcal{F}_{\rm WJ}$; note that $\mathcal{F}_{\rm mig}=\mathcal{F}_{\rm dis}+\mathcal{F}_{\rm WJ}$) and HJ formation within 1 Gyr ($\mathcal{F}_{\rm HJ}$), are provided in Table 2. For an $M_{p}=1M_{J},R_{p}=1R_{J}$ planet, LK chaotic tidal migration produces more HJs than standard LK migration ($\mathcal{F}_{\rm HJ}$ increases from $2.4\%$ to $3.9\%$; see Table 3 in ASL16). More importantly, chaotic tides can produce hot giant planets with a broad range of masses and radii. This is in stark contrast with the standard LK migration, which produces very few hot Saturns because of the severe tidal disruption experienced by such low-mass giants (ASL16; Munoz et al. 2016).

The orbital properties of planets from our population synthesis that survived chaotic tidal migration are shown in Figs. 17 and 19. The WJs that are produced directly by chaotic tidal dissipation have pericentre distributions that peak near $r_{\rm p,WJ}\sim 3r_{\rm tide}=3R_{p}(M_{\star}/M_{p})^{1/3}$. These pericentre values are smaller than those of observed high-eccentricity WJs such as HD80606 b (Hébrard et al., 2010). This is expected as the eccentric WJs formed by chaotic tides are “transient” and should move quickly through the high-eccentricity phase. The HJs that form via LK chaotic tidal migration exhibit a pile-up around a 3 day orbital period, depending on the planet’s mass-radius relation and the assumed tidal disruption criterion. This is a feature of all tidal migration mechanisms. However, our predicted HJ period distribution for LK chaotic tidal migration differs from that for the standard LK migration, particularly for low-mass giant planets (compare Fig. 19 with the middle row of Fig. 23 from ASL16). Chaotic tidal migration can produce Saturn-mass planets at periods longer than 5 days and generates a wider period distribution for such planets. Finally, we found that chaotic tidal migration yields similar spin-orbit misalignments as the standard LK migration (compare Fig. 19 with Fig. 24 of ASL16). For all three planet models we considered, the distribution of final spin-orbit misalignments is bimodal with peaks at $\theta_{\rm sl}^{\rm F}\sim 30^{\circ}$ and $\theta_{\rm sl}^{\rm F}\sim 130^{\circ}$.

In Section 6, we used an analytical method, developed in Muñoz et al. (2016), to understand how the HJ formation rate varies with planet mass and radius. We calculated a lower-bound on the eccentric WJ formation rate (not accounting for planets that are spared from tidal disruption by chaotic tides). The predicted (analytical) rates are in agreement with the results from the population synthesis study (see Section 5) and are shown in Fig. 20. We confirmed that the WJ formation rate from LK chaotic tidal migration is constant over a reasonable range of giant planet masses and radii, as suggested by, e.g., Fig. 14 and Table 2.

7.2 Discussion

The results presented in this paper show that chaotic tides endow the LK migration scenario with a number of “favorable” features (see also Wu, 2018). These not only reduce the theoretical uncertainties regarding tidal dissipation that are inherent in the theory, but also may help reconcile some of the discrepancies between observations and predictions of LK migration. Chaotic tides drastically reduce the amount of time that a gas giant spends at high eccentricity. This could explain the lack of observations of super-eccentric gas giants (Dawson et al., 2015). In addition, chaotic tides quickly decouple a gas giant experiencing LK oscillations from the stellar perturber. This allows planets at larger pericentre distances to migrate within the lifetime of their host stars. In this way, chaotic tidal migration naturally produces a period distribution with a longer tail. Indeed, HJs beyond the 3-day pile-up are observed but difficult to explain with standard high-e migration scenarios.

Although chaotic tides increase the HJ yield from LK migration in stellar binaries, particularly for low-mass planets, we should not expect this particular formation channel to account for all HJs. The occurrence rate of HJs produced by this channel can be computed from

[TABLE]

where $\mathcal{F}_{b}$ is the fraction of stars with a binary companion and $\mathcal{F}_{p}$ is the fraction of solar-type stars with a giant planet at a few au. Assuming $\mathcal{F}_{b}\sim 50\%$ (Raghavan et al., 2010; Ngo et al., 2015), $\mathcal{F}_{p}\sim 10\%$ (as in ASL16), and $\mathcal{F}_{\rm HJ}\sim 4\%$ from our population synthesis calculation, we obtain an estimate of $\mathcal{R}_{\rm HJ}\sim 0.2\%$, which is nearly an order of magnitude smaller than the observed occurrence rate of $1\%$ (Marcy et al., 2005; Wright et al., 2012; Fressin et al., 2013). Thus, LK chaotic tidal migration in stellar binaries can roughly account for $20-30\%$ of the observed HJ population.

However, we expect that many of the “nice” features of chaotic tides may also apply to other flavours of high-$e$ migration scenarios, such as LK migration induced by planetary companions and secular chaos in multi-planet systems. The eccentricity and pericentre ranges of planets that are susceptible to chaotic tides are set by the boundary for chaotic f-mode behaviour (see Section 2), and do not depend on a specific high-$e$ migration scenario. Chaotic tides can save a planet from tidal disruption when the time-scale for energy transfer to the planet’s oscillation mode is shorter than the time-scale for driving and maintaining the planet’s high-eccentricity (see Section 4.2) – this condition can be satisfied by all secular eccentricity excitation mechanisms. Indeed, the recent work by Teyssandier et al. (2018) showed that chaotic tides significantly increase the HJ formation fraction in the secular-chaos high-$e$ migration scenario. Overall, chaotic tides boost the importance of high-$e$ migration for the formation of HJs.

The story of chaotic tidal migration hinges upon the planet’s ability to survive rapid tidal heating. There is reason to expect that the planet interior could survive this process (Wu, 2018) if most of the dissipated energy goes into the outer layers of the planet, which can quickly radiate heat. However, if the tidal energy is deposited at a larger depth, the planet’s envelope may expand. Such changes in the structure of the planet would affect the f-mode frequency and could have a larger influence on the evolution of the planet’s orbit. A sudden expansion of the planet’s radius may also put the planet in danger of tidal disruption. Many of the planets that survive chaotic tidal migration come close to the tidal disruption radius. Expansion of the planet could easily lead to stripping of the outer layers and initiate mass transfer or mass loss. The effect of mode energy dissipation on the planetary structure is a very important problem for future study.

Acknowledgements

This work is supported in part by the NSF grant AST1715246 and NASA grant NNX14AP31G. MV is supported by a NASA Earth and Space Sciences Fellowship in Astrophysic. KRA is supported by a NASA Earth and Space Sciences Fellowship in Planetary Science.