Black Hole Mergers from Star Clusters with Top-Heavy Initial Mass Functions
Newlin C. Weatherford, Giacomo Fragione, Kyle Kremer, Sourav Chatterjee, Claire S. Ye, Carl L. Rodriguez, Frederic A. Rasio
DD RAFT VERSION J ANUARY
8, 2021
Preprint typeset using L A TEX style emulateapj v. 01/23/15
BLACK HOLE MERGERS FROM STAR CLUSTERS WITH TOP-HEAVY INITIAL MASS FUNCTIONS N EWLIN
C. W
EATHERFORD , G IACOMO F RAGIONE , K YLE K REMER , S OURAV C HATTERJEE , C LAIRE
S. Y E , C ARL
L.R
ODRIGUEZ , F REDERIC
A. R
ASIO
Center for Interdisciplinary Exploration & Research in Astrophysics (CIERA), Evanston, IL 60202, USA Department of Physics & Astronomy, Northwestern University, Evanston, IL 60202, USA TAPIR, California Institute of Technology, Pasadena, CA 91125, USA The Observatories of the Carnegie Institution for Science, Pasadena, CA 91101, USA Tata Institute of Fundamental Research, Homi Bhabha Road, Mumbai 400005, India and McWilliams Center for Cosmology, Department of Physics, Carnegie Mellon University, Pittsburgh, PA 15213, USA
Draft version January 8, 2021
ABSTRACTRecent observations of globular clusters (GCs) provide evidence that the stellar initial mass function (IMF) maynot be universal, suggesting specifically that the IMF grows increasingly top-heavy with decreasing metallicityand increasing gas density. Non-canonical IMFs can greatly affect the evolution of GCs, mainly because thehigh end determines how many black holes (BHs) form. Here we compute a new set of GC models, varyingthe IMF within observational uncertainties. We find that GCs with top-heavy IMFs lose most of their masswithin a few Gyr through stellar winds and tidal stripping. Heating of the cluster through BH mass segregationgreatly enhances this process. We show that, as they approach complete dissolution, GCs with top-heavy IMFscan evolve into ‘dark clusters’ consisting of mostly BHs by mass. In addition to producing more BHs, GCswith top-heavy IMFs also produce many more binary BH (BBH) mergers. Even though these clusters areshort-lived, mergers of ejected BBHs continue at a rate comparable to, or greater than, what is found for long-lived GCs with canonical IMFs. Therefore these clusters, although they are no longer visible today, could stillcontribute significantly to the local BBH merger rate detectable by LIGO/Virgo, especially for sources withhigher component masses well into the BH mass gap. We also report that one of our GC models with a top-heavy IMF produces dozens of intermediate-mass black holes (IMBHs) with masses M >
100 M (cid:12) , includingone with M >
500 M (cid:12) . Ultimately, additional gravitational wave observations will provide strong constraintson the stellar IMF in old GCs and the formation of IMBHs at high redshift. INTRODUCTIONAmong the densest environments in the Universe, globu-lar clusters (GCs) are ideal laboratories to investigate the im-portance of stellar dynamics in forming and evolving com-pact objects and compact binaries (see, e.g., Heggie & Hut2003). Frequent dynamical encounters between cluster mem-bers are fundamental in creating and explaining the existenceof a number of exotic populations, such as X-ray binaries(e.g., Clark 1975; Verbunt et al. 1984; Giesler et al. 2018;Kremer et al. 2018), radio pulsars (e.g., Lyne et al. 1987; Sig-urdsson & Phinney 1995; Ivanova et al. 2008; Ye et al. 2019),and gravitational wave (GW) sources (e.g., Rodriguez et al.2015; Askar et al. 2017; Banerjee 2017; Fragione & Kocsis2018; Kremer et al. 2019).Many results of GC modeling rely on the assumption thatthe stellar initial mass function (IMF) has the form of a canon-ical Kroupa (2001) IMF. Observations, however, suggest theIMF may not be universal. For example, Milky Way GCs withlow central densities appear deficient in low-mass stars (DeMarchi et al. 2007), while GCs in the Andromeda galaxy ex-hibit a trend between metallicity and mass-to-light ratio thatonly a non-canonical, top-heavy IMF could explain (Haghiet al. 2017). Ultra-compact dwarf galaxies have large dy-namical mass-to-light ratios and appear to contain an over-abundance of LMXB sources (Dabringhausen et al. 2009).Our Galactic center provides an extreme example of a non-canonical IMF; stars formed there in the past few Myr havemasses consistent with a top-heavy mass function (Bartkoet al. 2010). All these observations can be explained if thestellar IMF becomes increasingly top-heavy with decreasing metallicity and increasing gas density of the pre-GC cloud(Marks et al. 2012). This is theoretically expected giventhe Jeans mass instability in molecular clouds (Larson 1998)and self-regulation of accretion onto forming stars (Adams &Fatuzzo 1996; Adams & Laughlin 1996). Metallicity plays adecisive role by regulating line-emission cooling in the col-lapsing gas cloud and radiation pressure against stellar accre-tion (De Marchi et al. 2017).The IMF has long been known to strongly impact the dy-namical evolution and survival of GCs. Already Chernoff &Weinberg (1990) found enhanced mass loss rates in clustermodels with flatter IMFs, due to increased winds from high-mass stars and faster halo expansion and evaporation. Morerecent results show that GC models with top-heavy IMFs dis-solve much faster than models with canonical IMFs (Chatter-jee et al. 2017; Giersz et al. 2019). In particular, clusters withtop-heavy IMFs produce more numerous and more massiveblack holes (BHs). Crucially, this promotes the BH-burningprocess, in which strong dynamical encounters with BHs pro-vide energy to stellar populations, inflating the cluster halo(Mackey et al. 2007, 2008; Breen & Heggie 2013; Kremeret al. 2020a). Taken to an extreme by a top-heavy IMF, thismechanism will force rapid and unstable evaporation throughthe tidal boundary and early cluster dissolution (Giersz et al.2019).The impact on cluster evolution of varying BH abundancehas recently been further analyzed via a combination of ana-lytical calculations and N -body simulations (Breen & Heggie2013; Wang 2020). However, these studies only consider ide-alized star clusters with two components (stars and BHs) andno stellar evolution, which plays a crucial role in cluster dis- a r X i v : . [ a s t r o - ph . GA ] J a n ruption (e.g., Chernoff & Weinberg 1990). Apart from thefew models with non-canonical IMFs examined by Chatterjeeet al. (2017) and Giersz et al. (2019), there is no systematicstudy of the IMF’s influence on BH populations generated viastellar evolution. In this Letter, we extend the grid of clus-ter models in the CMC Cluster Catalog (Kremer et al.2020b), focusing on the role of a varying IMF. In particu-lar, we study how non-canonical IMFs shape the evolution ofGCs, their BH populations, and the number of dynamicallyproduced BBH mergers.The Letter is organized as follows. In Section 2, we de-scribe the parameters of the numerical models we evolve. InSection 3, we analyze the structural evolution of clusters anddescribe the properties of their black hole populations, includ-ing BBH mergers. Finally, in Section 4, we discuss the impli-cations of our findings and lay out our conclusions. CLUSTER MODELSTo evolve our cluster models, we use
CMC (Cluster MonteCarlo), a Hénon-type Monte Carlo code that includes thenewest prescriptions for wind-driven mass loss, compact ob-ject formation, and pulsational-pair instabilities (see Hénon1971a,b; Joshi et al. 2000, 2001; Fregeau et al. 2003; Chat-terjee et al. 2010, 2013; Pattabiraman et al. 2013; Rodriguezet al. 2015; Kremer et al. 2020b, and references therein).We generate 15 independent cluster simulations varying theinitial number of particles (single stars plus binaries) N i andinitial total mass M i , with uniform initial cluster virial radius( r v = 1 pc), metallicity ( Z = 0 . (cid:12) ), and Galactocentric dis-tance ( R gc = 8 kpc). These choices of r v , R gc , and Z are knownto result in models that closely match typical Milky Way GCswhen using a canonical IMF (Kremer et al. 2020b).We assume that all models are initially described by a Kingprofile with concentration W = 5 (King 1962). Primary stellarmasses are sampled from the Kroupa (2001) multi-componentpower-law IMF, ξ ( m ) ∝ m − . . ≤ m / M (cid:12) ≤ . m − . . ≤ m / M (cid:12) ≤ . m − α . ≤ m / M (cid:12) ≤ . . (1)We choose three different values for α = (1 . , . , . α = 2 . f b = 5% and draw secondary masses from auniform distribution in mass ratio (e.g., Duquennoy & Mayor1991). Binary orbital periods are sampled from a log-uniformdistribution, with orbital separations ranging from near con-tact to the hard/soft boundary, while binary eccentricities arethermal (e.g., Heggie 1975). For details on stellar evolution in CMC , see Kremer et al. (2020b). We compute GW recoil kicksfor BH merger products following the methods described inRodriguez et al. (2019) and references therein. We assumeall BHs have zero natal spin while BH merger products areassigned new spins of ∼ . T H ≈ α = 1 .
6, simulation terminatesonce there are fewer than 300 particles per CPU, the default Assuming a Milky Way-like potential (e.g., Dehnen & Binney 1998). particle count used to compute local densities throughout thecluster. This typically corresponds to a few × stars in theGC, beyond the point where the assumptions of the MonteCarlo method – namely spherical symmetry and a relaxationtimescale longer than the dynamical timescale – start to breakdown (Chatterjee et al. 2017). However, the exact cut-offpoint in the simulation does not affect the general evolutionarytrends of the dissolution process. Further accurate evolutionof the remaining ‘dark star cluster’ would require a switch todirect N -body methods (e.g., Banerjee & Kroupa 2011). RESULTSIn this Section, we discuss how different IMFs affect theevolution of the cluster mass, size, and BH population. We donot find significant qualitative differences between the caseswhere we keep the initial number of particles constant (mod-els 1-3, 7-9, and 10-12 in Table 1) or the initial total massconstant (models 4-9 and 13-15).3.1.
Cluster evolution
Variations in the high end of the IMF can dramatically af-fect cluster evolution and survival, especially by setting thenumber of BHs that are formed and retained, determiningthe degree to which BH-burning regulates the cluster’s en-ergy reservoir (Chatterjee et al. 2017). In this process, BHsquickly mass-segregate into the cluster core on a sub- Gyrtimescale, forming a central BH population that undergoesfrequent phases of core collapse and re-expansion (Morscheret al. 2015). During these events, the BHs mix with the restof the cluster and provide energy to passing stars in scatter-ing interactions (Breen & Heggie 2013). Cumulatively, theseinteractions inflate the cluster halo and force faster evapora-tion through the tidal boundary (Chatterjee et al. 2017; Gierszet al. 2019; Wang 2020).In Figure 1, we show the evolution of the total mass (top-left panel) for star clusters with different N i and α . Clustersevolve very differently depending on the choice of α . Themodels with α = 2 . T H , while the modelswith α = 1 . t ∼ N i . The different fates ofthe cluster models are understandable given the mechanismsthat power cluster mass loss at different epochs:1. Initially, mass is lost primarily via stellar winds frommassive stars (e.g., Chernoff & Weinberg 1990). Thisrelatively quiet first phase is visible early in the top-leftpanel of Figure 1, terminating around t ≈ (cid:46) t (cid:46)
60 Myr, is driven primarily by Type II su-pernovae (SNe) accompanying compact object forma-tion. This phase most clearly depends on the choiceof α . For α = 3, SNe-driven mass loss from few mas-sive stars is negligible compared to ongoing stellar windmass loss among the far more numerous low-mass stars.For α = 1 .
6, however, the SNe-driven mass loss phaseis nearly catastrophic, causing the cluster to lose closeto half its original mass by the time compact object for-mation slows ( t ≈
60 Myr). Ejections of lighter com-pact objects due to natal kicks also contribute to mass While these models may not survive with significant mass beyond a fewGyr, models with larger N i , larger Galactocentric distance, and/or smallervirial radius could feasibly survive much longer. Table 1
Model initial conditions and data on cumulative BH formation and BBH mergers from 0 ≤ t ≤
13 Gyr; upper-IMF slope α , initial number of particles N i ,initial mass M i , total BHs formed N BH , total BBH mergers N (merg)BBH , ejected BBH mergers N (merg)BBH , ejected , BBH mergers with a ‘mass-gap’ BH component( M BH > . (cid:12) ) formed via stellar collisions N (merg)BBH , gap , coll , multiple-generation BBH mergers with a mass-gap BH component N (merg)BBH , gap , + , and total BBHmergers at redshifts z < N (merg)BBH ( z < lower limits since we do not consider mergers resulting from dynamics in ‘dark clusters’ (see Section 3.3).BBH mergers α N i (10 ) M i (10 M (cid:12) ) N BH N (merg)BBH N (merg)BBH , ejected N (merg)BBH , gap , coll N (merg)BBH , gap , + N (merg)BBH ( z < Figure 1.
Evolution of the total mass (top-left), half-mass radius r h (bottom-left), number of BHs N BH (top-right), and fractional number of BHs N BH / N (bottom-right) for cluster models with different values of α (blue 1 .
6, black 2 .
3, red 3 . × , dashed 8 × , solid 1 . × ). The simulations terminate at T H = 14 Gyr, except for α = 1 .
6, where we cut off the plot just before r h starts todrop sharply (around this point, assumptions in our Monte Carlo methods stop being valid; see text). In the top-left panel, the four mass loss phases discussed inthe text are highlighted via shaded intervals and the dominant mass loss mechanisms indicated. The bifurcation of the third phase is due to decreased duration ofthis phase for α = 1 . loss in this phase starting at t ≈ N BH (top-right panel).3. A ‘third phase’ of mass loss begins at the end of rapidcompact object formation. Really just a return of thefirst phase, this period of mass loss is once again domi-nated by stellar winds, this time from lower-mass stars.For the models where α = 2 . α = 1.6 models – around t ≈
200 Myr.4. The fourth and final phase of mass loss is driven byevaporation through the tidal boundary. This phasedepends strongly on the IMF due to its influence onthe number of BHs formed, and thereby the degree towhich BH-burning inflates the halo. In the case of ex-treme BH-burning, the cluster quickly expands to fillits Roche lobe within the assumed Galactic potentialand starts to tidally disrupt, as seen in the models with α = 1 .
6. Note that similar accelerating mass loss is alsoseen at late times in the models with α = 3. In thesemodels, however, the runaway process of halo evapora-tion and tidal contraction is not stimulated by extremeBH-burning. Instead, these top-light, low-mass modelssimply start out with much smaller tidal radii. Hencethey see significant evaporation rates relatively quicklycompared to canonical models ( α = 2 . r h (bottom-left panel,Figure 1) also reflects this general picture. In all models, r h expands as a consequence of stellar mass loss. Enhanced BH-burning in models with α = 1 . t ∼ α = 3 do not exhibit rapid expansion, evolving moregradually as a result of reduced massive star formation andthe lack of a significant central BH population.3.2. Black hole population
In the top-right panel of Figure 1, we show the evolutionof the number of BHs N BH . Clusters with top-heavy IMFsproduce more high-mass stars and thereby more BHs. As de-scribed above, natal kicks start to eject newly-formed BHsafter t ≈ N BH to plateau for each model. Atthis point, we find that N BH in models with α = 1 . ∼ α = 2 . ∼
100 timeshigher than in models with α = 3 (see also the cumulativenumber of BHs formed in column 4 of Table 1). During thisphase of roughly constant N BH , dynamical friction causes BHsto segregate to the cluster center. This phase’s duration de-pends on α since clusters with top-heavy IMFs have higheraverage stellar mass and correspondingly longer segregationtimescales. Once the BHs have segregated to the core, theyare gradually ejected via strong dynamical encounters.In the bottom-right panel of Figure 1, we also show the evo-lution of the fractional number of BHs N BH / N . In modelswith α = 2 . N BH / N decreases in time as strong en-counters preferentially eject massive objects (e.g., BHs) fromthe cluster core. In models with α = 1 .
6, however, N BH / N increases in time due to earlier tidal evaporation. By the mo-ment these clusters retain only ∼
20% of their original mass,
Figure 2.
Evolution of the number of binary black holes N BBH bound to thestar clusters with different values of α (blue 1 .
6, black 2 .
3, red 3 .
0) and dif-ferent initial number of particles ( N i ). Different line styles represent different N i , as in Figure 1. Since BBHs continuously form and disrupt in central scat-tering interactions, we apply a rolling average over every 10 time-steps tosmooth the curves at low N BBH . N BH / N ≈ r h peaks). Note that thenumber and mass fraction of BHs at late times depends onthe cluster dissolution process beyond this point, which ourMonte Carlo methods are not designed to address. It is nev-ertheless clear that clusters with top-heavy IMFs can evolvethrough a ‘dark cluster’ stage during which their total mass isdominated by BHs (e.g., Banerjee & Kroupa 2011).Regardless of age, individual BHs are only directly de-tectable if they reside in binaries, either via detection of GWsfrom BBH mergers or through the presence of non-BH (nBH)companions. While dense star clusters are expected to be effi-cient factories of BBH and BH–nBH binaries, their numberstypically remain small as their dynamical assembly competeswith their disruption and ejection (e.g., Downing et al. 2010,2011; Morscher et al. 2015; Kremer et al. 2018). However,the choice of IMF greatly impacts BH binary formation.In Figure 2, we plot the evolution of the number of BBHs N BBH . Clusters with top-heavy IMFs produce more BBHsvia enhanced formation of high-mass stars and BHs. N BBH plateaus a little later than N BH (top-right panel, Figure 1), a de-lay reflecting the dynamical assembly of BBHs. An equilib-rium between N BBH and N BH is apparent during this plateau,with of order one BBH for every 100 bound BHs in the clus-ter. Subsequently, N BBH decreases more rapidly for clusterswith higher N BH , as evident by the slope past the plateau.This is unsurprising since BH-burning, enhanced for clusterswith higher N BH , is characterized by periodic collapse and re-expansion of the central BH population, a process that dis-rupts and ejects many BBHs (Morscher et al. 2015; Chatterjeeet al. 2017). 3.3. BBH merger rates
We report in Table 1, column 5, the cumulative numbersof BBH mergers N (merg)BBH in each model through 13 Gyr. Itis readily apparent that each reduction in α generally in-creases N (merg)BBH , with some nuance introduced depending onwhich values of α are compared and whether N i or M i is held Figure 3.
Cumulative binary black hole mergers with respect to redshift (andtime) for star clusters with different values of α (blue 1 .
6, black 2 .
3, red 3 . N i ). Different line styles representdifferent N i , as in Figure 1. All model GCs were assumed to be born 13 Gyrago and appropriate redshifts were then computed using Astropy ’s (As-tropy Collaboration et al. 2013) cosmology calculator under the flat Λ CDMmodel with H = 69 . / s / Mpc and Ω matter = 0 .
286 (see, e.g., Bennett et al.2014). Values for α = 1 . constant in the comparison. Most notably, N (merg)BBH in short-lived clusters with α = 1 . α = 3. This islargely due to the sheer number of merging BBHs they eject(Table 1, column 6).Even though clusters with top-heavy IMFs rapidly dissolve,the large number of merging BBHs they eject suggests theycould contribute significantly to the total BBH merger ratein the local Universe. Indeed, we show in Figure 3 that thedistribution of BBH mergers across redshift z is roughly com-parable between the models with α = 1 . .
3, even for z <
1. The number of BBH mergers in each model at z < α = 1 . α = 2 .
3. Thiscomparison is quite conservative since the merger counts forthe models with α = 1 . lower limits . First, we only con-sider BBHs ejected prior to the simulation’s end, ignoring thepotentially substantial contribution to the z < z < α = 1 . α = x at a given redshift z canbe approximated as Γ x ( z ) ≈ (cid:104) N (merg)BBH / M i (cid:105) ρ SF ( z ) f SF f x . Here, (cid:104) N (merg)BBH / M i (cid:105) is the mean number of BBH mergers per initial cluster mass, ρ SF is the cosmological density of the star for-mation rate, f SF is the fraction of the star formation rate as-sumed to occur in star clusters, and f x is the fraction of clus-ters born with α ≈ x . To mitigate uncertainties in the latterthree terms, we simply compute (cid:104) N (merg)BBH / M i (cid:105) for each valueof α and express the estimated merger rates from clusterswith non-canonical IMFs as ratios with respect to the better-studied rates from clusters with canonical IMFs.To compute (cid:104) N (merg)BBH / M i (cid:105) for each α , we extract the func-tional dependence of N (merg)BBH on M i from columns 3 and 5of Table 1. Though the number of BBH mergers per clus-ter scales roughly linearly with the present-day cluster mass,at least in clusters with canonical IMFs (e.g., Rodriguez et al.2015; Kremer et al. 2020b), we find a more general power lawof the form N (merg)BBH = a ( M i ) b fits the data slightly better given χ goodness-of-fit tests. Specifically, N (merg)BBH ≥ (9 . ± . (cid:16) M i M (cid:12) (cid:17) . ± . α = 1 .
6= (11 . ± . (cid:16) M i M (cid:12) (cid:17) . ± . α = 2 .
3= (1 . ± . (cid:16) M i M (cid:12) (cid:17) . ± . α = 3 . , (2)where the 1 σ uncertainties on the fit parameters are computedassuming Poisson uncertainties on N (merg)BBH and the ≥ symbolindicates a lower limit for top-heavy IMFs.Assuming the distribution of GC birth masses takes theform dN cluster / dM i ∝ M − i (Lada & Lada 2003), then (cid:28) N (merg)BBH M i (cid:29) = a (10 M (cid:12) ) b (cid:82) M H M L M b − i dM i (cid:82) M H M L M − i dM i = aM H ln( M H / M L )(10 M (cid:12) ) ( M H / M L − b = 2 aM b − L (2 − b )(10 M (cid:12) ) b (cid:104) − ( M L / M H ) − b − M L / M H (cid:105) b (cid:54) = 2 , (3)where M L and M H are the lower and upper limits of themass function for clusters capable of producing BBH merg-ers. Both bounds are somewhat arbitrary, but the lower boundhas a greater impact on the merger rate calculation since theintegrands in Equation 3 scale inversely with M i . To avoid ex-trapolating our fit to cluster masses that rarely produce BBHs,we set the lower bound for each α as the M i that producesan average of two stars with M >
25 M (cid:12) , given the assumedIMF. This is roughly the minimum cluster mass needed toproduce a single BBH merger since most progenitors with M >
25 M (cid:12) will collapse to a BH. Under this definition, M L / (100M (cid:12) ) ≈ (2 , ,
80) for α = (1 . , . , α inprinciple, observations do not well-constrain this value, so wenaïvely set the upper bound for all α to be M H = 10 M (cid:12) (Harris et al. 2014). With these assumptions, we find that (cid:104) N (merg)BBH / M i (cid:105) ≈ (59 , , / (10 M (cid:12) ), respectively. The ratesfor GCs with α = 1 . α = 2 .
3, are then roughly given by Γ . Γ . (cid:38) . f GC , . f GC , . , Γ . Γ . ≈ . f GC , . f GC , . . (4)Note again that the above ratio for clusters with α = 1 . M i and α nor any cumu-lative dependency of these ratios on the choices of initial pa- Figure 4.
Cumulative density functions (CDFs) of the component masses forbinary black hole mergers in star clusters with different values of α (blue 1 . .
3, red 3 .
0) and different initial number of particles ( N i ). Different linestyles represent different N i , as in Figure 1. rameters, such as virial radius, Galactocentric distance, metal-licity, and binary fraction. Furthermore, f x is notoriously un-certain and could depend on redshift, metallicity, and timeof birth. Nevertheless, these estimates qualitatively suggestthat if f . > f . (cid:38) f . , then the contribution to the cosmo-logical BBH merger rate from clusters could be a bit higherthan currently estimated (e.g., by Kremer et al. 2020b). Inturn, if f . > f . (cid:38) f . , then those rate estimates could betoo high already. Overall, uncertainties in the IMF may con-tribute about an order of magnitude to the uncertainty on thedynamical BBH merger rate (see also Chatterjee et al. 2017).Finally, it is worth commenting briefly on black hole-neutron star (BH-NS) and binary neutron star (BNS) mergers.In clusters that retain a BH population, NSs are pushed outof the dense, BH-filled core due to mass segregation. Hencethe merger likelihood is significantly smaller for BH-NS andBNS mergers than for BBH mergers in most GCs (Ye et al.2020). It is therefore unsurprising that all our models produceinsignificant numbers of BH-NS and BNS mergers. Whileour models with α = 3 – and therefore very few formed BHs– do exhibit about four times more frequent binary NS merg-ers than the models with α = 2 .
3, the numbers are still far toolow to account for the LIGO/Virgo estimated rates (Ye et al.2020; Fragione & Banerjee 2020).3.4.
BH and BBH merger masses
The average stellar mass is higher in clusters with top-heavyIMFs and lower in clusters with top-light IMFs. So, it isunsurprising that the average masses of BHs, BBH mergercomponents, and BBH merger products in our models all in-crease with decreasing α . Averaging across all our modelsfor each α value, we find a mean BH mass at formationof (16 , ,
9) M (cid:12) for α = (1 . , . , , ,
21) M (cid:12) and (74 , ,
40) M (cid:12) . Thistrend is further exhibited in Figure 4 showing the mass distri-bution of BBH merger components. The vertical jump in thedistributions at M = 40 . (cid:12) is due to the pile-up at the start of the ‘mass gap’ of BHs formed from massive stars via pul-sational pair-instability supernovae (PPISNe). CMC assumesthis gap begins at M > . (cid:12) (for details, see Belczynskiet al. 2016), extending to around 120 M (cid:12) . Hence, the first ver-tical jump indicates that clusters with top-heavy IMFs formsignificantly more BHs via PPISNe than clusters with canoni-cal IMFs. A second, smaller jump in the distributions for top-heavy IMFs at M ≈
77 M (cid:12) arises from 2nd-generation merg-ers with a component produced in an earlier merger betweentwo such PPISNe-generated BHs.The anti-correlation between the average mass of BBHmerger components and α arises for three reasons. First,clusters with top-heavy IMFs produce higher-mass stars, in-creasing not only the number of BHs formed, but also theiraverage mass. Second, for initial N held constant, clusterswith top-heavy IMFs also have greater total mass and highercollision rates in their densely populated cores. This collisionrate enhancement is obvious in Table 1, column 7, listing thenumber of mass gap BBH merger components formed fromthe stellar product of collisions (e.g., Di Carlo et al. 2020;Kremer et al. 2020c). While the column only lists collision-ally formed merger components in the mass gap, the overallcollision rate scales similarly with decreasing α . Finally, theincreased merger frequency in clusters with top-heavy IMFsitself increases the average mass of the components by in-creasing the chances that they will have already experienced amerger (e.g., Fragione et al. 2020a,b; Rodriguez et al. 2020).Indeed, the cumulative number of 2nd-generation or higher(2G+) BBH mergers is also enhanced in clusters with top-heavy IMFs, as seen in Table 1, column 8, listing the numberof 2G+ mergers with a component in the mass gap. For clus-ters with top-heavy IMFs, this is identical to the total numberof 2G+ BBH mergers, though 2G+ mergers with neither com-ponent in the mass gap are common for higher α . Overall,2G+ merger totals roughly double with each decrease in α from 3 to 2 . .
6, while the numbers of 2G+ mergers witha mass gap component roughly quadruple.Though we find no trend between initial N and merger com-ponent or product mass for α = 2 . α = 1 .
6. In this case, for N i = (4 , , × , theaverage BBH merger product mass is (54 , ,
79) M (cid:12) , respec-tively. This rise is in part due to 2G+ mergers, which accountfor (7% , , N are especially good at efficiently producingmany 2G+ mergers and many BHs in the mass gap.The increased collision and 2G+ merger rates in clusterswith top-heavy IMFs can also result in the formation of‘intermediate-mass’ BHs (IMBHs), which we define as BHswith masses exceeding 100 M (cid:12) . In particular, our highest- N model with α = 1 . M = 537 M (cid:12) . This particular IMBHformed in the merger of a 122 M (cid:12) IMBH with a 426 M (cid:12)
IMBH, which itself formed in the merger of two IMBHs with M ≈
200 M (cid:12) . While we intend to explore IMBH-formation inGCs with top-heavy IMFs more thoroughly in future work, wefor now direct the reader’s attention to our collaboration’s re-cent study on IMBH formation in GCs (González et al. 2020). DISCUSSION AND CONCLUSIONSObservations provide evidence that the stellar IMF may notbe universal. A non-canonical IMF can greatly affect a starcluster’s dynamical evolution, especially since its high end de-termines how many BHs form within, regulating the cluster’senergy budget and dynamical clock (Breen & Heggie 2013;Chatterjee et al. 2017; Kremer et al. 2020a; Wang 2020). Inthis Letter, we have extended the
CMC Cluster Catalog (Kremer et al. 2020b) to examine how a varying IMF affectsthe evolution of star clusters and their BH populations.We have shown that massive star clusters with top-heavyIMFs (low α ) are likely to lose most of their mass within afew Gyr, assuming they have low-to-average mass and Galac-tocentric distance for typical Milky Way GCs. The rapid massloss during dissolution occurs in stages, first driven by stel-lar winds and dynamical ejections, then by evaporation of thehalo through the tidal boundary. Extensive BH-burning en-hances the latter stage in clusters with top-heavy IMFs, whichproduce many BHs. Such clusters evolve through a pointwhere they consist mostly of BHs by mass (and up to at least3% by number). Further study with direct N -body methods(e.g., Banerjee & Kroupa 2011) is required to fully understandthe evolution of these clusters. Initializing direct N -body sim-ulations with the pre-dissolution states of CMC models couldbe especially useful in such an evaluation across the full clus-ter mass distribution.We note that the processes described above are also af-fected by the choices of initial cluster metallicity and natalkick distribution (Chatterjee et al. 2017). BHs in metal-richclusters have lower mass and do not inject as much energyinto the BH-burning process as BHs in metal-poor clusters.Thus, metal-rich clusters typically have higher densities anddispersion velocities, therefore processing the BH populationon shorter timescales and disrupting more binaries. Mean-while, high natal kicks will eject most BHs from the clus-ter during formation. In such a case, BHs and BBHs areexpected to be small in number regardless of the IMF. Top-heavy clusters with higher natal kicks also live longer; whilethey experience less tidal mass loss driven by BH-burning,these clusters lose even more mass due to kick-driven dynam-ical ejections of BHs (Chatterjee et al. 2017). For a recentstudy of top-heavy cluster evolution featuring different natalkick assumptions, see Haghi et al. (2020, published during re-view of this paper), who examine lower-mass clusters with di-rect N − body simulations incorporating gas expulsion physics.Notably, they find that top-heavy clusters, albeit drasticallyreduced in mass, may well survive to the present day if bornwith masses above M i (cid:38) × M (cid:12) . For times not too closeto disruption, where CMC ’s assumptions start to be challenged(see Section 2), our results are encouragingly compatible withthose of Haghi et al. (2020).Regardless of their long-term evolution and stability, wehave also shown that clusters with top-heavy IMFs – and cor-respondingly high BH production – may contribute signifi-cantly to the present-day binary BH merger rate. Even thoughthese clusters rapidly lose most of their mass within a few Gyr,mergers from ejected BBHs continue to contribute at latertimes (Fragione & Kocsis 2018), at rates comparable to orgreater than those for clusters with canonical Kroupa IMFs.The rate of 2nd-generation mergers with component massesin the mass gap may be especially enhanced in top-heavyGCs, motivating the existence of more GW190521-like merg-ers (Abbott et al. 2020a). In addition, the enhancement of col- Just like the fictional Maw cluster (Anderson 1994), this ‘dark cluster’stage of a GC born with a top-heavy IMF may be short-lived due to the rapidpace of tidal evaporation. lision rates and multiple-generation mergers in top-heavy GCsmay also lead to the formation of IMBHs and even IMBH-IMBH mergers, as demonstrated in one of our models.In general, we have shown that the high-mass slope of thecluster birth IMF may significantly impact the exact contri-bution to the cosmological BBH merger rate due to clusterdynamics. Specifically, if a large fraction of clusters wereborn with top-heavy IMFs, the cluster-dynamics merger ratemay be somewhat enhanced relative to recent estimates (e.g.,Kremer et al. 2020b). In turn, if a large fraction of clusterswere born with top-light IMFs, the cluster-dynamics mergerrate may be significantly reduced. With future observations ofgravitational waves providing unique information on the BBHmerger rate (Abbott et al. 2020b), it may be possible to lever-age this understanding to better constrain the IMFs of old starclusters. ACKNOWLEDGEMENTSThis work was supported by NSF grant AST-1716762 andthrough the computational resources and staff contributionsprovided for the Quest high-performance computing facilityat Northwestern University. NW acknowledges support fromthe CIERA Riedel Family Graduate Fellowship as well as theNSF GK-12 Fellowship Program under Grant DGE-0948017.GF acknowledges support from a CIERA Fellowship. KKis supported by an NSF Astronomy and Astrophysics Post-doctoral Fellowship under award AST-2001751. SC acknowl-edges support from the Department of Atomic Energy, Gov-ernment of India, under Project No. 12-R&D-TFR-5.02-0200.REFERENCES