Measuring the Properties of Active Galactic Nuclei Disks with Gravitational Waves

Gravitational waves from the mergers of binary black hole (BBH) systems have recently transformed astronomy. However, the astrophysical origins of these events are still uncertain. There are two main proposed astrophysical pathways to the mergers: (i) isolated binary evolution via mass transfer, including a common envelope phase in galactic fields; and (ii) dynamical formation in dense environments. Each pathway is associated with different distributions of black hole (BH) spin (Mandel & O'Shaughnessy 2010; Fishbach et al. 2017; Stevenson et al. 2017; Talbot & Thrane 2017; Wysocki et al. 2019) and binary eccentricity (Rodriguez et al. 2018; Samsing 2018; Lower et al. 2018; Romero-Shaw et al. 2019, 2021; Zevin et al. 2021; Gayathri et al. 2022). Measuring BBH spins and eccentricity with gravitational waves can therefore be used to determine how and where BBHs are assembled (Abbott et al. 2019, 2021).

Active galactic nuclei (AGNs) are expected to contain a dense population of stars and stellar remnants, such as stellar-origin BHs (Morris 1993; Miralda-Escudé and Gould 2000; Hailey et al. 2018; Generozov et al. 2018). BBH systems can form via close encounters in this dynamically "hot" environment, but are often rapidly "ionized" via tertiary encounters (Antonini & Rasio 2016; Fragione et al. 2019). The dense nuclear population and AGN gas disks (when present) can interact, resulting in an embedded population of stars and BHs within the disk. These embedded objects can weakly perturb the surface-density profile of the gas disk, resulting in gas torques within the AGN disk that allow for Type I (non-gap-opening) migration of the embedded objects (McKernan et al. 2012). The differential migration rates of the objects encourage binary formation, leading to compact binary mergers detectable with LIGO–Virgo (McKernan et al. 2014; Acernese et al. 2015; LIGO Scientific Collaboration et al. 2015; Bartos et al. 2017; Stone et al. 2017; Wang et al. 2021b). Kicked merger products are generally retained by the deep potential well, allowing for hierarchical BBH mergers. If the majority of the mergers observed by LIGO–Virgo are assembled in an AGN disk, it may therefore be possible to reverse engineer conditions beneath the AGN photosphere (McKernan et al. 2018; Wang et al. 2021b).

The gaseous disks in AGNs likely serve to align (to varying degrees) both BH spin vectors ⁷ ${\vec{\chi }}_{\mathrm{1,2}}$ and the binaries' orbital angular momentum vectors L —depending on the density and age of the disk (Bogdanović et al. 2007). On the other hand, tertiary encounters with the binaries in the disk tend to misalign ${\vec{\chi }}_{\mathrm{1,2}}$ relative to L —depending on the timescales of the encounters (Liu & Lai 2017; Tagawa et al. 2020a). The competing effects of the gaseous disk and the dynamical encounters on the BBHs in AGNs determine the distribution of the BBH spin orientations. On the other hand, binaries born in the field have ${\vec{\chi }}_{\mathrm{1,2}}$ nearly aligned to L , with a small spread due to supernova kicks (Kalogera 2000; Mandel & O'Shaughnessy 2010; Dominik et al. 2013; Giacobbo et al. 2017; Eldridge et al. 2017; Olejak et al. 2020). Finally, dynamically assembled binaries (e.g., in globular clusters) exhibit no correlation between L and ${\vec{\chi }}_{\mathrm{1,2}}$. In this paper, we propose a phenomenological model for the distribution of BH spins in AGNs, with a goal of capturing the salient features predicted from theoretical modeling.

The remainder of this paper is organized as follows. In Section 2, we review the spin orientations of the BBHs in AGNs. In Section 3.1, we present a phenomenological model describing AGN BBH spin orientations. In Section 4, we present the results of a simulated study, demonstrating the effectiveness of the model.

Modeling AGNs is challenging, due to the interplay between the dynamics of the gas, the scattering binaries, and the feedback from the central supermassive black hole (SMBH). Commonly used disk models span wide ranges of disk density and geometry (see, e.g., Sirko & Goodman 2003; Thompson et al. 2005) and provide broad estimates for merger rates in AGNs (e.g., McKernan et al. 2018; Gröbner et al. 2020; Tagawa et al. 2020a). However, we have qualitative predictions for the spin orientation population properties of merging BBH systems, using Monte Carlo and N-body simulations to identify the key features. Some of these predictions are tabulated in Table 1.

