Planetesimal accretion around binaries

Relative velocities among accreting planetesimals in binary systems: the circumbinary case

H. Scholl, F. Marzari, P. Thébault
Laboratoire Cassiopée, CNRS, Observatoire de la Côte d’Azur, B.P. 4229, F-06304 Nice, France
Department of Physics, University of Padova, Via Marzolo 8, 35131 Padova, Italy
Stockholm Observatory, Albanova Universitetcentrum, SE-10691 Stockholm, Sweden, and
Observatoire de Paris, Section de Meudon, F-92195 Meudon Principal Cedex, France
Accepted …. Received …..; in original form …

We numerically investigate the possibility of planetesimal accretion in circumbinary disks, under the coupled influence of both stars’ secular perturbations and friction due to the gaseous component of the protoplanetary disk. We focus on one crucial parameter: the distribution of encounter velocities between planetesimals in the 0.5 to 100km size range. An extended range of binary systems with differing orbital parameters is explored. The resulting encounter velocities are compared to the threshold velocities below which the net outcome of a collision is accumulation into a larger body instead of mass erosion. For each binary configuration, we derive the critical radial distance from the binary barycenter beyond which planetesimal accretion is possible. This critical radial distance is smallest for equal-mass binaries on almost circular orbits. It shifts to larger values for increasing eccentricities and decreasing mass ratio. The importance of the planetesimals’ orbital alignments of planetesimals due to gas drag effects is discussed.

extrasolar planets – accretion – planetesimals – binary stars

1 Introduction

At present more than 40 planets have been found in binary star systems, which corresponds to about 20 % of the whole sample of known extrasolar planets (Desidera & Barbieri, 2007). All of them are on so–called S–type or circumprimary orbits that encircle one component of the double star system. Planets may also revolve on dynamically stable orbits about both binary stars on so–called P–type or circumbinary orbits (e.g. Holman and Wiegert 1999). At present, there is only one planet assumed to revolve on a circumbinary orbit. The primary, HD 202206, is a metal rich main sequence star. Its companion, HD 202206b, is assumed to be a low–mass brown–dwarf with an estimated mass of 17.5 according to radial velocity measurements (Udry et al., 2002; Correia et al., 2005). This mass is beyond the widely accepted brown-dwarf limit of 10 Jupiter masses. The third body, HD 202206c, is assumed to be a Jupiter–like planet revolving about the binary. Adopting a mass of 2.41 for HD 202206c, Correia et al. (2005) show in a dynamical analysis that HD 202206c and HD 202206b might be in a 5/1 mean motion resonance. This resonant configuration is possibly a consequence of an inward migration of HD 202206c due to the forces exerted by a viscous circumbinary disk (Nelson 2003).

The signature of circumbinary planets in radial velocity measurements is in most cases difficult to be spotted because of the short–term large–amplitude velocity fluctuations on the primary induced by the companion star. However, planets might be very frequent among the population of close star couples. Mid–infrared emissions have in fact revealed substantial circumbinary material around PMS close binaries like DQ Tau, UZ Tau, GW Ori (Mathieu et al. 2000) and AK Sco (Jensen & Mathieu, 1997), and more detached binaries like GG Tau (Dutrey et al., 1994), UY Aur (Close et al., 1998). These disks are significantly more massive than the minimum–mass solar nebula which suggests abundant planetary formation. However, both concurrent mechanisms for extrasolar planet formation, namely the standard solid–core model forming planetary embryos from runaway and oligarchic accretion of planetesimals (e.g Greenberg et al., 1978; Wetherill & Stewart, 1989; Kokubo & Ida, 2000) possibly followed by gas infall onto the core for the formation of giant planets (Pollack et al. 1996; Bodenheimer and Lin 2002), and the alternative model of local gravitational collapse of disk material (Boss, 1997) might be affected by the secular perturbations of the binary.

