Free Access
Issue
A&A
Volume 530, June 2011
Article Number A103
Number of page(s) 10
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/201016375
Published online 19 May 2011

© ESO, 2011

1. Introduction

Although it is believed that approximately half of all stars belong to multiple stellar systems (e.g. Duquennoy & Mayor 1991),  ~90% of exoplanets are associated with single stars (Zsom et al. 2011). It is not yet clear whether this discrepancy is solely due to observational bias, or if the process of planetary formation may be seriously impaired even in very wide binary systems. Curiously, however, a few exoplanets have also been detected in very tight binary stellar systems, where the gravitational perturbations of the secondary component are so large that accretional collisions among small planetesimals are extremely difficult.

Perhaps the most extreme case is γ-Cephei, a binary stellar system whose most recent orbital determination (Neuhäuser et al. 2007) shows a secondary component of mass mB ~ 0.4   M orbiting a principal star mA ~ 1.4   M in an ellipse with semimajor axis aB ~ 20.2 AU and eccentricity eB ~ 0.41. Although both stars have a minimum mutual distance of only  ~12 AU, a giant planet has been detected at  ~2 AU from mA (Hatzes et al. 2003). So far, all attempts to understand the accretional history of this extrasolar planet have been unsuccessful, and it is difficult to visualize an scenario under which such a massive Jovian planet could form through accretional collisions from a primordial planetesimal swarm.

The γ-Cephei system then constitutes a paradigm. It may be argued that if we are able to comprehend planetary formation in such an extreme environment, we would have taken large steps towards a global understanding of planetary formation in any other system. It is then no surprise that this planetary system has caught the attention of several researchers over the past decade, and many dynamical and collisional studies may be found in the literature (e.g. Thébault et al. 2004, 2006; Haghighipour 2006; Verrier & Evans 2006; Tsukamoto & Makino 2007; Paardekooper et al. 2008; Xie & Zhou 2008; Kley & Nelson 2008; Paardekooper & Leinhardt 2010).

In Beaugé et al. (2010) we showed that the dynamics of the gas disk plays an important role in determining the evolution of small planetesimals. We found that, although an apsidal precession of the gas elements may play a disruptive role, especially in the inner parts of the disk, the combined effects of a non-zero precession rate plus a high forced eccentricity in the disk may in fact lower the relative velocity of solid bodies in the outer regions. We thus envisioned a possible scenario in which the growth from kilometer-size planetesimals to  ~100 km planetary embryos could initially take place near the truncation radius of the gas disk. As the embryos spiral down towards a lower semimajor axis, subsequent collisions could then lead to larger embryos, finally forming a giant planet core hopefully near the present location of the currently detected Jovian planet.

To test this idea, we performed preliminary N-body simulations of the dynamical and collisional evolution of planetesimal swarms in the outer regions of the gas disk. However, it was soon realized that the dynamics of the solid bodies in this region is extremely complex. First, we found that mean-motion resonances with the secondary star have significant effects, even though they are of very high order (e.g. 12/1, 11/1, 10/1). Second, the proximity to the secondary star also affects the secular dynamics, and the classical analytical models (e.g. Heppenheimer 1978) widely used in these problems fail to reproduce the correct orbital variations. The imprecision does not lie in the reduced expression adopted for disturbing potential, but in the averaging process used to eliminate short-period terms. Similar to the irregular satellites around the outer planets of our solar system (e.g. Ćuk & Burns 2004; Beaugé et al. 2006), higher order secular effects from the interaction of short-period terms (including evection) must be considered to represent the dynamics of planetesimals in close binary systems. However, contrary to the satellite problem, here the perturber lies in a high-eccentricity orbit, and classical high-order models cannot be applied directly (Correa Otto et al. 2010).

The purpose of this work is to present a general description of the secular dynamics of small planetesimals in circumstellar motion around the more massive star of a generic binary stellar system. We assume that all bodies share the same orbital plane. Although it is known that even moderate mutual inclinations can have significant effects in the accretional evolution of a planetesimal swarm (e.g. Marzari et al 2009a; Xie & Zhou 2009; Xie et al. 2010; Fragner et al. 2011), here we concentrate on the planar case and leave the extension to 3D for a future work. Finally, although our model will be generic, we apply the results to the particular case of γ-Cephei where we analyze how the uncertainties in the orbital fits of the secondary stellar companion may affect the evolution of a planetesimal swarm.

The paper is organized as follows. Section 2 presents our second-order perturbation model and analytical approximations for both the forced eccentricity and secular frequency. Since we focus our attention on γ-Cephei, Sect. 3 discusses the orbital parameters determined for the two stellar components and their precision. The secular dynamics of individual planetesimals, under the additional effects of a nonlinear gas drag, is analyzed in Sect. 4. We also present an example of resonance trapping obtained with divergent migration. Finally, discussions close the paper in Sect. 5.

2. Secular dynamics

Let us assume a small planetesimal of mass m in circumstellar motion around a star of mass mA, which is in turn part of a binary system with a smaller component of mass mB. Let aB be the mA-centric semimajor axis of mB and eB its orbital eccentricity. We further assume that all motion occurs in a plane.

