Trapping low-mass planets at the inner edge of the protostellar disc

Trapping low-mass planets at the inner edge of the protostellar disc

R. Brasser11affiliation: Earth Life Science Institute, Tokyo Institute of Technology, Tokyo, Japan , S. Matsumura22affiliation: Division of Physics, University of Dundee, Dundee, UK , T. Muto33affiliation: Division of Liberal Arts, Kogakuin University, Tokyo, Japan and S. Ida11affiliation: Earth Life Science Institute, Tokyo Institute of Technology, Tokyo, Japan

The formation of multiple close-in low-mass exoplanets is still a mystery. The challenge is to build a system wherein the outermost planet is beyond 0.2 AU from the star. Here we investigate how the prescription for type I planet migration affects the ability to trap multiple planets in a resonant chain near the inner edge of the protostellar disc. A sharp edge modelled as a hyperbolic tangent function coupled with supersonic corrections to the classical type I migration torques results in the innermost planets being pushed inside the cavity through resonant interaction with farther planets because migration is starward at slightly supersonic eccentricities. Planets below a few Earth masses are generally trapped in a resonant chain with the outermost planet near the disc edge, but long-term stability is not guaranteed. For more massive planets the migration is so fast that the eccentricity of the innermost resonant pair is excited to highly supersonic levels due to decreased damping on the innermost planet as it is pushed inside the cavity; collisions frequently occur and the system consists one or two intermediate-mass planets residing closer to the star than the disc’s inner edge. We found a neat pileup of resonant planets outside the disc edge only if the corotation torque does not rapidly diminish at high eccentricity. We call for detailed studies on planet migration near the disc’s inner edge, which is still uncertain, and for an improved understanding of eccentricity damping and disc torques in the supersonic regime.

celestial mechanics — planets and satellites: dynamical evolution and stability — planets and satellites: formation

1 Introduction

Preventing low-mass planets from migrating to their host star is a long-standing problem. Mergers with their host star can be prevented if the protostellar disc has a sharp inner edge at a few stellar radii (Masset et al., 2006; Ogihara et al., 2010). Terquem & Papaloizou (2007) showed that migrating protoplanets usually end up in resonances. Some of these resided inside the disc’s inner cavity. Ogihara & Ida (2009) also build a resonant chain of low-mass protoplanets that stalled near the disc’s inner edge when the migration was artificially slowed down and the reduction of the corotation torque was ignored. Recently Izidoro et al. (2017) trapped a high number of low-mass planets outside the disc edge in a resonant chain, that subsequently needed to break to account for the currently-observed exoplanet period distribution. On the other hand, Matsumura et al. (2017) had trouble trapping multiple planets near the disc’s inner edge, even though their migration prescription was very similar to that of Izidoro et al. (2017): both include supersonic corrections to the migration and eccentricity damping timescales. Terquem & Papaloizou (2007) also included such corrections, but they followed the prescription of Papaloizou & Larwood (2000) while Izidoro et al. (2017) and Matsumura et al. (2017) followed Coleman & Nelson (2014). The simulations of Matsumura et al. (2017) usually resulted in one or two hot Neptune planets rather than a multiplet of smaller planets. The disparity between all of these results warrants further study.

2 Disc model and planet migration

We employ the disc model of Ida et al. (2016), which is based on Garaud & Lin (2007) and Oka et al. (2011). Here we briefly summarise their model.

2.1 Disc parameters

We assume steady-state accretion onto the Sun. The gas accretion rate is


where is the gas surface density, is the disc scale height and is the orbital frequency. The -viscosity is assumed to be constant (Shakura & Sunyaev, 1973). The disc scale height is , where with , is the Boltzmann constant, the proton mass and is the mean atomic mass of the gas. Now evolves as (Hartmann et al., 1998)


The extra 0.1 Myr avoids the logarithmic singularity (Bitsch et al., 2015b).

For solar-type stars, the midplane temperature in the viscous part of the disc is


where is the distance to the star and we defined and . The reduced scale height is


Equations (1), (3) and (4) may be combined to compute the surface density of the gas.

2.2 Disc inner edge implementation

Near the star the surface density of the disc is assumed to smoothly decrease to zero. Cossou et al. (2014) suggested that


where is the planet trap location from the star; here the surface density is maximal. The trap is at  AU for solar-type stars; it is unclear how reliable the employed disc model is closer to the star where MHD effects become important. Trapping planets requires a sharp edge (Masset et al., 2006; Ogihara et al., 2010). Therefore we set the inner edge of the disc at , which is approximately 2 scale heights inside of . The surface density slope is computed as