We will focus here specifically on the “standard” solid–core scenario. in particular on the stage of planetary embryo accretion. This stage is just before the one investigated by Quintana and Lissauer (2006) who determined regions around close binaries where planets can form by accumulation of embryos. Within the frame of this model, the initial stage of planetesimal accretion is the runaway growth of isolated planetary embryos. This stage is very fast and efficient, provided that the encounter velocities among planetesimal are low, i.e. much smaller than the escape velocities of the growing embryos (Lissauer, 1993). As a consequence, runaway growth is particularly sensitive to any external perturbation able to stir up relative velocities in the planetesimal swarm. In the specific case of a circumbinary disk, relative velocities between planetesimals can be expected to be higher as compared to the single star case. In the vicinity of the binary, velocities could be increased to values for which runaway accretion is no longer possible, or even to higher values for which collisions may result in mass removal via cratering or fragmentation rather than accretion. The threshold velocity for the stop of runaway is of the order of the escape velocities of the growing bodies while the threshold velocity for the limit between accreting and eroding impacts is less straightforward to estimate and depends on the respective sizes and physical properties of the two colliding planetesimals. In any case, the distribution of encounter velocities is the key parameter controlling the fate of a collisionally interacting swarm of planetesimals. This is the issue which we address here, by numerically estimating the distribution of in a disk of planetesimals surrounding a binary system in order to derive assertive conclusions on the possibility of planet formation.

Only taking into account the gravitational potential of both stars might lead to incorrect estimates of relative velocities, in particular close to the star couple. Indeed, the planetesimal accretion phase is believed to take place while a large fraction of the primordial gas disk is still present. Gas drag is known to circularize orbits (Adachi et al., 1976) and thus reduce relative velocities in dynamically hot systems (Marzari & Scholl, 2000). Besides circularization, gas drag has the additional effect of aligning orbital apsides, an effect which is also called orbital synchronization (Marzari et al., 1997; Marzari & Scholl, 2000). The consequences of this synchronization on relative velocities depend on planetesimal sizes. Bodies of the same size all have their periapses aligned in the same direction. Smaller bodies have a much smaller dispersion of periapse directions than larger bodies. Furthermore, the mean direction of the periastron alignment is also different for objects of different sizes (Thébault et al., 2006). As a result, while relative velocities between small planetesimals of the same size may drop to almost zero, relative velocities among bigger bodies are more significant because of a larger periapse direction dispersion. The highest relative velocities can be expected between objects of different sizes due to the size–dependency of the periastron alignment. The relative velocity distribution depends, therefore, not only on the radial distance to the barycenter of the binary but also on the sizes of colliding planetesimals and on the gas density.

In a pioneering attempt to model planetary accretion in circumbinary disks, Moriwaki and Nakagawa (2004) used a purely gravitational model without gas drag. They found very stringent limits on the minimum distance from the binary beyond which planetesimal accumulation is possible. However, as pointed out above, neglecting the velocity damping effect and orbit alignment due to gas friction may yield misleading results for the efficiency of planetesimal accumulation, in particular in the vicinity of the binary. In the present work, we redress this issue by performing numerical simulations where the effects of gas drag are taken into account, although with unavoidable simplifying assumptions. Moreover, we compute the average encounter velocities between planetesimals from their trajectories and not from the values of forced eccentricities as Moriwaki & Nakagawa (2004) who used the relation , where denotes Keplerian velocity. This often used relation holds only under the condition that all orbits have fully randomized Keplerian angles. However, in circumbinary disks with gas, orbital phasing of periastra directions is of high importance because of the combined effects of secular perturbations and gas friction. Another limitation of the relation is that it works and it does not take into account the possibility of large radial excursions of bodies with high eccentricities (for more on this subject, see for instance the discussion in Thébault & Doressoundiram, 2003). Both these issues suggest that the relation used by Moriwaki and Nakagawa (2004) may not be appropriate in presence of large orbital eccentricities and possibly of gas drag. Our numerical approach, which takes into account all close particle encounters representing collisions, has the advantage to handle automatically orbital phasing (and dephasing) as well as radial mixing effects. This allows us to derive more reliable values for the limiting distance from the barycenter of the binary where planet formation is possible for different parameters of the binary. Not surprisingly, our results significantly differ from those of Moriwaki and Nakagawa (2004) in particular close to the binary where planetary accretion is more likely to form planets due to the higher density of material.

Due to obvious numerical constraints, mutual gravitational interactions among planetesimals, like for instance dynamical friction, are not included even if, in principle, they might be taken into account in N-Body simulations. However, the speed of current computers is still far too low to allow a systematic investigation of binary systems surrounded by a self-gravitating planetesimal disk. These simulations are still restricted to disks with massless planetesimals. Their size distribution, of course, does not evolve with time. To do so, particle–in–a–box computations, which can follow the evolution of a large number of massive bodies, are required, but the price to pay is a significant loss of accuracy in the computations of the orbital evolution of the planetesimals that, in the present scenario, is crucial. As a consequence, we decided to adopt a deterministic numerical approach, which, although neglecting important physical effects, gives a reliable estimate of encounter velocity distributions and thus enable to derive “maps” of the disk regions with relative velocities low enough to favor planetesimal accretion and planet formation.