Neglecting the gravitational effects of m on both stellar bodies, the orbit of mB will be a fixed ellipse, while the motion of the small planetesimal will be perturbed by the gravitational effects stemming from the secondary component. Thus, in our dynamical system, mB will play the role of the perturber, while m will be the perturbed mass.

2.1. The first-order secular model

Outside any significant mean-motion resonance, the orbital evolution of m will be dominated by the secular perturbations, and the short-period terms (associated to the mean longitudes) can be eliminated by a perturbation technique known as averaging. The expression for the secular disturbing function R usually employed for these studies was originally developed by Heppenheimer (1978) which, except for constant terms, is given by (1)(see Terquem & Papaloizou 2002), where is the gravitational constant, a is the mA-centric semimajor axis of the planetesimal, e its eccentricity, and ϖ its longitude of pericenter. The angle ϖB denotes the longitude of pericenter of the orbit of mB, assumed constant.

Expression (1) is constructed from Kaula’s (1962) expansion of the disturbing potential, truncated to second-order expansion in the eccentricity of the perturbed body, and performing a first-order “scissors” averaging (with respect to the masses) in the mean longitudes. We refer to the resulting expressions as a first-order model for the secular dynamics.

Since R does not depend explicitly on the mean longitude λ of the planetesimal, its semimajor axis is constant and equal to the proper value a. Consequently, the secular system is reduced to a single degree of freedom, and the differential equations governing the regular variables k = ecos(ϖ − ϖB) and h = esin(ϖ − ϖB) can be written as (2)where

Given arbitrary initial conditions (k0,h0), these equations admit periodic solutions of the form (5)where g is the secular frequency, is usually known as the proper (or free) eccentricity, and the phase angle is given by the expression tanφ0 = h0/(k0 − ef). The constant term ef is known as the forced eccentricity and is only present in systems with an eccentric perturber. Adopting fixed values for aB and eB, Eq. (4) implies that ef is a linear function of the proper semimajor axis (ef ~ a) while the secular frequency scales as g ~ a3/2.

2.2. Numerical simulations

Our first task is to assess the accuracy of the secular solutions (5) corresponding to the first-order model (1). Two quantities we particularly wish to test are ef and g. The forced eccentricity is crucial in determining the equilibrium eccentricity of planetesimals under the effects of gas drag from the protoplanetary nebula. Although any secular oscillatory motion is expected to be damped in a gas-rich scenario, the magnitude of g is important for establishing the validity of the averaging process of the disturbing function.

For our computations, we assume a generic binary system with mass ratio between the components of mB/mA = 0.25 and eccentricity eB = 0.36. This value is similar to the best-fit solution found by Hatzes et al. (2003) for γ-Cephei.

thumbnail Fig. 1

Forced eccentricity (top) and secular frequency (bottom), as function of the proper semimajor axis, calculated by three different methods: filtered exact numerical simulations (filled black circles), semi-analytical first-order averaging of the exact disturbing function (dashed lines), and the classical analytical first-order secular model (continuous lines).

Open with DEXTER

Figure 1 shows three different calculations of the forced eccentricity (top frame) and the secular frequency (bottom frame). The value of ef appears to grow linearly with the proper semimajor axis, reaching values of  ~ 0.1 for a ~ 0.2aB. The numerical calculations were obtained from a long-term integration of the exact equations of motion, after an online application of a low-pass FIR (finite impulse response) filter (e.g. Carpino et al. 1987).

A digital filter is a numerical tool that eliminates certain frequencies from an input signal. For example, given a certain time series (e.g. eccentricity as function of time) and a pass frequency νpass, applying a low-pass filter signal will yield an output that maintains all the periodic variations with frequencies ν < νpass while eliminating the rest. Digital filters are a common tool for constructing synthetic theories of long-term asteroid dynamics (e.g. Knezeviíc & Milani 2000) and planetary dynamics (e.g. Michtchenko & Ferraz-Mello 2001), and constitute a useful alternative to analytical perturbation theories when the Hamiltonian function is very complex.

For the present work, the parameters of the digital filter were chosen to eliminate all periodic variations with up to eight times the orbital period of the binary component. In a dynamical system displaying regular motion, the application of the filter is equivalent to a full (i.e. infinite-order) averaging of the Hamiltonian. The region located beyond a ~ 0.2aB shows dynamical instabilities that complicate the determination of the secular solution.

A comparison between the numerical and the analytical values shows significant differences. Although both methods yield similar results for low values of the semimajor axis, the exact secular frequencies are systematically underestimated by the analytical model, leading to differences of almost a factor of two for a/aB ~ 0.24. This limitation in the classical estimation of g was first noticed by Thébault et al. (2006), who presented an empirical functional correction term to expression (3). This correction reduced the discrepancy to values of around 5% in the same range of a.

Perhaps more important is that the value of the forced eccentricity also shows significant differences. While the analytical model predicts a monotonic increase in ef as function of a, the real value appears to reach a plateau around a/ab ≃ 0.17 (corresponding to ef ≃ 0.07) and decrease for larger radial distances. The scatter in the numerical values of both ef and g in this outer region stems from the action of high-order mean-motion resonances between the massless body and the binary star mB.