where . Now as and all planets in the type I regime cease migrating at or near .

To avoid the divergence of at , which could possibly cause numerical artefacts, we also tested a linear decrease in , for which . The resulting discontinuity of at was made smooth via a linear connection over length .

2.3 Planet migration

The gas disc exerts torques and tidal forces on the embedded planets which result in a combined effect of radial migration and the damping of the eccentricity and inclination. For low-mass planets the migration is of type I (Tanaka et al., 2002) while massive planets that are able to clear the gas in their vicinity experience type II migration (Lin & Papaloizou, 1986). Here we are only interested in the former.

We follow Coleman & Nelson (2014) for computing the torque and the direction of migration. Their formulae are based on Paardekooper et al. (2011) for the torque and on Fendyke & Nelson (2014) and Cresswell & Nelson (2008) for the eccentricity damping, including corrections to the damping timescale and corotation torque in the supersonic regime when the eccentricity . We restrict ourselves to planar orbits. The normalised torque is (Paardekooper et al., 2011)


where and are the corotation and Lindblad torques respectively and is a normalisation constant. The corotation torque in a non-isothermal disc with thermal diffusion becomes


Here and depend on , and (Paardekooper et al., 2011). The functions , and determine the amount of torque saturation and are dependent on the planet mass and disc scale height (Paardekooper et al., 2011). In steady state , where . The remaining contributions are then (Paardekooper et al., 2011)


where is the negative entropy gradient. The , and functions disallow writing the explicit dependence of on , but generally and .

The factors and are


where and . The eccentricity damping timescale is (Cresswell & Nelson, 2008)


where the wave timescale is (Tanaka & Ward, 2004)


The ’migration timescale’ (Cresswell & Nelson, 2008) is , and is


Despite its name, is not the actual migration timescale. By definition (Cresswell et al., 2007) and , so that when . Here is the timescale for semi-major axis evolution. In the supersonic regime hydrodynamical simulations show that the torque reverses direction when and the torque is maximal (and positive) when (Cresswell & Nelson, 2008). However, Cresswell & Nelson (2008) show that the planet always migrates inwards, despite the torque reversal at high eccentricity. The analytical approach of Muto et al. (2011) agrees with this result: inward migration persists at high eccentricity. Therefore, it is incorrect to use to compute the evolution of the semi-major axis; should be used instead (Matsumura et al., 2017). This is


Izidoro et al. (2017) adopted for the migration timescale, which led to the aforementioned difference. All timescales are positive for inward migration and negative for outward migration.

Figure 1: Plot of various timescales and normalised torques vs eccentricity for a planet. Left: at 1 AU. Right: at 0.099 AU. Dotted lines are for negative values.

2.4 Planet migration near the disc edge

The reciprocal migration timescale, , and its two components, and , are plotted in the top panels of Fig. 1 for the fiducial case (all supersonic corrections enabled) as a function of planetary eccentricity. The middle panels show the individual timescales and the bottom panels depict the normalised torques. The left column pertains to a 1  planet at 1 AU while the to the right the planet is at 0.099 AU. The outcome for the linear edge prescription is very similar.

At 1 AU when we generally have and angular momentum transfer between the disc and the planet usually leads to inward migration (); the example in Fig. 1, however, has outward migration when because for this mass and the strength of the corotation torque slightly exceeds that of the Lindblad torque (Paardekooper et al., 2011; Brasser et al., 2017); see bottom of Fig. 1. When we have and the inward migration is the result of eccentricity damping at constant angular momentum. These two cases are common to all the combinations of the torque formulae. However, the relations for , and are different when for the various configurations.

At 0.099 AU (near the disc’s inner edge), these two cases are also applicable, except that the difference between and when is smaller than at 1 AU. When at 0.099 AU, however, the relations for , and are different from that at 1 AU.

These fundamental differences are clearly seen between the two panels, assuming that equation (7) and the supersonic corrections are applicable at the disc edge.