Table 1. The Phenomenology of BH Spins for Binaries Merging in AGNs

	σ 12 ≤ 1		σ 12 > 1
	(a) Old, Dense AGNs		(b) Mid-aged, Dense AGNs
	${\overline{\tau }}_{\mathrm{AGN}}\gtrsim 5\,\mathrm{Myr}$	L ∥ J AGN	${\overline{\tau }}_{\mathrm{AGN}}\sim 0.5-5\,\mathrm{Myr}$	L ∥ J AGN
	$\bar{\rho }\gt {10}^{-11}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}$	${\vec{\chi }}_{1}\parallel {\boldsymbol{L}}$	$\bar{\rho }\gt {10}^{-11}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}$	${\vec{\chi }}_{1}\parallel {\boldsymbol{L}}$
σ 1 ≤ 1	tdamp < tenc	${\vec{\chi }}_{2}\parallel {\vec{\chi }}_{1}$	tdamp < tenc	${\vec{\chi }}_{2}$ isotropic
	Δm ≳ 0.01 mi		Δm ≪ 0.01 mi
	gas → 1		0 < gas < 1
	(c) Old, Dilute AGNs		(d) Young, Dilute AGNs
	${\overline{\tau }}_{\mathrm{AGN}}\gtrsim 5\,\mathrm{Myr}$	L isotropic	${\overline{\tau }}_{\mathrm{AGN}}\lesssim 1\,\mathrm{Myr}$	L isotropic
	$\bar{\rho }\lt {10}^{-11}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}$	${\vec{\chi }}_{1}\parallel {{\boldsymbol{J}}}_{\mathrm{AGN}}$	$\bar{\rho }\lt {10}^{-11}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}$	${\vec{\chi }}_{1}$ isotropic
σ 1 > 1	tdamp > tenc	${\vec{\chi }}_{2}\parallel {\vec{\chi }}_{1}$	tdamp > tenc	${\vec{\chi }}_{2}$ isotropic
	Δm ∼ 0.01 mi		Δm ≪ 0.01 mi
	0 < gas < 1		gas ≪ 1

Note. The two rows correspond to the different values of the population parameter σ₁, defined in Section 3.1, which controls the effective density of the AGN. The two columns correspond to the different values of the population parameter σ₁₂, defined in Section 3.1, which controls the effective age of the AGN. Each cell is divided into two. The left-hand side of the cell describes the properties of the AGN, while the right-hand side describes the distribution of BH spins at the time of formation.

Download table as:  ASCIITypeset image

In this section, we discuss three predictions for the spin distributions of merging BBHs at the time the binary is formed. The spin vectors subsequently evolve through general relativistic precession of the orbital plane. Nonetheless, the orientation of the spin vectors at the merger contains information about the orientation at formation. We describe the phenomenology of the spin vectors at formation, and then discuss the resulting phenomenology at the time of the merger.

• 1.  
  Gas accretion torques the BH spins to align with the disk. BHs embedded in the AGN disk early in the disk's lifetime should have isotropically distributed spins at formation (Mapelli & Gualandris 2016; Tagawa et al. 2020a). As the BHs of mass m migrate in the disk, they accrete Δm disk-gas mass, resulting in a torque pointing into the plane of the AGN disk. The magnitude of the torque on the BH spins depends on Δm (Bogdanović et al. 2007), and can be summarized as follows. In sufficiently long-lived, dense disks (Table 1(a)), fully embedded BHs accrete more than 1%–10% of their initial mass Δm ≳ 0.01 m_i . The resultant torque from gas accretion onto the embedded BHs reorients the BH spin vector $\vec{\chi }$ to align with the angular momentum vector for the AGN disk J _AGN (e.g., Bogdanović et al. 2007). Alternatively, for dilute disks ($\bar{\rho }\lt {10}^{-11}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}$, Table 1(c) and (d)), Δm ≪ 0.01 m_i implies that the BHs are not torqued into alignment with the disk in ${\overline{\tau }}_{\mathrm{AGN}}\sim 0.5-5\,\mathrm{Myr}$ (Bogdanović et al. 2007). Similarly, if the AGN disks are dense, but typically short-lived (${\overline{\tau }}_{\mathrm{AGN}}\lesssim 1\,\mathrm{Myr}$), this effect is weaker (Table 1(a) and (c)). Note that the details of accretion onto objects embedded in AGN disks are subject to much uncertainty. Feedback, turbulence, and interactions can alter the gas flow dynamics inside the BBH Hill sphere, possibly inhibiting a high rate of accretion onto the component BH, therefore limiting the average torque magnitude (Hankla et al. 2020). Lower accretion rates can result in longer timescales (${\overline{\tau }}_{\mathrm{AGN}}\gg 5\,\mathrm{Myr}$) for BHs to be torqued into alignment with the disk.
