On the mass-function of GWTC-2 binary black hole systems and their progenitors

The distribution of LIGO black hole binaries (BBH) shows an intermediate-mass range consistent with the Salpeter Initial Mass Function (IMF) in black hole formation by core-collapse supernovae, subject to preserving binary association. They are effectively parameterized by mean mass $\mu$ with Pearson correlation coefficient $r = 0.93\,\pm\,0.06$ of secondary to primary masses with mean mass-ratio $\bar{q}\simeq 0.67$, $q=M_2/M_1$, consistent with the paucity of intermediate-mass X-ray binaries. The mass-function of LIGO BBHs is well-approximated by a broken power-law with a tail $\mu\gtrsim 31.4M_\odot$ {in the mean binary mass $\mu=\left(M_1+M_2\right)/2$}. Its power-law index $\alpha_{B,true}=4.77\pm 0.73$ inferred from the tail of the observed mass-function is found to approach the upper bound $2\alpha_S=4.7$ of the uncorrelated binary initial mass-function, defined by the Salpeter index $\alpha_S=2.35$ of the Initial Mass Function of stars. The observed low scatter in BBH mass ratio $q$ evidences equalizing mass-transfer in binary evolution prior to BBH formation. At the progenitor redshift $z^\prime$, furthermore, the power-law index satisfies $\alpha_B^\prime>\alpha_B$ in a flat $\Lambda$CDM background cosmology. The bound $\alpha_{B,true}^\prime \lesssim 2\alpha_S$ hereby precludes early formation at arbitrarily high redshift $z^\prime \gg1$, that may be made more precise and robust with extended BBH surveys from upcoming LIGO O4-5 observations.


INTRODUCTION
It is well known that stellar black holes (BHs) are formed by gravitational collapse or fallback from corecollapse supernovae of high-mass stars (Oppenheimer & Snyder 1939;May & White 1966;Bethe 1990;Fryer 1999;Fryer & Kalogera 2001;Janka 2012). A key open question is their mass-function and its relation to the Salpeter stellar initial mass function (IMF) (Salpeter 1955), where M is the mass in units of M , dN is the number of stars in the corresponding mass range and α S = 2.35 (Salpeter 1955). In the case of M > 5 M , Eq.
(1) appears to be remarkably universal and robust against variations in metallicity following studies of star forming regions of the Large Magellanic Cloud with Z = 0.008 (Da Rio et al. 2009;Gouliermis et al. 2006), though some dependencies may exist more generally (Conroy & van Dokkum 2012;Geha et al. 2013;Kalirai et al. 2013). Through a complex process of massive stellar binary evolution, Eq. (1) is expected to have an imprint on the mass-function of their offspring in binary black holes (BBHs) recently observed by LIGO 1 -Virgo (Fig. 1). Starting from Eq. (1) normalized to ξ = M/M , the cumulative count ∝ M −α+1 and mass ∝ M −α+2 depends on the lower cut-off ξ c = M c /M , where M c denotes the minimum for hydrogen ignition. Slightly above the Jupiter mass, 0.01 < ξ c < 0.1, the precise value of which is uncertain. We assume ξ ∼ 20 for the low mass cut-off to produce BHs (e.g. Fryer 1999), whereby M BH /M tot = (ξ c /ξ) 1.35 ∼ 7 − 16 %, N BH /N tot = (ξ c /ξ) 0.35 ∼ 0.004 − 0.08 %. (2) Observational support for these percentages is found in population studies of the Milky Way (Elbert et al. 2018), showing a number fraction of BHs of about 0.09 % with aforementioned ξ c 0.1. This implies an expected mean BH-progenitor stellar massM BH = M BH /N BH = (M BH /M tot ) (N BH /N tot ) −1M * , whereM * = 0.1-− 0.5 M is the mean stellar mass. For the Milky Way, M * /N tot can be estimated from the stellar population by total mass 5 × 10 10 M (Licquia & Newman 2015) and count (1 − 5) × 10 11 . Based onM BH , the typical BH-progenitor masses are somewhere between 20 and 900 M , that we summarize by the geometric mean As a progenitor stellar mass, Eq. (3) is subject to mass-loss in the evolution of their progenitor binaries of massive stars. Preserving binary association limits mass-loss, in winds from the upper atmosphere of a star and in the core-collapse supernova producing the BHs. Regarding the first, simulations predict a positive correlation of mass loss with stellar mass, which tends to suppress the population of heavier stars (ameliorated by their relatively shorter lifetimes, Belczynski et al. 2012a) and restrict the maximum BH mass produced in the final core-collapse supernova (Fryer & Kalogera 2001). To illustrate the second, consider two (consecutive) spherically symmetric explosions in circular binaries, neglecting kick velocities and induced ellipticity (Wheeler et al. 1975). Mass-loss in each is bounded by one-half the total mass for the binary to survive. Let 0 < < 1 denote the fractional mass loss in each supernova. As the more massive star explodes (M 2 ≤ M 1 ), it blows off a stellar mass fraction leaving a BH mass M 1 ≤ (1 − ) M 1 and total mass of the binary where the single prime now refers to masses just prior to the supernova event, after mass-loss in winds. The condition for the binary to survive implies < (1 + M 2 /M 1 ). The explosion of the secondary similarly must satisfy M 2 < 1 2 M ≤ 1 2 (M 2 + (1 − ) M 1 ), giving the more restrictive condition < min 2 3 , 1 2 (1 + M 2 /M 1 ) . In the approximation that the mass-loss fraction in the progenitor stars suppresses the mass of their BH remnants, approximately equal-mass stellar binaries are expected to produceM To date, intermediate-mass black holes (IMBHs) such as Eq. (4) have remained elusive in electromagnetic (EM) observations for reasons not well understood. However, such may derive from a mass-correlation in binary stellar progenitors. If so, BBH systems form rather promptly with high-mass BHs with high-mass stellar companions practically undetectable (in the electromagnetic spectrum) by the short lifetime of the latter. Indeed, the observed black hole masses of X-ray binary systems are (3−18)×M (e.g., Remillard & McClintock 2006;Miller-Jones 2014) -much smaller than the geometric mean in Eq. (3) predicted by the Salpeter IMF.
Over the past few years, however, a new window opens to probe BBH systems by gravitational radiation (Fig. 1). The GW150914 event, the first LIGO-detected, is a merger of the BBH system with masses of 36 A key question is the cosmological origin in binary stellar progenitors of the BBHs, notably relative to the peak in the cosmic star formation rate around redshift z 1.9 (Madau & Dickinson 2014). To address whether the majority of the progenitor systems are found before or after, we study the binary mass-function of BHs of GWTC-2 in the context of Eq. (1) and possible effects of cosmic time-dilation.
To start, we consider some statistical properties and mass-function of the intermediate-mass BBHs in GWTC-2 ( §2), focused on power-law behavior in the tail of this population. Cosmic time-dilation effects on the power-law index is studied in §5. We summarize our results in §6.