Suppose there is a single planet at . When a second planet approaches the first one, it is likely to become trapped in resonance (Petrovich et al., 2013). The equilibrium eccentricity of both planets is the result of the balance of migration and damping, and is (Goldreich & Schlichting, 2014). At this eccentricity , but when , and the planets migrate inwards: the innermost planet has its eccentricity excited as the outer planet pushes the pair starward until the latter reaches ; the inner planet is now parked inside the disc cavity. Should a third, fourth and additional planets approach the inner pair the mechanism will likely repeat itself until the outermost planet is at the disc edge and all the other planets are inside of it, or until the outermost planets can no longer push the inner chain deeper into the cavity due to the chain’s inertia.

The above arguments assume that the planet is surrounded by a gas disk on both sides. Near the edge, the torque is one-sided: at the inner disk edge there is only the outer Lindblad torque that pushes the planet inwards, but no inner Lindblad torque. Liu et al. (2017) implemented simplified one-sided corotation and Lindblad torques and showed that small planets could be trapped near the edge, but their implementation only works for an infinitely sharp edge.

3 Numerical methods

To study the behaviour of planets near the disc edge we performed a set of numerical N-body simulations consisting of a solar-mass star and four equal-mass planets. These integrations used the symplectic -body code SyMBA (Duncan et al., 1998), which was heavily modified to include the effects of eccentricity and inclination damping as well as planet migration by the gas disc according to the formulation above (Matsumura et al., 2017).

We initially place the planets beyond their 2:1 mean-motion resonances. The planet’s masses are all either 0.1, 0.5, 2, 3, 5 or 8 Earth masses (). We compute the torques at each time step for each body; we apply eccentricity damping only when . There is a surface density maximum at  AU from the star. Closer to the star than the disc edge at there is no migration nor damping.

Simulations are run for 2 Myr with a time step of 0.146 d. Bodies are removed when they are closer than 0.02 AU or farther than 100 AU from the star, or when they collide. We assume perfect accretion during collisions. Initially and . For simplicity we keep fixed despite its potential to change close to the star where MRI effects (Bai & Stone, 2013) and disc ionisation are important (Gammie, 1996). The value of affects the torque in a complicated manner as described in equation (8) (Paardekooper et al., 2011; Brasser et al., 2017). As we show below, the most important factor in altering the torque is the supersonic corrections, in particular the corotation reduction, , which is independent of and .

We test the dependence of the migration on eccentricity by selectively enabling or disabling the supersonic corrections to the corotation torque, the Lindblad torque and the eccentricity damping timescale. Apart from Masset et al. (2006) we are not aware of any systemic hydrodynamical studies of embedded planets near the disc edge. Motivated by this, we enable or disable (’on’ or ’off’) the corotation reduction term ( in equation (10), ’CR’) as well as the supersonic corrections to the eccentricity damping timescale (the term inside the parentheses in equation (11), ’ED’) and to the Lindblad torque ( in equation (10), ’EL’). This simple approach should suffice in the absence of more detailed torque prescriptions near the edge.

4 Results

We show an example of the evolution of four planets of near the disc edge in Fig. 2. The evolution is very different depending on the migration prescription. In the broadest sense the planets become trapped in resonances outside the disc edge when the corotation reduction is not applied, and preferably when the supersonic correction to the Lindblad torque, , is applied. If either or both of these are applied, multiple planets will end up inside the disc cavity.

Figure 2: Migration of four 2 planets towards the disc edge. The panels all depict different migration prescriptions.

This outcome is in contrast with Izidoro et al. (2017) because they use to compute the semi-major axis evolution of the planets instead of . Since at the disc edge for any value of , their innermost planet can often stall any additional incoming planets even if it or the others are supersonic, though the exact evolution is mass dependent. For example, in their Figure 4, the innermost planet is 10 and all planets are trapped beyond 0.1 AU, while in their Figure 5 the innermost planet (initially) has a few and a chain of planets is pushed inward as more massive planets migrate in. In resonance the approximation breaks down because the planets are supersonic; equation (14) should be used instead.

Our results are also in disagreement with Ogihara et al. (2010), who found that a sharp edge was able to prevent the planets from falling into the cavity. Ogihara et al. (2010) concluded that the imbalance of increased drag on the planet at aphelion versus little to no drag at perihelion caused the planet to be stationary at the edge with a non-zero eccentricity. However, they did not include any supersonic corrections to their migration formulae.

In our approach the edge is fairly sharp, but the the tanh function quickly flattens beyond . Thus if the planet is near the drag at aphelion and perihelion is within 25% if the eccentricity is . Only when the planet is roughly halfway between and is the drag at aphelion (where the tanh function is 1) much stronger than at perihelion (where it is 0) and do we possibly recover the situation from Ogihara et al. (2010), but only if the supersonic reduction of the corotation torque is ignored.

