Excitation of Planetary Obliquities Through Planet-Disk Interactions

TL;DR

This paper explores a simple early-stage mechanism where planet-disk interactions can excite planetary obliquities through secular spin-orbit resonance, affecting planets within 0.3 to 2 AU, but likely not explaining Uranus and Neptune's high obliquities.

Contribution

It identifies and analyzes a fundamental mechanism for obliquity excitation via planet-disk interactions and spin-orbit resonance crossing during planet formation.

Findings

Obliquity can be excited when disk-induced orbital precession frequency crosses spin precession frequency.

Planets between 0.3 and 2 AU are most affected by this mechanism.

The mechanism is unlikely to explain Uranus and Neptune's high obliquities due to gravitational perturbations.

Abstract

The tilt of a planet's spin axis off its orbital axis ("obliquity") is a basic physical characteristic that plays a central role in determining the planet's global circulation and energy redistribution. Moreover, recent studies have also highlighted the importance of obliquities in sculpting not only the physical features of exoplanets but also their orbital architectures. It is therefore of key importance to identify and characterize the dominant processes of excitation of non-zero axial tilts. Here we highlight a simple mechanism that operates early on and is likely fundamental for many extrasolar planets and perhaps even Solar System planets. While planets are still forming in the protoplanetary disk, the gravitational potential of the disk induces nodal recession of the orbits. The frequency of this recession decreases as the disk dissipates, and when it crosses the frequency of a planet's spin axis precession, large planetary obliquities may be excited through capture into a secular spin-orbit resonance. We study the conditions for encountering this resonance and calculate the resulting obliquity excitation over a wide range of parameter space. Planets with semi-major axes in the range $0.3 \ \mathrm{AU} \lesssim a \lesssim 2 \ \mathrm{AU}$ are the most readily affected, but large-$a$ planets can also be impacted. We present a case study of Uranus and Neptune and show that this mechanism likely cannot help explain their high obliquities. While it could have played a role if finely tuned and envisioned to operate in isolation, large-scale obliquity excitation was likely inhibited by gravitational planet-planet perturbations.

Figures

Tables

Table 1. Table 1: Nodal secular frequencies for the giant planets in the Solar System (Laskar, 1988 ) . In the notation of Murray & Dermott ( 1999 ) , the subscripts i 𝑖 i are the eigenmode numbers, and f i subscript 𝑓 𝑖 f_{i} are the eigenfrequencies associated with the inclination/node solution. Planet i 𝑖 i is dominated by the mode of the same number.

$i$	$f_i$ (arcsec yr-1)
5	-
6	-26.34
7	-2.99
8	-0.69

Equations

$$T_{\alpha }=2\pi /(\alpha \cos ϵ),$$

$${\alpha }_0=\frac{k_2}{2C}n{f}_{\omega }{\left[\frac{M_{\star }}{M_p}{\left(\frac{R_p}{a}\right)}^3\right]}^{\frac{1}{2}}.$$

$$J_2=\frac{{\omega }^2{R_p}^3}{3G{M}_p}{k}_2.$$

$$\alpha =f_{\alpha }{\alpha }_0=\frac{3}{2}\frac{n}{f_{\omega }}\left[\frac{(k_2{f_{\omega }}^2/3)+q}{C+l}\right]{\left[\frac{M_{\star }}{M_p}{\left(\frac{R_p}{a}\right)}^3\right]}^{\frac{1}{2}}.$$

$$q$$

$$l$$

$$q$$

$$l$$

$$|{\nu }_s|\gg \alpha ,$$

$${\nu }_s={\dot{\Omega }}_s=-\frac{3}{2}{n}_s{J}_2{\left(\frac{R_p}{a_s}\right)}^2.$$

$$\frac{a_s}{R_p}\ll {\left[{\left(\frac{a}{R_{\star }}\right)}^3\frac{f_{\omega }C}{f_{\alpha }}\frac{{\rho }_p}{{\rho }_{\star }}\right]}^{\frac{2}{7}},$$

$$g=-\dot{\Omega }=\frac{3}{4}n\left(\frac{-\beta +2}{1+\beta }\right)\left(\frac{1-{\eta }^{-1-\beta }}{1-{\eta }^{-\beta +2}}\right)\left(\frac{M_d}{M_{\star }}\right){\left(\frac{a}{R_o}\right)}^3.$$

$$\Sigma (a,t)={\Sigma }_0(t){\left(\frac{a}{R_o}\right)}^{-\beta },$$

$$M_d(t)=\frac{2(1-{\eta }^{-\beta +2})\pi {R_o}^2}{-\beta +2}{\Sigma }_0(t)$$

$$\begin{split}\frac{\alpha}{g}&=0.004\bigg{(}\frac{k_{2}f_{\omega}f_{\alpha}}{C}\bigg{)}\bigg{(}\frac{\rho_{p}}{\mathrm{\footnotesize{g/cm^{3}}}}\bigg{)}^{-\frac{1}{2}}\bigg{(}\frac{M_{\star}}{M_{\odot}}\bigg{)}^{\frac{1}{2}}\\ &\times\bigg{(}\frac{\eta}{0.001}\bigg{)}^{\frac{5}{2}}\bigg{(}\frac{M_{d}/M_{\star}}{0.01}\bigg{)}^{-1}\bigg{(}\frac{R_{o}}{100\ \mathrm{AU}}\bigg{)}^{3}\bigg{(}\frac{a}{\mathrm{AU}}\bigg{)}^{-\frac{9}{2}}.\end{split}$$

$$\begin{split}a&>0.3\ \mathrm{AU}\bigg{[}\bigg{(}\frac{k_{2}f_{\omega}f_{\alpha}}{C}\bigg{)}\bigg{(}\frac{\rho_{p}}{\mathrm{\footnotesize{g/cm^{3}}}}\bigg{)}^{-\frac{1}{2}}\bigg{(}\frac{M_{\star}}{M_{\odot}}\bigg{)}^{\frac{1}{2}}\\ &\times\bigg{(}\frac{\eta}{0.001}\bigg{)}^{\frac{5}{2}}\bigg{(}\frac{M_{d}/M_{\star}}{0.01}\bigg{)}^{-1}\bigg{(}\frac{R_{o}}{100\ \mathrm{AU}}\bigg{)}^{3}\bigg{]}^{\frac{2}{9}}.\end{split}$$

