5 Smoothing issues and 3D expectations

The above simulations are performed in 2D, and the disk finite thickness is taken into account through the use of a smoothing coefficient $\epsilon$ in the potential, where $\epsilon$ is a sizable fraction of the disk vertical scale length  $H=r\frac{c_{\rm s}}{v_{\rm K}}$. One can wonder how accurate this description is, and what fraction $\eta =\epsilon /H$ of the disk thickness should be used in order to get as close as possible a result to the 3D expectations. Moreover, there is no reason that an adequate value of $\epsilon$ for the Lindblad torque (in the sense that this value of $\epsilon$ will provide a value for the one-sided and/or differential Lindblad torque in agreement with results obtained with three dimensional calculations, see e.g. Miyoshi et al. 1999) will also be an adequate value for the corotation torque. The following discussion is therefore divided in two parts: the effect of smoothing on the Lindblad torque is first investigated, then secondly on the corotation torque, and lastly it is discussed whether Lindblad and corotation torques can be described appropriately by the same smoothing coefficient.

5.1 Effect of smoothing on the Lindblad torques

Figure 18 shows the ratio of the one-sided Lindblad torque for a smoothed potential and the one-sided Lindblad torque for an unsmoothed potential, as a function of the ratio of the smoothing length to the disk thickness. The one-sided Lindblad torque is evaluated by the arithmetical average of the outer and inner Lindblad torque expressions, which are:

$$\Gamma_{\rm OLR}=\sum_{m=0}^{+\infty}\Gamma_+^m$$ (26)

and:

$$\Gamma_{\rm ILR}=\sum_{m=0}^{+\infty}\Gamma_-^m$$ (27)

where:

$$\Gamma_\varepsilon^m=-\varepsilon\pi m\frac{[\psi_m^\varepsilon(h,\epsilon)]^2h^ 3}{(rdD/dr)_{r_m^\varepsilon}}$$ (28)

in which $\varepsilon=\pm 1$ (for Outer/Inner Lindblad resonance respectively), and where:

$$\psi_m^\varepsilon(h,\epsilon) = \frac{\left(r\frac{\rm d}{{... ...\varepsilon}G_m^\epsilon(r_m^\varepsilon)}{(1 +4m^2h^2)^{1/2}}$$ (29)

in which $G_m^\epsilon(r)$ is a generalized Laplace coefficient (i.e. a Laplace coefficient for which a smoothing is introduced):

$$G_m^\epsilon(r)=\frac2\pi\int_0^\pi\frac{\cos(m\theta)}{(1-2r\cos\theta+r^2+\epsilon^2)^{1/2}} {\rm d}\theta$$ (30)

and where:

$$\Omega_m^\varepsilon=\Omega_{\rm p}[(1-h^2)^{1/2}+\varepsilon m^{-1}(1+m^2h^2)^{1/2}]^{-1}$$ (31)

is the Keplerian frequency at the effective location of the Lindblad resonances, and:

$$r_m^\varepsilon=\left(\frac{\Omega_m^\varepsilon}{\Omega_{\rm p}}\right)^{-2/3}$$ (32)

and we use:

$$\left(r\frac{{\rm d}D}{{\rm d}r}\right)_{r_m^\varepsilon}=-3m\varepsilon\Omega_m^\varepsilon \Omega_{\rm p}(1+m^2h^2)^{1/2}.$$ (33)

This calculation is identical to the calculation in Ward (1997) except for the introduction of the smoothing of the perturbing potential.

Figure 18: One-sided Lindblad torque cut-off as a function of the relative smoothing $\epsilon /H$. The long dashed line represents the approximate asymptotic solution for $h\rightarrow 0$ deduced from Eq. (34).

The cut-off function is shown for three different aspect ratios, and found to be only weakly dependent on h. At low aspect ratio, the functional dependence of the torque cut-off tends towards a fixed function, the graph of which should roughly coincide with the solid curve of Fig. 18. An approximation of the limit cut-off function when $h\rightarrow 0$ can be obtained by developing Eq. (30) to first order in mh, and approximating the $G_m^\epsilon$ coefficients with a standard technique as a function of the Bessel K₀ and K₁ functions. The summation over m can then be approximated as an integral. Using the variables $\xi=mh$ and $\eta =\epsilon /H$, one can finally write the following approximate formula for the one-sided Lindblad torque as a function of the relative smoothing $\eta$ only:

$$\gamma_{\rm LR}(\eta)=K\int_0^\infty\xi^2\frac{(1+\xi^2)^{1/... ...4\xi^2} F\left[\frac 49(1+\xi^2)+\xi^2\eta^2\right]{\rm d}\xi~$$ (34)

where K is a constant and:

$$F(x)=\frac{1}{3\sqrt{x}}{K}_1(\sqrt{x})+{K}_0(\sqrt{x})$$ (35)

where K₀and K₁ are the modified Bessel functions. The long dashed line of Fig. 18 is the graph of the cut-off $\gamma_{\rm LR}(\eta)/\gamma_{\rm LR}(0)$. Figure 19 shows the behavior of the differential Lindblad torque cut-off for the same set of disk thicknesses. This cut-off exhibits quite similar a behavior to the one-sided Lindblad torque cut-off.

Figure 19: Differential Lindblad torque cut-off as a function of the relative smoothing $\epsilon /H$.

In Figs. 18 and 19 the dot-dashed line shows the value of the smoothing for which the torque cut-off amounts to 0.43, which is the value quoted by Miyoshi et al. (1999) for the ratio of the Lindblad torque in a vertically resolved disk and in an infinitesimally thin disk, in the linear regime. The corresponding relative smoothing is $\eta=0.76$, both on the one-sided Lindblad torque and on the differential Lindblad torque. Therefore a potential smoothing length of 76% of the disk thickness should be used in 2D simulations in order to give correct results for the (differential) Lindblad torque. The fact that this value is independent of the disk aspect ratio and planet mass (at least in the linear regime) comes from the fact that most of the torque between the disk and the planet comes from a zone which is at a distance $\pm H$ from the orbit (the Lindblad resonances pile up for $m\rightarrow\infty$at a distance $\pm \frac 23H$ from the orbit). Such an argument cannot be used for the co-orbital corotation torque which needs a different treatment.

5.2 Effect of smoothing on the corotation torque

The corotation torque comes from the exchange of angular momentum between the planet and nearby, orbit crossing fluid elements (which can be either fluid elements librating on closed, horseshoe like streamlines, or fluid elements participating in the global viscous accretion of the disk, and which pass by the planet as they go from the outer disk to the inner disk). As a fluid element participating in the corotation torque, with impact parameter $x=r-r_{\rm p}$, will be located at a distance -x from the orbit after a back-scattering by the planet, the specific angular momentum that it gives to the planet is $4Bxr_{\rm p}$, and therefore this value does not depend on smoothing. The total value of the corotation torque does depend on it however, since the radial width of the libration region depends on the smoothing: the smaller the smoothing, the further the separatrices lie from the orbit, and the larger the corotation torque.

Figure 20: Total torque acting on the planet, as a function of the reduced viscosity, for runs $\mbox{S20R}17_5^i$, $\mbox{S40R}17_5^i$ and R175i.

This behavior can be checked at Fig. 20. The runs for a planet of mass $m_{\rm p}=16.7$ $M_\oplus$have been performed with three different smoothing values. At low viscosity (i.e. for a saturated corotation torque) the total torque is negative and has a larger absolute value for a smaller smoothing. This corresponds to expectations since the torque in this region (whether it includes the Lindblad torque coupling term of the corotation torque or not) scales as the (one-sided or differential) Lindblad torque. On the other hand, the difference between the maximum value of the torque with the minimum value at low viscosity should roughly correspond to the maximum value of the corotation torque main term, which scales as $x_{\rm s}^4$. At Fig. 21 we plot this difference as a function of the separatrix distance to the orbit $x_{\rm s}$.

Figure 21: Corotation torque main term estimate as a function of the separatrix position. The dotted line correspond to a fourth power of $x_{\rm s}$, and the dashed one to a third power of $x_{\rm s}$.

Despite the small number of measurements, the results are in good agreement with a $x_{\rm s}^4$ dependency of the corotation torque. This also confirms the argument exposed above that the corotation torque depends on the smoothing only through the value of $x_{\rm s}$.

It is worth emphasizing the dramatic dependence of the corotation torque on smoothing. In Fig. 20, for a smoothing length that is 60% of the disk thickness, the torque is always negative, although it can reach a value one order of magnitude smaller than the differential Lindblad torque linear estimate for an unsmoothed potential. On the other hand, a smoothing length that is 20% of the disk thickness leads to a torque reversal over a significant range of viscosity ( $1.5\times10^{-5} < \hat\nu < 6\times 10^{-5}$, which corresponds to $6\times 10^{-3}<\alpha<0.024$).

An idea of the co-orbital dynamics in a three dimensional situation can be obtained simply if one assumes that the motion is purely horizontal. In that case the motion in the slice of altitude z and infinitely small thickness dz is equivalent to a 2D situation, in which the potential smoothing length is z. Furthermore, as for $\vert z/r_{\rm p}\vert \ll 1$ the specific angular momentum radial gradient is still $2Br=\Omega_{\rm K}r/2$, the torque contribution of the slice is the same as the torque exerted by an infinitesimally thin disk (with a potential smoothed over a length z) of surface density $\rho(z){\rm d}z$. Let $X_{\rm s}(\epsilon, q, h)$ be the separatrix distance to the orbit in an infinitesimally thin disk in the equatorial plane, in which the potential is smoothed over a length $\epsilon$, the sound speed is $c_{\rm s}=v_{\rm K}h$, and the planet to primary mass ratio is q. The fully unsaturated corotation torque main term in a thick disk can therefore be written as:

$$\Gamma_{\rm c}=\int_{-\infty}^{+\infty}\frac 98\Omega_{\rm p... ...m s}^4\left[z,q,\frac{c_{\rm s}(z)}{v_{\rm K}}\right]{\rm d}z.$$ (36)

As this analysis is restricted to the case of isothermal Keplerian disks, Eq. (36) can be recast as:

$$\Gamma_{\rm c}=\int_{-\infty}^{+\infty}\frac 98\Omega_{\rm p}^2\rho_0{\rm e}^{-\frac 12(z/H)^2}X_{\rm s}^4(z,q,h){\rm d}z$$ (37)

where $\rho_0$ is the disk density in the equatorial plane. Eq. (36) can also be written as:

$$\Gamma_{\rm c}=\frac 98\Omega_{\rm p}^2\Sigma_0\overline{x}_{\rm s}^4$$ (38)

in which $\overline x_{\rm s}^4$ is defined as:

$$\overline x_{\rm s}^4=\int_{\infty}^{+\infty}\frac{\rho(z)X_{\rm s}^4(z,q,h)}{\Sigma_0}{\rm d}z$$ (39)

which can be recast as:

$$\overline x_{\rm s}^4=\frac{1}{\sqrt{2\pi}h}\int_{\infty}^{+\infty}{\rm e}^{-\frac 12(z/H)^2}X_{\rm s}^4(z,q,h){\rm d}z.$$ (40)

Therefore, if one wants a 2D simulation to give the correct amplitude for the fully unsaturated corotation torque main term, one needs to use the smoothing length which endows the co-orbital motion with a separatrix position $\overline x_{\rm s}$, where $\overline x_{\rm s}$ is given by Eq. (40). The integrand of the right hand side contains the function $X_{\rm s}$, which has to be tabulated a priori with a set of runs with different smoothing lengths, for a given planet mass and disk aspect ratio. The same tabulation is used to infer the correct value for $\epsilon$ once $\overline x_{\rm s}$ is known.

Figure 22: Separatrix distance to the orbit as a function of the relative smoothing $\eta =\epsilon /H$, for different aspect ratios and planet masses ( $q=1.67\times 10^{-5}$: dashed lines, $q=3.33\times 10^{-5}$: dotted lines, and $q=5\times 10^{-5}$: solid lines). As expected, $x_{\rm s}$ decreases as the smoothing increases, decreases as the aspect ratio increases (everything else being fixed), and increases as the mass increases (everything else being fixed).

Figure 22 shows such tabulations for different planet masses and disk aspect ratios. They are measured with a dichotomic search of the separatrix position on the results of low viscosity runs ( $\hat \nu=10^{-6}$), over 150 orbits.

Figure 23: Corotation torque adapted smoothing as a function of $R_{\rm H}/H$, inferred from the tabulations of Fig. 22.

Figure 23 shows the results of this integration for the nine situations of Fig. 22. It can be seen that on the contrary to the case of the Lindblad torque, the value of the correct smoothing depends on $R_{\rm H}/H$ even when this value is below 1/2. Furthermore, in this case, which corresponds to the linear regime for the Lindblad torque, the correct smoothing length is clearly smaller (50-60% of the disk thickness) than the correct smoothing length for the Lindblad torque (76% of the disk thickness). The choice of any value between 60 and 76% of the disk thickness will therefore underestimate the corotation torque and overestimate the Lindblad torque, and therefore in any case will underestimate the total torque, since the differential Lindblad torque is negative and the corotation torque is positive. Furthermore, as the corotation torque is a very sensitive function of the smoothing, there is little hope of finding a correct smoothing prescription which predicts a correct total torque value in a 2D simulation, especially in the cases of interest where the corotation torque and differential Lindblad torque almost cancel out. One can however make conservative assumptions to infer the direction of migration of a protoplanet by choosing a lower or upper limit for the smoothing coefficient. This is the object of the next section which presents four series of selected high resolution runs.