Termination of an inward migration of a gap-opening planet triggered by dust feedback

TL;DR

This study uses hydrodynamic simulations to show that dust feedback can halt inward migration of gap-opening planets, especially Jupiter-sized ones, potentially explaining observed ring and gap structures in protoplanetary disks.

Contribution

It introduces the effect of dust feedback on planet migration, demonstrating its role in stopping inward migration in a way not previously considered.

Findings

Dust feedback reduces outer disk surface density.

Outer disk torque decreases due to dust feedback.

Inward migration can be halted or slowed significantly.

Abstract

The planet migration due to the disk--planet interaction is one of the most important processes to determine the architecture of planetary systems. A sufficiently massive planet forms a density gap and migrates together with the gap. By carrying out two-dimensional and two-fluid (gas and dust grains) hydrodynamic simulations, we investigated the effects of the dust feedback on the migration of the gap-opening planet, which was not considered in previous studies. We found that the gas surface density at the outer edge of the gap becomes smaller due to the dust feedback, and thus the torque exerted from the outer disk decreases. This mechanism becomes effective as the gap becomes wider and deeper. In particular, when the mass of the planet is Jupiter-size and turbulent viscosity is $\alpha = 3\times 10^{-4}$, the planet can migrate outward due to the reduction of the torque exerted from the outer disk. Even for a smaller planet, the migration becomes significantly slow down. This termination of the inward migration triggered by the dust feedback may explain why ring and gap structures can be frequently observed within the protoplanetary disks.

Equations

$$\frac{1}{a_{\mathrm{p}}}\frac{\mathrm{d}{a}_{\mathrm{p}}}{\mathrm{d}t}$$

$$\Gamma$$

$${\Gamma }_{\mathrm{gas}}$$

$${\Gamma }_{\mathrm{dust}}$$

$${\Gamma }_0(R)$$

$$T_{\mathrm{gas}}(R)$$

$$T_{\mathrm{dust}}(R)$$

$$s_{\mathrm{d},0}$$

$$=57.1\left(\frac{{\Sigma }_0}{8.9\text{mbox}{g/cm}^2}\right){\left(\frac{{\rho }_p}{1\text{mbox}{g/cm}^3}\right)}^{-1}\text{mbox}cm,$$

$$t_{\mathrm{stop}}$$

$$St$$

Full text

planet-disk interactions — accretion, accretion disks — protoplanetary disks — planets and satellites: formation

††journal: ApJL

1 Introduction

A planet within a protoplanetary disk interacts with a surrounding disk gas and migrates in radial direction (e.g., Goldreich & Tremaine, 1980; Lin & Papaloizou, 1979). A sufficiently massive planet forms a density gap as a consequence of disk–planet interaction, and it migrates together with the gap (e.g., Lin & Papaloizou, 1986; Ward, 1997). This planet migration continues until the gaseous disk is depleted, which determine the final architecture of planetary systems. Therefore, the planet migration is one of the most important processes in the formation of the planetary system.

The radial migration of the gap-opening planet is called type II migration. In the ideal case, the gap induced by the massive planet halts the gas flow across it and hence the planet migrates in the same velocity of the viscous drift of the gas (e.g., Lin & Papaloizou, 1986; Armitage, 2007). However, recent numerical hydrodynamic simulations have shown that the planetary gap hardly halts the gas flow across it, even if the planet is as massive as Jupiter (e.g., Duffell et al., 2014; Dürmann & Kley, 2015, 2017; Kanagawa et al., 2018b; Robert et al., 2018). By carrying out two-dimensional hydrodynamic simulations in a broad parameter range, Kanagawa et al. (2018b) have found that the torque exerted on the planet is roughly proportional to the surface density at the bottom of the gap, though strictly speaking, the torque is also exerted from the bottom and the edge of the gap. This result indicates that the migration speed of the gap-opening planet is sensitive to the gas structure of the gap.

