Monte Carlo Simulations of Globular Cluster Evolution. III. Primordial Binary Interactions

APreprinttypesetusingLTEXstyleemulateapjv.11/12/01

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.PRIMORDIAL

BINARYINTERACTIONS

¨rkan2,K.J.Joshi3,&F.A.Rasio4J.M.Fregeau1,M.A.Gu

AcceptedforpublicationinApJ

arXiv:astro-ph/0301521v2 9 May 2003ABSTRACT

Westudythedynamicalevolutionofglobularclustersusingour2DMonteCarlocodewiththeinclusionofprimordialbinaryinteractionsforequal-massstars.Weuseapproximateanalyticalcrosssectionsforenergygenerationfrombinary–binaryandbinary–singleinteractions.Afterabriefperiodofslightcontractionorexpansionofthecoreoverthefirstfewrelaxationtimes,allclustersenteramuchlongerphaseofstable“binaryburning”lastingmanytensofrelaxationtimes.Thestructuralparametersofourmodelsduringthisphasematchwellthoseofmostobservedglobularclusters.Attheendofthisphase,clustersthathavesurvivedtidaldisruptionundergodeepcorecollapse,followedbygravothermaloscillations.Ourresultsclearlyshowthatthepresenceofevenasmallfractionofbinariesinaclusterissufficienttosupportthecoreagainstcollapsesignificantlybeyondthenormalcorecollapsetimepredictedwithoutthepresenceofbinaries.Fortidallytruncatedsystems,collapseiseasilydelayedsufficientlythattheclusterwillundergocompletetidaldisruptionbeforecorecollapse.Asafirststeptowardtheeventualgoalofcomputingallinteractionsexactlyusingdynamicalthree-andfour-bodyintegration,wehaveincorporatedanexacttreatmentofbinary–singleinteractionsinourcode.Weshowthatresultsusinganalyticalcrosssectionsareingoodagreementwiththoseusingexactthree-bodyintegration,evenforsmallbinaryfractionswherebinary–singleinteractionsareenergeticallymostimportant.

Subjectheadings:celestialmechanics,stellardynamics—globularclusters:general—methods:

numerical

1.introduction

Therealizationabout10yearsagothatprimordialbi-nariesarepresentinglobularclustersindynamicallysig-nificantnumbershascompletelychangedourtheoreticalperspectiveonthesesystems(see,e.g.,thereviewbyHutetal.1992a).Mostimportantly,dynamicalinteractionsbetweenhardprimordialbinariesandothersinglestarsorbinariesarethoughttobetheprimaryenergygenerationmechanismresponsibleforsupportingaglobularclusteragainstcorecollapse(Goodman&Hut19;McMillanetal.1990,1991;Gaoetal.1991).Theterm“binaryburn-ing”isnowoftenusedbyanalogywithhydrogenburningforstars.InthesamewaythathydrogenburningallowsastarliketheSuntoremaininthermalequilibriumonthemainsequenceforatimemuchlongerthantheKelvin-Helmholtztimescale,primordialbinaryburningallowsaglobularclustertomaintainitselfinquasi-thermalequilib-riumandavoidcorecollapseforatimemuchlongerthanthetwo-bodyrelaxationtimescale.

Inaddition,strongdynamicalinteractionsinvolvingbi-nariescanexplainverynaturallythelargenumbersofex-oticobjectsfoundindensestarclusters.Exchangein-teractionsbetweenhardprimordialbinariesandneutronstarsinevitablyproducelargenumbersofX-raybinariesandrecycledpulsarsinglobularclusters(Hutetal.1991;Sigurdsson&Phinney1995;Davies&Hansen1998;Rasioetal.2000).Resonantinteractionsofprimordialbinariesresultindramaticallyincreasedcollisionratesformain-sequencestarsinglobularclustersandevenopenclus-1234

ters(Baconetal.1996;Cheungetal.2003;Leonard19;Leonard&Linnell1992).Directobservationalevidenceforstellarcollisionsandmergersofmain-sequencestarsinglobularclusterscomesfromthedetectionoflargenum-bersofbrightbluestragglersconcentratedinthedenseclustercores(Bailyn1995;Bellazzinietal.2002;Ferraroetal.2001).Previouslyitwasthoughtthatprimordialbinarieswereessentiallynonexistentinglobularclusters,andsoothermechanismssuchastidalcaptureandthree-bodyencountershadtobeinvokedinordertoformbi-nariesdynamicallyduringdeepcorecollapse.However,theseothermechanismshavesomeseriousproblems,andaremuchmorelikelytoresultinmergersthanintheformationoflong-livedbinaries(Chernoff&Huang1996;Kochanek1992;Kumar&Goodman1996;PortegiesZwart&McMillan2002;Kim&Lee1999;Kimetal.1998;Lee&Ostriker1993).Multiplemergersofmain-sequencestarsandrunawaycollisionsinyoungstarclusterscouldleadtotheformationofamassivecentralblackholeinsomesystems(Lee1993;Gebhardtetal.2002;PortegiesZwart&McMillan2002).

Theprimordialbinaryfractionisthereforeakeyinputparameterforanyrealisticstudyofdensestarclusterdy-namics(Hutetal.1992a).Earlydeterminationsofbi-naryfractionsinglobularclusterscamefromobservationsofspectroscopicbinarieswithredgiantprimaries(Pryoretal.1988,see,e.g.,Coteetal.1996foramorerecentstudy)aswellaseclipsingbinaries(Mateoetal.1990;Yan&Mateo1994).HubbleSpaceTelescopeobservationshave

DepartmentofPhysics,MIT37-624A,77MassachusettsAve,Cambridge,MA02139;fregeau@mit.eduDepartmentofPhysicsandAstronomy,NorthwesternUniversity;ato@northwestern.eduIBMCorporation,404WymanStreet,Waltham,MA024;kjoshi@alum.mit.edu

DepartmentofPhysicsandAstronomy,NorthwesternUniversity;rasio@northwestern.edu

1

2

¨FREGEAU,GURKAN,JOSHI,&RASIO

provideddirectconstraintsonprimordialbinaryfractionsinthecentralregionsofmanyglobularclusters,wherebi-nariesareexpectedtoconcentratebecauseofmasssegre-gation.Rubenstein&Bailyn(1997)usedobservationsofabroadenedmainsequenceinNGC6752toderiveabinaryfractionintherange15%–38%fortheinnerclustercore.Theirmethodhasnowbeenappliedtomanyotherclus-ters.Forexample,Bellazzinietal.(2002)deriveasimilarbinaryfraction,intherange0.08–0.38,inthecentralre-gionofNGC288.Addingpropermotioninformationcanleadtomuchtighterconstraints,asinthecaseofNGC6397,whereCool&Bolton(2002)deriveabinaryfraction󰀁5−7%nearthecenter.

Despitethecrucialroleofprimordialbinariesinthedy-namicalevolutionofadensestarcluster,theoverallevo-lutionofthebinarypopulationwithinacluster,anditsdirectimplicationsfortheformationrateofobservablesystemssuchasrecycledpulsarsandbluestragglers,re-mainspoorlyunderstoodtheoretically.Onereasonisthattherelativeimportanceofbinaryinteractionsinacluster,likemanyotherdynamicalprocesses,dependsinacomplexmanneronthenumberofstarsinthesystem.ThismakesitdifficulttoextendresultsobtainedfromsmalldirectN-bodysimulationstorealisticglobularclustermodels.Inparticular,therateatwhichbinariesare“burned”and,ultimately,destroyedorejectedfromtheclusterdependsonthesizeofthecluster.Whentheinitialprimordialbi-naryfractionisbelowacertaincriticalvalue,aglobularclustercorecanrunoutofbinariesbeforetheendofitslife-time,i.e.,beforedisruptioninthetidalfieldoftheGalaxy(McMillan&Hut1994).Withoutthesupportofbinaries,theclusterwillthenundergoamuchdeepercorecollapse,perhapsfollowedbygravothermaloscillations(Sugimoto&Bettwieser1983;Breedenetal.1994;Makino1996).Atmaximumcontraction,thecoredensitymayincreasebymanyordersofmagnitude,leadingtogreatlyenhancedinteractionrates.

Detailednumericalstudiesofglobularclusterevolutionwithprimordialbinariesarestilllacking,forseveralrea-sons.First,theinclusionofevenamodestfractionofprimordialbinariesaddsaverysignificantcomputationaloverheadtoN-bodysimulations.Thisismainlyduetotheextracomputationsrequiredtotreatbinaryinterac-tions,butalsobecausethelifetimeofaclustercanbesignificantlyextended(byuptomanyordersofmagni-tude)throughbinaryburning.Inaddition,indirectN-bodysimulations,theextremelylargeratiooftheoverallclusterdynamicaltimetotheorbitalperiodofclosebina-ries(aslargeas∼1010inaglobularcluster!)introducesmanycomputationaldifficulties.ThismakesN-bodysim-ulationswith4primordialbinariesprohibitivelyexpensiveforN󰀂10stars,althoughspecial-purposesupercom-puterssuchasthenewGRAPE-6mayincreasethislimitinthenearfuture(Makino2001).Orbit-averagedcalcu-lations,likedirectFokker-PlanckintegrationsandMonteCarlosimulations,getaroundthisproblembytreatingbi-nariesjustlikesinglestars,exceptduringbriefperiodsofstronginteractions.Unfortunately,thisrequiresthatcrosssectionsforstronginteractionsinvolvingbinariesbeknownaccurately,forawiderangeofbinaryparameters(masses,semi-majoraxes,andeccentricities).Thesecrosssectionsaredifficulttodetermineingeneral,andreliable