At first hand, it seems natural to believe that the limitations of expressions (3), (4) are due to the truncation of the disturbing function to second-order terms in the eccentricities and/or third-order terms in the ratio a/aB. However, this is not the case. In Fig. 1 we have also plotted the same quantities determined using a semi-analytical model for the disturbing function. This expression is calculated directly as (6)where r and are the position vectors of m and mB, respectively, r and rB are their absolute values, λ and λB are the mean longitudes, and φ is the instantaneous angular distance between both bodies. The integrand is the exact expression for the disturbing function with no approximations, and the double integral is performed numerically. From this expression the value of ef can be estimated from the minimum value attained by  ⟨ R ⟩  in the line segment (ϖ − ϖ0) = 0, while the secular frequency is given by (7)at the same point. Here is the modified Delaunay canonical momenta conjugate to the longitude of the pericenter. Expression (7) may also be evaluated numerically for any initial condition.

This type of semi-analytical model has proved to be a powerful tool for mapping the phase space of complex dynamical systems, especially in the high-eccentricity regime where analytical approximations for the Hamiltonian are not available (e.g. Michtchenko & Malhotra 2004; Michtchenko et al. 2006). Formally, it is equivalent to a first-order averaging (in the masses) of the exact Hamiltonian function.

Figure 1 shows the values of both ef and g determined with this semi-analytical approach. Although the value of the secular frequency shows a significant improvement over the analytical estimation, there is still a discrepancy with the exact value. This is even more noticeable in the forced eccentricity, where there is practically no difference with the value determined from Eq. (3). Consequently, it appears that the limitations of the analytical model are not due primarily to the truncation of the disturbing function.

2.3. A second-order secular model

Since the errors in the estimation of both the forced eccentricity ef and the secular frequency g do not come from the limitations in the adopted disturbing function, their origin must lie in the construction of the averaged solution itself. As mentioned previously, Heppenheimer’s (1978) expressions are a first-order model with respect to the perturbing mass. Here, we extend the calculations to the second order.

One of the most widely used perturbation techniques is the so-called Hori’s averaging process (Hori 1966; see also Ferraz-Mello 2007), which employs Lie-type canonical transformations to eliminate the dependence of the Hamiltonian with respect to a given set of variables. The new Hamiltonian function is given by a power series in the small parameter (e.g. perturbing mass).

Since we adopt a Hamiltonian formulation, we first need to introduce canonical variables. We have chosen the modified Delaunay variables (L,Λ,G − L,λ,λB), where the canonical momenta are given in terms of the orbital elements, by (8)and Λ is the canonical conjugate of the mean longitude of the perturbing mass (i.e. λB). This third degree of freedom appears when passing to the extended phase space to eliminate the non-autonomous character of the perturbation.

The full Hamiltonian function governing the dynamics of the planetesimal m is given by (9)where nB is the mean-motion of the perturbing mass mB and R the disturbing function. We can express this Hamiltonian in a form adequate for perturbation theory: F = F0 + εF1, where (10)and is a small parameter that serves as a guide of the relative magnitudes between the perturbation term F1 and the unperturbed integrable Hamiltonian F0.

For the disturbing function, we adopt a Legendre expansion, truncated to fourth order in the ratio a/aB; in other words, we approximate the perturbation by (11)where Pi(cosφ) is the Legendre polynomial of degree i. Switching from a power series in cosφ to a harmonic decomposition in φ and transforming them to orbital elements, we can obtain a truncated expansion of the disturbing function leading to (12)where M and MB are the mean longitudes of both bodies, and Di,j,k,l may be obtained in terms of the Hansen coefficients (see Beaugé & Michtchenko 2003).

Having an explicit expression for F1 in mean variables, we may now apply Hori’s method. The idea is to search for a Lie-type canonical transformation B = εB1 + ε2B2 + ... to a new set of variables such that the transformed Hamiltonian F is independent of λ and . Up to second order in the small parameter, the new Hamiltonian function may be written as (13)where Δϖ = ϖ − ϖB. The different orders in expression (13) are given by (14)where  {  }  is the Poisson bracket,  ⟨  ⟩ λ,λB denotes the averaging with respect to both mean longitudes (keeping all other variables fixed), and B1 is the first-order generating function of Hori’s method. In terms of the adopted expansion for the disturbing function (12), it is given by (15)where the function must be evaluated in the new variables.

The construction of the new secular Hamiltonian F((G − L)ϖ;L) is cumbersome, although fairly straightforward when using an algebraic manipulator. Fortunately, it will not be necessary to write an explicit expression here. Let it suffice to say that F constitutes a second-order model of the secular system and a single degree of freedom system in variables ((G − L)ϖ). Employing the inverse transformation from Delaunay variables to orbital elements, we can also obtain an expression for F(eϖ;a) in terms of the mean eccentricity e and the proper semimajor axis a. Since the latter orbital element is constant, it appears in the Hamiltonian as an external parameter.

Finally, after solving the secular system and obtaining both e and Δϖ as functions of time, we may invoke the inverse Hori transformation to obtain the short-period variations of the original osculating variables. For the eccentricity, this yields (16)Because B1 explicitly depends on the mean longitudes, the second term models the short-period variations in the eccentricity, while the first term (e2(t)) gives the main secular contributions. Since the eccentricity is a positively defined function, the magnitude of the second term also specifies the minimum mean eccentricity e of the secular system for any given proper semimajor axis a. At the same time, it also gives the averaged semi-amplitude of the short-period variations Δe in the same orbital element.

