Subparsec supermassive black hole SMBH binaries are expected to be common in active galactic nuclei AGN as a result of the hierarchical buildup of galaxies via mergers. While direct evidence for these compact binaries is lacking a few hundred candidates have been identified most based on the apparent periodicities of their optical lightcurves. Since these signatures can be mimicked by AGN rednoise additional evidence is needed to confirm their binary nature. Recurring selflensing flares SLF occurring whenever the two BHs are aligned with the line of sight within their Einstein radii have been suggested as additional binary signatures. Furthermore in many cases lensing flares are also predicted to contain a dip whenever the lensed SMBHs shadow is comparable in angular size to the binarys Einstein radius. This feature would unambiguously confirm binaries and additionally identify SMBH shadows that are spatially unresolvable by highresolution VLBI. Here we estimate the number of quasars for which these dips may be detectable by LSST by extrapolating the quasar luminosity function to faint magnitudes and assuming that SMBH binaries are randomly oriented and have massratios following those in the Illustris simulations. Under plausible assumptions about quasar lifetimes binary fractions and Eddington ratios we expect tens of thousands of detectable flares of which several dozen contain measurable dips.
['Astrophysics - High Energy Astrophysical Phenomena']
Selflensing flares from black hole binaries IV the number of detectable shadows
{'Self-lensing flares from black hole binaries IV: the number of detectable shadows': "Kevin Park 1 , ∗ Chengcheng Xin 2 , Jordy Davelaar 3 , 4 , 5 , 2 , 6 , and Zolt´an Haiman 2 , 1  \n1 Department of Physics, Columbia University, New York, NY 10027, USA  \n2 Department of Astronomy, Columbia University, New York, NY 10027, USA  \n3 Department of Astrophysical Sciences, Peyton Hall,  \nPrinceton University, Princeton, NJ 08544, USA  \n4 NASA Hubble Fellowship Program, Einstein Fellow  \n5 Center for Computational Astrophysics, Flatiron Institute,  \n162 Fifth Avenue, New York, NY 10010, USA and  \n6 Astrophysics Laboratory, Columbia University, 550 W 120th St, New York, NY 10027, USA  \nSub-parsec supermassive black hole (SMBH) binaries are expected to be common in active galactic nuclei (AGN), as a result of the hierarchical build-up of galaxies via mergers. While direct evidence for these compact binaries is lacking, a few hundred candidates have been identified, most based on the apparent periodicities of their optical light-curves. Since these signatures can be mimicked by AGN red-noise, additional evidence is needed to confirm their binary nature. Recurring self-lensing flares (SLF), occurring whenever the two BHs are aligned with the line of sight within their Einstein radii, have been suggested as additional binary signatures. Furthermore, in many cases, lensing flares are also predicted to contain a 'dip', whenever the lensed SMBH's shadow is comparable in angular size to the binary's Einstein radius. This feature would unambiguously confirm binaries and additionally identify SMBH shadows that are spatially unresolvable by high-resolution VLBI. Here we estimate the number of quasars for which these dips may be detectable by LSST, by extrapolating the quasar luminosity function to faint magnitudes, and assuming that SMBH binaries are randomly oriented and have mass-ratios following those in the Illustris simulations. Under plausible assumptions about quasar lifetimes, binary fractions, and Eddington ratios, we expect tens of thousands of detectable flares, of which several dozen contain measurable dips.", 'I. INTRODUCTION': "Supermassive black holes (SMBHs) with masses between M ≈ 10 6 -10 9 M ⊙ are present in the nuclei of most nearby galaxies [1]. In hierarchical cosmologies, galaxies grow by frequent mergers, which deliver the nuclear SMBHs [2, 3], along with significant quantities of gas [4], to the central regions of the post-merger galaxy. The natural conclusion is that after the two SMBH's separation decreases, SMBH binaries (SMBHBs) should form frequently in galactic nuclei [5]. Hydrodynamical simulations have shown that electromagnetic (EM) emission from these SMBHBs, provided they are surrounded by circumbinary gas, should be detectable starting well before the merger, and should persist all the way to the merger [6-9].  \nThese compact SMBHBs are a fundamental ingredient of galaxy formation and are also prime targets to be observed in gravitational waves (GWs) by LISA [10, 11], and by pulsar timing arrays (PTAs). Indeed, PTAs have recently discovered a stochastic GW background (GWB) in the nHz bands, which is consistent with the cosmological population of coalescing SMBHBs [12-15]. Combining the EM and GW signals from the same source or even the same or overlapping populations of SMBHs would open windows to especially novel science, including astrophysics, cosmology, particle physics, and the nature of gravity [16].  \nThese so-called 'multimessenger' opportunities have stimulated strong interest in finding wider SMBH binaries in EM data. Approximately 300 SMBH binary candidates have been identified among bright active galactic nuclei (AGN) in large time-domain optical surveys [1719] based on their apparent periodicities, and a handful of additional SMBH binary candidates have been identified serendipitously, or through other tentative signatures involving double-peaked or offset emission lines, or spatial structures of radio jets and lobes [20, 21]. These candidates remain controversial because of the lack of a 'smoking gun' binary signature, and in the cases of the periodic candidates, stochastic red-noise AGN variability can mimic periodicities for a few periods [22].  \nApotential SMBHB signature that in some cases could help lift this degeneracy is a 'self-lensing flare' (SLF). If the two SMBHs are aligned within the line-of-sight to within the system's Einstein radius, then whenever one of the SMBHs passes behind the other, its emission will be strongly magnified. These lensing flares occur once or twice per orbit (depending on whether one or both BHs are active). Depending on the binary masses and separations, the flares can last from hours to weeks [23-26]. An AGN identified in the Kepler catalog, dubbed Spikey, has a light-curve consistent with relativistic Doppler modulation from an eccentric binary, with a narrow ∼ 10% spike at the expected orbital phase whose symmetric shape is well fit by a microlensing model [27].  \nIn a fraction of these self-lensing binaries, the size of the SMBH's shadow ( ∼ 5 times its gravitational radius) is commensurate with the (angular) Einstein radius, and an additional feature is imprinted on the light-curves, in  \nthe shape of a 'dip' near the peak of each flare. This feature could unambiguously confirm binaries and additionally identify SMBH shadows that are spatially unresolvable by high-resolution VLBI, as was done by the Event Horizon Telescope Collaboration [28, hereafter Paper I]. Toy models for the BH emission show that the precise shapes and sizes of these dips depend on the binary system's parameters [29, Paper II], and recent hydrodynamical simulations find that these 'dip' features exist even in the strongly distorted and fluctuating circumbinary gas, and can be recovered via phase-folding in the face of stochastic noise [30, III]. [31] perform a similar study of the ray-traced emission from SMBH binaries, but using a boosted binary metric and following the binary's inspiral with 3.5PN equations of motion. They recover the self-lensing flares and dips closely matching those in Paper I.  \nSelf-lensing flares require nearly edge-on viewing angles, and a natural question is how rare these lensed configurations are. The Vera Rubin Observatory's LSST [32] is expected to contain between 20-100 million bright quasars [33], and combined with its high cadence (with photometric points every few days), is an ideal dataset in which to search for rare sources with periodically recurring flares. Kelley et al. [34] have recently used SMBH populations from the Illustris simulation [35], combined with toy models for binary emission and lensing, and concluded that LSST could detect several hundred selflensed SMBHs with flares lasting for 30 days or longer.  \nIn this paper, we follow up on the above studies, to assess the number of self-lensing flares which additionally have detectable 'dips' in their light-curves. In principle, these dips require more stringent alignment, but their incidence rate rises steeply for shorter-duration flares. Here we estimate the number of quasars for which these dips may be detectable by LSST, by extrapolating the quasar luminosity function to faint magnitudes, and assuming that SMBH binaries are randomly oriented and have mass ratios following those in the Illustris simulations. Under plausible assumptions about quasar lifetimes, binary fractions, and Eddington ratios, we predict tens of thousands of detectable flares with durations down to ten days, of which several dozen contain measurable dips.  \nThe rest of this paper is organized as follows. In § II, we describe our methodology, including models for the SMBH binary populations, their emission and lensing, the criteria for detectable flares and dips, and the impact of finite source sizes. In § IV, we present our results, in terms of the number of flares and dips detectable in LSST as a function of binary parameters and observational thresholds (magnitude and flare duration). In § IV, we summarize our main results and their implications.", 'II. METHODS': 'In this section, we calculate the number of detectable self-lensing dips in LSST and its dependence on several  \nFigure 1. The number of quasars N QSO ( M,z ) in the redshift and mass range z ∈ [0 , 6] and M/ M ⊙ ∈ [10 5 , 10 9 . 5 ]. Integrating over the full mass and redshift ranges gives 20 million, 44 million, 100 million quasars above the assumed LSST magnitude detection limits of m i = 24 , 25 , 26, in agreement with the results of [33].  \n<!-- image -->  \nbinary parameters. Results from [36] demonstrate that the evolution of quasars can be reproduced in a model in which they are activated in galaxy mergers. As mentioned above, these galaxy mergers are expected to deliver the two SMBHs to the new galactic nucleus where they form a bound binary [5]. Given these results, our main assumption is that galaxy mergers are responsible for both quasar activity and for producing SMBH binaries [37]. First, we find the expected number of quasars ( N QSO ) in LSST using the quasar luminosity function ( § II A). We then modify N QSO to find the number of binary quasars as a function of the binary orbital period, using estimates of the quasar lifetime at each orbital period ( § II B) allowing the overall binary fraction to be a free parameter. We find the mass ratio distribution for these binaries in different mass and redshift bins using the Illustris simulations ( § II C) and we compute the probability of a detectable self-lensing dip given binary inclination ( § II D) as a function of the mass, mass ratio and orbital period. Finally, we compute similar probabilities for self-lensing flares in the point-source limit ( § II E) and accounting for the finite sizes of the emitting regions ( § II F).', 'A. Number of quasars and binaries': "We follow the calculations in [33] to obtain the number of quasars above LSST's detection threshold, based on the extrapolated quasar luminosity function [QLF; 38]. In Figure 1, we present the number of quasars in LSST's 20,000 deg 2 survey, N QSO , for three magnitude limits, m i = 24, 25 and 26, where m i = 24 corresponds to the  \nsingle-exposure magnitude limit of LSST in the i band, and m i = 26 corresponds to the co-added magnitude limit over the whole survey [39]. The BH mass and i -band magnitudes are related using Eq. 1:  \nm i = 24+2 . 5 log [ ( f Edd 0 . 3 ) -1 ( M 3 · 10 6 M ⊙ ) -1 ( d L ( z ) d L ( z = 2) ) 2 (1)  \nwhere f Edd is defined by the bolometric quasar luminosity L = f Edd L Edd and L Edd is the Eddington luminosity for total mass M . Initially, we assume f Edd = 0 . 3 , but in § IV below we will discuss the dependence of our results to varying f Edd . Quasars in the hatched area are discarded since their magnitudes are below m i = 26.  \nUnder the assumption that all quasars are triggered by mergers and are associated with SMBH binaries [40, 41], and that the bright quasar phase lasts for a typical lifetime of τ Q (=say 10 8 years), the number of quasars powered by SMBH binaries with τ m years left to merger scale linearly with the fraction τ m /τ Q ,  \nN ( M,z ) ≡ ( τ m τ q ) f bin N QSO ( M,z ) . (2)  \nN QSO ( M,z ) is the number of quasars shown in Figure 1, which assumes that all quasars are associated with binaries. In this analysis, we vary the fraction of quasars that are SMBH binaries, f bin , from 0.2 to 1. We use residence times τ m = τ res ( q, M, T orb , z ) appropriate to GWor gas-driven binary inspirals, depending on the mass M and orbital period T orb of the binary (see next section). By setting a maximum period T max , N ( M,z ) represents the number of binaries with an orbital period of T max or less, equivalent to the number of binaries with τ m years left to merger.", 'B. Residence time vs. orbital period': "We obtain the residence time of a binary quasar at a given orbital period (or separation) τ res as follows. For short periods, the orbital decay is primarily GW-driven, given by the quadrupole formula [42]:  \nτ GW = 1 . 11 × 10 7 q -1 s M -5 / 3 7 ( T orb yr ) 8 / 3 ( 1 + z 4 ) yr , (3)  \nwhere q s ≡ 4 q/ (1 + q ) 2 is the symmetric mass ratio, q ≡ M 2 /M 1 ≤ 1 is the mass ratio, M 7 ≡ M/ (10 7 M ⊙ ) is the total mass in units of 10 7 solar masses, and T orb is the orbital period in units of years.  \nFor long periods, we assume that the circumbinary gas dominates the binary's inspiral, and adopt the simple scaling from [43]. This assumes that the gas-driven inspiral timescale is proportional to the mass accretion timescale, roughly in agreement with the results of hydrodynamical simulations [44-49]. First, the mass accretion timescale is referenced in units of the Eddington  \n]  \nFigure 2. The solid blue curve shows the fiducial model for the residence time τ res ( q, M, T orb , z ) of binaries with total mass M = 10 7 M ⊙ , mass ratio q = 1 and redshift z = 1 as a function of observed orbital period T orb . Below the threshold given in Eq. 5 the binaries evolve via GW emission and above the threshold the binaries' inspiral timescale is assumed to be driven by circumbinary gas and increase linearly with T orb . In § IV we also present results for a more conservative model where a constant maximum lifetime of a/ ˙ a gas = 1 . 7 × 10 7 years is imposed (Eq.6), shown by the solid orange curve. Purely GW-driven inspiral, the dashed blue curve, is shown for reference.  \n<!-- image -->  \nrate:  \n˙ m Edd = 4 πGmm p ησ T c = 2 . 26 · 10 -2 × ( η 0 . 1 ) -1 × ( m 10 6 M ⊙ ) M ⊙ yr -1 , (4)  \nwhere G is the gravitational constant, m is the mass of the accreting BH, m p is the proton mass, c is the speed of light, σ T is the Thomson scattering cross section and η is the radiative efficiency, assumed to be 0 . 1. Assuming an Eddington accretion rate ˙ m = ˙ m Edd , then  \nm ˙ m ∼ 4 . 4 × 10 7 yr , (5)  \nand we follow [43] and adopt the residence time at binary separation a  \nτ gas = a ˙ a gas = 1 2 . 68 m ˙ m ∼ 1 . 7 × 10 7 yr , (6)  \nwhich is constant (independent of a ). At large separations for which the GW-inspiral timescale exceeds this threshold, we assume as our fiducial model that the binary gas-driven residence time increases linearly with the orbital period, as shown by the solid blue curve in Figure 2. As an alternative model, the residence time is fixed at a maximum of 1 . 7 × 10 7 years, shown by the solid orange curve. The illustrations in Figure 2 are for a  \nFigure 3. The left panel shows the distribution of the mass ratios of SMBH binaries in the Illustris simulations. In the middle and right panels, we illustrate the redshift and total mass distributions of Illustris binaries, in the four different mass-ratio ranges shown by colors in the left panel. Most binaries are symmetric in mass ( q > 0 . 1) and have relatively low redshift ( z < 3) and total mass ( M < 10 7 M ⊙ ).  \n<!-- image -->  \nbinary with mass M = 10 7 M ⊙ , redshift z = 1 and mass ratio q = 1. Finally, we constrain the maximum observed orbital period to T max = 5 years, given that for LSST's survey of 10 years, we want to observe a few recurrences of the dip.", 'C. Mass-ratio distributions': "The GW-driven inspiral timescale in Eq. 3 depends on the binary's mass ratio q , which is not constrained by the QLF or other observations. To estimate the mass-ratio distribution we instead use Illustris-3 [35], a cosmological hydrodynamical simulation that self-consistently follows the evolution and mergers of galaxies and their central SMBHs. Illustris provides the merger tree of their SMBH mergers, and we use this to compute to distribution of q as a function of total mass M and redshift z , shown in Figure 3. Illustris embeds MBHs of seed mass ∼ 10 5 M ⊙ which accrete and evolve dynamically. Due to the limitations to spatially resolve closely separated low-mass binaries ( M < 10 6 M ⊙ ), we initially implement a conservative mass cut of 10 6 M ⊙ of Illustris BH binaries and extrapolate the count distribution in ( M,z ) down to the binary total mass of 10 5 M ⊙ . We then divide the Illustris binaries by their mass ratios into 8 logarithmic bins between 10 -4 ≤ q ≤ 1. Figure 3 depicts the fraction of Illustris binaries in each of these 8 bins (left panel). The colors blue, yellow, green, and red each correspond to the binaries with mass ratios in the ranges of [10 -1 , 1] , [10 -2 , 10 -1 ] , [10 -3 , 10 -2 ] , [10 -4 , 10 -3 ]. Using Illustris, we calculate the probabilities that a BH binary with given M and z in one of 50 bins has a mass ratio q in one of the 8 logarithmic bins-see the middle and right panel of Figure 3. Using these mass ratio distributions we  \ncan evaluate the GW-driven inspiral timescale for binary quasars in the QLF.", 'D. Self-lensing dips': "As mentioned in § I, for SMBHBs with nearly edgeon orbital planes, periodic self-lensing flares occur as a result of gravitational lensing. Compact binaries close to merger, with an orbital period of 5 years or less in our fiducial model (Table I), are expected to have four or more self-lensing flares detectable within a full ten-year LSST survey. Additionally, GR ray tracing simulations of [28] reveal that self-lensing flares have observable dips, caused by the black hole shadow. In this work, we estimate the number of LSST binaries with detectable selflensing dips for different binary parameters. For this, we recap Eqs. 1-3 of [28], which give an analytical expression for the probability that a binary system with given binary mass ratio q , total mass M , orbital period T orb and redshift z has a detectable self-lensing dip.  \nFirst, assuming a circular binary with nearly edge-on binary inclination, the expected phase spacing between the two peaks before and after the dip is the ratio of the diameter of the BH shadow and the circumference of the orbit:  \n∆ ϕ = d shadow 2 πa orb , (7)  \nwhere a orb is the orbital radius, d shadow = 2 √ 27 GM source /c 2 is the BH shadow diameter of the source, and M source = qM/ (1 + q ) is the mass of the lensed BH (assumed here to be the lower-mass secondary). Expressing the orbital radius in terms of  \nTable I. Parameters in our fiducial model (top row) and their ranges considered (bottom row). f bin is the fraction of quasars associated with binaries, f Edd is their Eddington ratio, T max is the maximum orbital period of interest, τ Q is the average total bright quasar lifetime, τ res is the residence time, i.e. the duration a binary quasar spends at each orbital period (shown in Fig. 2), and τ f, min is the minimum required lensing flare duration (not imposed in the fiducial model).  \nbinary parameters via a orb = ( GMT 2 orb / 4 π 2 ) 1 / 3 gives  \n∆ ϕ = 5 . 63 qM 2 / 3 (1 + q ) T 2 / 3 orb radians , (8)  \nwhere G = c = 1 units were used. T orb has units R g /c , where R g = GM source /c 2 is the gravitational radius of M source . The dip in the self-lensing flare can be observed for a range of binary inclinations, ∆ i . When ∆ i is smaller than the angular size of the BH shadow on the sky, i.e. when Eq. 9 is satisfied, the focal point of the lens will be eclipsed by the BH shadow, causing dips at the local maxima of the flares. The inclination window for which this dip is visible is given by  \n∆ i ≤ sin -1 ( π ∆ ϕ ) ≡ ∆ i dip . (9)  \nAssuming that BH binary orbital inclinations are randomly distributed on a unit sphere, the minimum angular separation between the BHs on the sky occurs when the two BHs are at either end of the semi-minor axes of the orbit's projected ellipse. Then, for a binary chosen at random, the probability that it has an observable self-lensing dip during its orbit is given by  \nP dip = ∆ i dip ( q, T orb , M, z ) / 90 · , (10)  \ni.e. the probability that ∆ i , which can range from 0 · to 90 · , satisfies Eq. 9. The number of detectable selflensing dips is then N dips = P dip × N ( M,z ), integrated appropriately over the distribution of all binary parameters of mass ratio, observed orbital period, total mass and redshift.", 'E. Self-lensing flares - point source limit': 'We wish to compare our calculations with [34], who have calculated the number of detectable self-lensing flares, irrespective of whether dips from the BH shadow  \nare measurable, based on the abundance of binaries in Illustris. To do this, we construct a self-lensing probability similar to the self-lensing dip probability in Eq. 10. We adopt the point-source (PS) magnification [50]  \nM PS = u 2 +2 u √ u 2 +4 , (11)  \nwhere u = Re( u 1 + u 2 ) is the projected separation in units of the Einstein radius and  \nu j = [ ac 2 (cos 2 ϕ j +sin 2 i sin 2 ϕ i ) 4 G ( M -m j ) cos i sin ϕ j ] 1 / 2 , (12)  \nwhere i is again the binary inclination relative to the line of sight and ϕ j is the orbital phase of each BH. At the peak of the self-lensing flare ( ϕ 2 = π/ 2), we assume the secondary BH to be the source, and find the maximum orbital inclination ∆ i PS in which a 10% magnification occurs, i.e. M PS ( i ) > 1 . 1. The corresponding self-lensing probability is  \nP PS = ∆ i PS ( q, T orb , M, z ) / 90 · (13)  \nand the number of self-lensing flares in the point-source limit is then N PS = P PS × N ( M,z ), again integrated appropriately over the distribution of all binary parameters.', 'F. Self-lensing flares - finite source': "According to Figure 2 of [51], for binaries of total mass ≲ 10 8 M ⊙ and orbital periods of a few years, it is necessary to incorporate the accretion disc size of the secondary because the angular size of the accretion disc becomes comparable to the Einstein radius for lower total masses. We adopt their model to account for finitesource effects, based on a multi-color accretion disc extending from the innermost stable circular orbit (ISCO) r ISCO = 6 GM source /c 2 , to the secondary's tidal truncation radius r tidal = 0 . 27 q 0 . 3 a [52] The temperature profile of the disc is given by  \nσT 4 ( r ) = 3 GM source ˙ M source 8 πr 3 [ 1 -( r ISCO r ) 1 / 2 ] , (14)  \nwhere r is the distance from the central SMBH and ˙ M source is the accretion rate of the lensed BH, for which we assume an Eddington accretion rate. For r > r tidal and r < r ISCO , we set T ( r ) = 0 . The resulting flux from the disc is given by  \nF ν ( r ) = πB ν [ T ( r )] , (15)  \nwhere B ν is the Planck function. We calculate the lensing magnification by evaluating  \nM FS ν = ∫ 2 π 0 ∫ ∞ 0 F ν ( u ' , v ' ) M PS ( u ' ) u ' du ' dv ' ∫ 2 π 0 ∫ ∞ 0 F ν ( u ' , v ' ) u ' du ' dv ' , (16)  \nwhere M PS ( u ' ) is given in Eq. 11, evaluated in the lenscentered polar coordinates ( u, v ):  \nr 2 ∗ = ( u 2 0 + u 2 -2 u 0 u cos( v -v 0 )) r 2 E (17)  \nr = r ∗ √ cos 2 ϕ + sin 2 ϕ cos 2 ( π/ 2 -J ) (18)  \nsin ϕ = u sin v -u 0 sin v 0 √ ( u sin v -u 0 sin v 0 ) 2 +( u cos v -u 0 cos v 0 ) 2 . (19)  \nHere ( u 0 , v 0 ) is the position of the secondary, and J is the inclination of the source disc relative to the line of sight. We define ∆ i FS as the maximum orbital inclination such that M FS ν ( i ) > 1 . 1 at alignment and the probability that a binary has a detectable self-lensing flare in the finite source limit is  \nP FS = ∆ i FS ( q, T orb , M, z, ν, J ) / 90 · . (20)  \nThe number of self-lensing flares in the finite-source limit is N FS = P FS × N ( M,z ), integrated over binary parameters.  \nIn this paper, we consider the representative observed wavelengths near the center of the six LSST filters u, g, r, i, z, y at 380, 476, 622, 755, 870, 1015 nm, which are chosen by taking the average of the FWHM transmission points of each filter [39]. For example, for a binary at redshift z ∈ [0 , 6], the rest-frame wavelength of 380 / (1+ z ) nm will contribute flux to the u-band, and the rest-frame wavelength is used to evaluate the wavelengthdependent flux in Eq. 15.  \nAbsorption by hydrogen clouds in the intergalactic medium could attenuate quasar light depending on the observing wavelength and the quasar redshift. However, according to [53], for observed wavelengths of ν ∼ 400 nm, the mean cosmic transmission is close to 1 for z ≲ 3, where the majority of quasars are. For longer observed wavelengths of ν > 700 nm, the mean cosmic transmission is close to 1 for z ≲ 5, which means intergalactic absorption is mostly negligible for our purposes. In all calculations, we fix the source-disc inclination at a representative value of J = π/ 4 for simplicity.", 'A. Self-lensing dips and their mass- and redshift-dependence': 'We first present the number of detectable dips, for total SMBH binary masses and redshifts in the ranges of M/ M ⊙ ∈ [10 5 , 10 9 . 5 ] and z ∈ [0 , 6]. These mass- and redshift-distributions in the fiducial model are shown in Figure 4, and are also presented numerically in Tables II and III. Overall, we find 41-60 detectable dips, depending  \nFigure 4. The distribution of the number of detectable selflensing dips in the M -z plane. The contours for N dips are in units of per unit redshift per unit log mass bin, where the unit bin sizes are ∆ z = 0 . 12 and ∆log 10 M = 0 . 08. Three LSST sensitivity limits, covering the range from a single exposure to the fully co-added survey detection threshold, are shown for reference. Most detectable dips are from quasars at 1 . 5 ≲ z ≲ 3 with total (binary) SMBH masses of 10 6 M ⊙ ≲ M ≲ 10 8 M ⊙ .  \n<!-- image -->  \nTable II. Dependence of self-lensing dips on binary mass, integrated over all redshifts z ∈ [0 , 6] and mass ratios q ∈ [10 -4 , 1]. This result also assumes the fiducial parameters in Table I.  \non the magnitude threshold. Table II shows that ∼ 90% of these binaries with detectable dips have masses between 10 6 M ⊙ to 10 8 M ⊙ . Qualitatively, while the distribution of quasars (Fig. 1) is concentrated at low masses (below 10 6 M ⊙ ), the self-lensing probability in Eq. 10 increases with mass, resulting in the distribution clustering in the intermediate mass regime, see Table II.  \nTable III shows that approximately 80% of these binaries are between redshifts z = 1 -3, which is also expected considering that the lensing probability depends only weakly on redshift and the QLF peaks in this range.  \nFurthermore, the LSST limits on apparent magnitude, m i < 24 (or, optimistically, m i = 26) constrain the M -z parameter space, which is visualized in the contours of Figure 1. Therefore, we present all of our results below for three different LSST magnitude limits, m i = 24, 25 and 26, covering the range from a single exposure limit to  \nTable III. Redshift-dependence of self-lensing dips, integrated over all masses M/M ⊙ ∈ [10 5 , 10 9 . 5 ] and mass ratios q ∈ [10 -4 , 1]. The conditions for f bin and lensing duration were applied as for the total mass dependence in Table I.  \nTable IV. Mass-ratio dependence of self-lensing dips, integrated over all redshifts z ∈ [0 , 6] and total masses M ∈ [10 5 , 10 9 . 5 ] M ⊙ . This result also assumes the fiducial parameters in Table I.  \nthe fully co-added survey detection threshold. Figure 4 also shows these three magnitude cuts (blue curves), with quasars below m i = 26 are discarded as too faint to be detected (hatched region).', 'B. Mass-ratio dependence': 'There are no ab initio constraints on the Illustris mass ratios, but in practice we find all of them to be within q = 10 -4 -1, with almost all ( > 98%) binaries in the mass-ratio range of 0.1 to 1, as shown in Table IV. This follows from the mass ratio distribution of SMBH binaries from Illustris-the left panel of Figure 3. We note that these data from Illustris are quite uncertain, especially since all of the SMBH binary merger physics is sub-grid. In particular, mergers are assumed to take place instantly when their separation (modeled semi-analytically) decreases below a certain smoothing length. In reality, there can be delays and disruptions by a third SMBH, which can change the mass-ratio distribution [54], and in turn the number of self-lensing dips in each q -bin. We leave a full exploration of this dependence to future work.  \nTable V. Number of detectable dips varying with model parameters and LSST limits: f bin is the assumed fraction of quasars that are binaries, f Edd is the Eddington ratio, T max is the maximum allowed orbital period, and τ Q is the typical bright quasar lifetime. The first row corresponds to our fiducial model, while each subsequent row varies one of the parameters, shown in bold. Various maximum allowed lensing durations τ f, min were additionally considered and are shown separately in Table VI.', 'C. Dependence on model parameters': 'We next examine how the number of detectable dips depends on the orbital period and other model parameters. The dependence on T max comes from the probability of an observable self-lensing dip being proportional to ∆ ϕ ∝ T -2 / 3 orb , and the residence time as a function of period τ res ∝ T 8 / 3 orb at shorter orbital periods and τ res ∝ T orb at longer orbital periods. Doubling the fiducial maximum period to 10 years would yield an increase in the number of self-lensing dips by a factor between 4 and 1.26. This is consistent with the increase by a factor of ∼ 1 . 5 as shown in Table V when T max = 10 years.  \nOur results also depend on the assumed Eddington ratio f Edd , average quasar lifetime τ Q and the fraction of quasars that are binaries, f bin , as presented in Table V.  \nAccording to Eq. 1, changing f Edd from 0 . 3 to 0 . 1 would increase the total mass of our BHs by a factor of 3. The total mass dependence in our model also comes from the self-lensing dip probability ∝ M 4 / 3 , and the residence time, which is either ∝ M -5 / 3 for GW-driven inspiral or independent of total mass ( ∝ M 0 ) for gas-driven inspiral. Therefore, we expect decreasing f Edd from 0 . 3 to 0 . 1 to increase the number of detectable dips by a factor between 0.7 and 4.3, which is consistent with our results, where we find a factor of ∼ 2 increase.  \nFinally, our results are linearly proportional to τ -1 Q and f bin , which causes the number of dips to vary within a few orders of magnitude. As mentioned, initially we assume that all quasars are binaries ( f bin = 1). In reality, the fraction of quasars that are binaries is a function of  \nFigure 5. Distribution of self-lensing flares in the pointsource approximation, without requiring dips. Blue curves show LSST sensitivity limits, as usual, and the black contour lines show cutoffs resulting from requiring that lensing durations exceed 30, 20, and 10 days for our typical binaries with q = 1 , T orb = 5 years. The gray contour lines correspond to the same minimum lensing durations but for binaries with q = 0 . 1, which are the typical binary parameters when requiring a minimum lensing duration of 30 days. Changing the mass ratio from q = 1 to 0 . 1 increases M lens by a factor of ∼ 2 and lower total masses are needed to produce the same lensing duration contour. We further demonstrate the change in the mass-ratio distribution before and after the lensing requirement in Figure 6.  \n<!-- image -->  \nbinary parameters M and z . In Table V, we list the numbers for varying f bin . The numbers are also sensitive to the model for the residence time τ res . We find that if we conservatively set a maximum residence time of 1 . 7 × 10 7 years (orange curve in Figure 2), then there are only 5 , 6 , 7 self-lensing dips for the limits m i = 24 , 25 , 26. We also note that the expected phase spacing of the dip, Eq. 7, is valid only for circular orbits and would not hold for non-circular orbits.', 'D. Self-lensing flares and comparison to K21': "Although our main focus in this paper is on the number of detectable dips, we next discuss the number of detectable lensing flares, without considering whether or not a measurable dip may be present. As mentioned in the Introduction, the number of detectable flares has been previously estimated by [34]. We therefore make some adjustments to our fiducial model, to mimic the assumptions in K21 as closely as we can, and then compare our results to theirs in detail. The results, analogous to the number of dips, are shown in Figure 5.  \nTo make a direct comparison with the results of [34], we need to apply similar constraints on our SMBH population from the QLF, in terms of lensing amplifica-  \nTable VI. The number of detectable self-lensing dips and flares as a function of the required minimum lensing duration τ f , min and LSST magnitude limit. The point-source (PS) approximation is assumed for the flares. We show results for various LSST sensitivity limits. Likewise, the first row shows the results in our fiducial model, with the parameters that are varied shown in bold in subsequent rows. The number of detectable flares is a steep function of the minimum required flare duration. Figures 5 and 8 (right panel) visualize the effects of requiring these lensing durations.  \ntion, sensitivity limits, lensing durations and orbital periods. We therefore adopt K21's sensitivity limit of 3 × 10 -30 erg s -1 cm -2 Hz -1 , which corresponds to m i = 25 . 2. We also require that self-lensing durations exceed 30 days (10 intra-flare data points assuming three-day LSST cadence) using the Einstein radius at alignment  \nr E = √ 2 ar s cos I orb , (21)  \nwhere r s = 2 GM lens /c 2 is the Schwarzschild radius of the lens, M lens = M/ (1 + q ). The approximation for the lensing duration is given by  \nτ f = T orb π sin -1 ( r E a ) , (22)  \nwhere T orb is the observed orbital period. Finally, we consider only binaries with an observed orbital period of T max ≤ 5 years. Evaluating the number of detectable self-lensing flares in our BH population gives 737 detectable self-lensing flares, compared to their 450 +29 -40 flares, where the superscript and subscript denote the interquartile ranges for their expected number. See Table VII for a detailed summary of the two models, where we compare the parameter ranges of the BH populations, the SLF requirements, and the medians of various binary parameters of the SLFs.  \nIn both models, it is notable that as a result of requiring a minimum lensing duration of 30 days, the median mass of the binaries with SLFs are above ∼ 10 8 M ⊙ , significantly higher than the typical mass of both parent SMBH populations, 10 6 -7 M ⊙ . This occurs because the lensing duration increases with total mass. Furthermore, the lensing duration is longer when, for a given total mass M , the mass of the primary M lens = M/ (1 + q ) is larger. As a result, both models select relatively asymmetric mass ratios ( q = 10 -2 ∼ 10 -1 ) from the parent SMBH population, which is predominantly symmetric ( q = 1).  \nFigure 7. Peak magnification for a binary with q = 1 , M = 10 8 M ⊙ , z = 2 , T orb = 5 years. In the finite-source approximation, peak magnifications in the r-band (622 nm) and y-band (1015 nm) are shown as a function of orbital inclination in degrees. In this sample binary, the largest orbital inclination which admits a 10% flare is ∆ i FS1 / ∆ i PS = 1 . 24 times larger for the finite-source flare in the r-band compared to pointsource flare. This occurs for most of the binaries in our BH population. Between the 622 nm flare and the 1015 nm flare, there is only a ∆ i FS2 / ∆ i FS1 = 1 . 05 factor difference, which accounts for the mild increase in the number of detectable flares for longer wavelengths. For reference, for this binary ∆ i dip = 0 . 1 · , which is 50 -70 times smaller than ∆ i PS or ∆ i FS .  \n<!-- image -->  \nFigure 6. The mass-ratio distribution of dips and flares in the point-source approximation, where a minimum lensing duration of 0 days or 30 days is required. The two blue N dips distributions correspond to the mass ratio distributions of the 51 and 9 dips given in the m i = 25 column of Table VI. The two red N flares distributions correspond to the 97118 and 737 flares of the same m i = 25 column in Table VI.  \n<!-- image -->  \nIncreasing the orbital period increases the lensing duration and so binaries with orbital periods of several years tend to exhibit detectable self-lensing flares. Finally, the redshift dependence of the various lensing requirements is weak and the median redshift is mostly unchanged relative to the parent SMBH binary population.  \nWhile K21 requires a minimum flare duration of 30 days, equivalent to 10 intra-flare data points, it could be advantageous to include the much more numerous shorter flares (e.g. 10 or 20 days for 3 or 6 intra-flare data points) and search for additional data in higher-cadence surveys or perform targeted high-cadence follow-up observations. We find that there is a steep dependence on the lensing duration, where ∼ 10 4 -5 flares are reduced to just 734737 after requiring the lensing duration of 30 days. We show the dependence of the number of flares on their duration in Table VI. We fix all other parameters to their fiducial values (first row of Table I). In Figure 5, we further show the cut on the SLF binaries due to minimum duration requirements for two different mass ratios, q = 1 (black curves) and q = 0 . 1 (gray curves).  \nFinally, we demonstrate the mass-ratio distribution of the self-lensing flares and dips in our model in Figure 6. The mass-ratio distributions of flares and dips without requiring a minimum lensing duration are shown in red and blue. We assume m i = 25 and the other fiducial parameters given in Table I. In this case, most flares and dips have q = 1, following the mass-ratio distribution of Illustris binaries, previously shown in Figure 3. The mass-ratio distributions of flares and dips after requiring a minimum 30-day lensing duration are shown in dark red and dark blue. As a result of the minimum lensing-  \nduration requirements, the peak of the distribution shifts to q = 0 . 1, while the qualitative trend of N dips is relatively unchanged because there is only 1 dip for q < 0 . 1.", 'E. Self-lensing flares - finite source sizes': "To predict the number of detectable self-lensing flares in the previous section, we used the analytic methods in § II E, which assume that the background (lensed) SMBH is a point source. However, in some cases, the size of the source, assuming the emission arises from a minidisk modeled as a standard α -disk, becomes comparable to the Einstein radius, and must be taken into account [51].  \nTo see how a finite source size impacts our results, we re-compute the number of detectable lensing flares using the equations from § II F. Overall, we find a ∼ 25% increase in the number of detectable self-lensing flares relative to the point-source approximation, across the various wavelengths and LSST sensitivity limits, which we summarize in Table VIII.  \nTable VII. Comparisons between our model vs. Kelley et al. [34]. The top part of the table compares our assumed parent SMBH populations, where we rely on the QLFs and use Illustris only for mass ratios. In contrast, [34] self-consistently evolve ∼ 10 4 SMBH binaries extracted from Illustris using a semi-analytical model and a re-sampling scheme. Our QLF-based population extends to 10 times lower masses and has roughly 30 times more BHs. We also fix the Eddington ratio to 0.3, whereas [34]'s binaries have varying Eddington ratios typically at ∼ 0 . 01. For a binary with given total mass M and redshift z , our binaries are ∼ 30 times brighter than K21's binaries, which could account for why we predict 2 -10 times more detectable flares for minimum lensing durations of 10 -30 days. In the middle part , we compare the requirements imposed on BH binaries for detectable self-lensing flares. The only significant difference is in the minimum magnification, where we require 10%, whereas [34] require a minimum 5% magnification in addition to the inverse of signal-to-noise ratio SNR , which is dependent on their damped random walk model for intrinsic AGN variability. Finally, in the bottom part , we compare the median parameters of binaries that exhibit detectable self-lensing flares and also the final number of self-lensing flares. Despite our different approaches, the predicted number of flares are comparable (737 vs. 450).  \nTable VIII. The number of detectable self-lensing flares in the finite-source approximation for different required minimum lensing durations τ f , min in different LSST filters/wavelengths. The three right-most columns show results for various assumed detection limits.  \nThis result is somewhat counter-intuitive, since naively, one would expect point sources to be more highly magnified. To understand why there are more detectable flares from finite-sized sources, we consider ∆ i PS and ∆ i FS , which are defined by the largest orbital inclination that produces a self-lensing flare of 10%. These are proportional to the probabilities P PS , P FS that a binary exhibits a detectable self-lensing flare. We visualize the differences between ∆ i PS and ∆ i FS in Figure 7. This figure shows that for orbital inclinations close to zero, the point-source magnification at alignment is much greater than the finite-source magnification, but the opposite happens for larger orbital inclinations. For the typical binaries that constitute most of our BH population with q ∼ 1 , M/M ⊙ ∈ [10 6 , 10 8 ] , z ∈ [1 , 3] , T orb ∼ 5 years, we find that ∆ i FS / ∆ i PS = 1 . 1 -1 . 3, which accounts for the difference between the number of point-source flares and the number of finite-source flares. Notably, for less-common binaries with q < 0 . 1 or M > 10 9 M ⊙ , we find that the point-source and finite-source magnifications converge, i.e. ∆ i PS ∼ ∆ i FS , and these sources do not affect the overall difference in the number of flares.  \nIn Figure 7, we also compare the peak magnifications between different wavelengths as a function of orbital inclination. As expected, shorter wavelengths dominate at orbital inclinations close to zero but fall off at longer wavelengths. Thus, longer wavelengths lead to a wider  \n<!-- image -->  \n<!-- image -->  \nFigure 8. Illustration of going from the parent SMBH population to the detectable lensing dips. Left panel : Mass- and redshiftdistribution of all quasars (same as Fig. 1). For the LSST sensitivity limits of m i = 24 , 25 , 26, we find 20 million, 44 million, and 100 million quasars. Middle panel : Number of binaries with an orbital period of 5 years or less. For the LSST sensitivity limits of m i = 24 , 25 , 26, we find 7 million, 17 million, and 45 million binaries. Right panel : Number of detectable self-lensing dips in binaries (same as Fig. 4). For the LSST sensitivity limits of m i = 24 , 25 , 26, we find 41, 51, and 60 binaries in our fiducial model. The visualizations of the lensing duration cutoffs are shown for q = 1 (black) and q = 0 . 1 (gray).  \n<!-- image -->  \nrange of orbital inclinations that permit a detectable flare. Between the chosen LSST filter wavelengths, we typically find that discrepancies between ∆ i FS are a few degrees or less, which amount to few % or less differences in the self-lensing probabilities P PS , P FS and the number of detectable self-lensing flares.  \nThe steep dependence on the minimum lensing duration which occurred for the point-source flares also occurs in the finite-source case. Furthermore, for both pointlike and finite sources, the number of detectable flares depends sensitively on the LSST sensitivity limits when no lensing duration was required, but when a minimum lensing duration of 30 days is required, the dependence on the LSST sensitivity limits diminishes.  \nFinally, we note that in general, in addition to the finite radial extent of the 'minidisks' around the individual SMBHs, one must take into account the finite thickness of these minidisks, as well as the circumbinary gas. This raises the concern that a thick circumbinary disk for a nearly edge-on binary may obscure the lensing phenomena discussed in our paper. However, ray-tracing the emission through a hydrodynamical simulation [30] has shown that the circumbinary disk in the foreground typically does not obscure the self-lensing flare and dip. This is because lensing itself allows us to see into the cavity along bent photon paths.", 'IV. SUMMARY AND CONCLUSIONS': "In this paper, we estimated the number of binary quasars that are sufficiently edge-on to display detectable self-lensing flares and dips. Our approach relies on the assumption that the same galaxy mergers produce SMBH  \nbinaries and activate bright quasars. Combined with further assumptions about quasar lifetimes, binary fractions, and Eddington ratios, as well as the mass-ratio distributions of binaries extracted from the Illustris simulations, we computed the number of detectable flares, as well as the number of flares with detectable 'dips' due to the background SMBH's shadow, in LSST's expected catalog of tens of millions of AGN light-curves.  \nWe recover earlier results by K21 and find that several hundred lensing flares may be detected by LSST. However, this is based a requiring a minimum flare duration of 30 days. We demonstrated a steep dependence of this number on the minimum required flare duration and found that if much shorter flares (say 10 days) were recoverable, then there would be many more (tens of thousands) of these.  \nOur main novel result is that we estimate the number of self-lensing dips from SMBHBs. Under plausible assumptions, we find that several dozen of these should be present and detectable by LSST. Finding these dips, a feature of the BHs shadow, would indisputably prove the existence of SMBHBs in distant galaxies, that are currently out of reach for high-resolution VLBI experiments.", 'ACKNOWLEDGMENTS': "We thank Luke Krauth for useful discussions. We acknowledge support by NSF grant AST2006176 and NASA grants 80NSSC22K0822 and 80NSSC24K0440 (ZH). JD acknowledges support by a joint Columbia/Flatiron Postdoctoral Fellowship. Research at the Flatiron Institute is supported by the  \nSimons Foundation. This research was supported in part by the National Science Foundation under Grant No. NSF PHY-1748958 to the Kavli Institute for Theoretical Physics (KITP). This research has made use of NASA's Astrophysics Data System.  \n- [1] J. Kormendy and L. C. Ho, ARA&A 51 , 511 (2013), arXiv:1304.7762 [astro-ph.CO].\n- [2] V. Springel, T. Di Matteo, and L. Hernquist, ApJL 620 , L79 (2005), arXiv:astro-ph/0409436.\n- [3] B. Robertson, T. J. Cox, L. Hernquist, M. Franx, P. F. Hopkins, P. Martini, and V. Springel, ApJ 641 , 21 (2006), arXiv:astro-ph/0511053.\n- [4] J. E. Barnes and L. Hernquist, ARA&A 30 , 705 (1992).\n- [5] M. C. Begelman, R. D. Blandford, and M. J. Rees, Nature 287 , 307 (1980).\n- [6] B. D. Farris, P. Duffell, A. I. MacFadyen, and Z. Haiman, MNRASL 447 , L80 (2015), arxiv:1409.5124 [astro-ph].\n- [7] Y. Tang, Z. Haiman, and A. MacFadyen, MNRAS 476 , 2249 (2018).\n- [8] L. Major Krauth, J. Davelaar, Z. Haiman, J. R. Westernacher-Schneider, J. Zrake, and A. MacFadyen, arXiv e-prints , arXiv:2304.02575 (2023), arXiv:2304.02575 [astro-ph.HE].\n- [9] A. J. Dittmann, G. Ryan, and M. C. Miller, ApJL 949 , L30 (2023), arXiv:2303.16204 [astro-ph.HE].\n- [10] P. Amaro-Seoane, J. Andrews, M. Arca Sedda, et al. , Living Reviews in Relativity 26 , 2 (2023), arXiv:2203.06016 [gr-qc].\n- [11] M. Colpi, K. Danzmann, M. Hewitson, et al. , e-print arXiv:2402.07571 (2024).\n- [12] G. Agazie, A. Anumarlapudi, and A. M. Archibald et al. (The Nanograv Collaboration), ApJL 951 , L8 (2023), arXiv:2306.16213 [astro-ph.HE].\n- [13] J. Antoniadis, P. Arumugam, and S. Arumugam et al. (The EPTA and InPTA Collaborations), A&A 678 , A50 (2023), arXiv:2306.16214 [astro-ph.HE].\n- [14] D. J. Reardon, A. Zic, R. M. Shannon, G. B. Hobbs, M. Bailes, V. Di Marco, A. Kapur, A. F. Rogers, and E. Thrane et al., ApJL 951 , L6 (2023), arXiv:2306.16215 [astro-ph.HE].\n- [15] H. Xu, S. Chen, Y. Guo, J. Jiang, B. Wang, J. Xu, Z. Xue, R. Nicolas Caballero, J. Yuan, Y. Xu, J. Wang, L. Hao, J. Luo, K. Lee, J. Han, P. Jiang, Z. Shen, M. Wang, N. Wang, R. Xu, X. Wu, R. Manchester, L. Qian, X. Guan, M. Huang, C. Sun, and Y. Zhu, Research in Astronomy and Astrophysics 23 , 075024 (2023), arXiv:2306.16216 [astro-ph.HE].\n- [16] J. Baker, Z. Haiman, E. M. Rossi, E. Berger, N. Brandt, E. Breedt, K. Breivik, M. Charisi, A. Derdzinski, D. J. D'Orazio, S. Ford, J. E. Greene, J. C. Hill, K. Holley-Bockelmann, J. S. Key, B. Kocsis, T. Kupfer, S. Larson, P. Madau, T. Marsh, B. McKernan, S. T. McWilliams, P. Natarajan, S. Nissanke, S. Noble, E. S. Phinney, G. Ramsay, J. Schnittman, A. Sesana, D. Shoemaker, N. Stone, S. Toonen, B. Trakhtenbrot, A. Vikhlinin, and M. Volonteri, 'Multimessenger science opportunities with mHz gravitational waves,' (2019),  \nSoftware: python [55, 56], numpy [57], matplotlib [58]", 'DATA AVAILABILITY': 'The data underlying this article will be shared on reasonable request to the corresponding author.', 'arxiv:arXiv:1903.04417.': "- [17] M. J. Graham, S. G. Djorgovski, D. Stern, A. J. Drake, A. A. Mahabal, C. Donalek, E. Glikman, S. Larson, and E. Christensen, MNRAS 453 , 1562 (2015), arXiv:1507.07603.\n- [18] M. Charisi, I. Bartos, Z. Haiman, A. M. Price-Whelan, M. J. Graham, E. C. Bellm, R. R. Laher, and S. M'arka, MNRAS 463 , 2145 (2016), arXiv:1604.01020.\n- [19] Y.-J. Chen, S. Zhai, J.-R. Liu, W.-J. Guo, Y.-C. Peng, Y.-R. Li, Y.-Y. Songsheng, P. Du, C. Hu, and J.M. Wang, MNRAS 527 , 12154 (2024), arXiv:2206.11497 [astro-ph.GA].\n- [20] A. De Rosa, C. Vignali, T. Bogdanovi'c, P. R. Capelo, M. Charisi, M. Dotti, B. Husemann, E. Lusso, L. Mayer, Z. Paragi, J. Runnoe, A. Sesana, L. Steinborn, S. Bianchi, M. Colpi, L. del Valle, S. Frey, K. ' E. Gab'anyi, M. Giustini, M. Guainazzi, Z. Haiman, N. Herrera Ruiz, R. Herrero-Illana, K. Iwasawa, S. Komossa, D. Lena, N. Loiseau, M. Perez-Torres, E. Piconcelli, and M. Volonteri, New. Ast. Rev. 86 , 101525 (2019).\n- [21] T. Bogdanovi'c, M. C. Miller, and L. Blecha, Living Reviews in Relativity 25 , 3 (2022), arXiv:2109.03262 [astroph.HE].\n- [22] S. Vaughan, P. Uttley, A. G. Markowitz, D. Huppenkothen, M. J. Middleton, W. N. Alston, J. D. Scargle, and W. M. Farr, MNRAS 461 , 3145 (2016), arXiv:1606.02620 [astro-ph.IM].\n- [23] Z. Haiman, Phys. Rev. D 96 , 023004 (2017), arXiv:1705.06765 [astro-ph.HE].\n- [24] D. J. D'Orazio and R. Di Stefano, MNRAS 474 , 2975 (2018).\n- [25] P. Pihajoki, M. Mannerkoski, J. Nattila, and P. H. Johansson, The Astrophysical Journal 863 , 8 (2018).\n- [26] A. Ingram, S. E. Motta, S. Aigrain, and A. Karastergiou, MNRAS 503 , 1703 (2021), arXiv:2103.00017 [astroph.HE].\n- [27] B. X. Hu, D. J. D'Orazio, Z. Haiman, K. L. Smith, B. Snios, M. Charisi, and R. Di Stefano, MNRAS 495 , 4061 (2020), arXiv:1910.05348 [astro-ph.HE].\n- [28] J. Davelaar and Z. Haiman, Phys. Rev. Lett. 128 , 191101 (2022), arXiv:2112.05829 [astro-ph.HE].\n- [29] J. Davelaar and Z. Haiman, PhRD 105 , 103010 (2022).\n- [30] L. M. Krauth, J. Davelaar, Z. Haiman, J. R. Westernacher-Schneider, J. Zrake, and A. MacFadyen, Physical Review D 109 (2024), 10.1103/physrevd.109.103014.\n- [31] K. Porter, S. C. Noble, E. M. Gutierrez, J. Pelle, M. Campanelli, J. Schnittman, and B. J. Kelly, (2024), arXiv:2407.04089 [astro-ph.HE].\n- [32] Vera C. Rubin Observatory LSST Solar System Science Collaboration, R. L. Jones, M. T. Bannister, B. T. Bolin, C. O. Chandler, S. R. Chesley, S. Eggl, S. Greenstreet,\n- T. R. Holt, H. H. Hsieh, Z. Ivezic, M. Juric, M. S. P. Kelley, M. M. Knight, R. Malhotra, W. J. Oldroyd, G. Sarid, M. E. Schwamb, C. Snodgrass, M. Solontoi, and D. E. Trilling, in Bulletin of the American Astronomical Society , Vol. 53 (2021) p. 236.\n- [33] C. Xin and Z. Haiman, Monthly Notices of the Royal Astronomical Society 506 , 2408-2417 (2021).\n- [34] L. Z. Kelley, D. J. D'Orazio, and R. Di Stefano, Monthly Notices of the Royal Astronomical Society 508 , 2524-2536 (2021).\n- [35] D. Sijacki, M. Vogelsberger, S. Genel, V. Springel, P. Torrey, G. F. Snyder, D. Nelson, and L. Hernquist, Monthly Notices of the Royal Astronomical Society 452 , 575-596 (2015).\n- [36] G. Kauffmann and M. Haehnelt, Monthly Notices of the Royal Astronomical Society 311 , 576-588 (2000).\n- [37] Z. Haiman, B. Kocsis, and K. Menou, The Astrophysical Journal 700 , 1952-1969 (2009).\n- [38] G. Kulkarni, G. Worseck, and J. F. Hennawi, Monthly Notices of the Royal Astronomical Society 488 , 1035-1065 (2019).\n- [39] LSST Science Collaboration, P. A. Abell, J. Allison, S. F. Anderson, J. R. Andrew, et al. , arXiv e-prints , arXiv:0912.0201 (2009), arXiv:0912.0201 [astro-ph.IM].\n- [40] D. B. Sanders, B. T. Soifer, and G. Neugebauer, Dynamics and Interactions of Galaxies , 459-464 (1990).\n- [41] P. F. Hopkins, L. Hernquist, T. J. Cox, and D. Kereˇs, The Astrophysical Journal Supplement Series 175 , 356-389 (2008).\n- [42] P. C. Peters, PhR 136 , B1224 (1964).\n- [43] E. Bortolas, A. Franchini, M. Bonetti, and A. Sesana, The Astrophysical Journal Letters 918 (2021), 10.3847/2041-8213/ac1c0c.\n- [44] Y. Tang, A. MacFadyen, and Z. Haiman, MNRAS 469 , 4258 (2017), arXiv:1703.03913 [astro-ph.HE].  \n- [45] D. J. Mu˜noz, R. Miranda, and D. Lai, The Astrophysical Journal 871 , 84 (2019).\n- [46] M. S. L. Moody, J.-M. Shi, and J. M. Stone, ApJ 875 , 66 (2019), arXiv:1903.00008 [astro-ph.HE].\n- [47] D. J. Mu˜noz, D. Lai, K. Kratter, and R. Miranda, The Astrophysical Journal 889 , 114 (2020).\n- [48] P. C. Duffell, D. D'Orazio, A. Derdzinski, Z. Haiman, A. MacFadyen, A. L. Rosen, and J. Zrake, ApJ 901 , 25 (2020).\n- [49] C. Tiede, J. Zrake, A. MacFadyen, and Z. Haiman, ApJ 900 , 43 (2020), arXiv:2005.09555 [astro-ph.GA].\n- [50] B. Paczynski, The Astrophysical Journal 304 , 1 (1986).\n- [51] D. J. D'Orazio and R. Di Stefano, Monthly Notices of the Royal Astronomical Society 474 , 2975-2986 (2017).\n- [52] C. Roedig, J. H. Krolik, and M. C. Miller, The Astrophysical Journal 785 , 115 (2014).\n- [53] P. Madau, The Astrophysical Journal 441 , 18 (1995).\n- [54] M. Sayeb, L. Blecha, and L. Z. Kelley, Monthly Notices of the Royal Astronomical Society 527 , 7424-7437 (2023).\n- [55] T. E. Oliphant, Computing in Science & Engineering 9 , 10 (2007), https://aip.scitation.org/doi/pdf/10.1109/MCSE.2007.58.\n- [56] K. J. Millman and M. Aivazis, Computing in Science & Engineering 13 , 9 (2011), https://aip.scitation.org/doi/pdf/10.1109/MCSE.2011.36.\n- [57] S. van der Walt, S. C. Colbert, and G. Varoquaux, Computing in Science and Engineering 13 , 22 (2011), arXiv:1102.1523 ['cs.MS'].\n- [58] J. D. Hunter, Computing in Science and Engineering 9 , 90 (2007)."}