The dynamical evolution of massive black hole binaries - I. Hardening in a fixed stellar background

16 Pages·0.5 MB·English

by  Gerald D. Quinlan

The dynamical evolution of massive black hole binaries -
I. Hardening in a fixed stellar background.
Gerald D. Quinlan
LickObservatory, Universityof California, Santa Cruz CA 95060
Dept. of Physicsand Astronomy, Rutgers University,PO Box 849, Piscataway NJ 08855∗
∗Present address
6
9
17January1996
9
1

n ABSTRACT
a
The stellar ejection rate and the rates of change of the binary semimajor axis and
J
eccentricity are derived from scattering experiments for the restricted three-body
7
problem.They areusedto study the evolutionofbinariesinsimple models for galac-
1
tic nuclei, starting soon after the black holes become bound and continuing until the

evolutionisdominatedbythe emissionofgravitationalradiation,oruntil the ejected
1
mass is too large for the galaxy to be considered fixed. The eccentricity growth is
v
2 found to be unimportant unless the binary forms with a large eccentricity. The scat-
9 teringresultsarecomparedwithpredictionsfromChandrasekhar’sdynamical-friction
0 formulaandwithpreviousworkonthecaptureandscatteringofcometsbyplanetary
1 systems. They suggest that a binary with masses m1 ≥m2 should not be considered
0 hard until its orbital velocity exceeds the background velocity dispersion by a factor
6 that scales as (1+m1/m2)1/2.
9
/ 1 INTRODUCTION galacticnuclei,firstproposedbyBBR.Andtheywouldtell
h
uswhat toexpect ifthreeor moreBHs enterthecore of a
p
- The existence of massive black hole (BH) binaries follows galaxy, which can happen if the BHs are dragged in from
o fromtwowidely-acceptedassumptions:thatgalaxiesmerge thegalaxy’shaloorifthegalaxyundergoesmultiplemerg-
r with other galaxies, and that many galaxies contain mas- ers with other galaxies containing BHs. If the first binary
t
s siveBHs.ForiftwoBHsenterthecoreofamergedgalaxy, mergesfastitcanformabinarywithathirdBH,andonce
a
dynamical friction drags them to the center where they that merges it canform abinarywith afourth,andso on,
:
v form abinary.Thesubsequentevolutionwas firstoutlined leadingtoamassivecentralBH;butifthefirstbinarystill
i byBegelman, Blandford, and Rees(1980, hereafter BBR). exists when a third BH enters then one or all three of the
X
Initiallythebinaryhardens(i.e.itsseparationshrinks)be- BHs can be ejected in a sling-shot interaction. Arguments
r
a cause of the interaction between theBHs and all thestars like these can set limits on massive BHs as dark-matter
inthegalaxycore.ButthatisineffectiveoncetheBHsbe- candidates for galactic halos (see Hut and Rees 1992, Xu
come close because distant stars perturb the binary’s cen- and Ostriker1994 for conflicting limits for our Galaxy).
ter of mass but not its semimajor axis. The binary then Anotherquestion is what a binarymerger does tothe
hardens by giving kinetic energy to stars that pass in its surroundinggalaxy,i.e.whatobservablesignatureitleaves.
immediate vicinity; a hard binary can eject stars out of Massejectionduringtheevolutionshouldreduceagalaxy’s
the core at high velocity. If there are enough stars for the central density and expand its core (BBR). Ebisuzaki,
hardeningtocontinue(andgasaccretionontotheBHscan Makino,andOkumura(1991)haveproposedthisasanex-
help), eventually the BHs merge through the emission of planationforwhylargeellipticalgalaxieshavelowercentral
gravitational radiation; otherwise thehardeningstalls and densities and weaker density cusps than small ellipticals
thebinary survives for thelifetime of thegalaxy. (e.g. Kormendy et al. 1994).
Whether a binary merges or survives and how long it Wearefarfromhavingpreciseanswerstoanyofthese
spends in each stage of the evolution are questions rele- questions.BBRgavearangeofmergertimesforonetypical
vant to a number of problems in extragalactic astronomy. examplethatspannedthreeordersofmagnitudebecauseof
Their answers would help us predict the total BH merger theuncertaininfluencethatmassejectionhasonthehard-
rate and whether it is high enough for us to detect the ening rate. Fukushige, Ebisuzaki and Makino (1992) have
resulting gravitational waves (e.g. Thorne 1992, Haehnelt argued that dynamical friction causes a binary to become
1994).TheywouldhelpusassessBH-binarymodelsforthe highly eccentric and that this greatly reduces the merger
bending and apparent precession of radio jets from active time because gravitational radiation then becomes impor-
2 Gerald D. Quinlan
tant early in the evolution. Although their arguments are are used tostudycometary dynamics(see Fernandez1993
not convincing, they have called attention to the eccen- for a review) they are not of much help for our questions
tricity growth and ourignorance of its correct description. aboutBHbinaries,partlybecausetheyoftenconsideronly
Thereareuncertaintiesinhowthehardeningratedepends one mass ratio (for the Sun-Jupiter system), but mostly
ontheratioofthetwoBHmasses,inwhenabinarymakes because they are used to answer different questions, such
thetransitionfromsofttohard,andeveninwhatthewords as the cross section for the capture of interstellar comets,
softandhardshouldmeaninthiscontext.Andourknowl- orfortheconversionoflong-periodcometstoshort-period
edge of how a binary merger changes a galaxy is based on comets, the survival probability of comets once they are
back-of-the-envelope estimates and simple N-body experi- captured,andhowallthesedependonthecomet’sinclina-
mentswith unrealistic galaxy models. tion. There is nevertheless some overlap between the two
The first step towards resolving these questions is to problems.
understand how a massive binary evolves in fixed stellar Thegoal ofthispaperistopresentaccuratemeasure-
background. Consider a binary with masses m1 m2 and mentsofH,J,andKovertherangeofparametersofinter-
≥
semimajor axis a in an isotropic background of stars of est for the BH-binary problem, including the dependence
mass m∗ m2, density ρ, and one-dimensional velocity on the mass ratio, eccentricity, and degree of hardness of
≪
dispersion σ; let M12 and µ denote the total and reduced thebinary.Otherquantitiestobestudiedincludethecross
binary mass: section for a binary to capture stars into bound orbits,
for close encounters between stars and the binary mem-
M12 =m1+m2, µ=m1m2/M12. (1) bers, and the distribution of velocities with which stars
The binary evolution and its effect on the galaxy are de- are expelled from the binary. These add to our under-
termined by three dimensionless quantities: the hardening standing of H, J, and K, and are needed by themselves
rate for some applications. The results will be presented in a
model-independentwaysothattheycanbeappliedtoany
σ d 1
H = Gρdt a , (2) problem with m∗ ≪m2 ≤m1.
Once H, J, and K have been measured they can be
(cid:16) (cid:17)
the mass ejection rate (where Mej is the stellar mass that usedtostudytheevolutionofbinariesinfixedgalaxymod-
thebinary has ejected from thegalaxy core) els. That will be done here for some simple models. If the
BHsarelargetheywillofcourseejecttoomuchmassfrom
J = 1 dMej , (3) thegalaxycoreforittobeconsideredfixed.Buttheresults
M12 dln(1/a) will still be valid during the early stages of the evolution.
and theeccentricity growth rate And they will be helpful even in the later stages, because
wecanimagineatanyinstantthatthebinaryisembedded
de
K = . (4) in a fixed background whose properties are those of the
dln(1/a)
galaxy at that instant. The self-consistent evolution of a
The quantities H, J, and K can be found from scattering massivebinaryinarealisticgalaxymodelandthechanges
experimentsthattreatthestar-binaryencountersoneat a thisinducesinthemodelarebeststudiedbylargeN-body
time.Analyticapproximationssuchastheimpulseapprox- experiments. That will bedeferred topaper II,along with
imationarehelpfulduringtheearlystagesoftheevolution a discussion of what both papers imply for theastronomi-
(Gould 1991), but not once thebinary becomes hard. cal questionsmentioned above(Quinlan and Hernquist,in
Most publishedscatteringexperimentsassumethebi- preparation).
naries and stars to have equal or nearly equal masses (see
Heggie1988forareview).Therearesomeexceptions.Roos
(1981) performed scattering experimentsfor therestricted
three-bodyproblemtostudytheevolutionofhard,massive 2 COMPUTATIONAL METHOD
BHbinaries,andusedthemtocorrectamisplacedfactorof
2.1 Derivation of results from the restricted
m1/m2 intheBBRhardeningrate.HetriedtomeasureK three-body problem
buthisstatisticsweretoopoortogivedefiniteconclusions
(only 500 orbits per measurement). Hills (1983a) used the Wetreat thestar asa massless test particle movingin the
general three-body problem to study interactions between potential of the two BHs. From the changes in the star’s
a massive binary and low-mass intruders(m∗/m2 =0.01). energy and angular momentum per unit mass, ∆E∗ and
He gave results for the hardening rate for a wide range of ∆L∗,weinferthecorrespondingchanges∆E and∆Lthat
massratios(m1/m2 =1–300),butlikeRoosheconsidered the binary would have suffered if the star had been given
onlyveryhardbinaries.MikkolaandValtonen(1992)used a small but nonzero mass. The threebodies are treated as
therestrictedthree-bodyproblemtomeasureH andK for point masses and gravitational radiation is ignored.
equal-mass BH binaries with varying degrees of hardness. Inarealgalaxystarsapproachthebinarywithawide
Theirmeasurementsareaccurateforhardbinariesbuthave distribution of velocities at any given time. But the scat-
large error bars for binaries that are not hard. tering experiments are easiest to perform if the stars all
Ifm1 m2 thentheinteraction betweenastaranda start from the same velocity v at a large separation from
≫
BHbinaryissimilartotheinteractionbetweenacometand the binary (the initial velocity, or the velocity at infinity).
a planet orbiting a star. Although scattering experiments In that case we write
Massive Black Hole Binaries 3
d 1 Gρ 2.2 The scattering experiments
= H1, (5)
dt a v
(cid:16) (cid:17) Each scattering experiment requires five uniformly-
where the “1” indicates that the stars all have a single
distributedrandomnumbers(fourifthebinaryiscircular):
velocity v. The hardening rate for a Maxwellian velocity
oneforthesquareoftheimpact parameter(insomerange
distribution is then [0,b2 ]),andfourtofixthebinary’sorientationandphase:
max
∞ σ the cosine of the inclination ([-1,1]), the longitude of as-
H(σ)= dv 4πv2f(v,σ) v H1(v), (6) cending node ([0,2π]), the argument of pericenter ([0,2π]),
Z0 and the mean anomaly at some fixed time ([0,2π]). The
where numbersare chosen with thequasi-random numbergener-
1 ator sobseq of Press et al. (1992).
f(v,σ)= exp( v2/2σ2). (7)
(2πσ2)3/2 − The range [0,bmax] for impact parameters is split into
five intervals corresponding to ranges in scaled pericenter
Thehardeningrateisderivedfromtheaverageenergy
distancer /aof[0,1],[1,2],[2,4],[4,8],and[8,16].Eachout-
p
change for stars that scatter off the binary. We define a
putquantityismeasuredinanumberofsteps.Onthefirst
dimensionless energy change C by (Hills 1983a)
step the program spends short but equal amounts of cpu
C = M12 ∆E = a∆E∗. (8) time picking orbits from the five intervals. On each suc-
2m∗ E Gµ cessivestep theprogram doublesthecputimeandadjusts
its strategy so that the time it spends on each interval is
Thismustbeaveragedoverallangularvariablesdescribing
proportionaltotheuncertaintythatintervalcontributesto
the binary’s orientation and phase, to give C , and then
h i the quantity being measured. Once the uncertainty is re-
integrated over all impact parameters. The averaging and
ducedtoanacceptablelevel,orthecputimeexceedssome
integrating are done in a Monte Carlo fashion by picking
maximum allowed value, the results from all five intervals
orbits from suitable distributions. We sometimes describe
are combined with appropriate weights for a distribution
anorbitbyitsimpactparameterb,thedistanceatwhichit uniform in b2. For H, J, and K the last three intervals
wouldpassthebinaryifitfeltnoattraction,andsometimes
contributelittlebecausethechangesinenergyandangular
by its pericenter distance r , thedistance if it is attracted
p momentumfalloffrapidlywithincreasingimpact parame-
byapoint massM12.Thetwoarerelated bygravitational ter.
focusing:
The coordinates are chosen so that the binary’s cen-
ter of mass is at the origin and the star starts at infinity
b2 =rp2 1+ 2GrpMv212 . (9) with (x,y,z) = (b,0,∞) and (vx,vy,vz) = (0,0,−v). The
(cid:18) (cid:19) star is moved from r = to r = 50a along a Keplerian
∞
Wealso definea dimensionless impact parameter x by orbitaboutapointmassM12 attheorigin.Thenumerical
integration starts at r=50a.
x=b/b0, b20 =2GM12a/v2; (10) The orbits are integrated in double precision with an
b0 is approximately the impact parameter corresponding explicit,embeddedRunge-Kuttamethodoforder(7)8:the
to r =a if gravitational focusing is important. With this program dopri8 of Hairer, Norsett, and Wanner (1987).
p
notation we can write The program adjusts the integration stepsize to keep the
fractional errorperstepintheposition andvelocitybelow
∞
H1=8πIx(C)=8π dx x C , (11) some level ǫ, which was set to 10−9. With this choice the
Z0 h i change in a star’s Jacobi constant for a circular binary is
where thesecond equality definesthe operator I . at most 10−6GM12/a and often much smaller. The forces
x
from theBHs are not softened.
The derivation of the eccentricity growth rate is sim-
Someintegrationsaretimeconsumingbecausethestar
ilar. The change to the binary’s eccentricity from a single
gets captured into a weakly-bound orbit and makes many
scattering event is, if the change is small and the binary’s
revolutionsbeforeitisexpelled.Theintegrationstepsizeis
orbital angular momentum points in the z direction,
a small fraction of the binary’s period even if the period
∆e = (1−e2) 2 ∆Lz E 1 . (12) of thestar about thebinary is much longer. Thefollowing
∆ln(1/a) 2e − Lz ∆E − approximation expedites those experiments. If a captured
h (cid:16) (cid:17)(cid:16) (cid:17) i star moves further than rk = a(1010m2/m1)1/4 from the
Wedefineadimensionlessangular-momentumchangeBby
binary,thebinaryisreplacedbyapointmassM12 andthe
B= M12∆Lz = M12 ∆L∗,z , (13) starismovedalongaKeplerianellipseuntilitreturnsinside
− m∗ Lz µ [GM12a(1−e2)]1/2 rthk,ewcohrernecttheorfboirtcaelspfhroamse.tTheheBH(ms2a/rme1r)e1i/n4trmodauscsedscawliinthg
and can then write makes the quadrupole force from the binary at r about
k
(1 e2) I (B C) 10−10GM12/a2, independentof m2/m1.
K1 = − x − , (14) Orbits that get captured for long times tend to be
2e I (C)
x
(cid:20) (cid:21) highly chaotic. The integration for any particular orbit of
wherethe“1”hasthesamemeaningasbefore.Thederiva- that type is difficult to justify because a small change to
tion of K from K1 will bedescribed later. theintegrationprocedurecanmakeabigchangetotheout-
4 Gerald D. Quinlan
come.Buttheaverageresultsderivedfromalargenumber
m1/m2 H0 w/Vbin
of integrations can be correct even if the individual inte-
grations are not; that is suggested by shadowing lemmas 1 17.97 0.5675
that havebeen provedfor simple chaotic systems.The av- 4 20.54 0.4263
erage results presented here do not change noticeably if ǫ 16 21.87 0.2228
64 22.78 0.1043
is raised or lowered by a factor of 100, even though some
256 22.57 0.0573
orbits undergobig changes.
An integration is stopped when the star leaves the
sphere r = 50a with positive energy. The average results Table1.ParametersforfitstoH1(eq.16)foracircularbinary.
arenotsensitivetothelocationofthissphereprovidedthat
it is at least 10–15 times larger than a. Once the integra-
tion stops the program records the changes to the star’s
3 RESULTS FROM THE SCATTERING
energy and angular momentum, the minimum separations
EXPERIMENTS
between the star and the two BHs, and between the star
and the binary’s center of mass, the integration time, the 3.1 Hardening
number of integration steps, and the number of times the
3.1.1 Hardening rate
star’s radius passed through a minimum (which, if greater
than one, gives the number of revolutions that a captured ThehardeningrateH1(5)hasbeenmeasuredasafunction
orbitmade).Anorbitthatdoesnotget capturedtypically of the binary’s eccentricity and hardness for a wide range
takes a few hundred to a few thousand integration steps. ofmassratios.ItisplottedinFigure1versusthehardness
Captured orbits can take much longer. If an integration asgivenbytheratiooftheinitialstellarvelocityvandthe
lastsformorethan106 stepsitisabandoned.Thefraction binary’s orbital velocity V (the relative velocity of the
bin
of abandoned integrations is the largest for hard binaries two BHsif thebinary is circular):
with m1/m2 1, but even for those it is less than 0.1%.
≫
Vbin = GM12/a. (15)
The error in any average quantity has a systematic
component and a statistical component. Systematic errors The errpor bars show the statistical-error estimates; if not
arise, for example, because of errors in the numerical in- visible then theyare smaller than thesize of the points.
tegration, because integrations are abandoned if theytake ThevelocitydependenceofH1isfitbyafunctionwith
toolong,becausetheprogramimposesamaximumimpact two free parameters (whose values are given in Table 1):
parameter bmax, and because the orbits start and end at
tro=tal50syasitnesmteaatdicoefrraotrris=u∞sua.lBlyutmnuocnhesomfatlhleersethisanlatrhgee:stthae- H1 = [1+(vH/w0)4]1/2. (16)
tisticalerror.Thestatisticalerrorsareestimatedbytaking
thedifferencebetweenresultsfoundwithN orbits(thefinal Ttohedefcurnecatsieonashvasaappcroonascthaenstwv,alaunedHd0ecarteavses≪asw1,/svt2aratst
number)andN/2orbits,orsometimes—ifthatdifference
v w. This fits the data well at high and low veloc-
looks suspiciously small — onehalf thedifferencebetween ≫
N orbits and N/4 orbits. That gives a rough estimate of ities. It does not fit so well when H1 first starts to de-
crease;thosedatapointsweregivenlittleweightinthefit-
theerrorlevel.ThestatisticalerrorsdecreaseatlargeN as
(lnN)d/N when quasi-random numbers are used, where d tingprocedure.Athree-parameterfunction wasalso tried,
isthenumberofnumberspickedforeachexperiment(5in H1 = H0/[1+(v/w)κ]2/κ, but the exponent κ never dif-
fered bymuch from 4.
general, 4ifthebinaryis circular).That is fasterthan the
N−1/2 decrease that occurs with random numbers (Press Foran equal-mass binary theconstant hardeningrate
et al. 1992). at low velocity, H0 ≃ 18.0, agrees well with Mikkola
and Valtonen’s (1992) πR 18.2, and also with the
a
≃
The number of orbits needed to reduce the statistical results of Hills (1983a, 1992). At high velocity H1
errors to an acceptable level varies widely with the binary H0(w/v)2 5.8(Vbin/v)2, which can be compared wit≈h
≃
andthemeasurement.Measurementsaremoredifficultfor Gould’s (1991) analysis using the impulse approximation,
K1 than for H1 because of cancellation. Cancellation is a H1,Gould =(8π/3)(Vbin/v)2. The agreement is satisfactory
problemforH1tooathighvvalues.Andregardlessofwhat considering that the error bars are large at high veloc-
is being measured, binaries with m1 m2 require many ity and that the impulse approximation is justified only
≫
more orbits than equal-mass binaries because of the rare if v/Vbin 1. In fact it is surprising how well theimpulse
≫
close encounters with the mass m2. When the results are approximation works when v/Vbin 1.
≃
presented the number of orbits used will not be given be- Panel (a) shows that the hardening rate for an equal-
causeitisdifferentforeachmeasurement;someestimateof mass binary does not vary much with the eccentricity.
thestatistical errorwillbegiveninstead.Thenumbersare MikkolaandValtonen(1992)reachedthesameconclusion,
about 104–105 orbits per measurement for H1 (sometimes which remains true for all mass ratios. For later applica-
more at high v values) and 105–106 for K1. The complete tions the hardening rate for a circular binary will be used
set of experiments took about four monthsof cpu time on foralleccentricitiesbecausethevariationofH0andwwith
an IBM 580 RISCworkstation. eccentricity is too small to matter.
Massive Black Hole Binaries 5
Figure 1.Binaryhardening rateversus theinitial stellarvelocity v: (a) shows threeeccentricities form1/m2 =1;(b) and (c) show
fivemassratiosfore=0,onlinearandlogarithmicscales.Thelinesarethefitstoeq.(16).
Roos (1981) and Hills (1983a) showed that the low- For an equal-mass binary the coefficients of the log terms
velocity hardening rate H0 does not vary much with the differ by about 30%, which is satisfactory considering
massratio. Butthevelocitywdoes.Thevariationisfitby the uncertainties mentioned above. For non-equal masses
Gould’s log term does not have the correct mass depen-
w 0.85 m2 V =0.85 Gm2. (17) dence (Gould did not attempt to compute the log term
≃ rM12 bin r a accurately).
The hardening rate for a massive BH binary is some-
The physical significance of this mass dependence is the
timesderivedfromChandrasekhar’sdynamical-frictionfor-
following: if v<w thebinary can easily capturestars into
mula(e.g.Fukushigeetal.1992).Theerrorinthathasbeen
boundorbits; if v>w it cannot.
known for many years (Chandrasekhar 1944, Hills 1983a):
Theintegral(6)foraMaxwelliandistributionwaseval-
the distant encounters included in the friction formula do
uated numerically. The relation between H and H1 is fit not perturb the binary’s semimajor axis — they only per-
closely by theformula
turb its center of mass. It is an accident that the deriva-
tion gives a result like Gould’s for a Maxwellian distribu-
H(σ) 2 σ β
+ln 1+α , (18) tion if a suitable choice is made for the log term: if the
H1(√3σ) ≃rπ (cid:20) (cid:16)w(cid:17) (cid:21) same derivation is used for H1 it gives thenonsense result
that (for m1 = m2) H1 is zero at v = 0, rises as H1 v
with α = 1.16 and β = 2.40. In the limit of high velocity ∼
for v < V /2, and then drops abruptly back to zero at
bin
thisgives
v=V /2(becauseonlystarsmovingslowerthantheBHs
bin
H βH0 w2 ln σ2 . (19) cfoorntfurirbtuhteeridnisCcuhsasniodnr.asekhar’s formula). See Gould (1991)
≈ 6 σ2 w2
(cid:18) (cid:19) (cid:18) (cid:19) The velocity dependence of the hardening rate sug-
The log term looks like a familiar Coulomb logarithm but gests a new convention for the use of the word hard. A
comes from an integral over the velocity distribution, not hard binary is usually defined in one of three ways. The
over a range of impact parameters. The limit (19) can be first says that a binary with binding energy Eb is hard if
compared with Gould’s (1991) hardening rate, Eb m∗σ2 and soft if Eb m∗σ2 (p. 534 of Binney
≫ ≪
and Tremaine 1987). The second says a binary is hard if
16√2π Gµ σ2 it grows harder through interactions with stars and soft if
H = ln . (20)
Gould 3 aσ2 V2 it grows softer (Hut 1983). And the third, which is often
(cid:16) (cid:17) (cid:18) bin(cid:19)
6 Gerald D. Quinlan
stated as a corollary of the first or second rather than as
m1/m2 c1 c2 c3 c4 c5
an independent definition, says a hard binary is one that
“hardens at a constant rate,” i.e. at a rate H that is inde- 1 17.97 1.0066 3.5745 2.0865 0.6100
pendent of the hardness. The equivalence of the first two 4 20.54 0.7929 4.5326 1.2675 1.2377
definitionsiscalled Heggie’slaw.Butneitherofthosedefi- 16 21.87 0.4122 3.6588 1.2324 0.9754
64 22.78 0.1800 6.1855 0.5562 1.0087
nitionsis useful for amassive BH binary because both are
256 22.57 0.0846 8.1992 0.3856 0.9782
satisfied by almost any pair of massive BHs that is close
enough to be called a binary.
The third definition gives almost the same result as Table 2. Parameters for fits to Σcap (eq. 22) for a circular bi-
the first two when the masses are equal (see Hut 1983, or nary.
Figure6.3ofSpitzer1987) andisfarmoreusefulwhenthe
binaryismassive.ABHbinarywillthereforebecalledhard
literature). They say their cross sections are accurate to
ifithardensataconstantrate,i.e.ifσ worequivalently
ifVbin/σ (1+m1/m2)1/2.Itistemp≪tingtocallabinary within a factor of 2–3, although that is not clear because
≫ theyadjusttheirformulasinanadhocway—usingHills’s
soft if it is not hard, but that is confusing for massive bi-
(1983a)resultstoguidethem—forhardbinariesforwhich
nariesbecausethereisawidegapforthembetweenahard
the impulse approximation does not work.
binaryinthesenseusedhereandasoftbinaryinthefamil-
The measurements made here improve upon those of
iar sense that “soft binaries grow softer.” In later figures
Hills byusingareproduciblecapturedefinitionandbyex-
the properties of hard binaries are studied with scattering
ploringthedependenceonthebinary’sdegreeofhardness.
experiments using the lowest initial velocity v in Figure 1
Panels (a) and (b) of Figure 2 show the capture cross sec-
foreachmassratio;thosevelocitiesarelog (v)= 2.025,
10 − tion for a circular binary in units of the binary’s geomet-
2.25, 2.6, 2.8,and 3.0form1/m2 =1,4,16,64,and
− − − − rical cross section Σ , which includes the correction for
256. bin
gravitational focusing:
The scattering results and Gould’s (1991) analysis re-
futeHills’s (1990) statement that a binary grows harderif Σ =πa2 1+ 2GM12 . (21)
V > σ and softer if V < σ regardless of the values of bin av2
bin bin
m∗,m1,andm2.Althoughthemeanenergychange C at The velocity(cid:16)dependence(cid:17)is fit by thefunction
h i
zero impact parameter does change from positive to neg-
ative when the stellar velocity is raised from v < Vbin to Σcap =c1 1+ v c3 −c4ln 1+ c2Vbin c5 ; (22)
v>V ,thatsignchangedisappearswhen C isaveraged Σbin c2Vbin v
bin h i h (cid:16) (cid:17) i h (cid:16) (cid:17) i
overimpact parameter. thefiveparametersarelistedinTable2(thefitsshouldnot
The reason for the hard/not-hard transition at σ=w beextrapolatedtovelocitiesmuchhigherthanshowninthe
is best explained after we have examined the cross sec- figure). At low velocity Σcap/Σbin rises as ln(1/v) because
tion for stars to be captured by a binary, to have close the energy change C decreases exponentially with impact
encounters with the binary members, and the distribution parameter.AthighvelocityΣcap/Σbindecreasesasapower
of velocities with which stars are expelled from a binary. of v, which is clearer for the binaries with m1/m2 1.
≫
Thevelocityatthetransitionbetweenthelogarithmicand
power-law behavior,approximately c2Vbin,dependsonthe
mass ratio in the same way as w (eq.17).
3.1.2 Capture cross section
For most velocities and mass ratios the cross sections
We say that a binary captures an incoming star if the inFigure2agreetowithinafactoroftwowiththeapprox-
star’s orbital radius passes through more than one mini- imate crosssections ofPineault andDuquet(1993). There
mum.Almostallcapturedorbitsareeventuallyexpelledin are some larger differences at high velocity for an equal-
the three-body problem (there might be a set of measure mass binary, and at v/Vbin 0.01–0.1 when m1/m2 1.
zero that remain bound forever), but the star can survive The(v/Vbin)−4behaviorseen≃whenm1/m2 1alsoag≫rees
≫
for many revolutions before that happens. with that found by Pineault and Duquet.
Previousworkhasusedscatteringexperimentsandap- The capturecross section rises with the binary eccen-
proximate methods to derive capture cross sections. Hills tricity,butthedependenceis weak.Thedifferencein Σcap
has used scattering experiments to study the capture of for circular and highly-eccentric binaries is only 20–30%,
orbits by very hard, massive binaries (Hills 1983a, 1983b, too small to matter for most applications.
1992).Heunfortunatelydefinescapture—orwhathecalls Panels (c) and (d) of Figure 2 show the cross section
long-termcapture—inawaythatdependsonhisprogram for a captured orbit to survive for at least N revolutions,
(hesaysalong-termcaptureoccursiftheintegrationtakes i.e. for the radius to pass through at least N +1 minima
more than 150,000 steps). But he gives helpful informa- before the star is expelled. The results for binaries with
tiononhowthecaptureprobabilitydependsontheimpact m1/m2 1 are fit well by Σcap(N)/Σcap N−1/2. Ev-
parameter, eccentricity, and binary mass ratio. Pineault erhart (1≫976) noticed this N−1/2 scaling in∼the survival
and Duquet (1993) have used the impulse approximation of comets scattered by the Sun-Jupiter system, and inter-
to derive approximate capture cross sections for massive, preted it as resulting from a random-walk in the comet’s
circular binaries, for arbitrary mass ratios and degrees of energy, as in the gambler’s ruin problem from probability
hardness(they give many relevant references tothe comet theory(seeYabushita1979,Quinn,Tremaine,andDuncan
Massive Black Hole Binaries 7
Figure 2. The cross-section Σcap fora circular binaryto capture an orbit withinitial velocity v, plotted for five mass ratios on (a)
linear and (b) logarithmic scales. The solid lines in (a) and (b) are the fits to eq. (22); the dashed line varies as (v/Vbin)−4. Panels
(c) and (d) show the cross section foran orbitto be captured for at leastN revolutions for the samefive mass ratios, plotted in(c)
atthelowestv foreachmassratio,andin(d)atv=w.Thedottedlinesin(c)and(d)varyasN−1/2.
1990 for further discussion). The N−1/2 scaling does not were captured; the second allowed as many encounters as
work as well if thebinary is not hard or if m1 m2. necessary for thestars to beexpelled.
≃
The cross sections in Figure 2 place no limit on the ThecrosssectionforthelargerBHscalesasΣ/Σ
bin
∼
apocenter of the captured orbit. Some of the stars con- r for the single-encounter experiments because of gravi-
tributing to Σcap are captured into weakly-bound orbits tational focusing. For the multiple-encounter experiments
withapocentersmanyordersofmagnitudelargerthanthe Σ is larger but the increase is important only for r/a <
binary’ssemimajoraxis.Inarealgalaxythoseorbitswillbe m2/m1. The reason the increase is unimportant for r/a>
perturbedbypassingstarsandthegalacticpotentialbefore m2/m1 is that when a captured star orbits a binary with
they return to the binary. But that should not change the m1 m2 thestar’s distance of closest approach tom1 re-
≫
hardeningratemuchbecausethecontributionfromweakly- mainsnearly constantwhileitsenergyundergoes(approx-
boundcaptures is small, even when m1 m2. imately) a random walk. This is well known in cometary
≫
dynamics, where comets diffuse in energy at nearly con-
stantperihelion(e.g.Duncan,Quinn,andTremaine1987).
3.1.3 Close-encounter cross section
Thecrosssectionforacloseencounterwiththesmaller
Thecrosssectionforcloseencounterswiththebinarymem- BHisdifferent.Forthesingle-encounterexperimentsgrav-
bers is needed for applications to real problems where the itational focusing is important for r/a < m2/m1 but not
bodies are not point masses. For a massive BH binary we forr/a>m2/m1,soΣ/Σbin scalesasr orasr2 depending
need it to compute the rate at which stars are tidally dis- on whether r/a is smaller or larger than m2/m1. For the
ruptedbytheBHs, andto estimate howthose disruptions multiple-encounter experiments the cross section is larger
might change the hardeningrate. for all valuesof r,not just for r/a<m2/m1.
Figure 3(a) shows the cross section Σ for a star to The distance r/a = m2/m1 has a special importance
approach within a distance r of either of the BHs, for a for m2 if m1 m2 because the velocity of a star orbit-
≤ ≫
hard, circular binary with m1/m2 =64. The cross section ing m2 at that distanceequalsVbin.Figure3(b) showsthe
isplotted for twosetsof experiments:in thefirst thestars cross section Σ for such encounters as a function of the
wereallowedtoencounterthebinaryonlyonce,evenifthey massratio.Forthesingle-encounterexperimentsΣ/Σ
bin
∼
8 Gerald D. Quinlan
Figure 3. Close-encounter cross sections fora hard, circularbinary (v equals the lowest value inFig.1 foreach mass ratio). Panel
(a)showsthecrosssectionforanorbittoapproachwithin≤rofm1 orm2.Thelinesboundingtheuppershadedregionareform1:
thelowerlineresultswhentheorbitsencounter thebinaryonlyonce;theupperwhentheyhaveasmanyencounters asnecessaryfor
themtoleavewithpositiveenergy.Thelowershadedregionhasthesamemeaning,butform2.Thedottedlinesvaryasr−1andr−2;
thedashedlineisatr/a=m2/m1.Panel(b)showsthecrosssectionforanorbittoapproach withinr≤(m2/m1)aofm2,forboth
singleencounters (open circles)and multipleencounters (filledcircles).Panels (c) and(d) show thedifferential hardeningrates with
respecttothedistances ofclosestapproachtom1 andm2;thefivelinesareform1/m2=1(dotted), 4,16,64,and256(dashed).
(m2/m1)2; for the multiple-encounter experiments cap- at all, which would reduce the cross sections for both m1
turesraise that scaling almost toΣ/Σbin m2/m1. and m2. The two complications tend to cancel for m1.
∼
Panels (c) and (d) of Figure 3 show the differential
hardeningrateswithrespecttothedistancesofclosest ap-
3.1.4 Distribution of final velocities
proach to m1 and m2, normalized so that the area under
the curves is unity. The largest contribution to the hard- Thefinalvelocityisthevelocityofastaratinfinityafterit
ening comes from orbits that pass both BHs at a distance hasbeenexpelledbythebinary.Weneedtheirdistribution
notmuchsmallerthanthesemimajoraxis.Whenm1 m2 to compute themass ejection rate.
≫
thereisawidetailintheleftofpanel(d),butthecontribu- Everhart’s (1968, 1969) work on the scattering of
tion from close encounterswith m2 is still asmall fraction cometsbyplanetarysystemsisrelevanttothedistributions
of thetotal hardeningrate. tobeconsideredhere.Everhartusedanapproximateconic-
matching procedure to derive the probability h(U)dU for
In a real galaxy there will be two complications that the energy change U =∆E∗ to lie in the intervaldU after
can change these results. If weakly-bound captured stars a single encounter between the comet and the planet. The
areperturbedbynearbystarsorthegalacticpotentialthey distribution has three parts, which Everhart called A, B,
will not return to the binary in such a way as to keep and C. Parts A and B are for the small and intermediate
their distance of closest approach to m1 nearly constant. energy changesandarefitwell by(A)aGaussian and(B)
Thatwouldincreasethedifferencebetweenthesingle-and h(U) 1/U 3. Part C is for the large energy changes re-
∼ | |
multiple-encounter cross sections for m1. But if the cap- sultingfromrare,closeencounterswiththeplanet.Inparts
tured stars are perturbed too much they might not return A andB, h(U) dependsonly on U , butthat symmetryis
| |
Massive Black Hole Binaries 9
brokeninpartCwhereenergygainsaremorefrequentthan with a gain to its kinetic energy. The energy gain results
energy losses. mainly from the interaction with the smaller BH if m1
≫
Panel(a)ofFigure4showsthefinal-velocitydistribu- m2 because the larger BH acts as a fixed potential. The
tionforahardbinaryfromscatteringexperimentsdonein interactionforceofmagnitudeF Gm2/a2 actsforatime
amannersimilartoEverhart’s,sothat“final”meansafter ∆t (a3/GM12)1/2 to produce∼a velocity change ∆v
∼ ∼
a single encounter with the binary. If a star was captured F∆t (m2/M12)Vbin and a corresponding energy change
its final velocity was set to − 2|E∗|, with E∗ measured ∆E∗ ∼∼Vbin∆v ∼(m2/M12)Vb2in. That gives C ∼ m2/µ∼
when the star began returning to the binary for a second 1, which is sufficient to give a hardening rate H1 with no
p
encounter; otherwise the final velocity was set to √2E∗ at dependenceon thehardnessandalmost nodependenceon
the end of the integration. The figure shows the cross sec- the mass ratio.
tions Σ for the final velocity to be greater than v or less If this same derivation is repeated for a high-velocity
f
than v for some positive v . star (v > V ) it gives a hardening rate that rises as
f f bin
F−or the binaries with m1/m2 1 there is a range of H1 v2 instead of falling as H1 1/v2 as it should.
≫ ∼ ∼
velocities for which Σ is symmetric (depends only on v ) That is because the derivation ignores the orbits that lose
f
and varies as Σ 1/v4. This corresponds to Everha|rt’|s energy, which tend to cancel the ones that gain energy.
∼ f
part B. The hardeningrate would benearly zero for these The cancellation removes four powers of v and is the rea-
binaries if multiple encounters were not allowed because son the hardening rate is so difficult to measure at large
the positive and negative contributions would nearly can- v by the Monte Carlo method. For a hard binary there is
cel. The symmetry is not as good for binaries with equal no cancellation because the orbits that lose energy in the
massesforwhichthestarismorelikelytogainenergy.That first encounter are captured and eventually expelled with
is why the N−1/2 scaling in Figure 2 did not work so well anenergygain.Itisnotsurprisingthenthatthehard/not-
forequal-massbinaries.TheasymmetrythatEverhartpre- hard transition occurs at the velocity w where the binary
dictedforpartCisnotclearinthefigure,perhapsbecause begins capturingstars effectively.
the statistics are poor when m1 m2 for the rare, close
≫
encounters with the mass m2 (the m1/m2 = 256 results
come from 106 orbits, but that is still not enough). 3.2 Mass ejection
Panel (b) shows the final-velocity distribution when
thestarsareallowedtoencounterthebinaryasmanytimes Tomeasuretheejection rateweneedanejection criterion,
as necessary for them to be expelled. The Σ 1/v4 scal- i.e. a velocity vej such that a star with initial velocity v is
ing from the single-encounter experiments is r∼aisedfto ap- counted as ejected if it is expelled with final velocity vf >
proximately Σ 1/v3, but the probability of a star being vej. The conventional escape-velocity choice, vej = 2√3σ,
expelled with v∼f Vbfin is still small if m1 m2. leads to a problem for a Maxwellian distribution because
Panel (c) sh≃ows the differential harde≫ning rate with 0.7% of the stars have initial velocities v > vej and will
respect to the final velocity, normalized so that the area be counted as ejected if they receive any energy from the
under the curves is unity. The velocity at the maximum binary, no matter how little. We therefore choose vej =
scales with the mass ratio in the same way as w (eq. 17) max 1.5v,2√3σ ;theresultsdonot dependsensitively on
{ }
and is approximately 1.75w. There is a wide tail to the the numbers 1.5 and 2√3. Let Fej(x,v,σ) be the fraction
right of the maximum if m1 m2 but the high velocities ofstarsincidentuponthebinarywithimpactparameterx
contribute little to the hard≫ening. In fact the hardening andinitialvelocityvthatsatisfythiscriterion.Theejection
rate for m1 m2 can be computed quite accurately by rate is then
≫
considering just the positive velocities from panel (a) and 1 ∞ σ ∞
multiplying the result by two, i.e. by assuming that the J = dv 4πv2f(v,σ) 4π dx xFej(x,v,σ). (23)
H v
captured orbits eventually get expelled with thesame dis- Z0 Z0
tribution of final velocities as for the orbits that are not Theintegraloverthevelocity distributionisevaluatednu-
captured (thisworks only for veryhard binaries). merically after the inner integral is determined from final-
velocity distributions like those shown in Figure 4.
The ejection rate is plotted as a function of σ/V
bin
for five mass ratios in Figure 5. At low velocity J rises
3.1.5 Discussion
as ln(1/σ); at high velocity J falls, first gradually, then
We can now explain why a hard binary hardens at a con- precipitously.Thevelocitydependenceisfitbythefunction
stant rate that is independent of its mass ratio. The ex-
ptwhliaetnhdamtoimo1ningainvmten2cocbnoytmrRiebosuoftsrioo(m1n98too1r)btiihtsseinthhcaoarrtdrheencativn.egItclfioomsrepaleienbscinothuaranyt- J =j1"1+(cid:18)j2Vvbin(cid:19)j3#−j4ln(cid:20)1+(cid:16)j2Vvbin(cid:17)j5(cid:21); (24)
≫
ters with m2 and are expelled with high velocity. But the the five parameters are listed in Table 3. The parameters
scattering experiments show that the dominant contribu- giveaclosefittothedatainthefigurebutareerratic.Note
tion comes from orbits that do not have close encounters that the velocity at the bend in the curves in panel (b) is
withtheBHsandthatareexpelledwithavelocityvf w. not fit well by j2Vbin.
≃
Consideratypicalorbitthatstartswithalowvelocity The ejection rate for a hard binary can be estimated
v,passesatadistancer afrom thetwoBHs,andleaves bynotingthatcloseencounterswiththebinarygiveamean
≃
10 Gerald D. Quinlan
Figure 4.Distributionoffinalvelocitiesforahard,circularbinary(v equalsthelowestvalueinFig.1foreachmassratio).Σisthe
crosssectionforanorbittoleavewithvelocity ≥vf (solidlines),ortoremainboundwithenergy ≤−vf2/2(dotted lines).In(a)the
orbitsencounter the binaryonlyonce; in(b) as manytimes asnecessary forthem toleave withpositive energy. Panel (c) shows the
differentialhardeningratefortheexperimentsin(b).Themassratiosm1/m2 in(a)and(b)areasshownin(c),increasingfromright
toleft.
3.3 Eccentricity growth
m1/m2 j1 j2 j3 j4 j5
1 0.3779 0.9200 2.2572 22.415 0.3437 K1ismoredifficulttomeasurethanH1becauseofthecan-
4 0.1148 0.8815 1.5224 10.521 1.4162 cellationthatoccursinthenumeratorofequation(14).The
16 0.0284 0.6608 0.9404 6.7223 5.6247 B C distribution is wide and nearly centered on the ori-
64 0.0665 0.4438 0.8480 8.1901 2.1824 gin−withamean B C thatis10–100timessmallerthan
256 0.2800 0.0214 3.1294 0.5284 0.8108 h − i
thedeviationaboutthemean.Consequentlyweknowmuch
lessabouteccentricitygrowththanweknowaboutharden-
Table 3.ParametersforfitstoJ (eq. 24)foracircularbinary. ing. Roos (1981) tried tomeasure K1 with only 500 orbits
per measurement. He found K1 = 0.2 0.2 for a hard bi-
±
narywithe=0.6andconcludedthattheeccentricitycould
increase. Mikkola and Valtonen (1992) used 104 orbits per
energy change of C 1. It then follows from the defini-
tion (8) of C thaht aib≃inary must interact with about its measurementtostudythedependenceofK1 ontheeccen-
tricity andhardnessofan equal-massbinary.Theirresults
own mass in stars to shrinkby afactor of e.But “interact
are accurate for hard binaries, for which they found pos-
with”doesnotmean thesameas“eject.” Abinarythat is
nothardinteractswithmanystarsbutejects fewofthem. itive growth rates with a maximum of K1 = 0.19±0.04,
but havelarge error barsfor v V .
And even a hard binary need not eject its own mass to ≃ bin
shrink by a factor of e if it gives some stars much more The results derived here use 10 – 100 times more or-
energy than others. bitspermeasurementthanMikkolaandValtonenused(105
Figure (4)shows that ahardbinarywith m1/m2 1 per measurement at low stellar velocity,106 at high veloc-
≫
expelsfewstarswithv >V .SowhydoesJ notdecrease ity) and use quasi-random numbers rather than random
f bin
rapidlywithm1/m2 intheleft halfofFigure(5)?Because numbers to further reduce the statistical errors. The large
although thefraction of stars expelled with v >V does number of orbits required has two practical consequences.
f bin
decrease rapidly with m1/m2, the fraction expelled with First,onlyasmallnumberofvelocitiesandmassratioscan
v >wdoesnot,anditisthatfractionthatdeterminesthe beexamined.Second,theresultsshouldbeappliedonlyto
f
ejection rate when thebinary first becomes hard. problems where m∗ is much smaller than m1 and m2, for