Our numerical model is described in detail in Section 2. Section 3 is devoted to the computations performed in the gas-free case, a scenario similar to that of Moriwaki & Nakagawa (2004). In Section 4 we use the results obtained from models including gas drag in order to map the regions where we expect planet formation by accretion. Section 5 is dedicated to the dicussion of our results.

2 The numerical model

2.1 Main characterstics

Our deterministic code is an upgraded version from the one used in several earlier studies (Marzari & Scholl, 1998, 2000), adapted to the specificities of the circumbinary configuration. We consider a system of planetesimals, treated as non–gravitationally interacting test particles, submitted to the gravitational potentials of both stars and to a gas drag force following Adachi et al. (1976) and described in the following section. Planetesimal orbits are integrated using a fourth-order Runge-Kutta scheme with fixed step size. We like to point out that for the circumprimary problem (Thébault et al., 2006) we applied a different numerical scheme which is less convenient for the circumbinary case.

The binary stars move on an elliptic orbit of eccentricity and semi–major axis . For sake of simplicity and to limit the number of free parameters to a manageable value, we assume for all runs AU. This is the standard value also taken by Moriwaki & Nakagawa (2004) in their numerical explorations. We are thus here implicitly focusing on planetesimal accretion in the regions typically beyond 3–4 AU from the binary’s centre of mass. This region is sufficiently away from the critical instability region found by Holman and Wiegert (1999) to guarantee orbital stability for almost planar and circular motion. Both stars’ total mass is fixed and equal to 1 solar mass. The primary parameters of our model are thus the mass ratio and eccentricity of the binary star system. We use five different mass ratios for the binary, = 0.5, 0.4, 0.3, 0.2 and 0.1. The stellar masses of the primary and of the companion are then and , in solar mass units. For each mass ratio, models with six different binary eccentricities are simulated with = 0.0, 0.1, 0.2, 0.3, 0.4 and 0.5. This gives a total of 30 models. Each model is simulated with and without gas drag. In models taking into account gas drag, we consider planetesimal of physical radii ranging from 0.5km to 100km. The initial eccentricities and orbital inclinations are chosen with a uniform random distribution between 0 and . The initial particle semi–major axis are randomly distributed between AU and AU. The value for being approximately equal to the inner truncation of a circumbinary disk by the tidal influence of the binary (Artymowicz et al., 1994; Gunther & Kley, 2002). Whenever the semimajor axis of a planetesimal decreases below the inner limit because of drag friction, the body is removed.

The typical timescale for one run is 100 000 binary revolutions, which is approximately the time necessary to reach a steady state, with fully developed secular perturbations and the gas drag induced orbital phasing. After such a steady state is reached, we determine relative velocities among the planetesimal population. For this purpose, we use the positions and velocities of all planetesimals at a fixed time. For each planetesimal, we search its nearest neighbouring planetesimal and compute the relative velocity corrected for Kepler shear due the different radial distances. Relative velocities are sampled in bins in radial direction. A bin size of 0.5 AU is used. The median value for each bin is then taken as the representative value for the collisional or impact velocity of planetesimals at the corresponding radial distance. Since, as pointed out, the orbit alignment is size dependent, we compute the median relative velocities between planetesimal populations for each pairs of radii , in each bin. We typically use 20 to 30 close encounters in a bin to determine the corresponding median relative velocity.

2.2 The gas drag force

The drag force is modelled in laminar gas approximation by


where is the drag parameter given by (Kary et al., 1993) as


where is the planetesimal radius, its mass density, the gas density of the protoplanetary disk and a dimensionless drag coefficient related to the shape of the body ( for spherical bodies). Exploring the gas drag density profile as a free parameter would be too CPU time consuming and we shall restrict ourselves to one distribution. We assume the standard Minimum Mass Solar Nebula (MMSN) of Hayashi (1981), with and We take a typical value

