Binary Capture Rates for Massive Protostars
Draft version November 24, 2018
Preprint typeset using LATEX style emulateapj v. 08/22/09
BINARY CAPTURE RATES FOR MASSIVE PROTOSTARS
Nickolas Moeckel, John Bally
Center for Astrophysics and Space Astronomy, and
Department of Astrophysical and Planetary Sciences
University of Colorado, Boulder, CO
Draft version November 24, 2018
ABSTRACT
The high multiplicity of massive stars in dense, young clusters is established early in their evolu-
tion. The mechanism behind this remains unresolved. Recent results suggest that massive protostars
may capture companions through disk interactions with much higher efficiency than their solar mass
counterparts. However, this conclusion is based on analytic determinations of capture rates and esti-
mates of the robustness of the resulting binaries. We present the results of coupled n-body and SPH
simulations of star-disk encounters to further test the idea that disk-captured binaries contribute to
the observed multiplicity of massive stars.
Subject headings: binaries: general — circumstellar matter — stars: formation
1. INTRODUCTION
Massive stars such as those of the Trapezium in
the Orion Nebula cluster (ONC) have a high mul-
tiplicity compared to low mass stars (Mason et al.
1998; Preibisch et al. 1999; Stahler et al. 2000;
Garćıa & Mermilliod 2001). The formation of a
multiple system in the early life of a cluster can occur
through fragmentation of the prestellar material, or by
capture. Capture can occur via multi-body interactions
with other cluster members, or through disk interactions
in the protostellar phase. There is growing observational
evidence for massive, embedded disks surrounding
massive protostars (Cesaroni et al. 2007, for a recent
review). Simulations of massive star formation suggest
that the masses of the disks may build up to values ∼
30% of the central star’s mass before global instabilities
trigger a sudden accretion event (Krumholz et al. 2007).
Fragmentation of such massive disks into companions
is possible (Kratter & Matzner 2006), but in this work
we consider capture of a lower mass cluster member by
a massive star-disk system, continuing the analysis pre-
sented in Moeckel & Bally (2006) and Moeckel & Bally
(2007, hereafter MB07).
Disk assisted capture in low mass systems has been
studied using analytic methods (Clarke & Pringle 1991)
and smoothed particle hydrodynamics (SPH) codes
(Heller 1995; Boffin et al. 1998). The capture rates from
these studies are too low to be a significant contributer to
the binarity of roughly solar mass stars; for 20 M⊙ stars
with a 2 M⊙ disk in a Trapezium-like environment, MB07
show that the capture rate is high enough to account
for 50% binarity after 1 Myr. MB07 estimated analyti-
cally the likelihood of survival for these binaries, and con-
cluded that higher mass captured companions are more
likely to survive, while lower mass companions would be
preferentially ionized by encounters with other cluster
members. However, in calculating capture rates and es-
timating survival fractions, many averages and integrals
are taken. It is desirable to further test these results in
an unaveraged, more ‘realistic’ setting.
Electronic address: moeckel@colorado.edu
In this paper we describe the results of n-body simula-
tions of a cluster similar to the ONC, in which we include
the dissipative effects of passages through a circumstellar
disk surrounding a central 20 M⊙ star. This test demon-
strates that the rates derived by MB07 are reasonable in
a cluster setting, and that the most likely outcome of en-
counters between the captured binary and other cluster
members leaves the central star in a binary.
2. SIMULATIONS
In this study we combine the SPH results of MB07
with n-body simulations of a cluster similar to the ONC.
This is similar in spirit to the work of Scally & Clarke
(2001) and Pfalzner et al. (Pfalzner et al. 2006; Pfalzner
2006; Pfalzner & Olczak 2007), in which n-body simula-
tions of a cluster are used to determine encounter fre-
quency and parameters, which are then used to study
the effect of encounters on disks. This work is differ-
ent in that the results of the close encounters are in-
cluded in the simulation as it is running, similar to
McDonald & Clarke (1995), who used an analytic pre-
scription for disk-mediated binary formation in simula-
tions of a cluster with 10 stars.
The simulations are performed using Aarseth’s code
NBODY6 (Aarseth 2000). The cluster is set up as a
King model (e.g. Binney & Tremaine 1987) with W0 =
9, core radius rc = 0.2 pc, central density n0 = 2.0× 104
pc−3, and central velocity dispersion σ0 = 2.19 km
s−1. These parameters are similar to those of the ONC
(Hillenbrand & Hartmann 1998), and with a cut-off ra-
dius of 2.0 pc yield 4725 stars in the simulation. The
IMF used is that of Kroupa (2002) in the range 0.3 -
9.0 M⊙ . Each star is randomly placed according to the
King model density distribution, with a random veloc-
ity appropriate to its radius. In addition, we place a 20
M⊙ star at the center of the simulation with a random
velocity. 1000 sets of initial conditions were generated;
each is run with and without the effects of disk dissipa-
tion. In the dissipation runs, the only star with a disk
is the central star. This is an artificial situation; while
the most massive disk will dominate an interaction, the
presence of disks around all stars in the cluster would
http://arxiv.org/abs/0704.1162v1
mailto:moeckel@colorado.edu
Fig. 1.— Energy change as a function of remnant disk mass for
three of the binaries simulated in Moeckel & Bally (2006). The
outlying point for i = 180◦ shows that the true effect of disk dis-
ruption on repeated encounters is more complicated than a simple
mass scaling. However, for most encounters, a scaling between the
two bold lines is reasonable.
increase the effect of disk encounters.
MB07 used a modified version of the publicly released
SPH code GADGET-2 (Springel 2005) to simulate en-
counters between the massive star-disk system and im-
pactors with masses in the range 0.3 - 9.0 M⊙, perias-
tra in the range 50 - 550 AU, and inclinations from 0◦
- 180◦. The disk radius is 500 AU. Interpolation of the
data from these simulations gives, for any impactor mass,
periastron, and inclination angle, a change in orbital en-
ergy and change in disk mass associated with the en-
counter. We have modified NBODY6 to detect encoun-
ters between the central star and other cluster members.
At the encounter periastron, we determine the orbital pa-
rameters, and then modify the velocity of both encounter
partners in a momentum conserving fashion so that the
change in orbital energy is equal to that found via in-
terpolation of the data from MB07. The change in disk
mass due to the encounter is also tracked. Because we
only have data for periastra ≥ 50 AU, encounters with
smaller periastra are assumed to occur at 50 AU.
When the orbit of a relatively massive impactor is
coplanar to the disk, MB07 find that accretion and disk
capture can increase the impactor’s mass by up to 10%.
However, the change in orbital energy during these pas-
sages is an order of magnitude greater than can be ac-
counted for by accretion drag; the dominant contribution
to the energy change is the change in velocity due to disk
interactions. We make the simplification that the change
in orbital energy is due entirely to a change in the rel-
ative velocity of the two stars. Under this assumption,
at periastron the velocity kick δv on the impactor for a
given orbital energy change ∆E is given by
−2v +
4v2 +
m(1 +m/M)
v̂, (1)
where v is the pre-kick velocity, m is the impactor mass,
and M is the primary mass. The change in velocity for
Fig. 2.— The distribution of encounter frequency for runs with no
disk dissipation. Also plotted are the best fit Poisson distribution
to the realizations in which the central star remains in the cluster
core (solid), and the expected value from equation 3 (dotted). The
error bars are single sided 1-σ confidence limits (Gehrels 1986)
the primary is δV = −(m/M)δv. The change in the
total energy due to these kicks is tracked and included
in energy checks.
In order to account for the change in disk mass, we
scale the change in energy for an encounter according to
∆E(md, i,mi, rp) = ∆E0(i,mi, rp)
. (2)
Here ∆E is the orbital energy change from an encounter
with disk massmd, impactor massmi, inclination angle i,
and periastron rp. ∆E0 is the energy change for the same
parameters with the disk at its original mass md0, and is
found from interpolation of the data in MB07. The index
n we take to be 1 or 2, and scales the proportionality
of the the energy change with the remaining disk mass
fraction.
The scaling of energy change with disk mass is not
completely straightforward. Heller (1995), simulating en-
counters with disks of different masses, found that the
change in energy scales roughly linearly with the disk
mass, in which case n = 1. Moeckel & Bally (2006) sim-
ulated repeated encounters between a captured impactor
and the same remnant disk. Analysis of that data (Fig.
1) shows that the scaling of energy change and disk mass
is more like n = 2. However, this also takes into ac-
count the changing radius of the disk, an effect which
we do not include in this study. Since the proper scaling
depends on details of the encounters that would require
full SPH simulations of each close passage, we instead
run all simulations with both scalings and compare the
results.
Our scheme preserves the binary periastron separation,
which is shown by Moeckel & Bally (2006) to decrease
during repeated retrograde encounters, and increase with
prograde encounters. In the most extreme (and rarest)
case, an in-plane retrograde passage, the periastron sepa-
ration decreases by ∼ 25% after the first encounter. The
data in MB07 show that the change in energy scales com-
parably to the periastron radius; thus we would expect
an error . 25% due to our artifically fixed periastron.
This is similar in magnitude to the uncertainty in our
mass scaling, which as shown below has a negligible im-
pact on our results.
Each case is run for 0.5 Myr; we limit the analysis
to this time for three reasons. After this time the disks
are mostly destroyed. In our simulations this destruction
is by encounters, while in reality the additional effect of
photo-evaporation will contribute. Mass segregation also
begins to take effect after this time, which is an added
complication in the analysis of the results. Finally, the
multiplicity of stars is established early in their evolution
(Mathieu 1994), and a mechanism for binary formation
should work on short enough timescales to reflect that.
3. RESULTS
3.1. Encounter Rates
The calculation of binary capture rates (Clarke &
Pringle 1991; Heller 1995; Boffin et al. 1998; MB07) is
largely similar to the estimation of collision or encounter
rates. We begin by comparing the standard encounter
timescale calculation to the simulations.
The encounter rate for a star of mass M in a cluster
with number density n, velocity dispersion σ, and en-
counter radius r is given by
γ = 4
πnσr2
. (3)
Here M is M + m̄, with m̄ the average stellar mass in
the cluster (Binney & Tremaine 1987; MB07 for general
masses). For the parameters of our simulated clusters,
the encounter rate is γ = 1.02 × 10−5 yr−1, with 5.1
encounters expected over 0.5 Myr.
Plotted in figure 2 is the distribution of the number
of unique cluster members encountered by the central
star, without disk dissipation. Multiple encounters with
the same star, for instance in a binary, are counted only
once. Because some of the central stars have a high ini-
tial velocity and escape to the cluster outskirts, there is
an excess of low encounter-number runs. Therefore the
distribution for all runs is plotted, as well as only those in
which the central star remained within the cluster core.
One would expect that the distribution for stars in the
same environment would be Poisson; the data shows oth-
erwise. Because the number density and velocity disper-
sion are radially dependent, central stars that move to
larger radii are exposed to a much different environment
than those that remain near the cluster center. By limit-
ing our analysis to the cluster core, where the properties
are closest to those used in equation 3, the best fit (with
mean value 5.78) is in reasonable agreement with the
theoretical value of 5.1.
Equation 3 is averaged over the mass function. Because
of the dependence of the encounter rate on the mass of
the stars, encounters with more massive cluster members
are more frequent. The mass function of encounter part-
ners is well fit by a single power-law mass function of the
form ξ(m) ∝ m−α with α = 0.92, while a Salpeter mass
function has α = 2.35.
3.2. Binary Fraction and Mass Function
Of greatest interest is the fraction of the massive stars
at a given time that are in a multiple system. MB07
Fig. 3.— The number of runs with the central star in a binary,
NB(t), as a fraction of the total number of runs NT for the three
simulation series. Cases with both n = 1 and n = 2 in equation 2
are at ∼ 30% after 0.5 Myr.
calculate that for a cluster with the parameters of our
models, the binary formation rate is Γ = 0.6 Myr−1
(Figure 3 in MB07). We compare this to our simula-
tions as follows. For each simulation, a list of stars that
encounter the central star is generated, with the times of
their encounters. If the same star is involved in consec-
utive encounters, a binary has formed, and we track the
following events.
Ionization: The next two encounters are with differ-
ent stars. Even if the intruder and the original binary
partner form a binary, the central star is no longer in a
multiple system. Exchange: The next two encounters are
with the same star, which is different from the initial bi-
nary. Flyby: The next encounter is with a different star,
but the one after that is with the initial binary partner.
In the latter two cases, the central star is considered to be
in a binary throughout. In this scheme the possibility ex-
ists that random, consecutive encounters with the same
star could contaminate the binary statistics. In practice,
the number of binaries with semi-major axes larger than
the average inter-stellar distance in the cluster core is on
the order of 1%, and we consider the detected binaries
to be true binaries.
Plotted in figure 3 are the number of central stars in
binaries as a function of time, for each of the three series
of simulations. Considering first the series with no disks,
we see an initial rise in the binary fraction, followed by
a leveling off to ∼ 12% at 0.5 Myr, as binary creation
through random dynamical processes is balanced by ion-
ization. For the two series with disk dissipation, there is a
steady rise up to a value of approximately 30%, in good
agreement with the calculations in MB07. The binary
fraction does not appear to depend on the specifics of
the energy-change disk-mass scaling; the early passages,
when the disk still has nearly its original mass, account
for the increased binarity.
Since massive impactors are preferentially captured by
the disk and more likely to remain bound during ex-
changes and flybys, the mass function of binary partners
at 0.5 Myr is flatter than that of the encounter partners.
TABLE 1
Binary fate statistics at 0.5 Myr
Series NB NE NI NF NS (NE +NF )/NI
no disk 387 44 224 431 119 2.12
n = 1 708 186 217 1619 305 8.32
n = 2 698 169 238 1482 291 6.94
With no disks, the best fit single power-law mass func-
tion has α = 0.67. The disk case with n = 1 has α =
0.16, and for n = 2 we have α = 0.12.
3.3. Binary Robustness
MB07 estimated the survival probabilities of capture-
formed binaries, finding that massive companions are
more likely to survive encounters with other cluster mem-
bers, and that lower mass companions are likely to be ion-
ized. The binary fractions found here are in agreement
with calculations that don’t include ionization effects; in
order to explain this we turn to the question of ionization
versus exchange and flybys.
Shown in table 1 are NB, the number of binaries
formed, NE , the number of exchanges, NI , the number
of ionizations, and NF , the number of flybys that occur
in each series of 1000 simulations. NB includes all unique
binary pairings, so that an exchange contributes twice.
Thus the number of surviving binaries NS is given by
NB − NE − NI . The ratio of exchanges and flybys to
ionizations, (NE +NF )/NI , is indicative of the survival
chances of a binary in each series. For the diskless case,
n = 1, and n = 2 this ratio is 2.12, 8.32, and 6.94 respec-
tively. Encounters between binaries and other cluster
members are much more likely to end with the central
star in a binary for the simulations with dissipation com-
pared to the diskless case, and suggest that the binaries
formed via this capture mechanism are more robust than
indicated by the simple estimations of MB07. The in-
creased disk dissipation in the n = 1 series yields slightly
harder binaries, but the total binary fraction is not sig-
nificantly affected by the scaling.
4. DISCUSSION
The simulations presented here are intended to test and
verify the capture rates calculated in MB07. The conclu-
sions of that work are largely upheld; encounters occur
at approximately the expected frequency, and binaries
are captured at a rate consistent with the analytical es-
timates. In addition, once a binary is formed it is less
likely to be destroyed by ionization than the estimates
in MB07 suggest. Encounters between the captured bi-
naries and intruding cluster stars are far more likely to
result in a flyby or exchange than in an ionization, leav-
ing the central star in a binary.
The capture rates for the central star-disk system are
not high enough to fully account for the high multiplicity
observed in massive, cluster-bound stars. Our simula-
tions produce ∼ 30% binarity at 0.5 Myr, the time when
the disks are mostly destroyed by encounters or photo-
evaporation, an effect not modeled here. Recent work
(Krumholz et al. 2007) shows that during the formation
of a massive star, material moves through the protostel-
lar disk in sporadic, massive accretion events, between
which the disk builds up to large masses. A disk that is
∼ 30−50% of the central star’s mass, instead of the 10%
used here, could increase the capture rates by a factor of
several. Additionally, continued accretion onto the sys-
tem or the presence of disks around all the stars could
increase the capture rates. It is worth noting that our
simulations here are tailored to ONC-like systems. Since
the capture rate is linear with stellar density and drops
with higher velocity dispersion, changing the cluster pa-
rameters will affect the results.
The effect of encounters on accretion processes in mas-
sive star formation is unclear. It is possible that the de-
struction of the disk by repeated passages could truncate
accretion at the time of binary formation. Alternatively,
accretion could resume onto the binary and tighten the
orbit, leading to a massive, close binary system. The
frequency of the encounters modeled here is high enough
that such a situation warrants further investigation.
As concluded in MB07, this rate can not be ignored.
However, additional binary formation mechanisms must
be employed to explain the observed multiplicity of mas-
sive stars. In the Trapezium the massive stars have, on
average, 1.5 companions (Zinnecker & Bate 2002). Since
disks are effectively destroyed during binary capture, this
is a process that can only account for a single compan-
ion, unless further accretion creates a new circumstellar
or circumbinary disk.
This work was supported by NASA grant
NNA04CC11A to the CU Center for Astrobiology.
REFERENCES
Aarseth, S. J. 2000, in The Chaotic Universe, Proceedings of the
Second ICRA Network Workshop, Advanced Series in
Astrophysics and Cosmology, vol.10, Edited by V. G.
Gurzadyan and R. Ruffini, World Scientific, 2000, p.286, ed.
V. G. Gurzadyan & R. Ruffini
Binney, J., & Tremaine, S. 1987, Galactic Dynamics (Princeton,
NJ, Princeton University Press, 1987, 747 p.)
Boffin, H. M. J., Watkins, S. J., Bhattal, A. S., Francis, N., &
Whitworth, A. P. 1998, MNRAS, 300, 1189
Cesaroni, R., Galli, D., Lodato, G., Walmsley, C. M., & Zhang,
Q. 2007, in Protostars and Planets V, ed. B. Reipurth,
D. Jewitt, & K. Keil, 197–212
Clarke, C. J., & Pringle, J. E. 1991, MNRAS, 249, 584
Garćıa, B., & Mermilliod, J. C. 2001, A&A, 368, 122
Gehrels, N. 1986, ApJ, 303, 336
Heller, C. H. 1995, ApJ, 455, 252
Hillenbrand, L. A., & Hartmann, L. W. 1998, ApJ, 492, 540
Kratter, K. M., & Matzner, C. D. 2006, MNRAS, 373, 1563
Kroupa, P. 2002, Science, 295, 82
Krumholz, M. R., Klein, R. I., & McKee, C. F. 2007, ApJ, 656,
Mason, B. D., Gies, D. R., Hartkopf, W. I., Bagnuolo, W. G.,
Brummelaar, T. T., & McAlister, H. A. 1998, AJ, 115, 821
Mathieu, R. D. 1994, ARA&A, 32, 465
McDonald, J. M., & Clarke, C. J. 1995, MNRAS, 275, 671
Moeckel, N., & Bally, J. 2006, ApJ, 653, 437
—. 2007, ApJ, 656, 275
Pfalzner, S. 2006, ApJ, 652, L129
Pfalzner, S., & Olczak, C. 2007, A&A, 462, 193
Pfalzner, S., Olczak, C., & Eckart, A. 2006, A&A, 454, 811
Preibisch, T., Balega, Y., Hofmann, K.-H., Weigelt, G., &
Zinnecker, H. 1999, New Astronomy, 4, 531
Scally, A., & Clarke, C. 2001, MNRAS, 325, 449
Springel, V. 2005, MNRAS, 364, 1105
Stahler, S. W., Palla, F., & Ho, P. T. P. 2000, in Protostars and
Planets IV, ed. V. Mannings, A. Boss, & S. Russell (Tucson:
Univ. Arizona Press), 327
Zinnecker, H., & Bate, M. R. 2002, in ASP Conf. Ser. 267: Hot
Star Workshop III: The Earliest Phases of Massive Star Birth,
ed. P. Crowther, 209