$$\dot{\alpha }-\dot{g}\lesssim \alpha g\sin {ϵ}_0\sin i,$$

$$\frac{\alpha\ \tau_{\scriptscriptstyle d}\ i\sin 2\epsilon_{0}}{2}\bigg{(}\frac{1}{\cos\epsilon_{0}-f_{\tau}}\bigg{)}>1.$$

$$\alpha\ \tau_{\scriptscriptstyle d}\ i\epsilon_{0}\bigg{(}\frac{1}{1-f_{\tau}}\bigg{)}>1.$$

$$\begin{split}a&<1.93\ \mathrm{AU}\bigg{[}\frac{k_{2}f_{\omega}f_{\alpha}}{C}\bigg{(}\frac{\rho_{p}}{\mathrm{\footnotesize{g/cm^{3}}}}\bigg{)}^{-\frac{1}{2}}\bigg{(}\frac{1}{1-f_{\tau}}\bigg{)}\\ &\times\bigg{(}\frac{\tau_{\scriptscriptstyle d}}{\mathrm{Myr}}\bigg{)}\bigg{(}\frac{M_{\star}}{M_{\odot}}\bigg{)}\bigg{(}\frac{i}{1^{\circ}}\bigg{)}\bigg{(}\frac{\epsilon_{0}}{20^{\circ}}\bigg{)}\bigg{]}^{\frac{1}{3}}.\end{split}$$

$$\mathcal{H}=\frac{\alpha }{2}(1-e^2{)}^{-3/2}{X}^2+\dot{\Omega }i\sqrt{1-X^2}\cos (\Omega +\psi ).$$

$$\Phi =1-X\varphi =-\psi .$$

$$\begin{array}{rl}\mathcal{H}& =\frac{\alpha }{2}(1-e^2{)}^{-3/2}(1-\Phi {)}^2\\ & +\dot{\Omega }i\sqrt{1-(1-\Phi {)}^2}\cos (\dot{\Omega }t-\varphi )+\mathcal{T}.\end{array}$$

$${\mathcal{F}}_2=(\dot{\Omega }t-\varphi )\Phi +(t)\Xi .$$

$$\begin{split}\Theta&=-\Phi&\theta&=\dot{\Omega}t-\phi\\ \Xi&=\mathcal{T}-\dot{\Omega}\Theta&\xi&=t.\end{split}$$

$$\mathcal{H}=\frac{\alpha }{2}(1+\Theta {)}^2-g\Theta -gi\sqrt{1-(1+\Theta {)}^2}\cos \theta .$$

$$x=\sqrt{-2\Theta }\cos \theta y=\sqrt{-2\Theta }\sin \theta .$$

$$\mathcal{H}=\frac{\alpha}{2}\bigg{(}1-\frac{1}{2}x^{2}-\frac{1}{2}y^{2}\bigg{)}^{2}+\frac{1}{2}g(x^{2}+y^{2})-\frac{i}{2}gx\sqrt{4-x^{2}-y^{2}}.$$

$$M_d(a,t)=M_d(t=0)\exp (-t/{\tau }_d).$$

$$g_{pp}=\frac{n}{4}{b}_{3/2}^{(1)}(\alpha )\alpha \bar{\alpha }\frac{M^{\prime }}{M_{\star }}$$

$$b_{3/2}^{(1)}(\alpha )=\frac{1}{\pi }{\int }_0^{2\pi }\frac{\cos \psi }{(1-2\alpha \cos \psi +{\alpha }^2{)}^{3/2}}d\psi .$$

Full text

1 Introduction

The spin dynamics of planetary bodies influence many of their salient features, including atmospheric circulation, climate, and tidal evolution. The two most important features of the spin state are the rotation period and the obliquity, defined as the angle between the planet’s spin and orbital axes. The present-day obliquities of the Solar System planets are wide-ranging, from near-zero spin-orbit misalignment for Mercury, to $23^{\circ}$ for Earth, and $98^{\circ}$ for Uranus. Notably, five of the eight planets have obliquities larger than $20^{\circ}$.

Obliquities of extrasolar planets have not yet been measured directly.111Note these obliquities are distinct from the stellar obliquity, which has been constrained for many systems (Winn & Fabrycky, 2015). There are several prospects, however, for measuring obliquities in coming years. Ground-based, high-resolution spectroscopy permits measurements of the $v\sin i$ of the planet, offering a degenerate constraint on the spin rate and component of obliquity along the line-of-sight (Snellen et al., 2014; Bryan et al., 2018). Alternatively, measurements of planetary asphericity from transit photometry offer potential constraint on the spin rate, obliquity, and the coefficient of the gravitational quadrupole moment (Hui & Seager, 2002; Seager & Hui, 2002; Barnes & Fortney, 2003; Ragozzine & Wolf, 2009; Carter & Winn, 2010a, b; Zhu et al., 2014; Biersteker & Schlichting, 2017).

For close-in planets, features in full-phase, optical or infrared light curve observations may point toward non-zero obliquities (Schwartz et al., 2016; Rauscher, 2017; Ohno & Zhang, 2019a, b; Adams et al., 2019). Similarly, the detection of seasonal variability from infrared photometry or direct imaging has also been proposed as a means of constraining obliquity (Gaidos & Williams, 2004; Kawahara, 2016; Kane & Torres, 2017). A large obliquity may also be inferred indirectly by considering its effect on other planetary properties, such as its tidal dissipation rate (Millholland & Laughlin, 2018, 2019).

Despite the lack of direct detection thus far, there are theoretical expectations for the frequent occurrence of large planetary obliquities. Several means exist for exciting obliquities from their primordial values, and the Solar System planets are prime examples. First, collisions with massive proto-planetary cores, such as the Moon-forming giant impact, can produce large axial tilts. Similarly, giant impacts are routinely invoked in the origin scenarios of the extreme Uranian tilt (Safronov, 1966; Morbidelli et al., 2012). It is also possible that the obliquities of the giant planets were excited by an early process that twisted the total angular momentum vector of the Solar System, such as asymmetric infall during the system’s initial collapse and formation (Tremaine, 1991).