At the outer edge of the gap, relatively large dust grains, so-called “pebbles”, are highly pilled-up, because these dust grains cannot cross the gap region (e.g., Paardekooper & Mellema, 2004; Muto & Inutsuka, 2009; Zhu et al., 2012; Dong et al., 2015; Pinilla et al., 2015; Weber et al., 2018; Meru et al., 2019; Kanagawa et al., 2018a). The dust ring formed at the outer edge of the gap might be associated with bright rings in the protoplanetary disks observed by e.g., ALMA (e.g., ALMA Partnership et al., 2015; Nomura et al., 2016; Tsukagoshi et al., 2016; Akiyama et al., 2016; van der Plas et al., 2017; Fedele et al., 2017; Dong et al., 2018; Z. C. Long et al., 2018; Long et al., 2018; Andrews et al., 2018). The dust grains lose angular momentum due to friction with the surronding disk gas and fall down to the central star. Simultaneously the disk gas receives the same amount of angular momentum that the dust grains lose, which is referred to as dust feedback (Fu et al., 2014; Gonzalez et al., 2015, 2017; Taki et al., 2016; Dipierro & Laibe, 2017; Kanagawa et al., 2017; Weber et al., 2018; Dipierro et al., 2018; Kanagawa et al., 2018a). Recent gas–dust two fluid hydrodynamic simulations (Weber et al., 2018; Kanagawa et al., 2018a) have found that though it is not significant, the dust feedback can modify the outer edge of the gap. This slight modification induced by the dust feedback may significantly affect the migration speed, because it could break the balance between torques exerted from the inner and outer disks.

We investigate the effects of the dust feedback on the migration speed of the gap-opening planet, by carrying out the gas–dust two fluid and two-dimensional hydrodynamic simulations. In Section 2, our numerical method and setup are described. The results of our hydrodynamic simulations are shown in Section 3. In Section 4, we summarize and discuss our results.

2 Numerical method

We implemented the dust component into the hydrodynamic code FARGO (Masset, 2000), and numerically solved the equations of motion and the continuity equations for the gas and the dust grains (Kanagawa et al., 2017, 2018a). We include the gas-dust drag force in the equation of motion of the gas, as well as that of the dust grains. Due to this drag force in the equation of the motion of the gas (dust feedback), the gas is affected by the dust grains, similar to that the dust grains drift because of the drag force with the gas. In the present study, we take into account the planet migration, which was not considered by the previous studies. We choose a two-dimensional cylindrical coordinate system ($R,\phi$), and its origin is located at the position of the central star.

2.1 Planetary migration

Here we briefly describe the method to implement the migration of the planet. We assume that the eccentricity of the planet is much smaller than unity (this assumption is valid in our simulations shown in this letter). In this case, the equation of motion of the planet is written by

[TABLE]

where $M_{\rm p}$, $a_{\rm p}$, and $\Omega_{\rm K,p}$ denote the mass, semi-major axis, and the Keplerian angular velocity of the planet, respectively. In the following, the subscript “${\rm p}$” indicates a physical quantity at the orbital radius of the planet $R=R_{\rm p}$. Since we consider a small eccentricity, $R_{\rm p}\simeq a_{\rm p}$. The torque exerted on the planet $\Gamma$ is given by the sum of these exerted from the gas disk and the dust disk, as

[TABLE]

where $\Gamma_{\rm gas}$ and $\Gamma_{\rm dust}$ are given by (Kley & Nelson, 2012)

[TABLE]

where $\psi$, $R_{\rm in}$ and $R_{\rm out}$ denote the gravitational potential, and the radii of the inner and outer edge of the disk, respectively. The surface densities of the gas and the dust grains are represented by $\Sigma_{\rm gas}$ and $\Sigma_{\rm dust}$, respectively. We introduce a characteristic torque $\Gamma_{0}$ as

[TABLE]

where $\Omega_{\rm K}$ and $\Sigma_{\rm gas}{}_{\rm,un}$ denote the Keplerian angular velocity and the gas surface density of the unperturbed disk, respectively, which are the functions of $R$. Note that $\Gamma_{0}$ roughly corresponds to the absolute value of the torque expected by the Type I regime.

For convenience, we define a cumulative torque from the gas and the dust grains as

[TABLE]

Thus, the cumulative torque exerted from the gas and the dust grains are given by $T(R)=T_{\rm gas}(R)+T_{\rm dust}(R)$.

2.2 Setup of numerical simulations