thumbnail Fig. 2

Forced eccentricity as function of the proper semimajor axis, calculated by three different methods: filtered exact numerical simulations (filled black circles), first-order analytical model (dashed lines), and the new second-order secular model (continuous lines).

Open with DEXTER

Figure 2 shows an application of our second-order model to the same generic binary system as was discussed in Fig. 1. The plot shows the forced eccentricity, as a function of the ratio a/aB, calculated with three different methods. Recall that the first-order model predicts a linear increase of ef with the semimajor axis. Finally, the value of the forced eccentricity determined with our second-order Hamiltonian F is shown as a continuous curve. The agreement with the numerical data is very good, and the saturation in the value of ef is reproduced quite well. Since we have avoided all small denominators in the generating function B1, the model curve is smooth and shows no indication of the effects of mean-motion resonances.

2.4. Extending the Thébault et al. (2006) approximation

As mentioned before, although the second-order model leads to significant improvement in the secular solution, as well as allowing the magnitude of the short-period orbital variations to be modeled, it is much too complex to constitute a workable model. For this reason, we wondered whether the empiric correction term introduced by Thébault et al. (2006) for the secular frequency could be extended to reproduce both the forced eccentricity and the short-period variations. Of course it is not expected to yield the exact same results, but if the errors are not significant, such an empirical second-order approximation could constitute a simple quantitative analytical model.

Following the same approach as Thébault et al. (2006), we use ef0 and g0 to denote the first-order expressions for the forced eccentricity and secular frequency, and reserve ef and g for the second-order values. The idea then is to write ef = ef0(1 + εδef) (and a similar equation for g), and attempt to model the correction terms δef and δg. After several tests and multivariate linear regressions, we find that the expressions (17)agree with the complete second-order model very closely. There are some slight differences in g with respect to the original formula introduced by Thébault et al. (2006) but they are minor and not very significant. Finally, the expressions for ef0 and g0 are those given in (3) and (4).

In terms of (17) the secular Hamiltonian may be approximated well by (18)where k = ecos(Δϖ) and h = esin(Δϖ) are the new regular secular variables.

thumbnail Fig. 3

Variation in the forced eccentricity (top) and secular frequency (bottom), in terms of the proper semimajor axis, for three values of the binary eccentricity eB. As before, filled circles present results from filtered exact numerical simulations, dashed lines correspond to the first-order analytical model, while the empirical solutions (17) are shown in continuous lines.

Open with DEXTER

Finally, the semi-amplitude of the short-period variations in eccentricity can also be empirically modeled according to the expression (19)

Figure 3 once again compares the estimated values of ef and g, this time for three different values of the eccentricity eB of the binary component. Given the simplicity of these equations, the agreement with the numerical results is surprisingly good.

3. The γ-Cephei binary system

After specifying the basic ingredients of our second-order dynamical model, we attempt to apply it to γ-Cephei. As mentioned in the introduction, this is probably the best-studied tight binary system with a known planetary body. Since the main secular parameters g and ef strongly depend on the stellar masses and orbital elements of the secondary star, we begin our discussion by reviewing the accuracy of these parameters.

Throughout this work we refer to the more massive stellar component by γ-Cephei-A, while γ-Cephei-B is used to identify the less massive star. The giant planet orbiting γ-Cephei-A is called γ-Cephei-b. The masses of each body are denoted by mA, mB, and mp, in that order.

3.1. History and radial velocity data

Several years before the discovery of the first planetary body around a main sequence star (Mayor & Queloz 1995), Campbell et al. (1988) suggested the presence of a Jupiter-mass object in a 2.7 yr orbit around γ-Cephei-A. The authors, however, remained cautious about claiming a true planetary detection, since the observed periodic variations in radial velocity (RV) were at the very limit of the instrumental resolution. To complicate the problem even further, the variations in RV attributed to the Jovian planet, with a semi-amplitude of only about 25 m/s, were superimposed on a much larger variation caused by a previously unnoticed stellar companion with a much longer orbital period.

The planetary interpretation was questioned later when changes in the chromospheric activity were observed with similar period (Walker et al. 1992). Thus, it was proposed that the observed changes in RV were spurious and probably only due to changes in the spectral line profiles caused by surface inhomogeneities (spots).

The existence of a binary component (i.e. γ-Cephei-B) was only reevaluated several years later, when Griffin et al. (2002) combined several historical sources of radial velocities that include epochs from 1896 to 1980. This data set consisted of 88 RV observations, although many of them did not contain proper uncertainties, and a gap of some 50 years was present in the data set. Even so, the authors proposed a secondary stellar mass in the system with an orbital period of P ~ 66 yrs.

The presence of a third body, this time a planet around γ-Cephei-A, was only confirmed by Hatzes et al. (2003) after incorporating new high-precision velocity observations from the McDonald Observatory. They show convincingly that the 2.5 yr variation in radial velocity was coherent in phase and amplitude throughout the entire 20 yr interval, as would be expected for Keplerian motion, and that no changes were observed in the spectral-line bisectors.