Another important obliquity-exciting mechanism is spin-orbit resonance – a natural syncing of the frequencies of a planet’s spin axis precession and orbital precession. Orbital precession is inevitably driven by planet-planet interactions, and spin-orbit resonances of this type likely explain Jupiter’s $3^{\circ}$ and Saturn’s $27^{\circ}$ obliquities through interactions with Uranus and Neptune, respectively (Ward & Hamilton, 2004; Hamilton & Ward, 2004; Ward & Canup, 2006; Vokrouhlický & Nesvorný, 2015; Brasser & Lee, 2015; Nesvorný, 2018). Furthermore, in the case of extrasolar planets, Millholland & Laughlin (2019) showed that enhanced tidal dissipation generated by obliquity excitations in spin-orbit resonances may be a key factor in sculpting the orbital period distribution of systems with multiple transiting planets.

We note further that when several spin-orbit resonances overlap, large-amplitude, chaotic excursions of the obliquity can result. This has likely been important for all of the Solar System terrestrial planets (Laskar & Robutel, 1993), especially Venus (Correia et al., 2003; Correia & Laskar, 2003) and Mars (Ward, 1973; Touma & Wisdom, 1993), where the obliquity varies chaotically between $\sim 10^{\circ}$ and $50^{\circ}$ on timescales of hundreds to thousands of years.

Although planet-planet interactions are often considered as a driver of orbital precession, they are not unique. In fact, any deviation from a purely $\propto 1/r$ potential will give rise to nodal recession. Natural sources that act as the system is still forming include the gravitational field of a protoplanetary disk and the quadrupole potential of a rapidly-rotating young star. The precession frequencies associated with these torques decrease in magnitude as the nebular gas dissipates and the young star contracts. If these evolutions result in adiabatic spin-orbit resonance crossings where the spin and orbital precession frequencies are equal, then obliquities will inevitably be excited.

In this work, we focus on the first of these two torque generators: the protoplanetary disk. We present an explication of the role of disk-induced orbital precession in exciting planetary obliquities. A schematic representation of this process is shown in Figure 1. As we show below, the mechanism is extensive and robust due to the large magnitude of the initial disk-induced precession rate and the wide range of commensurabilities that are swept. However, it has not yet been widely appreciated in the literature222Winn & Holman (2005) discussed disk-induced secular spin-orbit resonance in the context of hot Jupiters. In this same context, Fabrycky et al. (2007) examined the role of stellar quadrupole-induced resonances. Millholland & Laughlin (2019) also touched on the influence of a young stellar quadrupole in short-period, compact “Kepler-like” systems..

It is enlightening to point out that this mechanism is closely analogous to other secular resonances that might occur early on in a system’s lifetime. One theory for the excitation of stellar spin-orbit misalignments in extrasolar systems involves a resonant commensurability between stellar spin axis precession and nodal regression of a circumstellar disk in a star-disk-binary system (Batygin, 2012; Batygin & Adams, 2013; Lai, 2014; Spalding & Batygin, 2014; Zanazzi & Lai, 2018). Specifically, the stellar spin precesses about the disk angular momentum vector. Meanwhile, the disk angular momentum vector precesses due to a binary star perturber. As the disk dissipates, a resonance ensues when the spin precession and disk precession frequencies match, resulting in adiabatic excitation of the stellar obliquity. This scenario is mathematically similar to the one at hand, except we are now considering the precession of the planetary (rather than stellar) spin about the orbital (rather than disk) angular momentum vector, and this vector itself precesses due to the disk (as opposed to the binary perturber).

This paper is organized as follows. In Section 2, we define the frequencies of spin axis precession and disk-induced orbital precession. We then use these expressions in Section 3 to define analytic criteria for spin-orbit resonance crossing and capture. In Section 4, we develop a simple perturbative model to study the evolution of a planet’s obliquity when this resonance is encountered, and map the resulting obliquity excitation over a wide range of parameter space. Finally, we present a special case study of the Uranian and Neptunian obliquities in Section 5 and conclude in Section 6.

2 Spin and Orbital Precession

In order to study the evolution of the spin and orbital precession into resonant commensurability, we must first define their characteristic frequencies.

2.1 Spin axis precession

The torque from the host star on a rotationally-flattened planet causes the planet’s spin axis to precess about its orbit normal. The period of precession is

[TABLE]

where $\epsilon$ is the planet’s obliquity (angle between the spin and orbital axes) and $\alpha$ is the precessional constant. The magnitude of $\alpha$ may be strongly enhanced by the presence of close-in satellites or a primordial circumplanetary accretion disk. This more rapid precession is due to the adiabatic gravitational coupling between the oblate planet and surrounding disk/satellites, which makes the system precess as a unit (Goldreich, 1965).

In the absence of a disk or satellites, $\alpha=\alpha_{0}$ is given by (Ward & Hamilton, 2004)

[TABLE]

Here $M_{\star}$ is the host star mass, $M_{p}$ and $R_{p}$ the planet mass and radius, $a$ the semi-major axis of the planet, and $n=2\pi/P$ its mean motion. The quantity $k_{2}$ is the planet’s Love number, a dimensionless value related to the planet’s central concentration of the density profile and its deformation response to tidal disturbance. The quantity $C$ is the planet’s moment of inertia normalized by $M_{p}{R_{p}}^{2}$. Finally, $f_{\omega}=\omega/\omega_{b}$ is defined as the planet’s spin rate as a fraction of break-up speed, $\omega_{b}=\sqrt{GM_{p}/{R_{p}}^{3}}$. The purpose of this parameterization will become clear later. The definition of $\alpha_{0}$ has incorporated the following form of $J_{2}$, the coefficient of the quadrupole moment of the planet’s gravitational field (Ragozzine & Wolf, 2009),

[TABLE]

As already mentioned above, if the planet is encircled by a circumplanetary disk or a massive satellite, the effective value of $\alpha$ will be larger. If we define $f_{\alpha}=\alpha/\alpha_{0}$ as the enhancement factor of the precessional constant over its value when there is no disk/satellite, then (Ward & Hamilton, 2004)

[TABLE]

The quantity $q$ is the effective quadrupole coefficient of the satellite system or disk, and $l$ is the angular momentum of the satellites/disk normalized by $M_{p}{R_{p}}^{2}\omega$. In the case of a single equatorial satellite with mass, $m_{s}$, and semi-major axis, $a_{s}$, $q$ and $l$ are defined by

[TABLE]

Similarly, for a circumplanetary disk with mass, $m_{cp}$, radius, $r_{cp}$, and surface density profile, $\Sigma_{cp}=\Sigma_{cp,0}(r/R_{p})^{-\gamma}$ with $\gamma<2$, $q$ and $l$ are defined,