MASS-FUNCTION INDEX IN GWTC-2
We consider LIGO compact binary mergers excluding events involving neutron stars (GW170817, GW190425 and GW190814) based on estimated masses, and three more by limited False Alarm Rates (Abbott et al. 2019;Abbott et al. 2021a).
In the BBHs in GWTC-2 ( Fig. 1), the primary (M 1 ) and secondary masses (M 2 ) are strongly correlated with Pearson correlation coefficient r and slope s, representing a meanq of their mass-ratio M 2 /M 1 . The minimum and maximum of the mean binary masses are 6.9 M and 79.1 M , respectively. Fig. 2 shows the normalized mass difference δM/µ, refers to the mean binary mass. Shown is a distinctly narrow skewed Gaussian mass-distribution, evidencing a relatively small mass-difference between primary and secondary in these binaries. The relatively modest standard deviation of 22% in the signed mass-differences ( Fig. 2) shows the BBHs are composed of roughly equalmass BHs, leaving the mean mass as the primary parameter in their mass-distribution. Figs. 1-2 highlight a rather tight correction between the two masses, see also Fishbach & Holz (2020). This permits parameterizing BBHs by their mean mass µ. After all, their observation is determined by their luminosity in gravitational-wave emission, primarily as a function of chirp mass, satisfying where x = δM/µ. The distribution of δM/µ shown in Fig. 2 shows x 2 0.04. By (7), the chirp mass M of our BBH binaries in Fig. 1 is tracked by the mean mass µ to better than 10%. In light of this and aforementioned slope in Fig. 1, µ provides a statistically equivalent parameterization to primary mass M 1 (Abbott et al. 2021b), satisfying where κ = 1 2 (s + 1) = 0.8350, ignoring any potential mass-dependence in this correlation.
We next turn to the observed and true mass-functions of the BBH mergers (Abbott et al. 2021b).

