Broad search for gravitational waves from subsolar-mass binaries through LIGO and Virgo's third observing run

We present a search for gravitational waves from the coalescence of binaries which contain at least one subsolar mass component using data from the LIGO and Virgo observatories through the completion of their third observing run. The observation of a merger with a component below $1\,M_{\odot}$ would be a clear sign of either new physics or the existence of a primordial black hole population; these black holes could also contribute to the dark matter distribution. Our search targets binaries where the primary has mass $M_1$ between 0.1-100$\,M_{\odot}$ and the secondary has mass $M_2$ from 0.1-1$\,M_{\odot}$ for $M_1<20\,M_\odot$ and 0.01-1$\,M_{\odot}$ for $M_1 \ge 20\,M_\odot$. Sources with $M_1<7\,M_\odot, M_2>0.5\,M_\odot$ are also allowed to have orbital eccentricity up to $e_{10} \sim 0.3$. This search region covers from comparable to extreme mass ratio sources up to $10^4:1$. We find no statistically convincing candidates and so place new upper limits on the rate of mergers; our analysis sets the first limits for most subsolar sources with $7\,M_{\odot}20\,M_{\odot}$ ($M_1<7\,M_{\odot}$). Using these limits, we constrain the dark matter fraction to below 0.3 (0.7)$\%$ for 1 (0.5)$\,M_{\odot}$ black holes assuming a monochromatic mass function. Due to the high merger rate of primordial black holes beyond the individual source horizon distance, we also use the lack of an observed stochastic background as a complementary probe to limit the dark matter fraction. We find that although the limits are in general weaker than those from the direct search they become comparable at $0.1 \,M_{\odot}$.


I. INTRODUCTION
Gravitational waves are now regularly observed by the ground-based observatories Advanced LIGO [1] and Advanced Virgo [2]. Their accomplishments include more than 90 observed binary black hole (BBH) mergers [3][4][5] and a handful of binary neutron star [6,7] and neutron star -black hole mergers [8]. These observations provide a wealth of knowledge for understanding the population of stellar-mass black holes [9], which may have arisen through standard stellar evolution [10]. Field evolution [11][12][13][14] and dynamical formation channels [15][16][17][18][19][20] for compact binaries have been proposed to describe this process. Several observations though have challenged current understanding of stellar formation; these include GW190521 [21,22], the observation of a merger that includes a black hole that may be in the "upper mass gap" caused by pair-instability supernovae [23][24][25]. Observations also confirm the existence of compact objects with secondary mass 1 − 3 M (e.g. GW190814) [9,26]. It has been proposed that such events may be composed of primordial black holes (PBHs) [27][28][29][30][31]. Studies have also shown that the current population of observed binary black holes is compatible with and may include contributions from a population of merging PBHs [32][33][34].
PBHs are hypothesized to form by direct collapse of overdensity in the very early Universe [35,36], and have implications for several astrophysical and cosmological scenarios [37]; these include seeding the first galaxies and the formation of supermassive black holes [38,39], explaining the recent excess power detected by pulsar tim- * alex.nitz@aei.mpg.de ing arrays [40][41][42], and most intriguingly, as a candidate for dark matter [36]. A variety of astrophysical observations have put constraints on PBH abundance (for a review see e.g. [43]). So far, dark matter has evaded all direct searches based on experiments on the earth, including the hunt for weakly interacting massive particles (WIMPs) [44] and axion dark matter [45]. The interest in macroscopic dark matter candidates such as PBHs has been revived in light of the discovery of black hole mergers by LIGO and Virgo [46][47][48].
Among the gravitational-wave catalog, GW190814 is a source of particular interest due to its mass ratio (q∼10) and low spin [26]. Its spin is consistent with a merger of primordial origin [28,49], however more mundane explanations are also possible [50]. Like GW190814, GW190425 [51], GW191219, GW200105, GW20015 [8], and GW200210 all contain secondary components with mass less than 3 M [9], which is lighter than the lower limit observed from X-ray binaries [52][53][54]. If there exists a distribution of binaries composed of PBHs, the population could be convincingly demonstrated by the observation of a similar source with a subsolar-mass black hole secondary. The detection of a subsolar-mass black hole would confidently establish the existence of PBHs [43,55,56] or other exotic physics able to produce subsolar-mass black holes [57][58][59] due to their inability to form through standard stellar evolution [60,61].
In this letter, we conduct a search to answer if LIGO and Virgo have observed any mergers that include a subsolar-mass secondary using the most recent data from the three completed observation runs (O1-O3) [62,63]. Previous work has directly searched for comparable mass subsolar mass sources up to the first half of the third observation run (O3a) [64][65][66][67][68][69], and high-mass-ratio sources through the second observation run (O2) [70]. We perform a broad analysis designed to be sensitive to sources where the primary mass M 1 is 0.1-100 M while the secondary can range from 0.01 − 1 M for M 1 > 20 M and 0.1 − 1 M for M 1 < 20 M . Fig. 1 shows the boundaries of this search and how it compares to previous analyses; our search encompasses the mass ranges of prior searches and for the first time includes the full range where the primary mass is 7 − 20M . The most significant candidate in our search had a false alarm rate of ∼ 1 per 1.8 years. Due to the time searched, we consider this a null observation. With the inclusion of the most recent data, our limits on the rate of mergers are ∼ 8× more stringent than our prior high-mass-ratio search [70] and 1.6× the comparable mass search [64,65] which only included up to O2 and O3a, respectively. We find the merger rate of 0.5-0.5 (1-1) M binaries is < 4400(700) Gpc −3 yr −1 and for 1 -20 M sources is < 65 Gpc −3 yr −1 .
Our non-detection constrains astrophysical population models which predict the binary PBH merger rate. For a fiducial monochromatic mass population, we find that the dark matter mass fraction of PBHs is ≤ 6%(0.3%) for component mass 0.1(1) M . For a two-point mass function where the primary mass is fixed to be 37 M (approximately the average value for all detections), the mass fraction is ≤ 3%(0.03%) for secondary mass 0.01(1) M . Because the stochastic gravitational-wave background [71][72][73] provides information about the population of unresolved, high-redshift sources and the PBH merger rate increases with redshift beyond the detection horizon for individual sources, we also use the nondetection of a stochastic background in O3 data [74] to infer additional constraints.

