A&A 506, L9L12 (2009)
LETTER TO THE EDITOR
Discplanet interactions in subkeplerian discs
S.J. Paardekooper
Department of Applied Mathematics and Theoretical Physics, University of Cambridge, Centre for Mathematical Sciences, Wilberforce Road, Cambridge CB3 0WA, UK
Received 26 August 2009 / Accepted 8 September 2009
Abstract
Context. One class of protoplanetary disc models, the Xwind
model, predicts strongly subkeplerian orbital gas velocities, a
configuration that can be sustained by magnetic tension.
Aims. We investigate discplanet interactions in these
subkeplerian discs, focusing on orbital migration for lowmass planets
and gap formation for highmass planets.
Methods. We use linear calculations and nonlinear hydrodynamical
simulations to measure the torque and look at gap formation. In both
cases, the subkeplerian nature of the disc is treated as a fixed
external constraint.
Results. We show that, depending on the degree to which the disc
is subkeplerian, the torque on lowmass planets varies between the
usual type I torque and the onesided outer Lindblad torque, which
is also negative but an order of magnitude stronger. In strongly
subkeplerian discs, corotation effects can be ignored, making migration
fast and inward. Gap formation near the planet's orbit is more
difficult in such discs, since there are no resonances close to the
planet accommodating angular momentum transport. The location of the
gap is shifted inwards with respect to the planet, leaving the planet
on the outside of a surface density depression.
Conclusions. Depending on the degree to which a protoplanetary
disc is subkeplerian, discplanet interactions can be very different
from the usual Keplerian picture, making these discs in general more
hazardous for young planets.
Key words: planetary systems: protoplanetary discs  planets and satellites: formation
1 Introduction
Planets form in circumstellar discs, and the gravitational interaction between planet and disc plays a major role in shaping planetary systems. Tidal waves excited by the planet lead to orbital migration (Goldreich & Tremaine 1979), which comes in several flavours. type I migration (Ward 1997) is thought to hold for lowmass planets, up to a few times the mass of the Earth ( ). Highmass planets, comparable to Jupiter (with mass ), can tidally truncate the disc and open up a deep annular gap around their orbit. The resulting type II migration is driven by the viscous evolution of the disc (Lin & Papaloizou 1986). A strongly dynamical form of migration, type III (Masset & Papaloizou 2003), can be achieved for intermediatemass planets that partially open up a gap in massive discs, or by planets placed on a strong density gradient (Peplinski et al. 2008). For a recent overview of discplanet interactions, see Papaloizou et al. (2007).
Most studies of discplanet interactions have focused on discs that have a Keplerian velocity profile, possibly with a slight correction for a radial pressure gradient. One class of protoplanetary disc models, however, predicts strongly subkeplerian motion due to magnetic tension (Shu et al. 2007). In an interesting study, Adams et al. (2009) point out that in these Xwind discs, migration due to the strong head wind experienced by the planet, which is moving with Keplerian velocity, can overcome type I migration for Earthsized planets. On the other hand, gap formation should not be affected that much, with the planet located slightly towards the inner edge of the gap.
In this Letter, we study discplanet interactions in subkeplerian discs, under the simplifying assumption that the subkeplerian nature of the disc is an external and fixed constraint. This allows us to use purely hydrodynamical models as a first approximation. Future models should selfconsistently include the evolution of the magnetic field. We outline the numerical and physical setup in Sect. 2, present the results in Sect. 3 and conclude in Sect. 4.
2 Disc model
The basic equations are conservation of mass, momentum, and energy in a 3D cylindrical geometry . We make two simplifications at this point: by considering an isothermal equation of state, we remove the need to solve the energy equation, and by integrating the resulting system vertically, we end up with a twodimensional problem. Although it has been pointed out that the isothermal approximation is not valid in the inner regions of protoplanetary discs (Paardekooper & Mellema 2006a,2008), we see below that in strongly subkeplerian discs these effects play essentially no role. The sound speed is given by , where H is the scale height of the disc and the equilibrium angular velocity (see below). We varied , where is the orbital radius of the planet, between 0.05 and 0.2. The higher values of H are appropriate if the disc close to the star is ``puffed up'' by direct irradiation. A kinematic viscosity was used, with such that the initial, constant surface density is a stationary solution with no accretion. In the text, we quote the corresponding parameter at the location of the planet, . The resulting set of equations can be found elsewhere (see Paardekooper & Papaloizou 2009). The value of the surface density at , , can be chosen arbitrarily in nonselfgravitating discs; we normalise the torque accordingly (see Sect. 3.1).
The gravitational potential entering the equations contains terms due to the central star and the direct and indirect component of the planet potential. We use a softened point mass potential for the planet, with a softening parameter b = 0.6 H (see Paardekooper & Papaloizou 2009). This value of b is appropriate to account for 3D effects. The planet is assumed to be on a fixed circular orbit. For the central star potential we use , where f is a dimensionless parameter governing the degree to which the disc is subkeplerian. It is easy to see that the resulting equilibrium angular velocity, neglecting pressure effects, is given by , where is the Keplerian angular velocity. For f=1, we have a purely Keplerian disc, while for f<1, the disc is subkeplerian. Xwind models suggest , or even during FU Orionis outbursts (Shu et al. 2007).
We combine a linear code (Paardekooper & Papaloizou 2009) with fully nonlinear hydrodynamic calculations with RODEO (Paardekooper & Mellema 2006b) to solve the resulting set of equations. The computational domain extends from to , and covers the full in azimuth. This domain is covered with a regular grid with 256 cells in the radial and 768 cells in the azimuthal direction. Nonreflecting boundary conditions were used. We vary the planet mass, , and f to study different regimes of discplanet interactions in subkeplerian discs.
3 Results
The most important difference between subkeplerian discs and Keplerian discs, as far as discplanet interactions are concerned, is the change in position of the Lindblad and corotation resonances. It was shown in Goldreich & Tremaine (1979) that an embedded object can exchange angular momentum with the disc only at the location of a resonance. Corotation resonances occur where the disc corotates with the planet, and inner (outer) Lindblad resonances are located inside (outside) corotation. The strongest Lindblad resonances occur at a distance H of the planet Ward (1997). In Keplerian discs, a planet on a circular orbit is located close to corotation, with Lindblad resonances on either side of its orbit. In strongly subkeplerian discs, the distance of the planet to corotation is large enough so that we can safely neglect any nonbarotropic effects associated with the corotation torque (Paardekooper & Papaloizou 2008). In this sense, subkeplerian discs are easier to handle, since we only have to deal with Lindblad resonances.
The approximate position of the Lindblad resonances up to 10th order are shown in Fig. 1 for different values of f. Pressure effects have been ignored in calculating these positions. In a purely Keplerian disc (f = 1), the planet interacts with both inner and outer resonances. For f < 1, the resonance positions move inward, with the result that the planet now mainly interacts with the outer Lindblad resonances. This has consequences for both migration and gap formation.
Figure 1: Position of the inner (open triangles) and outer (solid circles) Lindblad resonances in discs with varying f. Pressure effects have been ignored in calculating the resonance positions. 

