THE GRAVITATIONAL INTERACTION BETWEEN PLANETS ON INCLINED ORBITS AND PROTOPLANETARY DISKS AS THE ORIGIN OF PRIMORDIAL SPIN–ORBIT MISALIGNMENTS

A major open question in the study of exoplanets is the origin of their apparent obliquity properties—the distribution of the angle λ between the stellar spin and the planet's orbital angular momentum vectors as projected on the sky (see, e.g., the review by Winn & Fabrycky 2015). Measurements of the Rossiter–McLaughlin effect in hot Jupiters (HJs, defined here as planets with masses of ${M}_{{\rm{p}}}\gtrsim 0.3\,{M}_{{\rm{J}}}$ that have orbital periods of ${P}_{\mathrm{orb}}\lesssim 10\,$ days) have indicated that λ spans the entire range from 0° to 180°, in stark contrast with the situation in the solar system (where the angle between the planets' total angular momentum vector and that of the Sun is only ∼6°). In addition, there is a marked difference in the distribution of λ between G stars, where ∼1/2 of systems are well aligned ($\lambda \lt 20^\circ$) and the rest are spread out roughly uniformly over the remainder of the λ range, and F stars of effective temperature ${T}_{\mathrm{eff}}\gtrsim 6250\,{\rm{K}}$, which exhibit only a weak excess of well-aligned systems. There is, however, also evidence for a dependence of the obliquity distribution on the properties of the planets and not just on those of the host star; in particular, only planets with ${M}_{{\rm{p}}}\lt 3\,{M}_{{\rm{J}}}$ have apparent retrograde orbits ($\lambda \gt 90^\circ$).

Various explanations have been proposed to account for the broad range of observed obliquities, but the inferred dependences on T_eff and M_p provide strong constraints on a viable model. In one scenario (Winn et al. 2010; Albrecht et al. 2012), HJs arrive in the vicinity of the host star on a misaligned orbit and subsequently act to realign the host through a tidal interaction, which is more effective in cool stars than in hot ones. In this picture, HJs form at large radii and either migrate inward through their natal disk while maintaining nearly circular orbits or are placed on a high-eccentricity orbit after the gaseous disk dissipates—which enables them to approach the center and become tidally trapped by the star (with their orbits getting circularized by tidal friction; e.g., Ford & Rasio 2006).¹ The processes that initiate high-eccentricity migration (HEM), which can be either planet–planet scattering (e.g., Chatterjee et al. 2008; Jurić & Tremaine 2008; Beaugé & Nesvorný 2012) or secular interactions that involve a stellar binary companion or one or more planetary companions (such as Kozai–Lidov oscillations—e.g., Wu & Murray 2003; Fabrycky & Tremaine 2007; Naoz et al. 2011; Petrovich 2015b—and secular chaos—e.g., Wu & Lithwick 2011; Lithwick & Wu 2014; Petrovich 2015a; Hamers et al. 2017), all give rise to HJs with a distribution of misaligned orbits. In the case of classical disk migration, the observed obliquities can be attributed to a primordial misalignment of the natal disk that occurred during its initial assembly from a turbulent interstellar gas (e.g., Bate et al. 2010; Fielding et al. 2015) or as a result of magnetic and/or gravitational torques induced, respectively, by a tilted stellar dipolar field and a misaligned companion (e.g., Lai et al. 2011; Batygin 2012; Batygin & Adams 2013; Lai 2014; Spalding & Batygin 2014).

The tidal realignment hypothesis that underlies the above modeling framework was challenged by the results of Mazeh et al. (2015), who examined the rotational photometric modulations of a large number of Kepler sources. Their analysis indicated that the common occurrence of aligned systems around cool stars characterizes the general population of planets and not just HJs, and, moreover, that this property extends to orbital periods as long as ∼50 days, about an order of magnitude larger than the maximum value of P_orb for which tidal interaction with the star remains important. To reconcile this finding with the above scenario, Matsakos & Königl (2015) appealed to the results of planet formation and evolution models, which predict that giant planets form efficiently in protoplanetary disks and that most of them migrate rapidly to the disk's inner edge, where, if the arriving planet's mass is not too high ($\lesssim 1\,{M}_{{\rm{J}}}$), it could remain stranded near that radius for up to ∼1 Gyr—until it gets tidally ingested by the host star. They proposed that the ingestion of a stranded HJ (SHJ)—which is accompanied by the transfer of its orbital angular momentum to the star—is the dominant spin-realignment mechanism. In this picture, the dichotomy in the obliquity properties between cool and hot stars is a direct consequence of the higher efficiency of magnetic braking and lower moment of inertia of the former in comparison with the latter. By applying a simple dynamical model to the observed HJ distributions in G and F stars, Matsakos & Königl (2015) inferred that ∼50% of planetary systems harbor an SHJ with a typical mass of ∼$0.6\,{M}_{{\rm{J}}}$. In this picture, the obliquity properties of currently observed HJs—and the fact that they are consistent with those of lower-mass and more distant planets—are most naturally explained if most of the planets in a given system—including any SHJ that may have been present—are formed in, and migrate along the plane of, a primordially misaligned disk.² This interpretation is compatible with the properties of systems like Kepler-56, in which two close-in planets have $\lambda \approx 45^\circ$ and yet are nearly coplanar (Huber et al. 2013), and 55 Cnc, a coplanar five-planet system with $\lambda \approx 72^\circ$ (e.g., Kaib et al. 2011; Bourrier & Hébrard 2014).³ It is also consistent with the apparent lack of a correlation between the obliquity properties of observed HJs and the presence of a massive companion (e.g., Knutson et al. 2014; Ngo et al. 2015; Piskorz et al. 2015).

In this paper, we explore a variant of the primordial disk misalignment model first proposed by Batygin (2012), in which, instead of the tilting of the entire disk by a distant (∼500 au) stellar companion on an inclined orbit, we consider the gravitational torque exerted by a much closer (∼5 au) planetary companion on such an orbit, which acts to misalign only the inner region of the protoplanetary disk. This model is motivated by the inferences from radial velocity surveys and adaptive-optics imaging data (Bryan et al. 2016; see also Knutson et al. 2014) that ∼70% of planetary systems harboring a transiting HJ have a companion with mass in the range of 1–13 M_J and a semimajor axis in the range of 1–20 au, and that ∼50% of systems harboring one or two planets detected by the radial velocity method have a companion with mass in the range of 1–$20\,{M}_{{\rm{J}}}$ and semimajor axis in the range of 5–20 au. Further motivation is provided by the work of Li & Winn (2016), who re-examined the photometric data analyzed by Mazeh et al. (2015) and found indications that the good-alignment property of planets around cool stars does not hold for large orbital periods, with the obliquities of planets with ${P}_{\mathrm{orb}}\gtrsim {10}^{2}\,$ days appearing to tend toward a random distribution.

