text:
# Cosmological Equations for a Thick Brane

## I Introduction

The idea of extra-dimensions where ordinary matter is confined within a lower dimensional submanifold has received an enormous amount of attention during the last few years (see e.g.  for a recent review). This idea has been explored, in particular, in the cosmological context.

####

In the case of a single extra-dimension, a simplifying assumption has been to consider our brane-universe as infinitely thin along this extra-dimension. It has been shown that this approximation is quite valid as long as the energy scales one wishes to explore are much smaller that the energy scale related to the inverse thickness of the brane.

####

The purpose of the present paper is to explore the modification of the equations yielding the cosmological evolution within the brane when the finite thickness of the brane is taken into account. We consider this question in a purely phenomenological approach by assuming that the matter energy-momentum tensor describing the brane has some distribution over a finite length interval along the fifth dimension. Our approach (similar to that of ) is thus more crude, but also more general, than the domain wall description of a brane, for which the brane is embodied by a scalar field, such as in  (see also  for a recent study in 4d gravity and  for a review on domain walls, in particular in supergravity).

####

Two limiting cases are already known. The case of an infinitely thin brane gives the unconventional Friedmann equation of , which can be made compatible with ordinary cosmology at late times if one introduces both a negative cosmological constant in the bulk spacetime and a related tension (unmeasurable by cosmological observers) in the brane . The other and opposite limit is that of an infinitely thick brane, in a sense defined below. This limit effectively corresponds to a Kaluza-Klein picture where matter is homogeneously distributed over the extra-dimension(s) and, as it has been known for a long time, ordinary Friedmann equations are recovered. Our analysis gives access to a description of the intermediate cases.

####

The paper is organized as follows. Our effective Friedmann equation is established in Section II and written in terms of four-dimensional effective quantities. Section III is devoted to the cosmology for a brane with a tension, the cosmological effect of which is cancelled by a negative bulk cosmological constant. We consider the limits of a low brane energy (for all brane thickness) and of a small brane thickness (for all brane energy), and we propose uniform expressions for the cosmological equations valid over the whole parameter domain covered by these two limits. In Section IV, we determine the higher order corrections (in the normalized brane thickness) to the effective cosmological equations in the simpler case of a Minkowski bulk. In particular, we show that no effective four-dimensional reduction is possible in the regimes of order one normalized brane thickness.

## II Effective Friedmann equation for a thick brane

We consider a five-dimensional spacetime, which is homogeneous and isotropic along three spatial dimensions, and which contains a thick brane. For simplicity, we will restrict ourselves to the case of a mirror symmetric brane (and thus spacetime) along the extra-dimension. It is then always possible to find a Gaussian coordinate system, starting from the hypersurface representing the “center” of the brane, in which the line element takes the simple diagonal form

$$ds^2g_{AB}dx^Adx^B=n^2(t,y)dt^2+a^2(t,y)\delta _{ij}dx^idx^j+r_b^2dy^2,$$
(1)
where $`y`$ is the coordinate of the fifth dimension and $`r_b`$ is the brane thickness which will be assumed to be time independent. The brane is localized between $`y=1/2`$ and $`y=1/2`$. To be compatible with the spacetime symmetries, the energy-momentum tensor of the matter content in the brane is necessarily of the form

$$T_{}^{A}{}_{B}{}^{}=\mathrm{diag}[\frac{1}{r_b}(\rho _b,p_b,p_b,p_b),P_T],$$
(2)
where $`\rho _b`$, $`p_b`$, and $`P_T`$ are functions of $`t`$ and $`y`$. The presence of a negative cosmological constant in the bulk, $`\mathrm{\Lambda }`$, is accounted for by an energy-momentum tensor of the form

$$T_{}^{A}{}_{B}{}^{}|_{\mathrm{Bulk}}=\mathrm{\Lambda }g_{}^{A}{}_{B}{}^{}.$$
(3)
Note that, strictly speaking, $`\mathrm{\Lambda }`$ as defined here does not have the dimension of a cosmological constant but rather that of a five-dimensional energy density.

####

The purpose of this paper is to establish effective cosmological equations for an observer living in the brane. Due to the finite thickness of the brane, there is some arbitrariness in the definition of what the effective four-dimensional quantities should be. The simplest prescription one can think of consists in defining the four-dimensional effective quantity associated to a five-dimensional quantity $`Q(t,y)`$ as its spatial average over the brane thickness

$$Q(t)=_{1/2}^{1/2}Q(t,y)𝑑y.$$
(4)
In this paper, we will adopt this simple prescription. The four-dimensional “observable” counterparts of $`a`$, $`\rho _b`$, and $`p_b`$ are thus $`a=_{1/2}^{1/2}a𝑑y`$, $`\rho _b=_{1/2}^{1/2}\rho _b𝑑y`$, and $`p_b=_{1/2}^{1/2}p_b𝑑y`$. Note that, in contrast with the usual dimensional reduction in field theory, we do not integrate over the whole (compact) extra-dimension but only over the extension of the brane in the extra-dimension.

####

The next step is now to obtain evolution equations for the four-dimensional effective quantities from the five-dimensional Einstein’s equations

$$G_{AB}R_{AB}\frac{1}{2}Rg_{AB}=\kappa _{(5)}^2T_{AB},$$
(5)
where $`R_{AB}`$ is the Ricci tensor and $`R`$ its trace. The components of the Einstein equations, with the metric (1) and the energy-momentum tensors (2) and (3) can be obtained from the components of the Einstein tensor given in the appendix. In particular, the 05 component of the Einstein equations yields

$$n(t,y)=\xi (t)\dot{a}(t,y),$$
(6)
where the dot stands for a partial derivative with respect to time $`t`$ and $`\xi `$ depends on the normalization prescription for $`n`$. In the following, we will take the normalization $`n=1`$ which gives $`\xi =\dot{a}^1`$. The 00 component of the Einstein equations then reads

$$(a^2)^{\prime \prime }=\frac{2\kappa _{(5)}^2}{3}\left(r_b\rho _b+r_b^2\mathrm{\Lambda }\right)a^2+2r_b^2\dot{a}^2,$$
(7)
where a prime denotes a partial derivative with respect to the coordinate $`y`$. Integrating this equation over the brane, one obtains

$`a(1/2)^{}`$ $`=`$ $`{\displaystyle \frac{1}{a(1/2)}}\left[{\displaystyle \frac{\kappa _{(5)}^2}{6}}r_b\rho _ba^2{\displaystyle \frac{\kappa _{(5)}^2}{6}}r_b^2\mathrm{\Lambda }a^2+{\displaystyle \frac{1}{2}}r_b^2\dot{a}^2\right]`$ (8)
$`=`$ $`{\displaystyle \frac{a}{\alpha }}\left[\epsilon \eta {\displaystyle \frac{r_b^2}{\mathrm{}_\mathrm{\Lambda }^2}}\stackrel{~}{\eta }+{\displaystyle \frac{1}{2}}r_b^2H_a^2\right],`$ (9)
where we have introduced the effective four-dimensional Hubble parameter

$$H_a\dot{a}/a,$$
(10)
the Anti-de Sitter lengthscale associated to the (negative) cosmological constant in the bulk

$$\mathrm{}_\mathrm{\Lambda }\sqrt{\frac{6}{\kappa _{(5)}^2\mathrm{\Lambda }}},$$
(11)
and the dimensionless quantities $`\epsilon `$, $`\alpha `$, $`\eta `$, and $`\stackrel{~}{\eta }`$, respectively defined by

$$\epsilon \frac{\kappa _{(5)}^2}{6}r_b\rho _b,$$
(12)
and

$$\alpha \frac{a(1/2)}{a},\eta \frac{\rho _ba^2}{\rho _ba^2},\stackrel{~}{\eta }\frac{a^2}{a^2}.$$
(13)
Whereas $`\epsilon `$ characterizes the thickness of the brane, the quantities $`\alpha `$, $`\eta `$, and $`\stackrel{~}{\eta }`$ characterize the inhomogeneity of the brane along the fifth dimension (in the case of a homogeneous brane, one has $`\alpha =\eta =\stackrel{~}{\eta }=1`$). The 55 component of the Einstein equations can be written as

$$\dot{F}=\frac{2}{3}\kappa _{(5)}^2\dot{a}a^3P_T,$$
(14)
with $`F=a^2[(a^{}/r_b)^2\dot{a}^2+\kappa _{(5)}^2a^2\mathrm{\Lambda }/6]`$. Imposing the boundary condition $`P_T(\pm 1/2)=0`$, one has $`\dot{F}(\pm 1/2)=0`$ which gives, after time-integration,

$$r_b^2H_a^2=\left[\frac{a^{}(1/2)}{a}\right]^2+\frac{𝒞r_b^2}{\alpha ^2a^4}+\frac{\kappa _{(5)}^2}{6}\alpha ^2r_b^2\mathrm{\Lambda },$$
(15)
where $`𝒞`$ is a constant of integration. Injecting Eq. (8) into Eq. (15), one obtains

$$\alpha ^2r_b^2H_a^2=\left[\epsilon \eta +\frac{r_b^2}{\mathrm{}_\mathrm{\Lambda }^2}\stackrel{~}{\eta }+\frac{1}{2}r_b^2H_a^2\right]^2+\frac{𝒞r_b^2}{a^4}\frac{r_b^2}{\mathrm{}_\mathrm{\Lambda }^2}\alpha ^4.$$
(16)
This equation is one of the main results of this paper. It provides the generalization of the (first) Friedmann equation relating the effective four-dimensional Hubble parameter to the effective brane energy density and the brane thickness. It of course also includes some dependence on the exact profile of the energy-momentum tensor but this dependence has been conveniently reduced to the three dimensionless quantities (13), which in typical cases are of order unity. From Eqs. (7) and (16), one finds that a necessary condition for the brane to be homogeneous along the fifth dimension reads

$$\frac{𝒞}{a^4}=\frac{2\epsilon }{r_b^2}\frac{1}{\mathrm{}_\mathrm{\Lambda }^2},$$
(17)
which, in general, cannot be fulfilled uniformly in time unless the energy-momentum content of the brane reduces to a mere brane tension (as e.g. in the original static Randall-Sundrum model ). It follows that, in the cosmological context, a brane with a finite thickness cannot remain homogeneous along the fifth dimension.

####

Equation (16) can be rewritten in the form

$`H_a^2={\displaystyle \frac{2}{r_b^2}}\left(\alpha ^2+\epsilon \eta {\displaystyle \frac{r_b^2}{\mathrm{}_\mathrm{\Lambda }^2}}\stackrel{~}{\eta }\right)`$ (18)
$`\times \left[1\pm \sqrt{1{\displaystyle \frac{(\epsilon \eta \stackrel{~}{\eta }r_b^2/\mathrm{}_\mathrm{\Lambda }^2)^2+𝒞r_b^2/a^4\alpha ^4r_b^2/\mathrm{}_\mathrm{\Lambda }^2}{(\alpha ^2+\epsilon \eta \stackrel{~}{\eta }r_b^2/\mathrm{}_\mathrm{\Lambda }^2)^2}}}\right],`$ (19)
which admits real solutions if the arbitrary constant $`𝒞`$ fulfills the condition

$$\frac{𝒞r_b^2}{a^4}2\alpha ^2\left(\epsilon \eta \frac{r_b^2}{\mathrm{}_\mathrm{\Lambda }^2}\stackrel{~}{\eta }\right)+\alpha ^4\left(1+\frac{r_b^2}{\mathrm{}_\mathrm{\Lambda }^2}\right).$$
(20)
In the following, we will consider Expression (18) with the minus sign only, which corresponds to a well defined cosmology in the $`r_b0`$ limit. It is important to notice that this equation should be regarded as an implicit equation for $`H_a`$ because $`\alpha `$, $`\eta `$, and $`\stackrel{~}{\eta }`$ can depend on $`H_a`$.

####

To conclude this section, it can be instructive to consider the two opposite limits of small and large $`r_b`$. In the limit $`r_b0`$, Equation (18) reduces to

$`H_a^2`$ $``$ $`{\displaystyle \frac{1}{\alpha ^2}}\left({\displaystyle \frac{\epsilon ^2\eta ^2}{r_b^2}}+{\displaystyle \frac{𝒞}{a^4}}{\displaystyle \frac{\alpha ^4}{\mathrm{}_\mathrm{\Lambda }^2}}\right)`$
$`=`$ $`{\displaystyle \frac{1}{\alpha ^2}}\left({\displaystyle \frac{\kappa _{(5)}^2}{36}}\rho _b^2\eta ^2+{\displaystyle \frac{𝒞}{a^4}}+{\displaystyle \frac{\kappa _{(5)}^2}{6}}\mathrm{\Lambda }\alpha ^4\right).`$
It will be shown in the following that, in this limit, both $`\alpha `$ and $`\eta `$ tend to $`1`$. Thus, one exactly recovers the unconventional Friedmann equation of thin brane cosmology, where the brane energy density enters quadratically . In the opposite limit $`r_b+\mathrm{}`$, Equation (18) reduces to