Open with DEXTER 
Figure 2: Absolute value of the torque on a planet embedded in discs of constant initial surface density for different values of h and f. Solid curves denote the results of linear calculations, while the symbols denote results from hydrodynamical simulations. 

Open with DEXTER 
3.1 Migration of lowmass planets
We first consider the torque on lowmass planets, not massive enough to open up a gap. We compare results from linear calculations to hydrodynamical simulations for different values of f in Fig. 2. The torque is normalised by , where and is the planet's angular velocity. We plot the absolute value of the torque; it is negative for all values of f.
It is immediately clear from Fig. 2 that the torque very strongly depends on f. This can be understood in terms of a change in the position of the resonances with respect to the planet. For f=1 we find the classical type I torque, modified by nonlinear effects at corotation (Paardekooper & Papaloizou 2009). Around , the planet mainly interacts with the strongest outer Lindblad resonances, yielding a very strong negative torque. The total torque approaches the onesided Lindblad torque (Ward 1997), which would be obtained when considering only outer Lindblad resonances and which is a factor of 1/h stronger than the differential Lindblad torque. The maximum torque depends on the adopted smoothing, with higher values for lower values of b. For , the torque decreases, since the resonances it interacts with become weaker as the planet moves further away from corotation. For h = 0.05 and f = 0.6, the resulting torque is comparable in magnitude to the classical type I torque. Linear calculations and nonlinear simulations agree very well, which again indicates that nonlinear effects at corotation do not play an important role.
The results depicted in Fig. 2 indicate that the crucial parameter determining the total torque on the planet is p=(1f)/h. The maxima of the curves occur factors of 2 apart in Fig. 2, confirming this picture. The scaling of the torque with ensures that for f=1, the curves fall on top of each other. When , a sizeable fraction of the onesided Lindblad torque can be exerted on the planet. The onesided torque scales as h^{3}, while the differential Lindblad torque scales as h^{2}. This is the reason for the different maxima for different values of h. This scaling holds for all f < 1; only for f=1 do we get the h^{2} scaling of the differential Lindblad torque. Therefore, in subkeplerian discs, for a fixed value of p, the torque scales as h^{3}. For a fixed value of f, the situation can be different. For example, for 1f=0.4, the torque for h=0.2 is twice as strong as would be expected from an h^{3} scaling, which is due to the shift of the maximum torque to higher values of f.
The torque due to a headwind in a subkeplerian disc on a planet with radius
is given by Adams et al. (2009):
with a constant of order unity (Weidenschilling 1977). This means that the ratio of the torque due to a headwind and the classical type I torque is proportional to h, and is of order unity for an Earth mass planet located at 1 AU, in a disc with f=0.66 and h = 0.05. Since most resonances are located far inward with respect to the orbit of the planet, the actual type I torque is close to its classical value at f = 1 (see Fig. 2). It is however very difficult to generalise this result due to the complex behaviour of the type I torque when varying h and f. For example, while for f = 0.6, at 1 AU for a planet of 1 in a disc with h = 0.05, the type I torque will in fact dominate for higher values of h. The dimension of the parameter space is quite high (q, , h, etc.), so that it is difficult to make more general statements.
It is important to realise that subkeplerian discs are hazardous environments for lowmass planets. Since the planet mainly interacts with outer Lindblad resonances, and corotation torques play a minor role, migration will be directed inward. The time scale for migration can be even smaller than the classical type I time scale, in unfavourable cases more than an order of magnitude smaller.
Figure 3: Azimuthally averaged surface density, in units of the initial (constant) surface density , for a planet embedded in a h=0.05, disc for different values of f. The inset shows the averaged surface density at the orbital radius of the planet as a function of f. 