One possible origin for a giant planet on an inclined orbit with a semimajor axis a of a few astronomical units is planet–planet scattering in the natal disk. Current theories suggest that giant planets may form in tightly packed configurations that can become dynamically unstable and undergo orbit crossing (see, e.g., Davies et al. 2014 for a review). The instabilities start to develop before the gaseous disk component dissipates (e.g., Marzari et al. 2010; Matsumura et al. 2010), and it has been argued (Chatterjee et al. 2008) that the planet–planet scattering process may, in fact, peak before the disk is fully depleted of gas (see also Lega et al. 2013). A close encounter between two giant planets is likely to result in a collision if the ratio $({M}_{{\rm{p}}}/{M}_{* })(a/{R}_{{\rm{p}}})$ (the Safronov number) is $\lt 1$ (where M_* is the stellar mass and R_p is the planet's radius), and in a scattering if this ratio is $\gt 1$ (e.g., Ford & Rasio 2008). The scattering efficiency is thus maximized when a giant planet on a comparatively wide orbit is involved (see Petrovich et al. 2014). High inclinations might also be induced by resonant excitation in giant planets that become trapped in a mean-motion resonance through classical (Type II) disk migration (Thommes & Lissauer 2003; Libert & Tsiganis 2009), and this process could, moreover, provide an alternative pathway to planet–planet scattering (Libert & Tsiganis 2011). In these scenarios, the other giant planets that were originally present in the disk can be assumed to have either been ejected from the system in the course of their interaction with the remaining misaligned planet or else reached the star at some later time through disk migration. As we show in this paper, a planet on an inclined orbit can have a significant effect on the orientation of the disk region interior to its orbital radius when the mass of that region decreases to the point where the inner disk's angular momentum becomes comparable to that of the planet. For typical mass depletion rates in protoplanetary disks (e.g., Batygin & Adams 2013), this can be expected to happen when the system's age is ∼10⁶–${10}^{7}\,$ years, which is comparable to the estimated formation time of Jupiter-mass planets at $\gtrsim 5\,$ au. In the proposed scenario, a planet of mass ${M}_{{\rm{p}}}\gtrsim {M}_{{\rm{J}}}$ is placed on a high-inclination orbit at the time ${t}_{0}\gtrsim 1\,$ Myr that, on the one hand, is late enough for the disk mass interior to the planet's location to have decreased to a comparable value, but that, on the other hand, is early enough for the inner disk to retain sufficient mass after becoming misaligned to enforce the orbital misalignment of existing planets and/or form new planets in its reoriented orbital plane (including any Jupiter-mass planets destined to become an HJ or an SHJ).

The dynamical model adopted in this paper is informed by the smooth-particle-hydrodynamics simulations carried out by Xiang-Gruess & Papaloizou (2013). They considered the interaction between a massive (1–$6\,{M}_{{\rm{J}}}$) planet that is placed on an inclined, circular orbit of radius 5 au and a low-mass ($0.01\,{M}_{* }$) protoplanetary disk that extends to 25 au. A key finding of these simulations was that the disk develops a warped structure, with the regions interior and exterior to the planet's radial location behaving as separate, rigid disks with distinct inclinations; in particular, the inner disk was found to exhibit substantial misalignment with respect to its initial direction when the planet's mass was large enough and its initial inclination was intermediate between the limits of 0° and 90° at which no torque is exerted on the disk. Motivated by these results, we construct an analytic model for the gravitational interaction between the planet and the two separate parts of the disk. The general effect of an interaction of this type between a planet on an inclined orbit and a rigid disk is to induce a precession of the planet's orbit about the total angular momentum vector. In contrast with Xiang-Gruess & Papaloizou (2013), whose simulations only extended over a fraction of a precession period, we consider the long-term evolution of such systems. In particular, we use our analytic model to study how the ongoing depletion of the disk's mass affects the orbital orientations of the planet and of the disk's two parts. We describe the model in Section 2 and present our calculations in Section 3. We discuss the implications of these results to planet obliquity measurements and to the alignment properties of debris disks in Section 4 and summarize in Section 5.

2.1. Assumptions

The initial configuration that we adopt is sketched in Figure 1. We consider a young star (subscript s) that is surrounded by a Keplerian accretion disk, and a Jupiter-mass planet (subscript p) on a circular orbit. The disk consists of two parts: an inner disk (subscript d) that extends between an inner radius ${r}_{{\rm{d}},\mathrm{in}}$ and an outer radius ${r}_{{\rm{d}},\mathrm{out}}$, and an outer disk (subscript h) that extends between ${r}_{{\rm{h}},\mathrm{in}}$ and ${r}_{{\rm{h}},\mathrm{out}};$ they are separated by a narrow gap that is centered on the planet's orbital radius a. The two parts of the disk are initially coplanar, with their normals aligned with the stellar angular momentum vector ${\boldsymbol{S}}$, whereas the planet's orbital angular momentum vector ${\boldsymbol{P}}$ is initially inclined at an angle ${\psi }_{{\rm{p}}0}$ with respect to ${\boldsymbol{S}}$ (where the subscript 0 denotes the time $t={t}_{0}$ at which the planet is placed on the inclined orbit).

Zoom In Zoom Out Reset image size

We assume that, during the subsequent evolution, each part of the disk maintains a flat geometry and precesses as a rigid body. The rigidity approximation is commonly adopted in this context and is attributed to efficient communication across the disk through the propagation of bending waves or the action of a viscous stress (e.g., Larwood et al. 1996; see also Lai 2014 and references therein).⁴ Based on the simulation results presented in Xiang-Gruess & Papaloizou (2013), we conjecture that this communication is severed at the location of the planet. This outcome is evidently the result of the planet opening up a gap in the disk, though it appears that the gap need not be fully evacuated for this process to be effective. In fact, the most strongly warped simulated disk configurations correspond to comparatively high initial inclination angles, for which the planet spends a relatively small fraction of the orbital time inside the disk, resulting in gaps that are less deep and wide than in the fully embedded case. Our calculations indicate that, during the disk's subsequent evolution, its inner and outer parts may actually detach as a result of the precessional oscillation of the inner disk. This oscillation is particularly strong in the case of highly mass-depleted disks on which we focus attention in this paper: in the example shown in Figure 6 below, the initial amplitude of this oscillation is ∼40°.

The planet's orbital inclination is subject to damping by dynamical friction (Xiang-Gruess & Papaloizou 2013), though the damping rate is likely low for the high values of ${\psi }_{{\rm{p}}0}$ that are of particular interest to us (Bitsch et al. 2013). Furthermore, in cases where the precessional oscillation of the inner disk causes the disk to split at the orbital radius of the planet, one can plausibly expect the local gas density to become too low for dynamical friction to continue to play a significant role on timescales longer than the initial oscillation period (∼10⁴ years for the example shown in Figure 6). In light of these considerations, and in the interest of simplicity, we do not include the effects of dynamical friction in any of our presented models.

As a further simplification, we assume that the planet's orbit remains circular. The initial orbital eccentricity of a planet ejected from the disk by either of the two mechanisms mentioned in Section 1 may well have a nonnegligible eccentricity. However, the simulations performed by Bitsch et al. (2013) indicate that the dynamical friction process damps eccentricities much faster than inclinations, so that the orbit can potentially be circularized on a timescale that is shorter than the precession time (i.e., before the two parts of the disk can become fully separated). On the other hand, even if the initial eccentricity is zero, it may be pumped up by the planet's gravitational interaction with the outer disk if ${\psi }_{{\rm{p}}0}$ is high enough ($\gtrsim 20^\circ ;$ Teyssandier et al. 2013). This is essentially the Kozai–Lidov effect, wherein the eccentricity undergoes periodic oscillations in antiphase with the orbital inclination (Terquem & Ajmia 2010). These oscillations were noticed in the numerical simulations of Xiang-Gruess & Papaloizou (2013) and Bitsch et al. (2013). Their period can be approximated by ${\tau }_{\mathrm{KL}}\sim {({r}_{{\rm{h}},\mathrm{out}}/{r}_{{\rm{h}},\mathrm{in}})}^{2}(2\pi /| {{\rm{\Omega }}}_{\mathrm{ph}}| )$ (Terquem & Ajmia 2010), where we used the expression for the precession frequency ${{\rm{\Omega }}}_{\mathrm{ph}}$ (Equation (32)) that corresponds to the torque exerted by the outer disk on the misaligned planet. For the parameters of the representative mass-depleted disk model shown in Figure 6, ${\tau }_{\mathrm{KL}}\sim {10}^{6}\,\mathrm{years}$. This time is longer by a factor of ∼10² than the initial precession period of the inner disk in this example, implying that the Kozai–Lidov process will have little effect on the high-amplitude oscillations of ${\psi }_{{\rm{p}}}$. Kozai–Lidov oscillations might, however, modify the details of the long-term behavior of the inner disk, since ${\tau }_{\mathrm{KL}}$ is comparable to the mass depletion time τ (Equation (10)) that underlies the secular evolution of the system.

Our model takes into account the tidal interaction of the spinning star with the inner and outer disks and with the planet, which was not considered in the aforementioned simulations. The inclusion of this interaction is motivated by the finding (Batygin & Adams 2013; Lai 2014; Spalding & Batygin 2014) that an evolving protoplanetary disk with a binary companion on an inclined orbit can experience a resonance between the disk precession frequency (driven by the companion) and the stellar precession frequency (driven by the disk), and that this resonance crossing can generate a strong misalignment between the angular momentum vectors of the disk and the star. As it turns out (see Section 3), in the case that we consider—in which the companion is a Jupiter-mass planet with an orbital radius of a few astronomical units rather than a solar-mass star at a distance of a few hundred astronomical units—this resonance is not encountered. We also show that, even in the case of a binary companion, the misalignment effect associated with the resonance crossing is weaker than that inferred in the above works when one also takes into account the torque that the star exerts on the inner disk (see Appendix C).

2.2. Equations

We model the dynamics of the system by following the temporal evolution of the angular momenta (${\boldsymbol{S}}$, ${\boldsymbol{D}}$, ${\boldsymbol{P}}$, and ${\boldsymbol{H}}$) of the four constituents (the star, the inner disk, the planet, and the outer disk, respectively) due to their mutual gravitational torques. Given that the orbital period of the planet is much shorter than the characteristic precession timescales of the system, we approximate the planet as a ring of uniform density, with a total mass equal to that of the planet and a radius equal to its semimajor axis.

The evolution of the angular momentum ${{\boldsymbol{L}}}_{k}$ of an object k under the influence of a torque ${{\boldsymbol{T}}}_{{ik}}$ exerted by an object i is given by $d{{\boldsymbol{L}}}_{k}/{dt}={{\boldsymbol{T}}}_{{ik}}$. The set of equations that describes the temporal evolution of the four angular momenta is thus

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{S}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{ds}}+{{\boldsymbol{T}}}_{\mathrm{ps}}+{{\boldsymbol{T}}}_{\mathrm{hs}},\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{D}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{sd}}+{{\boldsymbol{T}}}_{\mathrm{pd}}+{{\boldsymbol{T}}}_{\mathrm{hd}},\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{P}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{sp}}+{{\boldsymbol{T}}}_{\mathrm{dp}}+{{\boldsymbol{T}}}_{\mathrm{hp}},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{H}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{sh}}+{{\boldsymbol{T}}}_{\mathrm{dh}}+{{\boldsymbol{T}}}_{\mathrm{ph}},\end{eqnarray} \tag{ 4 }$$

where ${{\boldsymbol{T}}}_{{ik}}=-{{\boldsymbol{T}}}_{{ki}}$. The above equations can also be expressed in terms of the precession frequencies ${{\rm{\Omega }}}_{{ik}}$:

$$\begin{eqnarray}&&\displaystyle \frac{d{{\boldsymbol{L}}}_{k}}{{dt}}=\displaystyle \sum _{i}\,{{\boldsymbol{T}}}_{{ik}}=\displaystyle \sum _{i}\,{{\rm{\Omega }}}_{{ik}}\displaystyle \frac{{{\boldsymbol{L}}}_{i}\times {{\boldsymbol{L}}}_{k}}{{J}_{{ik}}},\end{eqnarray} \tag{ 5 }$$

where ${J}_{{ik}}=| {{\boldsymbol{L}}}_{i}+{{\boldsymbol{L}}}_{k}| ={({L}_{i}^{2}+{L}_{k}^{2}+2{L}_{i}{L}_{k}\cos {\theta }_{{ik}})}^{1/2}$ and ${{\rm{\Omega }}}_{{ik}}={{\rm{\Omega }}}_{{ki}}$. In Appendix A, we derive analytic expressions for the torques ${{\boldsymbol{T}}}_{{ik}}$ and the corresponding precession frequencies ${{\rm{\Omega }}}_{{ik}}$.

2.3. Numerical Setup

The host is assumed to be a protostar of mass ${M}_{* }={M}_{\odot }$, radius ${R}_{* }=2\,{R}_{\odot }$, rotation rate ${{\rm{\Omega }}}_{* }=0.1{({{GM}}_{* }/{R}_{* }^{3})}^{1/2}$, and angular momentum

$$\begin{eqnarray}S & = & {k}_{* }{M}_{* }{R}_{* }^{2}{{\rm{\Omega }}}_{* }=1.71\times {10}^{50}\\ & & \times \left(\displaystyle \frac{{k}_{* }}{0.2}\right)\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right){\left(\displaystyle \frac{{R}_{* }}{2{R}_{\odot }}\right)}^{2}\left(\displaystyle \frac{{{\rm{\Omega }}}_{* }}{0.1\sqrt{{{GM}}_{\odot }/{(2{R}_{\odot })}^{3}}}\right)\,\mathrm{erg}\,{\rm{s}},\end{eqnarray} \tag{ 6 }$$

where ${k}_{* }\simeq 0.2$ for a fully convective star (modeled as a polytrope of index n = 1.5). The planet is taken to have Jupiter's mass and radius, ${M}_{{\rm{p}}}={M}_{{\rm{J}}}$ and ${R}_{{\rm{p}}}={R}_{{\rm{J}}}$, and a fixed semimajor axis, a = 5 au, so that its orbital angular momentum is

