Milky Way Globular Clusters: Nurseries for Dynamically-Formed Binary Black Holes

Abstract: We present a novel self-consistent theoretical framework to characterize the formation, evolution, and merger sites of dynamically-formed black hole binaries, with a focus on explaining the most massive events observed by the LIGO-Virgo-KAGRA Collaboration. Our approach couples the galaxy formation model GAMESH with cluster population synthesis codes to trace the cosmic evolution of globular clusters simultaneously with mergers of massive black holes. Our reference model, which includes prescriptions for both cluster formation and disruption depending on properties of specific galaxies, accurately reproduces the observed age-mass distribution of the Milky Way globular clusters. We find that approximately 30% of the globular clusters observed in our galaxy's halo may have originated from satellite galaxies of the Milky Way. We confirm that hierarchical black hole mergers provide a significant contribution to the formation of black holes in and above the pair-instability mass gap. However, quantifying their contribution is challenging, as different population synthesis codes yield divergent results in terms of black hole mass function and merger rates. Furthermore, we characterize the host galaxies where massive black holes form in terms of their dark matter, stellar mass, and metallicity. Ultimately, we demonstrate that the merger and birth rate densities of binary black holes increase with redshift till z = 5. This cosmic evolution is a crucial signature with significant implications for future detectors like the LISA, the Einstein Telescope and Cosmic Explorer, which will be capable to probe the high-redshift Universe.

• The paper presents a simulation framework that integrates galaxy formation and cluster evolution models to explain high-mass binary black hole mergers in Milky Way globular clusters.
• It employs calibrated CPS codes alongside cosmological models to analyze hierarchical mergers, cluster disruption, and their effects on BBH mass distributions and merger rate densities.
• Moreover, the study demonstrates that high-density, high-mass clusters drive BBH and intermediate-mass black hole formation, offering insights for next-generation gravitational-wave observations.

Milky Way Globular Clusters as Factories of Dynamically-Formed Binary Black Holes

Introduction

This paper develops a comprehensive, physically motivated simulation framework combining the GAMESH galaxy formation model with cluster population synthesis (CPS) codes (FASTCLUSTER and RAPSTER) to characterize the formation, assembly, and dynamical evolution of globular clusters (GCs) and their binary black hole (BBH) populations within a Local Group-like (LG-like) cosmological context. The primary aim is to address the origin and demographics of high-mass BBH mergers, particularly those systems challenging single-star evolutionary models, in light of LIGO-Virgo-KAGRA gravitational-wave detections.

Simulation Methodology

The study achieves an internally consistent modelling pipeline by integrating cosmological structure formation (GAMESH) with prescriptions for GC formation (with self-consistent gas, star formation, and feedback) and semi-analytic cluster population evolution. The approach leverages MW-like galaxy assembly histories to synthesize the observed MW GC system—reproducing properties such as the age-mass distribution and enabling analysis of accreted versus in situ GC demographics.

GCs are initialized based on galactic properties (gas surface density and star formation surface density), with cluster mass sampling drawn from a power-law initial cluster mass function. Detailed prescriptions for cluster disruption include both rapid, environment-dependent tidal dissolution (cruel cradle effect) within galactic disks and merger-induced GC migration into galactic halos. This is necessary to reproduce the observed MW halo GC population.

CPS codes (RAPSTER and FASTCLUSTER) are calibrated to consistent initial conditions and used to model dynamical BBH production, including the hierarchical growth of black holes by repeated mergers, mass segregation, and ejection mechanisms. The evolution is parameterized by relevant environmental (density, metallicity) and internal (binary fraction, mass function) quantities.

GC System Properties

In the model, GC formation is strongly peaked at high $z$ (main peak $z\sim 4$--5 for this LG-like, overdense volume), consistent with rapid early assembly of the MW and significant merger activity. While most GCs are initially formed at $z>2$, only a subset survives to $z=0$, shaped by the balance of rapid disruption in dense gas and preservation during mergers that relocate GCs into the MW halo.

Comparison with MW GC observations indicates fidelity in age distribution, particularly when accounting for uncertainties in observed GC ages. However, the model systematically overpredicts metallicity of GCs compared to canonical observations—this is attributed to differing metallicity definitions, possible overenrichment in the simulated interstellar medium, and lack of hydrodynamical mixing and cold accretion dilution in the GAMESH model.

