{'fraction_non_alphanumeric': 0.046561391690277586, 'fraction_numerical': 0.01878306354143964, 'mean_word_length': 4.40327588052671, 'pattern_counts': {'":': 0, '<': 12, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Applications of the close-limit approximation: horizonless compact objects and scalar fields

The ability to model the evolution of compact binaries from the inspiral to coalescence is central to gravitational wave astronomy. Current waveform catalogues are built from vacuum binary black hole models, by evolving Einstein equations numerically and complementing them with knowledge from slow-motion expansions. Much less is known about the coalescence process in the presence of matter, or in theories other than General Relativity. Here, we explore the Close Limit Approximation as a powerful tool to understand the coalescence process in general setups. In particular, we study the head-on collision of two equal-mass, compact but horizonless objects. Our results show the appearance of ``echoes'' and indicate that a significant fraction of the merger energy goes into these late-time repetitions. We also apply the Close Limit Approximation to investigate the effect of colliding black holes on surrounding scalar fields. Notably, our results indicate that observables obtained through perturbation theory may be extended to a significant segment of the merger phase, where in principle only a numerical approach is appropriate.

The ability to model the evolution of compact binaries from the inspiral to coalescence is central to gravitational wave astronomy. Current waveform catalogues are built from vacuum binary black hole models, by evolving Einstein equations numerically and complementing them with knowledge from slow-motion expansions. Much less is known about the coalescence process in the presence of matter, or in theories other than General Relativity. Here, we explore the Close Limit Approximation as a powerful tool to understand the coalescence process in general setups. In particular, we study the head-on collision of two equal-mass, compact but horizonless objects. Our results show the appearance of "echoes" and indicate that a significant fraction of the merger energy goes into these late-time repetitions. We also apply the Close Limit Approximation to investigate the effect of colliding black holes on surrounding scalar fields. Notably, our results indicate that observables obtained through perturbation theory may be extended to a significant segment of the merger phase, where in principle only a numerical approach is appropriate.

I. INTRODUCTION
A few years ago we witnessed the birth of a new science, gravitational wave (GW) astronomy [1][2][3][4][5]. Observatories like LIGO and Virgo are now able to see routinely the inspiral and coalescence of compact objects such as black holes (BHs) and neutron stars. Observations are so far consistent with all the predictions made by General Relativity (GR), even in regimes probing strong-field and dynamical processes. In the near future, similar observations will be pushed to exquisite precision with thirdgeneration detectors and the space-based LISA mission. With high quality data and low instrumental noise, these observatories will have access to virtually all the visible universe [6,7].
In this new era of GW astronomy, new precision tests of GR and search for new physics also become possible [8][9][10]. Of special relevance is the scrutiny of GR predictions regarding the gravitational interaction in the strong-field regime. A fundamental result of vacuum GR is that isolated BHs all belong to the same family of solutionsthe Kerr family [11] -fully described by two parameters alone, mass and angular momentum [12][13][14]. Thus, testing the Kerr nature of BHs means, to some extent, testing GR. In addition, foundational questions are associated with the presence of horizons, particularly issues concerning the breakdown of determinism associated with Cauchy horizons or the fate of singularities of the classical equations [10,[15][16][17][18]. While possible pathological behavior is conjectured to be hidden behind horizons, questions remain concerning the effect of quantum gravity on the near-horizon structure or even on horizons themselves: do horizons exist? Are the objects we observe really BHs, or are they extreme (and exotic) compact objects (ECOs) which mimic the BH behaviour? Gravitational wave astronomy can have an important role in this matter, by constraining the existence of "echoes" or assessing the tidal properties of compact objects [10,[19][20][21][22][23][24][25][26][27][28][29][30]. Precision GW astronomy can also inform us on the nature and distribution of dark matter. A nontrivial dark matter environment changes the inspiral of a compact binary, via accretion or dynamical friction. If dark matter consists on new fundamental light fields, then rotating BHs can become lighthouses of GWs [31][32][33][34][35].
The possibility of extracting meaningful results and information from GW observations depends on our ability to model the coalescence of compact objects. These processes are now well understood: the inspiral part is modeled using slow-motion expansions [36,37], the merger through fully non-linear numerical simulations (the socalled "numerical relativity" (NR)) [38][39][40][41], and the ringdown of the final object is modeled using spacetime perturbation theory [42][43][44][45][46]. Perturbation theory, extended to include self-force effects [47,48], can also be used to model extreme mass ratio binaries.
These models are built in the framework of GR, and rely on the hypothesis that all compact objects in the Universe (with masses (2.5 − 3) M ) are BHs. To test GR and the BH hypothesis beyond null or biased tests [49], one needs to develop models of the inspiral, merger and ringdown which do not assume GR and do not assume that the compact objects are BHs. In particular, to test GR and the BH hypothesis in the strongfield regime we need to extend the models of the merger phase. However, the extension of NR to modified gravity theories and to ECOs is not an easy task, since they require a well-posed and well-behaved formulation of the time-evolution problem [50][51][52][53][54][55].
A possible way to circumvent the above problem is to use an alternative approach, such as the close-limit approximation (CLAP) [56][57][58][59][60][61][62][63][64]. In this approach, the slice of spacetime describing a late stage of the merger is computed solving the appropriate constraints of the field equations, as is commonly done to find the initial data for NR simulations. These initial data, if the two bodies are close enough, can be interpreted as a deformation of the stationary spacetime describing the final BH. Therefore, the subsequent evolution of this geometry can be studied using perturbation theory: the perturbation parameter, in this case, is the separation between the two BHs. Thus, we do not need NR to evolve these initial data: it is sufficient to solve the perturbation equations (e.g. the Zerilli equation) to obtain the GWs emission from the merger and ringdown stages. Additionally, a recent method called the "Backwards One-Body (BOB) method" [65] assumes as well the end state of a binary coalescence to compute waveforms emitted by coalescing binaries; in this context, the dynamics of light rays around the final perturbed BH is the key to compute the full GW signal.
Since the CLAP showed a remarkable agreement with NR simulations [66,67], its extension beyond GR and beyond the BH hypothesis could be a valuable tool to model the merger stages of compact binary coalescence in an extended framework.
This article is the first step towards the extension of the CLAP to modified theories of gravity and to extreme compact objects. In the first part of the paper we shall consider head-on collisions of ECOs, finding the gravitational waveform emitted in this process. In this article we shall consider the simplest version of the CLAP, in which the initial data are assumed to be static. This is of course a strong approximation, which will be relaxed in subsequent work. Thus, our results should be interpreted only as an indication, showing the qualitative features of the emitted GW signal. We also remark that, since we consider the head-on collision of non-spinning objects, we do not capture the angular momentum emission which would be produced in a more realistic setup.
In the second part we shall consider a test scalar field in a binary BH (BBH) spacetime (which can be interpreted either as ultra-light fundamental fields, or as extra gravitational degrees of freedom in scalar-tensor gravity), finding an indication on how the dynamical evolution of the background during the merger -which we model using the CLAP -affects the quasi-normal mode (QNM) frequencies of the scalar field.
These are just the first steps of this project: in future works we shall implement non-head-on and non-static initial data, to study the emission from the merger of two ECOs, and the emission of a BBH merger in modified theories of gravity such as Einstein-dilaton Gauss-Bonnet gravity or dynamical Chern-Simons gravity.
In Sec. II we review the standard CLAP approach to describe the merger of BBHs in GR. In Sec. III we use the CLAP to describe the head-on collision of two ECOs. In Sec. IV we use the CLAP to show the effect of a BBH merger on the QNMs of a test scalar field. In the following we use geometric units G = c = 1. Greek indices refer to quantities defined on a four-dimensional (4D) space-time manifold, Latin indices to their three-dimensional spatial parts (3D).

II. THE CLOSE LIMIT APPROXIMATION
We shall here review the standard CLAP approach for BBH coalescences in GR. This approach is based on the assumption that in the final stages of the coalescence, when the two BHs are sufficiently close to each other, the spacetime is a small deformation of a single BH spacetime, and thus can be studied using the techniques of spacetime perturbation theory. To this aim, BBH initial data, originally developed in the context of NR, are recast as a perturbation of Kerr spacetime. In this article we shall consider, as a first step, the case of a head-on collision of non-rotating BHs, which leads to a non-rotating BH. Thus, the BBH initial data are recast as a perturbation of a Schwarzschild BH.

A. The 3 + 1 decomposition
We briefly recall the basic concepts of the 3 + 1 decomposition in NR. We refer the interested reader to one of the several excellent books and reviews on the subject [68][69][70][71] for further information and details.
In order to formulate GR (or any other gravitational theory) as a time evolution problem, we first decompose the 4D spacetime in a set of 3D spatial hypersurfaces Σ t , labeled by a time parameter t, each of them having a 3-metric γ ij given by the space components of Here, g µν is the 4D spacetime metric, and n µ is the unit timelike (n µ n µ = −1) vector normal to the hypersurfaces Σ t . Thus, we can write the spacetime metric as where α, β i are called lapse function and shift vector, respectively, and contain the information about how the coordinate system changes from a slice to another. The choice of the lapse and shift, thus, corresponds to the choice of the foliation of the spacetime. The embedding of the hypersurfaces in the 4D spacetime is described by the extrinsic curvature where L n = α −1 (∂ t − L β ) is the Lie derivative along the unit vector n µ . If a scalar field Φ is also present, we also define its momentum Then, all 4D quantities such as the Ricci scalar R, the Ricci tensor R µν , etc. can be decomposed in terms of the 3D metric γ ij . The 3D Ricci tensor and Ricci Scalar, 3 R ij , 3 R, the extrinsic curvature and K ij , its trace K = γ ij K ij , and the scalar field momentum Π follow from this 3D metric. We denote with ∇ µ the covariant derivative with respect to the 4D metric g µν , and with D i the covariant derivative with respect to the 3D metric γ ij . With this decomposition, Einstein's equations (possibly supplemented with the scalar field dynamical equation) lead to two sets of equations: (i) the evolution equations, which are (with a careful choice of variables) hyperbolic equations giving the time evolution of the 3D quantities; (ii) the constraint equations, which are ellyptic equations which have to be satisfied at any 3D slice Σ t . The constraint equations are the Hamiltonian constraint equation and the three momentum constraint equations In presence of matter, the constraint equations are where ρ is the energy density of the matter fields, and j i is its energy-momentum flux. For a scalar field they are where Π is the momentum of the scalar as defined in Eq. (4). In NR one first solves the constraint equations at the initial time t = t 0 , finding the initial data on Σ t0 of the system. Then, one solves the evolution equations, finding the spacetime at all times t ≥ t 0 . Within the CLAP, instead, the evolution of physical initial data is determined through the use of perturbation theory, therefore solving a linearized version of the evolution equations (see below).

B. Initial data for binary black holes in general relativity
Solving the constraint equations (5)-(6) and finding initial data appropriate to study a given system is a subject of study on its own [69][70][71][72]. The original formulation of the CLAP [56] used the Misner initial data [73] (a common choice for NR simulations decades ago) to describe BBH head-on collisions. Subsequently, it was found that the Brill-Lindquist (BL) initial data [74] are more appropriate for NR simulations, because they have a simpler form and are easier to be extended to the more realistic non-head-on case (the Bowen-York initial data [75]). Thus, later applications of the CLAP employ the BL (or Bowen-York) initial data [57,59,63,76]. In this article, we shall use BL initial data and their extensions.
In the case of a head-on collision of two non-rotating compact objects starting from rest, the extrinsic curvature identically vanishes. Additionaly, considering that in vacuum Φ = 0, the constraint equations (5)-(6) reduce to 3 R = 0. The BL three-metric γ ij is a conformally flat solution of this equation describing two BHs. It has the form where dΩ 2 = dθ 2 + dφ 2 sin 2 θ , 3 dη 2 is the flat threemetric, and ϕ BL (R, θ, φ) is the conformal factor. With this choice, the Hamiltonian constraint 3 R = 0 reduces to ∇ 2 ϕ BL = 0. The BL solution of this equation is where the vector R is the position vector in the flat threespace, and R i is the position of the i-th BH; M = m 1 + m 2 is the ADM mass of the entire spacetime. Note that m 1 and m 2 are not the ADM masses of the two BHs, which are [74] and d = |R 2 −R 1 |. We choose the origin of the coordinate system in the center of mass: In the region of spacetime where R > R 1 , R 2 , the conformal factor can be expanded in terms of Legendre polynomials: The coefficients ξ are given by [76] In the following we shall consider equal-mass binary systems; thus m 1 = m 2 = M/2. By choosing the Zaxis aligned with the motion, When the coefficients ξ vanish, ϕ BL = 1 + M/(2R) and the metric (11) describes the t = const. slices of the Schwarzschild's background. Indeed, by defining the Schwarzschild radial coordinate r in terms of the isotropic radial coordinate R by Eq. (11) gives where If ϕ BL = 1 + M/(2R), F BL = 1 and by defining the 4D spacetime metric as (corresponding to an appropriate choice of the shift vector and lapse function, i.e. to an appropriate gauge choice), it coincides with the Schwarzschild geometry. If, instead, the coefficients ξ are small but non-vanishing, the metric (11) is a perturbation of Schwarzschild's threemetric, and Eq. (21) is a perturbation of Schwarzschild's spacetime. Therefore the BL three-metric, in the coordinates (r, θ, φ), has the form If ξ 1 we can linearize in the parameters ξ , and Eq. (22) gives: The parameter Z 0 in Eq. (16) describes the (initial) distance between the BHs in the isotropic frame. This quantity does not have a direct physical interpretation; thus, several authors characterize the initial BH separation with the proper distance L between the apparent horizons of the two BHs, given by [76] The extrema Z 1 , Z 2 are the intersections of the apparent horizons of the two BHs with the Z-axis, and can be found by numerical integration of the equations describing the apparent horizons [77,78]; the procedure to compute these quantities is described in detail, for instance, in [63]. This computation shows that, for example, L = 3M for Z 0 0.5M , L = 3.5M for Z 0 0.7M , L = 4M for Z 0 0.85M .
A comparison with NR computations in the head-on case shows that the CLAP is accurate, in the equal-mass case, for L 4M (see [61] and the discussion in [63]), corresponding to Z 0 0.85M and thus ξ 2 0.36; as noted in [56], it is remarkable that the agreement extends far beyond the region ξ 1 in which the CLAP is expected to be applicable. Therefore, we shall apply Eq. (23) and, more generally, perturbation theory, also to initial separations for which the condition ξ 1 is only marginally satisfied.

C. Perturbations and their time evolution
The metric in Eq. (21) describes the Schwarzshild spacetime, which we consider as the background, with a perturbation with even parity, at a fixed "initial" time t = t 0 . Thus, it can be recast in the form where is the Schwarzschild spacetime (with f defined in Eq. (19)), and h µν (t, r, θ, φ) is the perturbation, whose evolution in the wave zone r M can be described in terms of a scalar function called the Zerilli function [43]. The original definition of the Zerilli function is appropriate for the study of oscillating solutions; in this case, we want to study the evolution of a given set of initial data at the time t = t 0 ; thus, it is more appropriate the definition of the Zerilli function in terms of gauge-invariant quantities given in [79,80].
We first note that the only non-vanishing perturbations in the metric (21)- (23) are where Using the notation of [80] for polar parity, axially symmetric perturbations, Eqs. (27) correspond to the perturbation functions H 2 (t, r), K (t, r) given by which implies that at the initial time t = t 0 , H 2 (t 0 , r) = K (t 0 , r) = g ξ , with = 2, 4, 6, . . . , (30) while the other perturbation functions (H 0 , H 1 ) and (G , h 0 , h 1 ) identically vanish at t = t 0 . Since the leading contribution comes from the quadrupole perturbations, we shall consider the = 2 contribution only. Then, following [80], we define the = 2 Zerilli function at t = t 0 , where λ = 1 + 3M/(2r) and (leaving implicit the index = 2 and the arguments of H 2 (t 0 , r) and K(t 0 , r)) The perturbation equations can be written as a wave equation for the Zerilli function, where r * = r+2M log r 2M − 1 is the tortoise coordinate, and is the Zerilli potential. Eq. (33) is called the Zerilli equation [43].
In most implementations of the CLAP, the Zerilli equation (33) with initial condition (31) is solved in the time domain. Alternatively, it can be solved in the frequency domain. Indeed, it can be shown (see e.g. [81]) that if the function ψ(t, r * ) is the solution of the Zerilli equation (33) with initial conditions (31) andψ(t 0 , r * ) = 0 (i.e., assuming stationarity of the initial data), then, choosing the time coordinate such that t 0 = 0, the Laplace transform is the solution of the Zerilli equation with source where and boundary conditions of ingoing wave at the horizon, outgoing wave at infinity (see Eq. (A4) in Appendix A 1). Finally, the total energy emitted in GWs in the collision is given by [56] We have solved Eq. (36) using two different approaches (see Appendix A 1 for details), reproducing in both cases the results in the literature for head-on BBH collisions with the CLAP approach [56]. The total energy obtained with our computation agrees with that of Ref. [56] within one percent.
In this article we only consider head-on collisions of non-rotating, equal mass compact objects starting from rest as in [56]; however, it is worth mentioning that the CLAP approach has been extended to more realistic setups, considering BHs with initial velocity [82], unequal masses [76] and non-head-on binary inspirals [59].

III. EXTREME COMPACT OBJECTS
We shall now apply the CLAP to describe the headon collision of two non-rotating, equal mass ECOs. We model an ECO (see e.g. [10] and references therein) as a spherically symmetric compact body with mass M and a surface at r = r 0 , with 1. The above parametrization for the ECO's surface well describes ECO models for which 1, i.e. gravastars, wormholes, fuzzballs etc. [10]. In the following we restrict to ultra-compact bodies whose surface can be described by Eq. (39). This means that the dynamical processes involving such ECOs resemble those of BHs (with due differences, as we will see shortly) rather than stellar objects such as neutron stars for which ∼ O (1)).
A fundamental difference between the surface of an ECO and the horizon of a BH is that an incoming wave is partially reflected by the ECO surface, while it is totally absorbed by a BH horizon. Thus, given a (scalar, gravitational, etc.) test field φ(t, r * ) ∼φ(r * )e iωt , its boundary condition near the surface r = r 0 is The parameter R, which in general depends on the frequency ω and on the spin of the field, is called reflectivity coefficient of the ECO.

A. Initial data for extreme compact objects
We shall first introduce the spacetime metric of a single, isolated, spherically symmetric ECO; then, using the procedure of Brill and Lindquist discussed in Sec. II B, we shall define initial data corresponding to two ECOs starting from rest (this approach can naturally be extended to the case of ECOs with initial momentum, using the procedure of Bowen and York [75]).
The GW signal produced by the collision of ECOs will consist of waves caused by excitations of the exterior spacetime [46,83], as well as contributions which probe the interior of the object [18][19][20][21][84][85][86]. For a wide class of asymptotically flat ECOs, the lapse in its interior is very small, leading to large delays in signals which have probed the interior geometry [10,18,87]. Thus, in the following we "freeze" the inner region, which means that in practice we cut it off from our domain, and all the information about the interior is replaced by boundary conditions at the surface.
Let us first consider a single, isolated (spherically symmetric) ECO. Due to Birkhoff's theorem, the exterior of the ECO r > r 0 is described by Schwarzschild's geometry (Eq. (26)) 1 . We define the isotropic radial coordinate as in Schwarzschild metric (17): In terms of the coordinate R, the location of the surface ). The interior of the ECO depends on the profile of its energy density ρ(R). As discussed above, we do not need to model the interior because, under our assumptions, the emitted GW signal does not depend on the structure of the interior; it only depends on the location of the surface and on its properties. We stress that "freezing" the internal degrees of freedom of the ECOs is a standard approximation in the literature of ECO modeling [10].
Let us now consider two spherically symmetric ECOs, with equal masses and total ADM mass M , initially at rest. We perform a 3 + 1 decomposition of the spacetime, as in Sec. II A. The three-metric of each slice is solution of the Hamiltonian constraint equation (7): Since the ECOs are initially at rest, the initial extrinsic curvature vanishes and Eq. (42) reduces to 3 R = 16πρ.
In the exterior of both ECOs, ρ = 0 and the metric can be written in the BL form, Eq. (23): with ξ = (Z 0 /M ) /2 as in Eq. (16). In the isotropic frame, the two ECOs move along the z-axis and are located at R = ±Z 0 , with 0 < Z 0 0.85M (see also Sec. II B). We call Eq. (43) the BL-ECO initial data.
We remark that although the BL initial data for BBHs in Eq. (23) and the BL-ECO initial data for binary ECOs in Eq. (43) are formally identical, the former are defined outside the BH horizon, while the latter are defined outside the ECO surface. This difference, as we shall show below, leads to a difference in the GW emission. Note also that, as discussed above, the region inside the surfaces of the ECOs does not contribute to the GW signal, by assumption. 1 Since 1, r 0 < 3M and thus the light ring r = 3M , where null circular geodesics are defined, lies in the exterior of the ECO. This plays a crucial role in the GW emission, as we discuss below.

B. Collision of extreme compact objects
We shall now recast the BL-ECO initial data in Eq. (43) as a perturbation of a single object, Then, we shall evolve them using the tools of perturbation theory.
Since we are not dealing with vacuum spacetimes, the outcome of the ECO collision depends, in principle, on the features of the ECOs, i.e. on their internal structure. Two outcomes are possible: either the collision leads to a single BH (which, in the head-on collision of non-rotating ECOs, is described by the Schwarzschild metric), or it leads to an ECO. In the former case the initial data should be recast as a perturbation of a Schwarzschild solution; in the latter, as a perturbation of an ECO. We shall treat these two cases separately.

Formation of a black hole
If the final object is a Schwarzschild BH, the background metric g (0) µν in Eq. (44) is Schwarzschild's metric. The procedure of recasting the BL-ECO initial data (43) as a perturbation of Schwarzschild metric is formally equivalent to the derivation in Sec. II C, leading to an = 2 Zerilli function at the initial time t = t 0 given by Eqs. (31)- (32).
Since the background is Schwarzschild's metric, the Zerilli equation (33) has to be solved with the same boundary condition as in the BBH case. Therefore, the gravitational waveform emitted in the collision of two ECOs into a Schwarzschild BH, computed with the CLAP (and neglecting the effect of perturbations inside the ECO's surface) coincides with the waveform emitted in a BBH collision. 2

Formation of an extreme compact object
If the final object is an ECO, the background metric g (0) µν is Schwarzchild metric only for r > r 0 .
It seems extraordinary that two compact objects coalesce and are prevented to collapse completely to a BH. However, we intend on mimicking the effects of a theory where BHs are absent (via, for example, high-energy phenomena). We may assume, for example, that a theory of quantum gravity forbids the existence of horizons via Planck-scale physics. We do not wish to build such a theory but merely to investigate some of its consequences. We do assume that departures from GR occur only close to the Schwarzschild radii.
We remark that, as discussed in the Introduction, this computation is based on strong assumptions (the bodies collide head-on, they have no spin, and most importantly the "initial data" are stationary, although they describe the late stage just before the merger). Therefore, our results should only be considered as a order-of-magnitude estimates. The computation itself should be considered as the first step in the modeling of coalescences with ECOs using the CLAP. These approximations will be relaxed in future publications.
The procedure of recasting the BL-ECO initial data (43) as a perturbation of the ECO metric follows the derivation in Sec. II C, leading to the Zerilli function at the initial time t = t 0 given in Eqs. (31) and (32). However, the Zerilli function ψ(t, r * ) is only defined for r * > r 0 * = r 0 + 2M log |r 0 /(2M ) − 1|, and it satisfies boundary conditions different from those of a BH.
Using the Laplace transform approach discussed in Sec. II C, one finds the Zerilli equation with source given in Eqs. (36) and (37). Although the equation is the same as for BBH collisions, the boundary conditions are different: since the surface is partially reflecting, they are (see Eq. (40)) where R is the reflectivity coefficient of the ECO 3 . We solved the Zerilli equation (36) with boundary conditions (45) using a shooting method (see Appendix A 2), for different values of and of R, finding the gravitational waveform emitted in the merger and ringdown phases of the collision. Since we only consider the leading quadrupolar contribution, the Zerilli function is proportional to ξ 2 , Eq. (31), and thus to Z 2 0 . Moreover, it has (in geometric units) the dimensions of length and thus, since the only dimensionful quantity characterizing the ECO is its ADM mass M , ψ(t, r * ) ∝ Z 2 0 /M . Figure 1 shows the Zerilli quadrupolar waveform generated by the head-on collision of two spherically symmetric compact objects starting from rest, evaluated at the extraction radius r extr * = 400 M . The waveforms are rescaled by 16Z 2 0 /M , and thus they do not depend neither on Z 0 nor on M . We compare processes leading to FIG. 1. Gravitational waveform for the head-on collision of two equal-mass, spherically symmetric BHs/ECOs starting from rest, evaluated at the extraction radius r extr * = 400 M , for two possible outcomes. If the final object is a BH, the waveform is identical to that resulting from the collision of two BHs: a sharp burst followed by ringdown (caused by the relaxation of the light ring). When the final product is an ECO, we observe a similar initial stage, followed at late times by echoes of the initial burst [10,20]. Here, the final object is taken to have a reflectivity of R = 0.1 and a surface r0 = 2M (1 + ) with = 10 −10 .

Results of the numerical integration
a Schwarzschild BH (dashed curve) and to an ECO with = 10 −10 and R = 0.1 (solid line). The signals are identical at early times while they differ at late times, since the ECO collision signal is characterized by a series of "echoes", which are a characteristic feature of GW signals from ECOs (see e.g. [10, 18-20, 84, 88]). We stress again that, in this approximation, the signal depends on the nature of the final object, but it does not depend on whether the colliding objects are BHs or ECOs.
A simple interpretation of the echo structure is the following. The light ring excitation by the initial data is followed by its "relaxation". Since both BHs and ECOs have equivalent geometries close to the light ring, they both relax in the same way [10,[18][19][20]. This relaxation produces an outgoing wave, which corresponds to the early part and ringdown of the signal. But there's also an "ingoing" wave which interacts with the ECO via boundary conditions. Due to a nonzero reflectivity, the ingoing pulse is partly reflected back, and interacts with the light ring again. The process continues and produces a sequence of distorted copies of the original burst, i.e. the echoes. Figure 2 displays the dependence of the Zerilli waveform on the ECO parameters and R. We can see that the echo delay increases with , while the amplitude of the echoes is larger for larger values of the reflectivity R.
We find that the time separation between echoes is well described (with a relative error smaller than 2%) by the following analytical fit, for 10 −10 < < 10 −6 . This logarithmic behaviour is consistent with our knowledge of echoes in the GW signals [10]. We only solved the Zerilli equation for ≥ 10 −10 ; for smaller values of , our numerical approach loses accuracy; however, since this limitation is of computational nature only, we expect the fit (46) to hold for smaller values of as well.
Finally, we compared the total energy E emitted in GWs during a compact object collision leading to an ECO, with that emitted in a BBH collision. It can be computed in terms of the Zerilli function from Eq. (38). We find that, for < 10 −3 , the energy radiated is very weakly dependent on , while it increases with the reflectivity R. This is due to the fact that the presence of echoes contributes to the energy loss of the system in GWs.
In Fig. 3 we show the total energy E emitted in GWs as a function of R, normalized by the energy emitted for R = 0 (which coincides with the result for a BBH collision [56]). This shows which part of the emitted energy is due to echoes. The function E(R) can be described by the following fit (accurate within 1.5% in the range FIG. 3. Total energy emitted in the head-on collision from rest of two equal-mass, spherically symmetric compact objects. Fixing the ECOs/BHs initial separation Z0, the ratio of energy emitted is independent of such quantity. The final object is assumed to be an ECO with = 10 −10 and the results are shown as a function of the reflectivity R. Inset: log-log scale of the same results, including R = 1. This also corresponds to the formation of a wormhole [19] (see main text). One can use the CLAP formalism to understand, in particular, the radiation given away during the formation of wormholes [89,90]. A class of these objects can be considered as ECOs with R = 1, when the "throat" is made of a rigid shell [19]. This corresponds to taking R = 1 in the above framework. The inset of Fig. 3 shows the total energy computed with the CLAP for head-on collisions forming one of these wormholes. The total radiated energy is slightly off the predictions of the fit (47) because points R 0.8 were not used in the fit. Notice that the total radiated energy is over one order of magnitude larger than the corresponding process forming a BH: substantial amount of energy is released in late-time echoes.
Furthermore, a spectral analysis shows that the ECO formation excites certain characteristic frequencies, which correspond to the quasinormal frequencies of the final object. In the particular case of a thin-shell wormhole described above, most of the radiation is in fact contained in such modes: the energy spectrum shows clear peaks at such quasinormal frequencies, which are in excellent agreement with the spectral findings of Ref. [19].

IV. BBHs AND SCALAR FIELDS
From the first observation of the Higgs boson by the ATLAS Collaboration [91], a growing recognition has been given in studying the effects of scalar particles, both at a cosmological and an astrophysical level [92][93][94][95][96][97][98][99]. Moreover, most of the modifications of GR which have been proposed so far can be reformulated in terms of couplings between gravity and extra fields, the simplest of which are scalar fields [8]. In this Section we study how the presence of scalar fields may affect a BBH collision, and whether they can leave observable imprints during the GWs generation.
We shall consider gravity minimally coupled with a complex scalar field. Since we are interested in BH solutions, we do not include matter fields in the model. Thus, the action is: The field equations obtained from this action are Einstein's equations coupled with the Klein-Gordon equations: where is the scalar field stress-energy tensor. A wide class of modified gravity theories in which gravity is nonminimally coupled with a scalar field -the so-called Bergmann-Wagoner scalar-tensor theories (see e.g. [8,100] and references therein), is formally equivalent to the theory in Eq. (48) 4 , if restricted to vacuum spacetimes. Thus, the scalar field Φ can be interpreted either as a fundamental "matter" field in GR, or as a gravitational degree of freedom in a modified gravity theory. The 3 + 1 decomposition of the action (48) has been discussed in Sec. II A. In particular, the constraint equations have the form where the energy density ρ and the energy-momentum flux j i of the scalar field are given in Eqs. (9)-(10). The theories described by the action (48) satisfy the no-hair theorem: stationary BH solutions are described by the Kerr metric, and thus they have vanishing scalar field [101][102][103] (see also [14] and references therein). Therefore, we know that the remnant of a BBH collision becomes -in the timescale of the QNM oscillations, i.e. of ∼ 1 − 10 M -a stationary BH solution, with vanishing scalar field.
The no-hair theorem does not tell us what happens before reaching the final stationary configuration: it does not constrain the dynamics of BH spacetimes. Still, a similar result applies to the inspiral part of a BBH coalescence: an analysis in the post-Newtonian (PN) approximation [36,37], which accurately describes the BBH inspiral, shows that the binary dynamics in the theory (48), up to 2.5 order in the PN expansion, is the same as in GR [104] 5 . However, we do not know whether the presence of a scalar field significantly affects the BBH dynamics during the merger and ringdown stages.
In order to address this problem, we shall study the QNMs of the scalar field during the merger and ringdown of a BBH (head-on) collision. We do not expect the scalar field to grow large before the collision, thus we shall treat it as a perturbation of the BBH spacetime; we define a perturbation parameter 1, such that Φ = O( ). Therefore, T µν = O( 2 ) and, to linear order in the perturbation, we can neglect the scalar field from Einstein's equations (49) (and in particular from the constraint equations (51)-(52)). We shall then study the linearized (O( )) field equation of the scalar field in the BBH spacetime, which is modeled using the CLAP approximation. For such configurations we shall compute the scalar field QNMs. Then, by comparing the scalar field QNMs in a BBH spacetime with the scalar field QNMs in a stationary BH spacetime, we will assess whether the binary dynamics significantly affects the scalar field dynamics.
We describe the BBH spacetime (neglecting the scalar field, as discussed above) using BL initial data. Therefore, we recast the BBH spacetime as a perturbation of the Schwarzschild metric. Including the leading-order quadrupolar contribution, as discussed in Sec. II C, the perturbed spacetime can be written as g µν = g (0) µν + h µν where g (0) µν = diag(−f, f −1 , r 2 , r 2 sin 2 θ) and h µν is given by Eqs. (29)- (30): where g = 4 (1 + M/(2R)) −1 M 3 /R 3 , the isotropic coordinate R is defined in Eq. (17), and ξ 2 = Z 2 0 /(2M 2 ) 1. We stress that we are assuming the background to be stationary, thus neglecting the motion of the BHs, in the timescale of the oscillation (t ∼ (1 − 10)M ). This is a crude approximation, since the BH separation changes, and their velocities become non-negligible in this timescale. Thus, the result of this computation should be considered as an order-of-magnitude estimate of the effect of the BH dynamics on the scalar QNMs.
Since h µν ∝ Z 2 0 , we can expand the D'Alembertian operator = 1 √ −g ∂ µ (g µν √ −g∂ ν ) in Eq. (50), for small separations, as where the explicit form of the operator (1) is given in Appendix B. Expanding the scalar field as we get a non-separable equation, where the harmonic component ψ m couples with the components ψ ±2 m . As discussed in Appendix B, we can follow the same approach used to study scalar field perturbations around rotating BHs (see e.g. [105,106]); indeed, the leadingorder rotational corrections are quadrupolar as well. Remarkably, the ↔ ± 2 couplings do not affect the QNM frequencies at leading order in the perturbations (see also [107,108]), and thus they can be neglected, leading to a decoupled, Schrödinger-like equation: where ≥ 1 6 . The derivation of Eq. (56) and the explicit form of the functions U A (r),Ũ A (r), W A (r) (A = 0, 1) are given in Appendix B.
In order to find the QNMs, we have solved Eq. (56) through direct integration with outgoing boundary conditions at infinity and with ingoing boundary conditions at the horizon. The boundary conditions have been determined as a perturbative, polynomial expansion at each boundary, whose coefficients have been found solving Eq. (56) order by order, as explained, e.g. in Ref. [109].
To validate our integration method, as a first step we computed the fundamental scalar QNMs of a Scwharzschild BH, setting Z 0 = 0. Our results agree with those in the literature (e.g. [110]) within 0.2%. Then, we computed the scalar QNMs in BBH collisions for different values of ≥ 1 and of Z 0 . Figure 4 shows the fractional percentage shift of real and imaginary parts of the QNMs with respect to those in Schwarzschild BHs: As we can see from Fig. 4, the QNMs shifts have a nontrivial dependence on the BH separation Z 0 , but each mode is shifted by just a few percent from the corresponding mode of an isolated BH. Despite the approximations discussed above, this result provides a strong indication that the BH dynamics does not significantly affect the behaviour of the scalar perturbations, at linear order. In particular, our results suggest that the scalar modes do not become unstable during the BH coalescence, and thus a small scalar field would remain small during the merger. For results concerning non-trivial effects of a self-interacting ambient scalar field on a BBH, we refer the reader to Refs. [111,112].
Thus, our results indicate that scalar fields in GR, or in modified gravity theories in which the no-hair theorem applies, do not significantly affect the coalescence of BBHs.
These results are also interesting in light of a completely different question: at what time the GW signal from a BBH coalescence can be described as a superposition of QNMs? This problem, i.e. the determination of the starting time of the ringdown, is widely debated [113][114][115][116][117][118]. Indeed, the procedure of constructing a GW template is based on joining different approximations, from the regime where the BHs separation is large (PN approximation) to the ringdown oscillations (perturbation theory), passing through the highly non-linear merger process that requires numerical approaches. Our results indicate that an observer measuring the scalar oscillations during a BBH collision may extend the validity of the ringdown treatment closer to the merger, where, in principle, only full numerical studies are accurate and reliable. The direct detection of scalar fields oscillations is not possible with current detectors, if one assumes GR as the correct theory of gravitation. However, we may take this as an indication that the same may hold true for gravitational waveforms, for which, instead, detections are possible through current GW interferometers.

V. CONCLUSIONS AND OUTLOOK
GW astronomy has the potential to answer crucial questions regarding the correct description of gravity. The full exploitation of such potential requires knowledge about the dynamics of compact objects in a generic theory of gravity. While NR is the tool of excellence for this, the evolution of a single binary within the context of a modified theory can take months to perform on supercomputers, and may require years of careful study of the relevant partial differential equations and associated well-posedness.
In this work, we explored the close-limit approximation as a "quick-and-dirty" tool to understand nonlinear coalescence processes. Its remarkable agreement with full nonlinear simulations is an important benchmark. In fact, albeit it is a perturbative scheme, it uses constraintsatisfying initial data, and their evolution works accurately even when the premises of the model are only partially satisfied. Consequently, this provides some confidence that this technique works well also when extending those studies beyond GR or BH spacetimes. The main requirement to use the CLAP consists in having solutions of the constraint equations. These can be solved, as we showed, also in the presence of fundamental fields (see also Ref. [119,120]).
This work is just a first step towards the implementation of the CLAP beyond the standard GR/BH scenario. Our computation is based on strong approximations; thus our results should be interpreted as qualitative indications. Further effort is required here. We expect the payoff to be significant: with much less computational time and effort one is able to investigate setups that describe better nonlinear geometries. We showed how the CLAP can work for the coalescence of compact, horizonless objects, and how it too predicts the existence of echoes in gravitational waves. This is a significant result in that it extends and complements other past perturbative calculations [10,19,20]. Moreover, we studied scalar fields minimally coupled with gravity in BH spacetimes (which are equivalent to non-minimally coupled scalar fields in a Bergmann-Wagoner scalar-tensor theory), estimating the scalar modes in the merger of a BH binary, and showing that they are very similar to those in a stationary BH spacetime.
In future works, we will relax some of the approximation used in this work, by considering non-head-on collisions of spinning objects, and non-stationary initial data. Moreover, we shall use the CLAP to investigate compact objects collisions in other -and perhaps more complicated, modified gravity theories such as Einsteinscalar-Gauss-Bonnet gravity. Acknowledgements.
L. A. acknowledges financial support provided by Fundaçao para a Ciência e a Tecnologia Grant number PD/BD/128232/2016 awarded in the framework of the Doctoral Programme IDPASC-Portugal. V. C. acknowledges financial support provided under the European Union's H2020 ERC Consolidator Grant "Matter and strong-field gravity: New frontiers in Einstein's theory" grant agreement no.
MaGRaTh-646597. This project has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No 101007855. We thank FCT for financial support through Project No. UIDB/00099/2020. We acknowledge financial support provided by FCT/Portugal through grants PTDC/MAT-APL/30043/2017 and PTDC/FIS-AST/7002/2020. We also acknowledge support from the Amaldi Research Center funded by the MIUR program "Dipartimento di Eccellenza" (CUP: B81I18001170001), and from the MIUR grant PRIN2017-MB8AEZ. The authors would like to acknowledge networking support by the GWverse COST Action CA16104, "Black holes, gravitational waves and fundamental physics." We shall discuss the numerical integration of the Zerilli equation with source in the frequency domain, Eq. (36) (see also [81,121]): whereψ(ω, r) is the Laplace transform (35) of the Zerilli function ψ(t, r), and S(ω, r) = iωψ(t = 0, r) = iω 4π 5 where Q(r) is given in Eq. (32) and ξ 2 = Z 2 0 /(2M 2 ) 1. We shall find the solution of Eq. (A1) satisfying ingoing boundary conditions at the horizon, outgoing boundary conditions at infinity, by employing two different approaches: the Green function approach and a shooting method, finding the same results.

Collision to a Schwarzschild BH
If the Zerilli equation (A1) describes perturbations of a Schwarzschild BH, it is defined in −∞ < r * < +∞. In this case, the source term (A3) does not vanish at the horizon, S(r * → −∞) = iωψ(t = 0, r * → ∞) =S = 0. The ingoing wave boundary conditions at the horizon r * → −∞, ∂ψ/∂r * = ∂ψ/∂t 7 translate in the Laplace transform space into ∂ψ/∂r * (ω) = −iωψ(ω) − ψ(t = 0) . Therefore, the boundary conditions ofψ(ω, r * ) are: with A H , A ∞ constants to be determined. The constant term is related to the fact that the BL initial data do not vanish at the horizon, and do not affect the GW emission at infinity. The Green function approach consists in finding two independent solutions of the homogeneous Zerilli equations:ψ H , satisfying ingoing wave conditions at the horizon, andψ ∞ , satisfying outgoing wave conditions at infinity, i.e. where W =ψ H ∂ψ ∞ /∂r * −ψ ∞ ∂ψ H /∂r * is the constant Wronskian of the homogeneous equation. By imposing the boundary conditions at infinity (A4) we find Apparently, the integral (A7) is not well defined at the lower bound, where the integrand reduces to the oscillating termSe −iωr * . This is due to the fact that, strictly speaking, the Laplace transform (A2) is well defined in the upper complex plane; the inverse Laplace transform can be computed along a path ω = ω R + i with 1 and −∞ < ω R < +∞. Thus, along this path e −iωr * → 0 as r * → −∞ and the oscillating term disappears. As 7 Note that these conditions are consistent withψ(t = 0) = 0 beacuse in the BL initial data ∂ψ/∂r * = 0 at the horizon.
suggested in [81], we can compute the integrals for real values of ω, as long as we subtract the ill-valued contribution at the horizon: wherer * is negative and very large. We have computed A ∞ (ω), by evaluating the integrals fromr * = −44 M to the extraction radius r extr * = 400 M . The time-domain Zerilli function at infinity can then be computed, as a function of the retarded time u = t − r * , as The shooting approach, instead, consists in the numerical integration of Eq. (A1), fromr * = −44 M to the extraction radius r extr * = 400 M , by imposing the boundary conditions (A4) and matching the solution at r extr * with an analytic expression obtained by an asymptotic expansion of Eq. (A1). For each value of the frequency ω, we performed the numerical integration of Eq. (A1) for different values of A H until we obtained an outgoing wave at infinity, as in Eq. (A4). In this way we computed the function A ∞ (ω) and then, by Eq. (A9), the Zerilli function at infinity.
The results of the two approaches perfectly agree with each other, and they also agree with the results of [56].

Collision to an Extreme Compact Object
If the outcome of the collision is an ECO, the Zerilli equation (A1) describes perturbations of the ECO spacetime which, as discussed in Sec. III, coincides with Schwarzschild's spacetime with the domain restricted to r 0 * < r * < +∞. Moreover, we impose at r * → r 0 * the partially reflecting boundary conditions (45): ψ = A H (e −iω(r * −r0 * ) + Re iω(r * −r0 * ) ) (r * → r 0 * ) +Ŝ ω 2 (A10) ψ = A ∞ e iωr * (r * → ∞) with R reflectivity coefficient which we assume, for simplicity, to be constant, and the constant term in the ingoing boundary conditions is due to the fact that S(r 0 * ) = iωψ(t = 0, r 0 * ) =Ŝ = 0. Note that since the tortoise coordinate does not extend to −∞, the Laplace transformψ(ω, r * ) is well-defined for real frequencies, and we do not need to worry about ill-defined contributions to the integrals. In this case we only perform the integration using the shooting approach, i.e. we integrate Eq. (A1), from r 0 * to the extraction radius r extr * = 400 M , by imposing the boundary conditions (A10), for different values of A H , until we obtained an outgoing wave at infinity. Then, using Eq. (A9) we compute the Zerilli function at infinity. therefore the = 0 equation is not affected by the O(Z 2 0 ) corrections. For 0 < ≤ 2, the quantities q (1,2)