• 2.  
  Gas torques dampen the BBH orbital angular momentum. BBH systems in AGN disks experience corotating (horseshoe) gas torques as they migrate through the disk (Tagawa et al. 2020b). Corotation torques dampen the binary's eccentricity and drive L into alignment with J _AGN on a characteristic timescale (Tanaka et al. 2002):

  $$\begin{eqnarray}&&{t}_{\mathrm{damp}}=\displaystyle \frac{{M}_{\mathrm{SMBH}}^{2}{h}^{4}}{{m}_{{\rm{b}}}{\rm{\Sigma }}{a}^{2}{\rm{\Omega }}},\end{eqnarray} \tag{ 1 }$$

  where M_SMBH is the SMBH mass, m_b is the BBH system's total mass, h = H/a is the disk aspect ratio, Σ is the disk surface density, and Ω is the Keplerian orbital frequency. In addition to aligning L with J _AGN, dynamical gas friction can promote binary hardening—the process of loosing orbital energy and tightening the orbit (Baruteau et al. 2011). If t_damp is smaller than the lifetime of the AGN (typically in long-lived AGNs, ${\overline{\tau }}_{\mathrm{AGN}}\gtrsim 5\,\mathrm{Myr}$, Table 1(a) and (c)), the BBH orbital angular momentum vector L will align with J _AGN by the time the binary merges.
• 3.  
  Dynamical encounters excite the BBH orbital angular momentum. The dense environment and high escape velocity in AGNs facilitate dynamical encounters (Yang et al. 2019a, 2019b; Samsing et al. 2022). Multiple migrators can quickly interact with each other, potentially leading to complex or chaotic dynamical encounters moderated by the disk gas (Wang et al. 2021a). Tertiary encounters of binaries with compact objects in the disk or in the spherical nuclear population component can harden or soften the BBH systems (Leigh et al. 2018; Yang et al. 2019a), increase the BBH orbital eccentricity (Samsing et al. 2022), and alter the orbital angular momentum of the BBH (Tagawa et al. 2020b, 2020a). Close encounters with a tertiary object on a disk-crossing orbit can perturb the orbital angular momentum of the BBH on the timescale of the encounter (t_enc; Leigh et al. 2018), which depends on the density of the nuclear star cluster (ρ_NSC), the BBH location in the disk, and the efficiency of the disk capture (a function of $\overline{\rho },{\overline{\tau }}_{\mathrm{AGN}}$). Small values of t_enc (Table 1(c) and (d)) lead to more binaries with L misaligned with J _AGN at the merger.

The relation between the dampening and the dynamical encounter timescales may provide further details about the AGN. For example, t_damp > t_enc could occur if (i) AGN disks are not long-lived ${\overline{\tau }}_{\mathrm{AGN}}\sim 0.5-5\,\mathrm{Myr}$, and the spherical population is not efficiently captured by the disk, resulting in lots of dynamic interactions (McKernan et al. 2018; Tagawa et al. 2020b); (ii) the BBHs are positioned in a short-lived inner disk, where the encounter rate with the spherical component is significantly higher (a is small; Leigh et al. 2018); (iii) ρ_NSC is large (e.g., the nuclear star cluster is cuspy, not cored; McKernan et al. 2018; Tagawa et al. 2020b); or (iv) the fraction of BBH systems hardened via gas torques _gas ⁸ is not particularly efficient relative to dynamical hardening (Stone et al. 2017).