More recently, Torres (2007) has again analyzed the historical sources of radial velocities using the extensive Harvard-Smithsonian Center for Astrophysics (CfA) database consisting of  ~250 000 spectra. Torres pointed out that some of the historical radial velocities showed large internal discrepancies when compared with other data taken at similar times and were consequently not reliable. The author constructed a reliable data set consisting of 30 RV observations. The complete sets of radial velocities (four sets by Hatzes et al. 2003; and one by Torres 2007) are shown in Fig. 4, where the errors bars indicate the uncertainty on each numerical value. The difference in precision is remarkable, showing how the incorporation of modern techniques in RV measurements lead to the detection of the planetary mass. This increase in precision also allowed the mean anomaly M and longitude of pericenter ϖ of the binary component to be accurately defined.

thumbnail Fig. 4

The five sets of RV data used for the orbital fit of γ-Cephei-B. The four datasets from Hatzes et al. (2003) are shown in black, while the dataset from CfA by Torres (2007) is shown in gray. The error bars correspond to the observational uncertainties given by the authors. Two orbital solutions are shown, one corresponding to a larger orbit for the binary (top) while the bottom panel represents a more compact configuration. Each plot also assumes a different value for mA.

Open with DEXTER

Independently and without attempting to identify any planetary body, Gontcharov et al. (2000) studied the mass and orbital parameters of the binary system using astrometric observations from several sources. They obtained an orbital period of  ~ 45 yr and a total mass of 3   M. Unfortunately, the individual masses were not specified. In his work, Torres (2007) also combined a total of 140 astrometric measurements obtained between the years 1898 and 1995. Because of the relatively short time span of these observations compared to the binary orbital period, Torres noted that there is a high probability that part of the orbital motion of the binary has been absorbed into the proper motion components reported by Hipparcos. The astrometric information from both these papers is quite different (compare Fig. 6 in Torres 2007 with Fig. 5 in Gontcharov et al. 2000) and the discrepancy in the data prior to 1979 prevent us from using any astrometric data in our analysis of the orbital determination.

Table 1

Published masses and orbital parameters for γ-Cephei.

Although the high-precision RV measurements from Hatzes et al. (2003) give very good definitions of the mass and orbital parameters of the planetary body, the orbit of binary itself is far from established. In Table 1 we summarize the results of four different best fits: not only are there noticeable differences in the semimajor axes, but the stellar masses also show large discrepancies.

The brightness of γ-Cephei has made it an easy target for spectroscopic studies to determine the effective temperature of the star. This parameter, together with the absolute magnitude, are used to determine its mass. However, the effective temperature for the star varies from 4300 to 5100 K (Torres 2007) yielding a wide variety of possible solutions.

3.2. Possible orbital solutions for γ-Cephei

Together with the RV data, Fig. 4 also shows two different orbital fits for the binary system, each constructed for different mA. The top frame adopts the value given by Torres (2007), which is significantly lower than the one employed by Hatzes et al. (2003), shown here in the bottom graph. We recall that the higher mass is usually used for dynamical studies in this system.

From the raw RV data, we redetermined the best fits for each value of mA. To this end we used the PISA code (Pikaia genetic algorithm + simulated annealing) that we developed for our studies of resonant exoplanetary systems (e.g. Giuppone et al. 2009). The minimization procedure implies a determination of ten free parameters: five for orbital parameters and five for the RV offsets. We neglected the presence of the planetary body, since it does not introduce any significant effect on the orbital calculation for the binary system. The values of mB, aB, and eB obtained for each solution are (20)As mentioned previously, the mean anomaly and the longitude of the pericenter are well specified, and show little change in both fits. We assumed an edge-on configuration. This is an important restriction, since any inclination of the orbital plane of the binary would lead to much higher values of mB (see Hatzes et al. 2003; Haghighipour 2006). Although both fits correspond to significantly different dynamical systems, they yield practically the same residuals: .

From the detection of exoplanets from RV data, we have learned (the hard way) that the best fit obtained with a limited data set does not necessarily correspond to the real configuration of the physical system. This is especially true for a small number N of data points, or when they only cover a fraction of the orbital period. In other words, even if we specify a value for mA, the real mass and orbit of mB does not necessarily have to coincide with the best fit.

A way to estimate the possible range of solutions compatible with the observational data is to analyze the shape of the residual function for a series of orbital fits around the minimum value . If this function shows a steep increase for small changes in the fitted parameters, then we may have a certain confidence that the best fit is probably very close to the real system. However, a shallow minimum could lead to a wide diversity of possible solutions with almost the same value of and, consequently, to different configurations that are statistically indistinct.

To test this idea for both choices of mA, we calculated a series of orbital solutions with predefined values of aB between 12 and 25 AU. Except for the pair (mA,aB), all the remaining parameters were free and determined with the minimization process. Results are shown in Fig. 5 for mA = 1.59   M and for mA = 1.18   M. The value of for each fit is shown in the top lefthand plot. Both families of solutions show similar minima, although displaced in semimajor axis.

thumbnail Fig. 5

Several multikeplerian orbital fits for γ-Cephei considering different values of mA and aB. Black lines show results for mA = 1.59   M, while red curves assume mA = 1.18   M. The dashed horizontal line in the upper-left panel corresponds to the 1σ confidence level around the global minimum of . The dashed line in the bottom-left panel corresponds to a constant pericenter distance equal to 11.48 AU.

Open with DEXTER