$$H_a^2\frac{2\epsilon \eta }{r_b^2}\frac{2\stackrel{~}{\eta }}{\mathrm{}_\mathrm{\Lambda }^2}+\frac{2}{r_b^2}\sqrt{\frac{\alpha ^2(\alpha ^22\stackrel{~}{\eta })}{\mathrm{}_\mathrm{\Lambda }^2}\frac{𝒞}{a^4}},$$

where terms up to $`O(r_b^1)`$ have been kept. In the case of a Minkowski bulk, $`\mathrm{}_\mathrm{\Lambda }^2=0`$, Eq. (20) yields $`𝒞/a^40`$ as $`r_b+\mathrm{}`$ and one obtains

$$H_a^2\frac{2\epsilon \eta }{r_b^2}=\frac{\kappa _{(5)}^2\eta }{3r_b}\rho _b,$$

which corresponds to the standard Friedmann equation, in which the energy density enters linearly with the effective four-dimensional gravitational constant $`\kappa _{(4)}^2\kappa _{(5)}^2\eta /r_b`$. Note that, in this limit, both sides of the necessary condition (17) vanish, which allows the possibility of a homogeneous brane along the fifth dimension ($`\eta =1`$), in agreement with the usual Kaluza-Klein picture.

## III Randall-Sundrum type cosmology of a thick brane

In thin brane cosmology, the unconventional Friedmann equation can be made compatible with the current cosmological observations if one assumes that the brane energy density consists of ordinary cosmological matter living in the brane and of a tension $`\lambda `$ adjusted so as to compensate the negative cosmological constant in the bulk . This adjustment corresponds to the fine-tuning condition of Randall-Sundrum , which can be reexpressed as the simple relation

$$\mathrm{}_\lambda =\mathrm{}_\mathrm{\Lambda },$$
(21)
where the lengthscale $`\mathrm{}_\lambda `$ is defined from the brane tension $`\lambda `$ as

$$\mathrm{}_\lambda \frac{6}{\kappa _{(5)}^2\lambda }.$$
(22)
In this section, we will also assume that the energy-momentum content of the brane can be decomposed into two parts,

$`\rho _b(t,y)=\lambda +\rho _m(t,y),`$ (24)
$`p_b(t,y)=\lambda +p_m(t,y),`$ (25)
where $`\lambda `$, which is strictly constant, corresponds to a tension of the brane. For the sake of simplicity, we assume here that $`\lambda `$ is independent of $`y`$. The case of a non-trivial profile of the brane tension along the fifth dimension could also be easily considered by writing $`\lambda =\lambda +\delta \lambda (y)`$ and by defining an effective “matter” energy-momentum content $`\stackrel{~}{\rho }_m(t,y)=\rho _m(t,y)+\delta \lambda (y)`$ and $`\stackrel{~}{p}_m(t,y)=p_m(t,y)\delta \lambda (y)`$. Since $`\delta \lambda (y)=0`$, it can be seen that the only modification induced by such a brane tension profile reduces to adding the geometrical term $`\delta \lambda a^2/(\lambda a^2)`$ to the parameter $`\epsilon _\rho `$ defined in Eq. (28). We will not analyse further on the effects of this extra term and restrict ourselves to the case of a $`y`$-independent brane tension. Injecting Eq. (24) into the effective Friedmann equation (18), it can be easily seen that the generalization of the Randall-Sundrum fine-tuning allowing for a finite brane thickness is given by

$$\mathrm{}_\lambda =\frac{\stackrel{~}{\eta }(1u)}{\alpha ^2}\mathrm{}_\mathrm{\Lambda },$$
(26)
where we have introduced the dimensionless ratio

$$u\frac{r_b\mathrm{}_\lambda }{\mathrm{}_\mathrm{\Lambda }^2}=\frac{r_b\mathrm{\Lambda }}{\lambda }.$$
(27)
In the limit where $`r_b`$ goes to zero, it can be checked that Eq. (26) reduces to the Randall-Sundrum condition (21), as it should be. From Eqs. (18), (24), and (26), and introducing the dimensionless parameters

$$\epsilon _\lambda \frac{r_b}{\mathrm{}_\lambda },\epsilon _\rho \frac{\eta \rho _m}{\lambda },$$
(28)
where $`\eta \rho _ma^2/(\rho _ma^2)`$, one finds that Eq. (18) (with the minus sign) can be rewritten as

$`H_a^2={\displaystyle \frac{2}{r_b^2}}\left[\alpha ^2(1+\epsilon _\lambda \mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })+\epsilon _\lambda \epsilon _\rho \right]`$ (29)
$`\times \left\{1\sqrt{1{\displaystyle \frac{\epsilon _\lambda ^2(2\epsilon _\rho \alpha ^2\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda }+\epsilon _\rho ^2+𝒞\mathrm{}_\lambda ^2/a^4)}{[\alpha ^2(1+\epsilon _\lambda \mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })+\epsilon _\lambda \epsilon _\rho ]^2}}}\right\}.`$ (30)
### A Cosmological equations in the low energy limit

In this subsection, we will consider the low energy regime in which the ordinary cosmological energy (and pressure) is small with respect to the brane tension. For the sake of simplicity, we will take $`𝒞=0`$. Consider the low energy limit in the regime defined by $`\epsilon _\rho \mathrm{min}(1,1/\epsilon _\lambda )`$. Then, at lowest order in $`\epsilon _\rho `$, Eq. (29) reduces to

$$r_b^2H_a^2=\frac{2\epsilon _\lambda ^2\epsilon _\rho }{\epsilon _\lambda +\mathrm{}_\mathrm{\Lambda }/\mathrm{}_\lambda },$$
(31)
and the quantities $`\alpha `$, $`\eta `$, and $`\stackrel{~}{\eta }`$ must be determined at zeroth order in $`\epsilon _\rho `$. In this limit, the 00 component of the Einstein equations \[see Eq. (7)\] reduces to

$$\left(\frac{a^2}{a^2}\right)^{\prime \prime }=4\epsilon _\lambda (1u)\frac{a^2}{a^2},$$
(32)
where we have neglected the term $`2r_b^2H_a^2`$ which is first order in $`\epsilon _\rho `$. Restricting ourselves to the case of physical interest $`(1u)>0`$ \[it can be checked a posteriori, cf. Eq. (42), that this inequality is always fulfilled\], and integrating Eq. (32) with the boundary condition $`a^{}(y=0)=0`$, one finds

$$a=a[𝒜\mathrm{cos}(ky)]^{1/2},$$
(33)
with

$$k^2=4\epsilon _\lambda (1u),$$
(34)
and

$$𝒜=\mathrm{cos}^{1/2}(ky)^2.$$
(35)
Note that here we necessarily have $`k<\pi `$, and the scale factor $`a`$ is always positive throughout the brane, reaching the extremity of the brane at $`y=1/2`$. It is possible to determine the metric outside the brane by matching, at $`y=1/2`$, the scale factor and its derivative with respect to $`y`$ to the general solution for the bulk with cosmological constant, obtained in , which is in fact the Schwarschild-AdS solution expressed in Gaussian normal coordinates.

####

These expressions can be used to determine the evolution of $`\rho _m`$ at lowest order in $`\epsilon _\rho `$. Since neither $`𝒜`$ nor $`k`$ depend on $`t`$, one simply has

$$\dot{\rho }_m+3H_a(\rho _m+p_m)=0,$$
(36)
from which it follows that, at this order, $`\rho _m`$ factorizes as $`\rho _m(y,t)=f(y)\rho _m(t)`$.

####

The effective four-dimensional brane cosmology in the limit $`\epsilon _\rho \mathrm{min}(1,1/\epsilon _\lambda )`$ is thus given by

$`H_a^2={\displaystyle \frac{\kappa _{(5)}^2}{3}}\left({\displaystyle \frac{\eta }{r_b+\mathrm{}_\mathrm{\Lambda }}}\right)\rho _m,`$ (38)
$`\dot{\rho _m}+3H_a(\rho _m+p_m)=0,`$ (39)
with

$$\eta =\frac{f(y)\mathrm{cos}(ky)}{f(y)\mathrm{cos}^{1/2}(ky)^2},$$

where $`k`$ and $`𝒜`$ are respectively given by Eqs. (34) and (35). Using the fine-tuning condition (26) with $`\alpha =[𝒜\mathrm{cos}(k/2)]^{1/2}`$ and $`\stackrel{~}{\eta }=(2𝒜/k)\mathrm{sin}(k/2)`$, one can reexpress the constant $`k`$ as a function of $`\epsilon _\lambda `$ in the parametric form

$`k=2\mathrm{tan}^1\sqrt{u/(1u)},`$ (41)
$`\epsilon _\lambda ={\displaystyle \frac{[\mathrm{tan}^1\sqrt{u/(1u)}]^2}{1u}},`$ (42)
with $`0u<1`$. Note that since, in this limit, $`\eta `$ does not depend on $`t`$, the only effect of the finite brane thickness is a modification of the relation between $`\kappa _{(5)}`$ and $`\kappa _{(4)}`$. Namely, one has

$$\kappa _{(4)}^2=\frac{\eta \kappa _{(5)}^2}{r_b+\mathrm{}_\mathrm{\Lambda }}.$$
(43)
Taking into account the constraint (26), one has the limits $`\kappa _{(4)}^2\kappa _{(5)}^3\sqrt{\mathrm{\Lambda }/6}\kappa _{(5)}^4\lambda /6`$ for $`r_b/\mathrm{}_\mathrm{\Lambda }1`$, and $`\kappa _{(4)}^2\eta \kappa _{(5)}^2/r_b`$ for $`r_b/\mathrm{}_\mathrm{\Lambda }1`$. In the simplest case $`f(y)=\mathrm{cst}.`$, $`\rho _m=\rho _m`$, and $`lim_{r_b+\mathrm{}}\eta 1.094`$. In Fig. 1 we have plotted $`k`$ as a function of $`u`$ (solid line), and $`\epsilon _\lambda `$ as a function of $`u`$ (dashed-dotted line). The infinitely thin limit corresponds to $`\epsilon _\lambda 0`$, in which case $`u0`$ and $`k0`$. The opposite limit $`\epsilon _\lambda +\mathrm{}`$, yields $`u1`$ and $`kk_{max}=\pi `$. Figure 2 shows the dimensionless coefficients, $`𝒜`$, $`\alpha `$ and $`\stackrel{~}{\eta }`$, characterizing the brane profile, as a function of $`k/k_{max}`$. In the infinitely thin brane limit ($`k0`$), these three coefficients tend to $`1`$. It can be seen that $`\stackrel{~}{\eta }`$ stays almost constant over the full range of $`k`$, with $`\stackrel{~}{\eta }(k_{max})1.094`$. As the brane thickness increases, $`\alpha `$ decreases down to zero for $`k=k_{max}`$ and $`𝒜`$ increases up to the value $`𝒜(k_{max})1.719`$.

####

All the results of this section have been obtained in the limit $`\epsilon _\rho \mathrm{min}(1,1/\epsilon _\lambda )`$ by keeping terms explicitely of lowest order in $`\epsilon _\rho `$ and of any order in $`\epsilon _\lambda `$. The actual small parameter appearing in the expansion of the square root on the right-hand side of Eq. (29) is $`\epsilon _\lambda ^2\epsilon _\rho (\alpha ^2\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })\epsilon _\lambda ^2\epsilon _\rho \mathrm{min}(1,1/\epsilon _\lambda )`$. These equations are thus correct, in the limit $`\epsilon _\rho \mathrm{min}(1,1/\epsilon _\lambda )`$, at lowest order in $`\epsilon _\lambda ^2\epsilon _\rho \mathrm{min}(1,1/\epsilon _\lambda )`$ and at any order in $`\epsilon _\lambda `$ and $`\epsilon _\mathrm{\Lambda }r_b/\mathrm{}_\mathrm{\Lambda }`$.

### B Cosmological equations for a small brane thickness

In the limit of a small brane thickness, $`\epsilon _\lambda 1`$, Eq. (29) with $`𝒞=0`$ reduces to

$$r_b^2H_a^2=\frac{\epsilon _\lambda ^2(2\epsilon _\rho \alpha ^2\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda }+\epsilon _\rho ^2)}{\alpha ^2(1+\epsilon _\lambda \mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })+\epsilon _\lambda \epsilon _\rho },$$
(44)
which is correct at third order in $`\epsilon _\lambda `$. At this order, $`\alpha `$ must be expressed at first order, and $`\eta `$ and $`\stackrel{~}{\eta }`$ at zeroth order, \[note that $`\eta `$ and $`\stackrel{~}{\eta }`$ behave as $`1+O(\epsilon _\lambda ^2)`$\]. To determine $`a`$ perturbatively in $`\epsilon _\lambda `$, we will use the following $`y`$-expansion

$`a(t,y)=a(t){\displaystyle \underset{n=0}{\overset{+\mathrm{}}{}}}\overline{a}_n(t)y^{2n}a(t)\mathrm{\Sigma }_a(t,y),`$ (46)
$`\rho _m(t,y)=\rho _m(t){\displaystyle \underset{n=0}{\overset{+\mathrm{}}{}}}\overline{\rho }_n(t)y^{2n}\rho _m(t)\mathrm{\Sigma }_\rho (t,y).`$ (47)
Inserting Eqs. (44) in Eq. (7), using the constraint (26), and neglecting the $`O(\epsilon _\lambda ^2)`$ terms, one obtains