We are here implicitly assuming that the gas disk is axisymmetric and pressure supported. We are aware that this axisymmetric prescription does probably not represent the real physical structure of the circumbinary gaseous component of the disk. This is a complex issue, which has only been addressed in a handful of elaborated numerical investigations (e.g Artymowicz et al., 1994; Gunther & Kley, 2002; Nelson, 2003; Pichardo et al., 2005). All these studies clearly show the inner tidal truncation of the disk, but are less detailed on the structure within the disk, which is particularly complex to be precisely computed, especially in the eccentric binary case. The only clear result is the onset of complex time–evolving spiral patterns extending from the inner edge of the disk towards the outer regions. At present it appears a very demanding task to perform a fully coupled numerical treatment of both the gas disk and an embedded planetesimal swarm. For this reason, we prefer in a first step a simplified approach where the gas drag force is just described by Equ.1. The onset of density waves in the disk induces radial and azimuthal variations both of the disk density and of the gas velocity where an additional radial component adds up to the Keplerian tangential component. These effects might alter the strength of the drag acceleration when a planetesimal enters and exits the higher or lower density region. However, it is reasonable to assume that, on the average, this effect would not strongly modify accelerations obtained by Equ.1 and thus the dynamical evolution of planetesimals, at least for planetesimals in the km–size range. In particular, we expect the occurrence of the periastron alignment.

As already mentioned, we start all planetesimals on almost circular orbits. This choice is from a dynamical point a view equivalent to assuming that, at , the proper eccentricity is almost equal to the forced one . This is a somewhat artificial and not a necessarily realistic starting configuration motivated by the fact that, at present, we have no generally accepted clues either on how planetesimals form within a disk with spiral waves nor on the initial orbital elements of the planetesimals as they detach from the gas. We integrate the model until a steady state is reached where periapses are aligned and the interaction with the gas becomes only a minor perturbation of the Keplerian orbit. This procedure is valid for small planetesimals, since they reach the same periapse alignment independently of their starting values with the same distribution of proper eccentrcities. Larger planetesimals, which are much less affected by gas drag, have their proper eccentricities increasing during their growth by collisions and close encounters with other large planetesimals. As a consequence, they too should lose memory of their initial orbits. As a consequence, we think that our initial choice of eccentricities does not introduce too significant effects on the calculation of relative impact velocities once the steady state is reached.

3 Impact velocities in the gas free case

In the gas free models we use 25000 particles with semimajor axes ranging from 3 to 50 AU. However, we focus here essentially on the region within 12AU, since, as it clearly appears in Fig.1, average encounter velocities drop down to almost zero beyond this limit. For the region, the we obtain might be fitted, for the range of stellar masses are orbital eccentricities explored, by the empirical formula


which is a second-order polynomial in and multiplied by an additional term equal to accounting for the radial dependence of on the local Keplerian velocity.

Since we neglect the gas drag force, the radius of the planetesimals is not a parameter of the problem. In Fig.1 we also compare our results to the computed with the Moriwaki & Nakagawa (2004) analytical approximation for the specific case with and . For reference, we also plot the numerical data obtained from the N-body integration. It can be clearly seen that our numerically derived values are always significantly smaller than the Moriwaki & Nakagawa (2004) estimates. This shows the weakness of the approximate relation for the present case, where orbital phasing and radial excursions of the planetesimals play a crucial role in determining the average impact velocity. As an example, while our simulations give low relative velocities, of the order of 50 m/s at about 11 AU from the barycenter, the analytical approximation of Moriwaki & Nakagawa (2004) predicts such low velocities only beyond 50 AU. Our numerical results give much lower values for impact velocities which allow planetary formation to occur much closer to the barycenter of the binary even in absence of gas friction. This is a first indication that planets on circumbinary orbits may be abundant due to the low relative velocities in the regions near the binary with higher surface density. In the following section, we show how this promising result is affected by gas drag.

Figure 1: Relative impact velocities between planetesimals computed with our equation 3 (continuous line) compared to the predictions of (Moriwaki & Nakagawa, 2004) (dashed line). The empty circles show the numerical results of the N-body integrations. The parameters of the binary are and .

4 Impact velocities in the presence of gas

We now turn to the more realistic case including gas drag. As previously mentioned, it is well known that the combination of the drag force and perturbations by the gravity field of the binary causes a strong paeriapse phasing of the orbits of the planetesimals. The planetesimals’ physical sizes are now a crucial parameter. Bodies of similar sizes have their periapses aligned towards approximately the same direction, which leads to very low impact velocities. Different size planetesimals are aligned towards different directions, increasing the non–tangential component of the impact velocities.