[TABLE]

We obtained these expressions by imagining the disk to be composed of a set of infinitesimal rings of mass and integrating.

To concretely examine typical values of $f_{\alpha}$, we must first define a set of fiducial system parameters. Here and in the remainder of the text unless otherwise noted, we use the parameters $k_{2}=0.125$ and $C=0.25$ (close to the values for Uranus and Neptune, Gavrilov & Zharkov 1977; Lainey 2016), $\rho_{\star}=1.41\ \mathrm{g/cm^{3}}$ (Solar density, which appears later), $\rho_{p}=1.27\ \mathrm{g/cm^{3}}$ (Uranus density), $\gamma=3/4$ (Canup & Ward, 2002), and $f_{\omega}=\omega/\omega_{b}=0.1$ (Bryan et al., 2018; Batygin, 2018). Additional parameters will be defined as necessary.

Figure 2 shows the enhancement of the spin axis precessional constant due to a satellite (solid lines) or circumplanetary disk (dashed lines). It is clear to see that the enhancement is significant, reaching several orders of magnitude for large mass ratios and orbital distances.

This $\alpha$-enhancement, however, is only active when the satellite/disk is close enough to the oblate planet to preserve adiabatic gravitational coupling. The adiabatic criterion requires that the rate of the satellite’s nodal recession due to the planet’s oblateness is much greater than the rate of the planet’s spin axis precession (Goldreich, 1965). This criterion may be expressed

[TABLE]

where $\nu_{s}$ is the frequency of the satellite’s nodal recession in the planet’s equatorial plane (Murray & Dermott, 1999),

[TABLE]

Using $\alpha=f_{\alpha}\alpha_{0}$ and equation 2 for $\alpha_{0}$, the adiabatic criterion reduces to

[TABLE]

where $\rho_{p}$ and $\rho_{\star}$ are the planet’s and star’s densities.

Equation 9 specifies the maximum satellite semi-major axis for which the gravitational coupling will be upheld for a given $\alpha$-enhancement and planet semi-major axis. An equivalent formula exists for a circumplanetary disk, where $a_{s}$ is simply replaced by $r_{cp}$. The black diagonal lines in Figure 2 show contours of the planet semi-major axis calculated by equality of the left-hand and right-hand sides of equation 9. At a given planet semi-major axis, satellite/disk orbital distances and $\alpha/\alpha_{0}$ enhancements below this line (toward the lower left) will maintain adiabatic gravitational coupling to the planet.

2.2 Orbit nodal recession

In addition to the spin precessional motion of the planet and satellite/disk system, there are also torques that drive orbital precession. The protoplanetary disk’s quadrupolar gravitational potential induces an orbit nodal recession for the planet with period, $T_{g}=2\pi/g$. For a planet on a circular orbit, $g$ is defined (Chen et al., 2013; Terquem, 2013)

[TABLE]

Here we have assumed the disk surface density profile follows a power law,

[TABLE]

where $R_{o}$ is the outer disk radius. The quantity $\eta=R_{i}/R_{o}$ is the ratio of inner and outer disk radii, and

[TABLE]

is the disk mass. We suppose that the disk profile follows the minimum-mass solar (or extrasolar) nebula (Hayashi, 1981; Chiang & Laughlin, 2013) and set the power law exponent equal to $\beta=3/2$.

Early evolution of the system can cause the spin precession and disk-induced orbital precession to reach a commensurable rate. The resulting resonant encounter can excite the planetary obliquity to large values, sometimes up to $90^{\circ}$. We will explore the conditions for this resonant crossing and capture in the next section.

3 Secular spin-orbit resonance

Capture into a secular spin-orbit resonance occurs when $T_{\alpha}=2\pi/(\alpha\cos\epsilon_{p})$ and $T_{g}=2\pi/g$ evolve such that $T_{\alpha}/T_{g}$ crosses through unity from above. As the protoplanetary disk decays, $\Sigma\rightarrow 0$ and $T_{g}\rightarrow\infty$. This crossing is therefore inevitable if $T_{g}<T_{\alpha}$ at some point during the system’s evolution. The $T_{g}<T_{\alpha}$ crossing condition is equivalent to $g\gtrsim\alpha$ when the initial obliquities are small. Using equations 4 and 10 above and setting $\eta\ll 1$, the ratio $\alpha/g$ simplifies to

[TABLE]

By setting the initial value of $\alpha/g$ to be less than 1, we find the minimum semi-major axis at which the resonance will be crossed:

[TABLE]

Resonant crossing is a necessary but insufficient criterion for resonant capture. Capture into the spin-orbit resonance also requires that the crossing is adiabatic, in other words, that the crossing timescale is slow in comparison to the resonant libration period (Hamilton & Ward, 2004). This criterion may be expressed

[TABLE]

where $\epsilon_{0}$ is the planet’s obliquity upon resonant crossing, and $i$ is the orbital inclination with respect to the invariable plane. The condition that $T_{\alpha}/T_{g}$ crosses unity from above (direction for capture) implies that $\dot{\alpha}-\dot{g}>0$.

Disk evolution leads $M_{d}(t)\rightarrow 0$. Though this decay may not be described by a simply-defined function of time, if it is adiabatic, the sole relevant timescale is the decay rate upon resonant crossing, $\tau_{\scriptscriptstyle d}$. With this definition, $g/\dot{g}=M_{d}/\dot{M}_{d}=-\tau_{\scriptscriptstyle d}$. Unlike the decrease in $g$, which is clearly linked to the disk decay, there are a variety of mechanisms that could drive evolution in $\alpha$. These include Kelvin-Helmholtz contraction as the planet cools, dissipation of the circumplanetary disk, and disk-driven migration. We define $\tau_{\scriptscriptstyle\alpha}=-\alpha/\dot{\alpha}$ and the ratio $f_{\tau}=\tau_{\scriptscriptstyle d}/\tau_{\scriptscriptstyle\alpha}$. Finally, noting that $g\sim\alpha\cos\epsilon_{0}$ at resonant crossing and assuming $i$ is small, the capture criterion may be written

[TABLE]

For small $\epsilon_{0}$, this reduces to

[TABLE]