$$(\mathrm{\Sigma }_a^2)^{\prime \prime }=4\epsilon _\lambda \left(\epsilon _\rho \mathrm{\Sigma }_\rho +\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda }\right)\mathrm{\Sigma }_a^2.$$
(48)
Since one needs $`a`$ at first order in $`\epsilon _\lambda `$, one can consider the constant (in $`y`$) component of Eq. (48) only. One finds $`\overline{a}_1=\epsilon _\lambda \overline{a}_0(\epsilon _\rho \overline{\rho }_0+\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })`$ and $`a=a\overline{a}_0[1\epsilon _\lambda (\epsilon _\rho \overline{\rho }_0+\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })y^2]`$. In this expression, $`\overline{\rho }_0`$ must be determined at zeroth order in $`\epsilon _\lambda `$, which yields $`\overline{\rho }_0=1`$ and $`\rho _m=\rho _m`$ solution to the usual 4D energy-momentum conservation equation. Thus, at first order in $`\epsilon _\lambda `$ one has $`a=a[1\epsilon _\lambda (\epsilon _\rho +\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })(y^21/12)]`$, $`\alpha =1(\epsilon _\lambda /6)(\epsilon _\rho +\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })`$ and $`\eta =\stackrel{~}{\eta }=1`$. Inserting these expressions in Eq. (26), one obtains (at first order in $`\epsilon _\lambda `$)

$$\frac{\mathrm{}_\lambda }{\mathrm{}_\mathrm{\Lambda }}=1\frac{\epsilon _\lambda (\epsilon _\rho 2)}{3},$$
(49)
and $`\alpha `$ reduces to

$$\alpha =1\frac{\epsilon _\lambda (\epsilon _\rho +1)}{6}.$$
(50)
####

From Eqs. (44), (49), and (50) it follows that the effective four-dimensional brane cosmology in the limit $`\epsilon _\lambda \mathrm{min}(1,1/\epsilon _\rho )`$ is given, at third order in $`\epsilon _\lambda `$, by

$`H_a^2={\displaystyle \frac{\kappa _{(5)}^2}{3}}\left[{\displaystyle \frac{(3+r_b/\mathrm{}_\lambda )+(34r_b/\mathrm{}_\lambda )\rho _m/(2\lambda )}{2r_b(1+\rho _m/\lambda )+3\mathrm{}_\lambda }}\right]\rho _m,`$ (52)
$`\dot{\rho _m}+3H_a(\rho _m+p_m)=0,`$ (53)
These equations have been obtained in the limit $`\epsilon _\lambda \mathrm{min}(1,1/\epsilon _\rho )`$, by keeping terms explicitely of third order in $`\epsilon _\lambda `$. The actual small parameter appearing in the expansion of the square root on the right-hand side of Eq. (29) is $`\epsilon _\lambda ^2\epsilon _\rho ^2\mathrm{max}(1,1/\epsilon _\rho )`$. It follows that these equations are correct, in the limit $`\epsilon _\lambda \mathrm{min}(1,1/\epsilon _\rho )`$, at third order in $`\epsilon _\lambda \epsilon _\rho `$ for $`\epsilon _\rho 1`$, and up to terms of order $`\epsilon _\lambda ^3\epsilon _\rho `$ for $`\epsilon _\rho 1`$. Note that in the domain $`\mathrm{min}(1,1/\epsilon _\rho )/\epsilon _\rho \epsilon _\lambda \mathrm{min}(1,1/\epsilon _\rho )`$, only the third order (in $`\epsilon _\lambda `$) correction corresponding to the term $`\epsilon _\lambda \epsilon _\rho `$ must be kept; the other corrections are smaller than (or of the same order as) terms already neglected. It can be checked that in the limit $`\epsilon _\lambda 1`$ and $`\epsilon _\rho 1`$, Eqs. (31) and (44) coincide at first order in $`\epsilon _\rho `$ and third order in $`\epsilon _\lambda `$, as it should be.

### C Uniform expressions in the regime: $`\epsilon _\lambda \epsilon _\rho 1`$