II. SUBSOLAR MASS SEARCH RESULTS
Our search is conducted using the compact-binary analysis included in the open-source PyCBC toolkit [77]. The analysis identifies candidates [78][79][80], checks for consistency between the data and astrophysical sources [81][82][83], and assesses each candidate's statistical significance. Our analysis targets the parameters space shown in Fig. 1 and is configured similarly to our previously conducted analyses for comparable mass sources in Refs. [64,65] and high-mass-ratio binaries in Ref. [70].
To detect sources, our matched-filter based analysis requires a model of the expected gravitational-wave signal. We model the signal using the TaylorF2 waveform template derived from the post-Newtonian expansion for noneccentric sources [84][85][86][87] where M 1 < 7 M and EOBNRv2 [88], which includes a model of the merger and ringdown phase of the signal for higher masses, elsewhere. We use TaylorF2e [89][90][91] [75,76], and the subsolar mass LVK search (blue dotted) [66]. The boundaries where our search also includes sources with eccentricity up to 0.3 (black dashed), along with the boundaries of our prior comparable mass (green dashed) [64,65] and high-mass ratio (orange dashed) searches [70] are shown. For comparison, we show the masses of reported gravitational-wave sources closest to our search region (yellow stars); if any of these is primordial in origin, the population may extend into our search region.
of the gravitational-wave signal. A discrete set of template waveforms is selected using a stochastic placement algorithm [92] to ensure the signal-to-noise ratio (SNR) loss due to template bank density is 3 − 5% on average. To control for computational cost, the search analyzes data from a minimum of 20 Hz, but this frequency is raised independently for each template to ensure its duration is no more than 512s where M 1 < 7 M and 60s everywhere else.
The public LIGO and Virgo dataset now contains data from all three observing runs through 2021 which amounts to 1.2 years of multi-detector time [62,63]. Our analysis finds no convincing gravitational-wave detections; the most significant candidate has a false alarm rate of 1 per 1.8 years which is consistent with a null observation given the observation time. Using a null detection, we can set an upper limit at 90% confidence on the rate of mergers (shown in Fig. 2) as a function of the binary parameters using where VT is the measured volume-time of the search [93]. We measure the VT of our analysis as a function of the component masses of a binary by empirically measuring the response of our analysis to simulated signal populations. Sources are assumed to have isotropic orientation and sky location in addition to a uniform   [66] are shown for reference. Our limits are up to 3x more constraining than the most recent LVK results, and have improved by ∼ 60% over our prior limits. For comparison, a 20:1 mass ratio system is also shown; high-mass-ratio sources will have reduced sensitivity in comparison to comparable mass at the same chirp mass due to the effect of higher-order gravitational-wave modes and the merger moving into the sensitive frequency band of the observatories. Right: The upper limit as function of the secondary mass for a selection of primary masses. Our new limits for M1 = 37 M (green solid) have improved by 8x over the prior state-of-the art (green dashed) [70]. The one sigma Monte-Carlo statistical uncertainty is shown with shading. distribution in volume. For comparable mass sources we use TaylorF2(e) as a source model. For M 1 > 7 where mass ratio can be up to 10 4 , we use the EOB-NRv2HM model which includes sub-dominant modes of the gravitational-wave signals [88,94]; this allows us to account for the effect of neglecting sub-dominant modes in our search. In addition, we assume that PBHs will have negligible spin which is consistent with the predicted spin distribution [95][96][97][98][99]. For large total mass binaries, M 1 + M 2 >∼ 35 − 80 M , however, we note that nonnegligible spin may be induced depending on the level of accretion [100].

