Publications

My full publication list is available on:


Paper highlights

Click any entry to expand it.

Broekgaarden 2026 — Lower Your Rates: On Claims of a Binary Black Hole Merger-Rate Crisis Jun 29, 2026

Floor S. Broekgaarden
arXiv:2606.28515

Explore the data and results yourself with interactive figures and tables: [interactive figures]

It has recently been suggested that isolated binary evolution models systematically overpredict the binary black hole (BBH) merger rate, and that the only way to bring simulations down to — or below — the rate observed by LIGO, Virgo, and KAGRA (GWTC-5) is to invoke high-metallicity star formation or high black hole natal kicks. This paper asks whether that is actually true. It is not.

To test this, I updated and greatly expanded the 2022 rates review, compiling metadata for each simulation and submodel from isolated binary evolution that reports a BBH merger rate — 1,490 simulated merger rates drawn from 57 studies. While most models do exceed the observed rate, a meaningful subset aligns with or falls below it, and crucially those low rates arise from many different, physically justified parameter variations — not just high metallicity or strong kicks. Simulated rates depend heavily on assumptions about natal kicks, common-envelope evolution, mass transfer, and star formation history, and low BBH rates can be reached along multiple independent routes through that parameter space.

Figure 1: Compilation of simulated local BBH merger rates from isolated binary evolution studies compared to the GWTC-5 observational constraint
Figure 1. Compilation of simulated local BBH merger rates from isolated binary evolution studies compared to the observational constraints from GWTC-5 (and earlier LVK catalogs). Each horizontal line is the BBH rate range from one study, with scatter points showing individual submodel rates. While many simulations find rates above the observed LVK range (gray), a substantial number reproduce or underpredict the observed local BBH rate (green). The simulated rates span several orders of magnitude, highlighting how strongly they depend on assumptions about massive-binary evolution and cosmic star formation history — and demonstrating that claims of a universal tension between isolated binary evolution and the observed BBH rate are highly model dependent.

A central finding is the appearance of simulation silos: individual computational frameworks show internal consistency that does not extend across different codes, so trends found within one framework do not straightforwardly generalize. The overall conclusion is that claims of a universal BBH merger-rate "crisis" are strongly model dependent, and that much broader parameter-space exploration is needed before declaring genuine tension between binary evolution theory and gravitational-wave data.

Figure 3: Simulated BBH merger rates grouped by population-synthesis framework, illustrating simulation silos
Figure 3. Simulated BBH merger rates grouped by population-synthesis framework rather than by parameter family, shown as cumulative distributions (top) and per-code distributions (bottom), with the LVK interval marked in gray. Different codes populate somewhat different regions of merger-rate space, reflecting differences in their underlying implementations of stellar and binary evolution, even as substantial overlap remains. This clustering by framework illustrates the simulation silos effect: conclusions drawn from a single population-synthesis code may not generalize across the broader landscape of isolated-binary evolution models.
Broekgaarden et al. 2026 — How Common Are Common Envelopes? Quantifying Their Role in Forming Gravitational-Wave Sources Jun 4, 2026

Excited that my paper is out! This was a team effort: Floor S. Broekgaarden, Ana Lam, Sasha Levina, Jakub Klencki, Kyle A. Rocha, Lieke van Son, Steffani M. Grondin, Monica Gallegos-Garcia, Brian D. Metzger, Enrico Ramirez-Ruiz, Angela Twum, Melanie Santiago, Julia Haynes, Tyler B. Smith, Amedeo Romagnolo, Edo Berger, and Lucas M. de Sá
arXiv:2606.05322

All data and results are publicly available with interactive figures and tables you can explore yourself: [interactive catalog]

One of the most consequential — and most uncertain — phases in the life of a binary star system is common-envelope (CE) evolution. When one star grows large and engulfs its companion, the two stars spiral together inside a shared gas envelope. If the system survives, the orbit shrinks dramatically, setting up the binary to eventually merge as a pair of compact objects detectable by LIGO, Virgo, and KAGRA. If it doesn't survive, the binary is destroyed. CE evolution is therefore a critical bottleneck in forming gravitational-wave sources — yet after decades of study, we still don't fully understand it.

This paper asks a deceptively simple question: is CE evolution actually required to form the binary black holes, neutron stars, and black hole–neutron star systems that gravitational-wave observatories detect? Rather than picking a single simulation and reporting its answer, we compiled and compared predictions from over 200 population-synthesis simulations spanning many different codes and physical assumptions, and asked what fraction of merging systems in each simulation formed through a CE phase.