To compute the relative velocity between any possible pair of different sized planetesimals in the swarm is far beyond the present computing capabilities. Therefore, we derive the mean relative velocities for two separate cases assumed to be representative for the accreting capability of a swarm:

  • between equal-size bodies of radius . In this case the periapse alignment is stronger.

  • for a target with radius and a projectile with radius . This size ratio has been chosen since it is the configuration for which, in the large majority of cases explored in our circumprimary study (Thébault et al., 2006), the kinetic energy delivered by the impactor is maximal. For the relative velocity decreases and approaches the low velocity value of equal–size impacting bodies. For , the increase in is too slow to compensate for the diminishing impactor mass.

Figure 2: Figure showing the relative velocity between planetesimals for different values of the binary eccentricity , mass ratio and planetesimal sizes and . The plots are organized as follows: top left (equal–mass stars case), km and km; top right , km and km; bottom left , km and km; bottom right , km and km.

For each of the 30 different (,) binary configurations considered, we perform 6 simulations exploring different target–impactor pairs (,) in the kmkm range (leading to a total of 180 runs). We compute for each pair the average impact velocity . In Fig.2 we show some examples of the relative velocity between different size planetesimal swarms. A general feature is that binaries with higher eccentricities always lead to higher . Conversely, the most perturbing configurations are the ones with lowest values, the equal–mass case () leading to the lowest . This result is logical when considering the structure of the secular term of the perturbing potential (see for example the Appendix of Moriwaki & Nakagawa, 2004): the forced eccentricity term is proportional to , i.e. it vanishes for the equal–mass case (where only short period perturbation terms remain). An alternative and simplified, but partially correct explanation is that the equal mass case is the one for which the radial excursion of both stars from the centre of mass is the more limited, i.e., they both stay at 0.5 AU from it, whereas for smaller values, star number 2 moves further away from the centre of mass and “enters” deeper into the circumbinary disk. Within the range of parameters which have been explored, the most perturbed case is thus logically the one with mass ratio and . In any case, even for different size impacting objects, one robust result is that the presence of gas drag always strongly reduces the relative velocity between planetesimals. This can be clearly seen in Table I, where we display, for one given set of binary parameters (the same as in Fig.1), relative velocity values for the case without gas drag, both from our simulations and from Moriwaki and Nakagawa (2004), and for the case with gas drag for a target planetesimal of size km. We consider both the case of equal-size planetesimals with km and different size planetesimals with km and km.

As discussed in the introduction, all values have to be compared to two critical values in order to determine the collisional evolution of a planetesimal swarm. One is the escape velocity of the target , which can easily be derived from the body’s radius and density. Only if can runaway growth occur. The other crucial threshold value is , the threshold velocity beyond which the net outcome of an impact is mass erosion instead of accumulation into a larger body. To estimate the value of for a given projectile and target is a complex task. It has to be derived from models of cratering and fragmentation which, at present, give uncertain results due to our poor knowledge of the size–strength scaling law. As discussed in Thébault et al. (2006), there are many different models in the literature giving the value of the critical specific shattering energy and the size distribution of fragments as a function of the projectile and target size, of composition and of impact velocity. To be conservative, we consider as in Thébault et al. (2006) three different estimates of computed according to Holsapple (1994), Marzari et al. (1995) and Benz & Asphaug (1999) prescriptions for and Thébault et al. (2003) prescription for the fragmented and cratered masses. Out of these we take the maximum value that gives an optimistic estimate of the region where accretion is possible, and the minimum and more conservative value.

r (AU) No gas No gas Gas drag Gas drag
(our simulations) (Moriwaki & Nakagawa, 2004) R1=R2=5km R1=5km R2=2.5km
4.25 849 m.s 2040 m.s 18 m.s 12 m.s
7.75 283 m.s 754 m.s 3 m.s 9 m.s
11.75 49 m.s 444 m.s 5 m.s 20 m.s
Table 1: Examples of impact velocity values obtained with different models, for the same example binary configuration as in Fig.1 (, ).

When we compare to and we can outline 3 different accretional modes:

  • 1) . In this case, runaway accretion can proceed as in the ’standard’ unperturbed scenario. The additional perturbations of the binary are not large enough to cancel out the gravitational focusing factor of the target, which is the source of this fast growth mode e.g. (Lissauer, 1993).

  • 2) For values of larger than but still smaller than , planetary accretion is still possible but runaway growth will either not occur at all or it will start only when large bodies have formed such that . In this case planetesimals will accumulate into planetary bodies but at a significantly slower rate.

  • 3) If exceeds , the majority of collisions end up with cratering and fragmentation that overcome accretion. The planetesimal population, instead of growing in size, will be slowly ground down to dust.