This criterion highlights the fact that capture into a secular spin-orbit resonance cannot take place in a purely axi-symmetric system. Due to the infinite precession period at zero inclination, small but nevertheless finite inclination and obliquity are required for secular resonant coupling to ensue.

We substitute equation 4 for $\alpha$ and solve for the maximum semi-major axis for which the crossing is adiabatic and capture is guaranteed,

[TABLE]

Figure 3 displays the lower and upper semi-major axis limits resulting from the resonant crossing (equation 14) and capture (equation 18) criteria, respectively. For a given upper limit, the region above the line will result in resonant crossing but not capture because the adiabatic criterion is not upheld. (In that case the spin axis precession period is too long compared to the crossing timescale.) We plot these relations as a function of $f_{\omega}f_{\alpha}=(\omega/\omega_{b})(\alpha/\alpha_{0})$ for different crossing timescales. Notably, both $\tau_{d}$ and $\tau_{\alpha}$ affect the crossing timescale. The capture domain between the upper and lower limits widens as $\tau_{\alpha}$ approaches $\tau_{d}$ because this results in slower evolution of $g/{\alpha}$. The capture domain is infinite (no upper semi-major axis) in the limit $\tau_{\alpha}\rightarrow\tau_{d}$.

Beyond the upper limit in semi-major axis, the resonance will not be captured, but the crossing will still induce an obliquity increase due to an impulsive encounter with the resonant separatrix. The magnitude of this depends on the proximity of the crossing to the adiabatic limit. These non-adiabatic cases will be investigated in the following section.

4 Analytical System Evolution

Here we consider an analytic perturbative model to describe the planet’s spin state evolution upon encountering the secular resonance. The planet’s spin evolution is best described in the non-inertial frame that precesses with the planet’s orbit plane. We begin by defining canonical action-angle variables, $X=\cos\epsilon$ and $\psi$, where $\epsilon$ is the obliquity and $\psi$ is the angle between the projection of the spin vector onto the orbit plane and the ascending node. When the orbit nodal precession is uniform ($\dot{\Omega}$ is constant), the Hamiltonian that governs the system is given by (Morbidelli, 2002)

[TABLE]

Here, $\alpha$ is the spin axis precessional constant, $e$ is the eccentricity, and $i$ is the inclination of the orbit with respect to the invariable plane. The Hamiltonian of this system is well-studied; its origin, phase space structure, and equilibrium points have been characterized in many previous works (e.g. Colombo, 1966; Peale, 1969, 1974; Ward, 1975; Henrard & Murigande, 1987).

To begin reducing the Hamiltonian to a time-independent, integrable form, we first apply a transformation of variables by introducing the canonically conjugate pair,

[TABLE]

Hamiltonian 19 now takes the form,

[TABLE]

Here we have also introduced a dummy action conjugate to time, $\mathcal{T}$, such that the Hamiltonian formally represents a two degree of freedom system.

Next, we remove the direct time-dependence by performing a canonical transformation arising from the following type-2 generating function:

[TABLE]

The transformation equations yield the following new canonically conjugate variables (Morbidelli, 2002):

[TABLE]

It is clear from direct substitution of the above equations into Hamiltonian 21 that the Hamiltonian is independent of $\xi$, such that $\Xi$ is a constant of motion and can be dropped. By applying the additional simplifications that $e=0$ and $g=-\dot{\Omega}$, we write the autonomous, one degree-of-freedom Hamiltonian in the form,

[TABLE]

An understanding of the dynamics encapsulated by Hamiltonian 24 may be obtained by examining its level curves in phase space. The spin vector’s motion is confined to reside on these level curves. Figure 4 presents phase space portraits for three values of the ratio $g/{\alpha}$, which decreases as the disk decays. The equilibrium points (extrema of the Hamiltonian) are advected as $g/{\alpha}$ evolves. These equilibria are called “Cassini states” and are special configurations in which the spin vector remains fixed in the reference frame that precesses with the planet’s orbit (Colombo, 1966; Peale, 1969; Correia, 2015). At a critical value of $g/{\alpha}$ that depends on the orbital inclination (Fabrycky et al., 2007) and which we define as $(g/{\alpha})_{\mathrm{crit}}$, a separatrix emerges that partitions the phase space area and produces additional equilibria. (See middle panel of Figure 4.)

If the evolution of $g/{\alpha}$ is sufficiently slow compared to the dynamical timescale of the resonant motion, then adiabatic invariance dictates that the phase space area enclosed by a level curve trajectory is a conserved quantity, that is, as long as the orbit does not encounter a separatrix. Accordingly, consider a spin state that, at large $g/{\alpha}$, is prograde with a small obliquity. We illustrate this case with an accentuated purple trajectory in Figure 4. The spin vector initially circulates on a trajectory with a small phase space area. As $g/{\alpha}$ decreases and passes through $(g/{\alpha})_{\mathrm{crit}}$, the conservation of the phase space area implies that the spin vector enters the resonant domain and starts librating about the equilibrium on a small, banana-shaped level curve. This resonant equilibrium is “Cassini state 2”, and it corresponds to a configuration in which the spin vector and orbital angular momentum vector coprecess on opposite sides of the vector normal to the invariable plane. Following resonant capture, the planet’s obliquity is forced to larger and larger values as $g/{\alpha}$ decreases. Specifically, the condition $g\approx\alpha\cos\epsilon$ is upheld when the orbital inclination is small. Moreover, any source of dissipation, such as tides, will tend to force the system toward the equilibria.

4.1 Outcome of resonant encounter

With a perturbative Hamiltonian now in hand, we use it to efficiently explore the outcome of resonant encounter for different system configurations. Formally speaking, Hamiltonian 24 is integrable only if $\alpha$ and $g$ are constant. If, however, $\alpha$ and $g$ evolve much slower than the libration period of $\theta$, then successive phase-space portraits entailed by integrable Hamiltonian 24 provide an excellent approximation to the real dynamics. Thus, it is instructive to solve the evolution of the system by an application of Hamilton’s equations. We do so numerically but first apply a final canonical transformation to convenient Cartesian coordinates,

[TABLE]

The Hamiltonian in these coordinates is

[TABLE]