BROKEN POWER LAW MASS DISTRIBUTION
LIGO BBHs mergers show a power law behavior in their true (astrophysical, per unit volume and observation time) mass distribution inferred from the observed distribution, here shown in Fig. 3 -a loglog plot of data from Figs. 3-4 of Abbott et al. (2021b).
Here, our focus is on the tail µ ≥ µ b2 = 31.4M , shown in Section C of Fig. 3, rather than the intermediate mass range of Section B (µ b1 = 4.64M < µ < µ b2 ). The location of the second break is consistent with µ break = 33.15M estimated by (8) ignoring systematic uncertainties in this conversion process. Revisiting uncertainty by scatter given the relatively modest number of events in the tail of the present survey, Fig. 4 repeats the estimate of α B = α B,obs using four bins in a loglog plot of event count versus µ. By (9), we infer α B = 4.08 ± 0.73, α B,true = 4.77 ± 0.73.
Our estimate is hereby consistent with but slightly less than α 2 = α B,true of Abbott et al. (2021b). In particular, α B,true satisfies the astrophysical bound discussed in the following section.

ASTROPHYSICAL BOUNDS ON POWER-LAWS
The observed steep index in Eq. (9) satisfies some priors derived from aforementioned Salpeter IMF as follows. From Eq. (1), the progenitor binary mass-function M 1 and M 2 of the BBH can be expressed as a function of mean µ and normalized mass difference x = ν /µ , ν = (M 1 − M 2 )/2. (Equivalently, consider µ and The two-parameter binary mass-function ψ B (µ , x) has two natural limits derived from essentially correlated or uncorrelated stellar masses.
On the first, we note that nearby Galactic open stellar clusters show a large fraction of O-type stars born as binaries (Sana et al. 2012). Radiation-hydrodynamic simulation suggests fragmentation of rotating gas disk by gravitational instabilities leads to the formation of massive binary stars with similar masses (Krumholz et al. 2009;Sana et al. 2012;Moe & Di Stefano 2017;Li et al. 2022). On the other hand, such mass-correlation is ex- pected to weaken as binary separation becomes large. Accordingly, therefore, the binary stellar mass-function is expected to be bounded by either of the two limits for the correlated and, respectively, uncorrelated case with M 1,2 = µ ± ν .
According to the above, we expect ψ B ∝ (µ ) The index Eq. (9) is significantly steeper than aforementioned Salpeter IMF with α S = 2.35 of their progenitor stellar mass-function. In fact, Eq. (9) shows an index essentially equal to twice the Salpeter value of a binary IMF with uncorrelated masses.
Some of the steepening in Eq. (9) might alternatively be attributed to cosmological time-dilation.

STEEPENING BY COSMIC TIME-DILATION
For a given cosmological background evolution, redshifts z of BBH progenitors can be traced back by times of coalescence (Peters 1964;Celoria et al. 2018;Mapelli 2018), where c is the speed of light, G is the gravitational constant, a is the initial separation, and M 1 and M 2 are mass of the primary BH and that of secondary BH. That is, where µ is the mean mass, and t H = 13.7 Gyr is the Hubble time. The PDF in the source frame, P (µ, z ) = P (µ, z) (∂z/∂z ) (1 + z ), derives from invariance of count N and µ, where dN = ψdµdT for an IMF ψ and time interval T , whereby P (µ, z)dµdzdt = P (µ, z )dµdz dt . Consequently, the transformation (z, µ) → (z , µ) with z = z (z, µ) gives P (µ, z )(∂(µ, z)/∂(µ, z ))dµdz = P (µ, z )Jdµdz with Jacobian J = ∂z/∂z , where z and z are related on given cosmological background evolution. Hence, P (µ, z ) = P (µ, z)(∂z/∂z )(dt/dt ), i.e.: as asserted, where dt/dt = 1+z , with the cosmic steepening factor For a concrete illustration, we consider the following Ansatz for progenitor redshift z : in the two parameters f > 0 and γ > 0. Here, µ 0 = 31.4 M is the minimum of the tail of the BBH mass-function (Fig. 3). For an initial separation a, the merging time t c satisfies with positive (γ > 0.75), negative (γ < 0.75) or neutral (γ = 0.75) correlation between t c and µ. For the PDF of the observed mergers, we consider where α B is the observed value Eq. (9), about the mean redshiftz = 0.51 in light of the modest standard deviation σ z = 0.2. Based on Eq. (14) and Eq. (17), we numerically solve for z of the progenitor given a merger event at z on a three-flat ΛCDM background. We highlight three model relationships between the initial separation a and coalescence time t c with mass parameterized by γ color-coded with blue, orange, and green, obtained with Python package cosmology with Hubble parameter H 0 = 67.8 km s −1 Mpc −1 and matter density Ω m = 0.307.
Figs. 6 shows the results for (14). Long merger times push the origin of high-mass binaries to high z . In treating the (tail of the) BBH of GWTC-2 as a uniform population, the condition z < ∞ sets a bound f 0.15. Shown is z versus µ for selected values γ = 0.30, 0.75, and 1.00 (Case 1-3), in blue, orange, and green, respectively. Table 1 summarizes the steepening in the PDF in the source frame compared to the observed PDF. While subtle (hard to discern by eye), the steepening present in (14) is consistent with the anticipated modest change due to J z shown in Fig. 5. Steepening appears for γ > 0.75 (Case 3), otherwise absent for γ ≤ 0.75 (Case 1-2).  Steepening by cosmic time-dilation is noticeable when progenitor and merger redshifts differ by order unity (z − z 1). Such may push α B across the upper bound 2α S in Eq. (12) in the binary IMF of progenitor stellar systems. Avoiding this limits the origin of the BBH progenitors to relatively low-z late-time cosmology, effectively posterior to the peak in the cosmic star formation rate (Hopkins & Beacom 2006;Madau & Dickinson 2014). The BBH mass-function of LIGO BBH mergers is identified with the Salpeter IMF of their stellar progenitors in the power-law tail µ ≥ 31.4 M (Fig. 3), effectively parameterized by mean mass µ due to the implied tight correlation with chirp mass Eq. (8). The mass scale Eq. (4) expected from the Salpter IMF, subject to preserving binary association in the final phase of binary evolution, is consistent with the observed masses in the tail shown in Figs. 3-4.

DISCUSSION AND CONCLUSIONS
Following a detailed consideration of GWTC-2 data, our main findings indicate 1. A tight correlation Eq. (5) between primary and secondary masses, consistent with the paucity in intermediate mass X-ray binaries; 2. A broken power-law mass-function with a tail beyond µ 31.4M . The power-law index α B,true 4.77 ± 0.73 is consistent with uncorrelated stellar progenitor masses at birth Eq. (11) by approaching the limit 2α S = 4.7 defined by the Salpeter index α S ; 3. A power-law index of µ is subject to steepening due to cosmological time-dilation (Fig. 5), e.g., when mass and orbital separation are positively correlated. The condition Eq. (12) hereby bounds the mean of progenitor redshiftz .
4. Subject Eq. 12, Table 1 suggests that a progenitor origin at z 1 is excluded, assuming the tail of the BBH population to derive from a uniform population. Progenitors hereby appears to be in the relatively recent epoch of cosmic star formation, about or posterior to the peak in the star formation rate.
Conceivably, the last finding can be made more rigorous with BBH surveys from upcoming O4-O5 observations. Such bounds hold promise to distinguish between an origin related to the peak in the cosmic star formation rate and an association with Pop III stars (e.g. Yajima & Khochfar 2015; Kulkarni et al. 2014;Kinugawa et al. 2014), previously considered for their relatively high mass of several tens of M (Kinugawa et al. 2014;Hosokawa et al. 2011;Kinugawa et al. 2016Kinugawa et al. , 2021. The tight correlation between primary and secondary BH masses and the uncorrelated binary progenitor masses in our findings is perhaps paradoxical. However, the pathway to black hole binaries from stellar progenitor systems is a complex process of binary stellar evolution including (uncertain) mass-losses in stellar winds (e.g. Krtička 2014;Chen et al. 2015;Belczynski et al. 2012b;Elbert et al. 2018), mass-transfer potentially equalizing masses (e.g. Kinugawa et al. 2014), terminating in two core-supernovae. These processes are subject to a stringent selection criterion of preserving binary association, i.e., a 50% mass-loss limit (in the idealized case of circular binary motion).
Equalizing masses in the progenitor systems offers some hints at their stellar evolution, perhaps with further contributions to steepening in the power-law index from wind mass-loss (cf. Vink et al. 2001);Vink et al. 2011;Chen et al. 2015), e.g., in binary association around cold red giants (Decin et al. 2020). While a detailed study of these complex radiation-hydrodynamical processes is outside the scope of this work, these processes and their down-selection effects in binary black hole formation are probably instrumental to understanding the detailed nature of the progenitor systems, here identified with ab initio uncorrelated progenitor stellar binary masses in view of α B 2α S in Eq. (5).
Data Availability. The data underlying this article were accessed from LIGO.