The two previous regimes cover the whole domain $`\epsilon _\lambda \epsilon _\rho \kappa _{(5)}^2r_b\rho _m/61`$ (with, of course, compatibility in the overlapping domain $`\epsilon _\lambda 1`$ and $`\epsilon _\rho 1`$ as mentioned at the end of the previous section. This suggests the following uniform expressions valid in the whole domain $`\epsilon _\lambda \epsilon _\rho 1`$:

$`H_a^2={\displaystyle \frac{\kappa _{(5)}^2}{3}}{\displaystyle \frac{[\alpha ^2(\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })+\eta \rho _m/(2\lambda )]}{r_b[\alpha ^2(\mathrm{}_\lambda /\mathrm{}_\mathrm{\Lambda })+\eta \rho _m/\lambda ]+\alpha ^2\mathrm{}_\lambda }}\eta \rho _m,`$ (55)
$`\dot{\rho _m}+3H_a(\rho _m+p_m)=0,`$ (56)
where

$`\alpha ^2={\displaystyle \frac{\mathrm{cos}(k/2)}{\mathrm{cos}^{1/2}(ky)^2}},`$ (58)
$`\eta ={\displaystyle \frac{f(y)\mathrm{cos}(ky)}{f(y)\mathrm{cos}^{1/2}(ky)^2}},`$ (59)
and

$`k=2\left(1+{\displaystyle \frac{\rho _m/\lambda }{1u}}\right)^{1/2}\mathrm{tan}^1\sqrt{{\displaystyle \frac{u}{1u}}},`$ (61)
$`{\displaystyle \frac{\mathrm{}_\lambda }{\mathrm{}_\mathrm{\Lambda }}}={\displaystyle \frac{\sqrt{u(1u)}}{\mathrm{tan}^1\sqrt{u/(1u)}}},`$ (62)
$`\epsilon _\lambda ={\displaystyle \frac{[\mathrm{tan}^1\sqrt{u/(1u)}]^2}{1u}},`$ (63)
with $`0u<1`$.

## IV Effective cosmological equations at higher orders in the case of a Minkowski bulk

In this section, we wish to extend the analysis of the previous section to higher orders. In order to simplify the calculations, we will restrict ourselves to the case where the bulk cosmological constant vanishes ($`\mathrm{\Lambda }=0`$), i.e. the spacetime outside the brane effectively is a five-dimensional Minkowski spacetime. In this case, one has $`\mathrm{}_\mathrm{\Lambda }^2=0`$ and Eq. (18) reads

$$H_a^2=\frac{2(\alpha ^2+\epsilon \eta )}{r_b^2}\left[1\sqrt{1\frac{\epsilon ^2\eta ^2+𝒞r_b^2/a^4}{(\alpha ^2+\epsilon \eta )^2}}\right].$$
(64)
As previously mentioned, in the limit of a vanishing brane thickness, Eq. (64) reduces to the unconventional Friedmann equation of thin brane cosmology

$$H_a^2=\frac{\kappa _{(5)}^4}{36}\rho _b^2+\frac{𝒞}{a^4},$$
(65)
(N.B. $`\alpha 1`$ and $`\eta 1`$ as $`r_b0`$). In the following, we will consider the limit $`\epsilon 1`$, $`r_b𝒞^{1/2}/a^2=O(\epsilon )`$, and determine the corrections to Eq. (65) due to a finite brane thickness (i.e. a finite $`\epsilon `$).

### A $`r_b^2H_a^2`$ at third order in $`\epsilon `$

At third order in $`\epsilon `$, Eq. (64) can be written as

$$r_b^2H_a^2=\frac{\epsilon ^2\eta ^2+𝒞r_b^2/a^4}{\alpha ^2+\epsilon \eta }.$$
(66)
At this order, $`\alpha `$ and $`\eta `$ must be expressed at first and zeroth order respectively \[note that $`\eta `$ behaves as $`1+O(\epsilon ^2)`$\]. To determine $`a`$ perturbatively in $`\epsilon `$ we will use the $`y`$-expansion (44) with $`\rho _b`$ replacing $`\rho _m`$. Inserting this expansion in Eq. (7) with $`\mathrm{\Lambda }=0`$, one obtains

$$(\mathrm{\Sigma }_a^2)^{\prime \prime }=4\epsilon \mathrm{\Sigma }_\rho \mathrm{\Sigma }_a^2+2r_b^2H_a^2.$$
(67)
Since one needs $`a`$ at first order in $`\epsilon `$, one can neglect the second term on the right-hand side of Eq. (67) \[which is $`O(\epsilon ^2)`$\], and consider the constant (in $`y`$) component of Eq. (67) only. One finds $`\overline{a}_1=\epsilon \overline{\rho }_0\overline{a}_0`$ and $`a=a\overline{a}_0(1\epsilon \overline{\rho }_0y^2)`$. In this expression, $`\overline{\rho }_0`$ must be determined at zeroth order in $`\epsilon `$, which yields $`\overline{\rho }_0=1`$ and $`\rho _b`$ solution to the usual 4D energy-momentum conservation equation. Thus, at first order in $`\epsilon `$ one has $`a=a[1\epsilon (y^21/12)]`$, $`\alpha =1\epsilon /6`$, and$`\eta =1`$. It follows that, at third order in $`\epsilon `$, the quantity $`r_b^2H_a^2`$ reads

$$r_b^2H_a^2=\frac{\epsilon ^2+𝒞r_b^2/a^4}{1+2\epsilon /3},$$

or, equivalently,

$$H_a^2=\frac{\kappa _{(5)}^4\rho _b^2/36+𝒞/a^4}{1+\kappa _{(5)}^2r_b\rho _b/9},$$
(68)
with

$$\dot{\rho _b}+3H_a(\rho _b+p_b)=0.$$
(69)
### B $`r_b^2H_a^2`$ at fourth order in $`\epsilon `$

At fourth order in $`\epsilon `$, Eq. (64) reads

$`r_b^2H_a^2`$ $`=`$ $`{\displaystyle \frac{\epsilon ^2\eta ^2+𝒞r_b^2/a^4}{\alpha ^2+\epsilon \eta }}\left[1+{\displaystyle \frac{\epsilon ^2\eta ^2+𝒞r_b^2/a^4}{4(\alpha ^2+\epsilon \eta )^2}}\right]`$ (70)
$`=`$ $`{\displaystyle \frac{\epsilon ^2\eta ^2+𝒞r_b^2/a^4}{\alpha ^2+\epsilon \eta }}\left[1+{\displaystyle \frac{1}{4}}r_b^2H_a^2+O(\epsilon ^3)\right],`$ (71)
which gives

$$r_b^2H_a^2=\frac{\epsilon ^2\eta ^2+𝒞r_b^2/a^4}{\alpha ^2+\epsilon \eta (\epsilon ^2\eta ^2+𝒞r_b^2/a^4)/4}.$$
(72)
In the numerator $`\eta `$ must be determined at second order in $`\epsilon `$. In the denominator $`\alpha `$ must be determined at second order in $`\epsilon `$ and $`\eta =1`$ \[N.B. $`\eta =1+O(\epsilon ^2)`$\]. Since one needs $`a`$ at second order in $`\epsilon `$, one must consider the components of Eq. (67) up to terms proportional to $`y^2`$. One finds $`\overline{a}_1=\epsilon \overline{\rho }_0\overline{a}_0+r_b^2H_a^2/(2\overline{a}_0)`$ and $`\overline{a}_2=(\epsilon /6)(2\overline{\rho }_0\overline{a}_1+\overline{\rho }_1\overline{a}_0)\overline{a}_1^2/(2\overline{a}_0)`$, which gives after some straightforward algebra

$`{\displaystyle \frac{a}{a}}`$ $`=`$ $`1\epsilon \left(y^2{\displaystyle \frac{1}{12}}\right)+{\displaystyle \frac{𝒞r_b^2}{2a^4}}\left(y^2{\displaystyle \frac{1}{12}}\right)`$
$`+`$ $`{\displaystyle \frac{\epsilon ^2}{12}}\left[5\left(y^2{\displaystyle \frac{1}{12}}\right)2\left(y^4{\displaystyle \frac{1}{80}}\right)\right]`$
$`+`$ $`{\displaystyle \frac{\epsilon \overline{\rho }_1}{12}}\left[\left(y^2{\displaystyle \frac{1}{12}}\right)2\left(y^4{\displaystyle \frac{1}{80}}\right)\right].`$
Thus, at second order in $`\epsilon `$ one has

$`\alpha =1{\displaystyle \frac{\epsilon }{6}}+{\displaystyle \frac{𝒞r_b^2}{12a^4}}+{\displaystyle \frac{\epsilon (11\epsilon +\overline{\rho }_1)}{180}},`$
$`\eta =1+{\displaystyle \frac{\epsilon (\epsilon 2\overline{\rho }_1)}{180}}.`$
To obtain $`\eta `$ one has used the expression $`\rho _b/\rho _b=1+\overline{\rho }_1(y^21/12)+O(\epsilon ^2)`$ and the fact that the $`O(\epsilon ^2)`$ term, the average of which must vanish, does not contribute to $`\rho _ba^2`$ at second order in $`\epsilon `$. Inserting these expressions of $`\alpha `$ and $`\eta `$ in (72), one obtains at fourth order in $`\epsilon `$

$$r_b^2H_a^2=\frac{\epsilon ^2[1+\epsilon (\epsilon 2\overline{\rho }_1)/90]+𝒞r_b^2/a^4}{1+2\epsilon /3\epsilon (9\epsilon \overline{\rho }_1)/90𝒞r_b^2/(12a^4)}.$$
(73)
The equation for $`\rho _b`$ reads

$$\dot{\rho _b}+3H_a(\rho _b+p_b)=3(\rho _b\mathrm{\Sigma }_\rho +p_b\mathrm{\Sigma }_p)\frac{\mathrm{log}\mathrm{\Sigma }_a}{\tau }.$$
(74)
Writing the right-hand side of this equation at second order in $`\epsilon `$, one obtains

$$\dot{\rho _b}=\frac{3H_a(\rho _b+p_b)}{1[\epsilon (\epsilon +\overline{\rho }_1)+p_b\epsilon (\epsilon +\overline{p}_1)/\rho _b]/60}.$$
(75)
To obtain this equation we have used $`\mathrm{\Sigma }_\rho =1+\overline{\rho }_1(y^21/12)+O(\epsilon ^2)`$, $`\mathrm{\Sigma }_p=1+\overline{p}_1(y^21/12)+O(\epsilon ^2)`$, and the fact that, at second order in $`\epsilon `$, there is no contribution of the $`O(\epsilon ^2)`$ terms. It remains to determine the evolution equation for $`\overline{\rho }_1`$ at lowest order in $`\epsilon `$, i.e. in which one keeps terms at most of first order. To this end let us assume an equation of state of the form $`p=w(\tau )\rho `$, (which implies $`\mathrm{\Sigma }_p=\mathrm{\Sigma }_\rho `$). In this case, the energy-momentum conservation reads

$$\rho _b\dot{\mathrm{\Sigma }}_\rho +[\dot{\rho _b}+3H_a(1+w)\rho _b]\mathrm{\Sigma }_\rho =3(1+w)\rho _b\mathrm{\Sigma }_\rho \frac{\mathrm{log}\mathrm{\Sigma }_a}{\tau },$$
(76)
which yields, at lowest order in $`\epsilon `$,

$$\dot{\overline{\rho }}_1=3(1+w)\dot{\overline{a}}_1=3(1+w)\dot{\epsilon }.$$
(77)
Note that the second term on the left-hand side of Eq. (76) is second order and must be neglected, as it can be seen from Eq. (75). If $`w`$ is constant over some time interval $`\tau _1\tau \tau _2`$, one can integrate Eq. (77), which gives

$$\overline{\rho }_1(\tau _1\tau \tau _2)=\overline{\rho }_1(\tau _1)+3(1+w)[\epsilon (\tau )\epsilon (\tau _1)].$$
(78)
To summarize, at fourth order in $`\epsilon `$, the quantity $`r_b^2H_a^2`$ reads \[Eq. (73)\]

$$r_b^2H_a^2=\frac{\epsilon ^2[1+\epsilon (\epsilon 2\overline{\rho }_1)/90]+𝒞r_b^2/a^4}{1+2\epsilon /3\epsilon (9\epsilon \overline{\rho }_1)/90𝒞r_b^2/(12a^4)},$$

with \[Eq. (75)\]

$$\dot{\epsilon }=\frac{3H_a(1+w)\epsilon }{1\epsilon [\epsilon (1+w)+\overline{\rho }_1(1+w^2)]/60},$$

and \[Eqs. (77) and (75)\]

$$\dot{\overline{\rho }}_1=9H_a(1+w)^2\epsilon .$$

This shows that at this order, and for the next orders as well, one must introduce auxiliary quantities, here $`\overline{\rho }_1`$, in order to describe adequately the effective cosmology from the brane point of view. This is a sign that one should not expect any effective four-dimensional reduction to be possible in the regimes of order one $`\epsilon `$.

## V Conclusion

In the present work, we have attempted to address in a phenomenological way the question of finite brane thickness and its influence on the effective four-dimensional cosmological equations. In order to simplify our analysis, we have made a number of assumptions, which it might be interesting to explore further or to relax in future investigations.

####

First we have defined in a heuristic way the notion of four-dimensional quantities by simply integrating over the fifth coordinate $`y`$, which is defined in a Gaussian normal coordinate system centered on the middle layer of the brane, supposed to be mirror symmetric. Although this averaging procedure corresponds to the standard one, one could maybe envisage, when having a warping factor, other types of averaging.

####

Another important assumption was to take the brane thickness constant in time. It would be worthwhile to generalize our results when the brane thickness is allowed to evolve in time.

####

After emphasizing a few limits of the present work, let us now summarize the main results. A very instructive result of our work is a generalized brane Friedmann equation, which interpolates between the familiar Kaluza-Klein picture, where the effective Friedmann equations after the averaging procedure yield exactly the usual four-dimensional Friedmann equations, and the initially surprising thin brane Friedmann equations, where the energy density enters quadratically. We hope that our generalized Friedmann equation will help in the intuitive understanding of the unconventional thin brane result, and clarify the delimitations of the various regimes where different Friedmann equations have to be applied.

####

Assuming the existence of a strict (in both patial and temporal senses) cosmological constant, we have also obtained a generalized version of the Randall-Sundrum condition of cancellation between the bulk cosmological constant and the brane tension term in the Friedmann equation. In our case, this cancellation condition depends explicitly on the brane thickness and yields back the familiar condition in the thin brane limit. Moreover, when this cancellation condition is assumed to hold and when the cosmological matter content of the brane is assumed to be small with respect to its tension, we were able to solve for the brane profile and obtain the effective cosmological equations.

####

Finally, we show that this is only in this case, as well as in the limit where the brane thickness is sufficiently small, that one can make sense of a four-dimensional effective cosmology. In the other cases, one can try to extend the regime of validity of the effective four-dimensional description but at the price of introducing auxiliary quantities, which is the sign that the full five-dimensional description is more adequate.

###### Acknowledgements.

We wish to thank P. Binétruy and N. Deruelle for very interesting discussions.

## A Five-dimensional Einstein equations

Inserting the ansatz (1) for the metric, the non-vanishing components of the Einstein tensor $`G_{AB}`$ corresponding to the antsatz (1) for the metric are:

$`G_{00}`$ $`=`$ $`3\left\{{\displaystyle \frac{\dot{a}}{a}}\left({\displaystyle \frac{\dot{a}}{a}}+{\displaystyle \frac{\dot{b}}{b}}\right){\displaystyle \frac{n^2}{r_b^2}}\left({\displaystyle \frac{a^{\prime \prime }}{a}}+{\displaystyle \frac{a^{}}{a}}\left({\displaystyle \frac{a^{}}{a}}{\displaystyle \frac{b^{}}{b}}\right)\right)\right\},`$ (A1)
$`G_{ij}`$ $`=`$ $`{\displaystyle \frac{a^2}{b^2}}\delta _{ij}\left\{{\displaystyle \frac{a^{}}{a}}\left({\displaystyle \frac{a^{}}{a}}+2{\displaystyle \frac{n^{}}{n}}\right)+2{\displaystyle \frac{a^{\prime \prime }}{a}}+{\displaystyle \frac{n^{\prime \prime }}{n}}\right\}`$ (A3)
$`+{\displaystyle \frac{a^2}{n^2}}\delta _{ij}\left\{{\displaystyle \frac{\dot{a}}{a}}\left({\displaystyle \frac{\dot{a}}{a}}+2{\displaystyle \frac{\dot{n}}{n}}\right)2{\displaystyle \frac{\ddot{a}}{a}}\right\},`$
$`G_{05}`$ $`=`$ $`3\left({\displaystyle \frac{n^{}}{n}}{\displaystyle \frac{\dot{a}}{a}}{\displaystyle \frac{\dot{a}^{}}{a}}\right),`$ (A4)
$`G_{55}`$ $`=`$ $`3\left\{{\displaystyle \frac{a^{}}{a}}\left({\displaystyle \frac{a^{}}{a}}+{\displaystyle \frac{n^{}}{n}}\right){\displaystyle \frac{r_b^2}{n^2}}\left({\displaystyle \frac{\dot{a}}{a}}\left({\displaystyle \frac{\dot{a}}{a}}{\displaystyle \frac{\dot{n}}{n}}\right)+{\displaystyle \frac{\ddot{a}}{a}}\right)\right\}.`$ (A5)
In the above expressions, a prime stands for a derivative with respect to $`y`$, and a dot for a derivative with respect to $`\tau `$.

Figure 1.

Figure 2.


format: text

text:
# Eccentricity Evolution for Planets in Gaseous Disks

## 1 Introduction

As of today, there are 77 extrasolar planets listed in the California & Carnegie Planet Search home page http://exoplanets.org/almanac.html. Almost all of the 57 planets with periods greater than 20 days move on orbits of substantial eccentricity, the mean and median being 0.35 and 0.37 respectively. By contrast, the 13 planets with periods less than 7 days have low eccentricity orbits, 0.03 mean and 0.02 median, apparently the result of eccentricity damping associated with tides raised in the planets by their central stars. The critical period separating small from large eccentricity orbits implies that extrasolar planets have tidal quality factors, $`Q10^5`$, similar to that of Jupiter.

Solar system planets, with the exception of Mercury, have low eccentricity orbits.111We consider Pluto to be a member of the Kuiper belt and not a bone fide planet. Thus the large orbital eccentricities of extrasolar planets came as a surprise. How might they have arisen? Interactions among planets have been the focus of most previous suggestions (Rasio & Ford, 1996; Lin & Ida, 1997; Weidenschilling & Marzari, 1996; Ford et al., 2001; Chiang et al., 2002). In systems with isolated planets, these involve speculative scenarios in which planets were either ejected from the system (Rasio & Ford, 1996) or merged with the remaining planet (Lin & Ida, 1997). Less attention has been paid to interactions between planets and the disks from which they formed. Analytic theory is not capable of resolving whether these interactions can produce eccentricity growth. However, it can identify crucial issues and thus provide guidance for targeted simulations. That is the goal of our paper.

Our investigation applies results derived in studies of satellite interactions with planetary rings by Goldreich & Tremaine (1979, 1980) and later extended by Ward (1986) and Artymowicz (1993b, a). Goldreich & Tremaine (1980, 1981) conclude that interactions at Lindblad resonances excite orbital eccentricity of both narrow rings and their shepherd satellites, whereas those at corotation resonances damp it. However, the balance between excitation and damping is a precarious one; damping exceeds driving by about 4.6% provided the corotation resonances are unsaturated. Ward (1986) and Artymowicz (1993a) point out that interactions at co-orbital Lindblad resonances damp its orbital eccentricity. The latter estimates that, in a gapless disk, damping by co-orbital Lindblad resonances is about a factor of three more rapid than excitation by external Lindblad resonances. In recent papers, Ward & Hahn (1998, 2000) claim that interactions at apsidal resonances, which damp eccentricity, are more potent than those at ordinary Lindblad resonances. On balance, current opinion favors the notion that planet-disk interactions damp eccentricity. We suggest that for massive planets the opposite may be true.

The structure of our paper is as follows. In §2, we argue following Rasio & Ford (1996) and Ford et al. (2001) that planet-planet interactions cannot account for the prevalence of eccentric orbits among isolated planets. §3 reviews the effects on eccentricity of planet-disk interactions at Lindblad, corotation, and apsidal resonances. We demonstrate in §4 that interactions at apsidal resonances are much less effective at damping eccentricity than previous estimates suggest. §5 is devoted to an evaluation of the balance between eccentricity damping at Lindblad resonances which overlap the planet’s orbit and eccentricity excitation at those which lie either inside or outside it. A planet which clears a sufficiently clean gap tilts the balance in favor of net eccentricity excitation by Lindblad resonances. However, linear theory predicts that corotation resonances enforce eccentricity damping. This leads to §6 in which we show that the torques at corotation resonances are weakened by saturation provided the planet’s orbital eccentricity is more than a few percent.

In assembling the case that orbital eccentricities of extrasolar planets may result from planet-disk interactions, we derive three new technical results.

* A correction to the standard Lindblad torque formula for apsidal resonances which applies when apsidal waves cannot propagate.
* An account of the balance between eccentricity excitation by Lindblad resonances and eccentricity damping by corotation resonances in a gap where the pressure gradient reduces the epicyclic frequency of the gas to a value that is much smaller than its orbital frequency.
* A treatment of the saturation of corotation resonances that includes the competition between the smoothing of the surface density gradient by the mass flux driven by the corotation torque and the action of viscous diffusion which acts to reestablish the large scale density gradient.

## 2 Planet-Planet Interactions

Is the observed eccentricity distribution of extrasolar planets the result of planet-planet interactions? In particular, can they account for the following properties:

* Typical eccentricities are of order a few tenths.
* Eccentricities smaller than one tenth are rare except for planet’s with orbital periods less than 20 days where tidal friction is likely to have damped eccentricity.
* Most of the planets discovered to date are not currently involved in a significant interaction with another planet of comparable or larger mass.

For planet-planet interactions to produce an eccentric orbit for an isolated planet requires at least one planet to disappear. The missing planet might have collided and merged with the remaining planet or effectively merged due to tidal capture. It might also have been ejected from the system or fallen into the central star. Numerical integrations of systems with two equal mass planets on initially circular orbits by Ford et al. (2001) produce a much greater fraction of isolated planets with low eccentricity orbits than is observed. These are a consequence of mergers which Ford et al. (2001) assume to occur whenever the separation between the planets drops below the sum of their radii. But the case against planet-planet interactions is even stronger than the results of Ford et al. (2001) indicate because they neglect tidal captures. Taking the tidal capture cross section for $`n=1`$ polytropes from Kim & Lee (1999), and the relative velocity at infinity as numerically calculated by Rasio & Ford (1996), we estimate a critical impact parameter for tidal capture between two and three times larger than the two radii assumed for merger by Rasio & Ford (1996) and Ford et al. (2001).

Next we consider some aspects of the merger process for a system whose initial state consists of two planets, each having mass $`M_p`$ and radius $`R_p`$, moving on coplaner circular orbits with radii $`r_1`$ and $`r_2`$ around a star of mass $`M_{}`$ and radius $`R_{}`$. We assume that the final state consists of a single or binary planet with mass $`2M_p`$ moving on an orbit with semimajor axis $`a`$ and eccentricity $`e`$.222Orbital elements for a binary refer to its center of mass. Applying conservation of energy and angular momentum, we relate the final orbit to the initial ones by

$$E=\frac{GM_{}M_p}{a}=\frac{GM_{}M_p}{2r_1}\frac{GM_{}M_p}{2r_2}\mathrm{\Delta }E,$$
(1)
and

$$H=2M_p\sqrt{GM_{}a(1e^2)}=M_p\left(\sqrt{GM_{}r_1}+\sqrt{GM_{}r_2}\right)\mathrm{\Delta }H.$$
(2)
Here $`\mathrm{\Delta }E`$ and $`\mathrm{\Delta }H`$ are the energy and angular momentum stored internally in either the merger product or in the relative orbit of the binary. Energy dissipated by impact or through tidal dissipation is accounted for by an increase of the binding energy.

We estimate $`\mathrm{\Delta }E`$ and $`\mathrm{\Delta }H`$ by noting that when the planets are separated by less than the Hill radius, $`r_{\mathrm{Hill}}(M_p/3M_{})^{1/3}a`$, they are effectively a two-body system. Thus333Ford et al. (2001) provide a similar derivation to ours. However they incorrectly assume that $`\mathrm{\Delta }E/E(M_p/M_{})^{1/2}`$ which leads to a higher eccentricity estimate.

$$\frac{\mathrm{\Delta }E}{E}\frac{M_pa}{M_{}r_{\mathrm{Hill}}}\left(\frac{M_p}{M_{}}\right)^{2/3}.$$
(3)
Moreover, since merging or tidal capture requires that the planet-planet periapse distance be no larger than a few times $`R_p`$,

$$\left|\frac{\mathrm{\Delta }H}{H}\right|\left(\frac{M_pR_p}{M_{}a}\right)^{1/2}\left(\frac{M_p}{M_{}}\right)^{2/3}\left(\frac{R_{}}{a}\right)^{1/2}.$$
(4)
The final expressions for $`\mathrm{\Delta }E/E`$ and $`\mathrm{\Delta }H/H`$ follow from the assumption that the planet and star have similar densities. In what follows we discard $`|\mathrm{\Delta }H/H|\mathrm{\Delta }E/E`$ since $`R_{}/a1`$.

Eccentricity is related to orbital energy and angular momentum by

$$e^2=1+\frac{2H^2E}{(2M_p)^3(GM_{})^2}.$$
(5)
Substituting for $`E`$ and $`H`$ using equations (1) and (2) yields

$$e^2=\frac{1}{4}\left(3\frac{r_1^2+r_2^2}{2r_1r_2}\sqrt{\frac{r_1}{r_2}}\sqrt{\frac{r_2}{r_1}}\right)\frac{(r_1+r_2+2\sqrt{r_1r_2})\mathrm{\Delta }E}{4GM_{}M_p}$$
(6)
Taking the leading contribution in $`|r_1r_2|/(r_1+r_2)`$ and writing

$$\frac{\mathrm{\Delta }E}{E}=A\left(\frac{M_p}{M_{}}\right)^{2/3},$$
(7)
where $`A1`$, we deduce that the eccentricity is bounded by

$$e^2=A\left(\frac{M_p}{M_{}}\right)^{2/3}\frac{3}{4}\left(\frac{r_1r_2}{r_1+r_2}\right)^2.$$
(8)
Maximal eccentricity is achieved if the initial orbits have the same semimajor axis. Mergers or captures are not possible for initial separations much larger than the Hill radius. We conclude that merging and tidal capture cannot produce eccentricities of more than a few percents for Jupiter like planets.

Our discussion of planet-planet interactions giving rise to orbital eccentricities is far from exhaustive. Here we briefly mention a scenario analyzed by Chiang et al. (2002). In it two planets undergo differential orbital migration which causes the ratio of their mean motions to diverge. Passage through mean motion resonances is shown to lead to eccentricity growth. A shortcoming of this work is that although planet-disk interactions are taken to be responsible for orbital migration, their direct effects on eccentricity evolution are not accounted for.

## 3 Planet-Disk Interaction - General

The gravitational potential of a planet can be expanded in a Fourier series in azimuthal angle $`\theta `$ and time $`t`$. Each term of this series is proportional to $`\mathrm{cos}[m(\theta \mathrm{\Omega }_{l,m}t)]`$ and has a radius dependent amplitude $`\varphi _{l,m}(r)`$. The mean motion $`\mathrm{\Omega }_p`$ is the unique pattern speed for a planet with a circular orbit. From here on, subscripts $`p`$ and $`d`$ will denote planet and disk, respectively. Pattern speeds for a planet which moves on an eccentric orbit may contain harmonics of the epicyclic frequency $`\kappa _p`$ and are denoted by $`\mathrm{\Omega }_{l,m}=\mathrm{\Omega }_p+(lm)\kappa _p/m`$. To first order in eccentricity $`e_p`$, each value of $`m`$ contributes three components, a principal one with pattern speed $`\mathrm{\Omega }_{m,m}=\mathrm{\Omega }_p`$ whose amplitude $`\varphi _{m,m}`$ is independent of $`e_p`$, and two first order components with pattern speeds444The first order $`m=0`$ terms are an exception. They should be combined into a single component proportional to $`\mathrm{cos}[\kappa _pt]`$, formally implying an infinite pattern speed. More on this is in the discussion section. $`\mathrm{\Omega }_{m\pm 1,m}=\mathrm{\Omega }_p\pm \kappa _p/m`$ whose amplitudes $`\varphi _{m\pm 1,m}`$ are proportional to $`e_p`$.555In this paper we consider potential components only up to first order in $`e_p`$.

Two kinds of resonance are associated with each potential component. Corotation resonances occur where the pattern speed matches the angular velocity of the disk material, $`\mathrm{\Omega }_{l,m}=\mathrm{\Omega }_d`$. A disk particle located at a corotation resonance experiences a constant torque which causes the radius of its orbit to change but does not excite its epicyclic motion. Lindblad resonances occur where the disk’s angular velocity differs from the pattern speed such that $`m(\mathrm{\Omega }_d\mathrm{\Omega }_{l,m})=\pm \kappa _d`$. The two Lindblad resonances associated with each potential component are distinguished by the adjectives inner and outer and are often denoted as $`ILR`$ and $`OLR`$. A disk particle located at a Lindblad resonance is subject to radial and azimuthal perturbation forces which vary at its epicyclic frequency. These excite its epicyclic motion and also change its semimajor axis.

Each $`m`$ has nine resonances associated with it: three potential components $`\varphi _{m,m}`$ and $`\varphi _{m\pm 1,m}`$, and three resonances for each potential component. Table 1 describes some properties of these resonances.

### 3.1 Understanding the sign of the torque

Each potential component is constant in a frame rotating with its pattern speed $`\mathrm{\Omega }_{l,m}`$, so its perturbations of the disk’s angular momentum $`H_d`$ and energy $`E_d`$ must preserve the Jacobi constant $`J_{l,m}=E_d\mathrm{\Omega }_{l,m}H_d`$. Thus

$$\frac{dE_d}{dt}=\mathrm{\Omega }_{l,m}\frac{dH_d}{dt}=\mathrm{\Omega }_{l,m}T_d,$$
(9)
where $`T_d`$ is the torque exerted on the disk. Provided the unperturbed disk is circular, $`e_d=0`$, equation (5) applied to the disk material on which the torque $`T_d`$ is applied yields

$$\left(\mathrm{\Omega }_{l,m}\mathrm{\Omega }_d\right)T_d=\frac{\mathrm{\Omega }_dH_d}{2}\frac{de_d^2}{dt}0.$$
(10)
Equality obtains if the epicyclic motion of the disk material is not excited. Thus a potential component may be viewed as transporting angular momentum from higher to lower angular velocity (Lynden-Bell & Kalnajs, 1972). Consequently, $`T_d>0`$ at an OLR and $`T_d<0`$ at an ILR.

Because $`\mathrm{\Omega }_{l,m}\mathrm{\Omega }_d=0`$ at a corotation resonance, the sign of the corotation torque depends upon whether the interaction is dominated by material inside or outside corotation. A simple understanding may be obtained by considering the behavior of collisionless particles. Disk particles close to the resonance experience slowly varying torques proportional to $`\mathrm{sin}[m(\theta \mathrm{\Omega }_{l,m}t)]`$ which change their orbital radii. Particles moving towards corotation drift at a decreasing rate with respect to the phase of the potential. The opposite pertains to particles moving away from corotation. As a result, there is a flux of particles towards corotation on both sides of the resonance. Material outside corotation loses angular momentum and that inside gains it. The corotation torque depends upon the difference and is proportional to $`d(\mathrm{\Sigma }/B)/dr`$. Oort’s constant $`B=(2r)^1d(r^2\mathrm{\Omega })/dr`$, which is equal to $`\mathrm{\Omega }/4`$ for Keplerian rotation, appears because the radial drift rate is equal to the torque divided by $`2B`$. Pressure effects, neglected here but accounted for later, spread the region of influence over a radial width which for a Keplerian disk is comparable to the disk’s vertical scale height. However, they do not affect the net corotation torque.

We note that the corotation torque is negative if material outside the resonance dominates the interaction. This might appear surprising in light of our finding, based on the Jacobi constant, that torques at Lindblad resonances transfer angular momentum from higher to lower angular velocity. Wouldn’t a similar conclusion apply to corotation resonances? The answer is no, but for a subtle reason. In evaluating the Jacobi constant near a Lindblad resonance we neglect the planet’s perturbation potential in the energy budget. That is fine there, but it is not appropriate near a corotation resonance.

### 3.2 Understanding the effects of resonances on the eccentricities.

Angular momentum and energy transferred to the disk at a resonance is removed from the planet’s orbit. Consequently,

$$\frac{dE_p}{dt}=\mathrm{\Omega }_{l,m}\frac{dH_p}{dt}=\mathrm{\Omega }_{l,m}T_d.$$
(11)
An alternate derivation proceeds from the fact that the linear density perturbation at a resonance rotates with the pattern speed. Thus the perturbed disk’s gravitational backreaction on the planet preserves the planet’s Jacobi constant. Combining equations (5) and (11), we find

$$e_p\frac{de_p}{dt}=\left(\left(1e_p^2\right)^{1/2}\mathrm{\Omega }_p\mathrm{\Omega }_{l,m}\right)\frac{H_p^2T_d}{(GM_{})^2M_p^3}.$$
(12)
At first order resonances $`\mathrm{\Omega }_{l,m}\mathrm{\Omega }_p`$, so to lowest order in $`e_p`$,

$$e_p\frac{de_p}{dt}=\left(\mathrm{\Omega }_p\mathrm{\Omega }_{l,m}\right)\frac{H_p^2T_d}{(GM_{})^2M_p^3}.$$
(13)
For first order Lindblad resonances, we make use of equation (10) to prove that

$$\mathrm{sgn}\left(\frac{de_p}{dt}\right)=\mathrm{sgn}\left[(\mathrm{\Omega }_p\mathrm{\Omega }_{l,m})(\mathrm{\Omega }_{l,m}\mathrm{\Omega }_d)\right].$$
(14)
Thus $`e_p`$ decreases if the resonance resides on the same side of corotation as the planet, and increases if it resides on the opposite side. Henceforth we adopt the terminology co-orbital666In a Keplerian disk these resonances overlap the planet’s semimajor axis. and external to distinguish these two classes of first order Lindblad resonances.

For first order corotation resonances the sign of the torque is opposite to the sign of $`d(\mathrm{\Sigma }/B)/dr`$. Therefore,

$$\mathrm{sgn}\left(\frac{de_p}{dt}\right)=\mathrm{sgn}\left[(\mathrm{\Omega }_{l,m}\mathrm{\Omega }_p)\frac{d}{dr}\left(\frac{\mathrm{\Sigma }}{B}\right)\right]$$
(15)
If the planet has cleared a gap, the density increases away from it so both first order corotation resonances cause $`e_p`$ to decay.

Next we examine how the principal Lindblad resonances affect the planet’s eccentricity. Since $`\mathrm{\Omega }_{m,m}=\mathrm{\Omega }_p`$, we retain the lowest order dependence on $`e_p`$ in equation (12) which reduces to

$$\frac{1}{e_p}\frac{de_p}{dt}=\frac{\mathrm{\Omega }_pH_p^2T_d}{2(GM_{})^2M_p^3}.$$
(16)
Goldreich & Tremaine (1980) compare the rates at which the principal and first order Lindblad resonances of the same $`m`$ change the planet’s eccentricity. They find that the former is smaller than the latter by a factor $`m`$, although both have the same dependence on $`e_p`$. We offer some additional comments. Because $`T_d`$ is negative at an ILR and positive at an OLR, the effect of the principal Lindblad resonances on $`de_p/dt`$ suffers from the same cancellation as their effect on $`da_p/dt`$. Rates of change of eccentricity and semimajor axis due to principal Lindblad resonances satisfy

$$\frac{4}{e_p}\frac{de_p}{dt}=\frac{1}{a_p}\frac{da_p}{dt}.$$
(17)
A planet which clears a gap migrates inward as the disk accretes, so its eccentricity increases. But since the timescale for eccentricity change is four times that for radial migration, principal Lindblad resonances are of negligible importance in eccentricity evolution.

Now we compare the timescale for eccentricity change, $`t_e`$, with that for semimajor axis migration, $`t_{vis}`$. The semimajor axis of a planet that opens a gap evolves with the disk on the viscous timescale

$$t_{vis}^1\frac{\nu }{r^2}\alpha \mathrm{\Omega }\left(\frac{h}{r}\right)^2.$$
(18)
Here $`\alpha `$ is the standard viscosity parameter for accretion disks, and $`hc_s/\mathrm{\Omega }`$ is the vertical thickness of the disk with $`c_s`$ the sound speed. The planet’s orbital eccentricity changes on the timescale

$$t_e^1\frac{1}{e_p}\left|\frac{de_p}{dt}\right|\left(\frac{r}{\mathrm{w}}\right)^4\left(\frac{M_p\mathrm{\Sigma }r^2}{M_{}^2}\right)\mathrm{\Omega },$$
(19)
where $`\mathrm{w}`$ is the gap’s width. Equation (19), adapted from Goldreich & Tremaine (1980), is obtained by summing contributions to $`de_p/dt`$ from first order resonances within a narrow ring separated from the planet’s position by an empty gap. Assuming Keplerian rotation, they find that $`e_p`$ decays. However damping by corotation resonances exceeds driving by Lindblad resonances by only a small, $`4.6\%`$, margin. The timescale quoted above for eccentricity change is that due to either corotation or Lindblad resonances acting separately. Conditions under which this, rather than the net contribution from both types of resonance, is the appropriate $`t_e`$ to compare with $`t_{vis}`$ are described in §6. To elucidate the comparison between $`t_e`$ and $`t_{vis}`$, we relate $`\mathrm{w}`$ to $`\alpha `$ and $`M_p/M_{}`$ by balancing the viscous torque777As written below, $`T_{vis}`$ is appropriate to a Keplerian disk. More generally, the factor of $`3`$ should be replaced by $`(2B)^1|rd\mathrm{\Omega }/dr|`$.

$$T_{vis}=3\pi \nu \mathrm{\Sigma }\mathrm{\Omega }r^2=3\pi \alpha \mathrm{\Sigma }r^2(\mathrm{\Omega }h)^2,$$
(20)
with the torque from the principal Lindblad resonances,

$$T_L\left(\frac{r}{\mathrm{w}}\right)^3\mathrm{\Sigma }r^2(r\mathrm{\Omega })^2\left(\frac{M_p}{M_{}}\right)^2,$$
(21)
to obtain

$$\frac{\mathrm{w}}{\mathrm{r}}(3\pi \alpha )^{1/3}\left(\frac{r}{h}\frac{M_p}{M_{}}\right)^{2/3}.$$
(22)
Our typical parameters, $`\alpha =10^3`$, $`h/r=0.04`$, $`M_p/M_{}=10^3`$, give $`\mathrm{w}0.4\mathrm{r}`$. Substituting this expression for $`\mathrm{w}`$ into equation (19) yields

$$t_e^1\frac{r}{\mathrm{w}}\frac{\mathrm{\Sigma }r^2}{M_p}t_{vis}^1.$$
(23)
For $`\mathrm{w}/\mathrm{r}1`$, eccentricity evolves faster than semimajor axis provided the disk is more massive than the planet.

Co-orbital Lindblad resonances can be ignored in investigations of the eccentricities of the orbits of shepherd satellites and narrow planetary rings because the surface density at the orbit of the shepherd satellite is negligible. But because the ratio of vertical thickness to orbital radius is so much larger in protoplanetary disks, $`h/r0.04`$, than in planetary rings, $`h/r10^6`$, even a massive planet may not clear a very clean gap so interactions at co-orbital Lindblad resonances could be significant. As shown by equation (14), these interactions damp the planet’s orbital eccentricity. Attempts to estimate the rate at which first order co-orbital Lindblad resonances damp eccentricity are described by Ward (1988) and Artymowicz (1993a). The latter concludes that, in a disk of uniform surface density, they damp eccentricity about three times faster than external first order Lindblad resonances excite it.

An apsidal resonance is a first order inner Lindblad resonance with $`l=0`$ and $`m=1`$. Its pattern speed is equal to the planet’s apsidal precession rate; $`\mathrm{\Omega }_{0,1}=\mathrm{\Omega }_p\kappa _p=\dot{\varpi }_p`$. As a co-orbital first order Lindblad resonance, it contributes to eccentricity damping. Application of the standard Lindblad resonance torque formula to the apsidal resonance led Ward & Hahn (1998, 2000) to conclude that an apsidal resonance damps eccentricity much faster than external first order Lindblad resonances excite it.

We have described three types of resonant planet-disk interactions that damp eccentricity and one that excites it. In what follows we argue that eccentricity can grow provided its initial value exceeds some minimal threshold and the planet clears a sufficiently clean gap. A minimal initial eccentricity is necessary for the nonlinear saturation of corotation torques. This involves a reduction of the surface density gradient at the resonance position. A clean gap renders ineffective damping by first order co-orbital Lindblad resonances. Although damping due to the apsidal resonance might appear to be the most potent of all, it is likely to be unimportant. Apsidal waves have such long wavelengths in protoplanetary disks that they do not propagate. Thus torques at apsidal resonances are much smaller than predicted by the standard Lindblad resonance torque formula.

## 4 Apsidal Resonances

Apsidal resonances are special because the doppler shifted forcing frequency detunes more gradually from the epicyclic frequency with distance away from an apsidal resonance than it does at other Lindblad resonances. This is why the standard Lindblad resonance torque formula predicts an apsidal torque

$$T_{0,1}^{ILR}e_p^2\left(\frac{M_p}{M_{}}\right)^2\frac{\mathrm{\Sigma }r^8\mathrm{\Omega }^2}{(\mathrm{\Delta }r)^4}\frac{\mathrm{\Omega }}{r|d\dot{\varpi }_d/dr|},$$
(24)
which is larger by $`\mathrm{\Omega }/\dot{\varpi }_d10^3`$ than that at an ordinary first order Lindblad resonance at a distance $`\mathrm{\Delta }rr`$ from the planet’s orbit (Goldreich & Tremaine, 1978; Ward & Hahn, 1998, 2000). For the same reason, apsidal waves have longer wavelengths than those at ordinary Lindblad resonances.

The long wavelengths of apsidal waves call into question the applicability of the standard torque formula. Of particular concern is the assumption that density waves propagate away from the resonance and ultimately dissipate without reflection. To estimate the first wavelength, $`\lambda _1`$, of a density wave at a Lindblad resonance, we start from the WKBJ dispersion relation, $`k^2c_s^2=m^2(\mathrm{\Omega }_d\mathrm{\Omega }_{l,m})^2\kappa _d^2D`$, which applies in a disk with negligible self-gravity.888Because of the long wavelengths of apsidal waves in protoplanetary disks, a more complete treatment should include the disk’s self-gravity. We define $`\lambda _1`$ by

$$_r_{}^{r_{}+\lambda _1}k(r)𝑑r=2\pi ,$$
(25)
where $`r_{}`$ is the resonance radius. Expanding $`D𝒟(rr_{})/r_{}`$ around the resonance position, we obtain

$$\lambda _1=(3\pi )^{2/3}(c_s^2/|𝒟|r_{}^2)^{1/3}.$$
(26)
For apsidal waves $`\mathrm{\Omega }_{0,1}=\dot{\varpi }_p`$, which implies $`D2\mathrm{\Omega }_d(\dot{\varpi }_p\dot{\varpi }_d)`$. Crudely approximating $`d\dot{\varpi }_d/dr\dot{\varpi }_d/r`$, we obtain $`𝒟2\mathrm{\Omega }_d\dot{\varpi }_d`$. Hence

$$\lambda _1/r(9\pi ^2/2)^{1/3}(\mathrm{\Omega }_d/\dot{\varpi }_d)^{1/3}(h/r)^{2/3}.$$
(27)
Substituting reasonable parameters for a protoplanetary disk, $`\mathrm{\Omega }_d/\dot{\varpi }_d10^3`$ and $`h/r0.04`$, we estimate that $`\lambda _14r`$, somewhat longer than the disk radius. How does the wavelength far away from the resonance compare to the local radius? In the absence of self gravity apsidal waves propagate toward smaller $`r`$. We expect that $`\dot{\varpi }_d`$ increases inward. Applying the dispersion relation to estimate $`\lambda =2\pi /k`$ where $`\dot{\varpi }_d\dot{\varpi }_p`$ yields

$$\lambda /r=\sqrt{2}\pi (\mathrm{\Omega }_d/\dot{\varpi }_d)^{1/2}(h/r)5,$$
(28)
provided we adopt the same values for $`\mathrm{\Omega }_d/\dot{\varpi }_d`$ and $`h/r`$ used in estimating $`\lambda _1`$. Our estimates for $`\lambda /r`$ imply that there are not many wavelengths, probably not even one, in an apsidal wave train.

Apsidal waves are not traveling waves subject to the nonlinear steepening and shock dissipation that is the leading form of damping for density waves generated at ordinary Lindblad resonances (Goodman & Rafikov, 2001). Instead, they are standing waves trapped in a cavity bounded by the resonance radius and the inner edge of the disk, with viscous dissipation as their principal source of damping. Consequently, the standard Lindblad resonance torque formula cannot be applied to determine the torque at an apsidal resonance in a protoplanetary disk. Fortunately, simple physical arguments enable us to estimate the factor by which the apsidal torque is reduced below the estimate provided by the standard formula. Suppose that the apsidal cavity contains $`N`$ wavelength. Comparing the viscous dissipation timescale, $`(\nu k^2)^1`$, with that for propagation across a wavelength, $`2\pi /(kv_g)=2\pi \mathrm{\Omega }_d/(kc_s)^2`$, we deduce that the apsidal torque is a factor $`2\pi N\alpha `$ smaller than the standard torque formula implies. Thus for plausible parameters, $`N1`$ and $`\alpha 10^3`$, it is even smaller than the torque at other first order Lindblad resonances.

We have deliberately glossed over several items that merit attention. These are briefly covered below.

* The location of the apsidal resonance is uncertain. Factors influencing the apsidal precession of the disk include its gravity, the planet’s gravity, and gas pressure. Only the disk’s gravity affects the precession of the planet’s orbit. Thus we might expect that $`\mathrm{\Delta }rr`$.
* Equation (24) implicitly assumes that $`\lambda _1r`$. For $`\mathrm{\Delta }r<\lambda _1<r`$, the prediction of the standard torque formula should be multiplied by the factor $`F(\mathrm{\Delta }r/\lambda _1)^2`$. Appropriate correction factors if either $`\mathrm{\Delta }r>r`$ or $`\lambda _1>r`$, or both apply, remain to be worked out.
* Then there is the issue of the cavity size measured in wavelengths. Our estimate of torque reduction assumes an amplitude for the standing cavity wave similar to that which a traveling wave would have in the absence of a reflecting boundary. This is an adequate approximation provided $`N1`$, but it would be an underestimate for a cavity tuned near resonance.
* As mentioned previously, a more careful treatment would include the disk’s self-gravity. However, we can safely state that if self-gravity dominates gas pressure in the determination of $`\lambda _1`$, then $`\lambda _1`$ would be larger than calculated by neglecting self-gravity which would further reduce the torque.

## 5 Lindblad and Corotation resonances

First order Lindblad resonances are of two kinds, external ones which excite the planet’s eccentricity, and co-orbital ones which damp it. In a disk of uniform surface density, damping by co-orbital resonances is about a factor of three faster than excitation by the external ones (Artymowicz, 1993a). Obviously, damping due to co-orbital resonances is weaker if the planet orbits within a gap. However, surface density gradients associated with a gap activate eccentricity damping by corotation resonances. Can a gap be clean enough to prevent eccentricity damping by co-orbital Lindblad resonances yet have sufficiently shallow density gradients so as to avoid eccentricity damping by corotation resonances? Although we have not been able to answer this question in general, the two examples described below suggest that this is not possible.999In this section we neglect possible saturation of corotation resonances.

Here we compare the rate of eccentricity damping by co-orbital Lindblad resonances plus corotation resonances with the rate of eccentricity excitation by external Lindblad resonances. Our analysis takes account of a gap in the surface density around the planet’s orbit but assumes that this does not significantly perturb the Keplerian rotation rate. For simplicity, we concentrate on the material external to the planet’s orbit. In this region the external and co-orbital Lindblad resonances are outer and inner ones respectively. We employ three constants, $`𝒞_{eLR}`$, $`𝒞_{cLR}`$, and $`𝒞_{CR}`$, to describe the relative strengths of the external Lindblad resonances, the co-orbital Lindblad resonances, and the corotation resonances. As described in §3, $`𝒞_{cLR}/𝒞_{eLR}3`$ and $`𝒞_{CR}/𝒞_{eLR}=1.046`$. We take the torque cutoff to occur at $`r_1a_ph`$ for each type of resonance (Goldreich & Tremaine, 1980; Ward, 1986; Artymowicz, 1993b).

Eccentricity excitation by first order external Lindblad resonances is evaluated from

$$\frac{1}{e_p}\frac{de_p}{dt}|_{eLR}=_{r_1}^{\mathrm{}}\frac{𝒞_{eLR}}{(ra_p)^5}\left(\frac{\mathrm{\Sigma }(r)}{B_d}\right)𝑑r.$$
(29)
A similar expression describes eccentricity damping by first order co-orbital Lindblad resonances, expect that $`\mathrm{\Sigma }`$ is evaluated at the planet’s semimajor axis,

$$\frac{1}{e_p}\frac{de_p}{dt}|_{cLR}=_{r_1}^{\mathrm{}}\frac{𝒞_{cLR}}{(ra_p)^5}\left(\frac{\mathrm{\Sigma }(a_p)}{B_d}\right)𝑑r=\frac{𝒞_{cLR}}{4(r_1a_p)^4}\frac{\mathrm{\Sigma }(a_p)}{B_d}.$$
(30)
The effect of corotation resonances on eccentricity evolution is given by

$$\frac{1}{e_p}\frac{de_p}{dt}|_{CR}=_{r_1}^{\mathrm{}}\frac{𝒞_{CR}}{4(ra_p)^4}\frac{d}{dr}\left(\frac{\mathrm{\Sigma }(r)}{B_d}\right)𝑑r.$$
(31)
Integrating by parts, we arrive at

$$\frac{1}{e_p}\frac{de_p}{dt}|_{CR}=_{r_1}^{\mathrm{}}\frac{𝒞_{CR}}{(ra_p)^5}\left(\frac{\mathrm{\Sigma }(r)}{B_d}\right)𝑑r+\frac{𝒞_{CR}}{4(r_1a)^4}\frac{\mathrm{\Sigma }(r_1)}{B_d}.$$
(32)
Except for the negative sign and the substitution of $`𝒞_{CR}`$ for $`𝒞_{eLR}`$, the first term is identical to that for external Lindblad resonances. The second term, which is positive, may be viewed as a correction due to disk material close to the planet’s orbit.

Equations (29), (30), and (32) combine to yield

$$\frac{1}{e_p}\frac{de_p}{dt}=_{r_1}^{\mathrm{}}\frac{\left(𝒞_{CR}𝒞_{eLR}\right)}{(ra_p)^5}\left(\frac{\mathrm{\Sigma }(r)}{B_d}\right)𝑑r\frac{\left(𝒞_{cLR}\mathrm{\Sigma }(a_p)𝒞_{CR}\mathrm{\Sigma }(r_1)\right)}{4(r_1a_p)^4B_d}.$$
(33)
If the planet clears a sufficiently clean gap, only the first term survives, and eccentricity damping by corotation resonances overcomes eccentricity driving by external Lindblad resonances, but only by a small margin since $`𝒞_{CR}/𝒞_{eLR}=1.046`$.101010This situation pertains to narrow planetary rings and their shepherd satellites (Goldreich & Tremaine, 1980). Residual material in the gap leads to additional eccentricity damping unless $`\mathrm{\Sigma }(a_p)/\mathrm{\Sigma }(r_1)<𝒞_{CR}/𝒞_{cLR}1/3`$. However, such a steep density gradient near the center of the gap is incompatible with our assumption that the disk maintains Keplerian rotation.

Next we examine how departures from Keplerian rotation might affect the balance between eccentricity driving and damping. Gas pressure perturbs the disk’s epicyclic frequency. It decreases $`\kappa _d`$ at the outskirts of a gap and increases it around the middle.111111For an isothermal gas the leading correction to $`\kappa _d^2`$ is proportional to the second radial derivative of the logarithm of the surface density. Provided the gap is sufficiently clean, only resonances in its outer parts need be considered. Numerical simulations of gap formation by planets show that this limit can be achieved, at least in two-dimensional disks. This requires $`1|d^2\mathrm{ln}\mathrm{\Sigma }/d\mathrm{ln}r^2|(r/h)^2`$. The upper limit is set by the requirements that principal Lindblad resonances are located away from the planet’s semimajor axis and by the Rayleigh stability criterion. It motivates us to investigate the eccentricity evolution of a planet that resides in a gap where the epicyclic frequency is much smaller than the orbital angular velocity.

As $`\kappa _d/\mathrm{\Omega }_d0`$, $`d\mathrm{\Omega }_d/dr2\mathrm{\Omega }_d/r`$, and the positions of the inner and outer Lindblad resonances associated with a given $`\varphi _{m1,m}`$ potential component approach, from opposite sides, the location of the corresponding corotation resonance. This simplifies comparison of the influences of the different resonances on eccentricity evolution. Inserting these approximations into the expressions for the torques at first order resonances given by Goldreich & Tremaine (1979), we obtain

$$T_{m1,m}^{CR}\frac{m\pi ^2\varphi _{m1,m}^2r}{\kappa _d^2}\frac{d\mathrm{\Sigma }}{dr},$$
(34)
and

$$T_{m1,m}^{OLR,ILR}\pm \frac{m^2\pi ^2\varphi _{m1,m}^2\mathrm{\Sigma }}{\kappa _d^3/\mathrm{\Omega }}.$$
(35)
Note that in the limit $`\kappa _d/\mathrm{\Omega }_d0`$, each Lindblad torque dominates the corotation torque. However, the Lindblad torques have opposite signs and when summed making use of their separations $`r\kappa _d/2m\mathrm{\Omega }_d`$, we discover that

$$T_{m1,m}^{ILR}+T_{m1,m}^{OLR}=\frac{m\pi ^2\varphi _{m1,m}^2r}{\kappa _d^2}\frac{d\mathrm{\Sigma }}{dr}.$$
(36)
So to leading order, the torques at the inner and outer Lindblad resonances just cancel that at the corotation resonance implying $`de_p/dt=0`$.

To summarize, neither of our examples yields a prediction of eccentricity growth. For Keplerian rotation eccentricity damping by corotation resonances is slightly faster than eccentricity excitation by external Lindblad resonances. In the extreme limit of vanishing epicyclic frequency eccentricity damping by corotation and co-orbital Lindblad resonances just balances eccentricity driving by external Lindblad resonances.

## 6 Saturation of Corotation Resonances

Linear perturbations of the disk are calculated near discrete resonances. Density wave fluxes and disk torques obtained from the linear perturbations are of second order in the planet’s potential. By transporting and depositing angular momentum, they cause a secular evolution of the disk. As the disk evolves, the linear perturbations, the fluxes, and the torques also change. In the particular case of a corotation resonance, the disk torque drives a mass flux which flattens the gradient of $`\mathrm{\Sigma }/B`$ thereby shutting down the torque. This is the phenomenon of resonance saturation.

### 6.1 collisionless particle disk

It is easy to understand the saturation of the corotation torque in a system of collisionless particles (see discussion in §3.1). Near a corotation resonance each particle experiences a slowly varying torque which causes a slow radial drift of its orbit. There is a net flux of particles towards corotation on each side of the resonance. Trapping of particles around peaks of the potential occurs within a region of width $`\mathrm{\Delta }r(\varphi /\mathrm{\Omega }^2)^{1/2}`$ on timescale $`\mathrm{\Delta }t(r^2/m\varphi )^{1/2}`$. Trapped particles circulate around potential peaks. Their angular momenta vary periodically but not secularly. After trapping occurs the corotation torque vanishes. A quantitative discussion is given by Goldreich & Tremaine (1981).

There is an alternate way to view corotation saturation in a collisionless particle disk that is more in keeping with what happens in a gas disk. Trapping kills the gradient $`d(\mathrm{\Sigma }/B)/dr`$ within a radial interval of width $`\mathrm{\Delta }r`$ around the resonance radius. Thus the action of the potential changes the disk so that the torque shuts off.

### 6.2 gas disk

The following discussion draws on material in §IV of Goldreich & Tremaine (1979). It also presents new results omitting details of their derivations. These will be the subject of a future paper.

The corotation torque density, that is, the torque per unit radius near corotation ($`r=r_{CR}`$), is given by

$$\frac{T_d}{r}=2\pi r\overline{\mathrm{\Sigma }_1\frac{\varphi _{l,m}}{\theta }}=\frac{T_{l,m}^{CR}}{2\mathrm{}}\mathrm{exp}(|rr_{CR}|/\mathrm{}),$$
(37)
where the overbar denotes an azimuthal average and the total corotation torque

$$T_{l,m}^{CR}=\frac{m\pi ^2\varphi _{l,m}^2}{2|d\mathrm{\Omega }/dr|}\frac{}{r}\left(\frac{\mathrm{\Sigma }}{B}\right)|_{CR}.$$
(38)
Note that as a result of gas pressure the torque density is spread over a radial distance $`\mathrm{}c_s/\kappa _d`$. In a Keplerian disk, $`\kappa _d=\mathrm{\Omega }_d`$ so $`\mathrm{}=h`$. A technical comment is in order here. Although the torque density is spread over an interval of width $`\mathrm{}`$ around corotation, the gradient $`d(\mathrm{\Sigma }/B)/dr`$ which appears in equation (38) is to be evaluated at the resonance radius. That is the meaning of the subscript $`CR`$ which we attach to the gradient in this section. In a viscous disk the value of the corotation torque is determined by the average value of this gradient within a boundary layer of width

$$\delta _\nu \left(\frac{\nu }{m|d\mathrm{\Omega }/dr|}\right)^{1/3},$$
(39)
where $`\nu `$ is the effective kinematic viscosity in the disk. Using the standard prescription $`\nu =\alpha h^2\mathrm{\Omega }`$,

$$\frac{\delta _\nu }{\mathrm{}}\left(\frac{\alpha h^2r}{m\mathrm{}^3}\right)^{1/3}.$$
(40)
Thus $`\delta _\nu \mathrm{}`$ provided $`\alpha m\mathrm{}^3/(h^2r)`$, which we assume to be true in what follows.

The corotation torque density drives an azimuthally averaged radial mass flux

$$F_\mathrm{\Sigma }\overline{\mathrm{\Sigma }v_r}=\frac{1}{2Br}\left[\frac{1}{2\pi r}\frac{T_d}{r}+\frac{1}{r}\frac{}{r}\left(r^3\nu \mathrm{\Sigma }\frac{d\mathrm{\Omega }}{dr}\right)\right].$$
(41)
Equation (41) has a simple physical interpretation. Each term in the square brackets is a torque per unit area, the first comes from the planet and the second from the viscous stress. Orbits of disk particles drift adiabatically in response to absorbing angular momentum at rate inversely proportional to $`2Brd(r^2\mathrm{\Omega })/dr`$, the radial gradient of specific angular momentum. However, a rigorous derivation is more complicated since the Reynold’s stress also transports angular momentum and angular velocity perturbations affect the angular momentum density. These two effects cancel each other leaving equation (41) intact.

Next we substitute the expression for the mass flux given by equation (41) into the equation of mass conservation and use equations (37) and (38) to eliminate the corotation torque density. Under the assumption that $`\mathrm{\Sigma }`$ varies on a shorter spatial scale than $`\nu `$, $`\mathrm{\Omega }`$, $`B`$, and $`d\mathrm{\Omega }/dr`$, we arrive at

$$\frac{\mathrm{\Sigma }}{t}=\frac{\nu r|d\mathrm{\Omega }/dr|}{2B}\frac{^2\mathrm{\Sigma }}{r^2}\mathrm{sgn}\left(rr_{CR}\right)\frac{\pi m\varphi _{l,m}^2}{16r^2\mathrm{}^2|d\mathrm{\Omega }/dr|B^2}\frac{\mathrm{\Sigma }}{r}|_{CR}\mathrm{exp}(|rr_{CR}|/\mathrm{}).$$
(42)
Specializing to near Keplerian rotation, we arrive at

$$\frac{\mathrm{\Sigma }}{t}=3\nu \frac{^2\mathrm{\Sigma }}{r^2}\mathrm{sgn}\left(rr_{CR}\right)\frac{2\pi m\varphi _{l,m}^2}{3rh^2\mathrm{\Omega }^3}\frac{\mathrm{\Sigma }}{r}|_{CR}\mathrm{exp}(|rr_{CR}|/h).$$
(43)
The large scale density gradient in a gap is maintained by a balance between torques from principal Lindblad resonances and the viscous stress. Understanding the action of the two terms on the right hand sides of equations (42) and (43) is straightforward. The mass flux driven by the differential corotation torque tends to decrease the density gradient within a region of width $`\mathrm{}`$ around corotation. Viscous diffusion acts to maintain the density gradient at the value it has on larger scales. For parameters of interest to us, the density gradient suffers only a small fractional reduction. By integrating the steady state version of equation (43) once with respect to $`r`$, we find

$$\frac{\mathrm{\Sigma }}{r}=\frac{\mathrm{\Sigma }}{r}|_{\mathrm{}}\frac{2\pi m\varphi _{l,m}^2}{9\alpha rh^3\mathrm{\Omega }^4}\frac{\mathrm{\Sigma }}{r}|_{CR}\mathrm{exp}(|rr_{CR}|/h).$$
(44)
A more revealing form for first order corotation resonances is obtained using $`\varphi _{m\pm 1,m}me(M_p/M_{})(r\mathrm{\Omega })^2`$, $`mr/\mathrm{w}`$, and substituting equation (22) into equation (44);

$$\frac{\mathrm{\Sigma }}{r}\frac{\mathrm{\Sigma }}{r}|_{\mathrm{}}e_p^2\frac{r}{h}\frac{\mathrm{\Sigma }}{r}|_{CR}\mathrm{exp}(|rr_{CR}|/h),$$
(45)
which implies

$$\frac{\mathrm{\Sigma }}{r}|_{CR}\frac{1}{1+e_p^2r/h}\frac{\mathrm{\Sigma }}{r}|_{\mathrm{}}$$
(46)
A more careful derivation establishes that equation (45) is not missing a significant numerical constant. With $`h/r0.04`$, $`e_p0.2`$ is required for complete saturation, but $`e_p0.05`$ suffices to tip the balance in favor of eccentricity driving by Lindblad resonances as opposed to eccentricity damping by corotation resonances.

Partial saturation of corotation resonances leading to eccentricity growth further requires that the initial value of $`e_p`$ be able to survive damping until the torque is adequately saturated. To assess the condition under which this is satisfied, we compare the timescale for eccentricity damping given by equation (19),

$$t_e\left(\frac{\mathrm{w}}{\mathrm{r}}\right)^4\frac{M_{}^2}{M_p\mathrm{\Sigma }r^2\mathrm{\Omega }},$$
(47)
with that during which $`d\mathrm{\Sigma }/dr|_{CR}`$ achieves a steady state,

$$t_{sat}\frac{h^2}{\nu }\frac{1}{\alpha \mathrm{\Omega }}.$$
(48)
Applying equation (22) yields

$$\frac{t_{sat}}{t_e}\frac{h^2}{r\mathrm{w}}\frac{\mathrm{\Sigma }r^2}{M_p}.$$
(49)
For plausible parameters, $`\mathrm{w}/\mathrm{r}1`$ and $`h/r0.04`$, this ratio is less than unity except for very large values of $`\mathrm{\Sigma }r^2/M_p`$. Under the dual assumptions of a steady state gap and steady state saturation, eccentricity growth is governed by

$$\frac{1}{e_p}\frac{de_p}{dt}=\frac{1}{t_e}\left(1\frac{𝒞_{CR}/𝒞_{eLR}}{1+re_p^2/h}\right)$$
(50)
Sample solutions of this equation are plotted in figure 2.

## 7 Discussion

Resonant interactions between a planet and a protostellar disk cause the planet’s orbital eccentricity to evolve on a short timescale. We have investigated whether these interactions might be responsible for the eccentric orbits that characterize recently discovered extrasolar planets. Our preliminary findings are as follows.

Co-orbital Lindblad resonances damp eccentricity. They dominate the eccentricity evolution of planets that are too small to clear gaps (Ward, 1988; Artymowicz, 1993a). External Linblad resonances excite eccentricity (Goldreich & Tremaine, 1980). Lindblad resonances, by themselves, lead to net eccentricity excitation for planets which open clean gaps. However, gradients of surface density are enhanced in the walls of gaps and these activate eccentricity damping by corotation resonances. Unsaturated, corotation resonances damp eccentricity slightly faster than external Lindblad resonances excite it (Goldreich & Tremaine, 1980). Corotation resonances can saturate because the planet’s torque drives a mass flux which tends to smooth the surface density gradient near corotation. However, saturation is only effective if the smoothing occurs on the timescale of viscous diffusion across a radial scale comparable to the disk thickness. This requires the planet to obtain an initial orbital eccentricity of a few percent. Thus eccentricity growth due to planet-disk interactions is a finite amplitude instability.

An apsidal resonance is a special kind of inner Lindblad resonance whose pattern speed is equal to the apsidal precession rate of the planet. Apsidal waves are longer than density waves excited at ordinary Lindblad resonances. The standard Lindblad resonance torque formula predicts that the apsidal torque is enhanced relative to that at other Lindblad resonances. Formally, this implies that eccentricity damping due to apsidal waves is more potent than other forms of eccentricity change due to planet-disk interactions. However, because of their long wavelengths, apsidal waves are standing rather than traveling disturbances in protoplanetary disks. Unlike traveling waves, they do not steepen into shocks and dissipate, but instead are weakly damped by viscosity. As a result, the apsidal torque is much smaller than the standard torque formula predicts. It is plausible that apsidal waves do not play a significant role in the eccentricity evolution of extrasolar planets.

Next we provide concrete example of a system in which eccentricity growth might occur. Suppose that a Jupiter mass planet with $`M_p/M_{}10^3`$, initial semimajor axis $`a_p`$AU, and initial eccentricity $`e_p0.06`$, is embedded in a protostellar disk with radius to scale high ratio $`r/h25`$ and viscosity parameter $`\alpha 10^3`$. The planet would open a wide gap with $`\mathrm{w}\mathrm{r}/2`$. If the gap were as clean as the Rayleigh stability condition, $`\kappa ^2>0`$, permits, the residual surface density at its center would be so small that eccentricity damping due to co-orbital Lindblad resonances would be negligible. Acting separately, either eccentricity damping by corotation resonances or eccentricity excitation by external Lindblad resonances would proceed on a timescale of order $`10^5`$ years. However, if the corotation resonances were unsaturated and they were acting together with the Lindblad resonances, eccentricity would damp on a significantly longer timescale. Given the initial eccentricity, partial saturation of the corotation resonances leading to eccentricity growth would occur on a timescale of order $`2\times 10^3`$ years, the viscous diffusion timescale across a radial interval comparable to the disk’s vertical scaleheight. As the eccentricity increased, so would the saturation of the corotation resonances, and thus the eccentricity growth rate. Figure 2 shows that the eccentricity could reach the observed mean value of $`e=0.35`$ within about $`10^6`$ years.

Investigations of eccentricity evolution due to planet-disk interactions are plagued by several major uncertainties.

* Angular momentum transport by internal stresses plays a major role in disk accretion, gap formation, and the saturation of corotation resonances. Presumably it arises from an instability, but we lack a basic understanding of its origin and character. At present, all we can do is assume that the disk possesses a kinematic $`\alpha `$ viscosity. Dissipation of disturbances near resonances occurs through shocks and turbulence. How these contribute to angular momentum transport is a subject for future investigation. Analogous phenomena are involved in maintaining sharp edges in planetary rings (Borderies et al., 1982a, b).
* Clean gaps are essential for the suppression of eccentricity damping by co-orbital Lindblad resonances. Two dimensional simulations produce gaps of the requisite cleanliness, but it is not yet known whether three dimensional simulations will confirm their reality (Bryden et al., 2000). The Rayleigh stability criterion, $`\kappa ^2>0`$, limits the magnitude of density gradients in gap walls. There are also weaker, nonaxisymmetric instabilities that come into play at positive $`\kappa ^2`$. Examining their role in determining gap shape is a task best suited for a coordinated attack by analytical calculations and numerical simulations.
* Jovian mass planets should be able to open gaps in protostellar disks. Those that were not accreted onto and consumed by their central stars must have been present while the disk dispersed. During the latter stages of disk dispersal, they would have been more massive than the disk. It is known from studies of narrow planetary rings that the growth and damping of a ring’s eccentricity is in most respects identical to that of the satellites which shepherd it (Goldreich & Tremaine, 1980). The ratio of the rates of eccentricity change of ring and satellite are inversely proportional to the ratio of ring to satellite mass. To maintain an eccentric shape a narrow ring must precess rigidly. This requires the wavelength of the apsidal wave with pattern speed equal to the precessional angular velocity to be longer than the ring’s radial width. Since apsidal waves in protostellar disks have wavelengths comparable to the local radius, large regions of these disks might maintain eccentric shapes.
* Our proposal for eccentricity growth due to planet-disk interactions involves a finite amplitude instability. However, we do not have an obvious candidate for giving a planet the requisite initial eccentricity of a few percent. Gravitational interactions with surface density perturbations created by instabilities in the walls of gaps is one possibility.
* Once an adequate core has been assembled, the accretion of gas to form a Jovian mass planet proceeds rapidly. As its mass increases, a planet begins to clear a gap which hinders its ability to accrete additional gas. Eccentricity growth could occur prior to gap attaining its equilibrium width provided it become sufficiently clean. This would have two positive effects on the planet’s potential for eccentricity growth. It would lessen the initial eccentricity required for adequate saturation of corotation resonances since $`e_{crit}\mathrm{w}^{3/2}`$ for fixed viscosity, and it would speed up the rate of eccentricity increase since $`de_p/dt\mathrm{w}^4`$.
* Saturation of corotation resonances involves the flattening of $`\mathrm{\Sigma }/B`$ in regions of width $`\mathrm{}`$ around each resonance. Neighboring resonances are separated by distances of order $`r/m^2`$, so these regions overlap for $`m>(r/\mathrm{})^{1/2}`$. Where overlap occurs, the interplay between gap maintenance and corotation saturation will be more complicated than we have described. Presumably $`\mathrm{\Sigma }/B`$ would flatten in a region narrower than $`\mathrm{}`$ and this would require a larger eccentricity to overcome viscous diffusion. Although this is an unresolved difficulty, it is unlikely to be of great significance since with our standard parameters the most important resonances lie well outside the region of overlap.
* The sites of first order corotation resonances coincide with those of principal Lindblad resonances. Linear perturbations associated with the latter are larger than those of the former by a factor of order $`(me_p)^1`$. Whether this affects saturation of corotation resonances is an open issue.
* We only consider principal and first order resonances in this paper. While this is adequate for examining the excitation and damping of small eccentricities, higher order resonances will have to be taken into account if we want to examine the growth of eccentricity to the large values characteristic of extrasolar planets.
* Hill radii of Jupiter mass planets are larger than the vertical scale heights of protostellar disks. Thus they might at least temporarily trap gas in quasi two dimensional horseshoe orbits. Then some of the angular momentum deposited at co-orbital Lindblad resonances could find its way back into the planet’s orbit. This would reduce the eccentricity damping rate of co-orbital Lindblad resonances.
* Our analysis is suitable for $`m1`$. However, for Jupiter mass planets, and for $`\alpha =10^3`$, the size of the gap $`\mathrm{w}\mathrm{r}`$ suggests $`m1`$. Significant corrections may apply in this limit.
* The first order axisymmetric, $`m=0`$, resonance deserves special attention. It can be thought of as a co-orbital Lindblad resonance, and as the arguments of §3 demonstrate, it leads to eccentricity damping. However, preliminary investigation shows that it is of little importance.

Our mechanism of eccentricity growth through saturation of corotation resonances can ultimately be tested by numerical simulations. Recently, Papaloizou et al. (2001) observed eccentricity growth for planets larger than about $`1020M_J`$ in an $`\alpha =4.5\times 10^3`$ disk. With these parameters the gap extends past the $`2:1`$ resonance making all first order corotation resonances impotent. Similar behaviour was observed earlier by Artymowicz et al. (1991) for binary stars surrounded by a disk. Both papers appear to support our claim that the apsidal resonance does not damp eccentricity as fast as external Lindblad resonaces excite it. In apparent contradiction with our conclusions, Papaloizou et al. (2001) did not obtain eccentricity growth for Jovian mass planets. Several explanations come to mind. Limitations of our analysis may have led us to predict eccentricity growth where it does not occur. Alternatively, the simulations may lack sufficient resolution to capture the partial saturation of corotation resonances. Indeed the grid spacing of their Jovian mass simulation was of order $`\mathrm{}`$. Finally, their initial $`e_p`$ may be below the limit required for adequate saturation.

This research was supported in part by an NSF grant awarded to PG and a Sherman Fairchild Fellowship held by RS. We thank Geoff Bryden for several illuminating conversations.