Tidal Disruptions of Stars by Binary Black Holes: Modifying the Spin Magnitudes and Directions of LIGO Sources in Dense Stellar Environments
Martin Lopez Jr., Aldo Batta, Enrico Ramirez-Ruiz, Irvin Martinez, Johan Samsing
II N PREPARATION FOR
ApJ . DRAFT OF D ECEMBER
5, 2018.
Preprint typeset using L A TEX style AASTeX6 v. 1.0
TIDAL DISRUPTIONS OF STARS BY BINARY BLACK HOLES: MODIFYING THE SPIN MAGNITUDES ANDDIRECTIONS OF LIGO SOURCES IN DENSE STELLAR ENVIRONMENTS M ARTIN L OPEZ J R . , A LDO B ATTA , E
NRICO R AMIREZ -R UIZ , I
RVIN M ARTINEZ , J OHAN S AMSING Department of Astronomy and Astrophysics, University of California, Santa Cruz, CA 95064, USA Niels Bohr Institute, University of Copenhagen, Blegdamsvej 17, 2100 Copenhagen, Denmark Department of Physics & Astronomy, Princeton University, Princeton, NJ, 08544
ABSTRACTBinary black holes (BBHs) appear to be widespread and are able to merge through the emission of gravitationalwaves, as recently illustrated by LIGO. The spin of the BBHs is one of the parameters that LIGO can inferfrom the gravitational wave signal and can be used to constrain their production site. If BBHs are assembledin stellar clusters they are likely to interact with stars, which could occasionally lead to a tidal disruptionevent (TDE). When a BBH tidally disrupts a star it can accrete a significant fraction of the debris, effectivelyaltering the spins of the BHs. Therefore, although dynamically formed BBHs are expected to have random spinorientations, tidal stellar interactions can significantly alter their birth spins both in direction and magnitude.Here we investigate how TDEs by BBHs can affect the properties of the BH members as well as exploringthe characteristics of the resulting electromagnetic signatures. We conduct hydrodynamic simulations with aLagrangian Smoothed Particle Hydrodynamics code of a wide range of representative tidal interactions. Wefind that both spin magnitude and orientation can be altered and temporarily aligned or anti-aligned throughaccretion of stellar debris, with a significant dependence on the mass ratio of the disrupted star and the BBHmembers. These tidal interactions feed material to the BBH at very high accretion rates, with the potential tolaunch a relativistic jet. The corresponding beamed emission is a beacon to an otherwise quiescent BBH.
Keywords: black holes, tidal disruptions, close binaries, dense stellar systems, LIGO INTRODUCTIONA watershed event occurred on September 14 2015, whenthe Laser Interferometer Gravitational-Wave Observatory(LIGO) succeeded in detecting the first gravitational wave(GW) signal (Abbott et al. 2016), GW150914, of a binaryblack hole (BBH) merger. This detection, followed by fiveothers, has unveiled a population of stellar mass BHs thatis significantly heavier than those inhabiting X-ray binaries(Farr et al. 2011).A large number of progenitor systems have been sug-gested, all designed to manufacture BHs in the observedmass range. The two most widely discussed scenarios en-compass dynamical assembly in dense star clusters (Sigurds-son & Hernquist 1993; Portegies Zwart & McMillan 2000;Downing et al. 2010, 2011; Ziosi et al. 2014; Samsing et al.2014; Rodriguez et al. 2015, 2016c,b; Samsing & Ramirez-Ruiz 2017; Samsing et al. 2018a) and isolated massive stellarfield binaries (Paczynski 1976; Iben & Livio 1993; Podsi-adlowski 2001; Voss & Tauris 2003; Kalogera et al. 2007;Taam & Ricker 2010; Dominik et al. 2012, 2013; Ivanovaet al. 2013; Postnov & Yungelson 2014; Belczynski et al.2016; Schrøder et al. 2018), including chemically homoge-neous stars (de Mink et al. 2009; Marchant et al. 2016; Man-del & de Mink 2016; de Mink & Mandel 2016).Other proposed scenarios include active galactic nuclei(AGN) discs (Bartos et al. 2017; Stone et al. 2017; McKernan et al. 2017), galactic nuclei (O’Leary et al. 2009; Hong & Lee2015; VanLandingham et al. 2016; Antonini & Rasio 2016;Stephan et al. 2016; Hoang et al. 2017), single-single GWcaptures of primordial BHs (Bird et al. 2016; Cholis et al.2016; Sasaki et al. 2016; Carr et al. 2016), and very mas-sive stellar mergers (Loeb 2016; Woosley 2016; Janiuk et al.2017; D’Orazio & Loeb 2017). Generally, these theoreti-cal predicted channels can be broadly tuned to be consistentwith the properties and rates of the BBH sources observedby LIGO so far, and the challenge remains to find reliableobservational tests.Recent work suggests that the key parameters that mighthelp discriminating between formation channels include theBH mass (e.g. Zevin et al. 2017), orbital eccentricity in LIGO(O’Leary et al. 2009; Kocsis & Levin 2012; Samsing et al.2014; O’Leary et al. 2016; Samsing & Ramirez-Ruiz 2017;Samsing & Ilan 2018; Samsing et al. 2018b; Samsing & Ilan2019; Samsing 2018; Samsing et al. 2018a; Zevin et al. 2018;Rodriguez et al. 2018a; Gondán et al. 2018) and LISA ( e.g. ,Samsing & D’Orazio 2018), and especially the dimension-less spin parameter χ eff (Farr et al. 2018, 2017; Rodriguezet al. 2016c; Schrøder et al. 2018). χ eff is the total massweighted BH spin components in the direction of the orbitalangular momentum, χ eff = M bh a bh + M bh a bh M bh + M bh · ˆL . (1) a r X i v : . [ a s t r o - ph . H E ] D ec L OPEZ J R . ET AL .Here a bh and a bh are the dimensionless spins of the BHsand ˆL is the direction of the orbital angular momentum.The spin measurements of BBHs arising from the isolatedmassive stellar field binary scenario roughly predicts align-ment of the BH spins and the orbital angular momentum(Kalogera 2000), while dynamically assembled BHs are ex-pected to have uncorrelated spins as they are formed andharden through a series of chaotic exchange interactions (Ro-driguez et al. 2016c).Here we will analyze the dynamical scenario and investi-gate whether the determination of χ eff allows for constraintsto be placed on the spin history of the BBH system betweenassembly and merger. Such a BBH becomes detectable onlythrough interactions with its gaseous environment. Gas thatis lost from nearby stars, or even stars plunging into such bi-naries, can produce detectable signatures. Through the useof Smoothed Particle Hydrodynamic (SPH) simulations, weshow how stellar material which is accreted following a tidaldisruption event (TDE) can alter the birth spin magnitudesand orientation of the individual BHs, possibly aligning ormisaligning them temporarily. Furthermore, the supply ofmaterial to the BBH is above the Eddington limit and couldlaunch a relativistically-beamed jet. The emerging class ofhigh energy transient bursts all have peak luminosities anddurations reminiscent of ultra-long γ -ray bursts. Tidal dis-ruptions of stars by BBHs thus uniquely probe the currently-debated existence of LIGO signals emanating from dense starclusters.The structure of the paper is as follows. Section 2 discussesthe dynamics of LIGO BBH (LBBH) TDEs in dense starclusters. Section 3 overviews the hydrodynamic formalismand presents the results as well as their significance for thespin magnitude and alignment of the individual BHs. Sec-tion 4 explores the implications of our results and possiblesources for upcoming high energy transient surveys. TIDAL DISRUPTION EVENTS BY LIGO BBHS2.1.
Single BH Dynamics
Canonical TDEs occur when a star with mass M ∗ andradius R ∗ gets disrupted when approaching a supermassiveblack hole (SMBH) with mass M bh at a pericenter distance R p = R τ = q − / R ∗ , where q = M ∗ /M bh (Rees 1988;Phinney 1989; Evans & Kochanek 1989). After the disrup-tion, about half of the star becomes unbound and ejected,while the other half becomes bound to the SMBH on ellip-tical orbits. 3D hydrodynamical simulations have quantifiedthe rate at which material falls back onto the SMBH (Guil-lochon & Ramirez-Ruiz 2013). A good fit to observed lightcurves of TDEs is obtained if one assumes that the accretionluminosity directly follows the fallback rate in the simula-tion (Mockler et al. 2018). However, it is not clear why thisshould be the case. Bound debris returns to the SMBH witha large range of eccentricities and orbital periods (Ramirez-Ruiz & Rosswog 2009) and it may take many Keplerian or-bits for fallback material to circularize and accrete (Guillo-chon & Ramirez-Ruiz 2015). Some mechanism is thereforerequired to quickly dissipate the kinetic energy of the fall-back material and circularize it into an accretion disk. In standard TDE discourse (Rees 1988), the disruptingSMBHs have masses M bh (cid:38) M (cid:12) yielding q (cid:28) , whichallows the semi-major axis of the most bound material to beapproximated as: a mb = (cid:18) M bh M (cid:63) (cid:19) / R τ = q − / R τ . (2)However for disrupting BHs within a LBBH, the mass ratiois near unity, making the extent of the star comparable to thetidal radius. In this case, the specific orbital energy of stellarmaterial varies significantly across the star: E ( r ) = − G M bh (cid:34) ∞ (cid:88) n =1 (cid:18) q / R (cid:63) (cid:19) n +1 r n (cid:35) , (3)where r is the distance from the star’s center of mass (CM).For material that is bound to the BH, this expression trans-lates into a range of semi-major axes given by: a ( r ) = − G M bh E ( r )= (cid:34) ∞ (cid:88) n =1 (cid:18) q / R (cid:63) (cid:19) n +1 r n (cid:35) − , (4)which for canonical TDEs ( q (cid:28) can be safely approx-imated to first order. As q approaches unity this approx-imation is no longer valid and the semi-major axis of themost bound material approaches the tidal radius and becomesequal to it at a critical mass ratio q crit = 0 . . The assump-tion that the circularization radius of the most bound materialis about twice the tidal radius (Cannizzo et al. 1990; Ulmer1998; Gezari et al. 2009; Lodato & Rossi 2011; Strubbe &Quataert 2011; Guillochon et al. 2014) also breaks down inthe LBBH regime. The circularization radius of the mostbound material R c , mb in this case is given by R c , mb = 2R τ (cid:104) − q / (cid:105) , (5)while the spread in circularization radii can be written as ∆R c R c , p = q / (cid:104) − q / (cid:105) , (6)where the circularization radius of the pericenter is R c , p =2R τ . In order for this material to circularize and form a disk,energy must be dissipated efficiently after disruption. Ma-terial falling to pericenter can be heated by hydrodynamicalshocks and Guillochon et al. (2014) show that the fractionalenergy dissipation per orbit, ν H , can be written as ν H = βq / (7)where β = R p / R τ . For disruptions in the LBBH regime,the energy dissipation via shocks at pericenter can be siz-able and lead to efficient circularization. This is in contrastto the standard case with q (cid:28) , for which hydrodynamicalshocks at pericenter are likely to be insufficient and rapid cir-cularization might only be achieved via general relativisticeffects (Shiokawa et al. 2015; Guillochon & Ramirez-Ruiz PIN E VOLUTION OF
LIGO S
OURCES x y z I x y II N d = 0.1 AU d = 1.0 AU d = 0.316 AU ε d = 0.316 AU d = 1.0 AU ε = -0.12 ε =0.035 Figure 1 . The CM energy distributions with respect to the disrupting BH and specific trajectories are shown for a sun-like star ( M ∗ = 1 M (cid:12) , R ∗ = R (cid:12) ) interacting with a (cid:12) equal mass BH binary with e = 0 . . The properties of the binary have been selected to reflect thosederived by Rodriguez et al. (2016a) for dynamically assembled LBBHs. Here we study the outcomes of TDE interactions and their associatedCM energy distributions by performing large set of numerical scattering experiments using the N -body code developed by Samsing et al.(2014). The panels shows the energy distribution for different binary separations ( d = 1 . . τ , d = 0 .
316 AU = 27 . τ , d = 0 . . τ ) and the trajectories of unbound (I) and bound (II) stellar orbits ( orange trajectories) . Here ε is the CM energy of thestar with respect to the disrupting BH at R τ in units of the binding energy of the star. Binary BH Dynamics
For BBH TDEs, the star does not necessarily follow aparabolic orbit and the orbital deviations before disruptiondepend strongly on the separation d and eccentricity e of thebinary. The CM energy distributions of a sun-like star withrespect to the disrupting BH part of a (cid:12) equal mass BBHwith e = 0 . are shown in Figure for three distinct binaryseparations. In this case R τ = 2 . (cid:63) = 0 .
01 AU . For R τ /d (cid:28) the CM energy is essentially parabolic while alarger fraction of unbound CM orbits are observed for tighterbinaries. This is partly due to the individual BHs evolvingfaster around their binary CM as BBHs get tighter, whichthen maps to a higher relative velocity at the time of disrup-tion and thereby a higher relative energy. Note here that stel-lar elements unbound with respect to the disrupting BH canstill be bound to the CM of the BBH, which then can lead tolater accretion.After disruption, the fate of the debris also depends sen-sitively on the ratio R τ /d . If R τ /d > , the disruption willtake place outside of the binary and the infalling material willform a circumbinary disk around the system. In what follows,we refer to this scenario as the circumbinary scenario (CS).When R τ /d (cid:46) , the star will be disrupted by one of the bi-nary members but the accretion history of the debris onto thesystem is determined by d . This is due to the debris orbitingaround the disrupting BH with a wide range of semi-majoraxes such that there is always some material that is able toreach the sphere of influence of the companion BH. In or-der to determine whether or not the non-disrupting BH canaccrete significant amounts of stellar debris we make use oftwo important characteristic scales. One is the semi-majoraxis a of the disrupted material whose orbit contains 90%of the stellar debris. In other words, a is the semi-major ofmaterial whose radius, measured from the most bound ma-terial of the star inwards, contains of the stellar mass. Therefore we classify a strong interaction as being one wherethe non-disrupting BH interacts with of stellar debris.The other scale is the Roche lobe radius R L , which deter-mines the gravitational sphere of influence of the disruptingBH. R L can be written (Eggleton 1983) as R L d min = 0 . q b2 / . q / + ln (cid:16) q / (cid:17) , (8)where q b is the mass ratio of the BBH and d min is theminimum separation of the binary. When a / R L < , asmall fraction of the debris is able to interact with the non-disrupting BH but most of the stellar debris will be accretedby the disrupting BH. In this case, the tidal interaction willresemble that caused by a single BH and we refer to this asthe single scenario (SS). On the other hand, disrupted mate-rial with a / R L (cid:38) will be influenced by the companionand a sizable fraction of debris can be accreted by the non-disrupting BH. A case we refer to as the overflow scenario (OS).In order to calculate the spin change due to accretion ofdisrupted material we use (Bardeen 1970) S (M bh , f ) = (cid:0) (cid:1) / bh M bh , f (cid:40) − (cid:20) (cid:16) M bh M bh , f (cid:17) − (cid:21) / (cid:41) , (9)which assumes an initially low or non-spinning BH. Here M bh , f = M bh + fM (cid:63) is the final mass of the BH after accret-ing a fraction f of the disrupted star. For a TDE of a star ina parabolic orbit ( f = 0 . ), the maximum mass that the BHcan accrete is . (cid:63) such that the maximum spin up, S max is given by S max ( q ) = (cid:0) (cid:1) / (cid:16) q (cid:17) (cid:40) − (cid:20) (cid:16) q (cid:17) − (cid:21) / (cid:41) , (10)The values of S max for a few characteristic q ’s are S max (cid:0) q = 1 × − (cid:1) = 1 . × − , S max ( q = 0 .
01) =0 . , and S max ( q = 0 .
5) = 0 . . This clearly illustratesthat for LBBHs, the digestion of stars during the lifetime of L OPEZ J R . ET AL .the binary could lead to noticeable spin changes. HYDRODYNAMICS3.1.
Set-Up
Our hydrodynamical simulations of LBBH TDEs usea modified version of the SPH code Stellar GADGET-3(Springel 2005; Pakmor et al. 2012). GADGET-3 allows oneto accurately follow the accretion of material into sink par-ticles and the compressibility of the gas is described with agamma-law equation of state P ∝ ρ γ . By solving the Lane-Emden equation and using the same method as in Batta et al.(2017), we created three-dimensional spherically symmetricdistributions of SPH particles by mapping polytropic stars inhydrodynamical equilibrium with a structural gamma Γ set toeither 5/3 or 4/3, representative of low and high-mass stars,respectively. During the simulation, the stars are evolved hy-drodynamically according to a γ = 5 / equation of state,with the difference between Γ and γ for higher-mass (or con-vective) stars being a consequence of radiation transfer in thestar’s interior. We ran test cases of the tidal disruption of a1 M (cid:12) star by an equal mass M bh = M bh = 15M (cid:12) LBBHwith varying resolutions between N = 10 and particles,which showed clear convergence for the accretion rates andmass bound to the system.3.2. Initial Conditions
All initial conditions (ICs) assume typical parameters forLBBHs and stars in globular clusters (GCs). We take e = 0 . for the LBBH’s eccentricity and assume that the individualspins of the BHs ( S and S ) to be initially zero, which isconsistent with the small spins observed for LIGO events sofar (The LIGO Scientific Collaboration & The Virgo Collab-oration 2018). By means of a three-body code, we obtainedthe dynamical properties of the LBBH and star prior to a tidaldisruption, tracing the trajectories for all three bodies back inthe time when the incoming star lies about six tidal radii awayfrom the disrupting BH. These dynamical properties were in-cluded in the GADGET-3 IC file.3.3. Simulation Results
In Section 2.2 we have outlined three representative sce-narios for LBBH TDEs: SS, CS and OS. In the SS case wehave R τ (cid:28) d and a < R L and the event resembles thatfrom a single BH TDE in which only one BH accretes. Inthe CS case we have R τ > d and the LBBH ends up be-ing embedded in a circumbinary disk. In the OS case wehave R τ (cid:46) d and a > R L and the accretion of the dis-rupted debris by both BHs is able to produce multiple TDEs.The simulation results for the various cases outlined here arepresented in Sections 3.3.1-3.3.4 and shown in Figure 2 andFigure 3. 3.3.1. The Single Scenario
The SS simulation is characterized here by R τ /d = 0 . and a / R L = 0 . . For these ICs, almost no significant in-teraction of the disrupted material is expected to occur withthe non-disrupting BH. The SS simulation shown here is con-sistent with the scenario shown in Figure 1 for an unbound stellar orbit. The top panels in Figure 2 show the gas columndensity in the orbital plane at three different times, which areshown in units of the dynamical timescale of the star. Thebound material is observed to circularize promptly and, as aresult, the mass accretion rate is observed to follow the stan-dard mass fallback rate. However, given that q = 0 . , theearly shape of the mass accretion rate curve differs from thatderived by Guillochon & Ramirez-Ruiz (2013), which wascalculated assuming q (cid:28) . By the end of the simulation,the disrupting BH accreted a total mass of . (cid:12) and hasan accretion disk with a leftover mass of about . (cid:12) andwhose angular momentum J disk is inclined about .
75 rad with respect to the orbital angular momentum of the binary J bin . This angle is consistent with that of the star’s angu-lar momentum at the moment of disruption. Assuming thatthe bound ≈ .
12 M (cid:12) of material is accreted by the BH, theresultant spin magnitude will be S ≈ . , resulting in ananti-aligned effective spin of χ eff ≈ − . .3.3.2. The Circumbinary Scenario
The CS simulation is parametrized by R τ /d = 2 . . Thetidal radii of each BH overlap and encompass the binary,resulting in a disruption where bound material forms a cir-cumbinary disk. At the moment of disruption the orientationof the angular momentum of the star’s CM with respect to J bin is approximately .
44 rad . As the most bound materialreturns to pericenter the binary exerts a torque on the streamand, as a result, alters the angle of J disk to ≈ ; seeSection 4.1 and Figure 5. The disk rapidly circularizes due tohydrodynamical dissipation at pericenter as well as collisionsbetween the returning stream caused by the time changingbinary potential ( middle panels in Figure 2). The material re-siding in the disk is slowly accreted onto both BHs throughviscous dissipation. We stopped the simulation at approx-imately ten percent of the time it would take to ingest theentire disk and found that each BH accreted about . (cid:12) and the accretion disk has . (cid:12) of gas leftover. If we as-sume that this material is evenly accreted by both BHs, theresultant spin magnitudes will be S = S ≈ . and,given that the spin angles of each BH are aligned with J disk , χ eff ≈ − . . 3.3.3. The Overflow Scenario
The OS simulation is characterized here by R τ /d = 0 . and a / R L = 5 . . This guarantees that after the disrup-tion, a significant amount of bound disrupted material willbe able to reach the sphere of influence of the non-disruptingBH. Within this scenario, accretion onto both BHs can oc-cur, which might result in temporary BH spin alignment oranti-alignment. The star survives after the initial disruptionleading to multiple resonant TDEs, as can be seen in the bot-tom panels of Figure 2. A total of four interactions take placewith the same BH in this scenario until the star is fully dis-rupted. The angular momentum of the star with respect to J bin changes in each disruption. By the end of the simula-tion, the mass accreted by the disrupting and non-disruptingBHs is .
19 M (cid:12) and .
02 M (cid:12) , respectively. The resultantangles are .
58 rad and .
24 rad with respect to J bin forthe disrupting and non-disrupting BHs, respectively. The PIN E VOLUTION OF
LIGO S
OURCES t * =11R τ t * =500R τ t * =1000R τ t * = 16070R τ t * = 53440R τ t * = 392500R τ t * =270R τ t * =570R τ t * =770R τ Figure 2 . Simulations of the tidal interaction of stars with a LBBH. Here t ∗ and R τ are the simulation times in dynamical time units and thecorresponding tidal radii. All panels are in the orbital plane of the LBBH. Top panels:
Simulation of the SS case at three different times, fromdisruption to the subsequent accretion onto the disrupting BH.
Middle panels:
Simulation of the CS case at three different times, from theinitial disruption occurring outside the LBBH to the assembly of the circumbinary disk.
Bottom panels:
Simulation of the OS case, from initialpartial disruption of the star, followed by a second and third disruption of the remaining stellar core. During this interaction, a total of fourdisruptions occur. The simulation parameters, listed as [SS, CS, OS], are: N = [10 , , ] , R (cid:63) = [1 , , (cid:12) , Γ = [4 / , / , / , d = [429 . , . , . (cid:12) , v ∞ = [30 , , / s . In all cases M (cid:63) = 1M (cid:12) , M bh = M bh = 15M (cid:12) , and e = 0 . . L OPEZ J R . ET AL . t * =570R Orbital Plane t * =770R Orbital Plane t * = 21.5R Orbital Plane t * = 53.7R Orbital Plane t * =570R Side View t * =770R Side View t * = 21.5R Side View t * = 53.7R Side View Figure 3 . A comparison between two OS simulations. Here t ∗ denotes the time in units of the star’s dynamical timescale while R τ shows thescale of the individual BH tidal radius. Left Panel:
Here we show the simulation snapshots for the OS shown in Figure 2. The additional sideviews plotted here clearly show how the orientation of the accretion disk changes between the multiple disruptions.
Right Panel:
Shown are thesimulation snapshots for the MOS discussed in Section 3.3.4. The orbital view shows the two interactions that take place until full disruption ofthe star. In the side view snapshots one can clearly see that the the orbital angular momentum of the binary is altered by the 3-body interaction.This change is significant in this case due to the higher mass ratio between the star and the BBH. The simulation parameters for the MOS are: N = 10 , M (cid:63) = 5 M (cid:12) , R (cid:63) = 6 R (cid:12) , Γ = 4 / , M bh = M bh = 10 M (cid:12) , d = 21 .
49 R (cid:12) , v ∞ = 30 km / s , e = 0 . . first disruption provides the majority of the accreted massfor the disrupting BH, while the the non-disrupting BH ac-cretes mass as it returns to the pericenter of the binary orbit.Therefore, the angle for the disrupting BH is similar to thatof the star’s angular momentum with respect to J bin at thetime of the first disruption, while the non-disrupting BH’sangle is aligned with J bin ; see Section 4.1 and Figure 5. Weobtain S ≈ . and S ≈ . which leads to a final χ eff ≈ . .3.3.4. The Massive Overflow Scenario
The changes in spin magnitude obtained in the scenar-ios discussed previously are expected to be small given that S max ( q = 0 . . . More sizable changes are ex-pected for larger values of q . Motivated by this, we run asimulation in which q = 0 . , which we refer to as the massiveoverflow scenario (MOS). The MOS simulation is character-ized by R τ /d = 0 . and a / R L = 16 . . A comparisonbetween the OS and MOS is shown in Figure 3.Both OS and MOS simulations lead to multiple disrup-tions and result in accretion onto both BHs. However, the ˙M curves shown in Figure 4 are significantly different. Inthe OS, accretion onto the disrupting BH proceeds like in acanonical TDEs, showing a fast rise and a subsequent power-law decay. Accretion onto the non-disrupting BH, which oc-curs as it plummets into the accretion disk around the disrupt-ing BH, is observed to be delayed and increases at a slowerrate. In the MOS panel, accretion onto both BHs occurs at asimilar time and the ˙M curves for both BHs are rather simi- lar yet differ from the canonical TDEs. In this case the firstdisruption was weaker and most of the material was madeavailable to the BHs until after the second disruption (Fig-ure 4). The star gets considerably closer to the BH duringthe second encounter and, as a result, the star is completelydisrupted. In what follows we refer to the disrupting BH asthe one responsible for the second disruption, which providesthe vast majority of the mass supply. The accretion disk thatforms after the second disruption can be seen in the right bot-tom panel of Figure 4 and is observed to be very extended,making it easy for the non-disrupting BH to accrete a sub-stantial amount of material, especially since the binary orbitis highly eccentric and the BH will eventually plunge into theaccretion disk.The mass accreted by the disrupting and non-disruptingBH at the end of the simulation is .
91 M (cid:12) and .
40 M (cid:12) respectively. This leads to S ≈ . at angle . withrespect to J bin for the disrupting BH and S ≈ . at an-gle .
14 rad with respect to J bin for the non-disrupting BH,leading to χ eff ≈ − . . The spin angle of the disruptingBH is consistent with the angle with respect to J bin of thestar at the time of the second disruption. The non-disruptingblack hole accretes the majority of the mass in the plane ofthe binary, as in the OS case. We note that the spin angle inthese interactions can change in due to multiple encounters,as can be clearly seen in Figure 3 for the OS scenario (seeSection 4.1 for further discussion). PIN E VOLUTION OF
LIGO S
OURCES time [days] −10 −9 −8 −7 −6 −5 ̇ M [ M ⊙ ⊙ s ] t = 0.24 days t = 14.28 days OSBH , R f = 0.3R ⊙ BH , R f = 0.3R ⊙ time [days] −10 −9 −8 −7 −6 −5 ̇ M [ M ⊙ ⊙ s ] t = 0.44 days t = 4.68 days MOSBH , R f = 0.3R ⊙ BH , R f = 0.3R ⊙ t =0.24 daysFallback Radius t = 14.28 daysFallback Radius t =0.44 daysFallback Radius t =4.68 daysFallback Radius Figure 4 . The mass accretion histories experienced by both BHs in the OS and MOS. Here R f is the Fallback Radius of the surface at which themass flux is calculated. All snapshot times are measured from the time the most bound material is accreted.
Left Panel:
The top panel showsaccretion curves for both BHs for the OS. A total of four tidal interactions take place. The first is responsible for feeding significant mass tothe disrupting, while the non-disrupting BH accretes material as it plunges into the disk around the disrupting BH. The bottom panels showsnapshots of the simulation at 0.24 and 14.3 days (shown as vertical lines in the top panel).
Right Panel:
The top panel shows the rate of massaccretion onto both BHs. In contrast to the OS, both BHs accrete promptly due to the larger amount of stellar material available. After thesecond disruption the accretion curves are almost identical. In this case, both BHs to accrete a notable percentage of their own mass during thedisruption. The bottom panels show snapshots of the simulation at two specific times: 0.16 and 4.42 days (shown as vertical lines in the
TopPanel ). 4.
DISCUSSIONThe detection of GW150914 and subsequent LBBHmerger GW observations have opened up many questionsabout LBBH formation history. Individual BH spins withinthe binary are often used to infer the specific formation chan-nel. In this paper we have explored the possibility and con-sequences of a LBBH experiencing a TDE during its life-time. The accretion that follows from a TDE can possiblyspin up each BH and align or anti-align their relative spins.The notion of temporary spin (mis)alignment contrasts withthe usual assumption that BH spins are non-evolving and re-main unaltered from BH formation to merger. The impli-cations of these tidal interactions are discussed as follows:Section 4.1 explores spin evolution from single and multipleTDEs; and Section 4.2 presents the possible observationalsignatures produced by these interactions.4.1.
Spin Evolution
Individual Disruptions
Section 2.2 outlines the possible scenarios for LBBHTDEs, while Section 3.3 shows how the spin magnitude andorientation of each scenario change as a result of these in-teractions. Following the disruption, accretion disks formaround either one or both BHs as shown in Figure 5. The an-gular momentum distribution of material is initially definedby the orbit of the star before disruption, yet the disk ori- entation can be tilted as the stream is torqued by the binary(Coughlin et al. 2017). The misalignment between J bin and J disk is expected to induce a precession of the accretion diskitself (Nixon & King 2016). The binary should, over longertimescales, induce a warped configuration in the disk with amagnitude depending on the local viscosity. If the accretiondisks are misaligned with respect to the rotation axis of aKerr BH, it will be also subject to Lense-Thirring precession(Bardeen & Petterson 1975). The reader is reminded herethat a particular LBBH experiencing a TDE might not neces-sarily merge and that these interactions are expected to onlytemporarily alter the spin orientation of the binary. WhileTDE interactions will undoubtably change the spin magni-tude of the the accreting BHs, subsequent interactions, ex-pected to take place preferentially with other BHs, will fur-ther modify χ eff .In Section 3.3 we discussed how the accreted spin can goalong J bin or J disk depending on the particular scenario. • For the SS, the disrupting BH is the only one that ac-cretes significant stellar debris. The accreted spin isobserved to be in the direction of J disk at approxi-mately .
75 rad , which is set by the angular momen-tum of the star at the time of disruption. • For the CS, the accreted spin of both BHs will bealigned with J disk . At the time of disruption, J disk hasan angle of about . with respect to J bin . As the L OPEZ J R . ET AL . t * = 393400R τ t * =770R τ t * = 21.5R τ Figure 5 . The structure of the accretion disks formed during the circumbinary, overflow, and massive overflow scenarios. Here t ∗ and R τ denotethe time in units of the star’s dynamical timescale and the individual tidal radius for each panel respectively. Left Panel:
Snapshot showing theaccretion disk structure at the end of the circumbinary scenario (CS) simulation. The angle of the angular momentum of the disk, J disk , relativeto the orbital angular momentum of the binary, J bin , is about . . Middle Panel:
Snapshot showing the accretion disk structure after thethird TDE (out of a total of four before full disruption) in the overflow scenario (OS). The angle of J disk relative to J bin is approximately .
69 rad . Right Panel:
Snapshot of the accretion disk after the initial TDE in the massive overflow scenario (MOS) case. The angle of J disk relative to J bin is ≈ .
85 rad . Accreted Spin (counter-clockwise)Accreted Spin (clockwise) J bin Single Scenario [SS]
Accreted Spin (clockwise)Accreted Spin (counter-clockwise) J bin Circumbinary Scenario [CS]
Figure 6 . Diagram illustrating how the accreted spin directions of the BHs are set by either J disk or J bin . The black and red arrows show thespin of the BHs expected from accretion of the stellar debris, whose angular momentum distribution can be clockwise or counter-clockwise. Left Panel:
In the SS, a single BH TDE occurs and only the disrupting BH accretes material. The resulting direction of the BH spin is expectedto be aligned with the angular momentum of the disk J disk . Right Panel:
In the CS disruption, a circumbinary disk is formed which allowsboth BHs to accrete material with similar specific angular momentum. stream of the most bound material returns to pericen-ter, the binary torques J disk to an angle of ≈ .The torqued stream is responsible for supplying thevast majority of the mass to the disk. As can be seenin the left panel of Figure 5, the initial stream remainsin the disruption plane. • For the OS, the accreted spin of the disrupting BH is inthe direction of J disk at the time of the initial disrup-tion (at .
58 rad ) while the accreted spin of the non-disrupting BH is aligned with J bin at angle of .
24 rad .The first disruption supplies the disrupting BH with themajority of the accreted mass. The middle panel ofFigure 5 shows the disk formed by the third disruption(out of a total of four) whose angle of J disk is .
69 rad with respect to J bin . • Contrary to the OS where a single BH is responsible for multiple disruptions, the MOS has disruptions oc-curring onto both BHs sequentially. Out of the twototal disruptions, the second and final disruption con-tributes the majority of mass accreted by the disruptingBH such that the accreted spin is aligned with J disk atan angle of . and the non-disrupting BH accretesspin in the direction of J bin at an angle .
14 rad . Theright panel of Figure 5 shows the disk arising from thefirst disruption at an angle for J disk of .
85 rad withrespect to J bin . • For the OS and MOS, where multiple disruptions arepossible, the angle of J disk in Figure 5 are differentfrom the final angular momentum distribution of thedisk. This is because the orientation of disk changesafter each disruption as a result of the chaotic nature ofthe three-body dynamics. The disruption resulting in PIN E VOLUTION OF
LIGO S
OURCES Aligned J bin Case OneCase Two
Anti-Aligned
Case OneCase Two J bin J bin J bin J bin J bin Accreted Spin (counter-clockwise)Non-Disrupting BH Accreted Spin J bin Accreted Spin (clockwise)
Overflow Scenario [OS]
Figure 7 . Two cases are depicted that produce aligned BH spins as well as two cases which give rise to anti-aligned BH spins, all of thembelonging to the OS. In all cases, the non-disrupting BH accretes spin in the direction of J bin , due to the density gradient it encounters when itenters the accretion disk. The vectors in the left and right panels indicate the velocity while the vectors in the bottom panels indicate the spinangular momentum unless noted otherwise. Left and Right Panels:
A view of the LBBH orbital plane before and after a star is disrupted. Ineach side panels two cases are depicted for the star’s orbital motion, which determines the orientation of J disk . Bottom Panel:
A side view ofthe LBBH with the final accreted spin directions from both the aligned and anti-aligned configurations. The spin of the non-disrupting BH isalways oriented in the direction of J bin . The alignment or anti-alignment of the BH spins is thus mainly determined by the motion of the starbefore it gets disrupted. OPEZ J R . ET AL .the most accretion will nonetheless determine the finalorientation of the BH spins.In general, for a subset of LBBH TDEs there is a possi-bility of relative alignment or anti-alignment between the in-dividual BH spins. Alignment or lack thereof is set by thespecific conditions of the stellar disruption as well as by theensuing orbital dynamics of the binary, as shown in Figures6 and 7. For the SS, the interaction is similar to a single BHTDE and only the disrupting BH accretes material and will,as a result, be spun up. Therefore, there will be no spin align-ment between the BHs at the end of the TDE. In this case, thethe spin direction of the accreting BH will be aligned with J disk (Figure 6). For the CS, the accretion disk is expectedto form outside of the binary such that the spin directionsof both accreting BHs will be similar and aligned with J disk (Figure 6). In the OS, accretion onto each BH is more com-plicated with the possibility of alignment or anti-alignment.In the case of a single passage disruption, the spin of thenon-disrupting BH will increase in the direction of J bin asmaterial is accreted. This is because a steep density gradientis encountered by the BH when it enters the disk region, asillustrated in Figure 7.The left panels of Figure 7 shows two cases that producealigned BH spins: • the star is disrupted outside of the LBBH in the direc-tion of the orbital motion, and • the star is disrupted inside the LBBH moving againstthe orbital velocity.The right panels of Figure 7 shows two cases that result inanti-alignment: • the star is disrupted outside the LBBH moving againstthe orbital motion, and • the star is disrupted inside the LBBH in the directionof the orbital velocity.We have discussed, in the context of LBBHs, the dynam-ics and subsequent accretion of stellar debris after a TDE. Inall the scenarios, we expect the direction of the star relativeto binary at the moment of disruption to be an essential pa-rameter in determining the resultant BH spins. To this end,we perform a large set of numerical scattering experimentsusing the N -body code developed by Samsing et al. (2014)in order to study the distribution of relative angles betweenthe star’s velocity and the binary orbital velocity upon dis-ruption. The relative angle distributions are plotted in Figure for a sun-like star disrupted by a (cid:12) equal mass BBHwith e = 0 . . From the scattering experiments we concludethat there is no preferred distribution and, as such, we predictequal probability for alignment and anti-alignment in the OS.It is expected that LBBHs will experience multiple interac-tions before merging (e.g., Rodriguez et al. 2016a) and assuch, any temporary alignment might be erased before co-alescence. TDE interactions from assembly to merge willnevertheless alter the spin magnitudes of the the LBBHs. Itis then tempting to try to constraint the spin properties of LBBHs experiencing multiple TDEs and it is to this issuethat we now turn our attention.4.1.2. Multiple TDEs and its Relevance to LBBH Growth
LIGO has uncovered a population of BHs that is more mas-sive than the population known to reside in accreting binaries(Remillard & McClintock 2006). One proposed model forthe formation of LIGO BHs is through hierarchical mergersof lighter BHs. In this case, repeated mergers are expected toleave a clear imprint on the spin of the final merger product(Fishbach et al. 2017; Gerosa & Berti 2017; Rodriguez et al.2018b; Samsing & Ilan 2019). For LBBHs forming hierar-chically, the distribution of spin magnitudes is universal andweighted towards high spins. Such a distribution appears tobe disfavored by current observations. This encourages us toinvestigate spin distributions emerging from LBBHs accret-ing from multiple TDEs.Three sets of simulations are explored here which areaimed at describing the evolution of LBBHs that undergomultiple TDEs before merging. Each simulation starts witha binary with M bh = M bh = 15M (cid:12) disrupting stars with M (cid:63) = 1M (cid:12) ( q = 0 . ). These binaries are assumed to dis-rupt stars isotropically with respect to J bin . Then for eachset of simulations we change the initial χ eff , which is pre-sumed to be set at BH formation or by the early disruption ofa more massive star when the cluster was younger. Figure 9shows our results. The top panel initializes the binary with χ eff = 0 , while the middle and bottom panels start the binarywith χ eff = 0 . and χ eff = 0 . , respectively. For simplicity,we assume the stars are on parabolic orbits and are fully dis-rupted in one passage. This results in a total mass accretedof about . (cid:63) per event, which is modified by an accretionefficiency that is dependent on the spin of the BH at the timeof disruption. This is done in order to account for the ra-diated energy required for a particle at the innermost stablecircular orbit to fall into the BH as described in Bardeen et al.(1972) and Misner et al. (2017). Figure 9 shows that if LIGOsources are built up through TDEs, | χ eff | (cid:46) . (see alsoMandel 2007). Furthermore, we show that an initial χ eff canbe significantly reduced if BH growth in the binary is furtherpromoted by TDEs.4.2. Observable Signatures
A primary source of interest of TDE interactions has beentheir prospects as transients sources. These tidal interactionsfeed material to the BH at rates that are orders of magni-tude above the Eddington photon limit (Figure 4). The to-tal energy, however, is similar from that of other phenom-ena encountered in astrophysics, and is in fact reminiscentof that released in gamma-ray bursts (GRBs; Gehrels et al.2009) and canonical TDE jets (e.g., De Colle et al. 2012).One attractive energy extraction mechanism in these sys-tems, which helps circumvent the Eddington restriction, isthe launching of a relativistic jet (Ramirez-Ruiz & Ross-wog 2009; Giannios & Metzger 2011). Such flows areable to carry both bulk kinetic energy and ordered Poynt-ing flux, which allows high energy radiation to be producedat large distances from the source, where the flow is opti-cally thin (e.g., MacLeod et al. 2014). The corresponding
PIN E VOLUTION OF
LIGO S
OURCES d = 0.1 AU d = 1.0 AU d = 0.316 AU N ε > ε < N θ θ Figure 8 . The distributions of relative angles θ between the stellar velocity vector and the binary orbital velocity vector upon disruption. Similarto Figure 1, the orbital trajectories are calculated for a sun-like star ( M ∗ = 1 M (cid:12) , R ∗ = R (cid:12) ) interacting with a (cid:12) equal mass LBBHwith e = 0 . . The panels show the distribution of θ for different binary separations ( d = 1 . . τ , d = 0 . . τ , d = 0 . . τ ). The trajectories for bound ( ε < ) and unbound ( ε < ) encounters are plotted separately. For comparison, anisotropic θ distribution is shown ( gray curve). Number of Disruptions e ff Mean Number of Disruptions e ff Number of Disruptions e ff Figure 9 . Three sets of simulations are shown, which are aimed atinvestigating the evolution of LBBHs that undergo multiple TDEs.In all cases we plot χ eff as a function of the number of disruptions.All simulations start with a LBBH with M bh = M bh = 15M (cid:12) disrupting M (cid:63) = 1M (cid:12) stars. The disruptions are assumed to beisotropically distributed. For each case, we change the initial χ eff .The top , middle and bottom panel start the binary with χ eff = 0 , χ eff = 0 . and χ eff = 0 . , respectively. beamed emission offers a promising observational signatureof LBBHs due to its expected high luminosity.Figure 10 shows the predicted luminosities for stars dis-rupted by LBBHs assuming that the jet power traces the mass supply to the BH: L j ∝ ˙M c . For comparison wealso plot the luminosities and durations of long γ -ray bursts(LGRBS), jetted TDEs from galactic nuclei as well as thosefrom the newly emerging class of ultra-long GRBs: GRB101225A, GRB 111209A, and GRB 121027A (Levan et al.2014). These ultra-long GRBs reach peak X-ray luminosi-ties of ≈ erg s − and show non-thermal spectra that isreminiscent of relativistically beamed emission. The derivedproperties of these LBBH TDEs appear to place them be-tween ultra-long GRBs and jetted TDEs from galactic nu-clei. Our ability to classify long duration transients as eventsemanating from LBBHs or massive BHs in galactic nuclei islikely to remain a challenge. One alternative in the near termis to search at the astrometric positions of these long tran-sients and see whether they are coincide with galactic cen-ters.Another idea is to look for interruptions in the observedlight curve caused by the binary companion, from whichone could extract the orbital time of the disrupting BBH andthereby its orbital parameters (e.g., Liu et al. 2014). The rel-ativistically beamed emission from these events is the onlycomponent that might be readily detectable since the diskemission is expected to be Eddington limited. We there-fore conclude that one avenue for constraining whether ornot LBBHs reside in star clusters is searching for their high-energy signatures. The possibility of collecting a sample ofsuch events in coming years with Swift appears promising,provided that the rate is similar to the LIGO merger rate ofLBBHs (for a detailed discussion on detectability the readeris refer to MacLeod et al. 2014).To get an estimate on the LBBH TDE rate from the GCpopulation we start by computing the rate per GC using Γ TDE ≈ N BBH × η s σ TDE v dis , where N BBH is the number ofBBHs per GC, η s is the number density of single stars, σ TDE is the TDE cross section, and v dis is the cluster velocity dis-persion. The cross section σ TDE can be written as a productof the binary-single interaction cross section and the proba-2 L
OPEZ J R . ET AL . t [s] l og L [ e r g / s ] LGRBsULGRBs TDEs
OS SS MOS
Figure 10 . The luminosity and duration of high energy transients,adapted from Levan et al. (2014). Shown are the predicted lu-minosities of three of the scenarios for LBBH TDEs discussed inthis paper, assuming L j ∝ ˙M and a 10% radiative efficiency. Forcomparison we plot the observed high-energy properties of GRBsand jetted TDEs. The timescales and durations of LBBH TDEs arewell removed from typical long GRBs, but lie between those of theemerging class of ultra-long GRBs and jetted TDEs. bility for an interaction to result in a TDE (e.g. Samsing et al.2017), i.e. σ TDE ≈ σ bs × P TDE . Assuming the gravitationalfocusing limit for σ bs and P TDE ≈ τ / a one finds, Γ gal . TDE ≈ − yr − (cid:18) η s pc − (cid:19) (cid:18) M bh (cid:12) (cid:19) / (cid:18) km/s v dis (cid:19) , where this rate is per galaxy ( LBBHs per GC, and
GCsper galaxy) derived for solar type stars ( (cid:12) , (cid:12) ) inter-acting with LBBHs of equal mass. This is about one orderof magnitude smaller than the LIGO merger rate, and will befurther observationally suppressed due to the expected beam-ing. Therefore, we expect observations of beamed LBBHTDEs to be lower than the inferred LBBH merger rate. How-ever, if one instead considers stellar tidal disruptions by sin-gle BHs in GCs the rate of beamed TDEs is higher by a factor roughly given by the number ratio of single BHs to the num-ber of LBBHs, Γ BHTDE ≈ Γ LBBHTDE × N BH N BBH , where Γ BHTDE ( Γ LBBHTDE ) is the rate from single (binary) BHstellar disruptions. Assuming the fraction of LBBHs to beat the percent level then this leads to that the rate of stellarsingle BH TDEs is ≈ − yr − per galaxy, which is muchcloser to observable limits. This scenario was recently stud-ied in Perets et al. (2016), and might also be used to con-strain the BH population that later forms LBBHs. We notethat our estimate might be at the optimistic side compared tothe rates derived in Perets et al. (2016), but any of these es-timates should be taken with caution and more sophisticated N -body methods must be used to explore this further.Irrespective of current uncertainties, the detection or non-detection of long duration transients from BH and LBBHstellar disruptions should offer strong constraints on the pop-ulation of LBBHs and the nature of the stellar clusters thathost them. In an upcoming paper we explore what the char-acteristic LBBH orbital parameters are for different clustertypes, as well as what we can learn about the dynamical for-mation of LBBH GW sources from observing the associatedpopulation of BH and LBBH TDEs.ACKNOWLEDGMENTSThe authors thank S. Schrøder, T. Fragos, B. Mockler, S. I.Mandel, W. Farr, C. Miller, D. J. D’Orazio, K. Hotokezaka,M. Gaspari and A. Askar for stimulating discussions. MLJRacknowledges that all praise and thanks belongs to Allah (anybenefit is due to God and any shortcomings are my own).ERR acknowledge support from the DNRF (Niels Bohr Pro-fessor) and NSF grant AST-1615881. JS acknowledges sup-port from the Lyman Spitzer Fellowship. The authors furtherthanks the Niels Bohr Institute for its hospitality while partof this work was completed, and the Kavli Foundation andthe DNRF for supporting the 2017 Kavli Summer Program.REFERENCES Abbott, B. P. et al. 2016, Physical Review Letters, 116, 061102, ADS,1602.03837Antonini, F., & Rasio, F. A. 2016, ApJ, 831, 187, ADS, 1606.04889Bardeen, J. M. 1970, Nature, 226, 64, ADSBardeen, J. M., & Petterson, J. A. 1975, ApJL, 195, L65, ADSBardeen, J. M., Press, W. H., & Teukolsky, S. A. 1972, ApJ, 178, 347, ADSBartos, I., Kocsis, B., Haiman, Z., & Márka, S. 2017, ApJ, 835, 165, ADS,1602.03831Batta, A., Ramirez-Ruiz, E., & Fryer, C. 2017, ApJ, 846, L15, ADSBelczynski, K., Holz, D. E., Bulik, T., & O’Shaughnessy, R. 2016, Nature,534, 512, ADS, 1602.04531Bird, S., Cholis, I., Muñoz, J. B., Ali-Haïmoud, Y., Kamionkowski, M.,Kovetz, E. D., Raccanelli, A., & Riess, A. G. 2016, Physical ReviewLetters, 116, 201301, ADS, 1603.00464Bonnerot, C., Rossi, E. M., Lodato, G., & Price, D. J. 2016, MNRAS, 455,2253, ADSCannizzo, J. K., Lee, H. M., & Goodman, J. 1990, ApJ, 351, 38, ADSCarr, B., Kühnel, F., & Sandstad, M. 2016, PhRvD, 94, 083504, ADS,1607.06077 Cholis, I., Kovetz, E. D., Ali-Haïmoud, Y., Bird, S., Kamionkowski, M.,Muñoz, J. B., & Raccanelli, A. 2016, PhRvD, 94, 084013, ADS,1606.07437Coughlin, E. R., Armitage, P. J., Nixon, C., & Begelman, M. C. 2017,MNRAS, 465, 3840, ADS, 1608.05711De Colle, F., Guillochon, J., Naiman, J., & Ramirez-Ruiz, E. 2012, ApJ,760, 103, ADS, 1205.1507de Mink, S. E., Cantiello, M., Langer, N., Pols, O. R., Brott, I., & Yoon,S.-C. 2009, A&A, 497, 243, ADS, 0902.1751de Mink, S. E., & Mandel, I. 2016, MNRAS, 460, 3545, ADS, 1603.02291Dominik, M., Belczynski, K., Fryer, C., Holz, D. E., Berti, E., Bulik, T.,Mandel, I., & O’Shaughnessy, R. 2012, ApJ, 759, 52, ADS, 1202.4901—. 2013, ApJ, 779, 72, ADS, 1308.1546D’Orazio, D. J., & Loeb, A. 2017, ArXiv e-prints, ADS, 1706.04211Downing, J. M. B., Benacquista, M. J., Giersz, M., & Spurzem, R. 2010,MNRAS, 407, 1946, ADS, 0910.0546—. 2011, MNRAS, 416, 133, ADS, 1008.5060Eggleton, P. P. 1983, ApJ, 268, 368, ADSEvans, C. R., & Kochanek, C. S. 1989, ApJ, 346, L13, ADS
PIN E VOLUTION OF
LIGO S
OURCES13