The key points from this section are as follows. Old, dense AGNs (Table 1(a)) produce a population of BBH systems with ${\vec{\chi }}_{1}$ preferentially aligned with the orbital angular momentum and ${\vec{\chi }}_{2}$ preferentially aligned with ${\vec{\chi }}_{1}$. Mid-aged, dense AGNs (Table 1(b)) produce BBH systems with ${\vec{\chi }}_{1}$ preferentially aligned with the orbital angular momentum, but ${\vec{\chi }}_{2}$ is not correlated with ${\vec{\chi }}_{1}$. Old, dilute AGNs (Table 1(c)) produce BBH systems where ${\vec{\chi }}_{1}$ is not preferentially aligned with the orbital angular momentum. However, ${\vec{\chi }}_{2}$ is preferentially aligned with ${\vec{\chi }}_{1}$. Young, dilute AGNs (Table 1(d)) produce BBH systems where ${\vec{\chi }}_{1}$ is uncorrelated with the orbital angular momentum and ${\vec{\chi }}_{2}$ is uncorrelated with ${\vec{\chi }}_{1}$.

3.1. Model Description

In this section, we construct a model for the spin orientation of BHs (at the time of formation) in merging binaries residing in the AGN disk. It will be useful to employ two coordinate systems. Vectors with no prime are measured with respect to the orbital angular momentum vector, such that

$$\begin{eqnarray}&&\hat{z}\propto {\boldsymbol{L}},\end{eqnarray} \tag{ 2 }$$

while primed vectors are measured with respect to the primary spin vector, such that

$$\begin{eqnarray}&&\hat{z}^{\prime} \propto {\vec{\chi }}_{1}.\end{eqnarray} \tag{ 3 }$$

Hence, ${\vec{\chi }}_{1}$ is given by

$$\begin{eqnarray}&&{\vec{\chi }}_{1}={\chi }_{1}\left(\begin{array}{c}\sin {\theta }_{1}\cos {\phi }_{1}\\ \sin {\theta }_{1}\sin {\phi }_{1}\\ \cos {\theta }_{1}\end{array}\right)\end{eqnarray} \tag{ 4 }$$

in the unprimed coordinate system, and

$$\begin{eqnarray}&&\vec{\chi }{{\prime} }_{1}={\chi }_{1}\left(\begin{array}{c}0\\ 0\\ 1\end{array}\right)\end{eqnarray} \tag{ 5 }$$

in the primed coordinate system. A schematic diagram depicting the primed coordinate system can be seen in Figure 1.

Zoom In Zoom Out Reset image size

Here, θ₁ is the zenith angle between L at formation and ${\hat{\chi }}_{1}$, and ϕ₁ is the azimuthal angle measured from $\hat{x}$, about $\hat{z}^{\prime}$ at formation. The secondary spin vector is

$$\begin{eqnarray}&&\vec{\chi }{{\prime} }_{2}={\chi }_{2}\left(\begin{array}{c}\sin {\theta }_{12}\cos {\zeta }_{12}\\ \sin {\theta }_{12}\sin {\zeta }_{12}\\ \cos {\theta }_{12}\end{array}\right),\end{eqnarray} \tag{ 6 }$$

where θ₁₂ is the zenith angle between the ${\vec{\chi }}_{1}$ and ${\vec{\chi }}_{2}$, and ζ₁₂ is the azimuthal angle measured from $\hat{x}^{\prime}$, about ${\vec{\chi }}_{1}$.

Our population model is framed in terms of θ₁ and θ₁₂, with at-formation angles between ${\boldsymbol{L}},{\hat{\chi }}_{1}$ and ${\hat{\chi }}_{1},{\hat{\chi }}_{2}$, respectively. The distributions of θ₁, θ₁₂, denoted with π(...), are conditional on the hyperparameters σ₁, σ₁₂, which encode the AGN physics:

$$\begin{eqnarray}&&\pi ({\theta }_{1},{\theta }_{12}| {\sigma }_{1},{\sigma }_{12}).\end{eqnarray} \tag{ 7 }$$

The hyperparameters determine the shape of the distribution.

Using the parameterization from Talbot & Thrane (2017), we assume that the cosine of the primary spin zenith angle $\cos {\theta }_{1}$ at formation is drawn from a truncated normal distribution (denoted ${{ \mathcal N }}_{{\rm{t}}}$) with mean = 1 and standard deviation σ₁:

$$\begin{eqnarray}&&\pi (\cos {\theta }_{1}| {\sigma }_{1})={{ \mathcal N }}_{{\rm{t}}}(\cos {\theta }_{1}| {\sigma }_{1}).\end{eqnarray} \tag{ 8 }$$

This distribution allows for preferred alignment between ${\hat{\chi }}_{1}$ and L , with a free parameter σ₁ controlling the typical misalignment angle. For small values of σ₁, the ${\hat{\chi }}_{1}$ distribution tends to be nearly aligned with the L . As σ₁ becomes large, the distribution becomes uniform, so that ${\hat{\chi }}_{1}$ is uncorrelated with L . We assume that the primary azimuthal spin angle ϕ₁ is drawn from a uniform distribution denoted U.

We assume that the secondary spin vector is preferentially aligned to the primary spin vector ⁹ at formation, such that

$$\begin{eqnarray}&&\pi (\cos {\theta }_{12}| {\sigma }_{12})={{ \mathcal N }}_{{\rm{t}}}(\cos {\theta }_{12}| {\sigma }_{12}).\end{eqnarray} \tag{ 9 }$$

Small values of σ₁₂ imply that ${\vec{\chi }}_{1}$ and ${\vec{\chi }}_{2}$ tend to point in nearly the same direction. As σ₁₂ becomes large, the directions of ${\vec{\chi }}_{1}$ and ${\vec{\chi }}_{2}$ become uncorrelated. We assume that ζ₁₂ is drawn from a uniform distribution.

Putting everything together, the prior for the spin vector orientations at formation is given by

$$\begin{eqnarray}\begin{array}{rcl}\pi ({\theta }_{1},{\theta }_{12}| {\sigma }_{1},{\sigma }_{12}) & = & {{ \mathcal N }}_{{\rm{t}}}(\cos {\theta }_{1}| {\sigma }_{1})\ {{ \mathcal N }}_{{\rm{t}}}(\cos {\theta }_{12}| {\sigma }_{12})\\ & & U({\phi }_{1})\ U({\zeta }_{12}),\end{array}\end{eqnarray} \tag{ 10 }$$

where U(ϕ₁) and U(ζ₁₂) are constants equal to 1/2π. Figure 2 displays distributions of $\cos {\theta }_{12}$ for different values of σ₁₂. In the following section, we delve into various configurations of σ₁ and σ₁₂, and their physical implications on AGNs.

Zoom In Zoom Out Reset image size

3.2. The Evolution of Spin Vectors with Time

As BBH systems evolve, the coupling of the component BH spin vectors and the orbital angular momentum vector results in Lense–Thirring precession, causing the spin and angular momenta to precess about the system's total angular momentum vector (Mashhoon et al. 1984). Thus, the angles θ₁, θ₁₂ evolve over time. While the distributions of θ₁, θ₁₂ modeled in Equation (10) describe the binary properties at formation (approximately when binary is infinitely separated; Johnson-McDaniel et al. 2021), these distributions are different by the time these binaries enter the observing band of audio-band gravitational-wave detectors (20–20,000 Hz).

Fortunately, the information present in the distribution of θ₁, θ₁₂ at formation is encoded in the distribution of the "effective inspiral spin parameter" χ_eff (Damour 2001) and the effective precession parameter χ_p (Schmidt et al. 2012), which are approximate constants of motion (Hannam et al. 2014; Gerosa et al. 2021). The parameter χ_eff measures the spin components aligned with the orbital angular momentum, while χ_p measures the spin components in the orbital plane.

Figure 3 shows joint distributions of χ_eff and χ_p , representing three different AGN populations, each with different values of σ₁, σ₁₂. The dashed orange distribution is created using σ₁ = 0.1 and σ₁₂ = 10 (a mid-aged, dense AGN). The solid purple distribution displays the distribution for σ₁ = 0.1 and σ₁₂ = 0.1 (an old, dense AGN). Finally, the solid green distribution is the distribution with σ₁ = 10 and σ₁₂ = 10 (a young, dilute AGN). The distinguishability of these distributions illustrates how the AGN properties my be imprinted on the distribution of quantities measured by LIGO–Virgo.

Zoom In Zoom Out Reset image size

In order to recast our model in terms of the quantities that are measured by LIGO–Virgo, it is necessary to compute