There are two principal parameters that govern whether the secular resonance will be crossed and captured: the semi-major axis, $a$, and disk decay timescale, $\tau_{\scriptscriptstyle d}$. We reduce the parameter space of potential system configurations to these two parameters and fix all others to the fiducial values assigned in Section 2 as well as the following: $M_{\star}=M_{\odot}$, $M_{d}(t=0)=0.01\ M_{\sun}$, $R_{o}=100\ \mathrm{AU}$, $\eta=0.001$, $i=5^{\circ}$, $f_{\alpha}=1$, and $f_{\omega}=0.1$ (Bryan et al., 2018; Batygin, 2018).

We keep the spin axis precession rate fixed ($\alpha=\mathrm{const}$, $f_{\tau}=0$), such that all evolution of the system is due to disk dissipation. The protoplanetary disk mass is set to evolve as a simple exponential decay,

[TABLE]

While this is certainly a simplification of the system’s evolution, it is sufficient for providing a sense of the behavior, which only depends on the evolution of the ratio $g/{\alpha}$. We note that in general, the spin rate, planetary radius, and possibly semi-major axis can also change during this time. Ultimately, however, $g/{\alpha}$ will cross through unity from above, and this is encapsulated by our model. With this parameterization in place, we apply Hamilton’s equations to Hamiltonian 26 and integrate numerically with a conventional Runge-Kutta ODE solver.

Five example system evolutions are displayed in Figure 5. They use the same parameters except for the semi-major axes. As shown in the second panel, all examples except for the $a=0.1$ AU case sweep through the resonance at some point during the disk decay. For the $a=0.3$, $0.5$, and $1.0$ AU cases, the resonance crossing is adiabatic and the planetary obliquities get excited all the way to $90^{\circ}$. The resonance crossings in the $2$ AU and $3$ AU cases are too fast for capture, so the planetary obliquities only receive impulsive, limited excitations.

Individual examples like these are helpful, but it is more informative to summarize the obliquity evolution over a wide range of system parameters. Figure 6 shows the final obliquity following resonant encounters for a grid of $a$ and $\tau_{\scriptscriptstyle d}$. The results are in excellent agreement with the analytic critical semi-major axis limits for resonant crossing and capture that we derived with equations 14 and 18. Below the upper limit, the obliquity is excited all the way to $90^{\circ}$ if the resonance is crossed. The heatmap also shows the magnitude of the obliquity excitation that occurs when the resonance is crossed too quickly for capture. This information is difficult to access through analytic considerations alone. Though it is possible to excite the obliquity without resonant capture, the excitation is not particularly significant unless it is relatively close to the adiabatic limit. However, excitations to $\sim 10-20^{\circ}$ still occur far past the limit.

4.2 Limiting effects from perturbing planets

In our model thus far, we have considered the interactions between a single planet, its host star and the protoplanetary disk in which it forms. In actuality, planets usually do not form in isolation, but rather in the presence of companion planets. Here we examine how the disk-induced secular spin-orbit resonance mechanism is altered as it pertains to multiple-planet systems.

For definitiveness, consider the simplified example where the perturbing planet’s orbit is taken to coincide with the plane of the disk. In this $i^{\prime}=0$ case, the frequency of a planet’s nodal recession, $g=-\dot{\Omega}$, never reaches zero. Rather, as the disk density diminishes, $g$ reaches a floor set by the secular frequencies of nodal recession induced by planet-planet perturbations. In the low-$i$, low-$e$ limit, these frequencies may be calculated using Laplace-Lagrange secular theory (Murray & Dermott, 1999). Specifically, the rate of nodal regression forced by a planetary companion reads:

[TABLE]

where $\alpha=a/a^{\prime}$, $\bar{\alpha}=\alpha$ if the perturbation is external and is unity otherwise, while the constant, $b_{3/2}^{(1)}(\alpha)$, is a Laplace coefficient defined by

[TABLE]

Although beyond the scope of our discussion, we note that when $i^{\prime}\neq 0$, the scenario is marginally more complicated due to multiple frequency components that are present (those associated with both the planet-disk interactions and the planet-planet interactions), implying that obliquity excitation in that case is not necessarily limited.

The floor on the nodal recession frequency can in some cases prevent resonant obliquity excitation, depending on the comparison between the planet’s spin axis precession frequency and the secular frequencies of nodal recession. If $g_{pp}>\alpha$, the nodal recession never reaches a commensurability with the spin axis precession, and obliquity excitation is blocked. Alternatively, if $g_{pp}<\alpha$, then limited excitation up to $\epsilon\approx\cos^{-1}(g_{pp}/{\alpha})$ will occur before the floor is reached. As will be shown below, it is the suppression of the outlined obliquity excitation mechanism by planet-planet interactions that renders this process inapplicable to the Solar System’s ice giants.

5 Application to Uranus and Neptune

Our analysis has indicated that disk-induced secular spin-orbit resonance is less likely for planets at distances greater than $\sim 5$ AU. This can be seen from the $a\sim 2$ AU fiducial upper limit at which the crossing is adiabatic and capture is ensured (see for example equation 18). Are there conditions for which the mechanism is relevant beyond $a\gtrsim 10$ AU? Figure 3 illustrates that planets with large $a$ are affected in the limit of significantly enhanced $\alpha$ (large $f_{\alpha}$), slow protoplanetary disk dissipation timescale $\tau_{d}$, and $f_{\tau}=\tau_{d}/{\tau_{\alpha}}$ close to unity. With all of these factors working together, the mechanism might be relevant to the ice giants, Uranus and Neptune.

The origins of Uranus’s $98^{\circ}$ and Neptune’s $30^{\circ}$ obliquities are not known with certainty. Their spin precession rates are much slower than any secular eigenfrequencies (Boué & Laskar, 2006; Nesvorný, 2018), implying that secular spin-orbit resonances induced by planet-planet interactions are not the cause of their large obliquities unless the spin precession was significantly enhanced by a massive satellite or circumplanetary disk (Boué & Laskar, 2010), or from the planets originally being closer to the Sun (Rogoszinski & Hamilton, 2018).

Collisional tilting from roughly Earth-mass giant impacts is the most well-accepted theory for the tilt of Uranus and Neptune (Safronov, 1966; Korycansky et al., 1990; Slattery et al., 1992; Morbidelli et al., 2012; Kegerreis et al., 2018). To explain Uranus’s prograde equatorial satellite system, Morbidelli et al. (2012) showed that the obliquity prior to the last tilting event must have been non-negligible, leading them to suggest that a multitude of such impacts were required. A particularly strong impact would evaporate ice from the ejected debris (Mousis, 2004), making the satellites devoid of ice, in contradiction with their $\sim 50\%$ water ice abundance (Szulágyi et al., 2018). Apart from collisions, the obliquities of Uranus and Neptune may also have been altered by some early process that twisted the Solar System’s total angular momentum vector (Tremaine, 1991).