Accreted GCs—those originally produced in MW satellites and delivered to the halo via mergers—comprise approximately 30% of the total surviving population, quantifying the contribution of ex situ assembly to MW GC demographics.

Figure 1: Mass-metallicity distribution of surviving MW GCs. The color bar measures central stellar number density; the red dashed histogram denotes observed MW GC masses.

Dynamical BBH Formation Pathways and GW Populations

The simulation tracks the BBH population produced in GCs throughout MW assembly. Two main features characterize the distribution of BBH masses:

• The primary black hole mass distribution is typically bimodal due to hierarchical mergers—in addition to a standard "first generation" peak around $35\,M_\odot$, dynamically-assembled systems produce a pronounced secondary peak at $70\!-\!90\,M_\odot$, extending into and above the expected pair-instability mass gap.

  Figure 2: Primary ($M_1$) and secondary ($M_2$) BH mass distributions for observable dynamically-formed BBHs; red points indicate the most massive GW events' medians and credible intervals.

• The merger rate density of BBHs increases with redshift up to $z\sim5$, with the high-mass bins (MBBHs, IMBH progenitors) peaking earlier than the low-mass population. Notably, a tail of systems with $M_1 > 1000\,M_\odot$ (IMBH regime) forms but remains largely undetectable at low $z$ with current-generation detectors.

BBH Birth Environments and Cluster Host Properties

GCs that form the most massive BBHs and IMBHs are highly biased toward high total mass and high central stellar density, with representative values $>10^6\,M_\odot$ and $>10^7\,\mathrm{pc}^{-3}$, respectively (see Figure 1). A key finding is the absence of a simple one-to-one correlation between metallicity and massive BBH formation; while more metal-poor clusters dominate the early Universe, massive, metal-enriched GCs are fully capable of producing high-mass hierarchical merger products, reflecting the physical dependence on cluster mass and density rather than metallicity alone.

Implications for Gravitational-Wave Astronomy

The results underscore the critical dependence of predicted BBH mass function and merger rate density on the physical prescriptions enforced in the CPS codes—particularly the treatment of initial binary fraction, hierarchical merger chain efficiency, and cluster disruption physics. The difference between the predictions of RAPSTER and FASTCLUSTER for the fraction of IMBHs and very high mass BBHs is substantial, even with consistent calibration, highlighting the imperative for more consistent treatments of cluster internal dynamics, massive binary retention, and ejection.

The model provides a physically motivated framework to interpret the apparent lack of a pair-instability mass gap “cutoff” in GW observations and predicts a non-negligible rate of dynamically-formed BBH mergers at $z\gtrsim3$ that will only be accessible to next-generation GW facilities (e.g., LISA, Einstein Telescope, Cosmic Explorer).

Conclusion

This work establishes that the formation of high-mass BBHs and IMBHs above $50\,M_\odot$ is strongly favoured in the densest, most massive GCs, which are in turn produced in the largest, most well-evolved dark matter halos of MW-like galaxies and their progenitors. The results demonstrate:

• Hierarchical merger processes in realistic GC populations can robustly account for GW events in and above the pair-instability mass gap, in contrast to canonical isolated binary evolution models.
• High-mass BBH/IMBH birth does not require low metallicity but does demand high-density, high-mass cluster environments.
• The shape of the local and cosmic BBH merger rate density (and the presence of observable IMBHs in the MW Halo) is acutely sensitive to details of cluster evolutionary physics and CPS methodology.
• Upcoming GW facilities targeting high $z$ and lower frequency bands will uniquely probe the early build-up of the GC and BBH population, constraining the full dynamical assembly of BBHs predicted by such models.

Future directions for this paradigm include integration with more sophisticated hydrodynamical simulations to resolve chemical inhomogeneity, cluster birth conditions, and further refinement of BBH evolutionary pathways in the hierarchical merging and galactic environment framework.

---

For direct comparisons, see GWTC catalogs [Abbott+2019, Abbott+2021, Abbott+2023, Abbott+2024]