We use the initial orbital radius of the planet $R_{0}$ and the mass of the central star $M_{\ast}$ as the unit of length and mass, respectively. Hence, the surface density is normalized by $M_{\ast}/R_{0}^{2}$. We also use the unit of time $t_{0}$ which is defined by $2\pi/\Omega_{\rm K}(R_{0})$. The computational domain runs from $R/R_{0}=0.3$ to $R/R_{0}=3.2$ with 512 radial zones (equally spaced in logarithmic space) and 1024 azimuthal zones (equal spaced). The initial surface density of the gas is set as $\Sigma_{\rm gas}{}_{\rm,un}=\Sigma_{0}\left(R/R_{0}\right)^{-1}$ with $\Sigma_{0}=10^{-4}$. For the dust grains, the initial surface density is set for the dust-to-gas mass ratio to be 0.01. The disk aspect ratio are assumed to be $h_{\rm g}/R=H_{0}(R/R_{0})^{1/4}$ which is assumed to be independent of time, where $H_{0}$ is set to be $0.05$ in this letter and $h_{\rm g}$ is the gas scale height. The initial angular velocity of the gas is given by $\Omega_{\rm K}\sqrt{1-\eta}$ and $\eta=-1/2(h_{\rm g}/R)^{2}\partial(\ln P)/\partial(\ln R)$. For the dust grains, the initial angular velocity is given by $\Omega_{\rm K}$. The initial radial velocities of the gas and the dust grains are set to be $-3\nu/R$ and zero, respectively, where $\nu$ is the kinetic viscosity. We use the $\alpha$-prescription of Shakura & Sunyaev (1973), and hence $\nu=\alpha h_{\rm g}^{2}\Omega_{\rm K}$. We also assume that $\alpha$ is constant throughout the computational domain. The smoothing parameter of the gravitational potential is set to be $0.6h_{\rm g,p}$. We set wave killing zones with the width of $0.1R_{0}$ for the gas at the inner and outer boundaries, as described in Kanagawa et al. (2018b).

We assume constant-size dust grains. The size of the dust grains $s_{\rm d}$ is set to be $0.1s_{\rm d,0}$, where $s_{\rm d,0}$ is defined by

[TABLE]

where $\rho_{p}$ denotes the internal density of the dust grains. Note that when $M_{\ast}=1M_{\odot}$ and $R_{0}=10\ \rm{AU}$, the dimensionless surface density $\Sigma=10^{-4}$ corresponds to 8.9 $\mbox{g/cm}^{2}$. Adopting the Epstein regime, one can write the stopping time of dust grains for the two-dimensional disk $t_{\rm stop}$ as

[TABLE]

The Stokes number of dust grains, $St$, is defined by

[TABLE]

3 Results

3.1 Halting inward planetary migration triggered by the dust feedback

First we present the result of our simulations in the case of a Jupiter-mass planet ($M_{\rm p}/M_{\ast}=1\times 10^{-3}$) with $\alpha=3\times 10^{-4}$ and $H_{0}=0.05$.

The top panel of Figure 1 compares the evolution of the semi-major axis of the planet in the case with and without the dust feedback. It is obvious that the migration rate is significantly lower when the dust feedback is included. In particular, after $t\simeq 12000\ t_{0}$, the inward migration is terminated and the planet migrates outward slowly.

The middle and bottom panels of Figure 1 show a total torque exerted on the planet and the torques exerted from the inner and outer disk, respectively. The absolute value of torque exerted on the planet in the case with the dust feedback is significantly smaller than that in the case without the dust feedback. After $t\simeq 12000t_{0}$, the total torque becomes positive and hence the planet moves outward in the case with the dust feedback. As can be seen from the bottom panel of Figure 1, only the torque exerted from the outer disk in the case with the dust feedback decreases as compared to that in the case without the dust feedback. Since the negative torque exerted from the outer disk descreases, the inward migration of the planet is halted and finally the planet migrates outward, in the case with the dust feedback.

The top and middle panels of Figure 2 show the azimuthal averaged surface densities of the gas and the dust grains at $t=12000\ t_{0}$, in the cases with and without the dust feedback. The dust feedback makes the gradient of the gas surface density small (in other words, the value of $\eta$ becomes small) in the outer edge of the gap, as shown by Kanagawa et al. (2018a); (see also Taki et al. (2016)). As a result, the surface density of the gas in the outer disk of $R>1.3R_{\rm p}$ becomes smaller than that without the dust feedback. Since the dust drift velocity is much faster than that of the planet, the dust grains drifting from the outer disk are trapped into the gas and they form the dust ring at the outer edge. In the inner disk, most of the dust grains drift to the central star. The distribution of the dust grains against the planet’s orbit is a consequence of the dust trap by the gap and the fast dust drift.