For sake of clarity we display our numerical predictions on planetesimal accretion in 4 summarizing graphs. In each of them we display, for all explored values of and , the radial distance beyond which planetesimal accretion is possible, for a population of objects having an initial size , all the way up to the maximal size of 100km considered in our simulations. Beyond that size, accretion dominates over erosion for all explored cases. In other words, is the limit beyond which for all () pairs in the range. We consider two cases for the starting size: one “standard” case with km, and one “small planetesimals” case for which km. Note that the value km is the commonly accepted minimum size for planetesimals: Independently of their formation process, when they reach this size they detach from the gas of the disk and move on independent Keplerian orbits, only perturbed by the gas. This value is however approximate and this is why we also explore a lower value of 0.5km for the initial planetesimal size. We also consider the two maximum and minimum values derived from the different prescriptions considered for .

Fig.3 displays the most favorable case for accretion, where we assume km and adopt the highest, and thus less restrictive value for . It can be seen that, in spite of the binary gravitational perturbations and of the periastra de–phasing of different size planetesimals, for mass ratios 0.5 and 0.4 planet formation proceeds in all the regions of the disk and for all binary eccentricities. The only region where accretion is inhibited is within the inner tidal gap of the disk. For smaller values of and high binary eccentricities , the inner border for accretion shifts to larger radial distances. For and planet formation can occur only beyond 10 AU: strong secular perturbations due to the large eccentricity and low mass ratio of the binary prevent planetesimal accumulation closer to the barycenter. Note however that for most cases where we found accretion to be possible the growth mode can probably not be runaway, or at least not runaway. This is because values are not for all explored values. The only exceptions, where unperturbed runaway accretion is possible, are obtained for the equal–mass case and low binary eccentricities (the white rectangles in Fig.3).

Figure 3: Map of the radial distance beyond which planetary formation is possible, i.e. where all () encounters lead to mass accretion () for all objects of size km, as a function of the binary mass ratio and binary eccentricity . The color coding used in Figures 3-6 is the following: -green: AU (the inner edge of our planetesimal disc) -pale blue: 4AUAU -dark blue: 6AUAU -red: 9AUAU -black: AU (the outer edge of our planetesimal disc) The radial distance given at the center of some rectangles is the minimum value beyond which runaway accretion is possible (i.e. ) for all km.

For a more conservative choice of the erosion velocity (Fig.4), the inner limit for planetesimal accretion moves significantly outward, as expected: only for , i.e. equal mass stars, is accretion possible at any radial distance, independently of .

The cases with smaller 0.5km initial planetesimals appear to be more critical (5 and Fig.6). While for mass ratios of 0.5 and 0.4 planet formation is (as in the km case) possible almost all over the disk, when becomes smaller than 0.4 in most cases accretion of small planetesimals is inhibited within 12 AU from the stars. The main reason for this accretion inhibiting behaviour is that for very small bodies (in the sub–kilometer range), orbital phasing due to gas drag is very efficient and, as a consequence, the differential phasing between objects of different sizes is very strong. This leads to high , especially when compared to the very low erosion threshold velocities corresponding to such small objects.

It is important to note that our results can be rescaled to different values of the reference gas density . Indeed, the behaviour of a population of planetesimals of size in a gas disk of density is similar to the behaviour of a population of planetesimals of size in a gas disk of density (see Equs.1&2). As an example, reducing by a factor 10 the size of the bodies is equivalent to model the evolution of larger bodies but increase the gas density by a factor 10. As a consequence, another way of reading Figs.6 and 5 is that in denser circumbinary disks, planetary formation may be easily inhibited even when starting from “standard” km initial planetesimals. Note however that the rescaling is not fully straightforward, because 5km and 0.5km objects do not have the same threshold velocity .

Figure 4: Same as in Fig.3 but assuming the minimum, and thus more restrictive value of .
Figure 5: Same as in Fig.3 (“optimistic” high value) but assuming a minimum value of equal to 0.5 km (small planetesimals).
Figure 6: Same as in Fig.4 (restrictive low value), but assuming a minimum value of equal to 0.5 km (small planetesimals).

5 Conclusions and Perspectives

With a full numerical approach we have investigated planetesimal accretion and the possibility of planet formation in circumbinary disks. The trajectories of 10000 planetesimals have been integrated over time. The evolution of one crucial parameter, i.e. the distribution of mutual encounter velocities within the system, is numerically explored for a wide range of binary parameters, specifically the mass ratio and orbital eccentricity , the semi-major axis being fixed at 1AU. From these data we map the minimum radial distance from the barycenter of the two stars beyond which accretion is always possible as a function of the binary mass ratio and eccentricity. Planets can form only in regions of the disk beyond this distance.