semi-analyticfitstonumericalscatteringexperimentsare

onlyavailableforsimpleconfigurationssuchasthosein-volvingequal-massstars.Forthesereasons,mostprevi-ousnumericalstudiesofglobularclusterswithprimordialbinarieshavebeenlimitedeithertoclusterswithequal-massstars(Gaoetal.1991;Spitzer3&Mathieu1980),ortoverysmallclusterswithN∼10−104stars(Heggie&Aarseth1992;Hurleyetal.2001;McMillanetal.1990,1991;McMillan&Hut1994).Simplifiedtreatmentshavealsobeenemployedinwhichthedynamicsofthebina-rieswasfollowedinastaticclusterbackground(Hutetal.1992b)orinabackgroundclustermodeledasanevolvinggassphere(Giersz&Spurzem2000).

TheresultsofGaoetal.(1991,hereafterGGCM91),basedondirectFokker-Planckintegrations,werethefirsttoclearlyillustratethedominanteffectofevenasmallfractionofprimordialbinariesontheevolutionofaglob-ularcluster.Inthispaper,wepresentthefirststudyofglobularclusterevolutionwithprimordialbinariesbasedonself-consistentMonteCarlosimulationswitharealis-ticallylargenumberofstars(N󰀂105).Partlyinor-dertoallowbettercomparisonofourresultswiththoseofGGCM91,weusesimilarinitialconditionsandcrosssectionsforbinary–binaryandbinary–singleinteractions,eventhoughourmethodforimplementingthesecrosssec-tionsintheMonteCarloschemeiscompletelydifferent.Inaddition,theresultsofGGCM91wereobtainedusinga1-DFokker-Planckmethod(inwhichisotropyinveloc-ityspaceisenforced).Morerealistic2-D(anisotropic)Fokker-Planckcalculationswithprimordialbinarieshaveneverbeenreportedintheliterature,tothebestofourknowledge.Evenforthe1-Dcalculations,andwithonlyasingleparameterrepresentingtheinternalstructureofbinaries(namely,theirbindingenergy),theinclusionofbinary–binaryinteractionssignificantlyincreasedtheover-allcomputationtime.SincetheFokker-Planckmethodusesdistributionfunctionstorepresentthesystem,everynewparameteraddsanewdimensiontothephasespace,makingtheFokker-Planckequationmoredifficulttosolvenumerically.Ithasalsobeenshownrecentlythatthe1-Dtreatmentisinadequateindealingwithsomeaspectsoftheevolution,suchastheescaperatefromtidallytrun-catedclusters(Takahashi&PortegiesZwart1998,2000).ManydifficultiesinthedirectFokker-Planckapproachcomefromthebasicrepresentationofthesystemintermsofsmoothdistributionfunctions.Neglectingthediscretenatureofthesystemmakesitimpossibletofollowthede-tailsofindividualinteractions,suchasbinary–singleorbinary–binaryinteractions.TheimplicitassumptionthatN→∞alsomakesitdifficulttoscaletheresultsofdirectFokker-Plancksimulationstofinitesystemswithdifferentnumbersofstars.

OurMonteCarlomethodprovidesanintermediateap-proach,whichcombinesmanyofthebenefitsofdirectN-bodysimulations(suchasthedescriptionoftheclusteronastar-by-starbasisandthepossibilitytotreateachin-dividualinteractionindetail)withthespeedofanorbit-averagedcalculation.Ourmethodisalso2-Dinvelocityspacebyconstruction,andhenceproperlyaccountsforanyvelocityanisotropyinthesystem.Anotherbenefitofthemethodisthatitallowsawiderangeofbinaryparam-eterstobeusedwithouthavingtomodifytheunderlying

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.

3

orbit-averagedcalculationoftherelaxationprocesses.Inprinciple,individualinteractionscanbetreatedinasmuchdetailasindirectN-bodysimulations,bycomputingallstrongencountersexactlyusingthree-bodyorfour-bodyintegrators.Asafirststep,forthispaper,wehaveincor-poratedathree-bodydynamicalintegratorintoourcode,whichallowsbinary–singleinteractionstobecomputedexactly(withoutreferencetoapproximate,pre-compiledcrosssections).Thisallowsustofollowtheoutcomesofinteractionsmoreprecisely,and,mostimportantly,willallowusinthefuturetoextendourcodetomulti-masssystems,forwhichanalyticcrosssectionsarenotavail-able.

2.treatmentofbinaryinteractions

WeusethebasicH´enon-typeMonteCarlomethodfor

modelingthedynamicalevolutionofclustersasasequenceofequilibriummodelssubjecttoregularvelocitypertur-bations(H´enon1971a,b);ourcodehasbeendescribedindetailbyJoshietal.(2000,2001,hereafterPapersIandII).Theregularvelocityperturbationsarecalculatedus-ingH´enon’smethodtorepresenttheaverageeffectofmanylong-rangesmall-anglegravitationalscatteringencountersusingonesuitablychosenencounterwithanearbystar(H´enon1971b).Ateachtimestep,wecalculatetheMonte-Carlorealizedradialpositionandvelocityofeachstar(as-sumingsphericalsymmetry),whichweusetocalculatewhethertwoobjects(binary-singleorbinary-binary)willinteractstrongly.Thesestronginteractionsareperformedusingeithersimplerecipesbasedoncrosssections,oradynamicalintegrator.Formostoftheworkreportedhere,weusecrosssectionsforthetreatmentofclosebinary–binaryandbinary–singleinteractions.Thesecrosssec-tionswerecompiledfromanalyticfitstotheresultsofnumericalscatteringexperimentsavailableinthelitera-ture.Giventheverylargeparameterspace,reliablecrosssectionsareavailableonlyforequal-massencounters,andsowestudyonlysingle-componentclustersinthispaper.Allsinglestarsareassumedtohavethesamemass,andallbinariescontaintwoidenticalstarswiththesamemassasthebackgroundsinglestars.Allstarsaretreatedaspointmasses,i.e.,weneglectphysicalcollisionsbetweenstarsduringinteractions(cf.Baconetal.1996;Cheungetal.2003).OurimplementationfollowscloselythatusedintheFokker-PlanckstudybyGGCM91,whichwillserveasthemaincomparisonforourwork.

2.1.UnitsandDefinitions

InourcodeweusethesystemofunitsdefinedbysettingG=M0=−4E0=1,whereM0istheinitialclustermass,andE0istheinitialclusterenergy(excludingthebindingenergyinbinaries).Thecorrespondingunitoftimeis

thentdyn(0)=GM5/2

0(−4E0)−3/2.However,thenaturaltimescaleforclusterevolutionistherelaxationtime

tNr(0)≡

0

r3h(0)

rh(0)

log10(γN0)

󰀁

ln(γN0)

󰀁

3

/2

4πGρ,(3)

0

󰀂1andthenumberofstarsinthecoreiscalculatedas

N4πr3c=

cρ0

3m¯

,(4)wherem¯istheaveragemassofastarinthecluster,andρc≃0.5ρ0(Spitzer1987).

2.2.Binary–SingleInteractions

Inasingletimestep,theprobabilitythatabinarywillstronglyinteractwithanotherobject(singleorbinary)isgivenby

P=σwn∆t,(5)whereσisthecrosssectionfortheinteraction,wisthe

relativevelocityatinfinity,nisthelocalnumberden-sityofstars(singleorbinary),and∆tisthetimestep.Forbinary–singleinteractions,σ=σn=nbs,thebinary–singleinteractioncrosssection,andwes,thelocalnumberdensityofsinglestars.Inourcodecalculatenalocalsamplingprocedure,andtakewtobetherelativesusingvelocitybetweenthenearestsinglestarandthebinary.Thetotalcrosssectionforclosebinary–singleinteractionsiscomputedasσbs=πb2max.Herebwhichgivesadistanceofclosestmaxistheimpactpa-rameterapproachbetween

4

¨FREGEAU,GURKAN,JOSHI,&RASIO

thebinaryandthesinglestarofrmin=3.5a,whereais

thebinarysemi-majoraxis.Forabinaryofmassmsinglestarofmassmwebandb2max=r2

min󰀁

have

1+

2G(m+mb)dy

=12.48πa

2

󰀁

w

/3

Gma2

mw2+0.04|ǫ2|

󰀂󰀄24

y󰀁

1+

7

2mw

2

,yielding

σbb=31.8

Gma2

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.

5