To estimate the limits of the 1σ confidence level, we employed the same procedure as employed in Giuppone et al. (2009). Given a best-fit algorithm with ν degrees of freedom, the value of associated to the 1σ confidence level is approximately given by (21)where is the minimum value. Since our problem contains 230 data points and 10 free parameters, we have ν = 220. For , Eq. (21) then gives .

Although both values only differ in  ~5%, the range of possible solutions with is surprisingly wide. Assuming the lower value for mA, the semimajor axis may be as small as 15 AU or as large as 22 AU. For mB = 1.59   M, the range is aB ∈  [17,23]  AU.

The top righthand plot of Fig. 5 presents the values of mB that give the best fit for each value of aB. As expected, when the components are more separated, the mass of the binary increases to produce same magnitude in radial velocity. The lower lefthand plot gives the binary eccentricity eB. Larger semimajor axes are accompanied by more elliptic orbits. The curves show a rough resemblance to the locus of constant pericentric distance pB, here drawn for one particular orbital solution. Since most of the RV data points are located near the pericenter of the binary’s orbit (see Fig. 4), the value of pB is much better defined than either aB or eB.

Finally, the lower-right panel shows that the offset from the CfA data varies significantly as a function of the semimajor axis of the binary. This behavior is not noted for the other offsets, which remain almost constant. All the best fit solutions showed almost no variation in the longitude of pericenter nor in the time of passage through the periastron.

Table 2 summarizes the minimum/maximum possible values of mB, aB, and eB leading to orbital fits with residuals within the 1σ confidence level. It is important to stress that statistically all are equally compatible with the observational data. In the next section we attempt to elucidate whether these uncertainties are important to the dynamics of small planetesimals orbiting the primary star and, consequently, whether they could have an effect on the planetary formation process.

4. Secular dynamics of planetesimals in γ-Cephei

4.1. The forced eccentricity

Planetary accretion requires low relative velocities which, in turn, implies similar orbits between colliding bodies. This condition may be satisfied if the orbital eccentricities are: (i) very small or (ii) very similar and the orbits are aligned. For small planetesimals where mutual perturbations are not crucial, the orbital eccentricities are determined by a complex interplay between several phenomena, including gas drag, collisional history, and the gravitational effects of the secondary star (e.g. Marzari & Scholl 2000; Thébault et al. 2006, 2008). In the secular approximation, these effects appear through the magnitude of the forced eccentricity ef.

Thus, one way to study planetary accretion under different orbital solutions for γ-Cephei would be to analyze the sensitivity of ef to the set (mA,mB,aB,eB) compatible with the RV data. Results are presented in Fig. 6 for two values of mA. Each panel shows level curves of constant forced eccentricity, as a function of the semimajor axis a of the planetesimal (abscissa), and for different semimajor axes aB for the binary pair (ordinate).

Contrary to our expectations, the range of eccentricities appears insensitive to the binary configuration. In all cases the values of ef extend from  ~0.03 for the small semimajor axis to  ~0.077 for a ≃ 4 AU. Adopting a lower value of mA seems to cause a slight reduction in the interval, but the change is not very significant. Consequently, and at least from this preliminary analysis, there appears to be no indication that different configurations for γ-Cephei could cause major changes in the secular dynamics of small bodies and, therefore, on the accretional possibilities of a planetesimal swarm.

Table 2

1σ Confidence limits for γ-Cephei-B.

thumbnail Fig. 6

Forced eccentricity ef, as a function of the semimajor axes of the binary aB and the planetesimal a, for all the orbital solutions of γ-Cephei with residuals within the 1σ confidence level of the best fit. Top and bottom frames assume two different values for mA. Each value of aB implies different values of both mB and eB, obtained from the families of orbital fits presented in Fig. 5. Values of aB for each best fit are shown with horizontal dashed lines.

Open with DEXTER

With hindsight, perhaps this result is not at all unexpected. Since all orbital solutions for γ-Cephei lead to practically the same amplitude in the RV signal, it is understandable that different values for the set (mA,mB,aB,eB) should also generate similar perturbative effects on other hypothetical bodies in the system; e.g. small planetesimals orbiting mA.

4.2. Simulations of individual particles with gas drag

A better test for the effects of different orbital fits on the secular dynamics is to compare the evolution of small planetesimals under the effects of a nonlinear drag force from a circumstellar gas disk centered on mA. We employ the same expression for the dissipative force as discussed in Beaugé et al. (2010). For the gaseous disk we assume a linear surface density profile with a total mass of 3 Jupiter masses, and an outer edge located at 5 AU. We consider planetesimals with a volumetric density of ρ = 3 g/cm3.

Figure 7 shows the result of numerical simulations of four different test planetesimals, with radius between s = 0.5 km (top) and s = 10 km (bottom). All were initially placed in circular orbits. The initial semimajor axis was equal to a = 4 AU for the first three integrations, and a = 3 AU for the largest planetesimal. We assumed a gas disk with constant eccentricity of eg = 0.05 and a rigid-body retrograde precession with a period of 1000 yrs.

We considered two different orbital fits for the stellar components. The mass and orbital parameters for the secondary star were taken from Eqs. (20) which are the best fit solutions for each value of mA.

thumbnail Fig. 7