Figure 2: Global overview of CE fraction diversity across 200+ population-synthesis simulations for BBH, BHNS, and BNS
Figure 2. Global overview of Level 1 formation-channel diversity across 200+ compiled population-synthesis simulations. Each point shows the fraction of systems forming with CE (x-axis) vs. without CE (y-axis) for binary black holes (BBH, 117 simulations), BH–neutron star systems (BHNS, 71 simulations), and binary neutron stars (BNS, 53 simulations). BBH and BHNS populations scatter across nearly the full allowed range — from CE-dominated to entirely CE-free — while BNS mergers cluster strongly toward the CE side. See the interactive catalog to explore individual simulations.

The headline result (Figure 2) is striking: for binary black holes and BH–neutron star systems, the simulations span nearly the full allowed range — from models where essentially every merger formed through a CE phase to models where none did — yet these same models all produce similar merger rates consistent with current LVK observations. This reveals a fundamental degeneracy: merger rates alone cannot tell us how these systems formed. We can match the data with completely opposite physical assumptions about CE evolution.

Binary neutron stars tell a very different story. Across essentially all models, BNS mergers form almost exclusively through channels involving at least one CE phase (≳90–100%). This makes BNS systems a much more powerful probe of CE physics — if we can measure their merger rates and mass distributions precisely enough, we can potentially break the degeneracy that plagues the BBH population.

Figure 3: Fractional CE contributions per simulation for the BBH population
Figure 3. Fractional contributions of formation channels (with CE in blue/orange, without CE in orange, other channels in lighter colors) for every BBH simulation in the compiled catalog, sorted by CE fraction. The right panel shows the corresponding intrinsic BBH merger rate R₀ for each simulation. The enormous spread in CE fractions — from 0% to ~100% — is immediately visible, and critically, this spread is not correlated with the merger rate: simulations at both extremes of CE fraction produce rates consistent with LVK observations.

Figure 3 makes the diversity explicit: sorted by CE fraction, the simulations fill essentially every possible value from 0 to 100%, and the accompanying merger rate (right panel) shows no strong trend. You can predict almost any CE fraction and still match the observed BBH rate. The diversity is driven by a handful of key assumptions — mass-transfer stability, angular-momentum loss, CE efficiency, and supernova natal kicks — but these parameters interact non-linearly and in correlated ways across different simulation codes and frameworks. This means that trends found by varying one parameter at a time within a single code don't straightforwardly generalize to the full landscape of models.

Figure 9: BBH CE fraction vs merger rate across simulation frameworks
Figure 9. Top panel: fraction of BBH mergers forming without a CE phase vs. total intrinsic BBH merger rate, for every simulation in the catalog colored by simulation framework. Within any single framework, varying parameters can shift both the CE fraction and the merger rate, but across the full compilation there is no universal correlation: CE-free and CE-dominated models alike span the full range of observed-compatible merger rates. The bottom panels show pairwise comparisons of key model parameters, illustrating the highly non-linear and framework-dependent nature of the relationships.

The bottom line is that we need to go beyond merger-rate measurements alone. Multi-messenger observations, mass and spin distributions, delay-time distributions, and eventually electromagnetic counterparts will all be needed to break the degeneracy and constrain what CE evolution actually looks like. The interactive catalog accompanying this paper is designed to make it easy to explore the full simulation landscape and identify which observables are most diagnostic. There is much more work to do here — but this paper is intended as a community resource to anchor those future efforts.

Smith et al. 2026 — Massquerade: Impacts of Mass Ratio Reversals on Binary Black Hole Merger Rates and Mass Distributions May 20, 2026

Tyler B. Smith, Floor Broekgaarden, Sasha Levina, Amedeo Romagnolo, Manasvini Komandur, Melanie Santiago, Kyle A. Rocha
arXiv:2605.21580

When two massive stars are born together in a binary system, the heavier star is expected to also produce the heavier black hole when it dies. But binary evolution can flip this expectation. Through a process called mass ratio reversal (MRR), the initially less massive star can end up forming the more massive black hole. This happens when mass is transferred between the stars during their lives — the lighter star receives material from its companion, grows a more massive core, and ultimately collapses into a heavier black hole than the one formed by the originally dominant star. The initially less massive star is then "massquerading" as the primary.

This paper investigates how common MRR is, and — crucially — how it shapes the black hole merger rate and mass distributions that gravitational-wave observatories like LIGO, Virgo, and KAGRA actually measure. The team uses two independent binary population synthesis codes, COMPAS and SEVN, to simulate large populations of binary stars from birth through to black hole merger, and compares the resulting distributions against current LVK observations.