Theenergyrequiredtodisruptthesofterbinary,aswellasthetotaltranslationalenergyreleasedinthecollision∆Earebothgeneratedattheexpenseofthesurvivingbinary.t,Thusthebindingenergyofthesurvivingpairincreasesbyanamountǫ2+y(ǫ1+ǫ2).AccordingtoMikkola(1983a),forcollisionsbetweenbinariesofequalbindingenergiesproducingabinaryandtwosinglestars,typicallyabout1/4ofthetranslationalenergyproducediscarriedawaybythebinary,andtheremainingisdistributedrandomlyamongthetwosinglestars.Forsimplicity,weassumethatthisprescriptionisapplicabletocollisionsbetweenbina-riesofunequalbindingenergiesaswell.Weselectthedirectionoftherecoilvelocitybetweenthebinaryandthesinglestarsrandomlyinthecenter-of-massframe.

Ifabinarydoesnotundergoastronginteractionwithasinglestaroranotherbinary,itisthentreatedasasinglestarintheusualtwo-bodyrelaxationstep(seePaperI),duringwhichitsinternalstructureisleftunchanged.

3.results

3.1.InitialConditionsandSummaryofModelResultsForourinitialclustermodelsweuseboththePlummermodel,assumedtobeisolated(i.e.,withnotidalboundaryenforced),andavarietyoftidallytruncatedKingmodels,assumedtobeonacircularorbitintheGalaxy(i.e.,withafixedexternaltidalpotential).MasslossthroughthetidalboundaryistreatedasinPaperII,usingacriterionbasedontheapocenterdistanceofeachstellarorbitinthecluster,andaniterativeproceduretodetermineboththemasslossandthenewpositionofthetidalboundaryaftereachrelaxationtimestep.Theinitialbinaryfractionfb(definedasthefractionofstars,bynumber,thatarebinaries)variesbetween0and30%.Inafewcases,forcalibration,wehavealsoperformedsimulationsinwhichthebinariesarepresent,butallinteractionsareturnedoff;thesemodelsareequivalenttotwo-componentmod-elsinwhichasmallfractionof(single)starshavetwicethemassofthebackgroundstars(seeWattersetal.2000andFregeauetal.2002forotherstudiesoftwo-componentclustersusingourcode).

Thebinariesaredistributedinitiallyintheclusterac-cordingtothesamedensityprofileasforsinglestars.Hencenoinitialmasssegregationisassumedforthebina-ries.Thedistributionoftheinternalbindingenergyofthebinariesisassumedtobeuniforminlogǫbetweenamin-imumvalueǫminandamaximumvalueǫmax.FollowingGGCM91,weconsideronlyhardbinaries,withthemin-imumbindingenergyǫmin=mσc(0)2,whereσvelocitydispersion.Softbinaries,c(0)istheinitialcentralifpresent,wouldbeassumedtobeionized(destroyed)assoonastheyparticipateinastronginteraction.Thereforetheywouldnotaffecttheoverallevolutionoftheclustersignificantly.Forthemaximumbindingenergywetakeǫwhichisapproximatelythebindingenergymax=133ǫofacontactmin,binaryfortwosolar-likestarsifσc(0)≃10kms−1.Theprecisevalueofǫveryhardbinariesmaxhaslittleinfluenceonourresults,sincebehaveessentiallyassinglemoremassivestars(withverysmallinteractioncrosssection).Table1liststheparametersofthemainmodelswecon-sidered,aswellasthemainresultsofoursimulationsforeachcluster.Thefirstcolumnidentifiestheinitialclus-termodel,PlummerorKing,andthevalueofthecon-

centrationparameterWmodels.Thesecond0(dimensionlesscentralpotential)forKingcolumngivestheinitialbi-naryfractionfN=3×105starsb.Allsimulationswereperformedwith(includingbinaries)initiallyintheclus-ter.Thefollowingcolumnssummarizethemainresultsofourdynamicalsimulations.Foreachmodelwefirstgivethetimeofcorecollapsettimetcc,inunitsoftheinitialhalf-massrelaxationrh(0),definedbyeq.(2).Herecorecollapseisdefinedasthemomentwhenthecoredensityreachesitsfirstmaximum.Thiscanbedeterminedtypi-callytowithinastatisticalerrorofatmostafewpercentinoursimulations.Wethengivethetotalmassoftheclusteratthemomentofcorecollapse(inunitsofitsini-tialtotalmass),andthefractionofbinariesthatremainatthatmoment.Forclustersthatdisruptcompletelybe-forereachingcorecollapse,welistthedisruptiontimetinsteadoftdiscc.

3.2.ComparisonwithDirectThree-BodyIntegrationAsasimpletestofourcodeandtheapproximatetreat-mentofinteractions,wecomparetheuseofcrosssectionswithdynamicalintegrationsofbinary–singleencounters.Infuturework,wewillalsoimplementdynamicalintegra-tionsofbinary–binaryinteractions,andwewillusemoredetailedcomparisonstore-calibratethevariousrecipesbasedoncrosssections.Hereourintentismerelytodemonstratethatthesesimplerecipesarereasonablyaccu-rate.Wehavenotchangedourprescriptionstotrytobet-termatchtheresultsofthedynamicalintegrations,sinceamaingoalinthisfirststudyistoprovidecomparisonswiththeFokker-PlancksimulationsofGGCM91thatusedthesamesimpleprescriptions.Forthistestbinary–binaryinteractionswereturnedoff.Inreality,theytendtodom-inatetheenergyproduction(seeSec.3.3).Thusthissim-pletestalsoallowsustostudyspecificallytheeffectsofthree-bodyinteractionsontheoverallclusterevolution.Figure1showstheevolutionofanisolatedclusterde-scribedinitiallybyaPlummermodelwithN=3×105starsand20%binaries.Solidlinescorrespondtothesim-ulationusingdirectthree-bodyintegrations,whiledashedlinesshowtheresultsusingoursimplecrosssections.Thetoppanelshowsthetotalmassinbinariesintheclus-ter,decreasingasbinaryburningproceeds.Sinceallbi-nariesinthemodelarehard,binary–singleinteractions(unlikebinary–binaryinteractions)cannotdestroyabi-nary,andthereforebinariescanonlybelostbyejectionfromfromthecluster(typicallyfollowingsignificanthard-eningthroughmultipleinteractions;seeHutetal.(1992b)andSec.3.5below).Therateofbinaryejectionacceler-atesabruptlyatt/trh(0)≃8−10nearcorecollapse.ThemiddlepanelofFigure1showstheenergygeneratedinbinary–singleinteractions,asafractionofthetotalinitialbindingenergyofthecluster.Bythetimeofcorecollapse,thisisonly∼0.1.Thisamountofenergyisnotsufficienttodelaycorecollapsesignificantly.Infactthebinaries,throughmasssegregation,acceleratecorecollapseinthis(artificial)simulation(recallthatthecorecollapsetimeofasingle-componentPlummermodelwithoutbinariesisgivenbytcc/trh(0)≃14).Thebottompanelshowsvari-ouscharacteristicradiiinthecluster:fromtoptobottom,thehalf-massradiusofsinglestars,thehalf-massradiusofbinaries,andthecoreradius,allinunitsoftheinitial

6

¨FREGEAU,GURKAN,JOSHI,&RASIO

half-massradius.

Theagreementbetweenthetwomethodsisstrong,al-thoughthetotalenergygeneratedinbinary–singleinterac-tionsisslightlysmallerwhencalculatedbydirectdynami-calintegrations.Consequently,themodelusingdynamicalintegrationsreachescorecollapsesoonerthanthemodelusingcrosssections,becauselessenergyisgeneratedtosupportthecoreagainstcollapse.Webelievethatthisdifferencecomesfromthedeterministictreatmentofbi-naryhardeningwithcrosssections,inwhicheverybinary–singleinteractionresultsinahardenedbinary.Inrealitythewidesthardbinariesinthesimulation,whicharerightaroundthehard/softboundary(andhavethelargestin-teractioncrosssection)haveroughlyequalprobabilitiesofhardeningandsofteninginaninteraction(Heggie1975).Topartlyrestoreconsistencybetweenthetwotreatments,weignoredynamicalintegrationoutcomesthatresultinionizationofthebinary.Weretheseincluded,thetotalenergygeneratedinbinary–singleinteractionswouldde-creasefurther,byroughly50%.Thiswouldcausethebi-narypopulationtobecomemorecentrallyconcentrated.Thus,inarealisticclustersimulation,wewouldexpecttheratiooftheenergygeneratedbybinary-singleinter-actionstobinary-binaryinteractionstodecreasebymorethan50%comparedtopredictionsofcross-sectionbasedrecipes.

3.3.IsolatedClusters

WeconsiderfirsttheevolutionofPlummermodelscon-tainingN=3×105starswitharangeofbinaryfractionsfb.Asafurthertestofourmethod,weshowinFig.2theevolutionofthevariousenergiesandthevirialratioofasystemwithfintoournumericalb=0.1.Sincedynamicalrelaxationisnotbuiltmethod,thedegreetowhichvirialequilibriumismaintainedduringasimulationisourmostimportantindicatorofnumericalaccuracy.Wemonitorthis,aswellasenergyconservation,inallourruns,andterminateacalculationwheneverthesequantitiesdeviatefromtheirexpectedvaluesbymorethanafewpercent(thistypicallyhappenswhenthenumberofbinarieshasbeenreducedtoaverysmallvalue,or,intidallytruncatedclusters,whenthetotalnumberofstarsremainingintheclusterbecomesverysmall;SeeSec.3.4below).