Our results significantly differ from those of previous studies (Moriwaki & Nakagawa, 2004) mainly for two reasons:

  • The use of the formula that relates the forced eccentricity of the binary companion to the average collisional velocity does not work properly when there is a significant radial mixing of planetesimals. The comparison of our full numerical results with those obtained from the forced eccentricity differ by almost a factor two. The estimate based on the forced eccentricity leads to artificially higher impact velocities.

  • We include in our simulations the crucial effect of the drag force due to the gaseous component of the disk. This significantly affects the relative velocities between planetesimals: the between similar-size bodies are significantly reduced by the alignment of the periapses, but different-size planetesimals are aligned towards different directions and they may experience higher impact velocities.

Summarizing all our results on equal and different-size planetesimals we outline the regions of the disk where accretion is always possible. We find that for a “standard” system where planetesimals have an initial size ofkm, accretion is always possible, in the considered AU region, for equal-mass binaries. For lower mass ratios, results significantly depend on the threshold velocity above which an impact by an impactor on a target leads to net erosion instead of net accretion. A general result is that lower and higher values result in higher and are thus less favorable to planetesimal accretion.

One aspect of the problem that needs further investigation is the initial orbital distribution of planetesimals. This parameter is linked to the difficult issue of when and how they detach from the gas. While fully formed planetesimals may be marginally affected by the spiral waves of the disk caused by the binary, smaller bodies would tend to follow the gas streamlines when they are small enough to be coupled to the gas. This is a complex problem since it is strongly tied to the planetesimal formation process. How is the dust aggregation process influenced by the gas waves and by the stellar secular perturbations? Is the gravitational instability still possible? These questions may be answered within a model capable to handle gas and dust at the same time. Such a complete model not being available yet, we chose to make in our simulations the simplifying assumption that at planetesimal orbits are circular, i.e., that their proper eccentricity is equal to the forced one. However, we think that this choice may not have too dramatic consequences on the estimation of average relative velocities. Indeed, we observe that the further evolution of the system, under the secular perturbations of the stars and of gas drag , leads to a steady state where the proper eccentricity is damped down and the periapses are aligned depending on body size. This steady state only weakly depends on the initial orbital parameters, especially for smaller planetesimals. In addition, since the major contribution to impact velocities is due to the different alignment of different-size planetesimals, rather than due to the absolute value of eccentricity, another choice of initial planetesimal orbits may not lead to significantly different results.

We considered a binary semimajor axis of 1 AU in our simulations, but equation 3 for the gas-free case can be easily scaled to larger (or smaller) values of . However, when we include gas drag, such an easy scaling is no longer possible. The dynamical evolution of planetesimals depends on the local gas density, which in turns depends on the radial distance to the star. As a consequence, the evolution of a planetesimal at a distance of the centre of mass of a binary of separation is not equivalent to the evolution of the same planetesimal at distance around a binary of separation (even when corrected by a timescale correcting factor), because the gas density is not the same at these 2 locations. Furthermore, there is no reason to assume that the global gas disk profile is the same around binaries with different orbital separations. As an example, if for AU it may be reasonable to adopt the minimum mass solar nebula density at 4 AU, in a disk around a binary with semimajor axis of AU, adopting the same density at 20 AU appears inconsistent. In other words, while the relative velocity in absence of gas scales with because of its dependence on the Keplerian velocity, its density does not necessarily scale with when gas drag is considered. As a consequence, our results may be applied to disks around binaries with small separations, but they cannot be scaled to large .

An additional effect that we neglect in our simulations is the physical outcome of mutual collisions (be it accretion or erosion). As noted in the introduction, a model capable to handle at the same time the collisional and dynamical evolution of a swarm of planetesimals is still far from the present computing capabilities.