Open with DEXTER 
3.2 Gap formation
We now turn to the issue of gap formation. The presence of a gap determines whether a planet moves in the type I or the type II regime, the latter being a much slower mode of migration (Ward 1997). There are two requirements that have to be fulfilled to open a gap: angular momentum transport due to the presence of the planet must exceed the viscous transport, otherwise the gap could be refilled by viscous evolution of the disc, and, in discs with low viscosity, the waves transporting angular momentum should be damped locally. For standard parameters, both criteria give similar minimum masses for gap opening, of the order of 100 . A unified criterion is derived in Crida et al. (2006).
Adams et al. (2009) derived a minimum mass for gapopening in subkeplerian discs, . They noted that for f=0.66 this mass is comparable to the classical estimates for Keplerian discs, although the planet would be located more towards the inner edge of the gap. In this section, we study gap formation for a (corresponding to 0.5 around a 1 star) planet embedded in a disc with h=0.05 and a viscosity corresponding to . For , this planet fulfils the criterion for gap formation as given above, if only slightly for f=0.4.
The results are depicted in Fig. 3, where we show the azimuthally averaged surface density for different values of f after 500 orbits of the planet. The system has reached a steady state by that time. For a Keplerian disc, the surface density close to the planet has dropped by more than an order of magnitude, confirming that this planet indeed opens up a gap. The main difference between Keplerian and subkeplerian discs is the position of the gap, which is shifting inward for f<1. This can again be understood in terms of resonance positions: the important resonances are located inside the orbit of the planet. The associated waves are damped locally through shocks, which results in a gap at . Therefore, contrary to the findings of Adams et al. (2009), the planet is located near the outside of the gap. The analysis of Adams et al. (2009) only applies to the case 1f < h, for which the location of the original gap (for f=1) is still densely populated by resonances. In this case, because the outer resonances are located closer to the planet (but still outside its orbit), the outer part of the gap forms more easily, leaving the planet closer to the inner edge of the gap.
The depth of the gap does not strongly depend on f. As long as 1f is not larger than the width of the gap for f=1, the angular momentum flux will be similar in both cases. Then, the ability of a planet to open up a gap is relatively independent of f. However, from the perspective of the planet the case f<1 can be radically different from the Keplerian case. For what accretion and migration are concerned, what matters is the surface density close to the planet. Since the planet is located more and more outside the gap for low values of f, it could accrete and migrate as if there were no gap. This is illustrated in the inset of Fig. 3, where we show the azimuthally averaged surface density at three different locations in the disc for different values of f. In Keplerian discs, planets usually interact with material within 2 scale heights from their orbit (Bate et al. 2003). This is the maximum distance that is considered in the inset of Fig. 3.
For a Keplerian disc, a gap forms with a half width that is approximately 2H, and in the inset of Fig. 3 all three curves fall below 0.1 for f=1. For f=0.8, there is still a density depression around the orbit of the planet, but for lower values of f there is no clear gap near , it has shifted inwards enough so that the planet is basically embedded in the disc again. However, since for the lowest values of f there are no more resonances located close to the planet, the type I torque will not be fully restored. We find that typically for gapopening planets. The influence of the head wind (see Eq. (1)) will be strong, however, since it does not rely on resonances. We comment that a gapopening planet is still tidally locked to the gap, just as in the Keplerian type II migration case.
The possibility of a highmass planet fully embedded in the disc may have some important consequences for gas accretion. In a Keplerian disc, accretion drops by an order of magnitude when a gap is formed (D'Angelo et al. 2003). In a subkeplerian disc, there a significant amount of mass remains near the orbit of the planet, making accretion potentially very efficient. It is not clear, however, if the planet is able to accept material that has such a high relative velocity.
Figure 4: Azimuthally averaged surface density, in units of the initial (constant) surface density , for a q = 10^{4} planet embedded in a h=0.05, disc for different values of f. 

Open with DEXTER 
The results for a more massive planet of q=0.001 (1 around a 1 star) are very similar to those presented in Fig. 3. The gap shifts inward, and for f<0.7 the planet is located on the outer edge of its own density depression. For smaller planets, which only open up a shallow density depression for f=1, remain fully embedded for f<1, as shown in Fig. 4. While for f=1, a q=10^{4} planet decreases the surface density around its orbit by a factor 0.7, for f<0.7 there is no more evidence for any density depression. For this lower planet mass, the important resonances become too weak to affect the surface density in strongly subkeplerian discs.
4 Discussion and conclusions
We have presented hydrodynamical simulations of planets embedded in subkeplerian discs. They represent the first step towards modelling fully 3D magnetised discs. It is known for Keplerian discs that magnetic fields can have a strong impact on planet migration: regular fields introduce magnetic resonances (Terquem 2003), while magnetic turbulence introduces stochastic migration (Nelson & Papaloizou 2004; Adams & Bloch 2009). It remains to be seen what impact a magnetic field configuration that gives rise to subkeplerian discs can have on the simple hydrodynamic picture presented here.
We have worked in the isothermal limit, but since corotation torques only play a minor role in strongly subkeplerian discs, results for more realistic discs should be similar. The twodimensional approximation, in combination with a gravitational softening parameter of order h, gives similar results to fully threedimensional simulations, again at least as far as the Lindblad torque is concerned (Paardekooper & Papaloizou 2009). We have found that the results do not depend strongly on the initial surface density profile.
We have considered migration of lowmass planets (the type I regime), finding that there is a strong dependence on f. For , the planet feels almost the full onesided Lindblad torque, which is a factor 1/h stronger than the classical type I torque. Such a disc would be very hazardous to lowmass planets, since inward migration is sped up by more than an order of magnitude compared to Keplerian discs. For 1f >h, the torque decreases because the resonances the planet interacts with become weaker. The dependence of the torque on h and f is quite complicated, and it is not easy to say for a disc of given f and h whether the type I torque will be stronger of weaker than the headwind torque.
Gap formation proceeds similar to that in Keplerian discs. However, because of the inward shift of the important resonances, the gap will be located inside the planet's orbit. This then leaves the planet on the outside of its own gap. For strongly subkeplerian discs, a gapopening planet can become fully embedded again. Since the torquegenerating resonances are located far away, this does not restore the full type I torque. Accretion time scales could be very short, if the planet is able to accept the available matter.
AcknowledgementsI acknowledge support from STFC in the form of a postdoctoral fellowship. I wish to thank Pawel Ivanov for his interest in retrograde orbiting embedded planets, which formed the start of this project, and the anonymous referee for an insightful report.
References
 Adams, F. C., & Bloch, A. M. 2009, ApJ, 701, 1381 [NASA ADS] [CrossRef]
 Adams, F. C., Cai, M. J., & Lizano, S. 2009, ApJ, 702, L182 [NASA ADS] [CrossRef]
 Bate, M. R., Lubow, S. H., Ogilvie, G. I., et al. 2003, MNRAS, 341, 213 [NASA ADS] [CrossRef]
 Crida, A., Morbidelli, A., & Masset, F. 2006, Icarus, 181, 587 [NASA ADS] [CrossRef]
 D'Angelo, G., Kley, W., & Henning, T. 2003, ApJ, 586, 540 [NASA ADS] [CrossRef]
 Goldreich, P., & Tremaine, S. 1979, ApJ, 233, 857 [NASA ADS] [CrossRef]
 Lin, D. N. C., & Papaloizou, J. 1986, ApJ, 309, 846 [NASA ADS] [CrossRef]
 Masset, F. S., & Papaloizou, J. C. B. 2003, ApJ, 588, 494 [NASA ADS] [CrossRef]
 Nelson, R. P., & Papaloizou, J. C. B. 2004, MNRAS, 350, 849 [NASA ADS] [CrossRef]
 Paardekooper, S.J., & Mellema, G. 2006a, A&A, 459, L17 [NASA ADS] [CrossRef] [EDP Sciences]
 Paardekooper, S.J., & Mellema, G. 2006b, A&A, 450, 1203 [NASA ADS] [CrossRef] [EDP Sciences]
 Paardekooper, S.J., & Mellema, G. 2008, A&A, 478, 245 [NASA ADS] [CrossRef] [EDP Sciences]
 Paardekooper, S.J., & Papaloizou, J. C. B. 2008, A&A, 485, 877 [NASA ADS] [CrossRef] [EDP Sciences]
 Paardekooper, S.J., & Papaloizou, J. C. B. 2009, MNRAS, 394, 2283 [NASA ADS] [CrossRef]
 Papaloizou, J. C. B., Nelson, R. P., Kley, W., et al. 2007, in Protostars and Planets V, ed. B. Reipurth, D. Jewitt, & K. Keil, 655
 Peplinski, A., Artymowicz, P., & Mellema, G. 2008, MNRAS, 386, 179 [NASA ADS] [CrossRef]
 Shu, F. H., Galli, D., Lizano, S., Glassgold, A. E., & Diamond, P. H. 2007, ApJ, 665, 535 [NASA ADS] [CrossRef]
 Terquem, C. E. J. M. L. J. 2003, MNRAS, 341, 1157 [NASA ADS] [CrossRef]
 Ward, W. R. 1997, Icarus, 126, 261 [NASA ADS] [CrossRef]
 Weidenschilling, S. J. 1977, MNRAS, 180, 57 [NASA ADS]
All Figures
Figure 1: Position of the inner (open triangles) and outer (solid circles) Lindblad resonances in discs with varying f. Pressure effects have been ignored in calculating the resonance positions. 

Open with DEXTER  
In the text 
Figure 2: Absolute value of the torque on a planet embedded in discs of constant initial surface density for different values of h and f. Solid curves denote the results of linear calculations, while the symbols denote results from hydrodynamical simulations. 

Open with DEXTER  
In the text 
Figure 3: Azimuthally averaged surface density, in units of the initial (constant) surface density , for a planet embedded in a h=0.05, disc for different values of f. The inset shows the averaged surface density at the orbital radius of the planet as a function of f. 

Open with DEXTER  
In the text 
Figure 4: Azimuthally averaged surface density, in units of the initial (constant) surface density , for a q = 10^{4} planet embedded in a h=0.05, disc for different values of f. 

Open with DEXTER  
In the text 
Copyright ESO 2009