Figures3,4and5showtheevolutionofmodelswithfb=0.02,0.1,and0.2,respectively.Themainimpactofintroducingbinariesinthemodelsisveryclear:corecol-lapseisdelayedconsiderably.Evenforaclusterwithonly2%binariesinitially(Fig.3),tccincreasesbymorethanafactor2.Clusterswithfb≃0.1−0.2canavoidcorecol-lapsefor∼100tclustersrh(0)(Figs.4and5).ForthevastmajorityofglobularinourGalaxy,wheretIfallrh(0)∼109yr,thistimescaleexceedsaHubbletime.globularclus-tersinourGalaxywerebornwithfb󰀂0.1,onlythosewithveryshortinitialrelaxationtimeswouldhavehadachancetoreachcorecollapse.However,forrealclusters,tidaltruncationandmassloss(Sec.3.4)aswellasstellarevolution(PaperII)complicatethispictureconsiderably.InFigure3,wealsoshowforcomparisontheevolu-tionofthecoreradiusforamodelinwhichbinariesarepresentbutallinteractionsareturnedoff(short-dashedlineinthebottompanel).Evenwithabinaryfractionassmallasfb=0.02inthiscase,corecollapseoccurs

significantlyearlierthaninasinglecomponentPlummer

model(attcc/trh(0)≃10insteadof14).Thisshowstheexpectedtendencyfortheheaviercomponentofbinariestoacceleratetheevolutiontocorecollapse,andtheresultisingoodagreementwithpreviousstudiesofcorecol-lapseintwo-componentclusters(see,e.g.,Wattersetal.2000).Notethatforsufficientlylargebinaryfractions,thesetwo-componentmodelsbecome“Spitzer-unstable,”i.e.,thecorecollapseisdrivenentirelybytheheaviercom-ponent.UsingthestabilitycriterionderivedbyWattersetal.,Λ≡(M2/Mm1)(m2/m1)2.4<0.32,herewithanindi-vidualmassratio2/m/(1−1=2andatotalcomponentmassratioM2/M1=2fbfb)weexpecttheSpitzerinsta-bilitytoappearwheneverthebinaryfractionfb󰀂0.03.Thusallourmodelswithbinaryfractionsaboveafewpercentshouldevolveonarelaxationtimescaletoastatewherethedynamicsoftheclustercoreislargelydominatedbythebinaries.Indeed,lookingatthemiddlepanelsofFigures4and5,weseethat,withfb=0.1−0.2,theenergygenerationislargelydominatedbybinary–binaryinteractions.Incontrast,forfb=0.02(middlepanelofFig.3),binary–binaryandbinary–singleinteractionscon-tributeroughlyequally.

Toquantifytheeffectofprimordialbinaryburningonthecorecollapsetime,wehaverepeatedcalculationswithbinariespresentbutallinteractionsturnedoffforsixdiffer-entmodelswithvaryingfb.Foreachmodel,wecanthenproperlycalculatetheratioofcorecollapsetimeswithandwithoutbinaryinteractions(butwithmasssegrega-tioneffectspresentinbothcases).TheresultsareshowninFigure6,wherethisratioisplottedasafunctionofthebinaryfraction.Asimplelinearfitgives

tcc≃tcc(fb=0)×(75fb+1),

andreproducesthenumericalresultstowithin∼30%intherangef“tb=0−0.3.Thenotationweusehere,cc(fb=0),”meansthecorecollapsetimeofaclusterwiththesamefractionof“inactive”binaries,ratherthanwithnobinaries.

AlsoshowninFigure3forcomparisonistheresultofasimulationinwhichallinteractionsareincluded,butbinary–singleinteractionsarecalculatedbydirectthree-bodyintegrationsasinSec.3.2(long-dashedlineinthebottompanel).Thiscomparisonisusefulagainasatestofthesimpletreatmentbasedoncrosssections,sincebinary–singleinteractionsplayanimportantroleasasourceofenergyinthismodelwithfb=0.02.OurconclusionisthesameasinSec.3.2:theagreementisverygooduntilt/trh(0)≃15,butthenthetwosimulationsdivergeandthemodelcomputedwithdirectthree-bodyintegrationscollapsesslightlyearlier(tcc/trh(0)≃19insteadof22).Inspiteofthisslightoffset,aftercorecollapsethecorere-expansionandgravothermaloscillationsalsolookverysimilarinthetwosimulations.

ThetoppanelsinFigures3–5showtheevolutionofthetotalclustermass,asafractionoftheinitialmass,andtheremainingmassinbinaries,asafractionoftheinitialmassinbinaries(thisisalsotheremainingfrac-tionbynumbersinceallbinarieshavethesamemass).Binary–binaryinteractionsarethemainprocessrespon-sibleforthedestructionofbinariesinthesesimulations(sincethesofterbinaryisassumedtobedisruptedineachinteraction).Intheabsenceofevaporationthroughatidal

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.

7

boundary,masslossfromtheclustercomesalmostentirelyfromstarsandbinariesejectedthroughrecoilfollowinganinteraction.Themasslossratethereforeincreaseswithin-creasingbinaryfraction.Atcorecollapse,thetotalmasslossfractionisabout5%,15%,and25%forfb=0.02,0.1,and0.2,respectively.However,whilethetotalnum-berofbinariesintheclustercontinuouslydecreases,theremainingfractionofbinariesatcorecollapseappearstoberoughlyconstant,around0.2,independentoffsufficienttopowermanycyclesofgravothermalb.Thisisoscil-lationsaftertheinitialcorecollapse,evenforinitialbi-naryfractionsassmallasafewpercent.Forfb=0.2(Fig.5),wewereabletoextend3ournumericalintegrationallthewaytoalmost∼10tarestillpresentrh(0),atwhichpointseveralthousandbinariesinthecentralregionofthecluster(thistimescalewouldofcoursevastlyexceedaHubbletimeformostGalacticglobularclusters!).

ThebottompanelsinFigures3–5showtheevolutionofseveralcharacteristicradii.Themostimportantisthecoreradiusrc(recallthat,byourdefinition,eq.(3),thecentraldensityscalesapproximatelyasρtralvelocitydispersionisapproximately0∝rc−2sincethecen-constant).Evenindeepcorecollapse,thecoreradiusofourmodelsneverdecreasesbymorethanafactor∼100(correspondingtoanincreaseinthecentraldensityby∼104).Modelswithhigherbinaryfractionscontractverylittle(seeFig.5:thefirstanddeepestcorecollapsecorrespondstoadecreaseinrcbylessthanafactor10).Alsoshownarethehalf-massradiiofthebinariesroftheh,singlebandofthesinglestarsrstarsalwaysincreasesh,s.Thehalf-massradiusmono-tonicallyfortheseisolatedclusters.Incontrast,thehalf-massradiusofthebinariestendstoincreaseonaveragebutshowsamuchmorecomplexbehaviorthatdependsstronglyonthebinaryfractionandontheparticulardy-namicalphaseintheevolutionofthecluster.Thetrendisforrh,btodecreaseduringnormalclusterevolution,asthebinariesmasssegregatetotheclustercore,andtoin-creasedramaticallyduringcorecollapse,asthedensityofbinariesinthecoregrowsandtherateofbinary–binaryinteractionsgrowsmorequicklythantherateofbinary-singleinteractions.Thiscausesmanysofterbinariesinthecoretobedisruptedandmanyharderbinariestobeejectedoutofthecorethroughrecoil.Forsufficientlylowbinaryfractions(Figs.3and4),rrcorecollapse.h,beventuallybecomeslargerthanh,safterForhighbinaryfrac-tions(Fig.5),thebinariesremainalwaysmuchclosertothecenterofthecluster.

WenowturntoamoredetaileddiscussionofthePlum-mermodelwith10%binaries(Fig.4),includingacom-parisonwiththeFokker-PlanckresultsofGGCM91(seetheirFigs.1–3),whoconsiderthistheir“standardmodel.”Qualitatively,ourresultsareinverygoodagreementuptocorecollapse.Afteraninitialphaseofcontractionlasting∼10trh,thecoreradiusbecomesnearlyconstantandtheclusterentersalongphaseofquasi-thermal-equilibrium.Thisisthestable“binaryburning”phase,analogoustothemainsequenceforastar.Duringthisphase,therateofenergyproductionthroughinteractionsinthecoreisbalancedbytherateatwhichenergyflowsoutintheouterhalo,whichcontinuouslyexpands(intheabsenceofatidalboundary).Corecollapseoccursrathersuddenlyattheendofthisphase.GGCM91findtcc/trh(0)≃50

forthismodel,whilewefindtisfollowedbycc/tgravothermalrh(0)≃70.Thisini-tialcorecollapseoscillations,whichareclearlystillpoweredbyprimordialbinaryburn-ing.Wewereabletofollowtheseoscillationsaccuratelyuntilt/trh(0)󰀂200,whileGGCM91terminatetheircal-culationatt/trh(0)≃90.