Figure 3: Left: Contour map of in Myr for a planet of 1 at 1 AU as a function of and . Middle: same at 0.099 AU. Right: in Myr at 0.099 AU. Bottom panels are the same but now i.e. there is no corotation reduction at high eccentricity.

4.1 Different torque prescriptions yield different outcomes

What aspects of the torque prescription are responsible for the different behaviour at the edge? At low eccentricity outward migration is caused by the corotation torque while inward migration is caused by the Lindblad torque; the latter only depends on the temperature and surface density gradients (Paardekooper et al., 2011). When supersonic corrections are considered the total torque also depends on . In general with our prescription is not a monotonically decreasing function of , but instead has a peak near and a trough at (Cresswell & Nelson, 2008). This behaviour is caused by how varies with eccentricity, and also because the maxima of and do not coincide (top-right panel of Fig. 1). The migration rate peaks for all planetary masses near , reaches a minimum when and a further maximum when . This non-monotonic behaviour is inconsistent with the analytical results of (Muto et al., 2011).

However, the behaviour is different when the planet is inside of the trap. At very low eccentricity, and for nearly all planet masses, is negative and low, implying slow outward migration. Migration is inward when and outward again at higher eccentricity. The top row of Fig. 3 is a contour map showing at 1 AU (left), near the disc’s inner edge (middle), together with near the edge (right).

The region of inward migration near moderate prevents the trap from stalling migrating planets in resonances, and any planet situated at is pushed deep inside the cavity, together with any planets interior to it. Therefore, the outcome of numerical simulations and the ability to trap planets near the disc edge depends on the exact migration prescription employed. When using the migration is always outwards.

Figure 2 suggests that eliminating the corotation reduction, , and weakening the Lindblad torque by applying provides the best prescription to trap multiple planets in resonance outside of , assuming the current torque formulae hold near the edge (cf. Liu et al. (2017)). The bottom row of Fig. 3 shows similar contour maps but now i.e. there is no corotation reduction. The behaviour is qualitatively different everywhere: at the edge migration is always outward, but the strength is a complicated function of both the planetary mass and the eccentricity.

Hydrodynamical simulations show that the corotation torque weakens as the eccentricity increases and mostly disappears at (Cresswell et al., 2007). It thus appears to be unphysical to remove the corotation reduction far from the disc edge, but it is unclear if this removal is applicable near the edge. The exponential reducion of Fendyke & Nelson (2014) does not appear to hold for low values of ; the corotation reduction also depends on how the torque is measured. Their Figs. 4 and 9 clearly show torque maxima near so that the exponential reduction may not be universally applicable. Clearly more work is needed, both on the reduction itself but also how it behaves near the disc edge.

Figure 4: Same as Fig. 2 but the planets are now 5.

In Fig. 2 the bottom-middle panel with all the supersonic effects enabled was able to temporarily trap the planets in a resonance even though the innermost planets were pushed into the cavity. Increasing the corotation torque kept the planets outside of the cavity (CR: OFF). This structure does not hold for higher-mass planets because they migrate faster and therefore excite themselves to higher eccentricities once the innermost planet is in the cavity and eccentricity damping is weak or non-existent. An example is shown in Fig. 4, which is the same as Fig. 2 but now the planets are 5.

5 Conclusions

The ability to trap multiple low-mass planets in a resonant chain outside the inner edge of the protostellar disc has been investigated. These low-mass planets execute type I migration which pulls them invariably towards the star. In the absence of a barrier these would all collide with the star. The disc’s inner edge could provide a trapping mechanism (Masset et al., 2006). We have tested two types of sharp inner edges of the disc: a hyperbolic tangent and a linear function, along with different migration prescriptions.

