[295, 254, 522, 298, 324, 606, 619, 194, 292]) to express the dynamics of GR as deviations from the Newtonian limit in terms of an expansion parameter . This parameter can be identified, for instance, with the typical velocities of the matter composing the source, or with the compactness of the source: [109*], but as a book-keeping parameter it is customary to consider . The spacetime metric and the stress-energy tensor are expanded in powers of and terms of order are commonly referred to as -PN corrections. The spacetime metric and the motion of the source are found by solving, order by order, Einstein’s equations.
Strictly speaking, the PN expansion can only be defined in the near zone, which is the region surrounding the source, with dimensions much smaller than the wavelength of the emitted GWs. Outside this region, and in particular in the wave zone (e.g., at a distance from the source), radiative processes make the PN expansion ill-defined, and different approaches have to be employed, such as the post-Minkowskian expansion, which assumes weak fields but not slow motion. In the post-Minkowskian expansion the gravitational field, described by the quantities (in harmonic coordinates, such that ) is formally expanded in powers of Newton’s constant . Using a variety of different tools (PN expansion in the near zone, post-Minkowskian expansion in the wave zone, multipolar expansions, regularization of point-like sources, etc.), it is possible to solve Einstein’s equations, and to determine both the motion of the source and its GW emission. Since each term of the post-Minkowskian expansion can itself be PN-expanded, the final output of this computation has the form of a PN expansion; therefore, these methods are commonly referred to as PN approximation schemes.
PN schemes are generally used to study the motion of -body systems in GR, and to compute the GW signal emitted by these systems. More specifically, most of the results obtained so far with PN schemes refer to the relativistic two-body problem, which can be applied to study compact binary systems formed by BHs and/or NSs (see Section 3.1.1). In the following we shall provide a brief summary of PN schemes, their main features and results as applied to the study of compact binary systems. For a more detailed description, we refer the reader to one of the many reviews that have been written on the subject; see e.g. [109*, 620, 676*, 454].
Two different but equivalent approaches have been developed to solve the relativistic two-body problem, finding the equations of motion of the source and the emitted gravitational waveform: the multipolar post-Minkowskian approach of Blanchet, Damour and Iyer [109*], and the direct integration of the relaxed Einstein’s equations, developed by Will and Wiseman [785*]. In these approaches, Einstein’s equations are solved iteratively in the near zone, employing a PN expansion, and in the wave zone, through a post-Minkowskian expansion. In both cases, multipolar expansions are performed. The two solutions, in the near and in the wave zone, are then matched. These approaches yield the equations of motion of the bodies, i.e., their accelerations as functions of their positions and velocities, and allow the energy balance equation of the system to be written as (see also [553*]), where is the (traceless) quadrupole moment of the source. The leading term in is then of 2.5-PN order (i.e., ), but since is computed in the Newtonian limit, it is often considered as a “Newtonian” term. A remarkable result of the multipolar post-Minkowskian approach and of the direct integration of relaxed Einstein’s equations, is that once the equations are solved at -th PN order both in the near zone and in the wave zone, is known at -PN order, and is known at -PN order with respect to its leading term, i.e., at -PN order. Once the energy and the GW flux are known with this accuracy, the gravitational waveform can be determined, in terms of them, at -PN order.
Presently, PN schemes determine the motion of a compact binary, and the emitted gravitational waveform, up to 3.5-PN order for non-spinning binaries in circular orbits [109*], but up to lower PN-orders for eccentric orbits and for spinning binaries [48*, 148*]. It is estimated that Advanced LIGO/Virgo data analysis requires 3.5-PN templates , and therefore some effort still has to go into the modeling of eccentric orbits and spinning binaries. It should also be remarked that the state-of-the-art PN waveforms have been compared with those obtained with NR simulations, showing a remarkable agreement in the inspiral phase (i.e., up to the late inspiral stage) [122*, 389].
An alternative to the schemes discussed above is the ADM-Hamiltonian approach [676*], in which using the ADM formulation of GR, the source is described as a canonical system in terms of its Hamiltonian. The ADM-Hamiltonian approach is equivalent to the multipolar post-Minkowskian approach and to the direct integration of relaxed Einstein’s equations, as long as the evolution of the source is concerned , but since Einstein’s equations are not solved in the wave zone, the radiative effects are only known with the same precision as the motion of the source. This framework has been extended to spinning binaries (see [726*] and references therein). Recently, an alternative way to compute the Hamiltonian of a post-Newtonian source has been developed, the effective field theory approach [358, 149, 627, 340], in which techniques originally derived in the framework of quantum field theory are employed. This approach was also extended to spinning binaries [626, 625]. ADM-Hamiltonian and effective field theory are probably the most promising approaches to extend the accuracy of PN computations for spinning binaries.
The effective one-body (EOB) approach developed at the end of the last century  and recently improved [247, 600*] (see, e.g., [240*, 249*] for a more detailed account) is an extension of PN schemes, in which the PN Taylor series is suitably resummed, in order to extend its validity up to the merger of the binary system. This approach maps the dynamics of the two compact objects into the dynamics of a single test particle in a deformed Kerr spacetime. It is a canonical approach, so the Hamiltonian of the system is computed, but the radiative part of the dynamics is also described. Since the mapping between the two-body system and the “dual” one-body system is not unique, the EOB Hamiltonian depends on a number of free parameters, which are fixed using results of PN schemes, of gravitational self-force computations, and of NR simulations. After this calibration, the waveforms reproduce with good accuracy those obtained in NR simulations (see, e.g., [240*, 249*, 600*, 61*]). In the same period, a different approach has been proposed to extend PN templates to the merger phase, matching PN waveforms describing the inspiral phase, with NR waveforms describing the merger [17*, 673*]. Both this “phenomenological waveform” approach and the EOB approach use results from approximation schemes and from NR simulations in order to describe the entire waveform of coalescing binaries, and are instrumental for data analysis .
To conclude this section, we mention that PN schemes originally treated compact objects as point-like, described by delta functions in the stress-energy tensor, and employing suitable regularization procedures. This is appropriate for BHs, and, as a first approximation, for NSs, too. Indeed, finite size effects are formally of 5-PN order (see, e.g., [239, 109*]). However, their contribution can be larger than what a naive counting of PN orders may suggest . Therefore, the PN schemes and the EOB approach have been extended to include the effects of tidal deformation of NSs in compact binary systems and on the emitted gravitational waveform using a set of parameters (the “Love numbers”) encoding the tidal deformability of the star [323*, 248*, 760*, 102*].
PN schemes are also powerful tools to study the nature of the gravitational interaction, i.e., to describe and design observational tests of GR. They have been applied either to build general parametrizations, or to determine observable signatures of specific theories (two kinds of approaches that have been dubbed top-down and bottom-up, respectively ).
Let us discuss top-down approaches first. Nearly fifty years ago, Will and Nordtvedt developed the PPN formalism [784*, 581*], in which the PN metric of an -body system is extended to a more general form, depending on a set of parameters describing possible deviations from GR. This approach (which is an extension of a similar approach by Eddington ) facilitates tests of the weak-field regime of GR. It is particularly well suited to perform tests in the solar system. All solar-system tests can be expressed in terms of constraints on the PPN parameters, which translates into constraints on alternative theories of gravity. For instance, the measurement of the Shapiro time-delay from the Cassini spacecraft  yields the strongest bound on one of the PPN parameters; this bound determines the strongest constraint to date on many modifications of GR, such as Brans–Dicke theory.
More recently, a different parametrized extension of the PN formalism has been proposed which, instead of the PN metric, expands the gravitational waveform emitted by a compact binary inspiral in a set of parameters describing deviations from GR [825*, 203*]. The advantage of this so-called “parametrized post-Einsteinian” approach – which is different in spirit from the PPN expansion, since it does not try to describe the spacetime metric – is its specific design to study the GW output of compact binary inspirals which are the most promising sources for GW detectors (see Section 3.1.1).
As mentioned above, PN approaches have also been applied bottom-up, i.e., in a manner that directly calculates the observational consequences of specific theories. For instance, the motion of binary pulsars has been studied, using PN schemes, in specific alternative theories of gravity, such as scalar-tensor theories [244*]. The most promising observational quantity to look for evidence of GR deviations is probably the gravitational waveform emitted in compact binary inspirals, as computed using PN approaches. In the case of theories with additional fundamental fields, the leading effect is the increase in the emitted gravitational flux arising from the additional degrees of freedom. This increase typically induces a faster inspiral, which affects the phase of the gravitational waveform (see, e.g., ). For instance, in the case of scalar-tensor theories a dipolar component of the radiation can appear [787*]. In other cases, as in massive graviton theories, the radiation has as in GR, but the flux is different. For further details, we refer the interested reader to [782*] and references therein.
[614, 613]. It took about thirty years to understand how to extend this computation at higher PN orders, consistently modeling the motion and the gravitational emission of a compact binary [109*, 785]. The state-of-the-art computations give the gravitational waveform emitted by a compact binary system, up to 3.5-PN order for non-spinning binaries in circular orbits , up to 3-PN order for eccentric orbits , and up to 2-PN order for spinning binaries . An alternative approach, based on the computation of the Hamiltonian , is currently being extended to higher PN orders [726, 457, 399]; however, in this approach the gravitational waveform is computed with less accuracy than the motion of the binary.
Recently, different approaches have been proposed to extend the validity of PN schemes up to the merger, using results from NR to fix some of the parameters of the model (as in the EOB approach [249, 600, 61, 240]), or matching NR with PN waveforms (as in the “phenomenological waveform” approach [17, 673]). PN and EOB approaches have also been extended to include the effects of tidal deformation of NSs [323, 248, 760, 102].
PN approaches have been extended to test GR against alternative theories of gravity. Some of these extensions are based on a parametrization of specific quantities, describing possible deviations from GR. This is the case in the PPN approach [784, 581], most suitable for solar-system tests (see [782, 783] for extensive reviews on the subject), and in the parametrized post-Einsteinian approach [825, 203], most suitable for the analysis of data from GW detectors. Other extensions, instead, start from specific alternative theories and compute – using PN schemes – their observational consequences. In particular, the motion of compact binaries and the corresponding gravitational radiation have been extensively studied in scalar-tensor theories [244*, 787*, 30].
[354, 32], or a BH undergoing small oscillations around a stationary configuration [487*, 316*, 95*].
In this approach, the spacetime is assumed to be, at any instant, a small deviation from the background geometry, which, in the cases mentioned above, is described by the Schwarzschild or the Kerr solution here denoted by . The deformed spacetime metric can then be decomposed as8 Einstein’s equations are linearized around the background solution, by keeping only first-order terms in (and in the other perturbation quantities, if present).
The simple expansion (15*) implies a deeper geometrical construction (see, e.g., ), in which one considers a family of spacetime manifolds , parametrized by a parameter ; their metrics satisfy Einstein’s equations, for each . The element of this family is the background spacetime, and the first term in the Taylor expansion in is the perturbation. Therefore, in the spacetime perturbation approach it is the spacetime manifold itself to be perturbed and expanded. However, once the perturbations are defined (and the gauge choice, i.e., the mapping between quantities in different manifolds, is fixed), perturbations can be treated as genuine fields living on the background spacetime . In particular, the linearized Einstein equations can be considered as linear equations on the background spacetime, and all the tools to solve linear differential equations on a curved manifold can be applied.
The real power of this procedure comes into play once one knows how to separate the angular dependence of the perturbations . This was first addressed by Regge and Wheeler in their seminal paper [641*], where they showed that in the case of a Schwarzschild background, the metric perturbations can be expanded in tensor spherical harmonics , in terms of a set of perturbation functions which only depend on the coordinates and . They also noted that the terms of this expansion belong to two classes (even and odd perturbations, sometimes also called polar and axial), with different behaviour under parity transformations (i.e., , ). The linearized Einstein equations, expanded in tensor harmonics, yield the dynamical equations for the perturbation functions. Furthermore, perturbations corresponding to different harmonic components or different parities decouple due to the fact that the background is spherically symmetric. After a Fourier transformation in time, the dynamical equations reduce to ordinary differential equations in .
Regge and Wheeler worked out the equations for axial perturbations of Schwarzschild BHs; later on, Zerilli derived the equations for polar perturbations [830*]. With their gauge choice (the “Regge–Wheeler gauge”, which allows us to set to zero some of the perturbation functions), the harmonic expansion of the metric perturbation is
It turns out to be possible to define a specific combination of the axial perturbation functions , and a combination of the polar perturbation functions which describe completely the propagation of GWs. These functions, called the Regge–Wheeler and the Zerilli function, satisfy Schrödinger-like wave equations of the form and represents nontrivial source terms. The energy flux emitted in GWs can be calculated straightforwardly from the solutions .
This approach was soon extended to general spherically symmetric BH backgrounds and a gauge-invariant formulation in terms of specific combinations of the perturbation functions that remain unchanged under perturbative coordinate transformations [555*, 346]. In the same period, an alternative spacetime perturbation approach was developed by Bardeen, Press and Teukolsky [75, 744*], based on the Newman–Penrose formalism [575*], in which the spacetime perturbation is not described by the metric perturbation , but by a set of gauge-invariant complex scalars, the Weyl scalars, obtained by projecting the Weyl tensor onto a complex null tetrad defined such that all their inner products vanish except . One of these scalars, , describes the (outgoing) gravitational radiation; it is defined as20*) with .
In this framework, the perturbation equations reduce to a wave equation for (the perturbation of) , which is called the Teukolsky equation . For a general account on the theory of BH perturbations (with both approaches) see Chandrasekhar’s book .
The main advantage of the Bardeen–Press–Teukolsky approach is that it is possible to separate the angular dependence of perturbations of the Kerr background, even though such background is not spherically symmetric. Its main drawback is that it is very difficult to extend it beyond its original setup, i.e., perturbations of Kerr BHs. The tensor harmonic approach is much more flexible. In particular, spacetime perturbation theory (with tensor harmonic decomposition) has been extended to spherically symmetric stars [753, 518, 266, 196*] (the extension to rotating stars is much more problematic ). As we discuss in Section 5.2.3, spacetime perturbation theory with tensor harmonic decomposition can be extended to higher-dimensional spacetimes. It is not clear whether such generalizations are possible with the Bardeen–Press–Teukolsky approach.
The sources describe the objects that excite the spacetime perturbations, and can arise either directly from a non-vanishing stress-energy tensor or by imposing suitable initial conditions on the spacetime. These two alternative forms of exciting BH spacetimes have branched into two distinct tools, which can perhaps be best classified as the “point particle” [250*, 179*, 569, 93*] and the “close limit” approximations [634*, 637*].
In the point particle limit the source term is a nontrivial perturbing stress-tensor, which describes for instance the infall of a small object along generic geodesics. The “small” object can be another BH, or a star, or even matter accreting into the BH. While the framework is restricted to objects of mass , it is generically expected that the extrapolation to yields at least a correct order of magnitude. Thus, the spacetime perturbation approach is in principle able to describe qualitatively, if not quantitatively, highly dynamic BHs under general conditions. The original approach treats the small test particle moving along a geodesic of the background spacetime. Gravitational back-reaction can be included by taking into account the energy and angular momentum loss of the particle due to GW emission [232, 445*, 548]. More sophisticated computations are required to take into account the conservative part of the “self-force”. For a general account on the self-force problem, we refer the interested reader to the Living Reviews article on the subject [623*]. In this approach is restricted to be a very small quantity. It has been observed by many authors [37, 718*] that promoting to the symmetric mass ratio describes surprisingly well the dynamics of generic BHs with masses .
In the close limit approximation the source term can be traced back to nontrivial initial conditions. In particular, the original approach tackles the problem of two colliding, equal-mass BHs, from an initial separation small enough that they are initially surrounded by a common horizon. Thus, this problem can be looked at as a single perturbed BH, for which some initial conditions are known [634*, 637].
A universal feature of the dynamics of BH spacetimes as given by either the point particle or the close limit approximation is that the waveform decays at late times as a universal, exponentially damped sinusoid called ringdown or QNM decay. Because at late times the forcing caused by the source term has died away, it is natural to describe this phase as the free oscillations of a BH, or in other words as solutions of the homogeneous version of Eq.(19*). Together with the corresponding boundary conditions, the Regge–Wheeler and Zerilli equations then describe a freely oscillating BH. In vacuum, such boundary conditions lead to an eigenvalue equation for the possible frequencies . Due to GW emission, these oscillations are damped, i.e., they have discrete, complex frequencies called quasi-normal mode frequencies of the BH [487*, 316*, 95*]. Such intuitive picture of BH ringdown can be given a formally rigorous meaning through contour integration techniques [506*, 95*].
The extension of the Regge–Wheeler–Zerilli approach to asymptotically dS or AdS spacetimes follows with the procedure outlined above and decomposition (16*); see also Ref. [176*]. It turns out that the Teukolsky procedure can also be generalized to these spacetimes [192*, 277*, 276*].
16*), but now augmented to include perturbations in matter fields, such as scalar or vector fields, or further polarizations for the gravitational field. Important examples where this formalism has been applied include scalar-tensor theories [668, 165*, 824*], Dynamical Chern–Simons theory [175, 554, 603], Einstein-Dilaton-Gauss–Bonnet [602*], Horndeski gravity [477, 478], and massive theories of gravity [135*].
[479*, 452] to GR in higher-dimensional spacetimes, with or without cosmological constant. Their approach generalizes the gauge-invariant formulation of the Regge–Wheeler-Zerilli construction to perturbations of Tangherlini’s solution describing spherically symmetric BHs.
Since many dynamical processes involving higher-dimensional BHs (in particular, the collisions of BHs starting from finite distance) can be described in the far field limit by a perturbed spherically symmetric BH spacetime, the Kodama and Ishibashi approach can be useful to study the GW emission in these processes. The relevance of this approach therefore extends well beyond the study of spherically symmetric solutions. For applications of this tool to the wave extraction of NR simulations see for instance [797*].
In the Kodama and Ishibashi approach, the -dimensional spacetime metric is assumed to have the form where is the Tangherlini solution and represents a small perturbation. Decomposing the -dimensional spherical coordinates into with angular coordinates , the perturbation can be expanded in spherical harmonics, as in the four-dimensional case (see Section 5.2.1). However, the expansion in is more complex than its four-dimensional counterpart: there are three classes of perturbations called the “scalar”, “vector” and “tensor” perturbations. The former two classes correspond, in , to polar and axial perturbations, respectively. These perturbations are decomposed into scalar (), vector () and tensor () harmonics on the -sphere and their gradients, as follows:
A set of gauge-invariant variables and the so-called “master functions”, generalizations of the Regge–Wheeler and Zerilli functions, can be constructed out of the metric perturbation functions and satisfy wave-like differential equations analogous to Eq. (19*). The GW amplitude and its energy and momentum fluxes can be expressed in terms of these master functions.
For illustration of this procedure, we consider here the special case of scalar perturbations. We define the gauge-invariant quantities
In summary, this approach can be used, in analogy with the Regge–Wheeler–Zerilli formalism in four dimensions, to determine the quasi-normal mode spectrum (see, e.g., the review [95*] and references therein), to determine the gravitational-wave emission due to a test source [98*, 94*], or to evaluate the flux of GWs emitted by a dynamical spacetime which tends asymptotically to a perturbed Tangherlini solution [797*].
The generalization of this setup to higher-dimensional rotating (Myers–Perry ) BHs is still an open issue, since the decoupling of the perturbation equations has so far only been obtained in specific cases and for a subset of the perturbations [564, 496*, 481].
Spacetime perturbation theory has also been used to study other types of higher-dimensional objects as for example black strings. Gregory and Laflamme [367, 368] considered a very specific sector of the possible gravitational perturbations of these objects, whereas Kudoh [495*] performed a complete analysis that builds on the Kodama–Ishibashi approach.
- Astrophysical systems. Perturbation theory has been applied extensively to the modelling of
BHs and compact stars, either without source terms, including in particular quasi-normal
modes [487, 316, 95], or with point particle sources. Note that wave emission from extended
matter distributions can be understood as interference of waves from point particles [400*,
693, 615]. Equations for BH perturbations have been derived for Schwarzschild [641*,
830*], RN , Kerr [744*] and slowly rotating Kerr–Newman BHs . Equations for
perturbations of stars have been derived for spherically symmetric [753, 518, 196] and slowly
rotating stars [197, 482].
Equations of BH perturbations with a point particle source have been studied as a tool to understand BH dynamics. This is a decades old topic, historically divided into investigations of circular and quasi-circular motion, and head-ons or scatters.
Circular and quasi-circular motion. Gravitational radiation from point particles in circular geodesics was studied in Refs. [551, 252*, 130*] for non-rotating BHs and in Ref.  for rotating BHs. This problem was reconsidered and thoroughly analyzed by Poisson, Cutler and collaborators, and by Tagoshi, Sasaki and Nakamura in a series of elegant works, where contact was also made with the PN expansion (see the Living Reviews article  and references therein). The emission of radiation, together with the self-gravity of the objects implies that particles do not follow geodesics of the background spacetime. Inclusion of dissipative effects is usually done by balance-type arguments [445, 446, 733, 338] but it can also be properly accounted for by computing the self-force effects of the particle motion (see the Living Reviews article  and references therein). EM waves from particles in circular motion around BHs were studied in Refs. [252, 130, 129].
Head-on or finite impact parameter collisions: non-rotating BHs. Seminal work by Davis et al. [250, 251] models the gravitational radiation from BH collisions by a point particle falling from rest at infinity into a Schwarzschild BH. This work has been generalized to include head-on collisions at non-relativistic velocities [660, 317, 524*, 93*], at exactly the speed of light [179, 93*], and to non-head-on collisions at non-relativistic velocities [269, 93*].
The infall of multiple point particles has been explored in Ref.  with particular emphasis on resonant excitation of QNMs. Shapiro and collaborators have investigated the infall or collapse of extended matter distributions through superpositions of point particle waveforms [400, 693, 615].
Electromagnetic radiation from high-energy collisions of charged particles with uncharged BHs was studied in Ref. [181*] including a comparison with zero-frequency limit (ZFL) predictions. Gravitational and EM radiation generated in collisions of charged BHs has been considered in Refs. [459, 460].
Head-on or finite impact parameter collisions: rotating BHs. Gravitational radiation from point particle collisions with Kerr BHs has been studied in Refs. [484, 483, 485, 486]. Suggestions that cosmic censorship might fail in high-energy collisions with near-extremal Kerr BHs, have recently inspired further scrutiny of these scenarios [71, 72] as well as the investigation of enhanced absorption effects in the ultra-relativistic regime [376*].
Close Limit approximation. The close limit approximation was first compared against nonlinear simulations of equal-mass, non-rotating BHs starting from rest . It has since been generalized to unequal-mass  or even the point particle limit , rotating BHs  and boosted BHs at second-order in perturbation theory . Recently the close limit approximation has also been applied to initial configurations constructed with PN methods .
- Beyond electrovacuum GR. The resurgence of scalar-tensor theories as a viable and important prototype of alternative theories of gravity, as well as the conjectured existence of a multitude of fundamental bosonic degrees of freedom, has revived interest in BH dynamics in the presence of fundamental fields. Radiation from collisions of scalar-charged particles with BHs was studied in Ref. [134*]. Radiation from massive scalar fields around rotating BHs was studied in Ref. [165*] and shown to lead to floating orbits. Similar effects do not occur for massless gravitons .
- Beyond four-dimensions and asymptotic flatness. The gauge/gravity duality and related
frameworks highlight the importance of (A)dS and higher-dimensional background spacetimes.
The formalism to handle gravitational perturbations of four-dimensional, spherically symmetric
asymptotically (A)dS BHs has been developed in Ref. , whereas perturbations of
rotating AdS BHs were recently tackled [192, 277, 276]. Gravitational perturbations of
higher-dimensional BHs can be handled through the elegant approach by Kodama and
Ishibashi [479*, 480], generalized in Ref.  to include perturbations of black strings.
Perturbations of higher-dimensional, rotating BHs can be expressed in terms of a single master
variable only in few special cases . The generic case has been handled by numerical
methods in the linear regime [270*, 395*].
Scalar radiation by particles around Schwarzschild-AdS BHs has been studied in Refs. [180, 178, 177]. We are not aware of any studies on gravitational or electromagnetic radiation emitted by particles in orbit about BHs in spacetimes with a cosmological constant.
The quadrupole formula was generalized to higher-dimensional spacetimes in Ref. [170*]. The first fully relativistic calculation of GWs generated by point particles falling from rest into a higher-dimensional asymptotically flat non-rotating BH was done in Ref. , and later generalized to arbitrary velocity in Ref. [94*]. The mass multipoles induced by an external gravitational field (i.e., the “Love numbers”) to a higher-dimensional BH, have been determined in Ref. .
[774, 14*, 707*, 93*, 489, 513*]. Furthermore, if the asymptotic states are an accurate description of the collision at all times, as for instance if the colliding particles are point-like, then one expects the ZFL to be an accurate description of the problem. [93*] and Appendix B2 of Ref. [513*]). Finally, one observes that the high-energy limit yields isotropic emission; when translated to a multipole dependence, it means that the energy in each multipole scales as in this limit.
The ZFL has been applied in a variety of contexts, including electromagnetism where it can be used to compute the electromagnetic radiation given away in -decay [181*, 455]; Wheeler used the ZFL to estimate the emission of gravitational and electromagnetic radiation from impulsive events ; the original treatment by Smarr considered only head-on collisions and computed only the spectrum and total emitted energy. These results have been generalized to include collisions with finite impact parameter and to a computation of the radiated momentum as well [513*, 93*]. Finally, recent nonlinear simulations of high-energy BH or star collisions yield impressive agreement with ZFL predictions [719*, 93*, 288*, 134*].
- Astrophysical systems. The zero-frequency limit for head-on collisions of particles was used by Smarr [707*] to understand gravitational radiation from BH collisions and in Ref.  to understand radiation from supernovae-like phenomena. It was later generalized to the nontrivial finite impact parameter case [513*], and compared extensively with fully nonlinear numerical simulations [93*]. Ref.  reports on collisions of an electromagnetic charge with a non-rotating BH in a spacetime perturbation approach and compares the results with a ZFL calculation.
- Beyond four-dimensional, electrovacuum GR. Recent work has started applying the ZFL to other spacetimes and theories. Brito  used the ZFL to understand head-on collisions of scalar charges with four-dimensional BHs. The ZFL has been extended to higher dimensions in Refs. [170, 513] and recently to specific AdS soliton spacetimes in Ref. .
In vacuum GR, a point-like particle is described by the Schwarzschild metric of mass . The gravitational field of a particle moving with velocity is then obtained by boosting the Schwarzschild metric. Of particular interest is the limiting case where the velocity approaches the speed of light . Taking simultaneously the limit so that the zeroth component of the 4-momentum, , is held fixed, , one observes an infinite Lorentz contraction of the curvature in the spatial direction of the motion. In this limit, the geometry becomes that of an impulsive or shock gravitational -wave, i.e., a plane-fronted gravitational wave with parallel rays, sourced by a null particle. This is the Aichelburg–Sexl geometry [16*] for which the curvature has support only on a null plane. In Brinkmann coordinates, the line element is:31*) is satisfied, where the Laplacian is in the flat 2-dimensional transverse space. Such a solution is given in closed analytic form by .
The usefulness of shock waves in modelling collisions of particles or BHs at very high energies relies on the following fact. Since the geometry of a single shock wave is flat outside a null plane, one can superimpose two shock wave solutions traveling in opposite directions and still obtain an exact solution of the Einstein equations, valid up to the moment when the two shock waves collide. The explicit metric is obtained by superimposing two copies of (31*), one with support at and another one with support at . But it is more convenient to write down the geometry in coordinates for which test particle trajectories vary continuously as they cross the shock. These are called Rosen coordinates, ; their relation with Brinkmann coordinates can be found in [420*] and the line element for the superposition becomes[766*]. This AH is the union of two surfaces,
for some functions to be determined. The relevant null normals to and are, respectively,
Instead of providing a bound on the inelasticity, a more ambitious program is to determine the exact inelasticity by solving the Einstein equations in the future of the collision. Whereas an analytic exact solution seems out of reach, a numerical solution of the fully nonlinear field equations might be achievable, but none has been reported. The approach that has produced the most interesting results, so far, is to solve the Einstein equations perturbatively in the future of the collision.
To justify the use of a perturbative technique and introduce a perturbation expansion parameter, D’Eath and Payne [257*, 258*, 259*] made the following argument. In a boosted frame, say in the negative direction, one of the shock waves will become blueshifted whereas the other will become redshifted. These are, respectively, the waves with support on and . The geometry is still given by (33*), but with the energy parameter multiplying terms ( terms) replaced by a new energy parameter (parameter ). For a large boost, , or in other words, in the boosted frame there are a strong shock (at ) and a weak shock (at ). The weak shock is regarded as a perturbation of the spacetime of the strong shock, and provides the expansion parameter to study this perturbation. Moreover, to set up initial conditions for the post-collision perturbative expansion, one recasts the exact solution on the immediate future of the strong shock, , in a perturbative form, even though it is an exact solution. It so happens that expressing the exact solution in such perturbative fashion only has terms up to second order:31*), since the latter are more intuitive than Rosen coordinates. The geometry to the future of the strong shock, on the other hand, will be of the form 36*). At higher orders, the problem can also be reduced to solving wave equations for , but now with sources provided by the perturbations of lower order [221*].
After obtaining the metric perturbations to a given order, one must still compute the emitted gravitational radiation, in order to obtain the inelasticity. In the original work [256*, 257*, 258*, 259*], the metric perturbations were computed to second order and the gravitational radiation was extracted using Bondi’s formalism and the Bondi mass loss formula. The first-order results can equivalently be obtained using the Landau–Lifshitz pseudo-tensor for GW extraction [420*]. The results in first and second order are, respectively:38*) are below the AH bound (35*), as they should. Secondly, and as we shall see in Section 7.6, the second-order result is in excellent agreement with results from NR simulations. Finally, as we comment in the next subsection, the generalisation to higher dimensions of the first-order result reveals a remarkably simple pattern.
The technique of superimposing two Aichelburg–Sexl shock waves [16*] was first used by Penrose in unpublished work but quoted, for instance, in Ref. [257*]. Penrose showed the existence of an AH for the case of a head-on collision, thus suggesting BH formation. Computing the area of the AH yields an upper bound on the fraction of the overall energy radiated away in GWs, i.e., the inelasticity. In the early 2000s, the method of superimposing shock waves and finding an AH was generalized to and non-zero impact parameter in Refs. [286*, 818*] and refined in Ref.  providing, in addition to a measure of the inelasticity, an estimate of the cross section for BH formation in a high-energy particle collision. A potential improvement to the AH based estimates was carried out in a series of papers by D’Eath and Payne [256*, 257, 258, 259]. They computed the metric in the future of the collision perturbatively to second order in the head-on case. This method was generalized to in first-order perturbation theory [420, 222] yielding a very simple result: . A formalism for higher order and the caveats of the method in the presence of electric charge were exhibited in . AH formation in shock wave collisions with generalized profiles and asymptotics has been studied in [19, 739, 31, 282].