Uponcloserexaminationandquantitativecomparisons,somemoresignificantdifferencesbecomeapparent.First,weseethattheinitialcontractionphaseappearsmuchdeeperinthemodelofGGCM91,withrcdecreasingbyalmostanorderofmagnitude,whileinourmodelthecorecontractsbyafactorofabout3.Duringthestablebinaryburningphase,theclusteralsoappearssomewhatmorecentrallyconcentratedinthemodelofGGCM91.Justbeforecorecollapse,theyfindrvalueis≃0.04.Ontheotherc/rhand,h(s)≃0.01,whileourtherateofbinaryburninganddestructionisnearlythesameinthetwomod-els.Compare,forexample,theevolutionofMb/Mb(0)inFigure4tothesamequantityplottedinFig.2aofGGCM91.Althoughthereareslightdifferencesintheshapesofthetwocurves,thereductionto0.8occursafterabout10trh(0)inbothcases,andthereductionto0.5afterabout28trh(0).Byt/trh(0)≃70thenumberofbinarieshasbeenreducedto0.2ofitsinitialvalueinbothmodels.Thisagreementisespeciallysurprisingsinceinourmodelthisisstill(just)beforecorecollapse,whileinGGCM91’smodelseveralcyclesofgravothermaloscillationshaveal-readyoccurred.

ThereareseveralreasonstoexpectdifferencesbetweenourresultsandthoseofGGCM91’sFokker-Plancksimu-lations,eventhoughourtreatmentsofindividualbinary–singleandbinary–binaryinteractionsareessentiallyiden-tical.

First,GGCM91’srepresentationofbinariesisintermsofaseparablecontinuousdistributionfunctioninE,theorbitalenergyinthecluster,andǫ,theinternalenergyofthebinary.Infact,thereisastrongandcomplexcorrela-tionbetweenabinary’sbindingenergyanditspositioninthecluster(orequivalentlyitsenergyinthecluster),withharderbinariesconcentratedneartheclustercore(seeHutetal.1992bandSec.3.5).WesuspectthatGGCM91’schoiceofaseparabledistributionfunctionhastheeffectofreducingtheenergygenerationrate,sincethenpro-portionatelymoresoftbinarieswillbechosenforbinaryinteractions—interactionsthatpredominantlyliberateaconstantfractionofthetotalbindingenergyavailable(seeSec.2andHeggie1975).

Second,1-DFokker-Planckresultsareknowntodifferfrom2-Dresultsingeneral(mostnotablyinthepredic-tionoftherateoftidalstripping;seePaperII).Enforcingisotropyinthestellarvelocitydistributionislikelytoaf-fectthedynamicsofthecorearoundthetimeofcollapse,whenthisdistributionmaybechangingrapidlyandtheincreasedinteractionratemaybecausinganisotropy.In-deed,Baumgardtetal.(2003)haverecentlyfoundwithN-bodysimulationsthatanisotropyneartheclustercen-terbecomessignificantduringcorecollapse.

Third,theonlyexplicitdependenceonNintheFokker-PlanckapproachisthroughtheCoulomblogarithm,so,eventhoughGGCM91setN=3×105fortheirtreatmentofinteractions,itisnotclearinwhatsensetheirresults,whichassumeasmooth,continuousdistributionfunction,

8

¨FREGEAU,GURKAN,JOSHI,&RASIO

malequilibrium.Thehighertheinitialbinaryfractionand

centraldensity(seebelow),thestrongerthetendencyforthecoretoexpandinitiallyinsteadofcontracting.Second,corecollapse5orcompletedisruptionalwaysoccursinlessthanabout45trh(0).Forfb≥0.1,thedisruptiontimetdisdecreaseswithincreasingbinaryfraction,andcompletedisruptionoccursbeforeanydeepcorecollapse.Thisisincontrasttothemodelswithlowerbinaryfractions(Figs.8and9),wherecorecollapsefollowedbygravothermalos-cillations(similartothoseobservedforisolatedclustersintheprevioussection)occurbeforedisruption.Thepos-sibilityforaclustertosuffercompletedisruptionbeforecorecollapseisaqualitativelynewbehaviorintroducedbyprimordialbinaries.Indeed,allKingmodelswithoutbinaries(andwithoutstellarevolution)reachcorecollapsebeforedisrupting(seePaperIIandQuinlan1996).

Figures12and13showtheevolutionofKingmodelswith10%binariesbutdifferentvaluesofW0.FortheverycentrallyconcentratedclusterwithW0=11(Fig.12),sig-nificantcoreexpansionoccursinthefirstfewrelaxationtimes(withrcincreasingbyaboutanorderofmagnitude).ThisisamoreextremeexampleofthebehavioralreadynotedinFig.11.Thefinalevolutionofthisclusterisalsopeculiar:thisisoneoffewexamples(seeTable1)weencounteredwherethebinariesarecompletelyexhaustedbeforetheclusterdisrupts.Att/trh(0)≃30,about20%oftheinitialclustermassremainsinsinglestars,andtheclusterundergoesdeepcorecollapse.Sincetherearenobi-nariesleft,andoursimulationsincludenoothersourceofenergy,nore-expansioncanoccurandwemustterminatethecalculation.

ForamodelwithW0=3,whichhasamuchmorenearly-uniformdensityprofileinitially,completedisrup-tionoccursbeforecorecollapseevenwithmuchlowerbi-naryfractions(seeTable1).Forthemodelwithfb=0.1showninFig.13,disruptionoccursatt/trh(0)≃15.Forcomparisonasingle-componentW0=3Kingmodelun-dergoescorecollapseatt/trh(0)≃12(PaperII).Alsonotehowthecorecontractsthroughouttheevolutionatanearly-constantratemuchfasterthaninmodelswithhighervaluesofW0(compare,e.g.,Fig.10).Justbeforefinaldisruption,thecoreradiushasdecreasedbyaboutanorderofmagnitudefromitsinitialvalue.

OnlyafewsmallN-bodysimulationsoftidallytrun-catedclusterswithprimordialbinarieshavebeenreportedpreviously(Heggie&Aarseth1992;McMillan&Hut1994).DetailedcomparisonsarenotpossiblebecausethesestudiesassumedratherdifferentinitialmodelsandtheN-bodyresults(forN∼1000−2000)areverynoisy.However,wedoseegoodqualitativeagreement,withmasslossrates∼10timeslargerthanforisolatedclusters,andcompletedisruptionalsoobservedafterafewtensofini-tialhalf-massrelaxationtimesinallN-bodysimulations.WearenotawareofanypreviousFokker-Plancksimula-tionsoftidallytruncatedclusterswithprimordialbinaries(GGCM91consideredonlyisolatedPlummermodels).

3.5.EvolutionoftheBinaryPopulation

correspondtothisparticularvalueofN.

Finally,wepointoutthatourresultsfortheinitialcon-tractionandcoresizeduringthebinary-burningphaseareinmuchbetteragreementwiththoseofdirectN-bodysim-ulationsincludingprimordialbinaries.Indeed,bothHeg-gie&Aarseth(1992)(seetheirFigs.5and18)andMcMil-lanetal.(1990)(seetheirFig.1)find,aswedo,thattheclustercorecontractstypicallybyafactorofaboutthreefromitsinitialsize.GGCM91,ontheotherhand,findcorecontractionbyanorderofmagnitude,adirectindicationthattheirmethodunderestimatestheenergygenerationrate.

Perhapsamoresignificantdifferencebetweenourre-sultsandthoseofGGCM91isinthepost-collapseevo-lution.GGCM91findmuchmorefrequent,erratic,anddeepergravothermaloscillationcycles.Ourmodelshowsalmostquasi-periodicoscillationswithperiod∼40trhandpeak-to-peakamplituderc,max/rc,min∼100.Instead,GGCM91find7oscillationsofwidelyvaryingperiodsbe-tweent/trh(0)=50and90,andrc,max/rc,min∼103.Webelievethismaypossiblybeduetodifferencesinthenu-mericalmethodofcalculatingNc(foruseineq.(1)—seediscussioninSec.2.1),althoughitisnotclearfromtheirpaperwhatmethodtheyactuallyuse.Tochecktheva-lidityofourresultsnearcorecollapse,wehaveexaminedmorecarefullythedynamicsgoverningcorere-expansionaftercollapse.InFigure7,weshowtheevolutionofthetemperatureprofileintheclusterasthesystemundergoesacorecollapseandrebound.Thetemperatureintheclus-ternormallydecreasesoutwardeverywhere,asinastarinthermalequilibrium.However,duringdeepcorecollapse,a“temperatureinversion”developsforashorttime.Thistemperatureinversionisresponsiblefordrivingtherapidre-expansionofthecore,asenergyisnowflowinginward.Thismechanismhasbeenpredictedtheoreticallyforalongtime(Sugimoto&Bettwieser1983;Heggie&Ramamani19),andhasbeenobserveddirectlyinrecentN-bodysimulations(Makino1996).However,toourknowledge,oursisthefirstnumericaldemonstrationofthiseffectforaclustercontainingprimordialbinaries.Inallpreviousstudies,thebinarieswereassumedtoformdynamicallyviathree-bodyinteractionsortwo-bodytidalcapturesduringdeepcorecollapse.Asnotedintheintroduction,thesemechanismarenowconsideredunrealistic,astheymostlikelyleadtostellarmergers(whichwerenottakenintoaccountinthepreviousstudies).