Disk-induced secular spin-orbit resonance is yet another option to consider. Here we outline the idealized scenario under which it could have been relevant. First, resonant capture would have required that $\alpha$ was enhanced by the presence of a circumplanetary disk. Szulágyi et al. (2018) used three-dimensional radiative hydrodynamic simulations to show that the ice giants likely did form gaseous circumplanetary disks during their accretion. The disk masses in their simulations were $m_{cp}\sim 10^{-3}M_{p}$. Moreover, the disks were capable of forming icy satellites matching the current-day moon system of Uranus, which itself matches the $\sim 10^{-4}M_{p}$ mass ratio of Jupiter’s and Saturn’s satellite systems (Canup & Ward, 2006). We therefore assume the existence of a circumplanetary disk in the calculations that follow in section 5.1.

5.1 Idealized conditions for capture

We begin with the condition 18 for resonant capture determined by the adiabatic criterion. Using known planetary parameters for Uranus and Neptune, we can use this to find the values of $f_{\omega}$, $f_{\alpha}$, $f_{\tau}$, and $\tau_{d}$ required for resonant capture. We use present-day values of the semi-major axes and densities, Love numbers equal to $k_{2U}=0.104$ and $k_{2N}=0.127$ (Lainey, 2016), and moment of inertia factors $C_{U}=C_{N}=0.225$. It is likely that these parameters were different early on, but because we have little constraint on their histories, it is best to use present-day values for this plausibility argument. Finally, we use inclinations, $i_{U}=1^{\circ}$ and $i_{N}=2^{\circ}$.

Figure 7 shows the resulting limits on the resonant capture domain. For a line with a given $f_{\tau}$, the space above the line represents the $f_{\omega}f_{\alpha}$ and $\tau_{d}$ required to guarantee capture by the adiabatic constraint. We also plot the maximum values of $f_{\omega}f_{\alpha}$ by assuming $f_{\omega}=1$, considering a circumplanetary disk with $m_{cp}/M_{p}=10^{-3}$, and applying equations 2, 4, and 6 to calculate $f_{\alpha}$. We use the largest $r_{cp}/R_{p}$ that is still below the limit for adiabatic gravitational coupling (equation 9).

Clearly, capture is difficult to attain. The resonant domain corresponds to large $f_{\omega}f_{\alpha}$, implying significant $\alpha$-enhancement from a circumplanetary disk or satellites. It also requires long timescales for $\tau_{d}$ unless $f_{\tau}$ is nearly unity. $f_{\tau}\sim 1$ can only occur if $\alpha$ is evolving at a similar rate as $g$. It turns out, however, that this is a good assumption, which we will now show.

Early on in the system’s evolution, $\rho_{p}$, $k_{2}$, $C$, and $f_{\omega}$ all evolve non-trivially, but $m_{cp}$ changes the most. Accordingly, let us consider the simplification that $m_{cp}$ is the only changing quantity. Then $\tau_{\alpha}$ becomes

[TABLE]

Equations 6 show that $q/{\dot{q}}=l/{\dot{l}}=m_{cp}/{\dot{m}_{cp}}$. If $r_{cp}/R_{p}\sim 100$, we can take the limits $q\gg k_{2}{f_{\omega}}^{2}/3$ and $l\ll C$. Equation 30 is then approximately

[TABLE]

We simplify further by assuming that m_{cp}\mathrel{\raisebox{3.01385pt}{\mbox{\underset{\vbox{\hbox{$\sim$}}}{\vbox{\hbox{$\propto$}}}}}}M_{d}. While this approximation is theoretically-motivated as a result of the scaling between the circumstellar and circumplanetary disk surface densities (Ward & Canup, 2010; Szulágyi et al., 2018), in reality both $m_{cp}$ and $M_{d}$ evolve non-trivially. The m_{cp}\mathrel{\raisebox{3.01385pt}{\mbox{\underset{\vbox{\hbox{$\sim$}}}{\vbox{\hbox{$\propto$}}}}}}M_{d} simplification implies however that $m_{cp}/\dot{m}_{cp}\approx M_{d}/\dot{M}_{d}$ and finally that $\tau_{\alpha}\approx\tau_{d}$ and $f_{\tau}\approx 1$. Therefore, this first-order argument suggests that it is reasonable that $f_{\tau}$ is close to unity, which greatly alleviates the constraints posed by Figure 7.

In addition to the adiabatic constraint, another limitation is the length of time required to tilt the planet while the protoplanetary disk is still present. The resonant tilting scenario is only plausible if it can be completed sufficiently quickly. Once the resonance is captured, $g\approx\alpha\cos\epsilon$ holds if the inclination is small, and the obliquity evolves as

[TABLE]

where we’ve assumed constant $\tau_{d}$ and $f_{\tau}$. Figure 8 shows $\epsilon(t)$ as a function of $t/{\tau_{d}}$ for several values of $f_{\tau}$. Obliquity excitations on the order of $\epsilon\sim 10^{\circ}-30^{\circ}$ – which notably includes the obliquity of Neptune – are easily attainable within the disk lifetime (a few times $\tau_{d}$) for a range of $f_{\tau}$’s. On the other hand, tilting a planet from near $0^{\circ}$ to near $90^{\circ}$ most likely cannot be accomplished while the disk is still present. This timescale constraint is also prohibitive in the mean-motion/spin precession resonant tilting scenario for Uranus (Quillen et al., 2018). The tilting timescale might be shortened, however, if $\tau_{d}$ and $f_{\tau}$ are not constant.

5.2 Suppression due to giant planet perturbations

The preceding discussion has outlined the fine-tuned conditions required for disk-induced obliquity excitation to be relevant for Uranus and Neptune. It is likely, however, that this scenario as envisioned cannot have occurred due to the suppression mechanism we discussed in section 4.2, which in this case involves secular interactions among the giant planets.