In the bottom panel of Figure 2, we show the cumulative torques. The torque is exerted from the bottom to the edge of the gap, as Kanagawa et al. (2018b) stated, and the torque in the case without the dust feedback is consistent to that given by the formula of Kanagawa et al. (2018b). In the case with the dust feedback, the slope of the cumulative torque is significantly less steeper than that in the case with the dust feedback because of the depletion of the outer disk due to the dust feedback, though they are similar to each other in other regions. Therefore, the torque exerted from the outer disk is smaller and finally the outward migration is caused in the case with the dust feedback, as seen in Figure 1.

3.2 Parameter study

In order to see how sensitive the effects of the dust feedback to planet masses ($M_{\rm p}$) and viscosities ($\alpha$), we performed a series of simulations.

First, in Figure 3, we show the case of $M_{\rm p}/M_{\ast}=5\times 10^{-4}$, but $H_{0}$ and $\alpha$ are the same as Figure 1. Similar to the case of Figure 1, the value of $|\Gamma/\Gamma_{0}|$ decreases because the torque exerted from the outer disk decreases, in the case with the dust feedback. As a result, the inward migration in the case with the dust feedback is much slower than that in the case without the dust feedback. However, in this case, since the gap is narrow and the dust ring is formed at closer from the planet, as compared to the case shown in Figure 1, the planet feels the negative torque exerted from the dust ring significantly. Becuase of it, the inward migration of the planet continues, though the torque exerted from the gas is positive after $t=12000\ t_{0}$.

The top panel of Figure 4 shows the dependence of the viscosity on the ratio of the torques with and without the feedback the dust feedback, $\Gamma_{\rm fb}$ and $\Gamma_{\rm nofb}$, for various $M_{\rm p}/M_{\ast}$ with $H_{0}=0.05$. Specifically, the value of torques in Figure 4 corresponds to the average over $11500\ t_{0}<t<12000\ t_{0}$ 111Only in the case of $M_{\rm p}/M_{\ast}=10^{-4}$ and $\alpha=10^{-3}$, we show the torque at $\sim 8000\ t_{0}$, because the planet reaches too close to the inner boundary at $\sim 12000\ t_{0}$.. Note that though we took the torque at $\sim 12000\ t_{0}$, the torque at that time may not exactly correspond to the stationary value, as in the case shown in Figure 1. Hence, the figure shows just a rough trend.

The top panel of Figure 4 indicates that $\Gamma_{\rm fb}/\Gamma_{\rm nofb}$ decreases as the decrease of $\alpha$, because the gap becomes wider with an decrease of $\alpha$. The bottom panel of Figure 4 shows this tendency more clear. The bottom panel is the same as the top panel, but the horizontal axis of the bottom panel is the dimensionless parameter $K^{\prime}=(M_{\rm p}/M_{\ast})^{2}(h_{\rm g,p}/R_{\rm p})^{-3}\alpha^{-1}$, which is related to the gap width of the gas in stationary state as $\Delta_{\rm gap}=0.41K^{\prime 1/4}$ (Kanagawa et al., 2016). Larger value of $K^{\prime}$ indicates wider gap. As can be seen in the bottom panel of Figure 4, the ratio of $\Gamma_{\rm fb}$ to $\Gamma_{\rm nofb}$ decreases as $K^{\prime}$ increases (the gap is wider), when $K^{\prime}\gtrsim 1$.

4 Summary and Discussion

In this letter, we have shown that the dust feedback reduces the surface density at the outer edge of the gap and the torque exerted on the planet from the outer disk. As a result, the radial migration of the planet can be significantly slow . In particular, when the viscosity is low, namely $\alpha=3\times 10^{-4}$, the inward migration is halted and the Jupiter-mass planet moves outward slowly.