3.4.TidallyTruncatedClusters

Wenowpresentourresultsformorerealistic,tidallytruncatedclusters.Figures8–11showtheevolutionofKingmodelswithW0=7andinitialbinaryfractionsfb=0.02,0.05,0.1,and0.2,respectively.Severalstrikingdifferenceswithisolatedclustersareimmediatelyappar-ent.First,weseethattheinitialcorecontractionphaseisabsent.Forfb≤0.1,thecoreradiusdecreasesslowlyandmonotonicallyallthewaytocollapse.Forhigherbi-naryfractions(Fig.11),thecoreradiusincreasesslightlyatfirst.Thisissimplybecausetheinitialbinaryburningrateinthismodelisclosetotherateneededtoreachther-5

Itshouldbenotedthat,inmuchoftheliteratureonglobularclusterdynamics,theterm“corecollapse”isusedtorefertotheinitialcorecontractionphase(whichweseeherecanactuallybeexpansioninstead),andwhatwecallthebinaryburningphaseisthencalledthe“post-collapse”phase.Clearlythisterminologynolongermakessense,andshouldbeabandoned.Whatwecall“corecollapse”inthispaperreferstothebriefepisodesofdeepcorecollapseattheonsetofandduringgravothermaloscillations.

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.

9

Inadditiontoaffectingtheglobalclusterevolution,bi-naryinteractionsalsostronglyaffectthepropertiesofthebinariesthemselves.Thestudyoftheevolutionofapri-mordialbinarypopulationdatesbacktotheseminalworkofHeggie(1975),butitisonlyrecentlythatdetailednu-mericalsimulationsoflargebinarypopulationsinglobularclustershavebeenperformed(Hutetal.1992b;Giersz&Spurzem2000).WecanusetheresultsofourMonteCarlosimulationstostudythedynamicalevolutionofthebinarypopulation.

Fig.14showstheevolutionofthebinaryfraction,fdifferentregionsofourevolvingKingmodels.Thereb,inisacleartrendforfbtoincreaseinthecoreanddecreaseinthehalowithtime,aswellasatrendforfWbtogrowmorewithsmallertendencyforbinaries0.Inspiteofmasssegregationandthetodominatethecentralregionofacluster(followingthedevelopmentoftheSpitzerin-stability;seeSec.3.3),corebinaryfractionsrarelyexceed0.5inourmodels.Thustherangeofinitialbinaryfrac-tionsweconsiderareatleastinroughagreementwiththemeasurementsofcorebinaryfractionsinglobularclusterstoday(∼0.1−0.4;seeSec.1).Note,however,thatthepresent-daybinaryfractioninthecoreofaclustercannotberelatedsimplytothecluster’sinitialbinaryfraction,asitmaydependinacomplicatedwayonseveralinitialparameters.Forexample,weseethataninitialWwithf.1has,duringmostofitsevolution,0=11modelaboutthesameb=0corebinaryfraction(≃0.15)asanini-tialW0=7modelwithfdefinitionoffb=0.05.Inaddition,recallthatourbinaries.Forbinthesesimulationsincludesonlythehardreasonabledistributionsofprimor-dialbinaryseparations,includingseveralmoredecadesonthesoftside,thetrueinitialbinaryfractionintheclustermighthavebeen∼2−3timeslargerthanourquotedvalueoffb(thisisthereasonwhywedidnotconsidervaluesoffb󰀂0.3,whichcouldnotberealistic,unlessdynamicsalreadyplaysanimportantroleduringtheprocessofstarformation;seeClarkeetal.2000).

Fig.15showstheevolutionoftheprimordialbinarypopulationinaW0=7Kingmodelwithfshowsthedistributionofbindingb=0.2.Each2-Dhistogramenergies(initiallyflatinlogǫ)andradialpositionsinthecluster.Inadditiontothecleartendencyformasssegregationandhardeningofthebinaries,wenotethedevelopmentofastrongcorrelationbetweenhardnessandradialdistribu-tion:harderbinariestendtoconcentrateintheclustercoremuchmorethansofterbinaries,inspiteofhavingallthesamemass(andincontrasttothefundamentalassumptionmadeintheFokker-PlanckcalculationsofGGCM91).Thegeneraltrendsobservedhereareingoodqualitativeagree-mentwiththeresultsofpreviousstudiesusingmoreide-alizedmodels(Hutetal.1992b;Giersz&Spurzem2000,see,e.g.,theirFig.25).

NeartheendoftheevolutionshowninFig.15(butalready∼10trh(0)beforecompletedisruption,whentheclusterstillretainsabout40%ofitsinitialmass),apartic-ularlystrikingsituationdevelopswhereallthesurvivingbinariesintheclustercoreareextremelyhard.Recallthatourinitial2upperlimitonthebindingenergyofabinary(∼10kT,where“kT”istheaveragekineticenergyofstarsinthecore)roughlycorrespondstocontactfortwosolar-likestars.Therefore,mostofthebinariesremaining

after∼30t103kTrh(0),withbindingenergiesnowintherange∼102−,wouldhavemergediftheycontainedsolar-likestars(perhapsformingbluestragglers).Ofcourse,inarealcluster,manyofthesebinariescouldcontaincom-pactobjects(mostlikelyheavywhitedwarfsandneutronstars)andwouldthenhavesurvived.Wecannotaddressanyoftheseissueshere,sinceoursimulationsareclearlytooidealized,butwepointoutthatglobularclustercoresareindeedobservedtocontainlargepopulationsofbluestragglers,WUMabinaries(eclipsingsystemscontainingtwomain-sequencestarsinacontactconfiguration;see,e.g.,Albrowetal.2001),andavarietyof“ultracompact”binariescontainingneutronstarsandwhitedwarfs(themostextremeexamplebeingperhapsthe“11-minute”X-raybinaryinNGC6624;see,e.g.,Deutschetal.2000).

4.summaryandcomparisonwithobservations

Wehaveperformed,forthefirsttime,discretesimula-tionsofglobularclusterswithrealisticnumbersofstarsandprimordialbinaries,usingour2DMonteCarlocodewithapproximateanalyticalcrosssectionsforprimordialbinaryinteractions.

Wehavecomparedtheuseofcrosssectionswithex-act,dynamicalintegrationsofbinary–singleencounters,andfindthattheagreementbetweenthetwomethodsisstrong,althoughourcurrentimplementationofthecrosssections,basedontheFokker-PlanckstudybyGaoetal.(1991),tendstooverestimateslightlytheenergygenera-tionrate.Consequently,modelsthatusecrosssectionstendtooverestimatecorecollapsetimesforclustersinwhichbinary–singleinteractionsdominate.However,wefindthatbinary–binaryinteractionsdominatetheenergygenerationforfb󰀂0.03,aresultthatisinquantitativeagreementwithasimpleSpitzer-typestabilitycriterionappliedtothecomponentofbinaries.