$$\begin{eqnarray}&&\pi ({\chi }_{\mathrm{eff}},{\chi }_{p}| {\sigma }_{1},{\sigma }_{12}).\end{eqnarray} \tag{ 11 }$$

In principle, an expression for this distribution may be obtained through a series of convolutional integrals, which most likely have to be evaluated numerically. An alternative to numerical integration is to estimate the probability density distribution at fixed values of σ₁, σ₁₂, using either histograms or kernel density estimators. These estimates can be used to interpolate the probability density for arbitrary values of σ₁, σ₁₂, by using a machine-learning algorithm, for example; see Hernandez Vivanco et al. (2020). Alternatively, a machine-learning algorithm can likely be trained to reproduce the results of the numerical integrals for different values of σ₁, σ₁₂.

However, as the focus of this work is to estimate how well LIGO–Virgo will be able to measure σ₁, σ₁₂, creating a machine-learning representation of Equation (11) goes beyond our present scope. Instead, in the demonstration that follows in Section 4, we pretend that the values of θ₁, θ₁₂ measured by LIGO–Virgo remain unchanged since the binary was formed. This is tantamount to assuming that the χ_eff, χ_p distribution encodes—without loss—all of the information in the distribution of θ₁, θ₁₂. Since there is likely to be some information loss, our results are overly optimistic. Nonetheless, judging by the distinguishability of the different populations shown in Figure 3, we believe this assumption still yields a ballpark estimate.

In Table 1, we summarize the implications of measuring various σ₁ and σ₁₂ on some parameters describing AGNs. Some of the implications in Table 1 are degenerate. For example, if AGN disks are relatively low density, on average ($\bar{\rho }\lt {10}^{-11}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}$), then their corresponding average lifetimes ${\overline{\tau }}_{\mathrm{AGN}}\gtrsim 5\,\mathrm{Myr}$ in order to torque ${\vec{\chi }}_{\mathrm{1,2}}$ into alignment with J _AGN. Moreover, the top left quadrant (small σ₁, σ₁₂ ≤ 1; Table 1(a)) is indistinguishable from a population of BBH systems assembled in the field binary (aligned spins). The bottom right quadrant (large σ₁, σ₁₂ > 1; Table 1(d)) is consistent with a dynamical assembly in a dense stellar environment such as a globular cluster (isotropic spins). Thus, depending on the nature of the BBH assembly in the AGN disks, this framework will not necessarily provide a useful means of testing the AGN hypothesis against a competing hypothesis. Rather, it is a means of probing AGN physics assuming the AGN hypothesis is true. Independent evidence, e.g., from electromagnetic counterparts, may be required to establish this premise.

In this section, we carry out a demonstration analysis to estimate the ability of Advanced LIGO–Virgo to measure AGN physics. We generate two simulated populations of 100 BBH signals, drawn from our population model described by Equation (10). The BBH systems are uniformly distributed in comoving volume. We employ the mass model from Talbot & Thrane (2018), with parameters consistent with the results from Abbott et al. (2021). For the sake of simplicity, we assume fixed values for the mass ratio q = 1 and dimensionless spin magnitude χ₁ = 0.6, χ₂ = 0.6. These simplifying assumptions roughly match the mass ratio and spin magnitudes of the 20% of merging BBH systems with nonzero spin (Roulet et al. 2021; Galaudage et al. 2021). We additionally use an arbitrary polarization angle ϕ = 0.1 to further simplify the analysis. We draw the BH spin tilts for each population from

$$\begin{eqnarray*}&&\begin{array}{l}\,\mathrm{Population}\ {\rm{A}}:\ {\pi }_{\mathrm{AGN}}({\vec{\chi }}_{1},{\vec{\chi }}_{2}\,| \ {\sigma }_{1}=0.5,{\sigma }_{12}=1),\\ \,\mathrm{Population}\ {\rm{B}}:\ {\pi }_{\mathrm{AGN}}({\vec{\chi }}_{1},{\vec{\chi }}_{2}\,| \ {\sigma }_{1}=3.0,{\sigma }_{12}=0.25),\end{array}\end{eqnarray*}$$