$$\begin{eqnarray}P & = & {M}_{{\rm{p}}}{({{GM}}_{* }a)}^{1/2}=1.89\times {10}^{50}\\ & & \times \left(\displaystyle \frac{{M}_{{\rm{p}}}}{{M}_{{\rm{J}}}}\right){\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{1/2}{\left(\displaystyle \frac{a}{5\mathrm{au}}\right)}^{1/2}\,\mathrm{erg}\,{\rm{s}}.\end{eqnarray} \tag{ 7 }$$

We consider two values for the total initial disk mass: (1) ${M}_{{\rm{t}}0}=0.1\,{M}_{* }$, corresponding to a comparatively massive disk, and (2) ${M}_{{\rm{t}}0}=0.02\,{M}_{* }$, corresponding to a highly evolved system that has entered the transition-disk phase. In both cases, we take the disk surface density to scale with radius as ${r}^{-1}$. The inner disk extends from ${r}_{{\rm{d}},\mathrm{in}}=4\,{R}_{\odot }$ to ${r}_{{\rm{d}},\mathrm{out}}=a$, and initially has 10% of the total mass. Its angular momentum is

$$\begin{eqnarray}D & = & \displaystyle \frac{2}{3}{M}_{{\rm{d}}}{({{GM}}_{* })}^{1/2}\displaystyle \frac{{r}_{{\rm{d}},\mathrm{out}}^{3/2}-{r}_{{\rm{d}},\mathrm{in}}^{3/2}}{{r}_{{\rm{d}},\mathrm{out}}-{r}_{{\rm{d}},\mathrm{in}}}\\ & \simeq & 1.32\times {10}^{51}\left(\displaystyle \frac{{M}_{{\rm{d}}}}{0.01\,{M}_{\odot }}\right){\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{1/2}{\left(\displaystyle \frac{a}{5\mathrm{au}}\right)}^{1/2}\,\mathrm{erg}\,{\rm{s}}.\end{eqnarray} \tag{ 8 }$$

The outer disk has edges at ${r}_{{\rm{h}},\mathrm{in}}=a$ and ${r}_{{\rm{h}},\mathrm{out}}=50\,\mathrm{au}$, and angular momentum

$$\begin{eqnarray}H & = & \displaystyle \frac{2}{3}{M}_{{\rm{h}}}{({{GM}}_{* })}^{1/2}\displaystyle \frac{{r}_{{\rm{h}},\mathrm{out}}^{3/2}-{r}_{{\rm{h}},\mathrm{in}}^{3/2}}{{r}_{{\rm{h}},\mathrm{out}}-{r}_{{\rm{h}},\mathrm{in}}}\\ & \simeq & 3.76\times {10}^{52}\left(\displaystyle \frac{{M}_{{\rm{h}}}}{0.09\,{M}_{\odot }}\right){\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{1/2}{\left(\displaystyle \frac{{r}_{{\rm{h}},\mathrm{out}}}{50\mathrm{au}}\right)}^{1/2}\,\mathrm{erg}\,{\rm{s}}.\end{eqnarray} \tag{ 9 }$$

We model mass depletion in the disk using the expression first employed in this context by Batygin & Adams (2013),

$$\begin{eqnarray}&&{M}_{{\rm{t}}}(t)=\displaystyle \frac{{M}_{{\rm{t}}}(t=0)}{1+t/\tau },\end{eqnarray} \tag{ 10 }$$

where we adopt ${M}_{{\rm{t}}}(t=0)=0.1\,{M}_{\odot }$ and $\tau =0.5\,\mathrm{Myr}$ as in Lai (2014). We assume that this expression can also be applied separately to the inner and outer parts of the disk. The time evolution of the inner disk's angular momentum due to mass depletion is thus given by

$$\begin{eqnarray}&&{\left(\displaystyle \frac{d{\boldsymbol{D}}}{{dt}}\right)}_{\mathrm{depl}}=-\displaystyle \frac{{D}_{0}}{\tau {(1+t/\tau )}^{2}}\hat{{\boldsymbol{D}}}=-\displaystyle \frac{{\boldsymbol{D}}}{\tau +t}.\end{eqnarray} \tag{ 11 }$$

For the outer disk, we assume that the presence of the planet inhibits efficient mass accretion, and we consider the following limits: (1) the outer disk's mass remains constant, and (2) the outer disk loses mass (e.g., through photoevaporation) at the rate given by Equation (10).⁵ We assume that any angular momentum lost by the disk is transported out of the system (for example, by a disk wind).

We adopt a Cartesian coordinate system ($x,y,z$) as the "lab" frame of reference (see Figure 1). Initially, the equatorial plane of the star and the planes of the inner and outer disks coincide with the x–y plane (i.e., ${\psi }_{{\rm{s}}0}={\psi }_{{\rm{d}}0}={\psi }_{{\rm{h}}0}=0$, where ${\psi }_{k}$ denotes the angle between ${{\boldsymbol{L}}}_{k}$ and the z axis), and only the orbital plane of the planet has a finite initial inclination (${\psi }_{{\rm{p}}0}$). The x axis is chosen to coincide with the initial line of nodes of the planet's orbital plane.

Table 1 presents the models we explore and summarizes the relevant parameters. Specifically, column 1 contains the models' designations (with the letters M and m denoting, respectively, high and low disk masses at time $t={t}_{0}$), columns 2–5 indicate which system components are being considered, columns 6–9 list the disk and planet masses (with the arrow indicating active mass depletion), and columns 10 and 11 give the planet's semimajor axis and initial misalignment angle, respectively. The last listed model (binary) does not correspond to a planet misaligning the inner disk but rather to a binary star tilting the entire disk. This case is considered for comparison with the corresponding model in Lai (2014).

Table 1.  Model Parameters

Model	${\boldsymbol{S}}$	${\boldsymbol{D}}$	${\boldsymbol{P}}$	${\boldsymbol{H}}$	${M}_{{\rm{d}}0}\ [{M}_{* }]$	${M}_{{\rm{h}}0}\ [{M}_{* }]$	${M}_{{\rm{t}}0}\ [{M}_{* }]$	Mp	a [au]	${\psi }_{{\rm{p}}0}\ [^\circ ]$
DP-M	⋯	√	√	⋯	$0.010\downarrow$	⋯	⋯	MJ	5	60
DP-m	⋯	√	√	⋯	$0.002\downarrow$	⋯	⋯	MJ	5	60
all-M	√	√	√	√	$0.010\downarrow$	$0.090\downarrow$	0.10	MJ	5	60
all-m	√	√	√	√	$0.002\downarrow$	$0.018\downarrow$	0.02	MJ	5	60
all-Mx	√	√	√	√	$0.010\downarrow$	0.090–	0.10	MJ	5	60
all-mx	√	√	√	√	$0.002\downarrow$	0.018–	0.02	MJ	5	60
retrograde	√	√	√	√	$0.002\downarrow$	$0.018\downarrow$	0.02	MJ	5	110
binary	√	√	√	⋯	⋯	⋯	$0.10\downarrow$	${M}_{\odot }$	300	10

Download table as:  ASCIITypeset image

The gravitational interactions among the different components of the system that we consider (star, inner disk, planet, and outer disk) can result in a highly nonlinear behavior. To gain insight into these interactions, we start by analyzing a much simpler system, one consisting only of the inner disk and the (initially misaligned) planet. The relevant timescales that characterize the evolution of this system are the precession period ${\tau }_{\mathrm{dp}}\equiv 2\pi /{{\rm{\Omega }}}_{\mathrm{dp}}$ (Equation (29)) and the mass depletion timescale $\tau =5\times {10}^{5}\,$ years (Equation (10)).

Figure 2 shows the evolution of such a system for the case (model DP-M) where a Jupiter-mass planet on a misaligned orbit (${\psi }_{{\rm{p}}0}=60^\circ$) torques an inner disk of initial mass ${M}_{{\rm{d}}0}=0.01\,{M}_{* }$ (corresponding to ${M}_{{\rm{t}}0}=0.1\,{M}_{* }$, i.e., to ${t}_{0}=0$ when ${M}_{* }={M}_{\odot };$ see Equation (10)). The top left panel exhibits the angles ${\psi }_{{\rm{d}}}$ and ${\psi }_{{\rm{p}}}$ (blue: inner disk; red: planet) as a function of time. In this and the subsequent figures, we show results for a total duration of 10 Myr. This is long enough in comparison with τ to capture the secular evolution of the system, which is driven by the mass depletion in the inner disk. To capture the details of the oscillatory behavior associated with the precession of the individual angular momentum vectors (${\boldsymbol{D}}$ and ${\boldsymbol{P}}$) about the total angular momentum vector ${{\boldsymbol{J}}}_{\mathrm{dp}}={\boldsymbol{D}}+{\boldsymbol{P}}$ (subscript j)—which takes place on the shorter timescale ${\tau }_{\mathrm{dp}}$ ($\simeq 9\times {10}^{3}\,\mathrm{years}$ at $t={t}_{0}$)—we display the initial 0.1 Myr in the top left panel using a higher time resolution and, in addition, show the projected trajectories of the unit vectors $\hat{{\boldsymbol{D}}}$, $\hat{{\boldsymbol{P}}}$, and ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$ in the x–y plane during this time interval in the top right panel. Given that $0.1\,\mathrm{Myr}\ll \tau$, the vectors $\hat{{\boldsymbol{D}}}$ and $\hat{{\boldsymbol{P}}}$ execute a circular motion about ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$ with virtually constant inclinations with respect to the latter vector (given by the angles ${\theta }_{\mathrm{jd}}$ and ${\theta }_{\mathrm{jp}}$, respectively), and the orientation of ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$ with respect to the z axis (given by the angle ${\psi }_{{\rm{j}}}$) also remains essentially unchanged. (The projection of ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$ on the x–y plane is displaced from the center along the y axis, reflecting the fact that the planet's initial line of nodes coincides with the x axis.) As the vectors $\hat{{\boldsymbol{D}}}$ and $\hat{{\boldsymbol{P}}}$ precess about ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$, the angles ${\psi }_{{\rm{d}}}$ and ${\psi }_{{\rm{p}}}$ oscillate in the ranges $| {\psi }_{{\rm{j}}}-{\theta }_{\mathrm{jd}}| \leqslant {\psi }_{{\rm{d}}}\leqslant {\psi }_{{\rm{j}}}+{\theta }_{\mathrm{jd}}$ and $| {\psi }_{{\rm{j}}}-{\theta }_{\mathrm{jp}}| \leqslant {\psi }_{{\rm{p}}}\leqslant {\psi }_{{\rm{j}}}+{\theta }_{\mathrm{jp}}$, respectively.

Zoom In Zoom Out Reset image size

A notable feature of the evolution of this system on a timescale $\gtrsim \tau$ is the increase in the angle ${\psi }_{{\rm{d}}}$ (blue line in the top left panel)—indicating progressive misalignment of the disk with respect to its initial orientation—as the magnitude of the angular momentum ${\boldsymbol{D}}$ decreases with the loss of mass from the disk (blue line in the bottom right panel). At the same time, the orbital plane of the planet (red line in the top left panel) tends toward alignment with ${{\boldsymbol{J}}}_{\mathrm{dp}}$. The magenta lines in the top left and bottom right panels indicate that the orientation of the vector ${{\boldsymbol{J}}}_{\mathrm{dp}}$ remains fixed even as its magnitude decreases (on a timescale $\gtrsim \tau$) on account of the decrease in the magnitude of ${\boldsymbol{D}}$. As we demonstrate analytically in Appendix B, the constancy of ${\psi }_{{\rm{j}}}$ is a consequence of the inequality ${\tau }_{\mathrm{dp}}\ll \tau$.

To better understand the evolution of the disk and planet orientations, we consider the (small) variations in ${\boldsymbol{D}}$ and ${{\boldsymbol{J}}}_{\mathrm{dp}}$ that are induced by mass depletion over a small fraction of the precession period. On the left-hand side of Figure 3, we show a schematic sketch of the orientations of the vectors ${\boldsymbol{D}}$, ${\boldsymbol{P}}$, and ${{\boldsymbol{J}}}_{\mathrm{dp}}$ at some given time (denoted by the subscript 1) and a short time later (subscript 2). During that time interval, the vector ${{\boldsymbol{J}}}_{\mathrm{dp}}$ tilts slightly to the left, and as a result it moves away from ${\boldsymbol{D}}$ and closer to ${\boldsymbol{P}}$. The sketch on the right-hand side of Figure 3 demonstrates that, if we were to consider the same evolution a half-cycle later, the same conclusion would be reached: in this case, the vector ${{\boldsymbol{J}}}_{\mathrm{dp}3}$ moves slightly to the right (to become ${{\boldsymbol{J}}}_{\mathrm{dp}4}$), with the angle between ${{\boldsymbol{J}}}_{\mathrm{dp}}$ and ${\boldsymbol{D}}$ again increasing even as the angle between ${{\boldsymbol{J}}}_{\mathrm{dp}}$ and ${\boldsymbol{P}}$ decreases. The angles between the total angular momentum vector and the vectors ${\boldsymbol{D}}$ and ${\boldsymbol{P}}$ are thus seen to undergo a systematic, secular variation. The sketch in Figure 3 also indicates that the vector ${{\boldsymbol{J}}}_{\mathrm{dp}}$ undergoes an oscillation over each precession cycle. However, when ${\tau }_{\mathrm{dp}}\ll \tau$ and the fractional decrease in M_d over a precession period remains $\ll 1$, the amplitude of the oscillation is very small and ${{\boldsymbol{J}}}_{\mathrm{dp}}$ practically maintains its initial direction (see Appendix B for a formal demonstration of this result). In the limit where the disk mass becomes highly depleted and $D\to 0$, ${{\boldsymbol{J}}}_{\mathrm{dp}}\to {\boldsymbol{P}}$, i.e., the planet aligns with the initial direction of ${{\boldsymbol{J}}}_{\mathrm{dp}}$ (${\theta }_{\mathrm{jp}}\to 0$ and ${\psi }_{{\rm{p}}}\to {\psi }_{{\rm{j}}}$). The disk angular momentum vector then precesses about ${\boldsymbol{P}}$, with its orientation angle ${\psi }_{{\rm{d}}}$ (blue line in top left panel of Figure 2) oscillating between $| {\psi }_{{\rm{p}}}-{\theta }_{\mathrm{dp}}|$ and ${\psi }_{{\rm{p}}}+{\theta }_{\mathrm{dp}}$.⁶ Note that the precession frequency is also affected by the disk's mass depletion and decreases with time (see Equation (29)); the time evolution of ${{\rm{\Omega }}}_{\mathrm{dp}}$ is shown in the bottom left panel of Figure 2.

Zoom In Zoom Out Reset image size

Figure 4 shows the evolution of a similar system—model DP-m—in which the inner disk has a lower initial mass, ${M}_{{\rm{d}}0}=0.002\,{M}_{* }$ (corresponding to ${M}_{{\rm{t}}0}=0.02\,{M}_{* }$, i.e., to ${t}_{0}=2\,\mathrm{Myr}$ when ${M}_{* }={M}_{\odot };$ see Equation (10)). The initial oscillation frequency in this case is lower than in model DP-M, as expected from Equation (29), but it attains the same asymptotic value (bottom left panel), corresponding to the limit ${J}_{\mathrm{dp}}\to P$ in which ${{\rm{\Omega }}}_{\mathrm{dp}}$ becomes independent of M_d. The initial value of ${J}_{\mathrm{dp}}/D$ is higher in the present model than in the model considered in Figure 2 ($\simeq 1.5$ versus $\simeq 1.1;$ see Equations (7) and (8)), which results in a higher value of ${\psi }_{{\rm{j}}}$ (and, correspondingly, a higher initial value of ${\theta }_{\mathrm{jd}}$ and lower initial value of ${\theta }_{\mathrm{jp}}$). The higher value of ${\psi }_{{\rm{j}}}$ is the reason why the oscillation amplitude of ${\psi }_{{\rm{d}}}$ and the initial oscillation amplitude of ${\psi }_{{\rm{p}}}$ (top left panel) are larger in this case. The higher value of ${J}_{\mathrm{dp}}/{D}_{0}$ in Figure 4 also accounts for the differences in the projection map shown in the top right panel (a larger y value for the projection of ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$, a larger area encircled by the projection of $\hat{{\boldsymbol{D}}}$, and a smaller area encircled by the projection of $\hat{{\boldsymbol{P}}}$).

Zoom In Zoom Out Reset image size

We now consider the full system for two values of the total disk mass: ${M}_{{\rm{t}}0}=0.1\,{M}_{* }$ (model all-M, corresponding to ${t}_{0}=0;$ Figure 5) and ${M}_{{\rm{t}}0}=0.02\,{M}_{* }$ (model all-m, corresponding to ${t}_{0}=2\,\mathrm{Myr};$ Figure 6), assuming that both parts of the disk lose mass according to the relation given by Equation (10). The inner disks in these two cases correspond, respectively, to the disk masses adopted in model DP-M (Figure 2) and model DP-m (Figure 4). The merit of first considering the simpler systems described by the latter models becomes apparent from a comparison between the respective figures. It is seen that the basic behavior of model all-M is similar to that of model DP-M, and that the main differences between model all-M and model all-m are captured by the way in which model DP-m is distinct from model DP-M. The physical basis for this correspondence is the centrality of the torque exerted on the inner disk by the planet. According to Equation (5), the relative magnitudes of the torques acting on the disk at sufficiently late times (after D becomes smaller than the angular momentum of each of the other system components) are reflected in the magnitudes of the corresponding precession frequencies. The dominance of the planet's contribution can thus be inferred from the plots in the bottom left panels of Figures 5 and 6, which show that, after the contribution of D becomes unimportant (bottom right panels), the precession frequency induced by the planet exceeds those induced by the outer disk and by the star.⁷

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

While the basic disk misalignment mechanism is the same as in the planet–inner-disk system, the detailed behavior of the full system is understandably more complex. One difference that is apparent from a comparison of the left-hand panels in Figures 5 and 2 is the higher oscillation frequency of ${\psi }_{{\rm{p}}}$ and ${\psi }_{{\rm{d}}}$ in the full model (with the same frequency also seen in the timeline of ${\psi }_{{\rm{s}}}$). In this case, the planet–outer-disk precession frequency ${{\rm{\Omega }}}_{\mathrm{ph}}$ (Equation (32)) and the inner-disk–outer-disk precession frequency ${{\rm{\Omega }}}_{\mathrm{dh}}$ (Equation (31)) are initially comparable and larger than ${{\rm{\Omega }}}_{\mathrm{dp}}$, and ${{\rm{\Omega }}}_{\mathrm{ph}}$ remains the dominant frequency throughout the system's evolution. The fact that the outer disk imposes a precession on both ${\boldsymbol{P}}$ and ${\boldsymbol{D}}$ has the effect of weakening the interaction between the planet and the inner disk, which slows down the disk misalignment process. Another difference is revealed by a comparison of the top right panels: in the full system, ${\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$ precesses on account of the torque induced by the outer disk, so it no longer corresponds to just a single point in the x–y plane. This, in turn, increases the sizes of the regions traced in this plane by $\hat{{\boldsymbol{D}}}$ and $\hat{{\boldsymbol{P}}}$. The behavior of the lower-${M}_{{\rm{t}}0}$ model shown in Figure 6 is also more involved. In this case, in addition to the strong oscillations of the angles ${\psi }_{i}$ already manifested in Figure 4, the different precession frequencies ${{\rm{\Omega }}}_{{ik}}$ also exhibit large-amplitude oscillations, reflecting their dependence on the angles ${\theta }_{{ik}}$ between the angular momentum vectors. In both of the full-system models, the strongest influence on the star is produced by its interaction with the inner disk, but the resulting precession frequency (${{\rm{\Omega }}}_{\mathrm{sd}}$) remains low. Therefore, the stellar angular momentum vector essentially retains its original orientation, which implies that the angle ${\psi }_{{\rm{d}}}$ is a good proxy for the angle between the primordial stellar spin and the orbit of any planet that eventually forms in the inner disk.

We repeated the calculations shown in Figures 5 and 6 under the assumption that only the inner disk loses mass while M_h remains constant (models all-Mx and all-mx; Figures 7 and 8, respectively). At the start of the evolution, the frequencies ${{\rm{\Omega }}}_{\mathrm{ph}}$ and ${{\rm{\Omega }}}_{\mathrm{dh}}$ are ∝M_h, whereas ${{\rm{\Omega }}}_{\mathrm{dp}}$ scales linearly (or, in the case of the lower-${M}_{{\rm{d}}0}$ model, close to linearly) with M_d (see Appendix A). In the cases considered in Figures 5 and 6 all these frequencies decrease with time, so the relative magnitude of ${{\rm{\Omega }}}_{\mathrm{dp}}$ remains comparatively large throughout the evolution. In contrast, in the cases shown in Figures 7 and 8 the frequencies ${{\rm{\Omega }}}_{\mathrm{ph}}$ and ${{\rm{\Omega }}}_{\mathrm{dh}}$ remain constant and only ${{\rm{\Omega }}}_{\mathrm{dp}}$ decreases with time. As the difference between ${{\rm{\Omega }}}_{\mathrm{dp}}$ and the other two frequencies starts to grow, the inner disk misalignment process is aborted, and thereafter the mean values of ${\psi }_{{\rm{d}}}$ and ${\psi }_{{\rm{p}}}$ remain constant. This behavior is consistent with our conclusion about the central role that the torque exerted by the planet plays in misaligning the inner disk: when the fast precession that the outer disk induces in the orbital motions of both the planet and the inner disk comes to dominate the system dynamics, the direct coupling between the planet and the inner disk is effectively broken and the misalignment process is halted. Note, however, from Figure 8 that, even in this case, the angle ${\psi }_{{\rm{d}}}$ can attain a high value (as part of a large-amplitude oscillation) when ${M}_{{\rm{t}}0}$ is small.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To determine whether the proposed misalignment mechanism can also account for disks (and, eventually, planets) on retrograde orbits, we consider a system in which the companion planet is placed on such an orbit (model retrograde, which is the same as model all-m except that ${\psi }_{{\rm{p}}0}$ is changed from 60° to 110°). As Figure 9 demonstrates, the disk in this case evolves to a retrograde configuration (${\psi }_{{\rm{d}}}\gt 90^\circ$) at late times even as the planet's orbit reverts to prograde motion. A noteworthy feature of the plotted orbital evolution (shown in the high-resolution portion of the figure) is the rapid increase in the value of ${\psi }_{{\rm{d}}}$ (which is an adequate proxy for ${\theta }_{\mathrm{sd}}$ also in this case)—and corresponding fast decrease in the value of ${\psi }_{{\rm{p}}}$—that occurs when the planet's orbit transitions from a retrograde to a prograde orientation. This behavior can be traced to the fact that $\cos {\theta }_{\mathrm{ph}}$ vanishes at essentially the same time that ${\psi }_{{\rm{p}}}$ crosses 90° because the outer disk (which dominates the total angular momentum) remains well aligned with the z axis. This, in turn, implies (see Equation (32)) that, at the time of the retrograde-to-prograde transition, the planet becomes dynamically decoupled from the outer disk and only retains a coupling to the inner disk. Its evolution is, however, different from that of a "reduced" system, in which only the planet and the inner disk interact, because the inner disk remains dynamically "tethered" to the outer disk (${\theta }_{\mathrm{dh}}\ne 90^\circ$). As we verified by an explicit calculation, the evolution of the reduced system remains smooth when ${\psi }_{{\rm{p}}}$ crosses 90°. The jump in ${\psi }_{{\rm{p}}}$ exhibited by the full system leads to a significant increase in the value of $\cos {\theta }_{\mathrm{ph}}$ and hence of ${{\rm{\Omega }}}_{\mathrm{ph}}$, which, in turn, restores (and even enhances) the planet's coupling to the outer disk after its transition to retrograde motion (see bottom panel of Figure 9). The maximum value attained by ${\theta }_{\mathrm{sd}}$ in this example is $\simeq 172^\circ$, which, just as in the prograde case shown in Figure 6, exceeds the initial misalignment angle of the planetary orbit (albeit to a much larger extent in this case). It is, however, worth noting that not all model systems in which the planet is initially on a retrograde orbit give rise to a retrograde inner disk at the end of the prescribed evolution time; in particular, we found that the outcome of the simulated evolution (which depends on whether ${\psi }_{{\rm{p}}}$ drops below 90°) is sensitive to the value of the initial planetary misalignment angle ${\psi }_{{\rm{p}}0}$ (keeping all other model parameters unchanged).

Zoom In Zoom Out Reset image size

In concluding this section it is instructive to compare the results obtained for our model with those found for the model originally proposed by Batygin (2012; see Section 1 for references to additional work on that model). We introduced our proposed scenario as a variant of the latter model, with a close-by giant planet taking the place of a distant stellar companion. In the original proposal, the disk misalignment was attributed to the precessional motion that is induced by the torque that the binary companion exerts on the disk. In this picture the spin–orbit angle oscillates (on a timescale of ∼1 Myr for typical parameters) between 0° and roughly twice the binary orbital inclination, so it can be large if observed at the "right" time. Our model retains this feature of the earlier proposal, particularly in cases where the companion planet is placed on a high-inclination orbit after the disk has already lost much of its initial mass, but it also exhibits a novel feature that gives rise to a secular (rather than oscillatory) change in the spin–orbit angle (which can potentially lead to a substantial increase in this angle). This new behavior represents an "exchange of orientations" between the planet and the inner disk that is driven by the mass loss from the inner disk and corresponds to a decrease of the inner disk's angular momentum from a value higher than that of the planet to a lower value (with the two remaining within an order of magnitude of each other for representative parameters). This behavior is not found in a binary system because of the large mismatch between the angular momenta of the companion and the disk in that case (and, in fact, it is also suppressed in the case of a planetary companion when the mass of the outer disk is not depleted).

As we already noted in Section 2.1, Batygin & Adams (2013) suggested that the disk misalignment in a binary system can be significantly increased due to a resonance between the star–disk and binary–disk precession frequencies. (We can use Equations (26) and (29), respectively, to evaluate these frequencies, plugging in values for the outer disk radius, companion orbital radius, and companion mass that are appropriate for the binary case.) Lai (2014) clarified the effect of this resonance and emphasized that, for plausible system parameters, it can be expected to be crossed as the disk becomes depleted of mass. However, for the planetary-companion systems considered in this paper the ratio $| {{\rm{\Omega }}}_{\mathrm{sd}}/{{\rm{\Omega }}}_{\mathrm{dp}}|$ remains $\lt 1$ throughout the evolution, so no such resonance is encountered in this case. In both of these systems ${{\rm{\Omega }}}_{\mathrm{sd}}$ is initially $\propto {M}_{{\rm{d}}}$, so it decreases during the early evolution. The same scaling also characterizes ${{\rm{\Omega }}}_{\mathrm{dp}}$ in the planetary case, which explains why the corresponding curves do not cross. In contrast, in the binary case (for which the sum of the disk and companion angular momenta is dominated by the companion's contribution) the frequency ${{\rm{\Omega }}}_{\mathrm{dp}}$ does not scale with the disk mass and it thus remains nearly constant, which makes it possible for the corresponding curves to cross (see Figure 12 in Appendix C). Since our formalism also encompasses the binary case, we examined one such system (model binary)—using the parameters adopted in Figure 3 of Lai (2014)—for comparison with the results of that work. Our findings are presented in Appendix C.

The model considered in this paper represents a variant of the primordial disk misalignment scenario of Batygin (2012) in which the companion is a nearby planet rather than a distant star and only the inner region of the protoplanetary disk (interior to the planet's orbit) becomes inclined. In this section, we assess whether this model provides a viable framework for interpreting the relevant observations.

The first—and most basic—question that needs to be addressed is whether the proposed misalignment mechanism is compatible with the broad range of apparent spin–orbit angles indicated by the data. In Section 3, we showed that the spin–orbit angle ${\theta }_{\mathrm{sd}}$ can deviate from its initial value of 0° either because of the precessional motion that is induced by the planet's torque on the disk or on account of the secular variation that is driven by the mass depletion process. In the "reduced" disk–planet model considered in Figures 2 and 4, for which the angle ${\psi }_{{\rm{d}}}$ is taken as a proxy for the intrinsic spin–orbit angle, the latter mechanism increases ${\theta }_{\mathrm{sd}}$ to ∼45°–50° on a timescale of 10 Myr for an initial planetary inclination ${\psi }_{{\rm{p}}0}=60^\circ$. The maximum disk misalignment is, however, increased above this value by the precessional oscillation, whose amplitude is higher the lower the initial mass of the disk. Based on the heuristic discussion given in connection with Figure 3, the maximum possible value of ${\psi }_{{\rm{d}}}$ (corresponding to the limit ${J}_{\mathrm{dp}}\to P$) is given by

$$\begin{eqnarray}&&{\psi }_{{\rm{d}},\max }=\arccos \displaystyle \frac{{D}_{0}+P\cos {\psi }_{{\rm{p}}0}}{{({D}_{0}^{2}+{P}^{2}+2{D}_{0}P\cos {\psi }_{{\rm{p}}0})}^{1/2}}+{\psi }_{{\rm{p}}0}.\end{eqnarray} \tag{ 12 }$$

For the parameters of Figure 4, ${\psi }_{{\rm{d}},\max }\approx 84\buildrel{\circ}\over{.} 5$, which can be compared to the actual maximum value ($\simeq 72^\circ$) attained over the course of the 10 Myr evolution depicted in this figure.⁸ Although the behavior of the full system (which also includes the outer disk and the star) is more complicated, we found (see Figures 5 and 6) that, if the outer disk also loses mass, the maximum value attained by ${\theta }_{\mathrm{sd}}$ ($\simeq 67^\circ$) is not much smaller than in the simplified model. Note that in the original primordial-misalignment scenario the maximum value of ${\theta }_{\mathrm{sd}}$ ($\simeq 2\,{\psi }_{{\rm{p}}0}$) would have been considerably higher ($\simeq 120^\circ$) for the parameters employed in our example. However, as indicated by Equation (12), the maximum value predicted by our model depends on the ratio $P/{D}_{0}$ and can in principle exceed the binary-companion limit if D₀ is small and P is sufficiently large.⁹ Repeating the calculations shown in Figure 6 for higher values of M_p, we found that the maximum value of ${\theta }_{\mathrm{sd}}$ is ∼89°, 104°, and 125° when ${M}_{{\rm{p}}}/{M}_{{\rm{J}}}$ increases from 1 to 2, 3, and 4, respectively. These results further demonstrate that the disk can be tilted to a retrograde configuration even when ${\psi }_{{\rm{p}}0}\lt 90^\circ$ if the planet is sufficiently massive, although a retrograde disk orientation can also be attained (including in the case of ${M}_{{\rm{p}}}\lesssim {M}_{{\rm{J}}}$) if the planet's orbit is initially retrograde (see Figure 9). A low initial value of the disk angular momentum D arises naturally in the leading scenarios for placing planets in inclined orbits, which favor comparatively low disk masses (see Section 1). The distribution of ${\psi }_{{\rm{p}}0}$ as well as those of the occurrence rate, mass, and orbital radius of planets on inclined orbits are required for determining the predicted distribution of primordial inner-disk misalignment angles in this scenario, for comparison with observations.¹⁰ However, this information, as well as data on the relevant values of ${M}_{{\rm{d}}0}$, are not yet available, so our results for ${\theta }_{\mathrm{sd}}$ are only a first step (a proof of concept) toward validating this interpretation of the measured planet obliquities.

Our proposed misalignment mechanism is most effective when the disk mass within the planetary orbit drops to ∼M_p. In the example demonstrating this fact (Figure 6), ${M}_{{\rm{d}}0}\approx 2\,{M}_{{\rm{J}}}$. In the primordial disk misalignment scenario, ${M}_{{\rm{d}}0}$ includes the mass that would eventually be detected in the form of an HJ (or a lower-mass planet) moving around the central star on a misaligned orbit. Furthermore, if the ingestion of an HJ on a misaligned orbit is as ubiquitous as inferred in the SHJ picture, that mass, too, must be included in the tally. These requirements are consistent with the fact that the typical disk misalignment time in our model (a few megayears) is comparable to the expected giant-planet formation time, but this similarity also raises the question of whether the torque exerted by the initially misaligned planet has the same effect on the gaseous inner disk and on a giant planet embedded within it. This question was considered by several authors in the context of a binary companion (e.g., Xiang-Gruess & Papaloizou 2014; Picogna & Marzari 2015; Martin et al. 2016). A useful gauge of the outcome of this dynamical interaction is the ratio of the precession frequency induced in the embedded planet (which we label ${{\rm{\Omega }}}_{\mathrm{pp}}$) to ${{\rm{\Omega }}}_{\mathrm{dp}}$ (Picogna & Marzari 2015). We derive an expression for ${{\rm{\Omega }}}_{\mathrm{pp}}$ by approximating the inclined and embedded planets as two rings with radii a and ${a}_{1}\lt a$, respectively (see Appendix A), and evaluate ${{\rm{\Omega }}}_{\mathrm{dp}}$ under the assumption that the disk mass has been sufficiently depleted for the planetary contribution (P) to dominate J_dp. This leads to ${{\rm{\Omega }}}_{\mathrm{pp}}/{{\rm{\Omega }}}_{\mathrm{dp}}\simeq 2{({a}_{1}/{r}_{{\rm{d}},\mathrm{out}})}^{3/2}$, which is the same as the estimate obtained by Picogna & Marzari (2015) for a binary system. In the latter case, this ratio is small ($\lesssim 0.1$) for typical parameters, implying that the embedded planet cannot keep up with the disk precession and hence that its orbit develops a significant tilt with respect to the disk's plane. However, when the companion is a planet, the above ratio equals ${({a}_{1}/a)}^{3/2}$ and may be considerably larger ($\lesssim 1$), which suggests that the embedded planet can remain coupled to the disk in this case.

A key prediction of our proposed scenario—which distinguishes it from the original Batygin (2012) proposal—is that there would in general be a difference in the obliquity properties of "nearby" and "distant" planets, corresponding to the different orientations attained, respectively, by the inner and outer disks. This prediction is qualitatively consistent with the finding of Li & Winn (2016) that the good spin–orbit alignment inferred in cool stars from an analysis of rotational photometric modulations in Kepler sources (Mazeh et al. 2015) becomes weaker (with the inferred orientations possibly tending toward a nearly random distribution) at large orbital periods (${P}_{\mathrm{orb}}\gtrsim {10}^{2}\,$days). The interpretation of these results in our picture is that the outer planets remain aligned with the original stellar-spin direction, whereas the inner planets—and, according to the SHJ model, also the stellar spin in ∼50% of sources—assume the orientation of the misaligned inner disk (which samples a broad range of angles with respect to the initial spin direction). Further observations and analysis are required to corroborate and refine these findings so that they can be used to place tighter constrains on the models.

The result reported by Li & Winn (2016) is seemingly at odds with another set of observational findings—the discovery that the orbital planes of debris disks (on scales $\gtrsim {10}^{2}\,\mathrm{au}$) are by and large well aligned with the spin axis of the central star (Watson et al. 2011; Greaves et al. 2014). This inferred alignment also seemingly rules out any interpretation of the obliquity properties of exoplanets (including the SHJ model) that appeals to a tidal realignment of the host star by a misaligned HJ. These apparent difficulties can, however, be alleviated in the context of the SHJ scenario and our present model. Specifically, in the SHJ picture the realignment of the host star occurs on a relatively long timescale (≲1 Gyr; see Matsakos & Königl 2015). This is much longer than the lifetime (∼1–10 Myr) of the gaseous disk that gives rise to both the misaligned "nearby" planets and the debris disk (which, in the scenario considered in this paper, are associated with the inner and outer parts of the disk, respectively). The inferred alignment properties of debris disks can be understood in this picture if these disks are not much older than ∼1 Gyr, so that the stellar spin axis still points roughly along its original direction (which coincides with the symmetry axis of the outer disk). We searched the literature for age estimates of the 11 uniformly observed debris disks tabulated in Greaves et al. (2014) and found that only two (10 CVn and 61 Vir) are definitely much older than 1 Gyr. Now, Matsakos & Königl (2015) estimated that ∼50% of systems ingest an SHJ and should exhibit spin–orbit alignment to within 20°, with the rest remaining misaligned. Thus, the probability of observing an aligned debris disk in an older system is $\sim 1/2$, implying that the chance of detecting two out of two such systems is $\sim 1/4$. It is, however, worth noting that the two aforementioned systems may not actually be well aligned: based on the formal measurement uncertainties quoted in Greaves et al. (2014), the misalignment angle could be as large as 36° in 10 CVn and 31° in 61 Vir. Further measurements that target old systems might be able to test the proposed explanation, though one should bear in mind that additional factors may affect the observational findings. For example, in the tidal-downsizing scenario of planet formation, debris disks are less likely to exist around stars that host giant planets (see Fletcher & Nayakshin 2016).

In this paper, we conduct a proof-of-concept study of a variant of the primordial disk misalignment model of Batygin (2012). In that model, a binary companion with an orbital radius of a few hundred astronomical units exerts a gravitational torque on a protoplanetary disk that causes its plane to precess and leads to a large-amplitude oscillation of the spin–orbit angle ${\theta }_{\mathrm{sd}}$ (the angle between the angular momentum vectors of the disk and the central star). Motivated by recent observations, we explore an alternative model in which the role of the distant binary is taken by a giant planet with an orbital radius of just a few astronomical units. Such a companion likely resided originally in the disk, and its orbit most probably became inclined away from the disk's plane through a gravitational interaction with other planets (involving either scattering or resonant excitation).

Our model setup is guided by indications from numerical simulations (Xiang-Gruess & Papaloizou 2013) that, in the presence of the misaligned planet, the disk separates at the planet's orbital radius into inner and outer parts that exhibit distinct dynamical behaviors even as each can still be well approximated as a rigid body. We integrate the secular dynamical evolution equations in the quadrupole approximation for a system consisting of the inclined planet, the two disk parts, and the spinning star, with the disk assumed to undergo continuous mass depletion. We show that this model can give rise to a broad range of values for the angle between the angular momentum vectors of the inner disk and the star (including values of ${\theta }_{\mathrm{sd}}$ in excess of 90°), but that the orientation of the outer disk remains virtually unchanged. We demonstrate that the misalignment is induced by the torque that the planet exerts on the inner disk and that it is suppressed when the mass depletion time in the outer disk is much longer than in the inner disk, so that the outer disk remains comparatively massive and the fast precession that it induces in the motions of the inner disk and the planet effectively breaks the dynamical coupling between the latter two. Our calculations reveal that the largest misalignments are attained when the initial disk mass is low (on the order of that of observed systems at the onset of the transition-disk phase). We argued that, when the misalignment angle is large, the inner and outer parts of the disk become fully detached and damping of the planet's orbital inclination by dynamical friction effectively ceases. This suggests a consistent primordial misalignment scenario: the inner region of a protoplanetary disk can be strongly misaligned by a giant planet on a high-inclination orbit if the disk's mass is low (i.e., late in the disk's evolution); in turn, the planet's orbital inclination is least susceptible to damping in a disk that undergoes a strong misalignment.

We find that, in addition to the precession-related oscillations seen in the binary-companion model, the spin–orbit angle also exhibits a secular growth in the planetary-companion case, corresponding to a monotonic increase in the angle between the inner disk's and the total (inner disk plus planet) angular momentum vectors (accompanied by a monotonic decrease in the angle between the planet's and the total angular momentum vectors). This behavior arises when the magnitude of the inner disk's angular momentum is initially comparable to that of the planet but drops below it as a result of mass depletion (on a timescale that is long in comparison with the precession period). This does not happen when the companion is a binary, since in that case the companion's angular momentum far exceeds that of the inner disk at all times. On the other hand, in the binary case the mass depletion process can drive the system to a resonance between the disk–planet and star–disk precession frequencies, which has the potential of significantly increasing the maximum value of ${\theta }_{\mathrm{sd}}$ (e.g., Batygin & Adams 2013; Lai 2014). We show that this resonance is not encountered when the companion is a nearby planet because—in contrast with the binary-companion case, in which the disk–binary precession frequency remains constant—both of these precession frequencies decrease with time in the planetary-companion case. However, we also show that when the torque that the star exerts on the disk is taken into account (and not just that exerted by the companion, as in previous treatments), the misalignment effect of the resonance crossing in the binary case is measurably weaker.

A key underlying assumption of the primordial disk misalignment model is that the planets embedded in the disk remain confined to its plane as the disk's orientation shifts, so that their orbits become misaligned to the same extent as that of the gaseous disk. However, the precession frequency that a binary companion induces in the disk can be significantly higher than the one induced by its direct interaction with an embedded planet, which would lead to the planet's orbital plane separating from that of the disk: this argument was used to critique the original version of the primordial misalignment model (e.g., Picogna & Marzari 2015). However, this potential difficulty is mitigated in the planetary-companion scenario, where the ratio of these two frequencies is typically substantially smaller.

The apparent difference in the obliquity properties of HJs around cool and hot stars can be attributed to the tidal realignment of a cool host star by an initially misaligned HJ (e.g., Albrecht et al. 2012). The finding (Mazeh et al. 2015) that this dichotomy is exhibited also by lower-mass planets and extends to orbital distances where tidal interactions with the star are very weak motivated the SHJ proposal (Matsakos & Königl 2015), which postulates that ∼50% of systems contain an HJ that arrives through migration in the protoplanetary disk and becomes stranded near its inner edge for a period of ≲1 Gyr—during which time the central star continues to lose angular momentum by magnetic braking—until the tidal interaction with the star finally causes it to be ingested (resulting in the transfer of the planet's orbital angular momentum to the star and in the realignment of the stellar spin in the case of cool stars). This picture fits naturally with the primordial misalignment model discussed in this paper. In this broader scenario, the alignment properties of currently observed planets (which do not include SHJs) can be explained if these planets largely remain confined to the plane of their primordial parent disk. In the case of cool stars, the planets exhibit strong alignment on account of the realignment action of a predecessor SHJ, whereas in the case of hot stars they exhibit a broad range of spin–orbit angles, reflecting the primordial range of disk misalignment angles that was preserved on account of the ineffectiveness of the tidal realignment process in these stars. A distinguishing prediction of the planetary-companion variant of the primordial misalignment model in the context of this scenario arises from the expected difference in the alignment properties of the inner and outer disks, which implies that the good alignment exhibited by planets around cool stars should give way to a broad range of apparent spin–orbit angles above a certain orbital period. There is already an observational indication of this trend (Li & Winn 2016), but additional data are needed to firm it up. A complementary prediction, which is potentially also testable, is that the range of obliquities exhibited by planets around hot stars would narrow toward $\lambda =0^\circ$ at large orbital periods. This scenario may also provide an explanation for another puzzling observational finding—that large-scale debris disks are by and large well aligned with the spin vector of the central star—which, on the face of it, seems inconsistent with the spin-realignment hypothesis. In this interpretation, debris disks are associated with the outer parts of protoplanetary disks and should therefore remain aligned with the central star—as a general rule for hot stars, but also in the case of cool hosts that harbor an SHJ, if they are observed before the SHJ realigns the star. This explanation is consistent with the fact that the great majority of observed debris disks have inferred ages of ≪1 Gyr, but the extent to which it addresses the above finding can be tested through its prediction that a sufficiently large sample of older systems should also contain misaligned disks.

We are grateful to Dan Fabrycky, Tsevi Mazeh, and Sean Mills for fruitful discussions. We also thank Gongjie Li and Josh Winn for helpful correspondence, and the referee for useful comments. This work was supported in part by NASA ATP grant NNX13AH56G and has made use of NASA's Astrophysics Data System Bibliographic Services and of matplotlib, an open-source plotting library for Python (Hunter 2007).

A.1. Torques

Figure 10 shows two concentric rings in a Cartesian coordinate system, oriented so that their mutual gravitational torques induce a rotation about the x axis. Because of the configuration's symmetry, the only nonzero component of the torque that ring i exerts on ring k is that along the x axis:

$$\begin{eqnarray}&&{[{{\boldsymbol{T}}}_{{ik}}]}_{x}=-{\int }_{k}{\int }_{i}{\left[{{\boldsymbol{y}}}_{k}\times \left(\displaystyle \frac{G{{dm}}_{k}{{dm}}_{i}}{{w}^{2}}\hat{{\boldsymbol{w}}}\right)\right]}_{x},\end{eqnarray} \tag{ 13 }$$

where

$$\begin{eqnarray}&&{{\boldsymbol{y}}}_{k}={r}_{k}\sin {\phi }_{k}\,\hat{{\boldsymbol{y}}},\end{eqnarray} \tag{ 14 }$$

$$\begin{eqnarray}&&\hat{{\boldsymbol{w}}}=\displaystyle \frac{{\boldsymbol{w}}}{w}=\displaystyle \frac{{{\boldsymbol{r}}}_{i}-{{\boldsymbol{r}}}_{k}}{w},\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}&&{{\boldsymbol{r}}}_{k}={r}_{k}\cos {\phi }_{k}\,\hat{{\boldsymbol{x}}}+{r}_{k}\sin {\phi }_{k}\,\hat{{\boldsymbol{y}}},\end{eqnarray} \tag{ 16 }$$

$$\begin{eqnarray}{{\boldsymbol{r}}}_{i} & = & {r}_{i}\cos {\phi }_{i}\,\hat{{\boldsymbol{x}}}+{r}_{i}\sin {\phi }_{i}\cos {\theta }_{{ik}}\,\hat{{\boldsymbol{y}}}\\ & & -{r}_{i}\sin {\phi }_{i}\sin {\theta }_{{ik}}\,\hat{{\boldsymbol{z}}},\end{eqnarray} \tag{ 17 }$$

$$\begin{eqnarray}w & = & [{r}_{k}^{2}+{r}_{i}^{2}-2{r}_{k}{r}_{i}(\cos {\phi }_{k}\cos {\phi }_{i}\\ & & +{\sin {\phi }_{k}\sin {\phi }_{i}\cos {\theta }_{{ik}})]}^{1/2},\end{eqnarray} \tag{ 18 }$$

$$\begin{eqnarray}&&{[{{\boldsymbol{y}}}_{k}\times {\boldsymbol{w}}]}_{x}={r}_{k}{r}_{i}\sin {\phi }_{k}\sin {\phi }_{i}\sin {\theta }_{{ik}}\,\hat{{\boldsymbol{x}}},\end{eqnarray} \tag{ 19 }$$

and ${\int }_{k},{\int }_{i}$ denote integrals over the masses m_k and m_i. These expressions can be readily generalized to a "continuum of rings"—i.e., a disk—with inner and outer radii of r_in and r_out, respectively. In the case of a ring ${dm}=\lambda r\,d\phi$, where $\lambda =m/2\pi r$ is the linear mass density, whereas in the case of a disk ${dm}={\rm{\Sigma }}r\,{drd}\phi$, where Σ is the surface density. Adopting ${\rm{\Sigma }}={{\rm{\Sigma }}}_{0}({r}_{0}/r)$ (as in Batygin 2012 and Lai 2014), where ${{\rm{\Sigma }}}_{0}$, r₀ are constants, gives $m=2\pi {{\rm{\Sigma }}}_{0}{r}_{0}({r}_{\mathrm{out}}-{r}_{\mathrm{in}})$. Therefore, dm can be expressed as

$$\begin{eqnarray}{dm}=\left\{\begin{array}{ll}\displaystyle \frac{m}{2\pi }\,d\phi & \quad {\rm{for}}\,{\rm{a}}\,{\rm{ring}},\\ \displaystyle \frac{m}{2\pi ({r}_{\mathrm{out}}-{r}_{\mathrm{in}})}\,{drd}\phi & \quad {\rm{for}}\,{\rm{a}}\,{\rm{disk}}.\end{array}\right.\end{eqnarray} \tag{ 20 }$$

Zoom In Zoom Out Reset image size

For ${r}_{i}\gg {r}_{k}$ one can approximate

$$\begin{eqnarray}\displaystyle \frac{1}{{w}^{3}} & \simeq & \displaystyle \frac{1}{{r}_{i}^{3}}+\displaystyle \frac{3{r}_{k}}{{r}_{i}^{4}}\cos {\phi }_{k}\cos {\phi }_{i}\\ & & +\displaystyle \frac{3{r}_{k}}{{r}_{i}^{4}}\sin {\phi }_{k}\sin {\phi }_{i}\cos {\theta }_{{ik}},\end{eqnarray} \tag{ 21 }$$

and thus the torque becomes

$$\begin{eqnarray}{[{T}_{{ik}}]}_{x} & \simeq & -\,{A}_{i}\,{B}_{k}\,\sin {\theta }_{{ik}}\cos {\theta }_{{ik}}\\ & = & -\,\left(3G{\displaystyle \int }_{i}\displaystyle \frac{{\sin }^{2}{\phi }_{i}}{{r}_{i}^{3}}\,{{dm}}_{i}\right)\\ & & \times \left({\displaystyle \int }_{k}{r}_{k}^{2}{\sin }^{2}{\phi }_{k}\,{{dm}}_{k}\right)\sin {\theta }_{{ik}}\cos {\theta }_{{ik}}\end{eqnarray} \tag{ 22 }$$

(the other terms integrate to zero), where

$$\begin{eqnarray}A=\left\{\begin{array}{ll}\displaystyle \frac{3{Gm}}{2{r}^{3}} & \quad \quad {\rm{for}}\,{\rm{a}}\,{\rm{ring}},\\ \displaystyle \frac{3{Gm}({r}_{\mathrm{out}}+{r}_{\mathrm{in}})}{4{r}_{\mathrm{out}}^{2}{r}_{\mathrm{in}}^{2}} & \quad \quad {\rm{for}}\,{\rm{a}}\,{\rm{disk}},\end{array}\right.\end{eqnarray} \tag{ 23 }$$

and

$$\begin{eqnarray}B=\left\{\begin{array}{ll}\displaystyle \frac{{{mr}}^{2}}{2} & \quad \quad {\rm{for}}\,{\rm{a}}\,{\rm{ring}},\\ \displaystyle \frac{m({r}_{\mathrm{out}}^{3}-{r}_{\mathrm{in}}^{3})}{6({r}_{\mathrm{out}}-{r}_{\mathrm{in}})} & \quad \quad {\rm{for}}\,{\rm{a}}\,{\rm{disk}}.\end{array}\right.\end{eqnarray} \tag{ 24 }$$

The torque that k exerts on i is ${[{T}_{{ki}}]}_{x}=-{[{T}_{{ik}}]}_{x}$.¹¹ Equation (22) can also be used when object k is a star by setting ${B}_{k}={k}_{q}{M}_{* }{R}_{* }^{2}{{\rm{\Omega }}}_{* }^{2}/({{GM}}_{* }/{R}_{* }^{3})$ and, in the case of a protostar, using ${k}_{q}\simeq 0.1$ (the value appropriate to fully convective stars; e.g., Lai 2014).

Figure 11 compares the torque that is calculated using Equation (13) (points) with the approximate expression of Equation (22) (lines) for the following systems: two rings (left-hand panel), a ring and a disk (middle panel), and two disks (right-hand panel). In all cases, the mass of each object is taken to be $m={M}_{{\rm{J}}}$, and we set ${\theta }_{{ik}}=30^\circ$. Each point in the figure (representing a superposed pair of + and × symbols) corresponds to a different system, characterized by its relevant parameters (the radius of the ring or the inner and outer radii of the disk). For the ring–ring system, ring 1 has radius ${r}_{1}=1\,\mathrm{au}$, and the different cases correspond to ${r}_{2}\in [0.1,10]$ au. For the ring–disk system, the ring has radius ${r}_{2}=2.8\,\mathrm{au}$, and the inner and outer edges of the disk lie in the ranges ${r}_{1\mathrm{in}}\in [0.013,133]$ au and ${r}_{1\mathrm{out}}\in [0.13,1,333]$ au, respectively. The same ranges are adopted for both disks in the case of the disk–disk system. The figure indicates that for the cases that are relevant to the present study—in particular, when a ring (representing a planet) is located at either the inner or the outer edge of the disk, or when the two disks are adjacent—Equation (22) provides a very good approximation to the torque.

Zoom In Zoom Out Reset image size

Our approach, which amounts to using only the lowest-order (quadrupole) term in the expansion of the interaction potential (e.g., Papaloizou & Terquem 1995), is less general than the Gaussian averaging method employed by Batygin (2012), but, as demonstrated in Figure 11, it is entirely adequate for our purposes. In our treatment of the disk, we make the further approximations that its constituent "rings" remain circular and that it behaves as a rigid body. In the case of a binary companion and ${M}_{{\rm{t}}}=0.01\,{M}_{* }$, Batygin (2012) verified by an explicit calculation that these approximations are well justified even if the binary moves on an eccentric orbit and self-gravity is the only mode of internal interaction in the disk.

A.2. Precession Frequencies

By combining Equations (5) and (22), we obtain an analytic expression for the precession frequencies:

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{{ik}}={{\rm{\Omega }}}_{{ki}}\simeq -{A}_{i}\,{B}_{k}\,\displaystyle \frac{{J}_{{ik}}}{{L}_{i}{L}_{k}}\cos {\theta }_{{ik}},\end{eqnarray} \tag{ 25 }$$

where, again, object i is taken to be "outside of" object k (${r}_{i}\gg {r}_{k}$). The six characteristic frequencies are

$$\begin{eqnarray}{{\rm{\Omega }}}_{\mathrm{sd}} & \simeq & -\,4.86\times {10}^{-5}\left(\displaystyle \frac{2{k}_{q}}{{k}_{* }}\right)\left(\displaystyle \frac{{M}_{{\rm{d}}}}{0.01{M}_{* }}\right)\\ & & \times \left(\displaystyle \frac{{R}_{* }}{2{R}_{\odot }}\right){\left(\displaystyle \frac{{r}_{{\rm{d}},\mathrm{in}}}{4{R}_{* }}\right)}^{-2}{\left(\displaystyle \frac{{r}_{{\rm{d}},\mathrm{out}}}{5\mathrm{au}}\right)}^{-1}\displaystyle \frac{{J}_{\mathrm{sd}}}{D}\\ & & \times \left(\displaystyle \frac{{{\rm{\Omega }}}_{* }}{0.1\sqrt{{{GM}}_{\odot }/{(2{R}_{\odot })}^{3}}}\right)\cos {\theta }_{\mathrm{sd}}\,\displaystyle \frac{2\pi }{\mathrm{years}},\end{eqnarray} \tag{ 26 }$$

$$\begin{eqnarray}{{\rm{\Omega }}}_{\mathrm{sp}} & \simeq & -\,5.14\times {10}^{-10}\left(\displaystyle \frac{2{k}_{q}}{{k}_{* }}\right)\left(\displaystyle \frac{{M}_{{\rm{p}}}}{{M}_{{\rm{J}}}}\right){\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{-1}\\ & & \times {\left(\displaystyle \frac{{R}_{* }}{2{R}_{\odot }}\right)}^{3}{\left(\displaystyle \frac{a}{5\mathrm{au}}\right)}^{-3}\displaystyle \frac{{J}_{\mathrm{sp}}}{P}\\ & & \times \left(\displaystyle \frac{{{\rm{\Omega }}}_{* }}{0.1\sqrt{{{GM}}_{\odot }/{(2{R}_{\odot })}^{3}}}\right)\cos {\theta }_{\mathrm{sp}}\,\displaystyle \frac{2\pi }{\mathrm{years}},\end{eqnarray} \tag{ 27 }$$

$$\begin{eqnarray}{{\rm{\Omega }}}_{\mathrm{sh}} & \simeq & -\,2.42\times {10}^{-9}\left(\displaystyle \frac{2{k}_{q}}{{k}_{* }}\right)\left(\displaystyle \frac{{M}_{{\rm{h}}}}{0.09{M}_{* }}\right)\\ & & \times {\left(\displaystyle \frac{{R}_{* }}{2{R}_{\odot }}\right)}^{3}{\left(\displaystyle \frac{{r}_{{\rm{h}},\mathrm{in}}}{5\mathrm{au}}\right)}^{-2}{\left(\displaystyle \frac{{r}_{{\rm{h}},\mathrm{out}}}{50\mathrm{au}}\right)}^{-1}\displaystyle \frac{{J}_{\mathrm{sh}}}{H}\\ & & \times \left(\displaystyle \frac{{{\rm{\Omega }}}_{* }}{0.1\sqrt{{{GM}}_{\odot }/{(2{R}_{\odot })}^{3}}}\right)\cos {\theta }_{\mathrm{sh}}\,\displaystyle \frac{2\pi }{\mathrm{years}},\end{eqnarray} \tag{ 28 }$$

$$\begin{eqnarray}{{\rm{\Omega }}}_{\mathrm{dp}} & \simeq & -\,2.23\times {10}^{-4}\left(\displaystyle \frac{{M}_{{\rm{d}}}}{0.01{M}_{* }}\right)\\ & & \times {\left(\displaystyle \frac{{r}_{{\rm{d}},\mathrm{out}}}{a}\right)}^{2}\displaystyle \frac{{J}_{\mathrm{dp}}}{D}\\ & & \times \left(\displaystyle \frac{{{\rm{\Omega }}}_{{\rm{p}}}}{\sqrt{{{GM}}_{\odot }/{(5\mathrm{au})}^{3}}}\right)\cos {\theta }_{\mathrm{dp}}\,\displaystyle \frac{2\pi }{\mathrm{years}}\end{eqnarray} \tag{ 29 }$$

$$\begin{eqnarray} & \simeq & -\,3.20\times {10}^{-5}\left(\displaystyle \frac{{M}_{{\rm{p}}}}{{M}_{{\rm{J}}}}\right){\left(\displaystyle \frac{{M}_{* }}{{M}_{\odot }}\right)}^{-1}\\ & & \times {\left(\displaystyle \frac{{r}_{{\rm{d}},\mathrm{out}}}{a}\right)}^{3/2}\displaystyle \frac{{J}_{\mathrm{dp}}}{P}\\ & & \times \left(\displaystyle \frac{{{\rm{\Omega }}}_{{\rm{p}}}}{\sqrt{{{GM}}_{\odot }/{(5\mathrm{au})}^{3}}}\right)\cos {\theta }_{\mathrm{dp}}\,\displaystyle \frac{2\pi }{\mathrm{years}},\end{eqnarray} \tag{ 30 }$$

$$\begin{eqnarray}{{\rm{\Omega }}}_{\mathrm{dh}} & \simeq & -\,1.51\times {10}^{-3}\left(\displaystyle \frac{{M}_{{\rm{h}}}}{0.09{M}_{* }}\right)\\ & & \times {\left(\displaystyle \frac{{r}_{{\rm{d}},\mathrm{out}}}{{r}_{{\rm{h}},\mathrm{in}}}\right)}^{2}\left(\displaystyle \frac{{r}_{{\rm{d}},\mathrm{out}}}{{r}_{{\rm{h}},\mathrm{out}}}\right)\displaystyle \frac{{J}_{\mathrm{dh}}}{H}\\ & & \times {\left(\displaystyle \frac{{{GM}}_{* }/{r}_{{\rm{d}},\mathrm{out}}^{3}}{{{GM}}_{\odot }/{(5\mathrm{au})}^{3}}\right)}^{1/2}\cos {\theta }_{\mathrm{dh}}\,\displaystyle \frac{2\pi }{\mathrm{years}},\end{eqnarray} \tag{ 31 }$$

and

$$\begin{eqnarray}{{\rm{\Omega }}}_{\mathrm{ph}} & \simeq & -\,3.02\times {10}^{-3}\left(\displaystyle \frac{{M}_{{\rm{h}}}}{0.09{M}_{* }}\right)\\ & & \times {\left(\displaystyle \frac{a}{{r}_{{\rm{h}},\mathrm{in}}}\right)}^{2}\left(\displaystyle \frac{a}{{r}_{{\rm{h}},\mathrm{out}}}\right)\displaystyle \frac{{J}_{\mathrm{ph}}}{H}\\ & & \times {\left(\displaystyle \frac{{{GM}}_{* }/{a}^{3}}{{{GM}}_{\odot }/{(5\mathrm{au})}^{3}}\right)}^{1/2}\cos {\theta }_{\mathrm{ph}}\,\displaystyle \frac{2\pi }{\mathrm{years}}.\end{eqnarray} \tag{ 32 }$$

To obtain an expression for the time evolution of ${{\boldsymbol{J}}}_{\mathrm{dp}}={\boldsymbol{D}}+{\boldsymbol{P}}$, we write

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{D}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{pd}}+{\left(\displaystyle \frac{d{\boldsymbol{D}}}{{dt}}\right)}_{\mathrm{depl}},\end{eqnarray} \tag{ 33 }$$

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{P}}}{{dt}}=-{{\boldsymbol{T}}}_{\mathrm{pd}},\end{eqnarray} \tag{ 34 }$$

and take their sum using Equation (11):

$$\begin{eqnarray}\displaystyle \frac{d{{\boldsymbol{J}}}_{\mathrm{dp}}}{{dt}} & = & {\left(\displaystyle \frac{{dD}}{{dt}}\right)}_{\mathrm{depl}}(\cos \phi ^{\prime} \sin \theta ^{\prime} \hat{x}^{\prime} +\sin \phi ^{\prime} \sin \theta ^{\prime} \hat{y}^{\prime} \\ & & +\cos \theta ^{\prime} \hat{z}^{\prime} ),\end{eqnarray} \tag{ 35 }$$

where we expressed $\hat{{\boldsymbol{D}}}$ in a cartesian coordinate system $(x^{\prime} ,y^{\prime} ,z^{\prime} )$ with $\hat{z}^{\prime} ={\hat{{\boldsymbol{J}}}}_{\mathrm{dp}}$ and with $\theta ^{\prime}$, $\phi ^{\prime}$ the spherical polar angles. Since the precession period is much shorter than the depletion time (for example, for the parameters that characterize model DP-M, the initial value of ${\tau }_{\mathrm{dp}}/\tau$ is $\simeq 0.017$), it is an excellent approximation to treat ${({dD}/{dt})}_{\mathrm{depl}}$ as a constant over one precession period. Averaging over ϕ therefore gives

$$\begin{eqnarray}&&\left\langle \displaystyle \frac{d{{\boldsymbol{J}}}_{\mathrm{dp}}}{{dt}}\right\rangle \simeq {\left(\displaystyle \frac{{dD}}{{dt}}\right)}_{\mathrm{depl}}{\hat{{\boldsymbol{J}}}}_{\mathrm{dp}},\end{eqnarray} \tag{ 36 }$$

where the angle brackets denote an average over a precession period. This implies that $\langle {{\boldsymbol{J}}}_{\mathrm{dp}}\rangle$ decreases in magnitude without changing its direction. The oscillation of ${\psi }_{{\rm{j}}}$ during a single precession period—described in Figure 3—is in practice so small (its amplitude is $\simeq {\tau }_{\mathrm{dp}}/\tau$ for $t\ll \tau$) that it cannot be picked out in Figures 2 and 4.

The formalism employed in this paper can also be used to treat the original variant of the primordial disk misalignment model, in which the companion is a distant star rather than a nearby giant planet. To validate our code, we consider one such system (model binary; see Table 1), which corresponds to the example presented in Figure 3 of Lai (2014). In that work, the evolution of a system consisting of a star, a binary star (subscript b), and a disk that undergoes mass depletion (according to the prescription given by Equation (10)) was studied by integrating the equations

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{S}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{ds}}\end{eqnarray} \tag{ 37 }$$

and

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{D}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{bd}}.\end{eqnarray} \tag{ 38 }$$

The results of solving these two equation with our numerical scheme are presented in the top two panels of Figure 12. These results are identical to those obtained by Lai (2014) and indicate that even a system with small initial misalignments between the star and the disk (${\theta }_{\mathrm{sd}}=5^\circ$) and between the disk and the binary (${\theta }_{\mathrm{db}}=5^\circ$, with ${\theta }_{\mathrm{sb}}=10^\circ$) can attain a large final spin–orbit angle if a resonance between the precessions frequencies ${{\rm{\Omega }}}_{\mathrm{sd}}$ and ${{\rm{\Omega }}}_{\mathrm{db}}$ is crossed. In this case, the precession frequency that the torque exerted by the disk induces in the stellar angular momentum vector is initially high enough (${{\rm{\Omega }}}_{\mathrm{sd}}\gt {{\rm{\Omega }}}_{\mathrm{db}}$) for the star–disk pair to remain coupled as the disk precesses under the influence of the binary. However, as the mass of the disk becomes depleted, ${{\rm{\Omega }}}_{\mathrm{sd}}$ decreases and eventually crosses ${{\rm{\Omega }}}_{\mathrm{db}}$. Beyond that point, the stellar angular momentum can no longer follow the precession of the disk's angular momentum, and the motion of these two vectors decouples. At resonance, the star–disk system may attain a large misalignment, which, in the absence of strong star–disk coupling, remains "frozen" during the ensuing evolution.

Zoom In Zoom Out Reset image size

In the model presented in Lai (2014) the torque that the star exerts on the disk is neglected. We now use the more general formulation employed in this work to extend that model by also including the torques exerted by the star. Thus, instead of Equations (37) and (38), we integrate

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{S}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{ds}}+{{\boldsymbol{T}}}_{\mathrm{bs}}\end{eqnarray} \tag{ 39 }$$

and

$$\begin{eqnarray}&&\displaystyle \frac{d{\boldsymbol{D}}}{{dt}}={{\boldsymbol{T}}}_{\mathrm{sd}}+{{\boldsymbol{T}}}_{\mathrm{bd}}.\end{eqnarray} \tag{ 40 }$$

In practice, only the torque that the star exerts on the disk plays a role, with the effect of ${{\boldsymbol{T}}}_{\mathrm{bs}}$ remaining negligible (see Equations (26) and (27)). The results of this integration are shown in the bottom two panels of Figure 12 and demonstrate that the back torque that the star exerts on the disk can significantly reduce the effectiveness of the resonance misalignment mechanism, so that its effect cannot in general be neglected.