Wehavestudiedtheevolutionofisolatedclusterswithvaryingbinaryfractions,andhavefoundthatthepres-enceofevenasmallfractionofbinariesissufficienttodelaysignificantlytheonsetofcorecollapse.Isolatedclusterswithabinaryfractiongreaterthanabout0.1–0.2canavoidcorecollapseforasmuchas∼102−103tcorecollapserh(0).Wefindasimplelinearrelationbetweenthetimeofacluster,twithbinarycc,andthecorecollapsetimeofthesameclusterinteractionsturnedoff(butmasssegregationstillpresent),denotedt1).Wehavecc(fb=0),givenbytcc≃tcc(fb=0)×(75fb+comparedourre-sultswiththoseofGaoetal.(1991),andfindreasonableagreement,withnearlyidenticalratesofbinaryburninganddestruction.Gaoetal.(1991),however,findashortercorecollapsetime,adeeperinitialcorecontraction,andsignificantlymoreerraticbehaviorduringthegravother-maloscillationphase.Weattributethedifferencesprimar-ilytotheirneglectofthestrongcorrelationbetweenbinaryhardnessandspatialdistributioninthecluster,aswellasfundamentaldifferencesbetweentheir1-DFokker-Planckmethodandour2-DMonteCarlomethod.OurresultsfortheinitialcorecontractionareinmuchbetteragreementwiththoseofdirectN-bodysimulations.Inaddition,wehavepresentedthefirstnumericaldemonstrationofthetheoreticallypredictedtemperatureinversionpoweringre-expansionaftercorecollapseandgravothermaloscillationsforaclusterwithprimordialbinaries.

10

¨FREGEAU,GURKAN,JOSHI,&RASIO

binaryfractions(fb󰀂0.05)toremainaround,oreven

convergeto,c≃1.5.ThisisinreasonableagreementwiththebottompanelofFig.17,whichshowstheobserveddis-tributionalsocenteredaroundc≃1.5.ThetoppanelinFig.16showsmostclustersinthebinaryburningphasewithrh/rc󰀁10,alsoinquitereasonableagreementwiththeobserveddistribution(Fig.17b),althoughtheob-servedpeakaroundrh/rc≃2wouldrequirethatmostini-tialmodelsbelesscentrallyconcentratedthanourW0=7Kingmodels.Wealsonotethatbothcandrh/rcincreasesignificantly,andsometimesdramatically,forclustersap-proachingtidaldisruption.Thusthesuggestionfromourresultsmightbethatclustersclassifiedobservationallyas“core-collapsed”arethoseinthelastfewrelaxationtimesbeforedestructionintheGalactictidalfield.MostofourKingmodelsappeartospendroughlythelast10–20%oftheirliveswithrh/rc󰀂10,orc󰀂2,againnottoodiffer-entfromtheobservedfractionof“core-collapsed”clustersinourGalaxy(seeFig.17a).Ofcourseamoreseriouscomparisonshouldtakeintoaccountrealclusteragesandthedistributionofinitialvaluesfortrh(0),whichisratheruncertain.

Wearecurrentlyintheprocessofimplementingexactdynamicalintegrationstohandlebinary–binaryinterac-tionsmoreaccuratelyinoursimulations,usingFewbody,anewsmall-Nintegratorwehavewritten.Thisintegra-torperformsautomaticclassificationofoutcomes,auto-maticstabilityanalysisofarbitrarilylargehierarchies,au-tonomousintegrationtermination,andstellarcollisions.WeareverygratefultoDouglasHeggie,SteveMcMil-lan,SimonPortegiesZwart,andSaulRappaportformanyhelpfuldiscussions.Thethree-bodyintegratorusedinthisworkispartoftheStarlabsoftwarepackagedevelopedbyPietHut,SteveMcMillanandSimonPortegiesZwart.ThisworkwassupportedbyNASAATPGrantNAG5-12044andNSFGrantAST-0206276.Someofournumeri-calsimulationswereperformedonparallelsupercomputersatBostonUniversityandNCSAunderNationalCompu-tationalScienceAllianceGrantAST980014N.

Wehavealsoconsideredmorerealistic,tidallytruncatedKingmodels.Wehavefoundthattheinitialcorecon-tractionphaseisabsentinthesesystems,orreplacedbyaninitialexpansionofthecore,andthatcorecollapseorcompletetidaldisruptionalwaysoccursinlessthan∼50trh(0).Forabinaryfraction󰀂0.1,thedisruptiontimedecreaseswithincreasingbinaryfraction,andcom-pletedisruptionoccursbeforeanydeepcorecollapse.Thepossibilityforaclustertosuffercompletedisruptionbeforecorecollapseisaqualitativelynewbehaviorintroducedbyprimordialbinaries.Ourresultsareingoodqualita-tiveagreementwithpreviousstudiesoftidallytruncatedclusterscontainingprimordialbinaries.

WehavealreadyarguedinSection3.5thatourre-sultsareingeneralagreementwithcurrentdetermina-tionsofbinaryfractionsinglobularclustercores,typically∼0.1−0.4.Wenowbrieflyconsiderourbasicpredictionsforthestructuralparametersofclustersduringthebinaryburningphaseandcomparethemtotheobservedstruc-turalparametersofglobularclusters.Whileourmodelsareclearlyfartooidealizedforanydetailedcomparison,itisusefultoexamineatleastthemostfundamentalstruc-turalparameters:thecoreradiusrc,thehalf-massradiusrh,and,fortidallytruncatedclusters,thetidalradiusrt.Sincetheoverallscaleislargelyirrelevant(althoughitcouldinprinciplebesetbyrelatingthemaximumbind-ingenergyofthebinariestoastellarradius),weconsideronlythetworatiosrh/rcandrt/rc(or,equivalently,theconcentrationparameterc≡log(rt/rc)oftenderivedbyobserversusingKingmodelfitstophotometricdata).Fig.16showstheevolutionofrh/rcandtheconcentra-tionparametercforseveralKingmodels(seeSec.3.4).Fig.17showsdistributionsofrh/rcandcforGalacticglobularclusters,withdatatakenfromthecompilationofHarris(1996).Thetoppanelshowsahistogramofrh/rcvaluesforallGalacticglobularclusters,includingthoseclassifiedobservationallyas“core-collapsed6.”Themiddleandbottompanelsshowthedistributionsofobservedval-uesforrh/rcandc,respectively,withthe“core-collapsed”clustersexcluded.First,inFig.16,notethetendencyfortheconcentrationparameterinclusterswithreasonable

REFERENCES

Albrow,M.D.,Gilliland,R.L.,Brown,T.M.,Edmonds,P.D.,Guhathakurta,P.,&Sarajedini,A.2001,ApJ,559,1060

Bacon,D.,Sigurdsson,S.,&Davies,M.B.1996,MNRAS,281,830Bailyn,C.D.1995,ARA&A,33,133

Baumgardt,H.,Heggie,D.C.,Hut,P.,&Makino,J.2003,MNRAS(inpress)

Bellazzini,M.,FusiPecci,F.,Montegriffo,P.,Messineo,M.,Monaco,L.,&Rood,R.T.2002,AJ,123,21

Breeden,J.L.,Cohn,H.N.,&Hut,P.1994,ApJ,421,195

Chernoff,D.F.&Huang,X.1996,inIAUSymp.174:DynamicalEvolutionofStarClusters:ConfrontationofTheoryandObservations,263–+

Cheung,P.,PortegiesZwart,S.,&Rasio,F.A.2003,MNRAS(submitted)

Clarke,C.J.,Bonnell,I.A.,&Hillenbrand,L.A.2000,ProtostarsandPlanetsIV,151

Cool,A.M.&Bolton,A.S.2002,inASPConf.Ser.263:StellarCollisions,MergersandtheirConsequences,163–+

Cote,P.,Pryor,C.,McClure,R.D.,Fletcher,J.M.,&Hesser,J.E.1996,AJ,112,574

Davies,M.B.&Hansen,B.M.S.1998,MNRAS,301,15

Deutsch,E.W.,Margon,B.,&Anderson,S.F.2000,ApJ,530,L21

6

Ferraro,F.R.,D’Amico,N.,Possenti,A.,Mignani,R.P.,&Paltrinieri,B.2001,ApJ,561,337

Ford,E.B.,Joshi,K.J.,Rasio,F.A.,&Zbarsky,B.2000,ApJ,528,336

Fregeau,J.M.,Joshi,K.J.,PortegiesZwart,S.F.,&Rasio,F.A.2002,ApJ,570,171

Gao,B.,Goodman,J.,Cohn,H.,&Murphy,B.1991,ApJ,370,567Gebhardt,K.,Rich,R.M.,&Ho,L.C.2002,ApJ,578,L41Giersz,M.&Spurzem,R.2000,MNRAS,317,581Goodman,J.&Hut,P.19,Nature,339,40H´enon,M.1971a,Ap&SS,13,284H´enon,M.H.1971b,Ap&SS,14,151Harris,W.E.1996,AJ,112,1487

Heggie,D.C.1975,MNRAS,173,729

Heggie,D.C.&Aarseth,S.J.1992,MNRAS,257,513Heggie,D.C.&Ramamani,N.19,MNRAS,237,757

Hurley,J.R.,Tout,C.A.,Aarseth,S.J.,&Pols,O.R.2001,MNRAS,323,630

Hut,P.1984,ApJS,55,301

Hut,P.,McMillan,S.,Goodman,J.,Mateo,M.,Phinney,E.S.,Pryor,C.,Richer,H.B.,Verbunt,F.,&Weinberg,M.1992a,PASP,104,981

Itisunclearwhattherelationis,ifany,betweenthisobservationalclassification—basedontheabsenceofagoodfittingKingmodel—andthevarioustheoreticaldefinitionsofcorecollapseusedintheliteratureondynamicalmodeling.

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.

Hut,P.,McMillan,S.,&Romani,R.W.1992b,ApJ,3,527Hut,P.,Murphy,B.W.,&Verbunt,F.1991,A&A,241,137Joshi,K.J.,Nave,C.P.,&Rasio,F.A.2001,ApJ,550,691

Joshi,K.J.,Rasio,F.A.,&PortegiesZwart,S.2000,ApJ,0,969Kim,S.S.&Lee,H.M.1999,A&A,347,123

Kim,S.S.,Lee,H.M.,&Goodman,J.1998,ApJ,495,786Kochanek,C.S.1992,ApJ,385,604

Kumar,P.&Goodman,J.1996,ApJ,466,946Lee,H.M.&Ostriker,J.P.1993,ApJ,409,617Lee,M.H.1993,ApJ,418,147Leonard,P.J.T.19,AJ,98,217

Leonard,P.J.T.&Linnell,A.P.1992,AJ,103,1928Makino,J.1996,ApJ,471,796

Makino,J.2001,inASPConf.Ser.228:DynamicsofStarClustersandtheMilkyWay,87–+

Mateo,M.,Harris,H.C.,Nemec,J.,&Olszewski,E.W.1990,AJ,100,469

McMillan,S.&Hut,P.1994,ApJ,427,793

McMillan,S.,Hut,P.,&Makino,J.1990,ApJ,362,522—.1991,ApJ,372,111

Mikkola,S.1983a,MNRAS,203,1107

11

—.1983b,MNRAS,205,733—.1984a,MNRAS,207,115—.1984b,MNRAS,208,75

PortegiesZwart,S.F.&McMillan,S.L.W.2002,ApJ,576,9PortegiesZwart,S.F.,McMillan,S.L.W.,Hut,P.,&Makino,J.2001,MNRAS,321,199

Pryor,C.P.,Latham,D.W.,&Hazen,M.L.1988,AJ,96,123Quinlan,G.D.1996,NewAstronomy,1,255

Rasio,F.A.,Pfahl,E.D.,&Rappaport,S.2000,ApJ,532,L47Rubenstein,E.P.&Bailyn,C.D.1997,ApJ,474,701Sigurdsson,S.&Phinney,E.S.1995,ApJS,99,609

Spitzer,L.1987,Dynamicalevolutionofglobularclusters(Princeton,NJ,PrincetonUniversityPress,1987,191p.)Spitzer,L.&Mathieu,R.D.1980,ApJ,241,618Spitzer,L.J.&Hart,M.H.1971,ApJ,1,399

Sugimoto,D.&Bettwieser,E.1983,MNRAS,204,19P

Takahashi,K.&PortegiesZwart,S.F.1998,ApJ,503,L49+—.2000,ApJ,535,759

Watters,W.A.,Joshi,K.J.,&Rasio,F.A.2000,ApJ,539,331Yan,L.&Mateo,M.1994,AJ,108,1810

Table1

Modelparametersandresults.

Model

fb

tcc[ortdis]/trh(0)

M(tcc)/M(0)

Mb(tcc)/Mb(0)

12

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.1.—Comparisonbetweentheuseofdirectthree-bodyintegrations(solidlines)andcrosssections(dashedlines)incalculatingthe

evolutionofaPlummermodelcontainingN=3×105starswith20%binariesinitially.Inbothcasesbinary–binaryinteractionswereturnedoff.Thetoppanelshowsthetotalmassinbinaries.Themiddlepanelshowstheenergygeneratedinbinary–singleinteractions,asafractionofthetotalinitialbindingenergyofthecluster.Thebottompanelshows(fromtoptobottom)thehalf-massradiusofsinglestars,thehalf-massradiusofbinaries,andthecoreradius,inunitsoftheinitialhalf-massradius.Timeisgiveninunitsoftheinitialhalf-massrelaxationtime.Theagreementbetweenthetwomethodsisstrong,althoughtheenergyproductionisslightlyoverestimatedinthetreatmentbasedoncrosssections,leadingtodivergentevolutionsnearcorecollapse(themodelcalculatedwithdirectthree-bodyintegrationscollapsesearlier).

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.13

Fig.2.—EvolutionofthekineticenergyT,potentialenergyWandtotalconservedenergyE,aswellasthevirialratio2T/|W|,foraPlummermodelwithN=3×105starsand10%binariesinitially.Theclusterremainsveryclosetovirialequilibriumthroughouttheintegration,andenergyconservationismaintainedtowithinafewpercent.NotethatEiscorrectedforboththeenergylostthroughevaporation,andtheenergygainedthroughbinary–binaryandbinary–singleinteractions,sothatE=T+W(exceptatt=0).Thetruetotalenergyofthecluster,T+W,increasessignificantlyovertimeasaresultoftheseinteractions.Hereweshowthequantitythatshouldbeconserved,whichwemonitor(inadditiontothevirialratio)forqualitycontrolpurposesinallourcalculations.

14

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.3.—EvolutionofanisolatedPlummermodelwithN=3×105starsand2%primordialbinariesinitially.Thetoppanelshowsthe

totalclustermassMandthetotalmassMbinbinaries.Themiddlepanelshowstheenergyreleasedthroughbinary–binaryandbinary–singleinteractions,inunitsoftheinitialbindingenergyofthecluster.Thelowerpanelshowsthecoreradiusrcofthecluster,thehalf-massradiusrh,sofsinglestars,andthehalf-massradiusrh,bofbinaries(solidlines).Forcomparison,thecoreradiusofanequivalentPlummermodelwith2%primordialbinariesbutwithallinteractionsturnedoffisalsoshown(short-dashedline).Weseethatevenaprimordialbinaryfractionassmallas2%cansignificantlydelaycorecollapse,withtccincreasingbymorethanafactor2inthiscase.Alsoshownforcomparisonandtestingisthecoreradiusofanequivalentmodelwherethebinary–singleinteractionswerecomputedwithdirectthree-bodyintegrationsinsteadofcrosssections(long-dashedline).Hereagainwenotethatthemodelbasedondirectintegrationscollapsesslightlyearlier.

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.15

Fig.4.—SameasFig.3,butforamodelwitha10%primordialbinaryfraction.Heretheenergygeneratedfrombinary–binaryinteractionsclearlydominatesthatfrombinary–singleinteractions.Weseethatanisolatedclusterwith10%binariescanbesupportedagainstcollapseforabout70trh.Severalcyclesofgravothermaloscillationspoweredbyprimordialbinariesareseenaftertheinitialcollapse.Theoscillationsinthiscaseappearquasi-periodicwithaperiodofroughly50trh(0).

16

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.5.—SameasFig.3,butforamodelwitha20%primordialbinaryfraction.Heretheenergygeneratedfrombinary–binaryinteractions

isevenmoreclearlydominant.Theclusterisinitiallysupportedagainstcollapseforabout125trh.Afterthefirstcorecollapse,gravothermaloscillationspoweredbyprimordialbinariescontinueupto∼103trh.Bythattimethetotalnumberofbinarieshasbeenreducedbyafactor∼10,buttheprimordialbinaryreservoirisstillnotexhausted.

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.17

Fig.6.—RatioofcorecollapsetimeswithandwithoutbinaryinteractionsasafunctionoftheinitialbinaryfractionfbforthePlummermodels.Thesolidlineshowsasimplelinearfit.

18

¨FREGEAU,GURKAN,JOSHI,&RASIO

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.19

Fig.8.—EvolutionofatidallytruncatedW0=7KingmodelwithN=3×105starsand2%primordialbinaries.ConventionsareasinFig.3.ComparedtoanisolatedPlummermodelwiththesamenumberofstarsandbinaries,theevolutionofthistidallytruncatedclustertocorecollapseisonlyslightlyfaster.

20

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.9.—SameasFig.8,butforamodelwitha5%primordialbinaryfractioninitially.

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.21

Fig.10.—SameasFig.8,butforamodelwitha10%primordialbinaryfractioninitially.Herecompletedisruptionoccursbeforecorecollapse.Notealsotheabsenceofanycorecontractioninitially.

22

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.11.—SameasFig.8,butforamodelwitha20%primordialbinaryfractioninitially.Notetheinitialexpansionofthecore.Here

alsocompletedisruptionoccursbeforecorecollapse.Theapparentre-expansionofthecoreradiusaftert/trh(0)≃42isanumericalartifactcausedbytheverysmallnumberofstarsleftinthecluster(ourmethodofcalculatingrcpicksupstarswellouttsidethetruecore).

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.23

Fig.12.—SameasFig.10,butforaW0=11Kingmodel(witha10%primordialbinaryfraction).Thisinitiallymuchmorecentrallyconcentratedmodelundergoesdeepcorecollapseasitrunsoutofbinariesbeforecompletedisruption.Notealsothesignificantinitialexpansionofthecoreneededtoreachquasi-equilibriuminafewtrh(0).

24

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.13.—SameasFig.10,butforaW0=3Kingmodel(witha10%primordialbinaryfraction).Thisclusterismuchlesscentrally

concentratedinitiallyandtherefore,asexpected,itisdisruptedbeforeundergoingdeepcorecollapse.Note,however,thatsignificantcorecontractionoccursthroughouttheevolution.

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.25

Fig.14.—EvolutionofthebinaryfractionfbindifferentregionsforvariousKingmodels.ThethinsolidlineisforaW0=11modelwith10%binariesandthethindottedlineisforaW0=3modelwith10%binaries.TheotherlinesareforW0=7modelswithincreasingbinaryfractionsfrombottomtotopasinFigs.8–11.Binaryfractionsashighas0.5–0.8(butmoretypically≃0.1−0.2)canbeexpectedinclustercores,whileintheouterhaloonehastypicallyfb󰀁0.1.

26

¨FREGEAU,GURKAN,JOSHI,&RASIO

MONTECARLOSIMULATIONSOFGLOBULARCLUSTEREVOLUTION.III.27

Fig.16.—Ratioofhalf-masstocoreradiusandconcentrationparameterc≡log(rt/rc)forvariousKingmodels.ThethinsolidlineisforaW0=11modelwith10%binariesandthethindottedlineisforaW0=3modelwith10%binaries.ThethicklinesareforW0=7modelswithincreasingbinaryfractionsfromtoptobottom,goingfrom2%,5%,10%,to20%.

28

¨FREGEAU,GURKAN,JOSHI,&RASIO

Fig.17.—Observeddistributionofrh/rcandconcentrationparametercforGalacticglobularclusters.Clustersclassifiedobservationally

as“core-collapsed”areincludedin(a)and(b),butnotin(c).Theobservedvaluesofthesebasicstructuralparametersforthenon-“core-collapsed”clustersareclearlyingeneralagreementwithvaluespredictedbyoursimplemodelsforclusterssupportedbyprimordialbinaryburning.