Orbital evolution of four different size planetesimals under the effects of a nonlinear gas drag in the γ-Cephei system. Black dots correspond to mA = 1.59   M and gray to mA = 1.18   M. In both cases the parameters of mB are those given by the best fits and detailed in Eqs. (20). The gaseous disk has a constant eccentricity of ef = 0.05 and a rigid retrograde precession rate with period 2π/|gg| = 1000 yr. In the lefthand panels, the orange curve shows the forced eccentricities as a function of the semimajor axis. In the righthand plots, the orange curves marks ϖ = 0.

Open with DEXTER

The lefthand panels show the orbital eccentricity as a function of the semimajor axis. The orange curves indicate the forced eccentricity ef, as obtained from our second-order model, for each value of mA. The righthand plots show the evolution of the longitude of pericenter ϖ. The origin was shifted so that ϖB = 0.

We note the existence of three distinct regions in the semimajor axis domain. For a ≲ 2 AU, the planetesimals show short-period oscillations around equilibrium values eeq and ϖeq. These equilibrium values depend on size. Small planetesimals show eeq < ef and ϖeq < 0. Larger bodies, however, appear coupled with the conservative secular solution (k,h) = (ef,0).

The second region lies roughly between 2 AU and 3 AU, and is characterized by very similar equilibrium values of both e and ϖ for all planetesimals with radius s ≳ 0.1 km. Moreover, the equilibrium eccentricities are virtually indistinguishable from the forced eccentricity ef. If planetary accretion is possible with these disk parameters, this region appears to be the most promising since orbital dispersion is kept at a minimum. The limit between both regions (here at a ~ 2 AU) depends on the surface density profile adopted for the disk. In our simulations we used a linear dependence with the radial distance (Beaugé et al. 2010), which translates to high densities close to the central star and low values beyond 3 AU (see Fig. 3 of Kley & Nelson 2008).

A third region is located beyond  ~3 AU. Although the secular dynamics also appear similar for all values of s, the simulations show large-amplitude oscillations. These are not only caused by short-period terms but also from several high-order mean-motion resonances between the particles and mB. In the plots these commensurabilities can be seen as spikes where the eccentricity is temporarily excited. Without a detailed analysis, it is not possible to establish whether these commensurabilities will inhibit or favor accretion. Although most of our simulations have shown scattering effects and significant orbital misalignment between resonant and non-resonant orbits, we have also found cases of resonant trapping. This appears to be a high-probability outcome for s ≳ 5 km.

thumbnail Fig. 8

Example of trapping of a s = 10 km planetesimal in a 10/1 mean-motion resonance with the secondary star of γ-Cephei, due to a nonlinear gas drag. Although the orbital migration is divergent, capture still occurs and leads to an apparently stable configuration. The resonant angle is θ = 10λB − λ − 9ϖB.

Open with DEXTER

An example is shown in Fig. 8 for a 10 km body placed in an initial circular orbit with a = 4 AU. After an initial decay in the semimajor axis, the body is captured in a 10/1 mean-motion resonance with the binary component. Both the resonant angle (22)and the difference in longitudes of pericenter librate around zero, although there seems to be a slow departure towards asymmetric librations at the end of the simulation. The resonant solution seems very stable, at least for timescales between 106 and 107 years. We have found similar outcomes in other commensurabilities, such as the 11/1 and 12/1, always for slow dissipative effects (i.e. large bodies). What is curious about this result is that trapping occurs during a divergent migration; in other words, the non-conservative force increases the separation of the bodies involved. Classical works (e.g. Neishtadt 1975; Henrard 1982; Beaugé & Ferraz-Mello 1993; Nelson & Papaloizou 2002) predict that capture is only possible in cases of convergent migration and, thus, the behavior shown in Fig. 8 should not occur.

The only reference we have been able to find describing similar findings is an abstract of a presentation in the 2007 DDA meeting (Hamilton & Zhang 2007). Although no details are available, it appears that divergent trapping is possible in high-order mean-motion resonances and with high-eccentricity perturbers. This may explain why the librating resonant angle includes the longitude of pericenter of the perturber instead of the planetesimal. However, a deeper analysis is necessary before we are able to understand this phenomena and establish its possible importance in planetary formation.

The results shown in Fig. 7 for both binary configurations show almost identical results. Although the forced eccentricity is slightly higher for mA = 1.18   M for most of the semimajor axis domain, the same three regions exist in both cases, and it is plausible to assume that the collisional evolution of a planetesimal swarm should be similar. Figure 9 shows a second series of integrations with the same four planetesimals as before, but this time we considered no precession for the gaseous disk. We also assumed that ϖg − ϖB = π, i.e. the disk is anti-aligned with the binary companion. This is consistent with the hydrosimulations of Marzari et al. (2009b) for disks with significant self-gravity. Although the magnitude of the resonant and short-period oscillations appear larger for mA = 1.18   M, the averaged secular dynamics is similar in both cases, and the same three regions discussed previously are also present.

thumbnail Fig. 9

Same as Fig. 7, but assuming a static gas disk anti-aligned with the secondary star mB.

Open with DEXTER

5. Conclusions

In this paper we have presented a second-order theory for the secular dynamics of massless particles orbiting a central star and perturbed by a secondary stellar component with high eccentricity. Only coplanar motion is considered. This dynamical problem is applicable to the motion of small planetesimals in tight binary systems such as γ-Cephei. Although the resulting expressions for the forced eccentricity ef and secular frequency g are complex, we were able to extend the empirical approximation originally introduced by Thébault et al. (2006) and deduce simple analytical formulas for both quantities.