We find that a neat pileup of resonant planets outside the disc edge is established if the corotation torque does not rapidly diminish at high eccentricity. The expectation is that if the resonant chain of the planets remains outside the inner disc edge they eventually start orbit crossing and instigate a phase of giant impacts. This may account for formation of similar-sized, regularly spaced, non-resonant low-mass planets that are found to be common in relatively close-in regions by Kepler observations. However, the eccentricity damping and disc torques in the supersonic regime remain uncertain near the disc’s inner edge. Due to resonant interactions, eccentricity is generally excited to values for which the migration is generally inward. Therefore we call for detailed studies on eccentricity damping and disc torques in the supersonic regime and near the disc edge. Such a study will play an important role in understanding the common architecture of compact systems.


  • Bai & Stone (2013) Bai, X.-N., & Stone, J. M. 2013, ApJ, 767, 30
  • Bitsch et al. (2013) Bitsch, B., Crida, A., Morbidelli, A., Kley, W., & Dobbs-Dixon, I. 2013, A&A, 549, A124
  • Bitsch et al. (2015a) Bitsch, B., Johansen, A., Lambrechts, M., & Morbidelli, A. 2015, A&A, 575, A28
  • Bitsch et al. (2015b) Bitsch, B., Lambrechts, M., & Johansen, A. 2015, A&A, 582, A112
  • Brasser et al. (2017) Brasser, R., Bitsch, B., & Matsumura, S. 2017, AJ, 153, 222
  • Coleman & Nelson (2014) Coleman, G. A. L., & Nelson, R. P. 2014, MNRAS, 445, 479
  • Cossou et al. (2014) Cossou, C., Raymond, S. N., Hersant, F., & Pierens, A. 2014, A&A, 569, A56
  • Cresswell et al. (2007) Cresswell, P., Dirksen, G., Kley, W., & Nelson, R. P. 2007, A&A, 473, 329
  • Cresswell & Nelson (2008) Cresswell, P., & Nelson, R. P. 2008, A&A, 482, 677
  • Duncan et al. (1998) Duncan, M. J., Levison, H. F., & Lee, M. H. 1998, AJ, 116, 2067
  • Fendyke & Nelson (2014) Fendyke, S. M., & Nelson, R. P. 2014, MNRAS, 437, 96
  • Gammie (1996) Gammie, C. F. 1996, ApJ, 457, 355
  • Garaud & Lin (2007) Garaud, P., & Lin, D. N. C. 2007, ApJ, 654, 606
  • Goldreich & Schlichting (2014) Goldreich, P., & Schlichting, H. E. 2014, AJ, 147, 32
  • Hartmann et al. (1998) Hartmann, L., Calvet, N., Gullbring, E., & D’Alessio, P. 1998, ApJ, 495, 385
  • Ida et al. (2016) Ida, S., Morbidelli, A. & Guillot, T. 2016. A& A 600, 154
  • Izidoro et al. (2017) Izidoro, A., Ogihara, M., Raymond, S. N., et al. 2017, MNRAS, 470, 1750
  • Lin & Papaloizou (1986) Lin, D. N. C., & Papaloizou, J. 1986, ApJ, 309, 846
  • Liu et al. (2017) Liu, B., Ormel, C. W., & Lin, D. N. C. 2017, A&A, 601, A15
  • Masset et al. (2006) Masset, F. S., Morbidelli, A., Crida, A., & Ferreira, J. 2006, ApJ, 642, 478
  • Matsumura et al. (2017) Matsumura, S., Brasser, R., & Ida, S. 2017, A&A, 607, A67
  • Muto et al. (2011) Muto, T., Takeuchi, T., & Ida, S. 2011, ApJ, 737, 37
  • Ogihara et al. (2010) Ogihara, M., Duncan, M. J., & Ida, S. 2010, ApJ, 721, 1184
  • Ogihara & Ida (2009) Ogihara, M., & Ida, S. 2009, ApJ, 699, 824
  • Ogihara et al. (2015) Ogihara, M., Morbidelli, A., & Guillot, T. 2015, A&A, 584, L1
  • Oka et al. (2011) Oka, A., Nakamoto, T., & Ida, S. 2011, ApJ, 738, 141
  • Paardekooper et al. (2011) Paardekooper, S.-J., Baruteau, C., & Kley, W. 2011, MNRAS, 410, 293
  • Papaloizou & Larwood (2000) Papaloizou, J. C. B., & Larwood, J. D. 2000, MNRAS, 315, 823
  • Petrovich et al. (2013) Petrovich, C., Malhotra, R., & Tremaine, S. 2013, ApJ, 770, 24
  • Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337
  • Tanaka et al. (2002) Tanaka, H., Takeuchi, T., & Ward, W. R. 2002, ApJ, 565, 1257
  • Tanaka & Ward (2004) Tanaka, H., & Ward, W. R. 2004, ApJ, 602, 388
  • Terquem & Papaloizou (2007) Terquem, C., & Papaloizou, J. C. B. 2007, ApJ, 654, 1110
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description