at a reference frequency of 20 Hz. ¹⁰ As discussed in Section 3.2, we ignore the evolution of the spin vectors from the time of formation to the moment they enter the LIGO–Virgo band, which is equivalent to assuming that the information in the θ₁, θ₁₂ distribution is encoded in the distribution of χ_eff, χ_p without any data loss. This assumption makes our results overly optimistic. Population A corresponds to mergers from AGNs that are old and dense (Table 1(a)). On the other hand, Population B corresponds to mergers from AGNs that are old and dilute (Table 1(c)).

We simulate the gravitational-wave signals from the BBH mergers using the waveform approximant IMRPhenomXPHM (Pratten et al. 2020, 2021; García-Quirós et al. 2020). We add the simulated signals into Gaussian noise, colored to the Advanced LIGO design sensitivity for the Hanford and Livingston detectors (LIGO Scientific Collaboration et al. 2015; Abbott et al. 2018). We ensure that the matched filter signal-to-noise ratio (S/N) ≥ 8. We perform Bayesian inference with parallel-bilby (Skilling 2004; Ashton et al. 2019b; Smith et al. 2020; Speagle 2020; Romero-Shaw et al. 2020) and GWPopulation (Talbot et al. 2019) to recover posterior probability densities for the parameters of the simulated signals using the same waveform approximant IMRPhenomXPHM, also evaluated at a reference frequency of 20 Hz. The prior and population probability distributions we use are documented in Table 2.

Table 2. Population and Prior Distributions for Demonstration

Population		Prior		Common
Parameter	Distribution	Parameter	Distribution	Parameter	Distribution
$| {\vec{\chi }}_{1}| ,| {\vec{\chi }}_{2}|$	δ(0.6)	$| {\vec{\chi }}_{1}| ,| {\vec{\chi }}_{2}|$	δ(0.6)	${ \mathcal M }/{M}_{\odot }$	U(15,60)
$\cos {\theta }_{1}$	${{ \mathcal N }}_{t}({\sigma }_{1})$	$\cos {\theta }_{1}$	U(−1, 1)	q, ϕ	δ(1), δ(0.1)
$\cos {\theta }_{12}$	${{ \mathcal N }}_{t}({\sigma }_{12})$	$\cos {\theta }_{2}$	U(−1, 1)	dL/Mpc	CM(200, 800)
η12, ϕ1	U(0,2π)	ϕ12, ϕJL	U(0,2π)	decl.	C(0, 2π)
ι	U(0, π)	θJN	S(0, 2π)	R.A., ψ	U(0, 2π)

Note. We use shorthand to represent the distributions: δ, delta; U, uniform; CM, uniform in comoving volume; C, cosine; and, finally, ${{ \mathcal N }}_{t}$, truncated normal (mean of 1). The definitions of the parameters are documented in Romero-Shaw et al. (2020; Table E1). The "Population" and "Prior" columns display the spin orientation distributions for the population and prior. The "Common" column displays some parameter distributions that are common between the population and prior.

Download table as:  ASCIITypeset image

In Figure 4, we plot cumulative distributions of $\cos {\theta }_{1},\cos {\theta }_{12}$ for Population A (old and dense; blue) and Population B (old and dilute; red). The true population distributions are the solid dark curves. The empirically observed distributions from the simulated catalogs are the shaded bands. The darker regions of the shaded bands indicate the 90% credibility ranges of the observed distributions, while the lighter regions indicate the 99% credibility ranges. Each row of the figure displays the observed distribution for a different population size. As the number of events in each population increases, the observed 99% credibility ranges shrink. Figure 4 demonstrates that the two populations can be visually distinguished once there are ${ \mathcal O }(10)$ gravitational-wave events. However, given recent results that suggest only ≈20% of BBH systems contain a BH with measurable spin (Miller et al. 2020; Galaudage et al. 2021; Roulet et al. 2021)—and given the loss of information moving from θ₁, θ₁₂ to χ_eff, χ_p —it is likely that $\gtrsim { \mathcal O }(50)$ events are required.

Zoom In Zoom Out Reset image size