The forced eccentricity, in particular, shows significant differences to the classical model (e.g. Heppenheimer 1978). While the first-order equations predict a linear dependence with the semimajor axis, numerical simulations show a quadratic functional form, and ef may actually reach a plateau for high values of a. Our model reproduces this behavior with good precision.

We also analyzed the reliability of the best fits presented in the literature for the stellar components of γ-Cephei. We found that the best solution depends on the adopted mass for the central star, and even for a fixed value of mA there may be many different configurations compatible with the observations. This is expected, since the radial velocity data covers less than one orbital period of the system. However, we have also found that the dynamics of small planetesimals appears to be only weakly dependent on the particular solution adopted for the binary. A comparative study of the evolution of planetesimals under the effects of gas drag from a circumstellar (mA centered) gas disk shows similar evolutionary paths. Although our integrations have not covered many initial conditions or disk parameters, we suspect that the accretional process of a planetesimal swarm should be practically equivalent in any case. Consequently, we believe that the difficulties in explaining planetary formation in tight binary systems cannot be attributed to uncertainties in the orbital fits. The solution must be found elsewhere, and the search is precisely what makes this problem intriguing.

Finally, we also presented a curious case of resonant trapping in divergent migration. There is practically no reference to this behavior in the literature and, although its importance in planetary formation is difficult to evaluate, we nevertheless believe the phenomena is intrinsically interesting and merits further analysis.

Acknowledgments

This work was supported by the Argentinian Research Council – CONICET – and by the Córdoba National University – UNC –. A.M.L. is a postdoctoral fellow of SECYT/UNC.

References

All Tables

Table 1

Published masses and orbital parameters for γ-Cephei.

Table 2

1σ Confidence limits for γ-Cephei-B.

All Figures

thumbnail Fig. 1

Forced eccentricity (top) and secular frequency (bottom), as function of the proper semimajor axis, calculated by three different methods: filtered exact numerical simulations (filled black circles), semi-analytical first-order averaging of the exact disturbing function (dashed lines), and the classical analytical first-order secular model (continuous lines).

Open with DEXTER
In the text
thumbnail Fig. 2

Forced eccentricity as function of the proper semimajor axis, calculated by three different methods: filtered exact numerical simulations (filled black circles), first-order analytical model (dashed lines), and the new second-order secular model (continuous lines).

Open with DEXTER
In the text
thumbnail Fig. 3

Variation in the forced eccentricity (top) and secular frequency (bottom), in terms of the proper semimajor axis, for three values of the binary eccentricity eB. As before, filled circles present results from filtered exact numerical simulations, dashed lines correspond to the first-order analytical model, while the empirical solutions (17) are shown in continuous lines.

Open with DEXTER
In the text
thumbnail Fig. 4

The five sets of RV data used for the orbital fit of γ-Cephei-B. The four datasets from Hatzes et al. (2003) are shown in black, while the dataset from CfA by Torres (2007) is shown in gray. The error bars correspond to the observational uncertainties given by the authors. Two orbital solutions are shown, one corresponding to a larger orbit for the binary (top) while the bottom panel represents a more compact configuration. Each plot also assumes a different value for mA.

Open with DEXTER
In the text
thumbnail Fig. 5

Several multikeplerian orbital fits for γ-Cephei considering different values of mA and aB. Black lines show results for mA = 1.59   M, while red curves assume mA = 1.18   M. The dashed horizontal line in the upper-left panel corresponds to the 1σ confidence level around the global minimum of . The dashed line in the bottom-left panel corresponds to a constant pericenter distance equal to 11.48 AU.

Open with DEXTER
In the text
thumbnail Fig. 6

Forced eccentricity ef, as a function of the semimajor axes of the binary aB and the planetesimal a, for all the orbital solutions of γ-Cephei with residuals within the 1σ confidence level of the best fit. Top and bottom frames assume two different values for mA. Each value of aB implies different values of both mB and eB, obtained from the families of orbital fits presented in Fig. 5. Values of aB for each best fit are shown with horizontal dashed lines.

Open with DEXTER
In the text
thumbnail Fig. 7

Orbital evolution of four different size planetesimals under the effects of a nonlinear gas drag in the γ-Cephei system. Black dots correspond to mA = 1.59   M and gray to mA = 1.18   M. In both cases the parameters of mB are those given by the best fits and detailed in Eqs. (20). The gaseous disk has a constant eccentricity of ef = 0.05 and a rigid retrograde precession rate with period 2π/|gg| = 1000 yr. In the lefthand panels, the orange curve shows the forced eccentricities as a function of the semimajor axis. In the righthand plots, the orange curves marks ϖ = 0.

Open with DEXTER
In the text
thumbnail Fig. 8

Example of trapping of a s = 10 km planetesimal in a 10/1 mean-motion resonance with the secondary star of γ-Cephei, due to a nonlinear gas drag. Although the orbital migration is divergent, capture still occurs and leads to an apparently stable configuration. The resonant angle is θ = 10λB − λ − 9ϖB.

Open with DEXTER
In the text
thumbnail Fig. 9

Same as Fig. 7, but assuming a static gas disk anti-aligned with the secondary star mB.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.