III. OBSERVATIONAL CONSTRAINTS ON PRIMORDIAL BLACK HOLE DARK MATTER CONTRIBUTION
Limits on the observed merger rate can constrain the dark matter mass fraction of PBHs given a population model for the binary merger rate. We use the model in Refs. [47,[101][102][103], which assumes a Poisson spatial distribution for PBHs when they initially form from large overdensity collapse in the early Universe. A nearby pair of black holes form a binary after decoupling from the cosmic expansion and then inspiral due to gravitational radiation. Given a general mass distribution P (m), the binary merger rate in units of Gpc −3 yr −1 is where m 1/2 are in unit of M , t is the cosmic time and t 0 is the age of the Universe, both in units of years, σ eq = 0.005 characterizes the variance of density inhomogeneity from dark matter at the mass density equality era. The effective mass fraction is defined to f 53/37 where f PBH is the true mass fraction and the suppression term S accounts for binary disruption after formation.
Refs. [104][105][106] showed the disruption by nearby PBH clusters can reduce the merger rate by orders of magnitude if PBHs account for nearly all dark matter but become negligible when f PBH O(0.1%).
In this letter we focus on a common fiducial mass distribution following [64,65,70]  FIG. 3. Left: 90% upper limits on the fraction of dark matter accounted for by PBHs assuming a delta function mass distribution. The orange lines are constraints from our direct searches with data through O2, O3a, and O3b. The blue solid and dashed lines are constraints from the stochastic background. Right: 90% upper limits on the abundance of PBHs for the secondary black hole assuming a two-point delta function for mass. To allow consistent comparison with previous work , we follow Ref. [64] to fix the abundance of the primary mass to be f primary PBH = 3 × 10 −3 which is obtained by fitting the observed black hole merger population [34]. The solid lines are constraints for different representative primary mass. The dash line is the prior constraint based on data only through O2.
distributions can also be derived from our limits (see e.g. [64,107,108]).

A. Constraints from our direct search
We first use our direct search limits to constrain the dark matter mass fraction from PBHs. Results are shown in Fig. 3. In the left panel the PBH mass function is assumed to be a delta distribution. For 1 and 0.1 M reference masses, the fraction of dark matter composed of PBHs can not exceed 0.3% and 6%, respectively. The right panel of Fig. 3 shows constraints where we choose the mass function to be a two-point delta distribution. We fixf PBH = 0.3% for the primary mass and constrain thef PBH of the secondary PBH. This choice is motivated by Ref. [34] which fitted a PBH mass function and abundance with the current black hole merger observations [3] assuming they are primordial in origin. We consider five representative primary component masses 1, 5, 20, 37 (approximately the average mass of all gravitational wave events) and 60 M , and secondary mass range in [0.01,1] M . The previous constraint for M 1 = 37 M , which included data only up to O2, are also shown for comparison. The O3b constraints improved over our previous by one order of magnitude due to upgrades of the advanced detectors for O3 and the longer duration of observation [3,109] B. Constraints from the stochastic background In addition to individually resolvable sources, the incoherent superposition of all binary PBH coalescences produces a stochastic background of gravitational waves [110,111]. While the astrophysical binary black hole merger rate approximately follows the star formation rate peaking at redshift ∼ 2 [112], the PBH merger rate gets higher when redshift increases [47,101]. The search for an isotropic stochastic gravitational-wave background by the LVK did not find any significant excess energy in their O3 analysis [74]. We investigate the implication of this non-detection of a stochastic background for the PBH abundance.
The stochastic background is characterized by where dρ is the gravitational-wave energy density in a frequency bin [ν, ν + dν], normalized by critical energy density for a flat Universe ρ c . Given the merger rate model for PBH binaries from Eq. 2, the stochastic background can be computed as where H(z) is the Hubble parameter at redshift z, dE s /df s is the energy spectrum from a single inspiraling source evaluated in the source frame. We use the waveform template IMRPhenomD [113,114] to compute dE s /df s and integrate Eq.4 up to z = 20; the contributions from higher redshifts can be neglected. We assume a delta mass distribution for PBH binaries and have verified its stochastic background spectrum follows a power law with slope index 2/3 within the LIGO/Virgo's most sensitive frequency band (∼ 20 − 100 Hz). Therefore, we apply the LVK O3 upper limit for a power law spectrum with index 2/3 that Ω GW (ν = 25Hz) < 3.4 × 10 −9 for a log-uniform prior on Ω GW and Ω GW (ν = 25Hz) < 1.2 × 10 −8 for a uniform prior [74] to constrain the PBH mass fraction. The result is shown in the left panel of Fig. 3 accompanied by the results from our direct search. As shown, the constraint is generally weaker than the direct search by one order of magnitude, but becomes comparable for 0.1 M where the individual source search is less sensitive. The stochastic background from PBH mergers at high redshift complements the direct search which only probes the local Universe.

IV. CONCLUSIONS
Subsolar mass black holes would establish the existence of a population of PBHs [43,55,56] or hint at new physics outside of the standard model of stellar evolution [57][58][59]. We have performed a broad search for compact-binary coalescences where at least one component is less massive than one solar mass; binary sources span from comparable mass to extreme mass ratios up to 10 4 . We use our null detection to set new upper limits on the rate of mergers and the implied fraction of dark matter composed of PBHs. The top candidates from our analysis along with the configuration files necessary to reproduce the search are available at https://github.com/gwastro/subsolar-O3-search.
We have investigated the constraints on PBH abundance from the non-detection of the direct search for individual sources and a search for a stochastic background from the superposition of unresolved sources [74]; these probe gravitational-wave sources from the local and high redshift Universe, respectively. The observation of binary black hole mergers at high redshift is a key way to distinguish primordial from stellar origin because the former can merge before star formation [115][116][117].
Assuming a non-detection, we expect the next observing run (O4), with expected data available at the end of 2024, will be able to improve constraints by another factor of 2−4× [118]. The future direct detection of subsolar mass black hole binaries or an excess in the merger rate at high redshift would give decisive evidence for PBHs.

ACKNOWLEDGMENTS
We acknowledge the Max Planck Gesellschaft and the Atlas cluster computing team at AEI Hannover for support. This research has made use of data, software and/or web tools obtained from the Gravitational Wave Open Science Center (https://www.gw-openscience.org), a service of LIGO Laboratory, the LIGO Scientific Collaboration and the Virgo Collaboration. LIGO is funded by the U.S. National Science Foundation. Virgo is funded by the French Centre National de Recherche Scientifique (CNRS), the Italian Istituto Nazionale della Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions by Polish and Hungarian institutes.