In the $i^{\prime}=0$ case where the perturbing planets (i.e. Jupiter and Saturn) are aligned with the disk, Uranus’s and Neptune’s nodal recession frequencies would have decreased during disk dispersal but reached a floor at the secular frequencies associated with interactions among the four giant planets. Obliquity excitation would not have been possible because the secular frequencies – shown in Table 1 – are faster than the spin axis precession constants of Uranus and Neptune. Specifically, when $\omega=\omega_{b}$ (as opposed to the present-day spin rates), $\alpha_{0}=0.014$ arcsec yr*-1* for Uranus $\alpha_{0}=0.0040$ arcsec yr*-1* for Neptune. Even when enhanced by circumplanetary disks up to a factor of 100-1000, the spin axis precession rates are not faster than the fastest secular frequency. Therefore, if $i^{\prime}=0$, the resonance capture mechanism would have been inhibited.

6 Summary and Discussion

Secular spin-orbit resonance is a robust mechanism for producing non-zero planetary axial tilts. While previous work has mainly focused on resonances involving orbital precession induced by planet-planet interactions, here we addressed primordial orbital precession driven by the gravitational field of a protoplanetary disk. Unlike planet-planet interactions, disk-driven orbital precession is transient. However, a wide range of precession frequencies are swept during the disk decay process, thereby producing many opportunities for obliquity-exciting commensurabilities between the planet’s spin and orbital precession frequencies.

We find that planets with semi-major axes, $0.3\ \mathrm{AU}\lesssim a\lesssim 2\ \mathrm{AU}$, can readily be tipped to $90^{\circ}$ if the interaction is adiabatic. The details of these limits depend on the exact parameter configuration, and the upper limit on $a$ increases significantly (up to $\sim 20\ \mathrm{AU}$) if the disk dissipates slowly or if the planet is accompanied by an evolving circumplanetary disk or massive satellite(s). Moreover, beyond the upper limit in $a$, the resonance will still be encountered impulsively and will excite the obliquity to a degree that depends on the proximity to adiabaticity. We employed a perturbative Hamiltonian model to show that the obliquity may be excited up to $\sim 10^{\circ}-20^{\circ}$ even far past the adiabatic limit. However, planet-planet perturbations can in some cases prevent obliquity excitation by providing a floor on the nodal recession frequency.

Uranus’s $98^{\circ}$ and Neptune’s $30^{\circ}$ obliquities are typically ascribed to giant impacts. Smooth particle hydrodynamic simulations show that a single impactor on Uranus would necessarily have been quite massive ($\gtrsim 2\ M_{\oplus}$, Kegerreis et al. 2018, 2019). If it was too large, however, the impact would have evaporated ices from the ejected debris (Mousis, 2004; Szulágyi et al., 2018), which presents tension with the icy composition of Uranus’s prograde equatorial satellites. Moreover, Uranus’s satellite configuration could have only formed if the pre-impact obliquity was non-negligible, comparable with that of Neptune’s present-day $30^{\circ}$ tilt (Morbidelli et al., 2012). This suggests an additional impact or some other tilting mechanism.

Disk-induced secular spin-orbit resonance may present a viable alternative to collisional tilting in an idealized setting. Uranus and Neptune likely formed gaseous circumplanetary disks during their formation (Szulágyi et al., 2018). These disks would have increased the rate of spin axis precession and caused the spin and orbital precession frequencies to evolve on similar timescales, allowing a potential resonant encounter to last far longer and the possibility of obliquity excitation up to $\sim 30^{\circ}$. This excitation would have been inhibited, however, due to secular interactions with Jupiter and Saturn if they were aligned with the disk.

If disk-induced spin-orbit resonance might have played a role in Uranus’s and/or Neptune’s elevated obliquities, the question is raised for the rest of the Solar System planets: were they affected? One important consideration is that the mechanism requires non-zero orbital inclinations. The inclinations of Jupiter and Saturn at the time of resonance crossing may have been much smaller than those of Uranus and Neptune, such that adiabatic excitation was not a possibility. Moreover, the primordial obliquities of the gas giants would have been perturbed after disk dispersal by their planet interaction-induced spin-orbit resonances (Ward & Hamilton, 2004; Hamilton & Ward, 2004; Ward & Canup, 2006). As for the terrestrial planets, they likely finished forming through collisional growth after nebular gas dispersal (see Izidoro & Raymond 2018 for a review), and such collisions – as well as chaotic obliquity variations arising from planet interactions (Laskar & Robutel, 1993) – would have altered their primordial obliquities.

Extrasolar planets are prevalent in the $0.3\ \mathrm{AU}\lesssim a\lesssim 2\ \mathrm{AU}$ range where the mechanism is most relevant. We therefore predict that it is important for many exoplanets, although their excited obliquities might be perturbed by subsequent giant impacts. Close-in planets on the short end of the semi-major axis range will also be susceptible to tidal interactions that gradually right any planetary tilt produced by the resonant encounter. During this damping process, however, spin-orbit resonances induced by planet-planet interactions can be captured. The disk-induced resonance thus plays an important role for close-in systems by increasing the likelihood that the planets later find a separate resonance that, even in the presence of tides, captures the obliquity in a long-term excited state.

A similar staging role can be played by spin-orbit resonance induced by the oblateness of a rapidly-rotating young star (Millholland & Laughlin, 2019). Apart from the disk, early on in the system’s lifetime there is an additional component to each planet’s orbit nodal recession from the quadrupole field of the young star (Batygin & Adams, 2013; Spalding & Batygin, 2017). Like the disk decay, this can produce obliquity-exciting resonant encounters as the star contracts and spins down. It is most relevant for close-in planets, since the nodal recession frequency is a steep inverse function of $a$.

The disk-induced secular spin-orbit resonance is a powerful mechanism for exciting non-zero planetary obliquities with many areas ripe for further exploration. For instance, it would be useful to investigate how warps and turbulence in the disk affects the first-order picture we have developed here. Gaps and rings in the disk, which we now know are quite common (e.g. ALMA Partnership et al., 2015; Andrews et al., 2016; Isella et al., 2016; Avenhaus et al., 2018; Andrews et al., 2018), may also impact the resonant encounters and therefore lend themselves as excellent topics for future study.

7 Acknowledgements

We wish to thank Chris Spalding, Alessandro Morbidelli, and Greg Laughlin for their helpful insights and comments. We also thank the referee for a useful report. S.M. is supported by the NSF Graduate Research Fellowship Program under Grant DGE-1122492. K.B. is grateful to the David and Lucile Packard Foundation and the Alfred P. Sloan Foundation for their support.