Gravitational waves can provide crucial insights about the environments in which black holes live. In this work we use numerical relativity simulations to study the behavior of selfinteracting scalar wavelike dark matter clouds accreting onto isolated and binary black holes. We find that repulsive selfinteractions smoothen the spike of an isolated black hole and saturate the density. Attractive selfinteractions enhance the growth and result in more cuspy profiles but can become unstable and undergo explosions akin to the superradiant bosenova that reduce the local cloud density. We quantify the impact of selfinteractions on an equalmass black hole merger by computing the dephasing of the gravitationalwave signal for a range of couplings. We find that repulsive selfinteractions saturate the density of the cloud thereby reducing the dephasing. For attractive selfinteractions the dephasing may be larger but if these interactions dominate prior to the merger the dark matter can undergo bosenova during the inspiral phase disrupting the cloud and subsequently reducing the dephasing.
['Astrophysics and astroparticle physics', 'General Relativity and Quantum Cosmology', 'Astrophysics - Cosmology and Nongalactic Astrophysics', 'High Energy Physics - Theory']
Selfinteracting scalar dark matter around binary black holes
{'Self-interacting scalar dark matter around binary black holes': "Josu C. Aurrekoetxea, 1, ∗ James Marsden, 1, † Katy Clough, 2, ‡ and Pedro G. Ferreira 1, §  \n1 Astrophysics, University of Oxford, DWB, Keble Road, Oxford OX1 3RH, United Kingdom 2 Geometry, Analysis and Gravitation, School of Mathematical Sciences, Queen Mary University of London, Mile End Road, London E1 4NS, United Kingdom  \nGravitational waves can provide crucial insights about the environments in which black holes live. In this work, we use numerical relativity simulations to study the behaviour of self-interacting scalar (wave-like) dark matter clouds accreting onto isolated and binary black holes. We find that repulsive self-interactions smoothen the 'spike' of an isolated black hole and saturate the density. Attractive self-interactions enhance the growth and result in more cuspy profiles, but can become unstable and undergo explosions akin to the superradiant bosenova that reduce the local cloud density. We quantify the impact of self-interactions on an equal-mass black hole merger by computing the dephasing of the gravitational-wave signal for a range of couplings. We find that repulsive self-interactions saturate the density of the cloud, thereby reducing the dephasing. For attractive selfinteractions, the dephasing may be larger, but if these interactions dominate prior to the merger, the dark matter can undergo bosenova during the inspiral phase, disrupting the cloud and subsequently reducing the dephasing.", 'I. INTRODUCTION': "In the absence of observations of weakly interacting massive particles (WIMPs) [1], light scalar fields like axions have undergone a resurgence as potential candidates for the dark matter (DM) of the Universe [2] (see [3-5] for reviews). Such light candidates exhibit distinctive signatures on astrophysical scales, potentially influencing the growth of cosmic structure, the cores of galaxies and the environments of compact objects like neutron stars and black holes (BHs). In the last case, which is the focus of this work, detecting the effects of the environment on the gravitational wave signal from merging sources offers one observational possibility [6-30].  \nAround spinning BHs, the density of scalar fields can grow through the superradiant instability - extracting energy and angular momentum from the BH [31-36] (see [37] for a review). However, in the non spinning Schwarzschild case these fields can also undergo more mundane processes like accretion, where scalar fields are drawn into the deep gravitational wells of the BH [38-43]. For a minimally coupled field obeying the Klein-Gordon equation, accretion leads to the accumulation of a density profile around the BH, with overdensities akin to the 'spikes' of heavy particle DM candidates [44-51]. Unlike in the particle case, if the Compton wavelength of the scalar field is close to the horizon radius of the BH, the wavelike nature of the scalar field manifests. In the case of solar mass or supermassive BHs, this corresponds to scalar field masses m ≈ 10 -11 eV and m ≈ 10 -17 eV, respectively.  \nThese differences change the impact of wave-like and particle-like candidates on gravitational waves emitted from binary black holes (BBHs). In the case of equal-mass mergers, previous work has indicated that the motion of the binary disrupts particle DM spikes, and thus has little or no effect on the gravitational wave signal near merger [52-55].  \nHence, the focus is on extreme/intermediate mass ratio inspirals (EMRIs/IMRIs), where the smaller secondary object does not disrupt the spike of the more massive primary BH. For scalar DM, on the other hand, we showed in [56, 57] (see also [58-61]) that scalar accretion onto equal-mass BHs creates an overdensity between them, which persists during the binary merger. This leads to a faster orbital decay (and hence dephasing of the gravitational-wave signal), which we found is maximized when the Compton wavelength of the scalar field λ c = 2 π / µ (where µ = mc / ¯ h ) is of the order of the binary's orbital separation.  \nMuch of our understanding of scalar field DM is based on the simplest model - a massive scalar field with potential µ 2 φ 2 . This approach mirrors the study of particle DM, where most research has focused on non-interacting, perfectly cold DM. As cosmological and astrophysical data have improved, interest has grown in more complex particle DM models that include multiple species, warm initial conditions, and selfinteractions [16]. Self-interactions can significantly affect galaxy formation, transforming a cusped density profile into a cored profile, and have recently been proposed as a solution to the final parsec problem in the particle case [62, 63]. Whilst the impact of interactions will be small at average galactic DM densities, they may play a larger role at higher densities if DM accumulates, e.g. around stars or BHs [46, 62, 64-82]. In the case of scalar DM, self-interactions enter as higher order terms in the potential, so a quartic term λφ 4 is a natural extension that preserves the parity symmetry. The coupling λ parameterises the strength of the self-interactions, with its sign indicating whether they are repulsive or attractive. This type of coupling provides the leading order term for DM with self-interactions, like those expected in the axion potential. Self-interacting scalar fields have been studied in the context of superradiance around isolated BHs, where it was proposed that attractive interactions could lead to rich and unique phenomenology, in particular bosenova explosions in which the cloud is radiated and absorbed into the BH [83-86]. However, more recent work has shown that in the scalar case the bosenova is unlikely to be reached due to superradiance alone, as the growth saturates when higher modes are excited [87-91]. Additional accretion may enhance the rate sufficiently to al-  \nhe instability to take hold [92], and it may play a role in dynamical capture around stars [64] or the higher densities of boson stars [93-95]. Vector bosons are more promising and it has been shown numerically that vector fields with a Higgslike mass mechanism can excite stringy bosenova explosions for sufficiently strong couplings [96, 97], and tensors should in principle be even more explosive, although the models suffer from some pathologies [98].  \nIn this paper we use numerical simulations to study the accretion of self-interacting scalar fields onto BHs. First, we describe our setup and chosen model in Sec. II. In Sec. III, we study the evolution of DM near an isolated BH, which gives us information about the effect self-interactions have on the spike of the primary BH of an EMRI system, and highlights the possibility of achieving a bosenova-like explosion via accretion. In Sec. IV we use numerical relativity simulations to study the growth of self-interacting scalar clouds during an equal-mass BH merger, and quantify the impact these clouds have on the BBH by extracting the gravitational-wave signal during the inspiral, merger and ringdown. We summarize our results and comment on the physical relevance of our findings in Sec. V.", 'II. SCALAR FIELD DARK MATTER MODEL': 'We model DM as a minimally coupled complex scalar field Φ described by the action  \nS = ∫ d 4 x √ -g ( R 16 π G -1 2 ( ∇ α Φ ) ∗ ( ∇ α Φ ) -V ( Φ ) ) , (1)  \nIn particular, we focus on the impact of self-interactions in the DMfield by adding a quartic self-interacting term  \nV ( | Φ | ) = µ 2 2 | Φ | 2 + λ 4 | Φ | 4 , (2)  \nwhere λ is a real, dimensionful coupling constant 1 that parameterises the strength of the self-interactions. The sign of λ specifies whether the self-interactions are repulsive ( λ > 0) or attractive ( λ < 0). If λ > 0, the potential becomes steeper for larger values of the field, which suppresses growth in the field. For λ < 0 on the other hand, the potential opens up asymptotically, making it easier to explore larger values. However, purely attractive couplings result in a potential that is unbounded from below for sufficiently large values of the field | Φ | . To avoid possible issues when studying attractive self-interactions, we instead simulate  \nV ( Φ ) = µ 4 2 λ ( exp [ λ µ 2 | Φ | 2 ] -1 ) (3) ≈ µ 2 2 | Φ | 2 + λ 4 | Φ | 4 + O ( | Φ | 6 ) ,  \nwhich approximates the quartic potential near the minimum but reaches a plateau rather than turning over after the maximum.  \nIn the flat-space, homogeneous limit, and in the absence of self-interactions, modelling the DM with a complex scalar field allows us to obtain a stationary density 2 . The Noether current associated with the U(1) symmetry of the system corresponds to particle number conservation, making the DM stable. On the complex plane, where the scalar field is defined, the behaviour is that of an ellipse with semi-axis and orientation set by the initial conditions. We choose homogeneous initial conditions Φ =( Re [ Φ ] , Im [ Φ ]) such that  \nΦ 0 =( φ 0 , 0 ) , ∂ t Φ 0 =( 0 , µφ 0 ) . (4)  \nThen, although the asymptotic value of the fields oscillates ( Φ ∼ exp [ i µ t ] ), the modulus | Φ ( t ) | remains constant throughout its evolution and the scalar field tracks circular orbits on the complex plane. This behaviour leads to a constant asymptotic energy density ρ ( t ) = ρ 0, with the magnitude determined by the initial conditions  \nρ 0 = µ 2 φ 2 0 + λ 4 φ 4 0 . (5)  \nWechoose the values φ 0 and λ so that the self-interaction term is initially small ( µ 2 ≫ λφ 2 0 ) and ρ 0 ≈ µ 2 φ 2 0 .  \nThe dynamics of the scalar field is given by the KleinGordon equation  \n∇ α ∇ α Φ -∂ V ( | Φ | ) ∂ Φ ∗ = 0 , (6)  \nwhere ∇ µ is the covariant derivative operator associated to the metric g αβ . In the first part of this paper we will focus on solving this equation on a fixed Schwarzschild background, so that the metric g αβ does not evolve. In the second part, we will use numerical relativity to solve the coupled system of the Einstein-Klein-Gordon equations to evolve both the scalar field and the spacetime metric for the case of an equal-mass binary merger.', 'III. DARK MATTER AROUND AN ISOLATED BLACK HOLE': "During an extreme-mass-ratio inspiral, a stellar-mass BH emits gravitational waves as it orbits around a supermassive BH. Given the separation of scales, the smaller BH can be treated as a perturber of the primary BH's spacetime. If the primary has a DM environment, then as the secondary object moves through it, its geodesic motion will deviate from the vacuum case, and this, combined with effects like dynamical friction and accretion, may leave distinctive imprints in the gravitational-wave signal. These signatures are one of  \nthe key targets of the Laser Interferometer Space Antenna [10, 11, 14, 18, 20, 23, 26, 67, 101]. However, characterizing these signatures requires a robust theoretical understanding of the expected environment. Motivated by this, we study the DMoverdensity that forms around an isolated BH for different self-interactions, neglecting the backreaction of the DM onto the spacetime of the BH. Specifically, we use the 3+1 opensource code GRDZHADZHA [102] to solve the Klein Gordon equation on a BH background in Cartesian Kerr-Schild coordinates [103]  \nds 2 = ( ηµν + 2 M r l µ l ν ) d x µ d x ν , (7)  \nwhere ηµν = diag ( -1 , 1 , 1 , 1 ) and M is the mass of the BH. Here l µ =( 1 , x / r , y / r , z / r ) is the ingoing null vector with respect to both g µν and ηµν . Asymptotically the spacetime is flat and in the absence of self-interactions the scalar field oscillates with frequency proportional to its mass Φ ∼ exp [ i µ t ] . Closer to the BH, the field is accreted and its amplitude and hence density grow [38-40]. After a short time, it reaches a quasi-stationary state with a density profile that (in the case of µ M ∼ 1) is proportional to r -3 / 2 . Notably, since our boundary conditions effectively feed in more scalar field over time, in the purely massive case, the overall amplitude of the cloud continues growing and does not stabilize, although the shape of the profile remains the same.  \nWe start from the homogeneous initial conditions in Eqn. (4) with φ 0 = 0 . 2 M , and study the evolution of the selfinteracting scalar field around a BH of mass M . We study cases where the mass of the scalar µ = M -1 and µ = 0 . 1 M -1 to explore lower and higher masses, although in both cases the Compton wavelength of the scalar field is roughly of the same order of the Schwarzschild radius. These values are partly set by the numerics but also cover the wave-like regime where the density near the BH is not strongly suppressed by the pressure effects arising from a large Compton wavelength (unlike relatively low mass scalars µ M ≪ 1). We use a grid resolution of N 3 = 128 3 with 6 levels of refinement, and a physical box size of L = 576 M , which gives a grid spacing of 0 . 07 M at the BH horizon.  \nWe use diagnostics that follow the approach of [56, 104107], and define a current J µ = ξ ν T µ ν in the direction ξ ν , as well as an associated charge and flux  \nQ = -n µ J µ , F = α NiJ i , (8)  \nwhere Ni is the outward normal direction to the surface that bounds the volume. If ξ ν is a Killing vector then ∇ µ J µ = 0 and the change in charge is balanced by a flux through a surface. For the fixed background case ξ ν t = ( 1 , 0 , 0 , 0 ) is a Killing vector, which gives us the density and fluxes across surfaces  \nρ E ≡ Qt = -αρ + β k S k , (9)  \nFE ≡ Ft = Ni ( β i ( αρ -β j Sj ) + α ( β k S i k -α S i ) , (10)  \nwhere Ni = si / √ ( γ jk sjs k ) is the normalised radial unit vector, with si = ( x , y , z ) / r . Here α = 1 / √ -g 00 and β i = g 0 i . The  \nquantities ρ and Si are the components of the stress-energy tensor measured by the normal observers  \nρ = n µ n ν T µν , Si = n µ T µ i , (11)  \nwhere n µ =( -α , 0 , 0 , 0 ) .  \nIn Fig. 1 we plot radial profiles of the density ρ E at different times. The top and bottom panels correspond to µ M = 1 and µ M = 0 . 1 respectively, and each column corresponds to attractive ( λ < 0), free ( λ = 0) and repulsive ( λ > 0) selfinteractions. For µ M = 1, the repulsive case ( λ > 0) exhibits two interesting properties compared to the density profiles of the free and attractive interaction cases: (i) The profile is flatter towards the centre, resulting in a more 'cored' structure. (ii) The density saturates, leading the profile to settle into a stable configuration. In contrast, the attractive ( λ < 0) and non-self-interacting ( λ = 0) cases show minimal differences from each other within the limited time range shown here. This limitation arises because the dynamics of scalar fields in concave potentials are more challenging to resolve numerically when the timescale of the oscillation is short. We also study the low-frequency regime by decreasing the mass of the scalar field to µ M = 0 . 1. We plot the results in the bottom panel of Fig. 1. In this limit, the power laws are approximately the same, but the results show that attractive self-interactions ( λ < 0) enhance the growth of the cloud, while repulsive selfinteractions ( λ > 0) stabilize and saturate the growth. It can be seen that the field density first stabilizes at the horizon, and as accretion continues, it extends out to larger radii. The accretion flux stabilises at a maximum value and the scalar field is able to reach a stationary state.  \nThis saturation and enhancement of the scalar cloud, depending on the sign of the self-interactions, is more apparent in the top panel of Fig. 2, where we plot the evolution of the density at r = 3 M . Initially, the scalar field evolves with small oscillations in the quadratic part of the potential for all three cases, so the growth due to accretion shows negligible differences at early times. However, as the scalar field amplitude increases, it begins to probe the quartic part of the potential, leading to variations in the density profile. The attractive case shows an enhancement in the growth of the cloud compared to the λ = 0 case, while repulsive self-interactions result in the saturation of the cloud's growth. We interpret these results using the bottom three panels of Fig. 2, where we extract the fluxes at two different radii, r = 3 M and r = 20 M . The difference between the fluxes, ∆ F , indicates the buildup of energy within the volume. The panels show that for the attractive case, ∆ F grows faster than in the λ = 0 case, leading to an enhancement in the cloud's density. In contrast, for the repulsive case, the outer and inner fluxes balance perfectly, causing the density profile to saturate and remain constant.  \nWe understand the saturation of the repulsive case as follows. The stability of the scalar cloud is related to the balance of the Newtonian gravitational potential, the gradient pressure of the scalar field, and the pressure induced by the self-interactions λ . When the scalar field probes the quadratic part of the potential, the gravitational potential dominates and results in the accretion and growth of the cloud. As we have seen, this implies that the scalar field increases its amplitude  \nFIG. 1. Growth of the dark matter density profile for different values of the mass µ and self-interaction λ . The top (bottom) panel shows the high (low) mass case as labelled. The left panels plot the evolution of the radial profiles for attractive self-interactions λ = -150 ( λ = -0 . 003), the middle for the purely massive cases and the right panel for repulsive cases λ = 1 . 6 ( λ = 0 . 015). The different lines show the field density profiles over time { t 1 , t 2 , t 3 } = { 100 , 200 , 220 } M ( { t 1 , t 2 , t 3 } = { 200 , 600 , 1000 } M ), normalised by the initial (homogeneous) density. The black dashed line shows a power law fit to the lines. We expect the power law profiles to continue to propagate outwards over time, but are limited in the length of our simulations by the domain size.  \n<!-- image -->  \nFIG. 2. Top panel: Density of the dark matter cloud at r = 3 M for different values of λ = -0 . 003, λ = 0 and λ = 0 . 015. Bottom panels: Fluxes at radii r = 3 M and r = 20 M . The growth of the density can be understood from the differences in the fluxes between these radii, which diverge, saturate or reach zero for each of the cases.  \n<!-- image -->  \nand starts probing the quartic part of the potential. In the λ > 0 case, the interaction potential creates a repulsive 'pressure' that allows the scalar cloud to reach hydrostatic equilibrium, with the repulsive self-interaction balancing the self-gravity.  \nThe interaction changes sign and becomes attractive when λ < 0, so the only contribution that can lead to hydrostatic equilibrium is the gradient pressure of the scalar field. Attractive self-interactions enhance the growth of the cloud, but the gradient pressure cannot compete with the attractive force due to self-interactions. When the density becomes too high, the system undergoes a rapid collapse - a bosenova-like event. This collapse increases the density and potential energy, leading to an explosive ejection of energy from the cloud and a subsequent reduction in density. The inset panel of Fig. 3, showing the evolution of density at r = 3 M , illustrates this process. The evolution of the density profile is shown in the top panel, with blue to green colours (and transparency) depicting time. The density starts homogeneous but grows near the BH. The density increases until t ≈ 1000 M , when it suddenly collapses and the cloud xplodes, radiating part of its mass outwards. In the bottom panel we plot the evolution of the radial fluxes. Initially, the fluxes are positive everywhere, so that there is an ingoing flux of DM towards the BH. The density grows, so the field explores flatter regions of the potential, so that the density and flux is enhanced. Eventually, it reaches a maximum at r ≈ 25 M , after which the flux becomes  \nFIG. 3. Evolution of the density and fluxes for the bosenova explosion for λ = -0 . 003. The top (bottom) panel shows the radial profiles of the density (radial flux). The inset panel tracks the density value at r = 3 M . Positive and negative fluxes correspond to ingoing and outgoing radiation, respectively. The density reaches a peak following which the cloud collapses, and the system undergoes a violent burst of energy that propagates outwards.  \n<!-- image -->  \nnegative, resulting in the explosion that reduces the density near the BH. This matter propagates outwards, shown by the the negative fluxes at late times in the bottom panel of the figure.", 'IV. DARK MATTER AROUND AN EQUAL MASS BLACK HOLE MERGER': 'As discussed in the introduction, unlike in the particle DM case, equal-mass BHs immersed in scalars with Compton wavelength similar to the binary separation may retain the environment up to merger. This motivates the study of equalmass binaries in such environments including backreaction, to quantify the potential impact on the GW signal. In such a case, the strong-field dynamics of the system needs to be resolved through numerical relativity. We decompose the line element in the general ADM form [108]  \nds 2 = -α 2 d t 2 + γ i j ( d x i + β i d t )( d x j + β j d t ) . (12)  \nHere γ i j is the 3-dimensional metric, and the quantities α and β i determine the choice of spatial hyperslicings and their coordinates, which follow the 1+log moving puncture gauge [109-112]. The extrinsic curvature tensor Kij = ( 2 D ( i β j ) -∂ t γ i j ) / 2 α is decomposed into a trace K and a traceless part Aij , i.e. Kij = Aij +( 1 / 3 ) K γ i j . We evolve the coupled system of Klein-Gordon and general relativity equations in the CCZ4 formulation [113] using the 3+1 open-source numerical relativity code GRCHOMBO [114-116].  \nWe study the tensor gravitational-wave modes emitted by the BBH. To do so, we extract the Newman-Penrose scalar Ψ 4 with tetrads proposed by [117], projected into spin-weight -2 spherical harmonics,  \nψ lm = ∮ S 2 Ψ 4 ∣ ∣ r = r ex [ -2 ¯ Y lm ] d Ω , (13)  \nwhere d Ω = sin θ d θ d ϕ is the area element on the S 2 unit sphere. The merger or coalescence time for the 10-orbit binary in the absence of a dark matter cloud is ¯ t c ≈ 2000 M , defined as when the amplitude of the dominant mode | ψ 22 | peaks. In Ref. [56], we showed that scalar DM around the binary induces a faster orbital decay, so that the gravitationalwave signal experiences a dephasing/delay of the coalescence time  \n∆ t c ≡ t c -¯ t c . (14)  \nIn particular, we found that the dephasing is maximized when the Compton wavelenght of the scalar particle matches with the initial orbital separation of the binary, 2 π / µ ≈ d 0. For our 10-orbit binary, where d 0 ≈ 12 M , this corresponds to a scalar particle with mass µ ≈ 0 . 5 M -1 .  \nIn this work we investigate the impact that repulsive and attractive self-interactions have in the growth of the DM cloud and the dephasing of the binary. For λ > 0, we evolve the system with the quartic potential in Eqn. (2). For λ < 0, we simulate the potential in Eqn. (3) to avoid unbounded values for large amplitudes of the field. We start from homogeneous initial conditions given by Eqn. (4) and choose BowenYork initial data for Aij for the boosted BHs [118]. We solve the Hamiltonian and momentum constraints using the CTTKhybrid method [119], so that the trace of the extrinsic curvature tensor is initially K 2 0 = 24 π G ρ 0. We fix ρ 0 ≈ 10 -9 M -2 , µ = 0 . 43 M -1 , and vary both λ and φ 0 such that ρ 0 is unchanged 3 via Eqn. (5). This approach ensures that we only need to solve the constraints once and the BHs have identical initial masses, velocities and directions for all our different runs. Then, any observed differences arise solely from their interaction with the DM cloud, and depend on the relative enhancements in the cloud density achieved for different scalar masses and interactions.  \nIn Fig. 4, we plot the evolution of the density (left to right) for three different self-interactions. The top and middle panels show the evolution of the DM cloud for attractive ( λ = -10 -4 ) and without self-interactions ( λ = 0). The profiles remain qualitatively unchanged, with both forming a similar overdensity between the BHs. In contrast, for repulsive self-interactions ( λ = 0 . 025, bottom panel), the density of the DM cloud saturates, and there are significant changes in the structure of the cloud. Unlike the other two cases, repulsive self-interactions suppress the formation of the central overdensity. In Fig. 5, we plot the relative dephasing  \nFIG. 5. Decrease in the coalescence time ∆ t c ≡ t c -¯ t c extracted from the gravitational-wave signal, for different self-interaction scenarios. For positive values of λ (repulsive interactions), the dephasing is rapidly lost as the self-interaction increases due to the saturation of the cloud density. For small negative λ (attractive interactions), the density of the cloud is enhanced and the binary experiences a larger dephasing. However, for sufficiently negative λ , the bosenova-like event decreases the density and the dephasing is again reduced.  \n<!-- image -->  \nFIG. 4. Evolution of the dark matter cloud for λ = -10 -4 (top), λ = 0 (middle) and λ = 0 . 025 (bottom). We plot 2D slices of the density. We see that in the repulsive case the overdensity around the binary is suppressed. In the attractive case shown here, the density does not become high enough to trigger a bosenova before the merger occurs, so its evolution is similar to the purely massive case, with only a slight enhancement in the density. A movie of the simulations can be found in [120].  \n<!-- image -->  \nof the gravitational-wave signal. Repulsive self-interactions saturate the growth of the cloud, resulting in smaller dephasings. For attractive (but small) self-interactions, the dephasing is slightly larger due to the increased density of the cloud.  \nWe can quantify more precisely the cause of this dephasing by computing the charges and fluxes associated with the currents in the time direction ξ ν t =( 1 , 0 , 0 , 0 ) , angular direction ξ ν φ = ( 0 , -y , x , 0 ) , and radial direction ξ ν r = -( 0 , x , y , z ) / r . The expressions for the time direction are in Eqns. (9-10), and the angular and radial directions  \nQ { φ , r } = Si ξ i { φ , r } , (15)  \nF { φ , r } = -Ni β i Sj ξ j { φ , r } + α NiS i j ξ j { φ , r } , (16)  \nwhere again Ni = si / √ ( γ jk sjs k ) and si =( x , y , z ) / r . We also track the flux through two inner surfaces that move together with the BHs, which introduces additional advection terms to the flux F BH = α N BH i J i -N BH i β i ( Q -S / 2 ) .  \nAs these directions do not correspond to Killing vectors, there must be an exchange of the charge between matter and curvature, which can be estimated by computing the source terms associated to these currents  \nS = α T µ ν ∇ µξ ν . (17)  \nIt is this quantity that corresponds to gravitational forces in the Newtonian limit, and that quantifies the way in which momentum is extracted from the binary by the matter 4 . The expres-  \nsions for the source terms for the energy, angular and linear momentum currents are  \nSt = -ρ∂ t α + Si ∂ t β i + α S i j ∂ t γ i j  \n2 (18) S { φ , r } = α S µ ν ∂µξ ν { φ , r } + α S µ ν ( 3 ) Γ ν µσ ξ σ { φ , r } -S νβ i ∂ i ξ ν { φ , r } + S νξ µ { φ , r } ∂µβ ν -ρξ µ { φ , r } ∂µα (19)  \nwhere ∂ t γ i j = -2 α Kij + Di β j + Dj β i , and both ∂ t α and ∂ t β i are given by our moving puncture gauge conditions. The quantities ∂µξ ν are all zero except ∂ x ξ y = 1 and ∂ y ξ x = -1.  \nIn Fig. 6 we plot the time integration of these quantities within a radius of 40 M from the centre of the binary. These quantify the extraction of mass ∆ M , angular momentum ∆ J and radial momentum ∆ Pr from the binary by the matter, which drive the dephasing and the change in the coalescence time ∆ t c since they affect the balance in the angular momemtum J = Mr 2 ( 2 π / T ) . We compare the attractive λ = -10 -4 and repulsive λ = 0 . 025 cases to the λ = 0 case. In these cases, the attractive self-interaction gives the greatest relative dephasing at around 8 . 2%, although it is very similar to the purely massive case. We attribute the slight difference in the coalescence time to the slightly higher values of ∆ J and ∆ Pr in Fig. 6. In the repulsive case, the dephasing is only around 1%: the exchanges in energy and angular momentum are highly  \nFIG. 6. Time integrals of the source terms for the energy, angular and radial momentum currents (top, middle and bottom panels, respectively). These quantify the extraction of mass ∆ M , angular momentum ∆ J and radial momentum ∆ Pr from the binary by the matter, illustrating the effect the dark matter cloud has on the binary. We compare the attractive λ = -10 -4 and repulsive λ = 0 . 025 cases to the λ = 0 case, which are those shown in Fig. 4.  \n<!-- image -->  \nsuppressed (top and middle panels in Fig. 6), whilst the radial exchange of momentum shown in the bottom panel is smaller.  \nIn principle one may naively expect that increasing the attractive self-interactions further enhances the density of the cloud, but Fig. 5 shows that for λ ⪅ -10 -4 , the dephasing ∆ t c decreases. This is because, as in the isolated BH case, when the field reaches sufficiently high values, the non-linear attractive interactions of the scalar field dominate and lead to a bosenova-like event. There is a sudden collapse at the overdensity between the binary that results in a burst of matter, disrupting the cloud, as illustrated in Fig. 7. This behaviour seems to be generic for any kind of attractive self-interaction with λ | Φ | 4 at leading order, but we expect that the exact form of the ejection and recovery to the steady growth phase after the bosenova will depend somewhat on the form of the potential beyond the λ | Φ | 4 regime. The critical value of λ at which it occurs will also depend on factors like the asymptotic density, the scalar mass and the parameters of the binary.', 'V. DISCUSSION': "We have studied the evolution of self-interacting dark matter around black holes using numerical simulations. Motivated by EMRIs, we have explored the accretion-driven growth of a scalar DM overdensity around an isolated Schwarzschild BH. We have found that the profile of the cloud is cored for repulsive self-interactions, and more cuspy when these are attractive. Repulsive self-interactions saturate the density of the cloud, while attractive self-interactions lead to  \na faster growth. However, when attractive self-interactions dominate the dynamics, this enhancement can be temporarily reversed by a bosenova-like event, as the cloud collapses and the system emits a violent burst of energy.  \nWe used numerical relativity simulations that include backreaction onto the metric to study the growth of a DM cloud during an equal-mass BH merger. We have extracted the gravitational waves radiated during the inspiral, merger and ringdown phases and quantified the impact that self-interactions have in the waveform by computing the dephasing of the gravitational-wave signal. We have found that repulsive selfinteractions reduce the amount of dephasing compared to the purely massive case because they suppress the DM density that accumulates around the binary. Attractive selfinteractions enhance it, but if the cloud grows significantly, the binary triggers a bosenova-like event as in the isolated BH case. The cloud is then disrupted by the motion of the BHs, and the dephasing is again reduced relative to the purely massive case.  \nIn our simulations, we are constrained to use larger values of the density in order to see the dephasing effects on the short timescales on which we are able to simulate the merger. Smaller values should result in similar effects but on longer timescales. To put the values we use into context - constraints exist on the total mass of dark matter around the supermassive BH at the centre of our galaxy, putting it at less than 0.1% of the mass of the central BH within 10 4 Schwarzschild radii using constraints from the S2 star [121, 122]. This corresponds to densities of roughly 10 -15 M -2 in our units, which is significantly lower than the value of 10 -9 M -2 that we use as the asymptotic density. However, the densities we use are lower than possible densities from superradiant clouds where 4-10% of the BH mass can be contained within a cloud out to around 1-10 Schwarzschild radii [123], which gives values up to ∼ 10 -5 M -2 . They are comparable to those expected from the isolated adiabatic evolution of dark matter spikes [124].  \nWhilst the densities that we study are large compared to galactic average values, the self-interactions that we study here are relatively small and not at a level that conflict with current observational constraints. The dimensionless coupling for the self-interaction ¯ λ for a scalar mass m = ¯ h µ / c is constrained by BBN and structure formation to be [125-127]  \n¯ λ ≲ 1 ( m eV ) 4 . (20)  \nIn our simulations we have fixed the scalar mass to be m ≈ 5 × 10 -17 ( M / 10 6 M ⊙ ) -1 eV, which corresponds to the constraint  \n¯ λ ≲ 6 × 10 -66 ( M 10 6 M ⊙ ) -4 . (21)  \nWe have chosen the parameter λ ≈ O ( 0 -0 . 1 ) M -2 , for which the upper end corresponds to the dimensionless value  \n¯ λ ∼ 10 -89 ( M 10 6 M ⊙ ) -2 , (22)  \nor equivalently a self-interaction energy scale f = √ λ / m ∼  \nFIG. 7. Bosenova explosion and disruption of the cloud during an equal-mass black hole merger for λ = -10 -3 . The density of the dark matter cloud is enhanced due to the attractive self-interactions (compare the middle panel of this figure and those in Fig. 4). However, attractive self-interactions dominate and result in the sudden collapse of the cloud, with a burst of radiation that reduces the cloud density, and hence the dephasing. A movie of the simulation can be found in [120].  \n<!-- image -->  \nM Pl . This interaction strength is below the bound above for all realistic BH masses 5 .  \nBased on our results, self-interactions appear to mainly reduce rather than enhance dephasing effects. Even for the relatively small values of self-interaction studied, the non-linear part of the potential played a significant role, and both attractive and repulsive values mostly resulted in lower dephasing than the purely massive case. We are only considering a short period close to merger, and so given the much longer timescales available for build-up of the cloud prior to this point, the effects could become relevant well in advance of the actual merger. In the case of attractive self-interactions, this could mean bosenova-like explosions reccuring during the inspiral. The absence or observation of such bursts can be used to put constraints on the existence of attractive scalar DM, as has been proposed in [64, 93-95].  \nAs this article was under review, a relevant paper [129] studied the effect of self-interactions in the context of superradiance.  \n- [1] J. Aalbers et al. (LZ), Phys. Rev. Lett. 131 , 041002 (2023), arXiv:2207.03764 [hep-ex].\n- [2] H.-Y. Schive, T. Chiueh, and T. Broadhurst, Nature Phys. 10 , 496 (2014), arXiv:1406.6586 [astro-ph.GA].\n- [3] L. Hui, Ann. Rev. Astron. Astrophys. 59 , 247 (2021), arXiv:2101.11735 [astro-ph.CO].\n- [4] L. A. Ure˜na L'opez, Front. Astron. Space Sci. 6 , 47 (2019).", 'VI. ACKNOWLEDGEMENTS': "We thank Joshua Eby, Ben Elder, Rodrigo Vicente and Miguel Zumalacarregui for useful comments and discussions, and the GRChombo collaboration (www.grchombo.org) for their support and code development work. JCA acknowledges funding from the Beecroft Trust and The Queen's College via an extraordinary Junior Research Fellowship (eJRF). JM acknowledges funding from STFC. KC acknowledges funding from the UKRI Ernest Rutherford Fellowship (grant number ST/V003240/1) and STFC Research Grant ST/X000931/1 (Astronomy at Queen Mary 2023-2026). PGF acknowledges support from STFC and the Beecroft Trust.  \nThis work used the DiRAC Memory Intensive services Cosma8 and Cosma7 at Durham University, managed by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The DiRAC service at Durham was funded by BEIS, UKRI and STFC capital funding, Durham University and STFC operations grants. DiRAC is part of the UKRI Digital Research Infrastructure.  \n- [5] J. C. Niemeyer, (2019), 10.1016/j.ppnp.2020.103787, arXiv:1912.07064 [astro-ph.CO].\n- [6] L. Barack et al. , Class. Quant. Grav. 36 , 143001 (2019), arXiv:1806.05195 [gr-qc].\n- [7] E. Barausse et al. , Gen. Rel. Grav. 52 , 81 (2020), arXiv:2001.09793 [gr-qc].\n- [8] N. Afshordi et al. (LISA Consortium Waveform Working Group), (2023), arXiv:2311.01300 [gr-qc].\n- [9] G. Caneva Santoro, S. Roy, R. Vicente, M. Haney, O. J. Piccinni, W. Del Pozzo, and M. Martinez, Phys. Rev. Lett. 132 , 251401 (2024), arXiv:2309.05061 [gr-qc].\n- [10] V. Cardoso and A. Maselli, Astron. Astrophys. 644 , A147\n- (2020), arXiv:1909.05870 [astro-ph.HE].\n- [11] E. Barausse, V. Cardoso, and P. Pani, Phys. Rev. D 89 , 104059 (2014), arXiv:1404.7149 [gr-qc].\n- [12] N. Yunes, B. Kocsis, A. Loeb, and Z. Haiman, Phys. Rev. Lett. 107 , 171103 (2011), arXiv:1103.4609 [astro-ph.CO].\n- [13] B. Kocsis, N. Yunes, and A. Loeb, Phys. Rev. D 84 , 024032 (2011), arXiv:1104.2322 [astro-ph.GA].\n- [14] C. F. B. Macedo, P. Pani, V. Cardoso, and L. C. B. Crispino, Astrophys. J. 774 , 48 (2013), arXiv:1302.2646 [gr-qc].\n- [15] V. Cardoso and C. F. B. Macedo, Mon. Not. Roy. Astron. Soc. 498 , 1963 (2020), arXiv:2008.01091 [astro-ph.HE].\n- [16] G. Bertone and T. Tait, M. P., Nature 562 , 51 (2018), arXiv:1810.01668 [astro-ph.CO].\n- [17] R. Alves Batista et al. , (2021), arXiv:2110.10074 [astroph.HE].\n- [18] L. Zwick, A. Derdzinski, M. Garg, P. R. Capelo, and L. Mayer, Mon. Not. Roy. Astron. Soc. 511 , 6143 (2022), arXiv:2110.09097 [astro-ph.HE].\n- [19] V. Cardoso and F. Duque, Phys. Rev. D 101 , 064028 (2020), arXiv:1912.07616 [gr-qc].\n- [20] A. Maselli, N. Franchini, L. Gualtieri, T. P. Sotiriou, S. Barsanti, and P. Pani, Nature Astron. 6 , 464 (2022), arXiv:2106.11325 [gr-qc].\n- [21] P. Amaro-Seoane, Living Rev. Rel. 21 , 4 (2018), arXiv:1205.5240 [astro-ph.CO].\n- [22] V. Cardoso, K. Destounis, F. Duque, R. Panosso Macedo, and A. Maselli, Phys. Rev. Lett. 129 , 241103 (2022), arXiv:2210.01133 [gr-qc].\n- [23] P. S. Cole, G. Bertone, A. Coogan, D. Gaggero, T. Karydas, B. J. Kavanagh, T. F. M. Spieksma, and G. M. Tomaselli, Nature Astron. 7 , 943 (2023), arXiv:2211.01362 [gr-qc].\n- [24] A. Krolak and B. F. Schutz, Gen. Rel. Grav. 19 , 1163 (1987).\n- [25] H. Koo, D. Bak, I. Park, S. E. Hong, and J.-W. Lee, Phys. Lett. B 856 , 138908 (2024), arXiv:2311.03412 [astro-ph.GA].\n- [26] M. Boˇskovi'c, M. Koschnitzke, and R. A. Porto, (2024), arXiv:2403.02415 [gr-qc].\n- [27] T. K. Karydas, B. J. Kavanagh, and G. Bertone, (2024), arXiv:2402.13053 [gr-qc].\n- [28] G. M. Tomaselli, T. F. M. Spieksma, and G. Bertone, (2024), arXiv:2403.03147 [gr-qc].\n- [29] G. M. Tomaselli, T. F. M. Spieksma, and G. Bertone, (2024), arXiv:2407.12908 [gr-qc].\n- [30] B. C. Bromley, P. Sandick, and B. Shams Es Haghi, Phys. Rev. D 110 , 023517 (2024), arXiv:2311.18013 [astro-ph.GA].\n- [31] Y. B. Zel'Dovich, Soviet Journal of Experimental and Theoretical Physics Letters 14 , 180 (1971).\n- [32] W. H. Press and S. A. Teukolsky, Nature 238 , 211 (1972).\n- [33] T. J. M. Zouros and D. M. Eardley, Annals Phys. 118 , 139 (1979).\n- [34] S. L. Detweiler, Phys. Rev. D 22 , 2323 (1980).\n- [35] V. Cardoso, O. J. C. Dias, J. P. S. Lemos, and S. Yoshida, Phys. Rev. D 70 , 044039 (2004), [Erratum: Phys.Rev.D 70, 049903 (2004)], arXiv:hep-th/0404096.\n- [36] W. E. East and F. Pretorius, Phys. Rev. Lett. 119 , 041101 (2017), arXiv:1704.04791 [gr-qc].\n- [37] R. Brito, V. Cardoso, and P. Pani, Lect. Notes Phys. 906 , pp.1 (2015), arXiv:1501.06570 [gr-qc].\n- [38] K. Clough, P. G. Ferreira, and M. Lagos, Phys. Rev. D 100 , 063014 (2019), arXiv:1904.12783 [gr-qc].\n- [39] L. Hui, D. Kabat, X. Li, L. Santoni, and S. S. C. Wong, JCAP 06 , 038 (2019), arXiv:1904.12803 [gr-qc].\n- [40] J. Bamber, K. Clough, P. G. Ferreira, L. Hui, and M. Lagos, Phys. Rev. D 103 , 044059 (2021), arXiv:2011.07870 [gr-qc].\n- [41] B. Bucciotti and E. Trincherini, JHEP 11 , 193 (2023),  \n- arXiv:2309.02482 [hep-th].\n- [42] M. de Cesare and R. Oliveri, Phys. Rev. D 108 , 044050 (2023), arXiv:2305.04970 [gr-qc].\n- [43] N. Sanchis-Gual, J. C. Degollado, P. Izquierdo, J. A. Font, and P. J. Montero, Phys. Rev. D 94 , 043004 (2016), arXiv:1606.05146 [gr-qc].\n- [44] P. Gondolo and J. Silk, Phys. Rev. Lett. 83 , 1719 (1999), arXiv:astro-ph/9906391.\n- [45] V. De Luca and J. Khoury, JCAP 04 , 048 (2023), arXiv:2302.10286 [astro-ph.CO].\n- [46] L. Berezhiani, G. Cintia, V. De Luca, and J. Khoury, JCAP 06 , 024 (2024), arXiv:2311.07672 [astro-ph.CO].\n- [47] L. Sadeghian, F. Ferrer, and C. M. Will, Phys. Rev. D 88 , 063522 (2013), arXiv:1305.2619 [astro-ph.GA].\n- [48] O. Y. Gnedin and J. R. Primack, Phys. Rev. Lett. 93 , 061302 (2004), arXiv:astro-ph/0308385.\n- [49] D. Merritt, S. Harfst, and G. Bertone, Phys. Rev. D 75 , 043517 (2007), arXiv:astro-ph/0610425.\n- [50] D. Merritt, Phys. Rev. Lett. 92 , 201304 (2004), arXiv:astroph/0311594.\n- [51] S. L. Shapiro and D. C. Heggie, Phys. Rev. D 106 , 043018 (2022), arXiv:2209.08105 [astro-ph.GA].\n- [52] D. Merritt and M. Milosavljevic, in 4th International Heidelberg Conference on Dark Matter in Astro and Particle Physics (2002) pp. 79-89, arXiv:astro-ph/0205140.\n- [53] G. Bertone and D. Merritt, Phys. Rev. D 72 , 103502 (2005), arXiv:astro-ph/0501555.\n- [54] B. J. Kavanagh, D. Gaggero, and G. Bertone, Phys. Rev. D 98 , 023536 (2018), arXiv:1805.09034 [astro-ph.CO].\n- [55] P. S. Cole, A. Coogan, B. J. Kavanagh, and G. Bertone, Phys. Rev. D 107 , 083006 (2023), arXiv:2207.07576 [astro-ph.CO].\n- [56] J. C. Aurrekoetxea, K. Clough, J. Bamber, and P. G. Ferreira, Phys. Rev. Lett. 132 , 211401 (2024), arXiv:2311.18156 [grqc].\n- [57] J. Bamber, J. C. Aurrekoetxea, K. Clough, and P. G. Ferreira, Phys. Rev. D 107 , 024035 (2023), arXiv:2210.09254 [gr-qc].\n- [58] G. Ficarra, in 55th Rencontres de Moriond on Gravitation (2021) arXiv:2105.05918 [gr-qc].\n- [59] Y.-P. Zhang, M. Gracia-Linares, P. Laguna, D. Shoemaker, and Y.-X. Liu, Phys. Rev. D 107 , 044039 (2023), arXiv:2209.11814 [gr-qc].\n- [60] S. Choudhary, N. Sanchis-Gual, A. Gupta, J. C. Degollado, S. Bose, and J. A. Font, Phys. Rev. D 103 , 044032 (2021), arXiv:2010.00935 [gr-qc].\n- [61] Q. Yang, L.-W. Ji, B. Hu, Z.-J. Cao, and R.-G. Cai, Res. Astron. Astrophys. 18 , 065 (2018), arXiv:1706.00678 [gr-qc].\n- [62] G. Alonso-' Alvarez, J. M. Cline, and C. Dewar, Phys. Rev. Lett. 133 , 021401 (2024), arXiv:2401.14450 [astro-ph.CO].\n- [63] M. S. Fischer and L. Sagunski, (2024), arXiv:2405.19392 [astro-ph.CO].\n- [64] D. Budker, J. Eby, M. Gorghetto, M. Jiang, and G. Perez, JCAP 12 , 021 (2023), arXiv:2306.12477 [hep-ph].\n- [65] S. L. Shapiro and V. Paschalidis, Phys. Rev. D 89 , 023506 (2014), arXiv:1402.0005 [astro-ph.CO].\n- [66] K. Kadota, J. H. Kim, P. Ko, and X.-Y. Yang, Phys. Rev. D 109 , 015022 (2024), arXiv:2306.10828 [hep-ph].\n- [67] A. Boudon, P. Brax, P. Valageas, and L. K. Wong, Phys. Rev. D 109 , 043504 (2024), arXiv:2305.18540 [astro-ph.CO].\n- [68] W.-X. Feng, A. Parisi, C.-S. Chen, and F.-L. Lin, JCAP 08 , 032 (2022), arXiv:2112.05160 [astro-ph.HE].\n- [69] G. Alvarez and H.-B. Yu, Phys. Rev. D 104 , 043013 (2021), arXiv:2012.15050 [hep-ph].\n- [70] H. S. Chia, C. Doorman, A. Wernersson, T. Hinderer, and S. Nissanke, JCAP 04 , 018 (2023), arXiv:2212.11948 [gr-qc].\n- [71] T. F. M. Spieksma, E. Cannizzaro, T. Ikeda, V. Cardoso, and Y. Chen, Phys. Rev. D 108 , 063013 (2023), arXiv:2306.16447 [gr-qc].\n- [72] N. Siemonsen, C. Mondino, D. Egana-Ugrinovic, J. Huang, M. Baryakhtar, and W. E. East, Phys. Rev. D 107 , 075025 (2023), arXiv:2212.09772 [astro-ph.HE].\n- [73] E. Cannizzaro, L. Sberna, A. Caputo, and P. Pani, Phys. Rev. D 106 , 083019 (2022), arXiv:2206.12367 [hep-ph].\n- [74] M. Boskovic, R. Brito, V. Cardoso, T. Ikeda, and H. Witek, Phys. Rev. D 99 , 035006 (2019), arXiv:1811.04945 [gr-qc].\n- [75] J. Marsden, J. C. Aurrekoetxea, K. Clough, and P. G. Ferreira, (2024), arXiv:2403.17595 [gr-qc].\n- [76] M. Y. Khlopov, B. A. Malomed, I. B. Zeldovich, and Y. B. Zeldovich, Mon. Not. Roy. Astron. Soc. 215 , 575 (1985).\n- [77] Y. Chen, X. Xue, and V. Cardoso, (2023), arXiv:2308.00741 [hep-ph].\n- [78] Y. Chen, J. Shu, X. Xue, Q. Yuan, and Y. Zhao, Phys. Rev. Lett. 124 , 061102 (2020), arXiv:1905.02213 [hep-ph].\n- [79] Y. Chen, Y. Liu, R.-S. Lu, Y. Mizuno, J. Shu, X. Xue, Q. Yuan, and Y. Zhao, Nature Astron. 6 , 592 (2022), arXiv:2105.04572 [hep-ph].\n- [80] Y. Chen, C. Li, Y. Mizuno, J. Shu, X. Xue, Q. Yuan, Y. Zhao, and Z. Zhou, JCAP 09 , 073 (2022), arXiv:2208.05724 [hepph].\n- [81] X. Gan, L.-T. Wang, and H. Xiao, (2023), arXiv:2311.02149 [hep-ph].\n- [82] J. H. Kim and X.-Y. Yang, (2024), arXiv:2407.14604 [astroph.CO].\n- [83] J. L. Roberts, N. R. Claussen, S. L. Cornish, E. A. Donley, E. A. Cornell, and C. E. Wieman, Phys. Rev. Lett. 86 , 4211 (2001).\n- [84] A. Arvanitaki and S. Dubovsky, Phys. Rev. D 83 , 044026 (2011), arXiv:1004.3558 [hep-th].\n- [85] H. Yoshino and H. Kodama, Prog. Theor. Phys. 128 , 153 (2012), arXiv:1203.5070 [gr-qc].\n- [86] H. Yoshino and H. Kodama, Class. Quant. Grav. 32 , 214001 (2015), arXiv:1505.00714 [gr-qc].\n- [87] H. Omiya, T. Takahashi, T. Tanaka, and H. Yoshino, Phys. Rev. D 110 , 044002 (2024), arXiv:2404.16265 [gr-qc].\n- [88] H. Omiya, T. Takahashi, T. Tanaka, and H. Yoshino, JCAP 06 , 016 (2023), arXiv:2211.01949 [gr-qc].\n- [89] H. Omiya, T. Takahashi, and T. Tanaka, PTEP 2022 , 043E03 (2022), arXiv:2201.04382 [gr-qc].\n- [90] H. Omiya, T. Takahashi, and T. Tanaka, PTEP 2021 , 043E02 (2021), arXiv:2012.03473 [gr-qc].\n- [91] M. Baryakhtar, M. Galanis, R. Lasenby, and O. Simon, Phys. Rev. D 103 , 095019 (2021), arXiv:2011.11646 [hep-ph].\n- [92] L. Hui, Y. T. A. Law, L. Santoni, G. Sun, G. M. Tomaselli, and E. Trincherini, Phys. Rev. D 107 , 104018 (2023), arXiv:2208.06408 [gr-qc].\n- [93] J. Arakawa, M. H. Zaheer, J. Eby, V. Takhistov, and M. S. Safronova, (2024), arXiv:2402.06736 [hep-ph].\n- [94] J. Arakawa, J. Eby, M. S. Safronova, V. Takhistov, and M. H. Zaheer, (2023), arXiv:2306.16468 [hep-ph].\n- [95] J. Eby and V. Takhistov, (2024), arXiv:2402.00100 [hep-ph].\n- [96] W. E. East and J. Huang, JHEP 12 , 089 (2022), arXiv:2206.12432 [hep-ph].\n- [97] W. E. East, Phys. Rev. Lett. 129 , 141103 (2022), arXiv:2205.03417 [hep-ph].\n- [98] W. E. East and N. Siemonsen, Phys. Rev. D 108 , 124048 (2023), arXiv:2309.05096 [gr-qc].\n- [99] L. Berezhiani and J. Khoury, Phys. Rev. D 92 , 103510 (2015), arXiv:1507.01019 [astro-ph.CO].\n- [100] L. Berezhiani, B. Elder, and J. Khoury, JCAP 10 , 074 (2019),  \n- arXiv:1905.09297 [hep-ph].\n- [101] M. Rahman, S. Kumar, and A. Bhattacharyya, JCAP 01 , 035 (2024), arXiv:2306.14971 [gr-qc].\n- [102] J. C. Aurrekoetxea, J. Bamber, S. E. Brady, K. Clough, T. Helfer, J. Marsden, M. Radia, D. Traykova, and Z. Wang, J. Open Source Softw. 9 , 5956 (2024), arXiv:2308.08299 [grqc].\n- [103] M. Visser, in Kerr Fest: Black Holes in Astrophysics, General Relativity and Quantum Gravity (2007) arXiv:0706.0622 [grqc].\n- [104] K. Clough, Class. Quant. Grav. 38 , 167001 (2021), arXiv:2104.13420 [gr-qc].\n- [105] R. Croft, Class. Quant. Grav. 40 , 105007 (2023), arXiv:2203.13845 [gr-qc].\n- [106] D. Traykova, K. Clough, T. Helfer, E. Berti, P. G. Ferreira, and L. Hui, Phys. Rev. D 104 , 103014 (2021), arXiv:2106.08280 [gr-qc].\n- [107] D. Traykova, R. Vicente, K. Clough, T. Helfer, E. Berti, P. G. Ferreira, and L. Hui, Phys. Rev. D 108 , L121502 (2023), arXiv:2305.10492 [gr-qc].\n- [108] R. L. Arnowitt, S. Deser, and C. W. Misner, Gen. Rel. Grav. 40 , 1997 (2008), arXiv:gr-qc/0405109.\n- [109] C. Bona, J. Masso, E. Seidel, and J. Stela, Phys. Rev. Lett. 75 , 600 (1995), arXiv:gr-qc/9412071.\n- [110] J. G. Baker, J. Centrella, D.-I. Choi, M. Koppitz, and J. van Meter, Phys. Rev. Lett. 96 , 111102 (2006), arXiv:grqc/0511103.\n- [111] M. Campanelli, C. O. Lousto, P. Marronetti, and Y. Zlochower, Phys. Rev. Lett. 96 , 111101 (2006), arXiv:grqc/0511048.\n- [112] J. R. van Meter, J. G. Baker, M. Koppitz, and D.-I. Choi, Phys. Rev. D 73 , 124011 (2006), arXiv:gr-qc/0605030.\n- [113] D. Alic, C. Bona-Casas, C. Bona, L. Rezzolla, and C. Palenzuela, Phys. Rev. D 85 , 064040 (2012), arXiv:1106.2254 [grqc].\n- [114] T. Andrade et al. , J. Open Source Softw. 6 , 3703 (2021), arXiv:2201.03458 [gr-qc].\n- [115] M. Radia, U. Sperhake, A. Drew, K. Clough, P. Figueras, E. A. Lim, J. L. Ripley, J. C. Aurrekoetxea, T. Franc¸a, and T. Helfer, Class. Quant. Grav. 39 , 135006 (2022), arXiv:2112.10567 [grqc].\n- [116] K. Clough, P. Figueras, H. Finkel, M. Kunesch, E. A. Lim, and S. Tunyasuvunakool, Class. Quant. Grav. 32 , 245011 (2015), arXiv:1503.03436 [gr-qc].\n- [117] J. G. Baker, M. Campanelli, and C. O. Lousto, Phys. Rev. D 65 , 044001 (2002), arXiv:gr-qc/0104063.\n- [118] J. M. Bowen and J. W. York, Jr., Phys. Rev. D 21 , 2047 (1980).\n- [119] J. C. Aurrekoetxea, K. Clough, and E. A. Lim, Class. Quant. Grav. 40 , 075003 (2023), arXiv:2207.03125 [gr-qc].\n- [120] 'Self-interacting scalar dark matter around binary black holes,' (2024), https://youtu.be/ZeoMlUHknUM .\n- [121] R. Abuter et al. (GRAVITY), Astron. Astrophys. 657 , L12 (2022), arXiv:2112.07478 [astro-ph.GA].\n- [122] A. Foschi et al. (GRAVITY), Mon. Not. Roy. Astron. Soc. 524 , 1075 (2023), arXiv:2306.17215 [astro-ph.GA].\n- [123] C. A. R. Herdeiro, E. Radu, and N. M. Santos, Phys. Lett. B 824 , 136835 (2022), arXiv:2111.03667 [gr-qc].\n- [124] G. Bertone, A. R. A. C. Wierda, D. Gaggero, B. J. Kavanagh, M. Volonteri, and N. Yoshida, (2024), arXiv:2404.08731 [astro-ph.CO].\n- [125] B. Li, T. Rindler-Daller, and P. R. Shapiro, Phys. Rev. D 89 , 083536 (2014), arXiv:1310.6061 [astro-ph.CO].\n- [126] A. Arbey, J. Lesgourgues, and P. Salati, Phys. Rev. D 68 , 023511 (2003), arXiv:astro-ph/0301533.\n- [127] A. Diez-Tejedor, A. X. Gonzalez-Morales, and S. Profumo, Phys. Rev. D 90 , 043517 (2014), arXiv:1404.1054 [astroph.GA].\n- [128] P. Sikivie and Q. Yang, Phys. Rev. Lett. 103 , 111301 (2009),", 'arXiv:0901.1106 [hep-ph].': '- [129] T. Takahashi, H. Omiya, and T. Tanaka, (2024), arXiv:2408.08349 [gr-qc].'}