In the rest of this letter, we discuss implications of our results for the planet formation and the observations. The planet migration within the disk continues until the disk gas removal, typically at a few million years, which determine the final position of the planet. Figure 5 shows the migration timescales, $\tau_{\rm mig}=-({\rm d}\ln a_{\rm p}/{\rm d}t)^{-1}$, calculated by our simulations against the dimensionless parameter $K^{\prime}$ (which is mentioned in the previous section). The estimated width of the gap from $K^{\prime}$ is indicated in the upper horizontal axis of the figure. This estimated width is just a reference because the simulations of Kanagawa et al. (2016) does not consider the dust feedback, and also the dust gap can be wider than the gap of the gas, depending on the size of the dust grains (e.g., Zhu et al., 2012; de Juan Ovelar et al., 2013). We can consider the estimated width as a rough lower limit of the gap width.

As can be seen in Figure 5, the migration timescale is longer as the increase of $K^{\prime}$. In particular, when $K^{\prime}\gtrsim 1$ ($\Delta_{\rm gap}/R_{\rm p}\gtrsim 0.5$), the timescale of the inward migration is about $10\mbox{ Myr}$, which is a order-of-manugnitude larger than the typical disk lifetime. Moreover, when the gap is sufficiently wide, namely, $K^{\prime}=20$, the planet migrates outward as shown in Section 3.1. Note that in the figure, the relative light disk is assumed as $M_{\rm disk}=0.006M_{\odot}$. For massive disks such as the disk of HL Tau ($M_{\rm disk}\simeq 0.1M_{\odot}$), the migration timescale is shorter ($\times\sim 0.1$) than that shown in Figure 5. Even for the massive disk, in the cases with the timescale of $\sim 10\mbox{ Myr}$ shown in Figure 5, the planet can survive until the removal of the disk gas.

Recent observations have revealed that a number of the protoplanetary disks have ring and gap structures. Here we discuss the observability of the gap formed by the planet. The timescales of the gap formation and planet migration are proportional to the Keplerian rotation period. That is, these timescales are longer in the region with larger radii. Hence, in the outer region, the gap shape may not reach stationary state yet and the width of the gap is narrower than the stationary width. By considering the long timescale of the gap formation and migration, it is understandable that we can observe the narrow gap in the outer region. However, in the inner region, the migration timescale becomes shorter as the Keplerian rotation is faster. Hence, the gap is difficult to be observed in the inner region due to the fast inward migration, without the mechanism which slows the migration down. Only when the planet forms a relatively wide gap, namely $\Delta_{\rm gap}/R_{\rm p}\gtrsim 0.5$, the planet could survive during the disk lifetime, due to the effect of the dust feedback. Hence, the width of the gap would be wide in the inner region, if the observed gap is formed by the planet. Actually, ALMA observations have found such a wide gaps in the relatively inner region, for instance, in the disks of HD143006 (Huang et al., 2018), HL Tau (ALMA Partnership et al., 2015; Kanagawa et al., 2016), and PSD 70 (Z. C. Long et al., 2018). These gaps could be formed by the planet. In the near future, ALMA will reveal the gaps within the region of $\lesssim 10\mbox{ AU}$. Our results imply that these gaps are as wide as $\Delta_{\rm gap}/R_{\rm p}\gtrsim 0.5$.

Within the dust ring in which the dust grains are highly accumulated, the planetesimals and protoplanets can be formed via e.g., the streaming instability (e.g., Youdin & Goodman, 2005). The planet migration is very slow until the dust grains are removed by the planetesimal formation. Once the planet formation starts within the dust ring, a number of planets can be formed within the dust ring as discussed in Kanagawa et al. (2018a). Otherwise, a remnant of the dust ring may survive. Such a path of the planet formation may explain the large diversity of the exoplanets revealed by observations.

The suppression and turnover of the inward migration of the gap-opening planet shown in this letter are critical for the planetary system formation and the origin of the observed gaps. Although we investigated only the cases with $H_{0}=0.05$, the further survey covering a wider parameter range should be done to determine the condition of the onset of the outward migration and quantitative relation of reduced migration time and planet mass.

KDK would like to thank Takayuki Muto for an insightful discussion and constructive comments, and Professor Yasushi Suto for helpful comments. KDK also thank the anonymous referee for thoughtful comments. KDK was supported by JSPS Core-to-Core Program “International Network of Planetary Sciences” and JSPS KAKENHI grant No. 19K14779. Numerical computations were carried out on the Cray XC50 at the Center for Computational Astrophysics, National Astronomical Observatory of Japan.