This paper introduces a physically motivated phenomenological population model describing the at-formation spin orientations of merging BBHs assembled in AGNs. By measuring the distribution of the spin orientations with gravitational waves, we may be able to learn about the AGN environment; whether it is old or young, whether it is dilute or dense. We demonstrate that ${ \mathcal O }(10)$ gravitational-wave events from BBH mergers with S/N ≥ 8 and spin magnitudes χ ∼ 0.6 are required to infer the population parameters describing the shape of the $\cos {\theta }_{1}$ and $\cos {\theta }_{12}$ distributions. However, since the majority of BBH mergers include BHs with negligible spin (Miller et al. 2020; Roulet et al. 2021; Galaudage et al. 2021)—and since we have optimistically assumed that the distribution of χ_eff, χ_p retain all the information of the spin orientations at formation—a total of $\gtrsim { \mathcal O }(50)$ BBH detections are necessary to begin to resolve the properties of AGNs.

To derive this estimate, we assume that all BBH mergers take place in AGNs. This is a reasonable starting point, since it is desirable to see whether all BBH detections can be understood within a single channel, and it would be somewhat surprising if it turns out that two different channels produce comparable merger rates. However, it may turn out that AGNs provide only some fraction of the observed population. If so, our effort to infer the properties of AGN disks will be complicated by contamination from other channels. One could use the model proposed here as a subpopulation in a mixture model, but this goes beyond our present scope.

Our model is cast in terms of the orientation of spin vectors at the time of BBH formation. While the spin vectors subsequently evolve due to precession, the initial orientation of the spin vectors is imprinted on the effective spin parameters χ_eff, χ_p ; see Figure 3. In order to further develop this model such that it can be used for Bayesian inference, it will be necessary to recast the model in terms of these effective parameters (see Equation (11)). Future work will focus on the development of a computationally efficient representation of this distribution, by using a machine-learning algorithm, for example. With a computationally efficient model, it will be possible to apply the model to current gravitational-wave catalogs that may now have enough events to begin resolving the properties of AGN physics.

We thank Shanika Galaudage and Colm Talbot for their technical help. We thank Ilya Mandel for helpful advice.

We gratefully acknowledge the Swinburne Supercomputing OzSTAR Facility for computational resources. All analyses (inclusive of test and failed analyses) performed for this study used 130k core hours on OzSTAR. This would have amounted to a carbon footprint of ∼ 7.6 t CO₂ (Australian Government—Department of the Environment and Energy 2021). However, as OzSTAR is powered by wind energy from Iberdrola Australia, the electricity for computations produces negligible carbon waste.

This material is based upon work supported by NSF's LIGO Laboratory, a major facility fully funded by the National Science Foundation. This research has used data, software, and web tools obtained from the Gravitational Wave Open Science Center (https://www.gw-openscience.org), a service of the LIGO Laboratory, the LIGO Scientific Collaboration, and the Virgo Collaboration. The LIGO Laboratory and Advanced LIGO are funded by the United States National Science Foundation (NSF), as well as the Science and Technology Facilities Council (STFC) of the United Kingdom, the Max Planck Society (MPS), and the State of Niedersachsen/Germany, to support the construction of Advanced LIGO and the construction and operation of the GEO600 detector. Additional support for Advanced LIGO was provided by the Australian Research Council. Virgo is funded through the European Gravitational Observatory (EGO), by the French Centre National de Recherche Scientifique (CNRS), the Italian Istituto Nazionale di Fisica Nucleare (INFN), and the Dutch Nikhef, with contributions by institutions from Belgium, Germany, Greece, Hungary, Ireland, Japan, Monaco, Poland, Portugal, and Spain.

A.V., R.S., and E.T. are supported by the Australian Research Council (ARC) Centre of Excellence CE170100004. B.M. and K.E.S.F. are supported by NSF AST-1831415 and Simons Foundation Grant 533845.

Software: parallel-bilby (Smith & Ashton 2021, v1.0.1), bilby (Ashton et al. 2019a, v1.1.4), bilby-pipe (Ashton et al. 2020, v1.0.4), GWPopulation (Talbot et al. 2019, v0.6.3), dynesty (Speagle 2020, v1.0.1), GWpy (Macleod et al. 2020, v2.1.0), LALSimulation (LIGO Scientific Collaboration 2020, v7.0), matplotlib (Hunter 2007, v3.2.0), NumPy (Harris et al. 2020, v1.8.1), SciPy (Virtanen et al. 2020, v1.4.1), pandas (Pandas Development Team 2020, v1.0.2), python (Oliphant 2007; Millman & Aivazis 2011, v3.7).