Our work must be considered a step forward of a better understanding of planet formation in circumbinary disks. There are certainly aspects of the problems that need further investigation in future papers like planetary accretion in habitable regions around binaries, planetary migration driven by disk interaction and by planetesimal encounters.


  • Adachi et al. (1976) Adachi, I., Hayashi, C., Nakazawa, K., 1976, Prog. Theor. Phys. 56, 1756-1771
  • Artymowicz et al. (1994) Artymowicz, P., Lubow, Stephen H., 1994, ApJ 421, 651-667
  • Benz & Asphaug (1999) Benz, W., Asphaug, E., 1999, Icarus 142, 5-20
  • Bodenheimer & Lin (2002) Bodenheimer, P., Lin, D. N. C., 2002, Annual Review of Earth and Planetary Sciences 30, 113-148
  • Boss (1997) Boss, A. P., 1997, Science 276, 1836-1839
  • Close et al. (1998) Close, L.M., Dutrey, A., Roddier, F., Guilloteau, S., Roddier, C., Northcott, M., Menard, F., Duvert, G., Graves, J. E., Potter, D., 1998, ApJ 499, 883-888
  • Correia et al. (2005) Correia, A. C. M., Udry, S., Mayor, M., Laskar, J., Naef, D., Pepe, F., Queloz, D., Santos, N. C., 2005, A&A 440, 751-758
  • Desidera & Barbieri (2007) Desidera, S., Barbieri, M., 2007, A&A 462, 345-353
  • Dutrey et al. (1994) Dutrey, A., Guilloteau, S., Simon, M., 1994, A&A 286, 149-159
  • Greenberg et al. (1978) Greenberg, R., Hartmann, W. K., Chapman, C. R., Wacker, J. F., 1978, Icarus 35, 1-26
  • Gunther & Kley (2002) Gunther, R., Kley, W., 2002, A&A 387, 550-559
  • Hayashi (1981) Hayashi, C., 1981, Prog. Theor. Phys. Suppl. 70, 35-53
  • Holman & Wiegert (1999) Holman, M. J., Wiegert P. A., 1999, AJ 117, 621-628
  • Holsapple (1994) Holsapple, K. A., 1994, Planetary and Space Science 42, 1067-1078
  • Jensen & Mathieu (1997) Jensen, E. L. N., Mathieu, R. D., 1997, AJ 114, 301-316
  • Kary et al. (1993) Kary, D. M., Lissauer, J. J., Greenzweig, Y., 1993, Icarus 106, 288-307
  • Kokubo & Ida (2000) Kokubo, E., Ida, S., 2000, Icarus 143, 15-27
  • Lissauer (1993) Lissauer, J. J., 1993, Ann. Rev. Astron. Astrophys. 31, 129-174
  • Marzari et al. (1995) Marzari, F., Davis, D., Vanzani, V., 1995, Icarus 113, 168-187
  • Marzari et al. (1997) Marzari, F., Scholl, H., Tomasella, L., Vanzani, V., 1997, Planetary and Space Science 45, 337-344
  • Marzari & Scholl (1998) Marzari F., Scholl H., 1998, Icarus 131, 41-51
  • Marzari & Scholl (2000) Marzari F., Scholl H., 2000, ApJ 543, 328-339
  • Mathieu et al. (2000) Mathieu, R. D., Ghez, A. M., Jensen, E. L. N., Simon, M., 2000, Protostars and Planets IV (Book - Tucson: University of Arizona Press; eds Mannings, V., Boss, A.P., Russell, S.S.), p. 703
  • Moriwaki & Nakagawa (2004) Moriwaki, K., Nakagawa, Y., 2004, ApJ 609, 1065-1070
  • Nelson (2003) Nelson, R. P., 2003, MNRAS 345, 233-242
  • Paardekooper (2007) Paardekooper, S., 2007, A&A 462, 355-369
  • Pichardo et al. (2005) Pichardo, B., Sparke, L. S., Aguilar, L. A., 2005, MNRAS 359, 521-530
  • Pollack et al. (1996) Pollack, J. B., Hubickyj, O., Bodenheimer, P., Lissauer, J. J., Podolak, M., Greenzweig, Y., 1996, Icarus 124, 62-85
  • Quintana & Lissauer (2006) Quintana, E. V., Lissauer, J. J., 2006, Icarus 185, 1-20
  • Thébault et al. (2003) Thébault P., Augereau, J.-C., Beust, H., 2003, A&A 408, 775-788
  • Thébault & Doressoundiram (2003) Thébault, P., Doressoundiram, A., 2003, Icarus 162, 27-37
  • Thébault et al. (2006) Thébault, P., Marzari, F., Scholl, H., 2006, Icarus 183, 193-206
  • Udry et al. (2002) Udry, S., Mayor, M., Naef, D., Pepe, F., Queloz, D., Santos, N. C., Burnet, M., 2002, A&A 390, 267-279
  • Wetherill & Stewart (1989) Wetherill, G. W., Stewart, G. R., 1989, Icarus 77, 330-357
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

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

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