A key finding is that the two codes make qualitatively different predictions. In COMPAS, MRR systems dominate the high-mass end of the distribution (primary masses above ~20 M☉, secondary masses above ~12 M☉), while in SEVN the MRR contribution is more diffuse and remains subdominant across the full mass range. Despite this difference, both codes agree that MRR systems preferentially populate the high mass-ratio regime (q ≳ 0.6, meaning the two black holes have similar masses). The upshot is that the observed mass distribution cannot simply be read as a direct map of the original stellar masses — MRR blurs that connection.

The paper also identifies three distinct evolutionary channels that produce MRR (Figure 6 below): (1) core growth, where stable mass transfer fattens the secondary's helium core until it collapses into the heavier black hole; (2) PPISN shrinkage, where the primary loses mass through violent pulsational pair-instability episodes and ends up lighter than the secondary; and (3) asymmetric core-collapse supernovae, where differential stripping leaves the secondary with a heavier remnant.

Figure 2: BBH merger rate density broken down by MRR and non-MRR contributions for COMPAS and SEVN
Figure 2. Intrinsic binary black hole merger rate density as a function of primary mass M₁ (top), secondary mass M₂ (middle), and mass ratio q (bottom) at redshift z ≈ 0.2, for COMPAS (left) and SEVN (right). The total population (black) is split into MRR (purple) and non-MRR (green) contributions and compared to LVK observational constraints (gray shaded regions). COMPAS predicts MRR systems dominate at high masses, while SEVN shows a more diffuse, subdominant contribution. Both models agree that MRR preferentially populates the high mass-ratio regime (q ≳ 0.6).
Figure 6: Representative evolutionary pathways leading to MRR
Figure 6. Representative evolutionary histories for each of the three MRR channels. Each panel shows the total stellar masses (black and red) and core masses (magenta and green) of both stars as a function of time, with vertical dotted lines marking the two supernova events. From left to right: the core-growth channel in COMPAS, the core-growth channel in SEVN (illustrating how the two codes treat this process differently), the PPISN-shrinkage channel, and the asymmetric core-collapse supernova channel.
Levina et al. 2026 — From cosmological simulations to binary black hole mergers: The impact of using analytical star formation history models on gravitational-wave source populations January 28, 2026

Sasha Levina, Floor Broekgaarden, Lieke van Son, Emanuele Berti, Amedeo Romagnolo, Ruediger Pakmor, Ana Lam
arXiv:2601.20202

To predict how many binary black hole mergers gravitational-wave detectors like LIGO and Virgo should see, theorists must combine a model of how stars form across cosmic time — the star formation history (SFH) — with a model of how binary stars evolve into merging black holes. In practice, most studies use simple analytical SFH models (e.g., power laws or fits to galaxy survey data), but these are approximations of a far more complex underlying reality. This paper asks: how much does the choice of analytical SFH model actually matter?

The team uses the IllustrisTNG cosmological simulation as a ground truth for the star formation history and metallicity evolution of the universe, and compares it against four widely-used analytical SFH prescriptions. By running binary population synthesis models with COMPAS on top of each SFH, they quantify how the choice of SFH propagates through to predictions for BBH merger rates and mass distributions detectable by current and future gravitational-wave observatories.

A key finding (Figure 4 below) is that the predicted merger rate can vary by up to a factor of ~2–3 depending on which analytical SFH model is used, with the differences being most pronounced at high redshifts. This has direct implications for next-generation detectors like the Einstein Telescope and Cosmic Explorer, which will be sensitive to mergers across cosmic history. The metallicity evolution built into each SFH model is the primary driver of the spread: models that assign lower metallicities to star-forming gas at high redshift predict more BBH mergers because low-metallicity stars retain more mass and more readily form heavy black holes.

Figure 6 shows how the mass distributions of detectable BBH mergers shift depending on the SFH choice, with the IllustrisTNG-based prediction bracketed by the analytical models. The results underscore that SFH uncertainty is a non-negligible systematic in population-level gravitational-wave analyses — and that using cosmological simulations as a calibration anchor can help identify and bound this uncertainty.

Figure 4: BBH merger rate as a function of redshift for different SFH models
Figure 4. Binary black hole merger rate density as a function of redshift for five star formation history models: IllustrisTNG (used as ground truth), and four analytical prescriptions. The choice of SFH leads to rate differences of up to a factor of ~2–3, with the spread growing at higher redshifts where metallicity evolution assumptions diverge most strongly.
Figure 6: Detectable BBH mass distributions for different SFH models
Figure 6. Distributions of primary mass M₁ for detectable binary black hole mergers under each SFH model, for current LIGO-Virgo sensitivity (left) and next-generation detector sensitivity (right). The IllustrisTNG prediction (black) is bracketed by the analytical models, illustrating how SFH choice shifts the relative contribution of high-mass mergers — an effect that will become increasingly important as detector sensitivity improves.