An Introduction to Gravitational Wave Theory

TL;DR

This paper provides an overview of the theoretical foundations of gravitational waves, covering general relativity, detection methods, and binary system waveforms, serving as educational lecture notes.

Contribution

It offers a comprehensive introduction to gravitational wave theory, integrating fundamental concepts with recent developments for educational purposes.

Findings

Summarizes key principles of gravitational wave physics

Describes detection techniques and waveform modeling

Serves as an educational resource for students and researchers

Abstract

Introduction to the theoretical foundations of gravitational waves: from general relativity to detection and binary system waveforms. Lecture notes prepared for the MaNiTou summer school on gravitational waves. Draft chapter for the CNRS contemporary Encyclopaedia Sciences to be published by ISTE.

Tables

Table 1. Table 1: Characteristics of different GW detectors and the corresponding GW wavelength.

	Characteristic detector	GW frequency	$f_{\mathrm{GW}}L$	$L$ vs ${\lambda }_{\mathrm{GW}}$
	size (km)	detectability range (Hz)
LVK	$\sim 1$	$\sim 10-{10}^4$	$f_{\mathrm{GW}}L\ll 1$	$L\ll {\lambda }_{\mathrm{GW}}$
LISA	$\sim {10}^6$	$\sim {10}^{-4}-1$	$f_{\mathrm{GW}}L\sim 1$	$L\sim {\lambda }_{\mathrm{GW}}$
PTA	$\sim {10}^{17}$	$\sim {10}^{-9}-{10}^{-7}$	$f_{\mathrm{GW}}L\gg 1$	$L\gg {\lambda }_{\mathrm{GW}}$

Equations

$$10\text{Hz}\lesssim f_{\mathrm{GW}}\lesssim 5\text{kHz}\text{(LVK)}.$$

$$10^{-4}\text{Hz}\lesssim f_{\mathrm{GW}}\lesssim 1\text{Hz}\text{(LISA)}.$$

$$\text{few Hz}\lesssim f_{\mathrm{GW}}\lesssim 10^4\text{Hz}\text{(ET, CE...)}.$$

$$10^{-7}\text{Hz}\lesssim f_{\mathrm{GW}}\lesssim 10^{-9}\text{Hz}\text{(PTA)}.$$

$$f_{\mathrm{GW}}=\frac{1}{\pi }{\left(\frac{G\mathcal{M}}{c^3}\right)}^{-5/8}{\left(\frac{5}{256\tau }\right)}^{3/8}$$

$$\mathcal{M}\equiv \frac{(m_1{m}_2{)}^{3/5}}{(m_1+m_2{)}^{1/5}}$$

$$\tau =t-t_c$$

$$f_{\text{merger}}=\frac{1}{6^{3/2}\pi }\left(\frac{c^3}{Gm}\right).$$

$$f_{\text{merger}}\simeq 1.5\text{kHz}(\text{BNS})$$

$$f_{\text{merger}}\simeq 60\text{Hz}(\text{stellar mass BBH}).$$

$$f_{\text{merger}}\simeq 10^{-3}\text{Hz}(\text{supermassive Binary BHs})$$

$$T\sim 10^{-3}{f}_{\mathrm{low}}^{-8/3}{\left(\frac{c^3}{G\mathcal{M}}\right)}^{5/3}$$

$$h\sim \frac{4}{R}{\left(\frac{G\mathcal{M}}{c^2}\right)}^{5/3}{\left(\frac{\pi f_{\mathrm{GW}}}{c}\right)}^{2/3}.$$

$$R\sim 400\text{Mpc}$$

$$r\sim {\left(\frac{Gm}{f_{\mathrm{GW}}^2}\right)}^{1/3}$$

$$G_{\mu \nu }+\Lambda g_{\mu \nu }=\frac{8\pi G}{c^4}{T}_{\mu \nu },$$

$$\frac{8\pi G}{c^4}\simeq 10^{-43}{\mathtt{kg}}^{-1}{\mathtt{m}}^{-1}{\mathtt{s}}^2.$$

$${\varphi }^{\prime }(x^{\prime })=\varphi (x),v^{\prime }{}^{\mu }(x^{\prime })=\frac{\partial x^{\prime }{}^{\mu }}{\partial x^{\nu }}{v}^{\nu }(x),g_{\mu \nu }^{\prime }(x^{\prime })=\frac{\partial x^{\rho }}{\partial x^{\prime }{}^{\mu }}\frac{\partial x^{\sigma }}{\partial x^{\prime }{}^{\nu }}{g}_{\rho \sigma }(x).$$

$${\delta }_{\xi }\varphi :=\varphi (x)-{\varphi }^{\prime }(x)={\xi }^{\mu }{\partial }_{\mu }\varphi \equiv {\pounds }_{\xi }\varphi ,$$

$${\delta }_{\xi }{v}^{\mu }:=v^{\mu }(x)-v^{\prime }{}^{\mu }(x)={\xi }^{\nu }{\partial }_{\nu }{v}^{\mu }-v^{\nu }{\partial }_{\nu }{\xi }^{\mu }\equiv {\pounds }_{\xi }{v}^{\mu },$$

$${\delta }_{\xi }{g}_{\mu \nu }:=g_{\mu \nu }(x)-g_{\mu \nu }^{\prime }(x)={\xi }^{\rho }{\partial }_{\rho }{g}_{\mu \nu }+2g_{\rho (\mu }{\partial }_{\nu )}{\xi }^{\rho }=2{\nabla }_{(\mu }{\xi }_{\nu )}\equiv {\pounds }_{\xi }{g}_{\mu \nu },$$

$${\mathcal{L}}_{\mathrm{EH}}=\frac{c^3}{16\pi G}(R-2\Lambda )\sqrt{-g}$$

$${\pounds }_{\xi }\mathcal{L}={\partial }_{\mu }({\xi }^{\mu }\mathcal{L}).$$

$${\delta }_{\xi }{\mathcal{L}}_{\mathrm{M}}$$

$$={\pounds }_{\xi }{\mathcal{L}}_{\mathrm{M}}+{\partial }_{\mu }{\tilde{\theta }}^{\mu }={\partial }_{\mu }({\xi }^{\mu }{\mathcal{L}}_{\mathrm{M}}+{\tilde{\theta }}^{\mu }),$$

$${\delta }_{\xi }{\mathcal{L}}_{\mathrm{M}}={\partial }_{\mu }({\xi }^{\mu }{\mathcal{L}}_{\mathrm{M}}+{\tilde{\theta }}^{\mu })-\frac{\delta {\mathcal{L}}_{\mathrm{M}}}{\delta g_{\mu \nu }}{\pounds }_{\xi }{g}_{\mu \nu }.$$

$$j_{\xi }^{\mu }=\frac{c^3}{8\pi G}((G^{\mu }{}_{\nu }+\Lambda {\delta }_{\nu }^{\mu }){\xi }^{\nu }-{\nabla }_{\nu }{\nabla }^{[\mu }{\xi }^{\nu ]}).$$

$${\nabla }_{\nu }{\nabla }^{[\mu }{\xi }^{\nu ]}=\frac{1}{2}(R^{\mu }{}_{\nu }{\xi }^{\nu }-\square {\xi }^{\mu }+{\nabla }^{\mu }{\nabla }_{\nu }{\xi }^{\nu }),$$

$$Q_{\xi }[S]={\oint }_S{\nabla }_{\nu }{\nabla }^{[\mu }{\xi }^{\nu ]}d{S}_{\mu },$$

$$Q_{\xi }[S_2]-Q_{\xi }[S_1]={\int }_V{R}^{\mu }{}_{\nu }{\xi }^{\nu }d{V}_{\mu }\widehat{=}\frac{8\pi G}{c^4}{\int }_V\left(T^{\mu \nu }{\xi }_{\nu }-(\Lambda +\frac{T}{2}){\xi }^{\mu }\right)d{V}_{\mu }.$$

$$g_{0\mu }=(-1,0,0,0),$$

$$\square x^{\mu }={\Gamma }_{\nu \rho }^{\mu }{g}^{\nu \rho }={\partial }_{\nu }(\sqrt{-g}{g}^{\mu \nu })=0,$$

Full text

1. I Overview and characteristic scales

2. I.1 Aims

3. I.2 On wave-like solutions and relativity

4. I.3 Detectors and GW frequencies

5. I.4 Compact binary systems: orders of magnitude and characteristic scales

6. I.4.1 The chirp signal

7. I.4.2 Merger frequency

8. I.4.3 Time to merger

9. I.4.4 Amplitude and distance

10. I.4.5 Distance between objects at merger

11. I.5 Roadmap

12. II Einstein’s equations: general covariance, Noether’s theorem and gauge transformations

13. II.1 Einstein’s equations and general covariance

14. II.2 Noether’s theorem and diffeomorphisms as gauge symmetries

15. II.3 Gauge and observers

16. III Perturbative treatment of Einstein’s equations

17. III.1 The idea: general background spacetime

18. III.2 Weak-field approximation

19. III.3 Spin and helicity

20. III.4 De Donder and TT gauges

21. III.5 Vacuum solutions

22. III.6 Gauge-invariant description: independent and constrained degrees of freedom

23. III.7 Gauge-fixing with sources

24. IV Detection of GWs

25. IV.1 Coordinate displacements versus physical displacements

26. IV.2 Interferometers

27. V Generation of GWs from sources

28. V.1 Introducing sources

29. V.2 Source multipoles

30. V.3 Solving the wave equation with sources

31. VI Dissipation by gravitational waves

32. VI.1 Energy of gravitational waves

33. VI.2 Dissipation equations

34. VI.3 On the validity of the quadrupole formula

35. VI.4 Back-reaction

36. VII GWs from binary systems: elliptical, circular and hyperbolic orbits

37. VII.1 Newtonian equations

38. VII.2 Energy and angular momentum fluxes

39. VII.3 Back-reaction and waveform: circular orbits

40. VII.4 Back-reaction: Elliptical orbits

41. VII.5 Frequency content of elliptical orbits

42. VII.6 Hyperbolic orbits

43. VII.6.1 Linear memory effect and low-frequency GWs

44. VII.6.2 Capture due to GW emission

45. VIII GWs in curved space-time, cosmology

46. VIII.1 General background metric

47. VIII.1.1 Linearised equations

48. VIII.1.2 WKB approximation

49. VIII.2 FLRW metric: background

50. VIII.3 FLRW metric: gravitational waveforms

51. VIII.4 Measuring cosmological parameters with GWs: outline of GW cosmology

52. IX Acknowledgements

53. A Second order action for perturbations around any background solution

54. B Green’s functions

55. C Landau-Lifshitz approach

I Overview and characteristic scales

I.1 Aims

Indirect evidence of gravitational waves (GWs) has been available for many years (see e.g. Will:2014kxa for a review). Direct detection on the other hand had to wait until 2015 LIGOScientific:2016aoc , that is nearly 100 years after GWs were predicted to exist in general relativity (GR). These detections by the LVK collaboration, consisting of the network of LIGO LIGOScientific:2014pky , Virgo VIRGO:2014yos and more recently KAGRA KAGRA:2020tym GW interferometers, are ongoing with new GWs signals being observed on a weekly basis LIGOGracedb . In the future more sensitive detectors on earth, together with ones working in different frequency bands such as the Laser Interferometer Space Antenna (LISA) LISA:2024hlh as well as Pulsar Timing Arrays (PTAs), will lead to new observations of the universe, potential new discoveries, and unprecendented tests of general relativity, cosmology and astrophysics.

The aim of these lectures The aim of these first two chapters is to present the basic introductory material required to understand GWs. We will address some of following questions:

• •

What are GWs? How do they emerge from GR? How does one deal with the symmetries (diffeomorphism invariance) of GR to fix gauges and coordinates, and what do they imply for the stress energy tensor of GWs?

• •

To what GW frequencies $f_{\text{GW}}$ are current and future GW detectors sensitive? Why are those detectors designed to be sensitive to particular GW frequency ranges?

• •

For a source consisting of two bound compact binaries objects (such as black holes) of masses $m_{1}$ and $m_{2}$ at some distance $R$ from an observer, what is the characteristic frequency, amplitude etc of the GWs emitted? Up to what distances $R$ can such sources be detected?

• •

Using the quadrupole formula (which we derive) what is the waveform of the emitted GWs and how does it depend for example on the ellipticity of the bound orbit?

• •

What sources correspond to the GW events detected by LVK? Are there other possible GW sources? We give an example of compact binary sources on unbound orbits and discuss the GW memory effect.

• •

If we consider sources on cosmological distance scales, how are their amplitude, frequency e.t.c. affected by the cosmological expansion?

I.2 On wave-like solutions and relativity

Gravitational waves are a natural expectation from GR, simply because it is a relativistic theory of gravity. To understand why, let us first take a step back to non-relativistic Newtonian gravity. When the famous apple drops on Newton’s head, the mass distribution of the Earth changes, and so does the gravitational field created. In Newton’s time, this variation, however negligible, was assumed to be the effect of some instantaneous “action at a distance”. After the discovery that the speed of light is finite, and that all effects in our universe appear to follow this causal limitation, it seems natural to expect that also the variations of the gravitational field will not be felt instantaneously in the whole universe, but will rather be propagated at the speed of light — or less.111We will see that Einstein’s relativistic theory of gravity, GR, predicts that — whatever their wavelength — these variations propagate at exactly the speed of light, and if some future experiment shows that they propagate at a lesser speed, then this would be an explicit violation of GR. The propagation of this perturbation of the gravitational field is intuitively what we call a gravitational wave.

Conceptually, a gravitational wave is similar to a water wave or an electromagnetic (EM) wave. However, while those propagate a modification in the depth of water or the intensities of the electromagnetic field, a GW propagates a modification of the structure of spacetime itself. As for producing one, it is natural to expect that the Earth emits GWs when orbiting the sun, thus carrying away energy and making the orbit decay, just like a charged particle emits EM waves when moving on an accelerating trajectory.

In practise, however, understanding GWs is very subtle for a number of reasons. First of all, the concept itself of propagation makes reference to a background spacetime, and in GR there is no fixed background structure. Splitting the dynamical spacetime into a reference background and a perturbation on top of it is a delicate process in which potential ambiguities have to be dealt with. In fact many decades passed before a consensus was reached, to the point that, famously, Einstein himself initially doubted the physical existence of GW, for reasons that we will briefly review and clarify below. Secondly, if we think of GWs as waves propagating in a medium, this medium is extraordinarily rigid: the waves go as fast as possible and have very tiny amplitudes. To give an idea, the power emitted by the Earth-Sun system in the form of GW is around 200W! This rigidity has to do with the weakness of the gravitational coupling constant. One may think that gravity is strong when for e.g. trying to beat a high jump record, or when skydiving, but this strength is ridiculously small compared to the much much stronger electro-magnetic force that dominates our daily life. These two points — background independence and weakness of the signal — are typical issues that one has to face when studying GWs.

Indeed only very massive and energetic objects can produce GWs of amplitudes that are actually detectable. Amongst the most massive and compact astrophysical objects known are black holes (BH), neutron stars (NS) and white dwarfs (WD). The GW sources detected to date by the LVK collaboration are all ‘compact binary systems’ made of a bound pair of BH and/or NS on closed orbits. As a result of the energy lost through GW emission, the two bodies making up the bound system approach closer to each other, inspiralling inwards, and eventually merging into one final object.In fact the GW signals detected by the LVK collaboration correspond to the last moments in the life of these systems including their merger — they are known as ‘compact binary coalescences’ (CBC’s). For comparison with the earth-sun system mentioned above, the energy emitted in GWs by the very first detected GW event GW150914 LIGOScientific:2016aoc , which was due to the coalescence of two BHs of masses $m_{1}\sim 36M_{\odot}$ and $m_{1}\sim 29M_{\odot}$, was almost $10^{48}$Joules in 0.2 seconds.

The direct detection of GWs can be used to test many aspects of gravity, for instance in the strong field regime, see e.g. Yunes:2013dva , as well as to probe cosmology, as will be discussed later. Indeed, the gravitational interaction is so weak that the universe is almost completely transparent to a gravitational wave. As a consequence, one can potentially collect pristine information about any cosmological era through GWs, and in particular through the detection and characterisation of a stochastic gravitational wave background. Sources relevant to cosmology include primordial GWs produced during inflation but there are also potential new sources to be discovered, such as primordial black holes, cosmic strings, and other exotic objects, see e.g. Caprini:2018mtu ; Caprini:2024ofd for reviews.

Our aim in these lectures is not to provide an introduction to the broad set of fascinating GW sources, confirmed or hypothetical, nor to the many creative ideas to detect them that have been proposed, investigated and realised in practise; but only to provide an introduction to the field, and to that end, we decided to focus on the most common type of sources, and most common type of detectors: CBCs and laser interferometers. In the rest of this overview section we review the characteristic properties of GWs emitted by CBCs and the relevant frequency bands of laser interferometers, in particular explaining why LVK detectors are sensitive to the merger of stellar mass BH, whilst LISA for example to that of supermassive BHs. The rest of the chapter will present the theoretical derivation of GWs from GR.

I.3 Detectors and GW frequencies

The LVK interferometers and future LISA detector are essentially Michelson-Morley interferometers, designed to be as sensitive as possible to time-varying changes in the separation between two freely falling test-masses — mirrors in the case of interferometers. The invariant distance between the test masses varies when a GW passes (see Section IV), leading to a change in the observed interference pattern in the detector.

• •

The LVK interformeters are on earth (in Livingston and Handford in the USA, in Pisa in Europe, and in Kamioka in Japan) and have a typical arm length $L\sim 3$km. They are sensitive to GWs with frequency of order

[TABLE]

• •

The LISA interferometer LISA:2024hlh was adopted by ESA on the 25th january 2024, and should be operational in 2037. The distance between the spacecraft which make up arms of LISA is $L\sim 2.5\cdot 10^{6}$km. LISA will be sensitive to GWs with frequencies in the range

[TABLE]

• •

There are plans to build new interferometers on earth beyond LVK. These include the Einstein Telescope in Europe Sathyaprakash:2012jk and Cosmic Explorer in the USA Reitze:2019iox , both of which should have $L\sim 10$km, and

[TABLE]

• •

An alternative to interferometers are PTAs which search for GWs by exploiting the variation in distance $L\sim 10^{17}$km between the earth and a typical distant galactic pulsar due to GWs. Pulsars emit EM pulses with extreme regularity $\Delta t$, typically of the order of milliseconds. If GWs are present, then as the EM pulses propagate from the pulsar to the earth, the observed $\Delta t_{obs}$ will be modulated. PTA experiments searching for these modulations are sensitive to GWs in the frequency band of the inverse year,

[TABLE]

In 2023 different PTA experiments presented strong evidence for the existence of a stochastic GW background EPTA:2023fyk ; NANOGrav:2023gor ; Reardon:2023gzh ; Xu:2023wog .

• •

Prior to the success of interferometers, there was an effort pioneered by Weber in the 60’s to use resonant bars to detect GWs, building material bars whose acoustic modes would resonate at a frequency as near as possible to that expected from the optimal sources Weber1960 . These experiments would typically have a narrow-band sensitivity around $10^{3}$Hz. In spite of constant experimental evolution throughout the 90’s, no observation has occurred in this way.

Table 1 summarises the different characteristics of the existing experiments, and in particular the ratio of their characteristic size $L$ to the GW wavelength $\lambda_{\rm{GW}}=c/f_{\rm{GW}}$.

I.4 Compact binary systems: orders of magnitude and characteristic scales

LVK and LISA were conceived in order to be sensitive to the particular range of frequencies that are not only within experimental reach, but also that are likely to constitute a rich source according to the known astrophysical data. Amongst those GW sources are compact binary systems. We now focus on orders of magnitude and characteristic scales for such compact binary system, consisting of two masses $m_{1,2}$ at a distance $R$ from the detectors, see figure 1. The expressions given here will be derived later in section VII. Furthermore, the expansion of the universe, neglected here, is considered in Section VIII.

As shown in figure 2, as a consequence of GW emission, the two masses $m_{1,2}$ approach each other — the inspiral phase — until they merge — the merger phase — and form a single object. This object will keep radiating GWs, in the so-called ringdown phase, until it settles down to an equilibrium state (which for BHs is expected to be represented by the Kerr or Schwarzschild solutions, according to theoretical and numerical evidence) after which no further emission occurs. The typical corresponding waveform, related to the GW amplitude, is shown in figure 2 as a function of time.

• •

The inspiral phase can be understood with perturbation theory (the “post-Newtonian (PN) expansion” of the Einstein equations) presented below, more details in Thorne:1980ru ; Blanchet:2006jqj ; poisson .

• •

The merger phase generally requires numerical relativity, or other techniques such as effective one-body techniques, see e.g. DeruelleUzan2018 for an introduction. These techniques can also applied to the late inspiral phase, in particular to understand accurately the GW signal as the merger is reached.

• •

The ringdown phase can also be approached with perturbative methods, namely BH perturbation theory, see e.g. Kokkotas:1999bd .

I.4.1 The chirp signal

During the inspiral phase the GW frequency increases with time according to the well-known chirp signal. Using the dominant quadrupolar mode contribution for point masses $m_{1}$ and $m_{2}$ (with spins set to zero), and assuming circular orbits, the time dependence of the frequency is given by

[TABLE]

see Eq. (212), section VII.3. Here the chirp mass is

[TABLE]

and

[TABLE]

is the time to coalescence, with $t_{c}$ the coalescence time. Clearly Eq. (5) will break down before $\tau=0$ where formally $f_{\rm GW}$ diverges. We thus define a “merger time” $t_{\rm merger}<t_{c}$ up to which Eq. (5) is assumed valid, and whose meaning we now discuss.

I.4.2 Merger frequency

We now assume that the two objects are Schwarzschild BHs, and that merger occurs at the innermost stable circular orbit (ISCO) namely a distance $a=6GM/c^{2}$ with $m=m_{1}+m_{2}$. It then follows from Keplers laws (see Sec. VII) together with Eq. (5) that

[TABLE]

(Note given a length scale $a$ and a mass $m$, $\sqrt{Gm/{a^{3}}}$ has dimensions of frequency. Setting $a=6Gm/c^{2}$ gives, modulo factors of $2\pi$, Eq. (8).)

• •

For a binary neutron stars (BNS) system, with say $m_{1,2}\sim 1.4M_{\odot}$ then Eq. (8) gives

[TABLE]

This is in the upper part of the LVK frequency band.

• •

For a stellar mass binary black hole (BBH) system with for instance $m_{1,2}\sim 35M_{\odot}$,

[TABLE]

This is right in the frequency band of LVK.

• •

For a massive black hole binary (MBHB) system with for instance $m_{1,2}\sim 10^{6}M_{\odot}$

[TABLE]

which is in the frequency band of LISA.

• •

Notice that PTA frequencies do not correspond to the merger frequency of any know astrophysical system. Rather, they correspond to the inspiral phase of super MBHB at times much before merger, as can be seen from Eq. (8). Hence these are on broad orbits, with periods of the order of years.

Detailed figures for the Ligo-Virgo, CE and ET sensitivities as a function of frequency can be found for instance in Maggiore:2024cwf . LISA sensitivities can be found for instance in LISA:2024hlh .

I.4.3 Time to merger

If GWs emitted during the inspiral enter the frequency band of a given detector at frequency $f_{\rm low}$, then it is straightforward to integrate Eq. (5) from $f_{\rm low}$ to $f_{\rm merger}$ to find the total duration of the GW signal as will be observed by the experiment. Assuming $f_{\rm merger}\gg f_{\rm low}$ for simplicity, one finds that the total duration of the signal is

[TABLE]

• •

For BNS entering the LVK band with $f_{\rm low}\sim 20$Hz, this gives $T\sim 4$ minutes.

• •

For BNS entering the ET band with $f_{\rm low}\sim 1$Hz, then $T\sim 5$ days.

(This implies for example that effects of the rotation of the earth cannot be neglected when calculating the GW properties in more detail, see e.g. Iacovelli:2022bbs and references within. Also one might expect other GW signals to be produced in such a long period, overlapping with the BNS one. This makes data analysis more complex Samajdar:2021egv .)

• •

For stellar mass BHs, with say $m_{1,2}\sim 35M_{\odot}$ entering the LVK band with $f_{\rm low}\sim 20$Hz, then $T\sim 0.1$ seconds.

• •

For stellar mass BHs, with say $m_{1,2}\sim 35M_{\odot}$ entering the ET band with $f_{\rm low}\sim 1$Hz, then $T\sim 300$ seconds.

• •

For MBHB with $m_{1,2}\sim 10^{6}M_{\odot}$ entering the LISA band with $f_{\rm low}\sim 10^{-4}$Hz, then $T\sim 1$ month.

(The orbital motion of LISA will thus also be non-negligible and e.g. Doppler effects must considered. Furthermore other LISA sources will overlap with the MBHB signal.)

I.4.4 Amplitude and distance

The dimensionless amplitude of the GW signal scales with distance $R$ to the source and GW frequency $f_{\rm GW}$ as

[TABLE]

As an example, consider say stellar mass BBH with $m_{1,2}\sim 35M_{\odot}$ for which $f_{\rm merger}\sim 60$Hz. In order to generate (at merger) a signal with amplitude $h\sim 10^{-21}$, which is accessible to LVK, requires

[TABLE]

which is of the order of cosmological scales (for comparison, the observable universe has a scale of $c/H_{0}\sim$Gpc, where $H_{0}$ is the Hubble constant).

Clearly from Eq. (13), given ${\cal{M}}$ and $f_{\rm GW}$, the more sensitive a detector, namely the smaller $h$ can be detected, the further one can detect a given GW source. The “detection volume” of LVK has been steadily increasing with the different observing runs of LVK obviously leading to increasing numbers of detected GW events.

Notice that if such a GW signal is detected, then from the time dependence of the GW frequency one can directly obtain chirp mass Eq. (5). With that, from the amplitude one can obtain the distance through Eq. (13). Distance measurements can thus directly be obtained with GW observations from binaries, hence their name standard sirens Schutz:1986gp ; Holz:2005df . This should be contrasted with the case of EM observations (standard candles) for which the determination of the distance is particularly difficult. See Sec. VIII for more information about GWs as distance indicators and their use in cosmology.

I.4.5 Distance between objects at merger

When GWs are emitted with frequency $f_{\rm GW}$, the two bodies in the compact binary are separated by a characteristic scale

[TABLE]

(see also the discussion after Eq. (8)). Since, from Eq. (5), the GW frequency increases during inspiral, the distance $r$ between the two bodies decreases. The minimum distance is at the merger frequency $f_{\rm merger}$. For example, for stellar mass BBH with $m_{1,2}\sim 35M_{\odot}$ and $f_{\rm merger}\sim 60$Hz then from Eq. (15) $r\sim{\cal{O}}(100)$km.

A distance $r\sim{\cal{O}}(100)$km is tiny compared to the characteristic size of a star. Some of the most dense stars in the universe — for instance WDs — have a size $\sim 10^{3}$km. Main sequence stars have a size which can go up to millions of km. Thus if GW signals are seen from objects which reach minimum approach distances $\sim{\cal{O}}(100)$km, those objects cannot be stars as they would already have collided. We must be dealing with BH (or possibly NS) for which the minimum distance will be determined by the Schwarzschild radius.

I.5 Roadmap

Having gone through the overview and discussed these orders of magnitude, the remainder of this chapter aims to derive formal results on GWs starting from Einstein’s equations. Many introductions and reviews on GWs already exist, see for instance poisson ; Maggiore1 ; Maggiore2 ; Blanchet:2006jqj ; AnderssonBook ; DeruelleUzan2018 to mention a few. It is a rich and intricate topic, and each of these reviews tends to have a different angle on it, whose mutual compatibility may not always be clear to somebody entering the field. We have strived at presenting the material in a way that allows one to understand how the different approaches relate to one another. We have also strived to spend time on some of the subtleties and delicate conceptual aspects of GR and GWs which are often left to the side in gravitational wave introductions, such as gauge dependencies, asymptotic charges and memory effects, and which are becoming more and more relevant as theoretical research and experiments advance into more accurate comparisons. In these lecture notes, a reader will therefore find discussions of questions of such as coordinate invariance, diffeomorphisms and gauge invariance, spin and helicity, the controversies about the stress energy tensor of GWs, GW memory effects, and Noether charges. The hope is that although most of the more advanced material is not needed for a first introduction, its inclusion here will stimulate the reader, and provide a useful reference for delving further into the topic.

II Einstein’s equations: general covariance, Noether’s theorem and gauge transformations

II.1 Einstein’s equations and general covariance

Einstein’s great discovery about gravitation was that it can be understood as the manifestation of the curvature of spacetime. In Wheeler’s words, spacetime tells matter how to move, matter tells spacetime how to bend.222While pictorially charming, this statement is not exactly true: spacetime can be extraordinarily bent even in the absence of matter, as black hole solutions show. Understanding gravity as a dynamical spacetime has changed profoundly our understanding of inertia. If we go back to Galilean relativity, an inertial observer is defined as one moving on a straight line at constant velocity. Special relativity introduces a non-trivial mixing of space and time, but leaves this notion unaffected: Inertial observers are still moving on a straight line, even though they are now related by Poincaré transformations as opposed to Galilean transformations, so to account for the experimental invariance of the speed of light. But in a curved spacetime, straight lines may no longer exist. The notion that encompasses them is the one of geodesics, which describe free-falling observers. Constant motion on a straight line is simply the flat-spacetime version of free falling. The understanding offered by general relativity thus has the merit of not only explaining gravity, but also explaining the origin of inertia. On the other hand, it changes the perspective on it radically: You reading these notes at your desk are inertial in Newton’s terms, but accelerated in Einstein’s, since you are being held by the ground against Earth’s gravitational attraction and not following a geodesic.

Another profound consequence of a dynamical spacetime metric is that the field equations of gravity and matter are covariant under general coordinate transformation, as we will review below, introducing a new paradigm that goes under the name of principle of general covariance. In a curved spacetime, there are no more preferred Cartesian coordinates, no more Poincaré transformations relating inertial observers, and familiar physical concepts such as time evolution and energy become surprisingly subtle. These aspects of GR are often glossed over in lectures aiming at introducing gravitational waves, where one can blissfully rely on the background spacetime introduced by the weak field approximation and ignore most of them. However we believe they are important in order to better appreciate some of the properties of gravitational waves, provide an understanding that is more conceptual and less application-driven, and mostly because frankly who’d need yet another introduction to GWs if we didn’t attempt something different? So we will briefly review these aspects below, and use them as benchmark to discuss some conceptual aspects of gravitational waves. For instance, the lack of preferred clocks in a curved spacetime is relieved in the weak field approximation, where one can use the flat Minkowski background to introduce a class of Cartesian observers, and select their proper time as preferred time. But the lack of well-defined notion of energy density is a subtlety that persists also in the weak-field approximation, and has to be dealt with.

Let us start by recalling Einstein’s equations

[TABLE]

where $G_{\mu\nu}:=R_{\mu\nu}-\frac{1}{2}Rg_{\mu\nu}$ is the Einstein tensor with $R_{\mu\nu}$ the Ricci tensor and $R$ the Ricci scalar, and $T_{\mu\nu}$ is the (symmetric) stress-energy tensor of matter. We use the definitions and conventions of poisson , in particular mostly-plus convention for the spacetime metric $g_{\mu\nu}$. The constants $G/{c^{4}}$ and $\Lambda$ are respectively the relativistic gravitational coupling constant and the cosmological constant. The first can be determined from local gravitational experiments to be

[TABLE]

This value is stupendously small, and it is the origin of the ‘rigidity’ of spacetime mentioned in the overview section. The smallness of this parameter has, on the other hand, a positive side: the gravitational force is so weak that many of the observed phenomena, and virtually all solar system experiments, can be studied using the weak field approximation, namely a perturbative expansion around the Minkowski metric. This is quite helpful because Einstein’s equations are non-linear and it is in general very difficult to find exact solutions. Strong gravity effects occur only near very compact objects, and to study them one has to resort to numerical techniques, or be able to push the perturbative treatment to high orders.

The cosmological constant $\Lambda$ can be determined from the observed acceleration of the expansion of the universe assuming homogeneity and isotropy on large scales, and turns then out to be $\Lambda\simeq 10^{-52}\,{\tt m}^{-2}$. This coupling constant can also be interpreted as a sort of averaged ‘vacuum energy’ density, often referred to as dark energy since it is not associated to visible matter, and whose value is $\rho_{DE}={\Lambda c^{2}}/G\sim 10^{-28}$kg/m3. The presence of $\Lambda$ affects the propagation of gravitational waves on cosmological distances, but it can be ignored for a first understanding of the perturbative treatment. We will set $\Lambda=0$ for now, and restore it below in Section VIII when discussing cosmological effects.

Analysis of the 10 field equations in Eq. (16) shows (see subsection III.6) that: four are redundant, because of the Bianchi identities; four are elliptic, hence describe gravitational degrees of freedom constrained by the sources; two are hyperbolic, hence contain independent degrees of freedom. This three-sided structure is a common feature to Maxwell and Yang-Mills theories, with the role of the Gauss constraint generating gauge transformation replaced by the so-called Hamiltonian and vector constraints generating diffeomorphisms, and it is our first indication that coordinate transformations are a gauge symmetry. A second indication comes from Noether’s theorem, but before talking about it, let us review how coordinate transformations act.

Recall that a tensor is a quantity that transforms homogeneously under general coordinate transformations $x^{\mu}\to x^{\prime}{}^{\mu}(x^{\nu})$. For instance a scalar field, a vector field and the metric transform respectively as

[TABLE]

The transformation law of a scalar is such that its value at one point $P$ is the same after the diffeomorphism, since both $x$ and $x^{\prime}$ identify the same point, just in different coordinates.333This is sometimes misstated by saying that scalars are invariant under coordinate transformations, which is not true. A quantity is invariant under coordinate transformations if it satisfies the stronger property that its value does not depend on the coordinates used, which for the scalar field would be the equation $\phi(x^{\prime})=\phi(x)$. This is not true in general, but only for isometries — more on this below. A typical example of coordinate invariance is the integral over the whole manifold of a scalar times the volume form, as experience from solving integrals via change of coordinates should show. The vector and metric do the same, but furthermore their indices are mixed up using the Jacobian of the coordinate transformation, or its inverse. A vector field is said to transform as a contravariant tensor of order one, and the metric as a covariant tensor of order two. Since coordinate transformations are typically restricted to be differentiable, namely continuous and connected to the identity, they are also invertible, and correspond to mathematical transformations called diffeomorphisms. In this language, (18) is a diffeomorphism of the metric.

If the coordinate transformation is infinitesimal, we can write it as $x^{\prime}{}^{\mu}=x^{\mu}+\xi^{\mu}(x)$, and approximate the transformation rules in (18) using the Taylor expansion (applied to both the field’s argument and the Jacobian). This defines the infinitesimal transformations

[TABLE]

where we introduced the use of round brackets for index symmetrization (and we will later on also use square brackets for index anti-symmetrization). Notice that in all cases we recover as infinitesimal transformation the Lie derivative. This is a general result valid for any tensor. The third equality in (21) is on the other hand special to the metric tensor, and follows from the expression of the connection in terms of the metric.

Being written in terms of tensors, Einstein’s equations are automatically covariant under general coordinate transformations. This is the principle of general covariance, that played a key role in guiding Einstein to formulate his theory. One immediate implication is that locally we can always find coordinates such that the metric takes the Minkowski expression at a point, which is one version of the principle of equivalence. A more subtle implication is that coordinate transformations must be symmetries of the theory, in other words a solution can be equivalently written in any coordinate system. To understand this point, let us consider the Lagrangian description of the dynamics. The field equations (16) are Euler-Lagrange equations of ${\cal L}={\cal L}_{\scriptscriptstyle\rm EH}+{\cal L}_{\scriptscriptstyle\rm M}$, where

[TABLE]

is the Einstein-Hilbert Lagrangian density, and ${\cal L}_{M}$ the matter contribution, left arbitrary for the moment. Here $g=\rm{det}(g_{\mu\nu})$, and the word density has a double meaning: in the physical sense, since $c{\cal L}$ has the dimensions J m*-3* of an energy density, but also in the mathematical sense, since $\sqrt{-g}$ makes it transform not as a scalar but as a scalar density of weight 1. This has the following consequence. Recall that the general formula for the variation of a determinant is $\delta g=gg^{\mu\nu}\delta g_{\mu\nu}$. This implies that $\pounds_{\xi}\sqrt{-g}=\frac{1}{2}\sqrt{-g}g^{\mu\nu}\pounds_{\xi}g_{\mu\nu}=\sqrt{-g}\nabla_{\mu}\xi^{\mu}$, whence $\pounds_{\xi}(\sqrt{-g}\phi)=\partial_{\mu}(\sqrt{-g}\xi^{\mu}\phi)$ for any scalar $\phi$. Thus, a Lagrangian density transforms as a total derivative under diffeomorphisms:

[TABLE]

Since total derivative do not affect the field equations, the transformed solutions are still solutions. We conclude that any diffeomorphism is a symmetry of a general covariant Lagrangian.

We stress that this result relies crucially on the fact that in a general covariant theory the metric is a dynamical field, and not a fixed background. To appreciate this point and the difference with non-general relativistic physics, let us consider the matter Lagrangian ${{\cal L}}_{\scriptscriptstyle\rm M}(g,\psi)$, which depends on both the metric $g_{\mu\nu}$ and the matter fields, which we denote collectively as $\psi$. Applying the chain rule, we find

[TABLE]

where the second equality follows from (19) and (21), and $\tilde{\theta}^{\mu}$ is the boundary term that arises because the Lagrangian depends on derivatives of the fields as well. Since the result is a total derivative, the field equations are unchanged, and this means that diffeomorphisms are symmetries also of the matter sector. But this conclusion relies crucially on treating the metric as a dynamical variable! In non-general relativistic physics the metric is a non-dynamical, ‘background’ field. Accordingly, there is no variation with respect to the metric, and no second term after the first equality of (24). Lacking this term the second equality breaks down, and $\delta_{\xi}{{\cal L}}_{\scriptscriptstyle\rm M}$ is no longer a boundary term. We find instead

[TABLE]

This means that only those diffeomorphisms which are isometries are symmetries of the non-general relativistic physics, e.g. the familiar Poincaré invariance of special relativity. Whereas an arbitrary diffeomorphism does not map a solution into a new solution in non-general relativistic physics, and there is no invariance under general coordinate transformations.

The discussion highlights why general covariance is often referred to as background independence, namely the absence of any fixed background metric in the theory, or as diffeomorphism invariance, since every physical observable should be independent of the coordinate used to describe it. General covariance, background independence, or diffeomorphism invariance, are thus different terms used to capture the same underlying property of general relativity.

The fact that every solution can be equivalently described in any coordinate system has a useful analogy with electromagnetism, where every solution can be described in any choice of gauge for the Maxwell potential. It is actually much more than an analogy, there is in fact a precise mathematical sense in which gauge transformations in Maxwell and Yang-Mills theories have the same property of coordinate transformations in general relativity, which we discuss next. Making this analogy precise is also useful in the context of gravitational waves as it will allow us to understand the origin of their gauge dependence and distinction between physical and unphysical modes.

II.2 Noether’s theorem and diffeomorphisms as gauge symmetries

Noether’s theorem proves that every differentiable444Namely, continuous and including the identity transformation. symmetry defines a current $j^{\mu}$ which is conserved on solutions, namely $\nabla_{\mu}j^{\mu}\,\hat{=}\,0$. Here the symbol $\,\hat{=}\,$ means an equality valid only for solutions, or ‘on-shell’, in theoretical physics jargon. Conserved currents are extremely useful to study the properties of the dynamics of the system, and to extract general physical predictions. We have seen above that diffeomorphisms are symmetries of a general covariant Lagrangian, continuous and connected to the identity. Therefore Noether’s theorem guarantees that there will be a conserved current associated to any diffeomorphism. However, the application of Noether’s theorem to general relativity is quite subtle. Let us first recall the difference between ‘proper symmetries’ and ‘gauge symmetries’. Both map solutions of the field equations into new solutions. If the new solution is physically distinguishable, we say that it is a proper symmetry, or a physical symmetry. If the new solution is on the other hand physically indistinguishable, we say that it is a gauge symmetry: This typically occurs when there is a redundancy of the field equations, which leaves some quantities undetermined but irrelevant for the physics. Noether’s theorem provides a simple test to distinguish the two cases: in the latter, the Noether current itself vanishes on shell, and not just its divergence. This is precisely the case with diffeomorphisms in general relativity. In fact, the conserved current associated with diffeomorphisms of the Einstein-Hilbert Lagrangian (22) is given by

[TABLE]

One can immediately verify using the Bianchi identities that $\nabla_{\mu}j^{\mu}_{\xi}\,\hat{=}\,0$. On the other hand, the first term above vanishes on-shell, and the second term is a total derivative. Therefore, the Noether current itself vanishes on-shell, as anticipated, and there are no conserved quantities (in the absence of boundaries). Trivial conserved quantities is a hallmark of gauge symmetries as opposed to physical symmetries, hence the result provides a precise mathematical sense in which coordinate transformations in general relativity have the same status as gauge transformations in Maxwell and Yang-Mills theories.555A more rigorous approach is to look at the symplectic 2-form, and show that it is degenerate along gauge transformations and diffeomorphisms.

For this reason, diffeomorphisms are also referred to as the gauge symmetry of general relativity, and fixing a coordinate choice as fixing the gauge in general relativity.

Having said so, there is a special situation that stands out: when the diffeomorphism corresponds to an isometry, namely a transformation that does not change the metric. This occurs when (21) vanishes, and the corresponding equation $\nabla_{(\mu}\xi_{\nu)}=0$ is called Killing equation, and $\xi$ a Killing vector. One should keep in mind that for a generic metric, this equation does not admit any solutions: Isometries occur only for very special metrics. These special metrics are, however, important for physical applications,666A cow is always a spherical object in the initial investigations of a theoretical physicist. hence Killing vectors play an important role. First of all, anyone who is familiar with the study of geodesics on spacetimes with isometries knows that there are conserved quantities associated with the Killing vectors, and which can be derived as Noether charges for the test particles’ dynamics.

More importantly for us, isometries play also an important role in the study of gravitational waves, because the dynamics of perturbations on a given background is such that the isometries of the background induce proper symmetries for the perturbations. We will see this in details in Section III below.

Before moving on, let us also mention another aspect in which isometries are important for the full theory. This is a more advanced topic, and will not be needed in the following, but it allows us to give a more complete picture, and also a first intuition of how boundaries introduce non-vanishing Noether charges for diffeomorphisms. If the spacetime has isometries, the Noether current (26) gives rise to a useful conservation law analogue to the Gauss law in electromagnism, which we recall states that the total charge in a region is equal to the flux of the electric field. To see this, we first observe that

[TABLE]

an identity which follows from the definition of the Riemann tensor as the commutator of two covariant derivatives. If $\xi^{\nu}$ is a Killing vector, the second term gives $-R^{\mu}{}_{\nu}\xi^{\nu}$ and the last term vanishes. Then integrating both sides of the equation over a 3d portion of space $V$ delimited by two boundaries $S_{1}$ and $S_{2}$, and using Stokes’ theorem, we find

[TABLE]

The Noether charge (28) obtained in this way is known as Komar charge. If the right-hand side of (29) vanishes, the Komar charge is conserved in the sense that it has the same value no matter which surface $S$ is used, and its value changes only when the deformations of $S$ include some source terms. If the right-hand side does not vanish, the Noether charge varies by an amount determined by the total quantity of energy-momentum in the enclosed region, see Fig 3. As an example, one can consider the Kerr solution, which possesses two Killing vectors corresponding to stationarity and axial symmetry. Evaluating (28) on an arbitrary 2-sphere $S$ encompassing the singularity gives respectively the mass and angular momentum (up to numerical coefficients to be fixed), independently of the coordinate used and independently of deformations of $S$.

While Komar charges are limited to isometries, it is possible to generalize the construction of Noether charges and canonical generators to arbitrary spacetimes, at least in so far as they admit boundaries with non-trivial residual diffeomorphisms.777There is also ongoing research on constructing charges on arbitrary regions in arbitrary spacetimes. Let us mention three important examples. First, spacetimes that are asymptotically flat at spatial infinity. The residual diffeomorphisms compatible with the boundary conditions are the Poincaré transformations of the flat boundary metric. One can construct Noether charges and canonical generators for these boundary diffeomorphisms (see e.g. Iyer:1994ys ), and the result coincides with the Arnowit-Deser-Misner (ADM) charges that were previously derived with canonical methods. Second, spacetimes that are asymptotically flat a null infinity. This case is particularly relevant to understand gravitational waves at the non-perturbative level. The residual transformations are a generalization of Poincaré transformations in which translations are angle-dependent, an infinite-dimensional extension known as Bondi-Van der Burg-Metzner-Sachs (BMS) transformations.888Intuitively, the extension comes about because the induced metric on a null hypersurface is degenerate, hence any deformation along that direction leaves the system invariant. Noether charges for the BMS symmetry were constructed in Ashtekar:1981bq ; Dray:1984rfa ; Wald:1999wa , and there is a large body of recent literature on the subject motivated by ongoing applications and developments. Among these, the application of the Noether approach to horizons and more generally null boundaries, see e.g. CFP ; AKKL .

II.3 Gauge and observers

Even though every physical phenomenon can be described in any coordinate system, some choices can stand out because they simplify the description, or because they are naturally associated with a class of observers of interest. For instance in flat spacetime, Cartesian coordinates make the Christoffel’s symbols vanish, thus simplifying many calculations, and can be associated to inertial observers, in the sense that they label their rods and clocks. In curved spacetimes the situation is more complicated because there are no Cartesian coordinates. To take a simple example, consider the Schwarzschild black hole. The most common coordinates used to describe this solution are the so-called ‘static ones’, which make time-independence manifest since the metric is independent of the time coordinate $t$. Or in better terms, since the time-translation Killing vector is simply $\partial_{t}$. The vector field $\partial_{t}$ describes a family of non-inertial observers static at a fixed distance outside the black hole, whose time delays are related by $\sqrt{-g_{tt}}$. The Schwarzschild static coordinates are naturally associated to these observers. An alternative choice are Gullstrand-Painlevé coordinates, whose $t$ describes the proper time of observers radially free-falling into the black hole. Another choice is the ‘temporal gauge’ defined by

[TABLE]

in which the radially free-falling observes are all synchronized,999This requires giving them non-zero energy, as opposed to the Gullstrand-Painlevé observers that have zero energy. whence the alternative name of ‘synchronous gauge’. For Schwarzschild, this gauge can be achieved using Lemaitre coordinates. The synchronization may look like a nice feature, but in these coordinates the spatial part of the Schwarzschild metric is explicitly dependent on the coordinate $t$! Therefore its staticity is hidden, and has to be verified by the existence of a time-translational Killing vector. In this example we can see the analogy between coordinate choices and gauge choices in electromagnetism very clearly. An electrostatic potential is more conveniently described in the Coulomb gauge because it makes the potential manifestly time-independent. But it can be described in any other gauge, and if we use the temporal gauge $A_{0}=0$, the potential acquires an inconvenient time dependence which is pure gauge and hides the staticity of the system.

The temporal gauge (30) can be chosen for any spacetime in a given coordinate chart. If it is done, it fixes completely the 4-dimensional diffeomorphism freedom in the chart. It is thus an example of complete gauge fixing,101010Complete here refers to the 4-dimensional picture. There remains the freedom of time-independent 3-dimensional diffeomorphisms, namely of choosing the coordinates on one – and one only – given hypersurface. and the resulting coordinates describe free-falling observers with synchronised clocks. It means that in the temporal gauge the coordinates are attached to free-falling test bodies, and this can be a very useful choice in many situations.

One shortcoming of the temporal gauge is that it may ‘squeeze’ physical information in field components one would not naturally expect, as the Schwarzschild example in Lemaitre coordinates shows. Another one is that it is non-symmetric among the components, in the sense that it relies on an initial choice of coordinate to be taken as the time. An example of symmetric gauge is the harmonic gauge, which is satisfied by coordinates such that

[TABLE]

where $\square=\nabla_{\nu}\nabla^{\nu}$ is the curved spacetime d’Alembertian. The symmetry makes this a convenient choice in many dynamical situations. For instance, it is the gauge in which it is easiest to see that the initial value problem is well-posed ChoquetBruhat1952 , and the one in which it is easiest to study gravitational waves, where it gives rise to a gauge preserving the Lorentz covariance of the Minkowski background. Notice that (31) is not a complete gauge fixing of the 4-dimensional diffeomorphism freedom, as there are infinitely many solution of the wave equation, hence infinitely many choices of harmonic coordinates for a given metric. To obtain a complete gauge fixing one has to specify a unique set of harmonic coordinates with additional conditions.

One version of the equivalence principle states that purely local experiments can not distinguish the presence of gravity. In the formalism of general relativity, this is embodied in the fact that at any given point in spacetime it is possible to find coordinates so that the metric is flat and its first derivative vanish. A coordinate system that achieves this is called a local inertial frame. Only an experiment that can probe second-order variations in the metric would be able to see the effect of gravity in a coordinate-independent way, and these variations are the tidal forces that show up in the geodesic deviation equation. An example of local inertial frame is provided by Riemann normal coordinates, which are constructed around a given point so that the Taylor expansion of the metric components around that point taken as the origin gives

[TABLE]

where $c_{\mu\nu}=-1/3$ $\forall\mu,\nu$.111111This formula is only valid in certain coordinate choices, hence the non-covariant notation with the $\mu\nu$ indices repeated but not summed over. This gauge fixing specifies coordinates only in the neighbourhood of a point and not in a full coordinate chart, but it is still perfectly sufficient to describe the physics of a local inertial frame around that point. It is also possible to find coordinates such that the Christoffel’s symbols vanish everywhere along a chosen time-like geodesic. Such coordinates are known as Fermi normal coordinates, and describe a free-falling local inertial frame. The metric takes the same form (32), but with $c_{00}=-1,c_{0a}=-2/3,c_{ab}=-1/3$. The difference between temporal gauge and Fermi normal coordinates is that in the first case every observer at constant spatial coordinates is free falling, whereas in the second case only the observer going through the origin along the $\partial_{t}$ geodesic is free falling.

III Perturbative treatment of Einstein’s equations

III.1 The idea: general background spacetime

Perturbation theory can be set up choosing a background metric $\bar{g}_{\mu\nu}$ and writing

[TABLE]

with the assumption that $|h_{\mu\nu}|\ll 1$ in some chosen coordinate system, so that it can be treated as a perturbation. One can then systematically Taylor-expand all metric functionals around the background, starting with the inverse metric $g^{\mu\nu}=\bar{g}^{\mu\nu}-h^{\mu\nu}+O(h^{2})$ and Levi-Civita connection

[TABLE]

and attempt to solve the field equations order by order:

[TABLE]

Here $\bar{G}_{\mu\nu}=G_{\mu\nu}(\bar{g})$, $\bar{T}=T(\bar{g})$, $G^{\scriptscriptstyle\rm(1)}_{\mu\nu}=G_{\mu\nu}(\bar{g};h)$, and so on. Explicit expressions will not be needed here, but are given in Appendix A for completeness. The lowest order of the procedure is straightforward, one just has to be consistent with the treatment of the matter energy-momentum tensor: if we want to expand around a background which is a vacuum solution, we need $\bar{G}_{\mu\nu}=0$ hence we should treat matter as a first order perturbation. Then the first order equation to solve is

[TABLE]

Solving this equation determines $h_{\mu\nu}$ in terms of its independent degrees of freedom, the background solution and the matter content. To go to second order, we add a second perturbation, writing

[TABLE]

Then $G^{\scriptscriptstyle\rm(2)}_{\mu\nu}=G^{\scriptscriptstyle\rm(2)}_{\mu\nu}(h)+G^{\scriptscriptstyle\rm(1)}_{\mu\nu}(h^{\scriptscriptstyle\rm(2)})$ has two contributions, and we solve for $h^{\scriptscriptstyle\rm(2)}$ using

[TABLE]

Notice that the first order solution feeds back as a source for the second order solution, a standard procedure from perturbatively solving non-linear equations. The same procedure applies also to the case when matter fields contribute to the background as well. Although we wrote the perturbative expansion using dimensionless quantities only, it is always possible to rescale the metric perturbation by $\sqrt{G}$, so that the kinetic term of the free-field Lagrangian is canonically normalized. Then the expansion is a power series in $G$, and it is also referred to as post-Minkowskian (PM) expansion, when the background is flat. In many cases the equations are still to hard to solve at each order, and one needs to look for additional approximations. A very common one is the non-relativistic approximation, also known as post-Newtonian (PN) expansion, in which one starts from a source that moves at small velocity with respect to the background flat metric, namely with $v\ll c$, and expands in powers of $v^{2}/c^{2}$. Alternative approximation schemes include the Bondi asymptotic expansion which is perturbative in the inverse distance from the source but valid at all orders in $G$ and $v^{2}/c^{2}$ Bondi:1962px ; Sachs:1962wk ; the extremal mass ratio inspiral (EMRI) and self-force expansion for a two-body system where the small parameter is the mass-ratio Barack:2018yvs ; the effective one-body approach based on resummed PN results Buonanno:1998gg , and more recently the effective field theory approach based on tools from quantum field theory Goldberger:2004jt .

On top of the technical difficulties, there is also a conceptual one: truncating at any order beyond the lowest breaks covariance. Expanding both sides of (21) with (33), we obtain

[TABLE]

We assume that the background metric is fixed once and for all and unaffected by diffeomorphisms, namely $\delta\bar{g}_{\mu\nu}=0$ and $\delta_{\xi}\bar{g}_{\mu\nu}=0$,121212It is also possible to interpret (39) as $\delta_{\xi}\bar{g}_{\mu\nu}=\pounds_{\xi}\bar{g}_{\mu\nu}$ and $\delta_{\xi}h_{\mu\nu}=\pounds_{\xi}h_{\mu\nu}$. The proof that this is a symmetry is identical. However this definition of perturbed transformations is less interesting physically, because it is more natural to compare perturbations when they are defined with respect to the same background. and define the effect of the linearized diffeomorphism on the perturbation to be

[TABLE]

The key point is that the right-hand side contains terms of different order in $h_{\mu\nu}$. This is why different orders of the perturbative expansion must be included in order for diffeomorphisms to be a symmetry, and conversely truncating at a fixed order breaks covariance. It is instructive to prove this in detail, because it highlights the special features that occur for the quadratic Lagrangian, namely the free theory, for which covariance is on the other hand possible.

To study the symmetries of the perturbative expansion, we write the perturbed Lagrangian as follows

[TABLE]

where barred quantities only depend on $\bar{g}$ and not on the perturbation. Using (40), the infinitesimal variation gives

[TABLE]

All terms can be collected into Lie derivatives:

[TABLE]

and so on. Each Lie derivative gives a boundary term through (23), hence (40) is indeed a symmetry of the full Lagrangian. However, in every case except the lowest one, getting a boundary term for $\bar{\cal L}^{\scriptscriptstyle\rm(n)}$ requires both $\bar{\cal L}^{\scriptscriptstyle\rm(n)}$ and $\bar{\cal L}^{\scriptscriptstyle\rm(n+1)}$, hence if we truncate the series at a fixed order, we lose covariance.

There are two special features that occur at the quadratic order. First, if we take the background to be a solution, the first term in (44) vanishes, hence

[TABLE]

is a symmetry of the quadratic Lagrangian. Second, if the background has isometries, then the term proportional to $\bar{\cal L}^{\scriptscriptstyle\rm(3)}$ in (45) drops out, and then

[TABLE]

is a symmetry of the quadratic Lagrangian. This shows that isometries play a special role in perturbation theory. For a generic on-shell background the symmetry of the quadratic Lagrangian is (46). But if the background has isometries, we have two different realization of the diffeomorphism symmetry in the quadratic Lagrangian: (46) for a generic diffeomorphism, and (47) for a Killing vector. To make this consistent with the perturbative expansion, we treat a Killing $\xi$ as zero-th order, and a non-Killing $\xi$ as first order.

We stress that (46) and (47) are symmetries only for the quadratic Lagrangian. From the cubic Lagrangian onwards, there is no symmetry at fixed order in $h_{\mu\nu}$. The only symmetry is the combined (40) and requires two different perturbative orders of the Lagrangian. This fact has immediate consequences for the expanded Einstein tensor (35). Since the leading order $G^{\scriptscriptstyle\rm(1)}$ comes from the quadratic Lagrangian, it is invariant under (46). But $G^{\scriptscriptstyle\rm(2)}_{\mu\nu}$ which is derived from the cubic Lagrangian is not invariant. As a consequence, the quantity $t^{\scriptscriptstyle\rm G}_{\mu\nu}$ appearing in (38) is not gauge invariant. This means that it does not provide a meaningful notion of gravitational energy-momentum tensor. We will come back to this important point below.

Even though both symmetries (46) and (47) descend from the same diffeomorphism invariance of the theory, which is a gauge symmetry, they have a different status at the perturbative level. Applying Noether’s theorem to the quadratic Lagrangian, we find that the generic diffeomorphisms still have vanishing conserved current, and therefore maintain their status of gauge symmetries. However the diffeomorphisms corresponding to isometries of the background have non-vanishing Noether charges, indeed just like a standard theory on flat spacetime. A related difference is that the field equations are invariant under (46), hence there are linear dependencies in the equations and some field components are left undetermined, and covariant under (47), and there are no constraints associated with them.

III.2 Weak-field approximation

The weak-field approximation is a special case of the perturbative expansion in which the background spacetime is Minkowski,

[TABLE]

For instance in the solar system spacetime is approximately flat, hence it can be described using this approximation.131313In first approximation, we can describe the solar system’s metric using the Schwarzschild solution with $M$ the mass of the sun. Then the curvature scale can be estimated writing the Kretschmann scalar in static coordinates, and this gives $R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma}=48G^{2}M^{2}c^{-4}r^{-6}\sim 10^{-27}{\tt m}^{-4}(R_{\odot}/r)^{6}$ which is extremely small. We further use Cartesian coordinates, so all covariant derivatives become partial derivatives. In this case the linearized Einstein equations (36) take the simple form

[TABLE]

where $\square$ is the d’Alembertian in flat spacetime.

Let us study the symmetries of (49). The Minkowski background has ten isometries, the Poincaré transformations. We parametrize them with vector fields

[TABLE]

where $a^{\mu}{}_{\nu}$ and $b^{\mu}$ are constants, and $a_{(\mu\nu)}=0$. According to the general discussion of the previous section, we expect two different types of symmetries in (49), both induced from diffeomorphism invariance of the full theory, and corresponding to those linearized diffeomorphisms that are isometries or not of the background. For diffeomorphisms corresponding to isometries, (47) gives

[TABLE]

This is the transformation of a rank-2 tensor in Minkowski under Poincaré transformations. For generic diffeomorphisms that are not isometries, (46) gives

[TABLE]

Both transformations are symmetries of (49), the difference being that the first changes the equations but in a covariant way, whereas the second leaves them invariant. This can be checked easily on the left-hand side. For the right-hand side, one has to pay attention to the behaviour of the energy-momentum tensor. Assuming the background $\eta$ to be a solution requires that we are treating the matter fields as first order in perturbation theory. Then it transforms as $\pounds_{\xi}T_{\mu\nu}$ for a Killing $\xi$, whereas its transformation under generic diffeomorphisms is second order, hence it does not affect the linearized equations.

The presence of the d’Alambertian in the linearized field equations suggests that wave solutions are indeed possible. However, there is an intricate tensorial structure that needs to be dealt with. Before doing so, let us discuss what it means for a wave to be ‘tensorial’. The waves that we are most familiar with, such as water waves or sound waves, are scalar waves, namely the quantity whose perturbation propagates following the wave equation is a scalar function, such as the height of the water or the pressure of the air. Electromagnetic waves on the other hand propagate changes in the electric and magnetic field which are described by vectors, and are thus ‘vectorial’ waves. The difference between scalar, vectorial and tensorial waves can be described in terms of spin. The reader already familiar with these concepts, or the reader interested in looking first at the explicit solutions without these details, can skip the next subsection.

III.3 Spin and helicity

Tensors in Minkowski spacetime belong to finite-dimensional representations of the Lorentz group. This property can be used to decompose each tensor into irreducible parts, namely parts that are not mixed with one another by a Lorentz transformation. The irreducible parts can be labelled by a pair of half-integers $(j_{1},j_{2})$, and contain $(2j_{1}+1)(2j_{2}+1)$ components. For instance, a 4-vector $v^{\mu}$ transforms under the irreducible Lorentz representation $\bf(\frac{1}{2},\frac{1}{2})$ with 4 components, and a symmetric tensor $h_{\mu\nu}$ transforms under the reducible Lorentz representation ${\bf(1,1)\oplus(0,0)}$ with $9+1=10$ components.

The representations can be further subdivided if we pick a time direction $\tau^{\mu}$, and restrict attention to the rotation subgroup of the Lorentz group that preserves it. The subsets of the tensor which are irreducible with respect to the rotation subgroup are called spin representations, and their allowed values given by the Clebsch-Gordan addition rule $(j_{1}+j_{2},j_{1}+j_{2}-1,\dots,|j_{1}-j_{2}|)$. For instance a vector contains the two spin representations $\bf 1\oplus 0$, and a symmetric tensor the four spin representations $\bf 2\oplus 1\oplus 0\oplus 0$. To be more explicit let us take $\tau^{\mu}=(1,0,0,0)$. The rotation subgroup that preserves it has the form

[TABLE]

It is then immediate to see that the spin-1 and spin-0 representation of the 4-vector are the spatial vector $v^{a}$ and the spatial scalar $v^{0}$ respectively. Similarly for a 1-form $v_{\mu}$, the spatial and time components $v_{a}$ and $v_{0}$. For a symmetric tensor, the four spin representations are

[TABLE]

The spin 0 representation has one component, the spin 1 has three, and the spin 2 has five. These different components can be classified choosing a reference spatial axis and looking at the eigenmodes of the rotation generator along that axis. To fix ideas let us choose the $z$ axis. The rotation matrix that preserves it is

[TABLE]

Inserting this in (53) and acting on a 4-vector $v^{\mu}$ we obtain $v^{\prime}{}^{\mu}=R^{\mu}{}_{\nu}v^{\nu}$, where

[TABLE]

The temporal and the $z$ component, which is longitudinal along the rotation axis, are invariant. The components $x$ and $y$, which are transverse to the rotation axis, transform among themselves, and we can diagonalize the action introducing

[TABLE]

In this way we have identified all 4 eigenmodes of the vector. We can equivalently describe the decomposition in terms of a basis of eigenvectors. If we start from the canonical basis $e^{\mu}_{I}:=\delta^{\mu}_{I}$, where $I=0,\ldots 3$, then the eigenvectors are

[TABLE]

and their corresponding eigenvalues

[TABLE]

The integers 0 and $\pm 1$ that appear in front of the rotation angle $\theta$ are the projection of the spin along the $z$ axis. Since our convention for vector rotations is anti-clockwise, the $+1$ (resp. $-1$) component can be also referred to as right-handed (resp. left-handed).

Acting on a symmetric tensor $h_{\mu\nu}$ we obtain $h^{\prime}_{\mu\nu}=R^{\rho}{}_{\mu}R^{\sigma}{}_{\nu}h_{\rho\sigma}$, where

[TABLE]

We have thus identified all ten eigenmodes of a symmetric rank-2 tensor under the $z$ rotations.

To write the eigenvectors, we introduce a canonical basis in the space of symmetric $4\times 4$ matrices,

[TABLE]

Then

[TABLE]

with eigenvalues

[TABLE]

The spin-2 components $e^{{\scriptscriptstyle\rm T}\pm}$, $e^{{\scriptscriptstyle\rm L}\pm}$ and $e^{\scriptscriptstyle\rm LL}$ carry respectively the $\pm 2,\pm 1$ and [math] modes of the $z$-projection. The spin-1 components $w^{\pm}$ and $w^{\scriptscriptstyle\rm L}$ carry $\pm 1$ and [math] modes as in (59), and $e^{\scriptscriptstyle\rm 0,s}$ are the remaining two spin-0 modes from (54).

The discussion so far concerned global vectors and tensors. In the case of electromagnetism and gravity linearized around Minkowski we work with vector and tensor fields, $A_{\mu}(x)$ and $h_{\mu\nu}(x)$. The mode decompositions described above are then applied locally in Fourier space. Each Fourier mode $\tilde{A}_{\mu}(k)$ and $\tilde{h}_{\mu\nu}(k)$ is characterized by a wave vector $k$, and we can use its spatial direction to classify the spin components. With some abuse of language, we will often refer to the wave vector as the momentum of the wave. The eigenvalues of the projection of the spin along the momentum are called helicities, hence when this basis is chosen the mode decomposition is called helicity decomposition. The notion of helicity is closely related to the notion of polarization. More precisely, modes of helicity $\pm 1$ describe waves of right-handed and left-handed circular polarization respectively, and different linear combinations can be taken to describe for instance linear or elliptic polarizations. The helicity-1 linear polarizations are simply the $x$ and $y$ components, and the helicity-2 linear polarizations are

[TABLE]

The labels stand for ‘plus’ (not to be confused with the plus used in the circular polarizations, and which will not be used in the following) and ‘cross’, and the reason will become clear below.

The helicity decomposition with an arbitrary momentum $\vec{k}$ can be conveniently described introducing the transverse and longitudinal projectors

[TABLE]

Using these, we can decompose a spin-1 vector as follows,

[TABLE]

where

[TABLE]

and a spin-2 tensor as follows,

[TABLE]

where

[TABLE]

We can then identify the helicities of the different projectors studying how they transform under a rotation with axis $\vec{k}$. For simplicity let us consider the case when $\vec{k}$ is in the $z$ direction, so that we can use the formulas already derived for the eigenvectors. In this case,

[TABLE]

Then (67) reduces to $A_{a}^{\scriptscriptstyle\rm T}=(A_{x},A_{y},0)$ and $A_{a}^{\scriptscriptstyle\rm T}=(0,0,A_{z})$. Comparing these to the earlier decomposition (59), we conclude that the transverse projector contains the $\pm 1$ helicity modes of a spin-1 field, and the longitudinal projector the [math] helicity mode. The spin-2 projectors based on (71) give

[TABLE]

Comparing these to the earlier decomposition (62), we conclude that the spin-2 part of the gravitational perturbation can be decomposed into 5 helicity modes $\pm 2,\pm 1,0$, which are carried respectively by the TT, L and LL components.

The spin-helicity interpretation of the components of $h_{\mu\nu}$ is a kinematical classification based on a choice of reference frame given by the time direction141414Which we have chosen to coincide with $t$ of the background Cartesian coordinates. However the whole spin-helicity description is Lorentz covariant. The spin projectors for an arbitrary time-direction $\tau^{\mu}$, $\tau^{2}=-1$, are

${}P^{\scriptscriptstyle\rm{\bf(1)}}{}^{\mu}_{\nu}=q^{\mu}_{\nu}:=\delta^{\mu}_{\nu}+\tau^{\mu}\tau_{\nu},\qquad P^{\scriptscriptstyle\rm{\bf(0)}}{}^{\mu}_{\nu}=-\tau^{\mu}\tau_{\nu},\qquad 1^{\bf(\frac{1}{2},\frac{1}{2})}=P^{\scriptscriptstyle\rm{\bf(1)}}+P^{\scriptscriptstyle\rm{\bf(0)}},$

for a vector, and

$\displaystyle P^{\scriptscriptstyle\rm{\bf(2)}}{}^{\mu\nu}_{\rho\sigma}=\delta^{\mu}_{(\rho}\delta^{\nu}_{\sigma)}+\tau^{\mu}\tau_{(\rho}\delta^{\nu}_{\sigma)}-\frac{1}{2}q^{\mu\nu}q_{\rho\sigma},\qquad P^{\scriptscriptstyle\rm{\bf(1)}}{}^{\mu\nu}_{\rho\sigma}=-\tau^{\mu}\tau_{(\rho}q^{\nu}_{\sigma)},$

$\displaystyle P^{\scriptscriptstyle\rm{\bf(0)}}{}^{\mu\nu}_{\rho\sigma}=\tau^{\mu}\tau_{\rho}\tau^{\nu}\tau_{\sigma},\qquad P_{s}^{\scriptscriptstyle\rm{\bf(0)}}{}^{\mu\nu}_{\rho\sigma}=\frac{1}{2}q^{\mu\nu}q_{\rho\sigma},$

(73)

$\displaystyle 1^{{\bf(1,1)\oplus(0,0)}}=P^{\scriptscriptstyle\rm{\bf(2)}}+P^{\scriptscriptstyle\rm{\bf(1)}}+P^{\scriptscriptstyle\rm{\bf(0)}}+P^{\scriptscriptstyle\rm{\bf(0)}}_{s},{}$

for a symmetric tensor. One can also write the helicity projectors as 4d covariant objects encoding $\tau$ in a choice of null vector transverse to $k^{\mu}$, see e.g. Weinberg . and the spatial direction $\vec{k}$. The next question is which of these components are dynamical. If the field satisfied the simple wave equation $\square h_{\mu\nu}=0$, then all components would be dynamical, and the field would carry 10 independent degrees of freedom corresponding to all the helicity states described above. But this is not the case of (49), because of the additional derivative operators present, and the gauge redundancy. One way to identify the degrees of freedom is to fix the gauge and study the resulting solutions. This is what we do next.

III.4 De Donder and TT gauges

The gauge symmetry can be exploited to simplify the linearized field equations. This can be done using (52) to put the metric perturbation in a form where it satisfies

[TABLE]

Doing so eliminates all terms containing divergences. The linearized field equations (49) are thus equivalent to

[TABLE]

The short-hand notation $\bar{h}_{\mu\nu}$ is introduced for convenience, and it is known as trace-reversed perturbation, since $\bar{h}=-h$. The condition (74) is called De Donder gauge in the literature, but also Lorenz gauge, in analogy with electromagnetism, and harmonic gauge, because it preserves harmonic coordinates, as can be seen linearizing (31). The compatibility of the coupled system is guaranteed by the fact that $\partial_{\mu}T^{\mu\nu}$ vanishes on solutions of the matter field equations, more on this in Section V.1 below.

The condition (74) does not fix completely the 4-dimensional diffeomorphism symmetry. Intuitively, this occurs because we are not fixing metric coefficients, but only their derivatives. To make the residual freedom explicit, observe that the diffeomorphism required to put an arbitrary metric perturbation in De Donder form is given by a solution of the equation

[TABLE]

But for a given initial metric, this equation admits infinitely many solutions, parametrized by the zero modes of the d’Alembertian operator. In other words, once the De Donder condition is satisfied, there remains a residual freedom of gauge transformations that satisfy $\square\xi^{\mu}=0$. Thus the De Donder gauge contains in fact an infinite family of inequivalent gauge fixings. For this reason, it may be more appropriate to refer to (74) not as a gauge fixing condition, but rather a family of gauge fixings. The simpler albeit vaguer term De Donder gauge is, however, the one in use in most of the literature.

If need be, a unique representative of the De Donder gauge can be specified fixing the solution of $\square\xi^{\mu}=0$ in terms of initial data at a reference $t=0$ hypersurface. Such initial data can always be chosen such that any four components of $h_{\mu\nu}$ and their time derivatives vanish there. For vacuum solutions, this means that the chosen components vanish everywhere. It is convenient to apply this procedure and set to zero the components $\bar{h}$, so that $\bar{h}_{\mu\nu}=h_{\mu\nu}$, and $h_{0a}$. Then the De Donder condition implies that also $\partial_{a}h^{ab}$ and $\partial_{0}h_{00}$ vanish everywhere. The vacuum equations then imply $\vec{\partial}^{2}h_{00}=0$, hence the only solution with vanishing boundary conditions is $h_{00}=0$. We conclude that in this gauge, vacuum solutions satisfy

[TABLE]

This shows explicitly that there are only two independent degrees of freedom. We call this choice the transverse-traceless gauge (TT gauge in short), since in this gauge, vacuum solutions coincide with the transverse-traceless perturbations. We stress that this property and the equations (77) are only valid for vacuum solutions. In other words, (77) is not the definition of a gauge condition, but rather the specific value that solutions take in a certain gauge. The analogue for solutions with sources will be discussed below. Notice also that $h_{0\mu}=0$ is the linearized approximation of (30). Therefore the TT gauge of vacuum solutions implies the temporal gauge, hence the TT coordinates describe free-falling observers.

In the presence of sources on the other hand, even if we set to zero 4 components and their time derivatives on the chosen initial hypersurface, they will no longer be zero in the future causal domain of the sources. In this case there are more convenient choices to specialize the De Donder gauge condition, as we will discuss below.

III.5 Vacuum solutions

Let us flesh out these considerations by looking at the explicit form of the vacuum solutions. This will make it clear that the identification of the degrees of freedom with the transverse-traceless modes is a gauge-invariant statement. The vacuum equations can be solved straightforwardly taking linear combination of plane waves via the Fourier transform

[TABLE]

Imposing the vacuum equations and the De Donder condition (75) requires

[TABLE]

The first equation is solved by $k^{0}=\pm|\vec{k}\,|$, and we denote $\omega:=k^{0}c>0$ for a future-pointing 4-momentum, thus

[TABLE]

The second equation gives a linear system of 4 conditions on the 10 components of the matrix. In the frame where $\vec{k}$ is along the $z$ axis, we take as independent components

[TABLE]

and then the system is solved by

[TABLE]

Therefore the general solution for a real, monochromatic wave with frequency $\omega$ and propagating along the $z$ axis is a 6-parameter family given by

[TABLE]

The solution, however, contains gauge redundancy. Recall in fact that there is residual gauge freedom in the form of diffeomorphisms $\xi^{\mu}$ that satisfy the vacuum wave equation. In Fourier space the gauge transformation (52) read $\delta_{\xi}\tilde{h}_{\mu\nu}=2ik_{(\mu}\tilde{\xi}_{\nu)}$, where $\tilde{\xi}^{\mu}(k)$ is the Fourier transform of $\xi^{\mu}(x)$ with $k^{2}=0$ in order to satisfy the vacuum wave equation and be an admissible residual gauge. Under the residual gauge transformation,

[TABLE]

The second property follows immediately from the fact that $k_{\mu}$ has no transverse components. We conclude that only the two components $h_{\scriptscriptstyle\rm+,\times}$ are gauge-invariant and thus physically relevant. These can be recognized as the TT components. The four coefficients $\tilde{h}_{0\mu}$ can be changed arbitrarily, and have no physical meaning. We can in particular set them to zero simultaneously, choosing $i\tilde{\xi}_{\mu}=\frac{c}{\omega}(\frac{1}{2}\tilde{h}_{00},\tilde{h}_{0x},\tilde{h}_{0y},\tilde{h}_{0z}-\frac{1}{2}\tilde{h}_{00})$. Then

[TABLE]

We have thus achieved the TT form (77). In this gauge, the only components of the vacuum solutions are the gauge-invariant transverse-traceless ones. Accordingly, the solution can also be written as a sum over the two TT polarization tensors (65), namely

[TABLE]

where $p=+,\times$ and $A_{p}=h_{p}$.151515The polarization decomposition can also be used to solve the wave equation. In this approach one writes (86) with the sum including all ten polarization modes (62) and ten arbitrary coefficients $A_{i}$. Then the relations (82) between Cartesian components are replaced by relations between the polarization coefficients. The algebra is slightly more involved but the end result is the same. We will see the use of the the polarization mode procedure on cosmological backgrounds in Sec VIII.

The resulting perturbed metric is

[TABLE]

Let us pause for a short historical digression, about some of the confusion that hindered historically the understanding of gravitational waves. Let us consider any member of the general solution (83) with $h_{\scriptscriptstyle\rm+}=h_{\scriptscriptstyle\rm\times}=0$. This looks like a genuine wave, and it is a solution of the linearized Einstein’s equations. However, it is a pure gauge solution, and can be set to vanish identically without loss of physical information. In other words, there are wave solutions which are in the end only coordinate artefacts. For the same reason, it is also possible to change coordinates so that the argument $ct-z$ in its cosine is replaced by $vt-z$ for an arbitrary constant $v$. Hence the pure gauge modes don’t really propagate, and if a gauge is chosen so that they look like they are propagating, well one can do this with an arbitrary speed, the speed is not constrained in any way by the dynamics. To use Eddington’s words, the non-physical gauge modes propagate at the “speed of thought”. This initial confusion was clarified by the identification of gauge-invariant components, and the fact that for the physical modes, the propagation speed is fixed to be the speed of light by Einstein’s equations. However, additional doubts persisted, because the metric (87) is not a solution of the exact Einstein’s equations; only of the linearized theory. In other words, there are vacuum solutions that take approximately the form (87) in some regions of spacetime, but none that has that exact form everywhere in spacetime. This raised the issue of whether gravitational waves existed in the full theory, or were only an artefact of the linearized approximation. This more complicated issue was solved only much later, when a non-perturbative identification of the wave degrees of freedom and their energy was made possible by the work of Bondi, Sachs and many others Bondi:1960jsa ; Bondi:1962px ; Sachs:1962wk ; Newman:1968uj ; Ashtekar:1981bq ; Dray:1984rfa ; Wald:1999wa . For a historical overview, see Kennefick:1997kb .

Coming back to the physical solution (85), we see that the resulting wave tensor is transverse to the direction of propagation. Comparing with (62), we conclude that the physical gauge-invariant modes have helicity $\pm 2$. All modes of helicity $\pm 1$ and [math] have dropped out, either by gauge fixing, or by solving the vacuum field equations. Because the physical modes have helicity $\pm 2$, it takes a rotation of an angle $\pi/4$ to turn one mode into the other. This can be compared with the electromagnetic case, whose modes have helicity $\pm 1$, and it takes a rotation by $\pi/2$ to turn one mode into the other, see Fig. 5.161616The general formula is that two modes of helicity $\pm j$ are related by $\pi/2j$. For the electromagnetic case, this angle $\pi/2$ corresponds to the orthogonality of the oscillations of the electric and magnetic fields. For the gravitational case, it corresponds to the antipodal symmetry in tidal forces. We will see below in Section IV how this intuition can be made precise by studying the effect of a gravitational wave on test particles, see in particular Fig. 5.

The solution (85) describes a monochromatic wave. The most general solution has $\infty^{3}$ Fourier components (one per choice of $\vec{k}$), and two independent degrees of freedom per component (the values of $h_{\scriptscriptstyle\rm+,\times}$), therefore it is described by $2\times\infty^{3}$ arbitrary numbers. These are the independent degrees of freedom of gravitational waves. It is the same number of the full theory, so the linearized approximation simplifies the dynamics but preserves the number of independent variations of the gravitational field that can occur in the full theory. Notice that the number of independent degrees of freedom is the same of electromagnetism, or of two scalar fields. Apart from the dynamical behaviour, what changes is also the behaviour of these degrees of freedom under Lorentz transformations, because of their different spins.

The analysis has also shown that only the transverse-traceless modes are gauge invariant. We did this using the partial gauge fixing provided by the De Donder condition, but it can be proved in full generality starting from the projector $P^{\scriptscriptstyle\rm{\bf(2)}TT}$ on transverse-traceless modes defined in (70), and observing that it annihilates gauge transformations. Since we are going to use this projector often, we drop the label ${\bf(2)}$ from now on. We also introduce the notation $\hat{k}=\vec{k}/k^{0}$, using which we can write the explicit form of the projector as

[TABLE]

Note that the symmetrization on the indices here and in (70) is omitted in some books poisson ; Maggiore1 , under the premises that one is applying it to symmetric tensors only anyway. In conclusion, we do not need to use the TT gauge, nor even the De Donder gauge, in order to identify the independent degrees of freedom. Whatever gauge we are using, we can always extract them via

[TABLE]

and the result of this projection is gauge invariant. If we align the frame so that $\hat{k}$ coincides with the $z$ axis, it is given by (72g). To treat a general direction, we parametrize it using polar coordinates $(\theta,\varphi)$ on the sphere, and write

[TABLE]

We then pick a basis in the plane orthogonal to $\hat{k}$, given by

[TABLE]

We use this basis to distinguish the two polarizations, and define

[TABLE]

We obtain in this way

[TABLE]

which coincides with (11.46) of poisson . If the wave is travelling along the $x$ axis for example, then $(\theta,\phi)=(\frac{\pi}{2},0)$, and $h_{+}=\frac{1}{2}(h_{zz}-h_{yy})$, $h_{\times}=-h_{yz}$.

We can also use this formula to relate the polarizations in an arbitrary direction $\hat{k}$ to the polarizations in the $\hat{z}$ direction. In this case $h_{xx}=-h_{yy}=h_{+}(\hat{z})$ and $h_{xy}=h_{\times}(\hat{z})$, while all other components vanish. The expressions then simplify to

[TABLE]

In physical applications, the coordinate system may be set from characteristics of the system, for instance for binaries, adapted so that the orbital plane coincides with the $(x,y)$ plane, and $x$ is in the direction of the pericenter. Then our rotation identifies $\theta$ with the inclination $\iota$ of the orbital plane, and the $y$ axis is along the ‘line of nodes’, namely the intersection between the orbital plane and the plane of sight, pointing in the direction of the descendent node. On the other hand, astronomers often choose the convention to align the line of nodes along the $x$ axis, pointing in the direction of the ascendent node. This means that in our reference system, the ‘longitude of ascending node’ $\Omega$ is $-\pi/2$, and the ‘longitude of pericenter’ $\omega$ taken from the ascending node is $\varphi-\pi/2$. That is,

[TABLE]

This maps our basis (90-91) to (3.44-3.45) of poisson . See Fig. 4.

III.6 Gauge-invariant description: independent and constrained degrees of freedom

The analysis of vacuum solutions has allowed us to identify the independent degrees of freedom of the gravitational field. There are also dependent degrees of freedom, namely components of the field that are gauge invariant, but uniquely determined from the sources. These can be exposed looking at the constraints, namely the $00$ and $0a$ components of (49) which give rise to elliptic equations, as opposed to hyperbolic ones. The intricacies of the tensorial structure make it however convenient to do first a kinematical analysis of gauge invariance using the transverse and longitudinal projectors. To that end, we first observe that the projectors can be described also in configuration space without doing the Fourier transform. In this case (66) is replaced by

[TABLE]

These expressions are somewhat implicit because one needs to specify boundary conditions in order to have a well-defined inverse of the Laplace operator $\vec{\partial}^{2}$. Requiring the fields to vanish at spatial infinity makes this definition equivalent to the one in momentum space. Notice that the transverse-longitudinal decomposition in configuration space can be recognized as the Helmholtz decomposition of a 3d vector field into solenoidal and irrotational parts. This representation also makes it clear that the projectors are local operators in Fourier space, but non-local in spacetime. It has the important consequence that it is not possible to identify exactly the gauge-invariant modes with local observations only.

It is also convenient to dispose of the projectors for the components with non-maximal helicities, by introducing auxiliary fields with smaller spin. For instance, we denote $W_{a}:=h_{0a}$ the spin-1 part of the gravitational perturbation, and write its longitudinal part as the gradient of a scalar:

[TABLE]

Similarly for the spin-2 part, we can write the mixed transverse-longitudinal and fully longitudinal modes introducing a transverse vector $B^{a}$ and a scalar $B$,

[TABLE]

Let us summarize. With this new notation at hand, the gravitational perturbation around Minkowski can be decomposed into

[TABLE]

with

[TABLE]

and

[TABLE]

Next, we look at the behaviour of these different helicities under gauge transformations. We have used the Poincaré symmetry of the Minkowski background to organize the ten components of $h_{\mu\nu}$ in terms of spin and helicity. But Minkowski is not invariant under the general diffeomorphism symmetry (52), hence there is no reason to expect that this decomposition be gauge-invariant. Consider for instance the spin-0 part $h_{00}$: this is just a metric component, and manifestly not invariant under diffeomorphisms. Only the TT component is gauge-invariant, as can be seen explicitly replacing $h_{\mu\nu}$ in (100a) with its gauge transformation (52) and observing that it vanishes thanks to the transverse projection on both indices. This is the gauge-invariance of the two TT components already observed in (the second equation in) (82). The remaining components transform as follows:

[TABLE]

It is possible to combine them to find a maximal number of 6 gauge-invariant quantities, given by

[TABLE]

The numerical factors and powers of $c$ have been chosen for later convenience. We will also use a dot for the derivative with respect to the time coordinate $t=x^{0}/c$, e.g. $\partial_{0}W=c^{-1}\dot{W}$. Our conventions for the physical dimensions are summarized as follows:

[TABLE]

Inserting the parametrization (98) in the linearized Einstein’s equation (49) one finds that gauge dependent quantities drop out and only the gauge-independent ones remain. This allows us to decouple the tensorial equations into two sets, an hyperbolic one featuring the d’Alambertian operator alone, and an elliptic one featuring the Laplace operator alone:

[TABLE]

The sources on the right-hand side of these equations are the components of $T_{\mu\nu}$ projected in the same way as (98), namely

[TABLE]

Rewriting the linearized Einstein equations (49) in the equivalent form (102) makes their three-sided structure, mentioned in subsection II.1, manifest. Four equations were redundant and have dropped out, hence the field has four undetermined components, which can be assigned arbitrarily choosing a specific gauge. Four equations are elliptic, and describe constrained degrees of freedom, namely dynamical components of the gravitational field which are uniquely determined by the sources. Finally, two equations are hyperbolic, hence they contain free data, and describe how these independent degrees of freedom propagate and react to the sources. This analysis therefore establishes that the gravitational field has two independent degrees of freedom, which are carried by the two $h_{ab}^{\scriptscriptstyle\rm TT}$ components of the metric. These are gauge-invariant, and describe the $\pm 2$ helicities of a spin-2 wave.

As for the remaining gauge-invariant components, we have seen that they satisfy Poisson equations, hence these are degrees of freedom that are entirely determined by the sources. For this reason they are sometimes called ‘Coulombic’ degrees of freedom. Their meaning can be elucidated looking at the post-Newtonian expansion, in which sources are moving slowly with respect to the speed of light. To begin with, let us first consider perfectly static sources.

For static sources in a given frame, $T_{00}=c^{2}\rho$ is the only non-vanishing component on the right-hand side of Einstein’s equations. Then the second and third equations in (102b) imply171717Assuming trivial boundary conditions, see App B. $\Phi_{a}=0$ and $\Psi=\Phi$, whilst the first gives Newton’s equation, and we can identify $\Phi$ with Newton’s potential. This is also the way in which one fixes the coupling constant of the full Einstein’s equations in terms of $G$ and $c$. How about the other potentials? They are sourced by moving bodies, and their existence is a consequence of relativistic invariance, akin to the electromagnetic occurrence of the vector potential next to the Coulomb potential. They describe effects that collectively go under the name of ‘gravito-magnetism’. These include additional contributions to the precessions of equinoxes, light bending, and frame dragging or Lense-Thirring effect. Their effect can be studied looking at the Lagrangian for a massive test particle, which gives

[TABLE]

where $v^{a}=dx^{a}/dt$. The term in round bracket is the first post-Newtonian correction. As we will see below, $h^{\scriptscriptstyle\rm TT}\sim c^{-4}$, hence the last term there is higher order: dissipative effects for massive particles only appear at 2PN.

III.7 Gauge-fixing with sources

Solving the decoupled equations (102) determines the gauge-invariant quantities (101) in terms of the sources and the initial conditions. These solutions do not determine a metric. Doing so requires the additional step of specifying the coordinates to be used. Only after the coordinates are given, or in other words only after a gauge is chosen, the physical degrees of freedom can be described in terms of a metric tensor.

The gauge-invariant approach is conceptually satisfying because it identifies the physical degrees of freedom and decouples the equations, making them easier to solve in principle. However, it is very limited in applicability. Firstly, the decoupling and simple identification of gauge-invariant quantities occur only for very special backgrounds, such as flat spacetime or homogeneous and isotropic.181818On a general background, one can still build the analogue of the spin-helicity decomposition replacing the partial derivatives with covariant derivatives, although care is needed to invert the Laplacian and handle its non-commutativity with the covariant derivative. However, the metric now enters explicitly the decomposition of the energy-momentum tensor, hence the decoupling will be lost in general. The non-commutativity of covariant derivatives also hinders the identification of gauge-invariant quantities. Secondly, even when the background is Minkowski, a general identification of gauge-invariant quantities at a fixed order in perturbation theory is only possible at the linear level, as explained earlier.

To go beyond these limitations, it is easier to put to the side the gauge-invariant description, and work instead in a fixed gauge. In the gauge-fixed approach, one chooses coordinates that impose restrictions on the metric, and then solves for individual metric components in that gauge, like we did in Section III.4. Namely we do not solve (102), but the original system (49), coupled to additional equations fixing the gauge. The additional equations remove the problem that the field equations are redundant and do not determine all metric components.

A simple example of gauge fixing is the temporal gauge

[TABLE]

This is the linearized version of the non-perturbative temporal gauge described in Section II.3, and provides a complete gauge fixing of the 4-dimensional diffeomorphism symmetry. In this gauge the only non-trivial components of the metric are the spatial ones, and (101) reduces to

[TABLE]

All gauge-invariant potentials are encoded in the components of the spin-2 mode (98d). A related choice would be to replace the condition on $h_{00}$ with the trace condition $h=0$. In this case

[TABLE]

This time it is the four metric components $B$, $B_{a}$ and $h_{00}$ that are fixed in terms of the sources, while $h_{0a}$ and $h$ vanish.

Another simple gauge fixing is the gravitational equivalent of the Coulomb gauge, defined so that both spin-2 and spin-1 parts are purely transverse, namely

[TABLE]

In this gauge the non-trivial components of the metric are

[TABLE]

and their relation to the gauge invariant potential is given by

[TABLE]

These examples have the advantage that the gauge-invariant potentials are identified with individual components of the metric. On the other hand, the radiative gauge-invariant modes $h^{\scriptscriptstyle\rm TT}_{ab}$ are complicated functions of the metric components. When fixing the coordinate gauge, one should keep in mind that some choices may be better than others, as we discussed with the Schwarzschild example in Section II.3.191919If $GM/c^{2}\ll r$ in the static spherical coordinates, we can treat it as a perturbative solution with $h_{\mu\nu}=2GM/(c^{2}r)(\delta^{t}_{\mu}\delta^{t}_{\nu}+\delta^{r}_{\mu}\delta^{r}_{\nu})$. We then see that these coordinates correspond to the Coulomb gauge fixing. Changing coordinates so to have temporal gauge with $h_{00}=0$ would push the Newton potential in the $B$ component, and introduce a dependence on the $t$ coordinate, as mentioned earlier. The examples above are non-covariant with respect to the Lorentz symmetry of the background, because they make reference to a given time foliation, and treat time and space components differently. In the presence of radiation, it is best to use a covariant gauge, because it simplifies the analysis of the solutions. We can, in fact, remark that while the above examples simplify the relation between the potentials and individual metric components, the remaining field equations for the propagating degrees of freedom are complicated. Whereas with the covariant De Donder gauge, all field equations took the simpler form (75). However while the De Donder condition (74) can be also imposed in the presence of sources, we can no longer select a unique representative satisfying the TT condition (77) globally. If need be, a unique representative of the De Donder gauge in the presence of sources can be specified as follows. The general solution of (75) is

[TABLE]

where $G(x,x^{\prime})$ is the d’Alembertian’s Green function (see Appendix B), $\bar{h}_{\mu\nu}^{\circ}$ any solution of the homogeneous equation, and the gauge condition is maintained via $\partial_{\mu}T^{\mu\nu}=\partial_{\mu}\bar{h}^{\circ\mu\nu}=0$. Since both $\bar{h}_{\mu\nu}^{\circ}$ and the residual diffeomorphism parameters $\xi^{\mu}$ satisfy the vacuum wave equation, we can use the residual freedom to set to zero any four components of $\bar{h}^{\circ}_{\mu\nu}$. For instance, we can choose $\bar{h}^{\circ}=0$, so that $\bar{h}^{\circ}_{\mu\nu}=\bar{h}^{\circ}_{\mu\nu}$, as well as $\bar{h}^{\circ}_{0a}=0$. Then the De Donder condition implies that four more components of the homogeneous solution also vanish everywhere, specifically $\bar{h}^{\circ}_{00}$ and $\partial_{b}\bar{h}^{\circ ab}$. At this stage the gauge is completely fixed, and the only components left in the homogeneous solution are the gauge-invariant ones $\bar{h}^{\circ\scriptscriptstyle\rm TT}_{\mu\nu}$, which carry the independent degrees of freedom of gravitational waves. The general solution is

[TABLE]

It is a complete gauge fixing that singles out a unique element of the De Donder family, and reduces to (77) for vacuum solutions. For generic sources, all components of the metric perturbation are non-zero in this gauge.202020More precisely, in the region causally connected to the sources, since the retarded Green function vanishes outside the light cone. So if the sources are present at all times, the metric perturbation is non-zero everywhere, whereas if the sources are ‘turned on’ at some initial time, then the perturbation vanishes outside the causal domain of the sources from that initial time.

Alternatively, it is also possible to select a unique representative of the De Donder gauge requiring that the metric components $h_{00}=h_{0a}$ and their derivatives vanish on a given initial value surface, as considered in WaldBook . However these metric components will remain zero only in the region outside the causal domain of the sources from the initial value surface. This prescription thus achieves a TT gauge valid also in the presence of sources, but only in a finite region of spacetime. A similar construction of TT gauge in a local region of finite time outside the sources is presented in Flanagan:2005yc . While interesting in principle, the ‘local TT gauge’ seems to be of minor practical use. The procedure described in the previous paragraph achieves instead a structure of the solutions valid globally and at all times, and it is the analogue of the complete gauge-fixing used for instance in electro-magnetism by the Lienard-Wiechert potentials of a moving charge.

A special situation occurs if the sources are static. In this case, it is possible to specialize the De Donder gauge to the Coulomb gauge (106). To prove this, we first observe that when rewritten in terms of the gauge-invariant quantities (101), the De Donder condition (74) implies

[TABLE]

For static sources $\Phi_{a}=\dot{\Phi}=0$ and $\Psi=\Phi$, hence all right-hand sides vanish. It is then possible to specialize the De Donder gauge choosing $W=B=B_{a}=0$ everywhere. Having done so, $h_{00}=-2c^{-2}\Phi=h_{s}/3$. We have thus achieved the Coulomb gauge described earlier, as opposed to the TT gauge, and again this is a complete gauge fixing singling out a member of the De Donder family. In other words, the De Donder gauge is compatible with the Coulomb gauge for static sources. Notice that this is what happens also in electromagnetism, where the Lorenz gauge is compatible with the Coulomb gauge for static sources.

We conclude with a word on another important gauge fixing that can be used to study radiation including sources, the Bondi, or Bondi-Sachs gauge. In this case the background Minkowski metric is not in Cartesian coordinates, but rather spherical coordinates, and time is replaced with retarded time $u=t-r$. In these coordinates the Minkowskian Christoffel symbols do not vanish. However a crucial advantage is that it is a gauge that can be defined also in the full theory, and allows a non-perturbative description of gravitational waves, at the price of introducing an asymptotic expansion away from the sources.

IV Detection of GWs

IV.1 Coordinate displacements versus physical displacements

The simplest way to study the effect of a gravitational wave, is to compute how it changes the physical distance between free-falling test masses. These follow time-like geodesics, whose tangent vector field $u^{\mu}$ satisfies the geodesic equation

[TABLE]

Here $\tau$ is the proper time $\tau$, and $u^{2}=-c^{2}$. If the masses are initially at rest, we have $u^{a}=0$ and

[TABLE]

For vacuum solutions in the TT gauge $h_{0\mu}=0$. Hence $\dot{u}^{\mu}$ vanishes at lowest order, and the coordinate distance (as well as the coordinate time delay) between two nearby time-like geodetics remains the same during the passage of the wave. This result provides us with an interpretation of the TT gauge: it is a choice of coordinates which are labelled by the position of free-falling test masses. It is thus the linearized version of the temporal gauge in the full theory discussed earlier. Now, even though the coordinate distance between two test masses remains the same in this gauge, their physical distance does not. It is given by

[TABLE]

Here $\hat{e}^{a}$ is the tangent to the curve connecting the test masses, and $\lambda$ an arbitrary parametrization thereof. Since we are working in perturbation theory, we can consider as curve the background geodesic connecting the test masses, and since we are working in Cartesian coordinates of the background, we can choose $\hat{e}^{a}$ to be constant along a coordinate axis, and take $\lambda$ as one of the coordinates. Then $L_{0}$ is the coordinate distance. It coincides with the physical distance in the background flat metric, but differs from it when the spacetime is perturbed. If we further assume that the wavelengths of the wave are much bigger than $L_{0}$, we can ignore the space dependence of $h^{\scriptscriptstyle\rm TT}_{ab}$ and write the result as

[TABLE]

The approximation becomes of course exact if the direction of propagation of the wave is orthogonal to the axis connecting the masses.

The discussion offers an example of one of the most important lessons of general relativity, namely the importance of distinguishing coordinate effects from physical results. General covariance guarantees that all calculations can be performed in any coordinate system. But one has to always make sure that the physical consequences derived from the calculations are coordinate independent. In this example, we found that the coordinate position of test particles are unaffected by the passage of the wave. This is true only in the TT gauge. The relevant coordinate-independent quantity is the physical distance, and we found that it changes. Furthermore, it would change in the same way in any coordinate system preserving the endpoints of the integral: the geodesic distance between two physical points is a gauge-invariant observable.

As we have already discussed, even though any coordinate system can be chosen, choosing a good coordinate system is important to simplify calculations. Even thought TT coordinates hide the passage of the wave, they are convenient because we can write the physical distance using fixed extrema in the integral at all times.212121The same effect is used when choosing the synchronous gauge in cosmology, it is often convenient to choose a coordinate system such that the values of the coordinate grid represent galaxies, so that their coordinate distance does not change, while the physical one does. Another convenient gauge is the one of Fermi normal coordinates, which works in the opposite way: the metric is trivial at lowest order near two geodesics, see (32), hence the physical distance is entirely captured by the coordinate distance. To see this, let us consider the geodesic deviation equation. It is given by

[TABLE]

where $\xi^{\mu}$ is a vector connecting neighbouring geodesics, chosen such that $\xi\cdot u=[\xi,u]=0$. At first order in $h$, and with the assumption of vanishing initial velocities, it reduces to

[TABLE]

In the TT gauge, the last term vanishes and the first one too if the initial velocity was zero. The coordinate distance between the geodesics stays constant, in agreement with the result already derived using the geodesic equation above.

The description changes completely if we use a gauge corresponding to a local inertial frame, such as the Fermi normal coordinates, which can be used to set to zero the Christoffel symbols all along a chosen geodesic. Doing so, the LHS of (116) becomes simply a second partial derivative, and the equation reads

[TABLE]

where in the second equality we have neglected any contribution from the potentials. This shows that in a local inertial frame, the deviation between nearby geodesics can be described directly in terms of their coordinate separation. This is consistent with what previously seen in (32), the metric only changes at quadratic order in the coordinate distance from the origin, hence coordinate distances coincide with physical distance at first order. The geodesic deviation equation (118) shows that an apparatus that can detect tidal effects will be sourced only by the physical components of the GW, and not by the gauge ones. The equation can be easily solved at first order in $h$, with

[TABLE]

Since in this gauge the coordinate distance coincides with the physical distance at first order, we recover the gauge-invariant result (115) but where this time the metric is unchanged, and it is the extremum of the integral that has moved.

The change of proper distance (115) also shows the meaning of the wave’s polarizations. Using the example (85) of a monochromatic wave propagating along the $z$ axis, and setting the test masses at $z=0$, we can write the relative change in physical distance as

[TABLE]

An $h_{\scriptscriptstyle\rm+}$ polarization would cause pairs of masses along the $x$ and $y$ axis to periodically approach and recede, hence drawing a $+$-like pulse in time, see Figure 5. An $h_{\times}$ polarization would cause the same effect but along the axis $\hat{e}^{a}=(1,1,0)/\sqrt{2}$, namely rotated by 45 degrees. This type of deformation is also called ‘shear’ of the congruence of time-like geodesics followed by the test masses. It can be visualized even more clearly if we use (119). We consider a circular distribution of test masses centered around the origin in the plane perpendicular to the direction of propagation of the wave, see Fig. 5. Then we can identify the displacement vector with the coordinate vector of each mass (labelled by $i$), and the effect of a monochromatic wave of frequency $\omega$ is

[TABLE]

The effect is shown in Fig. 5, where the period $T=2\pi/\omega$. The spin-2, quadrupolar nature of the gravitational force is evident from the shape of the shear deformation, and can be compared with the tidal deformation of earth’s oceans.

If external forces are present, on top of the gravitational one, then (112) and (116) acquire additional terms on the right hand side. So in particular (118) now reads

[TABLE]

as an example, consider a material bar. While the effect of the gravitational wave is to stretch spacetime changing the physical distance between the molecules of the bar, there are also electromagnetic forces that hold the bar together, and which are intrinsically much stronger. For instance, the Coulomb interaction between two electrons one angstrom apart is $\sim 10^{42}$ times stronger than its Newtonian counterpart. For this reason, one could in principle use simply a rigid ruler measuring the distance between two freely falling masses to detect gravitational waves. The problem with this idea is the weakness of the waves, and more sophisticated experiments are required to actually observe the waves.

A special case of external forces appearing on the right-hand side of (121) are the inertial forces. These, by the equivalence principle, can be reabsorbed in a coordinate transformation of the metric. For instance if the frame has both an acceleration $\vec{a}$ and an angular velocity $\vec{\Omega}$ with respect to a local inertial frame, then Ni:1978zz

[TABLE]

and (121) becomes

[TABLE]

All the gravitational effects as well as further non-inertial effects such as centrifugal acceleration are $O(x^{2})$. So in order to be capable of detecting gravitational waves, a detector must first of all be freed from all the external forces that would otherwise drown the signal in noise.

In realistic physical systems, the emission will not be a plane wave, but rather a wave packet with finite temporal extension. The effect on the circular distribution will then be a superposition of different frequencies and different helicities, each with their own (time-dependent) amplitude. The temporal finiteness of the signal can also lead to a new type of effect: after the wave has passed, the distribution will stop oscillating, but its shape will in general not be the same as before the wave’s arrival. This effect is called displacement memory, and we will see below in Sec.VII.6 an explicit example. The effect carries the memory of the wave, since it permits in principle to detect the passage of a gravitational wave even after the event. In practise though the detection is very difficult, because the external forces that make up the matter distribution will act and bring it back to its rest configuration. It is nonetheless one of the targets of future detectors Grant:2022bla .

IV.2 Interferometers

Let us briefly describe how the formulas above are used in the most common type of detectors, laser interferometers. Other chapters in this collection will cover more details as well as the types of detectors. The basic idea of a laser interferometer is to detect physical changes like (114) from the time-of-flight of monochromatic light signals. This can be done easily in the linear theory if we make the additional approximations that gravitational potentials can be neglected,222222The potentials generated by the source can be naturally neglected because they fall off faster than the radiative modes, so this approximation concerns mostly the potentials of the local gravitational field of the observer. and that the wavelength $\lambda$ of the signal is much longer that the arms of the interferometer. The first approximation guarantees that the only source of curvature comes from the wave, hence the Riemann tensor scales like $\lambda^{-2}$. We can then set up a free falling frame say in Fermi normal coordinates centered on the beam splitter’s geodesic. Thanks to the second approximation, the spatial projection of the null geodesics follows straight lines, hence the time of flight is directly related to the physical distance along the interferometer’s arms. The latter is given by (115) regardless of the direction of the wave, thanks again to the assumption that the wavelength is much larger than the arms’ length. Denoting $\hat{e}^{a}_{1,2}$ the two axis, we have

[TABLE]

For a typical signal $h\sim 10^{-21}$ (see overview Section I.4), hence the difference in arrival time would be $\Delta T=(L_{2}-L_{1})/c\sim 10^{-26}s$ which is way too small to be measurable. Two ingenious ideas come to the rescue. The first is Michelson-Morley’s idea to measure not time but phase interferences, and the second it to increase the effective path of light through Fabry-Perot cavities. If we set the lasers so that the phases at the beam splitter are identical, the phase shift after the travel to and back from the mirrors will be

[TABLE]

where $N_{p}=1$ for a Michelson device and up to 300 for the Fabry-Pérot one used by LIGO and Virgo, and the laser frequency $\nu$ can go as high as $10{\tt kHz}=10^{4}{\tt s^{-1}}$, thus improving the strength of a typical signal to $10^{-20}$, which is actually – and remarkably – an observable phase difference.

Combining the last two formulas we have

[TABLE]

where now $\hat{e}^{a}$ are the unit vectors giving the direction of each arm. It is possible to express the result in terms of the two wave polarizations by replacing $h_{ab}^{\scriptscriptstyle\rm TT}$ with (92). To compute the scalar products between the vectors corresponding to the detector’s arms and the vectors defining the polarizations in the plane orthogonal to the propagation, we introduce a rotation from the detector’s frame to the frame of propagation, plus a reflection to take into account the fact that the axis of propagation is opposite to the direction of the source acting also on $y$ to keep right-handed orientation of the frame. If we denote $(\theta,\varphi)$ the angles identifying the direction of propagation, the result is

[TABLE]

where the coefficients

[TABLE]

are called detector’s pattern functions.232323This formula assumes that the arms are perpendicular, and a fixed polarization basis. It can be generalized to include an additional rotation of the two polarization basis, as well as a non-perpendicular angle between the arms. This shows that while a single two-armed interferometer is sensitive to both polarizations, it cannot distinguish them. It also shows that the sensitivity depends on the relative orientation with respect to the sources. The dependence is very strong, to the point that there are directions in which the detector is completely blind, like $(\theta,\varphi)=(\frac{\pi}{2},\frac{\pi}{4})$. Hence the importance of multiple detectors in order to increase sensitivity in every direction and the possibility of distinguishing the polarizations. Multiple detectors also allow studying the localization of the source via triangulation.

If the approximation $\lambda\gg L_{0}$ is no longer valid, then one has to take into account the redshift changes during the time of flight, see for instance discussion in Cornish:2002rt ; AnderssonBook .

V Generation of GWs from sources

V.1 Introducing sources

Following the principle of general covariance, the matter Lagrangian should satisfy the property (24), namely be written solely in terms of the dynamical matter fields and spacetime metric, and no additional background fields. The simplest way to obtain a viable matter Lagrangian is then to start from the one used in the absence of gravity, and ‘covariantize’ it by the replacements

[TABLE]

Doing so introduces a minimal coupling of matter to the gravitational field. Additional interactions can be included if phenomenologically or theoretically motivated, provided they respect the principle of general covariance, embodied for instance by condition (24). Having done so, we define the matter energy-momentum tensor

[TABLE]

Inserting this definition in (24), and using (21), we obtain

[TABLE]

The first term on the right-hand side is proportional to the matter’s equations of motion, and vanishes on-shell. The rest is a total derivative and vanishes in the absence of boundaries or with appropriate boundary conditions. Since the identity holds for any $\xi$, we conclude that on-shell,

[TABLE]

This equation replaces the familiar conservation of the energy-momentum tensor guaranteed by Noether’s theorem in flat spacetime. More precisely, the Noether current of the total Lagrangian ${\cal L}_{\scriptscriptstyle\rm EH}+{\cal L}_{\scriptscriptstyle\rm M}$ is

[TABLE]

where $E$ are Einstein’s equations (16), namely (26) with the vacuum equations replaced by the equations in the presence of matter), and whose conservation requires to be on-shell of both the Einstein’s and matter’s field equations:

[TABLE]

Even though (132) is often referred to as the general covariant version of energy-momentum conservation, it is important to remark that it is not a conservation equation in the usual sense. To understand this point, let us follow the usual procedure to obtain Noether charges from the current, and apply Stokes’s theorem to a finite region $M$ with boundary $\partial M$. To do so we need a scalar, which we obtain by contracting the left-hand side of (132) with a vector $\xi^{\mu}$. After integrating by parts, we find

[TABLE]

for a boundary $\partial M$ with induced metric $q$, normal $n_{\mu}$ coordinates $y$. If $\nabla_{(\mu}\xi_{\nu)}$ vanishes, namely if the Killing equation is satisfied, then (132) can be turned into a conservation law. To do so, we consider the case in which $\partial M$ consists of two space-like hypersurfaces $\Sigma_{1}$ and $\Sigma_{2}$ connected by a time-like boundary $\cal T$, see Fig.7. If the fields satisfy conservative boundary conditions at $\cal T$ (typically $\cal T$ asymptotically far away and fall-off conditions on the fields), then

[TABLE]

for each Killing vector $\xi$. Therefore (132) gives as many conserved quantities as there are isometries in spacetime. For flat spacetime, these are the ten Poincaré charges. For a generic dynamical spacetime, there are none.

The validity of (132) implies the matter equations of motion, as we have seen from its derivation. In particular, if matter consists of test particles, namely free motion without self-interaction and ignoring the back-reaction on the metric, this equation implies the geodesics equation in curved spacetime. This is for instance how one can derive the relativistic corrections to the Kepler problem, by evaluating (132) on the Schwarzschild background. At lowest order in the weak-field expansion (48), (132) reduces to the energy-momentum conservation law in flat spacetime,

[TABLE]

This means that at lowest order the matter can interact with itself, but not with the gravitational field: the sources follow geodesics in flat spacetime (that is, straight lines). To include the effect of gravity on the sources we must go beyond the lowest order. In other words, the linearized theory still describes gravity in the Newtonian way, namely as a force acting in flat spacetime. Of course, it already contains departures from Newton’s theory, since it includes the special relativistic effects such as the gravito-magnetic interaction and radiation.

V.2 Source multipoles

Let us study the conserved quantities that arise on the Minkowski background. We choose $\Sigma$ to be a global hypersurface of constant time $t$, write its unit normal as $n_{\mu}=-\partial_{\mu}t$, and $\xi$ is one of the ten Poincaré Killing vectors (50). We can then use (136) to identify ten conserved quantities. Four are the energy and momentum

[TABLE]

corresponding to $a^{\mu}{}_{\nu}=0$ and unit values of $b^{\mu}$. The remaining six are the relativistic angular momentum

[TABLE]

corresponding to $b^{\mu}=0$ and unit values of $a^{a}{}_{\nu}$ and $a^{0}{}_{\nu}$ respectively. Here $\epsilon_{abc}$ is the completely anti-symmetric Levi-Civita symbol in flat spacetime. Their conservation can be easily checked. We start by separating (137) in time and space components,

[TABLE]

Then using Stokes’ theorem and vanishing boundary conditions we immediately see that

[TABLE]

The first of (138) is the total energy, but we followed the custom in the literature to denote it $M$ and refer to it as ‘mass’, using Newtonian language. The first of (139) is the angular momentum with respect to the frame defined by $n_{\mu}=-\partial_{\mu}t$. The second conserved quantity is the ‘boost charge’, and can be rewritten as $K^{a}=I^{a}-tP^{a}$, where we introduce the center-of-mass position

[TABLE]

Conservation of $K^{a}$ is thus the statement that the center of mass moves following the total momentum.

The conserved quantities can be used to fix a reference frame as follows. First, we can choose the rest frame, where $P^{a}=0$. This removes the freedom of Lorentz boosts. Then, we can fix the origin to be in the center-of-mass, where $I^{a}=0$. This removes the freedom of spatial translations. The rotation freedom can be fixed choosing the axis so that $L^{a}$ has only one component (say $z$), and the remaining SO(2) freedom is fixed choosing an axis in the plane perpendicular to $L^{a}$. Finally the time translation symmetry is fixed setting the zero value of the clock.

The quantity $I^{a}$ is also called mass-dipole moment. The terminology comes about if we see $\rho=c^{-2}T_{00}$ as a distribution, then $I^{a}$ is the first moment of that distribution. Following this logic, we introduce a multi-index notation for the higher multipole moments:

[TABLE]

The conservation laws (140) together with integration by parts in the absence of boundary terms provide relations between multipole moments and time variations of higher multipoles, such as

[TABLE]

and so on. The first one above is the conservation of $K^{a}$ already seen, and relates the momentum monopole to the mass dipole time variation. The second one allows one to determine the total effect of the stresses in the matter in terms of the second time derivative of the mass quadrupole. These relations are useful because it is typically easier to measure and interpret the multipole moments of the mass and momentum distributions, rather than the spatial stresses.

When working with multipoles, it is typically convenient to organize them into irreducible representations of the rotation group, which are label by a an integer number $l$ and have $2l+1$ components each, as recalled earlier. This can be achieved expanding the distribution in spherical harmonics, e.g. $\rho=\sum_{l,m}\rho_{l,m}Y_{l,m}$, then the integrals of the modes $\rho_{l,m}$ are the irreducible multipoles. It is possible although more cumbersome to do this composition directly in Cartesian coordinates without introducing spherical harmonics. One then gets

[TABLE]

and so on. Notice that while the multipole moments are useful at all orders in perturbation theory, the conservation laws (141) are only valid at lowest order.

V.3 Solving the wave equation with sources

We are interested in the emission of gravitational waves from matter sources, without incoming radiation. This can be imposed choosing the retarded Green function and setting to zero the independent degrees of freedom $h^{\circ\scriptscriptstyle\rm TT}_{\mu\nu}$. The general solution is then

[TABLE]

using (326) and the specialized De Donder gauge with no homogeneous solution. Even in the linearized approximation, the integral is in general very complicated and there is no analytic solution. So we resort to approximation schemes. In particular, we introduce two independent approximations:

• $(i)$

Wave-zone approximation: we assume to be very far away from the sources, that is $R:=|\vec{x}|\gg|\vec{x}\,^{\prime}|$. This allows us to expand the integrand in powers of $1/R\ll 1$. For the numerator, we have

[TABLE]

where $\vec{N}:=\vec{x}/R$, and

[TABLE]

where we introduce the retarded time242424Namely the time at which a signal travelling at the speed of light was sent in order to arrive at $t$.

[TABLE]

For the denominator, we have

[TABLE]

Furthermore, the direction of propagation of the wave coincides with the direction from the source, namely $-\vec{N}$ if we take the origin of the coordinates inside the source. Hence the TT projector can be written in terms of $\vec{N}$ instead of the wave vector.

• $(ii)$

Slow dynamics: We assume that the dynamics of the source is slow, so that time derivatives in (149) are small corrections. To understand why, consider that the integration coordinate $\vec{x}^{\prime}$ spans at most the size of the source, and if this has a typical frequency scale $\omega_{s}$ (for instance in a binary, the frequency of the orbit), then $v_{s}:=|\vec{x}^{\prime}|\omega_{s}$ is the velocity scale of the source. It follows that

[TABLE]

is suppressed by $v/c$. The Taylor expansion (149) is therefore controlled by the parameter $v/c\ll 1$, and it is called post-Newtonian expansion.

The approximated solution can thus be written as

[TABLE]

The first term is the leading order; the second term is the first of the PN corrections; the third term is the first of the $1/R$ corrections. Using the multipole definitions (143) and their conservation laws (144), we can rewrite the different components of the solution (153) as

[TABLE]

These are the first few terms of the double expansion in velocities and distance from the sources. We are not giving all metric components to the same higher order, this partial result is sufficient for our purposes. In the PN expansion this is the metric parametrization in the wave zone, with the no-incoming radiation condition, and completely gauge fixed. The lowest order of the time-time component reproduces the Newtonian result.252525In particular the lowest order of the Schwarzschild metric is obtained as the special case with constant $M$ and all the rest vanishing. Notice in fact that $h_{00}=\bar{h}_{00}+\bar{h}/2=2GM/c^{2}R$. The first PN correction is the movement of the source, and can always be set to zero by going to the rest frame. Doing so eliminates the lowest order of the $\bar{h}_{0a}$ component. The first corrections in that component contain the gravito-magnetic effects relevant to the Lense-Thirring effect, for instance. Notice also that the angular momentum is sub-leading in $R$, as one could have expected from a large distance expansion of Kerr’s metric.

The radiative degrees of freedom are in the spatial components (154c) and can be extracted acting with the projector (89). We have the mass quadrupole at leading order, and the first PN correction features the mass octupole and momentum quadrupole. We can immediately remark the absence of monopole and dipole contributions to the emission of waves. This is a direct consequence of the conservation laws, since they imply that the mass monopole and dipole have vanishing second time derivatives. As a consequence, an oscillating spherical distribution would not emit gravitational waves, in agreement with Birkhoff theorem in the full theory, nor would a distribution with axial symmetry rotating at constant velocity, in agreement with Kerr’s solution.

Applying the TT projector removes any trace, hence $h^{\scriptscriptstyle\rm TT}=\bar{h}^{\scriptscriptstyle\rm TT}$ and one can replace $I_{ab}$ with the irreducible quadrupole moment $Q_{ab}$, and obtain at lowest order

[TABLE]

This is the celebrated first quadrupole formula, derived by Einstein in 1918: The dominant radiation in the slow-motion approximation arises from the acceleration of the quadrupole moment of the mass distribution. From this we can also obtain the expressions for the two independent polarizations. If $\vec{k}=\hat{z}$, we can use (72g) and

[TABLE]

For a general $\vec{k}$ it is obtained replacing $h_{ab}\to(G/c^{4}R)\ddot{Q}_{ab}$ in (93). Notice also that $P^{\scriptscriptstyle\rm TT}(Q)=P^{\scriptscriptstyle\rm TT}(I)$ since the projector removes the trace, hence we can replace $Q_{ab}$ with $I_{ab}$ in these expressions.

Let us make some order-of-magnitude estimates. By dimensional analysis, the mass multipoles scale like $Mr^{l}$, where $r$ is the typical size of the source. If the dynamics of the system has a typical velocity scale $v$, then $Q\sim Mr^{2}$ and $\ddot{Q}\sim Mv^{2}$. This gives

[TABLE]

For example, if we extrapolate this formula to relativistic speeds $v\sim c$ for the merger of two 10-solar-masses black holes, we get a $10^{-18}$ amplitude at galactic distances, and $10^{-21}$ at 100 Mpc where the Virgo cluster is located.

This estimate is the lowest order of various approximations, which is useful to recap here: (1) weak-field, PM expansion; (2) long-distance, multipolar expansion; (3) small velocities, PN expansion. To obtain more accurate results, one has to include higher order corrections. Doing so is actually far from simple. Not only do we have three different expansion parameters with non-trivial hierarchies among them, we also have to face both technical and conceptual challenges. Let us list a few, and tools used to deal with them. The PN expansion is not a convergent series, but rather what is known as an asymptotic series. Its accuracy degrades as we increase $R$. Dealing with this mathematical problem requires techniques such as the matched asymptotic expansion. Related to this is also the more conceptual issue that the causal propagation determined by the Green’s function at lowest order follows the null cones of the background Minkowski metric. But null cones are bent by the gravitational interaction, hence higher order corrections have to also modify the retarded time to the correct one. For instance for the Schwarzschild metric the correct retarded time is

[TABLE]

Hence higher orders change the notion of retarded time.262626This issue can be solved non-perturbatively using the Bondi coordinates mentioned earlier, but then one is switching from a PN expansion to an asymptotic expansion of the full theory, and the approach is both conceptually and technically different. Another tricky effect comes in at higher orders: the waves backscatter and self-interact, causing a delay in part of the signal, which starts travelling inside the light-cone, similar to light slowing down in a medium due to interactions with the medium. Then the total signal includes a ‘tail’ that comes after the main part of the signal. To take this into account one has to include effects that arise from integration over time. Another problem is that divergences appear after the first iteration, because convolution of Poisson integrals diverge even if the initial source has compact support. To regularize this unphysical divergence one has to split the integrals into near-zone and far-zone integrations.

These and other types of difficulties plagued the theory throughout most of the seventies, and were addressed thanks to the work of many brilliant researchers, including pioneers like Thorne, Will, and Damour. On the phenomenological side, people thought for a while that the lowest quadrupole order would have been enough to match experiments, given the weakness of the waves. Later theoretical work, e.g. the seminal paper Cutler:1992tc , clarified the observational sensitivity to the PN corrections and justified the importance of the endeavour. The task is very challenging, and researchers have come up with different approaches. We refer to the specialized literature Thorne:1980ru ; Blanchet:2006jqj ; poisson ; Goldberger:2004jt for reviews of this more advanced topic. In the following we will content ourselves to stay at lowest order, which is enough to understand the basics of the physics, if not for a detailed match to observations.

VI Dissipation by gravitational waves

The first historical evidence of gravitational waves was the orbital decay of the Hulse-Taylor pulsar. Accurate measurements showed that the orbital decay was consistent with the prediction of general relativity. Indeed, general relativity predicts that gravitational waves carry energy away from a system that produces them. In this Section we describe how this prediction arises at the lowest order in perturbation theory, and in the next Section we show how it can be applied to predict the orbital decay.

VI.1 Energy of gravitational waves

Let us look at a gravitational wave as a spin-2 field propagating on the Minkowski background. Thanks to the Poincaré invariance of the background, we can apply Noether’s theorem and derive a conserved energy-momentum tensor for $h_{\mu\nu}$. An explicit calculation starting from the linearized Lagrangian gives

[TABLE]

where the label N stands for Noether, and we assumed here the De Donder condition to simplify the expression.This tensor is conserved, namely $\partial_{\mu}t^{{\scriptscriptstyle\rm N}\mu\nu}\,\hat{=}\,0$, but has no clear physical meaning, because it is not gauge-invariant: It changes under a linearized diffeomorphism (52), including those compatible with the De Donder condition, and consequently assigns a non-zero value of energy-momentum to pure gauge modes. Furthermore, we can make it vanish entirely at any point using Riemann normal coordinates, since in these coordinates the first derivatives of the metric vanish at that point. Since it is zero at one point in one coordinate system but not in others, it is not a tensor. It is usually referred to as ‘pseudo-tensor’.

The only gauge invariant quantities that can be extracted from (158), at least at lowest order, are global ones, such as the total energy and momentum on a space-like hypersurface $\Sigma$ of constant time,

[TABLE]

In fact, at lowest order in perturbation theory, a gauge transformation acts as (52) with $\xi=O(h)$, under which (158) transforms as a total derivative. Spatial total derivatives do not contribute to (159) under the standard condition that the perturbation goes to zero at infinity, and temporal ones can be converted to spatial derivatives on shell of the wave equation, hence the same argument applies. Similarly, one can define the total angular momentum as

[TABLE]

For our applications below, the total energy-momentum and angular momentum will be sufficient. They are not sufficient however for iterating the perturbation theory beyond lowest order, and we now briefly describe some of the conceptual aspects of what is done in the literature in this respect. First of all, the reader may recall that Noether currents are not unique, being defined only up to adding total derivatives whose conservation is trivial. In the case at hand, it means that any quantity of the type

[TABLE]

where $U$ has the same index symmetries as the Riemann tensor, is an equally valid Noether current. One may hope that there exists a representative in the equivalence class (161) that would be gauge invariant, but this is not the case: the lack of gauge-invariance is simply the linearized version of the fact that there cannot be any local tensor representing the energy of the gravitational field. Such quantity will have to be zero in a local free-falling frame where the effects of gravity are absent, and if it were a tensor, it would then be zero in any frame.272727A tensorial quantity capturing some aspects of gravitational energy can be constructed using the Bel-Robinson tensor, but it is fourth-order in derivatives, therefore does not have the right physical dimensions, and will capture only higher-order terms of the gravitational energy. This is therefore a direct consequence of the equivalence principle, and its mathematical implementation via general covariance.282828It is instructive to put this problem in perspective with what happens in the electromagnetic case. If one computes the canonical energy-momentum tensor of Maxwell’s theory using the Noether formula, one also finds a meaningless gauge-dependent expression. However, the Noether construction only defines the tensor up to total divergences, and it is possible to find one that gives a gauge-invariant expression, and which is furthermore symmetric and coincides with the one derived from the variation with respect to the metric. In gravity there is an analogue problem, but even adding total divergences it is not possible to find a local gauge invariant quantity.

The obstruction to identifying a local energy density for the gravitational field shows up very prominently in the full theory. We have already discussed how Noether charges for generic diffeomorphisms are all trivial in the bulk of spacetime. Similarly, the bulk Hamiltonian one finds from the Legendre transform of the Lagrangian is a sum of constraints, and thus identically zero when evaluated on solutions. Any attempt to work around these facts and define quasi-local observables representing the gravitational energy unavoidably run into trouble with ambiguities and dependence on coordinates or other unphysical background structures Szabados:2009eka . The clearest well-defined resolution to this problem is to consider not the energy density, but only the total energy. This is useful when describing isolated systems, namely spacetimes that are fully dynamical in a certain region, but become well approximated by flat spacetime at large distances from this region. In this case, one can introduce a physically meaningful notion of boundary to the spacetime, and exploit the fact that the Hamiltonian picks up a boundary contribution which is non-vanishing on solutions. The resulting surface charges can be used to characterise the total energy momentum and angular momentum of the system, and can be derived as Noether charges as well. Examples of this construction are the ADM charges at spatial infinity, and the BMS charges at future null infinity, as mentioned in Sec. II.2.

While looking at global quantities such as the total energy of a free GW (159), or the surface charges in asymptotically flat spacetimes mentioned at the beginning of the section, is the safest way to define energy in the full theory, the perturbative treatment offers an alternative, ‘quasi-local’ possibility. Perturbatively in fact, it is possible to construct gauge-invariant quantities by introducing a spacetime averaging procedure based on the properties of the background. We consider a region $L$ whose size is much larger than the typical wavelength $\lambda$ of the perturbation, but much smaller than the typical wavelength $\lambda_{B}$ of the background (which is infinite for a flat background), and we define the averaging of a functional $F$ as $\langle{F}\rangle:=\frac{1}{L}\int_{L}F$. If applied to an expression quadratic in the Fourier modes like (158), the procedure suppresses combinations with different frequencies or different phases, in a way completely similar to how the total energy in a standard background-dependent theory comes mainly from positive interference superposition of waves. The difference is that in background-dependent theories averaging the energy is a choice, since the local energy density is theoretically also well defined. In gravity it is not a choice but mandatory, since there is no meaningful local energy density, and furthermore care is needed to define correctly the procedure in a way to make it compatible with general covariance. Detailed analysis Isaacson:1968hbi ; Isaacson:1968zza ; Burnett shows that the result of the procedure is that expressions under the averaging sign can be freely integrated by parts in space and, upon going on-shell, also in time derivatives since a wave propagates on the light-cone. For instance,

[TABLE]

outside the sources. Under this procedure, we find

[TABLE]

in a region outside the sources. One can show that the averaging procedure makes the right-hand side gauge-invariant Isaacson:1968zza . This means that it can be expressed in terms of the TT projection and the gauge invariant potentials. The latter can be neglected if the sources variation (induced by the partial derivatives in the expression above) occurs over much longer time scales than the $h^{\scriptscriptstyle\rm TT}$ wavelengths.292929In poisson this step is called short-wave approximation, and it is performed without the averaging, in the context of the Landau-Lifshitz approach described in Appendix C. This motivates the definition of

[TABLE]

This quantity is actually gauge-invariant at lowest order, since the only non-invariant terms are the partial derivatives, and these transform linearly in $\xi\sim O(h)$. It follows that in so far as only lowest order results are needed, we can use the simpler expression (164) as a proxy for the averaging procedure. Notice also that it coincides with (158) in the TT gauge.

The expression (158) in the TT gauge is the one used by Einstein to determine the energy carried away by gravitational waves. It has, however, a limited applicability. We have already discussed its gauge dependence. Another issue is that it relies heavily on the specific background chosen, and had we worked with a non-isometric one, then there would be no Noether charge to begin with. Furthermore, it is not clear how to extend this construction to treat higher orders in perturbation theory. These shortcomings can be addressed if we look at a different definition for the gravitational energy-momentum pseudo-tensor, based on the actual back-reaction on the metric caused by the gravitational waves. In fact the actual “effective” source that determines the second-order metric perturbation is not (158), but rather the second order expansion of the Einstein tensor, that we denoted $t^{\scriptscriptstyle\rm G}_{\mu\nu}$ in (38). The candidate gravitational energy-momentum pseudo-tensor $t^{\scriptscriptstyle\rm G}_{\mu\nu}$ obtained in this way is also conserved. In fact, an explicit calculation shows that it differs from (158) precisely by a term like (161), with $U$ a certain quadratic expression in derivatives of $h_{\mu\nu}$. It has the improved property that it depends on second derivatives of the metric, so it cannot be made to vanish at any given point. However, it is still not gauge-invariant. Therefore one has to invoke again the averaging procedure in order to extract gauge-invariant information. Upon doing so, one finds that the two prescriptions give a consistent answer Isaacson:1968zza :

[TABLE]

This matching supports averaging as a viable way to extract unambiguous and gauge-independent quantities. The prescription $t^{\scriptscriptstyle\rm G}$ overcomes some limitations of the Einstein-Noether construction. It can be used in perturbation theory around an arbitrary background, and can be systematically extended to any order in perturbation theory, by computing higher order corrections $G^{\scriptscriptstyle\rm(n)}_{\mu\nu}$ and evaluating them on the perturbed solution. This procedure is however not very handy, and a better scheme is the one proposed by the Landau-Lifshitz reformulation of Einstein’s equations, see e.g. Blanchet:2006jqj ; poisson . There one changes variables from the metric to a densitized inverse metric ${\mathfrak{g}}^{\mu\nu}:=\sqrt{-g}g^{\mu\nu}$. One advantage of this reformulation is that it provides a full non-perturbative expression for a candidate energy-momentum, known as Landau-Lifshitz pseudo-tensor. This has the usual limitations (gauge-dependence and vanishing at any point in a local inertial frame) dictated by the equivalence principle, but has the merit of being set up in a way that makes it very natural to develop a systematic perturbative expansion, since the pseudo-tensor is defined already at non-perturbative level, and does not need to be determined order by order as in the previous approach. Furthermore, it provides a prescription for the energy, momentum and angular momentum as surface charges that, even though restricted in validity to Cartesian coordinates in the region far away from the sources, can be evaluated including higher orders, and bypasses the need for the spatial averaging of volume integrals. For these reasons, the Landau-Lifshitz formulation is widely used by the community working in the post-Newtonian expansion. We review it briefly in Appendix C. The lowest order in the weak-field approximation of the Landau-Lifshitz pseudo-tensor differs from the previous two options by a term like (161), and also matches the gauge-invariant result after averaging:

[TABLE]

Similar considerations apply also to define the angular momentum of gravitational waves. In this case, the result of the averaging procedure starting from any of the three prescriptions described above motivates the following definition dewitt2011gravitation ; Thorne:1980ru ; poisson

[TABLE]

As before, one can drop the TT projector and still get a gauge-invariant quantity after global integration or averaging.

We conclude that in so far as one is interested only in lowest order results, any of these three choices are equally good. For a systematic perturbative expansion, the latest is the better one. In fact even if the construction of gauge invariant quasi-local quantities via averaging is conceptually useful to clarify how gauge-invariant information could be extracted in principle, it is not very practical. To set up a systematic perturbative expansion, it is easier to work with gauge-fixed quantities at all intermediate steps, and then extract only at the end the physical predictions in terms of gauge-invariant observables. For instance, there is no problem in working with the non-averaged notions of energy-momentum and angular momentum pseudo-tensors, as long as one does not attempt to give them a direct physical interpretation. The idea is to use them to perform calculations, and at the end read off the physical dynamics not from their evolution but from that of gauge-invariant quantities such as the amplitude and frequency of TT modes, or the evolution of relative distances such as the periastron of an orbit. This is the logic used in the PN expansion, and based on the Landau-Lifshitz reformulation Blanchet:2006jqj ; poisson .

VI.2 Dissipation equations

The fact that (164) is conserved means that we can derive identities between time and spatial derivatives like those that led to the conservation laws (141) for the matter sources. The key difference however is that the matter sources had compact support, hence we could neglect boundary contributions when integrating by parts. This is no longer true for the gravitational contributions, since the waves have non-compact support. The non-vanishing of the boundary terms has the effect that the ‘charges’ corresponding to energy, momentum and angular momentum are no longer conserved. This dissipation is precisely the statement that gravitational waves carry energy and have a physical impact on the system.

At lowest order, we do not even need to use the conservation equation in order to study the dissipation, because of a special property of the explicit solution (154). Each metric component has functional dependence on coordinate of the form $f(t_{\scriptscriptstyle\rm R},N^{a})$. For such functions, it is easy to check that

[TABLE]

The leading order of this approximation plays an important role in simplifying many formulas in the wave zone, where $R\gg 1$.

Let us begin our analysis from the flux of gravitational energy, namely the emitted power. This is given by

[TABLE]

Stokes theorem choosing as boundary a 2-sphere of radius $R$ in the asymptotic region (hence the outgoing unit normal is simply $N_{a}$, and $dS=R^{2}d^{2}\Omega$ where $d^{2}\Omega=\sin\theta d\theta d\phi$) and (164) in the last equality. The spatial derivative can be replaced at lowest order with a time derivative using again (168), and we arrive at

[TABLE]

where in the last step we used the explicit form (154) of the solution, in particular (155). To evaluate the integral, we observe that the only angular dependence occurs in the TT projector (89). Using the following formula,

[TABLE]

we find

[TABLE]

This is the second famous quadrupole formula of Einstein Einstein1918 . It gives the instantaneous power radiated at a distance $R$ from the source and a time $t$, as a function of the quadrupole time variation at the retarded time $t-R/c$. Notice that the index contraction occurs over all indices of the (traceless) quadrupole moment. The effect of the $TT$ projection goes into a numerical factor, after the integration (171).

For the linear momentum loss,

[TABLE]

where we used twice (168) in the last equality. Since $N^{a}$ is an odd function on the sphere, the integral vanishes: there is no loss of momentum at lowest order (namely at order $G/c^{6}$, once we use the first quadrupole formula). A change in the total momentum of the system caused by the emission of GWs (‘kick’) occurs only at the next order $G/c^{7}$, when mixing of multipoles of different parity occurs.

For the angular momentum loss,

[TABLE]

Using the quadrupole formula (155) and performing the integrals using identities similar to (171), one arrives at

[TABLE]

Angular momentum loss occurs at the same order as energy loss, and involves one lesser time derivative.

VI.3 On the validity of the quadrupole formula

The quadrupole formula (155), and its application leading to the second quadrupole formula (172), were derived in the linear approximation. In this approximation, the energy-momentum tensor satisfies the flat spacetime conservation law (137), and the geodesics of matter are straight lines. It is thus valid only for systems whose gravitational interaction is negligible. It is not valid, in particular, for a binary system held together by gravity, even at the non-relativistic, Newtonian level. To treat a gravitational binary, one has to go beyond the linear approximation, using the iterative scheme described around (37). In the first step, one should include only the kinetic and non-gravitational pieces of the dynamics in the moment of inertia sourcing the quadrupole formula, consistently with the motion being along straight lines. The result is then used to source the effective energy-momentum tensor $t_{\mu\nu}$ needed for the second iteration of the field equations, which determines $h^{\scriptscriptstyle\rm(2)}_{\mu\nu}$ of (37). When this approach is applied to a gravitational binary, one finds that the second-order correction restores precisely the Newtonian contribution to the quadrupole moment. This is quite remarkable, and means that the same quadrupole formula is also valid if one includes the Newtonian interaction. See for instance poisson ; DeruelleUzan2018 for details.

We will make use of this fact below, and deduce the right results for gravitational binary systems using the quadrupole formula with the Newtonian potential included in the moment of inertia, without going into the technical details required to solve the second iteration, for which we refer to the cited literature. However, we should keep in mind that this formula can only be trusted because it has been derived including the second iteration, and not from the linearized approximation alone. We should also keep in mind that restoring the validity of the quadrupole formula after the second iteration in gravitational binaries is a special and remarkable fact, and not a general result. For instance, the linearized quadrupole formula gives the wrong answer for gravitational binaries in modified theories of gravity, see e.g. Taherasghari:2025mlf .

The fact that the linearized approximation provides potentially correct answers beyond its regime of validity contributed to the controversy that heated the debate around the quadrupole formula, and that were ultimately solved only with the systematic and rigorous developments of the early 80’s, see Damour:1982wm ; Kennefick:1997kb ; Will:2011nz for discussions.

VI.4 Back-reaction

Observable consequences of the emission of gravitational waves can be studied looking at how they impact the dynamics of the source. This can be done for instance evaluating the first PM correction to the source trajectories, by solving

[TABLE]

at various orders in the PN expansion. The details of this calculation can be found in poisson . In some cases however, it is possible to consider the following shortcut. At zeroth order, we have the Newtonian dynamics, and this comes with a clear identification of conserved quantities such as energy $E$ and angular momentum $J$. We then assume that the first order correction is obtained allowing these quantities to be not conserved, and equating their change to the dissipation caused by gravitational waves. That is, we posit

[TABLE]

insert the expressions (159) and (160) on the right-hand side, and solve the resulting equations. From the solutions we deduce how the source should change in time in order for its dynamics to be consistent with the dissipation caused by gravitational waves. This is what we do in the next Section, and for which the total expressions (159) and (160) are enough. Detailed calculations using the proper method (176), see e.g. poisson , confirm the validity of this shortcut, at least in so far as the lowest order in the PN expansion is concerned.

VII GWs from binary systems: elliptical, circular and hyperbolic orbits

We now apply the results of the previous section to determine the GW signal from binary systems. We will first consider the case of a bound system, with circular or elliptical orbits. These provide a simple yet realistic model of astrophysical sources that corresponds to the signals observed by LVK. We will see how one can express the two quadrupole formulas (and more generally the dissipation equations) in terms of the dynamics of the sources, compute the backreaction leading to orbital decay and increased wave emission, and produce analytic waveforms. We will also see explicitly the importance of the averaging procedure, which in the case of bound binary systems neatly separates the effects related to the two time-scales involved: the period of each orbit, and the ‘secular’ effects that cumulate over many orbits. We will then consider the case of unbounded, hyperbolic orbits, produce their waveforms. These orbits are interesting because they provide the simplest examples of displacement memory and gravitational capture. Throughout this Section, we will approximate the gravitational bodies with non-spinning point particles. This provides a good approximation at lowest order: detailed PN analysis shows that spinning and finite-size effect only enter at higher orders.

VII.1 Newtonian equations

We first recall the Newtonian equations of motion for two non-spinning point-particles of masses $m_{1,2}$, with relative position $\vec{r}=\vec{x}_{1}-\vec{x}_{2}$ and relative velocity $\vec{v}=\vec{v}_{1}-\vec{v}_{2}$. In the center-of-mass (CM) frame, the dynamics can be described by a single particle with position $\vec{r}$ and reduced mass

[TABLE]

The dimensionless quantity $\nu$ is introduced here for later convenience. The total energy $E=\frac{1}{2}\mu v^{2}-G\mu m/r$ and angular momentum $\vec{L}=\mu\vec{r}\times\vec{v}$ are conserved. The latter implies that the motion is confined to a plane, and we choose coordinates so that this is the $(x,y)$ plane. We then parametrize

[TABLE]

and introduce a second vector to form an orthogonal basis in the plane of the dynamics:

[TABLE]

The observer’s detector is at position $\vec{R}=R\vec{N}$ where in spherical coordinates the unit vector is $\vec{N}=(\sin\varphi\cos\theta,\sin\varphi\sin\theta,\cos\varphi)$, see figure 8.

Newton’s equations in the CM frame can be conveniently rewritten as

[TABLE]

where $e=\sqrt{1-(b/a)^{2}}$ is the eccentricity, and $p=a(1-e^{2})$ the semi-latus rectum, of an ellipse with semi-major axis $a$ and semi-minor axis $b$, see right panel of Fig. 8. These geometric quantities are related to the physical conserved quantities by

[TABLE]

and

[TABLE]

It follows from Eqs. (179)-(182) that

[TABLE]

where we have chosen the origin $\psi=0$ at periastron, and

[TABLE]

The Newtonian approximation requires $\vec{v}^{2}\ll{c^{2}}$, and thus the dimensionless ratio ${Gm}/{c^{2}p}\ll 1$.

Bound systems have eccentricity $e<1$, while unbound ones have $e>1$. The border case $e=1$ corresponds to parabolic orbits.

• •

Elliptical Orbits have $0<e<1$ with $-\pi\leq\psi<\pi$, $r_{\rm min}=p/(1+e)$ and $r_{\rm max}=p/(1-e)$. The orbital angular frequency $\omega_{0}$ and period $T$ satisfy Kepler’s laws

[TABLE]

• •

Circular orbits have $e=0$ and radius $r=p$, and orbital frequency

[TABLE]

• •

Hyperbolic orbits have $e>1$. Now $\psi_{-}(e)\leq\psi<\psi_{+}(e)$ where

[TABLE]

and correspondingly $\sin\psi_{\pm}=\pm e^{-1}{\sqrt{e^{2}-1}}$. These orbits are not periodic, but have a characteristic time-scale (a burst time scale) related to the characteristic frequency scale

[TABLE]

The closest distance of approach $r_{\rm min}=p/(1+e)$ at $\psi=0$. As $\psi\rightarrow\psi_{\pm}$, $v\rightarrow v_{\infty}$ with

[TABLE]

The GR modifications of the Newtonian equations of motion can be separated in two classes. First, the contributions of the potentials, of which the most famous is of course Einstein’s initial precession calculation: over an orbital period the perihelion of elliptical orbits advances by $\Delta_{h}=2\pi(3Gm/c^{2}p)$, while for hyperbolic orbits $\Delta_{h}=(\Delta_{r}/3)\left\{6\arccos(-1/e)+{e^{-2}}{\sqrt{e^{2}-1}}\left[2(2+e^{2})+5\nu(e^{2}-1)\right]\right\}$, see e.g. Damour1985 . Second, non-conservative or dissipative effects caused by the emission of gravitational waves. In the following we will ignore the first, and focus on the second. We will also ignore the precession introduced by dissipative effects, since this is higher order in the PN expansion.

VII.2 Energy and angular momentum fluxes

In the quadrupole approximation, the TT component of the waveform is given by Eq. (155), in terms of the traceless quadrupole tensor $Q_{ij}=I_{ij}-\frac{1}{3}I\delta_{ij}$. The mass quadrupole moment of the two point particles is

[TABLE]

This is where the discussion of Section VI.3 becomes crucial. If we are working in the linear approximation, we are only allowed to include in $I_{ab}$ the kinematical contribution from the velocities, and not the Newtonian gravitational potential, which would curve the geodesics away from the straight lines. However going through the second iteration has ultimately the neat effect of restoring the quadrupole formula with the Newtonian potential included. We will make use of this fact, and assume in the following that the quadrupole formula was derived already going through the second iteration of the field equations, and thus include the Newtonian potential in $I_{ab}$. Doing so, and using the equations of motion $d\vec{v}/dt=-Gm\vec{n}/r^{2}$, it follows that

[TABLE]

Thus from Eqs. (179) and (185) the non-zero components of $\ddot{I}_{ab}$ are

[TABLE]

and $\ddot{I}=2\nu mc^{2}\left(\frac{Gm}{c^{2}p}\right)e(e+\cos\psi)$. We have written the coupling constants as

[TABLE]

in order to highlight that $\ddot{I}_{ab}$ has the dimensions of energy. The third derivatives of the quadrupole tensor are straightforwardly obtained from Eq. (194) and (182) and read

[TABLE]

The GW perturbation is given by substituting these expressions into Eq. (155). In the direction $\vec{N}=\hat{z}$ the plus and cross polarisations are given by (see Eq. (156)),

[TABLE]

where the dimensionless amplitude is

[TABLE]

The time-dependence is determined from $\psi(t)$ which is a solution of Eq. (182). The polarizations $h_{+,\times}(t,\theta,\varphi)$ in an arbitrary direction $\vec{N}=(\sin\theta\cos\varphi,\sin\theta\sin\varphi,\cos\theta)$ can be obtained plugging (194) in (93), or written directly in terms of (197-198) using (94). Figure 9 shows surfaces of constant $h_{+,\times}$ as a function of $(\theta,\varphi)$, for a fixed value of $\psi=\pi$ and $e=0$. The quadrupolar nature is clearly visible.

For circular orbits ($e=0$), only the terms in $\cos(2\psi)$ and $\sin(2\psi)$ remain, and furthermore the angular velocity is constant, in particular $\psi=\omega_{0}t$ from (182). Thus for circular orbits, the GW angular frequency is twice the orbital frequency:

[TABLE]

Figure 10 shows the waveforms as a function of retarded time, in units of $T$, for 2.5 orbital periods. For elliptical orbits additional frequencies are present, both larger and smaller than $\omega_{0}$. In fact since the angular velocity is not constant, infinitely many harmonics are emitted. These are shown in Fig. 11.

In the quadrupole approximation, the energy and angular momentum fluxes (172) and (175) read

[TABLE]

(only the $z$-component of angular moment is relevant since the binary is in the $xy$-plane ). Substituting Eq. (196) gives

[TABLE]

where for circular orbits the constant rates of emission are given by

[TABLE]

in terms of the dimensionless coefficient $Gm/c^{2}p$. The above expressions are valid for all $e\geq 0$ provided $\vec{v}^{2}\ll c^{2}$. The power emitted is constant for circular orbits, and maximal at periastron for elliptical orbits. Figure 10 shows the waveforms Eq. (197)-(198) and power emitted Eq. (203) over 2.5 periods of a circular orbit with $e=0$, in unit of $t_{R}/T$ where $T$ is the orbital period. Figure 11 shows the same for an elliptical orbit with $e=0.3$. In both cases the periodic motion is clear. Over longer time-scales $t\gg T$, however, the emission of energy and angular momentum backreact on the orbital trajectories and must be considered. For a hyperbolic orbit, the corresponding plots are given in figure 12. The motion is obviously no-longer periodic and simply amounts to a fly-by: thus back-reaction effects do not accumulate over time and will be less significant (see subsection VII.6).

We now evaluate the effect of energy and angular momentum dissipation on the system’s dynamics in order to study the back-reaction.

VII.3 Back-reaction and waveform: circular orbits

The effect of the gravitational wave emission produces an effect that modifies the dynamics of the system, known as radiation-reaction force. As explained at the end of the previous Section, this effect can be deduced equating the non-conservation of the energy and angular momentum with the corresponding gravitational wave fluxes, see (177), and solving to find how the orbital element change in time. In order to do so we have to replace the constant orbital elements with functions of time, a method already used to compute perturbations to Keplerian orbits and known as ‘osculating orbits’. For circular orbits, there is only one independent quantity, it is thus sufficient to look at the energy equation. On the source side, we have (183)

[TABLE]

On the gravitational wave side, we have (205). Equating the two gives

[TABLE]

and this tells us how the orbit evolves. Using Kepler’s law (188), we can obtain the evolution of the angular frequency as

[TABLE]

where ${\cal{M}}$ is the chirp mass (6). This in turns tells us the change in the frequency of the emitted waves, which as seen in (200) is twice the orbital frequency:

[TABLE]

or equivalently in terms of $f=\omega/2\pi$,

[TABLE]

This shows that as gravitational waves are emitted, the orbit decays, the angular frequency increases, and so does the GW’s frequencies as well, leading to an even greater orbital decay. This run-away effect stops when the orbit decays completely and the two body coalesce. Clearly before that happens, higher order effects become important and one should improve the calculation. It is however instructive to get the lowest-order approximation and a qualitative overall picture to assume that the quadrupole approximation is valid throughout the evolution. We can then integrate between an initial time $t$ and the coalescing time $t_{c}$, where the frequency formally diverges. This gives

[TABLE]

Thus for a binary inspiral on a circular orbit $f^{-8/3}$ is linear in time, with a slope which determines directly the chirp mass. From this solution we can also immediately deduce the evolution of the orbital frequency $\omega_{0}$, the orbital radius $p$, and the GW amplitude (199).

To get the complete waveform, we only need the time dependence of $\psi$. This can be computed observing from (182) that $\dot{\psi}=\omega_{0}$, therefore

[TABLE]

where $\Phi_{c}$ is the phase at coalescing time. Combining these results, we get for the polarizations (197-198) in the direction perpendicular to the orbital plane

[TABLE]

For circular orbits, it is also easy to give the polarizations in an arbitrary direction $\vec{N}=(\sin\theta\cos\varphi,\sin\theta\sin\varphi,\cos\theta)$, using Eq. (94):

[TABLE]

As described after (94), the angle $\theta$ can be identified with the inclination $\iota$ of the source relative to the detector, see Fig. 1, whereas $2\omega=-2\varphi+\pi$ is the longitude of pericenter. Notice also that $\omega$ can be absorbed into a redefinition of the initial time. This is the analytic expression of the curve plotted in blue in figure 13.

VII.4 Back-reaction: Elliptical orbits

The case of elliptic orbits is more intricate, and offers a few interesting insight into the dynamics: there is a full spectrum of emission, and not a monochromatic one, which shows up in a non-symmetric waveform; there is loss of both energy and angular momentum, and at different rates, which shows up in the orbital back-reaction losing eccentricity faster than it decays; GW emission is not constant during the orbit, but stronger at periastron.

To analyse the system, it is convenient to distinguish the effects on two different time scales: short-time effects, namely the variations within a single orbit; and long-time, or secular, effects, namely the cumulative changes over many orbits. For instance, the dependence of the power emitted on the position $\psi$ illustrated in Fig.11 is a short-time effect. The orbital decay on the other hand is a secular effect. To consider secular effects, we introduce the average over one orbital period:

[TABLE]

where $\dot{\psi}$ is given in Eq. (182). We then replace (177) with their time-averages,

[TABLE]

The aim of this subsection is to solve these equations to determine the secular evolution of $e(t)$, $p(t)$, and thus $h_{+,\times}(t)$ with backreaction included.

Substituting (203) and integrating gives the Peter-Mathews formula Peters1963 :

[TABLE]

(This expression is only valid for $e<1$ as we are dealing with elliptical orbits.) Keeping $p$ constant, the radiation increases from $e=0$, to a maximum at $e\sim 0.5$ before decreasing and vanishing at $e=1$. The averaged angular momentum radiation is similarly determined using (204) and gives

[TABLE]

We now return to Eqs. (219), where on the left hand side the time-dependence is in $e(t)$ and $p(t)$. By definition, see (184), $L=\nu\sqrt{Gmp}$ from which

[TABLE]

This combined with Eqs. (219) and Eq. (221) gives

[TABLE]

The energy of the orbit is given in Eq. (183), from which

[TABLE]

Then plugging in (220) and (223) gives

[TABLE]

These coupled equations Eq. (223) and (225) can be solved using hypergeometric function to get $e(t)$ and $p(t)$, or alternatively combined to determine $p(e)$.

Observe that both $p(t)$ and $e(t)$ decrease with time. An elliptical orbit with initial eccentricity $e\neq 0$ will thus become more circular due to GW radiation. Whereas an initially circular orbit with $e=0$ remains circular for all times. This is the reason why often it is a good approximation to consider circular orbits, particularly when studying the last moments before the merger of the binary system. (This is the case of the events observed by LVK.) At first order in $e$, we can approximate (225) with

[TABLE]

where $T$ is the orbital period (187). This approximation provides a time-scale, $\tau_{R}$, of the radiative decay of $e$ and $p$. From Eq. (186), $\left(\frac{Gm}{c^{2}p}\right)\sim|\vec{v}|^{2}/c^{2}\ll 1$, and thus $\tau_{R}\sim(c/v)^{-5}T\gg T$.

The decrease of $p$ and $e$ also implies that $T$ decreases with time. Indeed from Eqs. (187), (225) and (223)

[TABLE]

where ${\cal{M}}$ is the chirp mass (6). Thus the orbital frequency $\omega_{0}$ increases, and GWs are emitted with increasing frequencies. The amplitude also increases, since has we can see from (199) it is inversely proportional to $p$.

To have the quantitative behaviour of the waveforms $h_{+,\times}$, we also need the time evolution of $\psi$. This is obtained solving (182) with $e(t)$ and $p(t)$ the solutions of Eqs.(225) and (223), and can be done numerically. The result is plotted in Fig. 13, for 4 different initial values of the eccentricity $e=0,0.3,0.5$ and $0.7$, for the plus polarization in the $z$ direction (197). The upper waveform in Fig. 13 is for circular orbits. The increasing amplitude and frequency of the GWs is clearly visibile and will be quantified in the discussion below. The waveform diverges when $p$ reaches zero, though clearly this is beyond the regime of applicability of the quadrupole approximation which assumes $|\vec{v}|/c\ll 1$. Since $|\vec{v}|\sim 1/\sqrt{p}$ this is clearly violated as $p\rightarrow 0$. That is the reason why, in Section I.4.2, we invoked the ISCO as a possible minimum distance, which then defined a merger frequency through Eq. (8). The elliptic orbits are not symmetric, and the amplitude is maximal at the periastron, where we also have the maximum emitted power. The different plot also show that the waveform diverges at earlier and earlier times as $e$ increases, in agreement with the fact that the emitted powers increase with $e$.

To get also a geometric intuition about the dynamics of the system during evolution, we plot in figure 14 the orbits for $e=0$ and $e=0.3$, showing the decrease in orbital radius and eccentricity.

VII.5 Frequency content of elliptical orbits

As a side result, let us come back to elliptic orbits, without dissipation, and discuss how one can study the frequency content of the waves emitted by computing the contributions at each frequency to the total power (201). To that end, we expand the dynamics into Fourier modes multiples of the fundamental harmonics given by the orbital frequency $\omega_{0}=\sqrt{Gm/p^{3}}$. So for instance we write

[TABLE]

with the dimension-full factor in front so to have dimension-less Fourier components

[TABLE]

Recall that $I_{ab}$ is given by (192) with $\vec{r}$ parametrized as in (179). To solve this integral one then needs the explicit time dependence of $r$ and $\psi$, namely the solution of the (unperturbed) Keplerian orbit of ellipticity $e$. Details can be found in e.g. Maggiore1 , and one finds

[TABLE]

and

[TABLE]

in terms of the Bessel’s functions $J$. Taking three time derivatives leads to factors of $n^{3}\omega_{0}^{3}$. Substituting into (201) and averaging over one period to calculate $\langle P_{\rm GW}\rangle$ leads to terms such as $\langle\sin n\omega_{0}t\sin m\omega_{0}t\rangle\sim\delta_{mn}$ meaning that the different harmonics do not interfere. In conclusion one finds

[TABLE]

where

[TABLE]

where $P_{\rm{GW}}^{e=0}$ is given in Eq. (205). The $\langle P_{n}\rangle$ are plotted in figure 15 for different values of $e$. For circular orbits all power is emitted via the $n=2$ mode. The power is the more and more distributed among other harmonics as $e$ increases, and the frequency at which the maximum power is radiated also increases with $e$.

VII.6 Hyperbolic orbits

While all GW detections to date are from bound elliptical/circular CBCs with $e<1$, many other potential GW sources exist, for instance non-spherical spinning NSs and asymmetric core collapse Supernovae. In this brief subsection we discuss another possible source, namely unbound binary systems on hyperbolic orbits. That is, we consider cases in which the eccentricity $e>1$ see Eq. (183) and Fig. 16.

Hyperbolic orbits are interesting not only because unbound orbits are expected to exist in nature (and hence the waveform for such events is and will be searched for by GW detectors NANOGrav:2023vfo ; Gasparotto:2023fcg ; Goncharov:2023woe ; Inchauspe:2024ibs ), but also because this simple system provides a first example of a gravitational wave memory effect. There are many different kinds of memory effects (see e.g. Favata:2010zu ), the simplest of which is the linear memory effect which occurs already at the lowest order in the PN expansion and which is illustrated by hyperbolic orbits. Memory effects occur when there is a permanent change $\Delta h_{ab}^{\rm{TT}}$ in the gravitational waveform, and thus leads to permanent displacement $\Delta L$ of the arms of GW detector for example, see Eq. (115).

VII.6.1 Linear memory effect and low-frequency GWs

The two GW polarisations for hyperbolic orbits are shown in Fig. 12 for $\vec{N}=(0,0,1)$. The waveform is not periodic but rather burst-like, and $h_{\times}$ has a non-zero variation between $t=\pm\infty$: this is the linear memory effect. More generally, the variation of the metric perturbation between $t=\pm\infty$ is given by

[TABLE]

where we introduced the shorthand notation $h_{cd}^{\text{TT},\pm}$ for the evaluation at $t=\pm\infty$. It follows from (155) that there will be a linear memory effect $\Delta h^{\text{TT}}_{ij}\neq 0$ when there is a net change in the second time derivatives of quadrupole moments of the system. This is precisely the case for hyperbolic orbits, since at initial and final times the accelerations are zero, but the velocity vectors are different. At lowest PN order,

[TABLE]

where $v_{\infty}$ is given in (191) and the unit-norm $n_{a}$ refers to incoming and outgoing directions. From this one can compute the polarizations, for instance in the $\hat{z}$ direction,

[TABLE]

See Braginsky:1987kwo for more details.

The time-scale over which the GW signal varies in Fig 12, namely the burst time-scale, is determined by the inverse of the characteristic frequency ${\omega_{c}^{-1}}$ given in Eq. (190). Up to factors of eccentricity, this also determines the characteristic frequency scale of the emitted GWs on hyperbolic orbits. Indeed, contrary to the case of periodic elliptical orbits discussed in section VII.5, GWs of all continous frequencies are emitted, and one can determine the GW power as a function of frequency by now Fourier transforming the power emitted. Using the convention $\tilde{I}_{ab}(\omega)=\int dtI_{ab}(t)e^{-i\omega t}$, as well as the quadrupole approximation, the total energy emitted in GWs is

[TABLE]

The emitted power in GWs is thus

[TABLE]

Direct calculation analogous to that of section VII.5 (see e.g. Brax:2024myc ) shows that ${\cal{P}}_{\rm GW}(\omega)$ is peaked at a value fixed by $\omega_{c}$ but which increases with $e$, see figure 17. Notice that ${\cal{P}}_{\rm GW}(0)\neq 0$. This is due to the linear memory effect: indeed Eq. (237), written in Fourier space reads

[TABLE]

thus a non-vanishing linear memory effect implies ${\cal{P}}_{\rm GW}(0)\neq 0$.

VII.6.2 Capture due to GW emission

One can estimate the energy emitted in GWs between $\psi_{-}<\psi<0$ by calculating

[TABLE]

where $P_{\rm{GW}}(\psi)$ is given in Eq. (203). Let us consider an orbit which is only slightly unbound, thus $e=1+\epsilon$ with $0<\epsilon\ll 1$ so that the orbital energy is $E\simeq\nu Gm^{2}\epsilon/p$, see Eq. (183). Then it is straightforward to determine $\Delta E_{\rm{GW}}$ to lowest order in $\epsilon$, leading to

[TABLE]

Note that here we assumed $p$ constant, that is we have neglected the backreaction of the emitted GR on $p$: this is a reasonable approximation for this individual burst process since $p$ changes on a time-scale $\sim(c/v)^{-5}\omega_{c}^{-1}$ as discussed previously. The orbital energy $E\sim\nu Gm^{2}\epsilon/p$ will thus be reduced by $\Delta E_{\rm{GW}}$, and if

[TABLE]

then the GW energy loss will convert the hyperbolic orbit into a bound orbit before the pericenter. This can be rewritten as the condition

[TABLE]

Instead of parametrising the orbit in terms of $(p,\epsilon)$, for such scattering trajectories it is more convenient to work with $(b,v_{\infty})$ where $b$ is the impact parameter and $v_{\infty}$ is the orbital velocity at infinite separation. From Eq. (191), $v^{2}_{\infty}\simeq 2(Gm/p)\epsilon$ and substituting $\epsilon$ from Eq. (245) gives

[TABLE]

Furthermore, basic trigonometry gives

[TABLE]

Saturating the bound in (246) identifies the impact parameter $b_{\rm capture}$ for which capture occurs. The corresponding capture cross-section $\sigma_{{\rm GW}}:=\pi b_{\rm capture}^{2}$ is thus given by

[TABLE]

Such a process could play an important role for instance in dense star clusters and galactic nuclei, see e.g. Capozziello:2008ra ; OLeary:2008myb ; Hoang:2020gsi .

VIII GWs in curved space-time, cosmology

In the previous sections we discussed GWs in Minkowski space; our aim is now to generalise the results presented there to a cosmological space-time. There is now a further scale of interest other than the characteristic size $d$ of the source and the GW wavelength $\lambda_{\rm{GW}}$, namely the cosmological horizon. The results of the previous sections are valid in the so-called the local wavezone of the source, namely at distances scales $R$ which are large compared to the GW wavelength but small compared to the cosmological horizon, $d\ll\lambda_{\rm GW}\ll R\ll{\rm horizon}$. We now aim to extend them to cosmological scales: we will see that the expansion of the universe dampen the GW amplitude, and redshift frequencies and masses.

VIII.1 General background metric

VIII.1.1 Linearised equations

We now consider the linearized equations around an arbitrary curved background $\bar{g}_{\mu\nu}$ satisfiying $\bar{G}_{{\mu\nu}}=(8\pi G/c^{4})\bar{T}_{{\mu\nu}}$. Thus as in (36) $g_{\mu\nu}=\bar{g}_{\mu\nu}+h_{\mu\nu}$, and similarly the source is perturbed as ${T}_{{\mu\nu}}=\bar{T}_{{\mu\nu}}+T^{\scriptscriptstyle\rm(1)}_{\mu\nu}$. The linearised equations are derived explicitly in Appendix A, and they can simplified using the de Donder gauge (74), generalised to curved spacetime as303030We are using a bar to indicate both the background metric, and the trace-reversed perturbation. We believe that no confusion should arise, because of the different places where these quantities enter.

[TABLE]

In this gauge, we have (see Eq. (321) for details)

[TABLE]

Here $\bar{R}_{\mu\alpha\nu\beta}$ is the Riemann curvature tensor of the background space-time, and $\bar{R}_{\mu\nu}$ the Ricci-tensor. The terms containing the background curvature depend on the ratio of $\lambda_{\rm{GW}}$ relative to the scale of variation of the background metric. In the following we will assume that they are negligible. The justification for this is that below we will restrict attention to the homogeneous and isotropic FLRW metric, where the scale of variation of the background is the cosmological horizon, which is much greater than $\lambda_{\rm{GW}}$.313131In a perturbed FLRW metric, there could be local bumps in the curvature of scale similar to $\lambda_{\rm{GW}}$: we do not consider this case here. See also Isaacson:1968hbi for more details. The right-hand side has generically two contributions, one from the metric perturbations, and one from the matter perturbation. In cosmology, it is often useful to approximate matter using a perfect fluid description, for which

[TABLE]

where $\rho$ and $p$ are respectively the energy and pressure densities of the fluid. In this case, $T^{\mu}{}_{\nu}$ is metric-independent, and this allows us to separate the metric and matter perturbations using the identity $T^{\scriptscriptstyle\rm(1)}_{\mu\nu}=\bar{g}_{\mu\rho}T^{{\scriptscriptstyle\rm(1)}\rho}{}_{\nu}+h_{\mu}{}^{\rho}\bar{T}_{\rho\nu}$. On-shell, the second term is proportional to the background curvature and therefore can be neglected again by the same argument as given above. We conclude that in the absence of matter perturbations, (250) reduces to

[TABLE]

VIII.1.2 WKB approximation

For a general spacetime, the solution of Eq. (252) can be obtained in the WKB approximation. The underlying physical assumption is that the amplitude of the wave is slowly varying with respect to the frequency of the wave, and hence we write

[TABLE]

where the sum is over polarisations $p$ with polarisation tensor $\epsilon_{\mu\nu}$ satisfying $\epsilon_{\mu\nu}\epsilon^{\mu\nu}=1$ (recall that to this leading order in $h$, indices are raised and lowered with the background metric) and $A$ is the corresponding amplitude. The parameter $\delta\rightarrow 0$, and we define

[TABLE]

Now substituting (253) the de Donder condition Eq. 249 becomes (dropping the $p$ sum for simplicity)

[TABLE]

which to leading order in $\delta$ implies

[TABLE]

and hence that $\epsilon_{\mu\nu}$ is the transverse polarisation tensor. Substitution into the equation of motion Eq. (252) leads to terms in $\delta^{-2}$, $\delta^{-1}$ which are, respectively

[TABLE]

The first equation, the Eikonal equation (geometric optics limit), implies that GWs are massless with dispersion relation $\omega^{2}=\vec{k}^{2}$ and propagate on null geodesics. Contracting the second equation (258) with $\epsilon^{{\alpha\beta}}$, and using that $\epsilon_{\mu\nu}\epsilon^{\mu\nu}=1$ so that $(\bar{\nabla}^{\alpha}\epsilon_{\mu\nu})\epsilon^{\mu\nu}=0$, leads to

[TABLE]

This gives the decay of the GW amplitude $A$ along the null geodesics. Finally, substituted back Eq. (259) into Eq. (258) gives

[TABLE]

which implies that the polarization tensor $\epsilon_{\alpha\beta}$ of the GW is parallel propagated along the null geodesics.

To summarise, the solution of

[TABLE]

in the WBK approximation, and in the Lorenz gauge, is $\bar{h}_{\mu\nu}(x)=\real\left[A(x)\epsilon_{\mu\nu}(x)e^{ik_{\mu}x^{\mu}}\right]$ with

[TABLE]

We now consider these equations in a FLRW metric.

VIII.2 FLRW metric: background

The flat Friedmann-Lemaitre-Robertson-Walker (FLRW) metric is

[TABLE]

where $\vec{x}$ are comoving coordinates and $a(t)$ is the scale factor (normalised such that today, at $t=t_{0}$, $a(t_{0})\equiv a_{0}=1$). For a perfect fluid source, Einsteins equations Eq. (16) reduce to the Friedmann equations

[TABLE]

where $\rho$ and $P$ are respectively the energy density and pressure of the perfect fluid, and $H=\frac{\dot{a}}{a}$ is the Hubble parameter, whose value today is the Hubble constant $H_{0}$. These two equations imply the conservation equation $\dot{\rho}+3H(\rho+P)=0=\nabla_{\nu}T^{\mu\nu}$. In terms of conformal time $\eta$ defined by $d\eta=dt/a(t)$, the metric in Eq. (266) is conformally related to the Minkowski metric

[TABLE]

Consider now a source (of photons or GWs) at fixed radial position $r=0$, and an observer at $r_{o}$. On a $(\eta,r)$ space-time diagram null radial geodesics propagate at 45 degrees. If two null geodesics are emitted at a conformal time interval $\delta\eta_{s}$ by the source, then they arrive at the observer with the same conformal time interval $\delta\eta_{o}=\delta\eta_{s}$. This implies the standard time-dilation relation

[TABLE]

where $z$ is the redshift of the source (the observer is at $t_{o}=t_{0}$), and equivalently that the emitted (or ‘source’) frequency $f_{s}$ is related to the observed frequency $f_{o}$ by

[TABLE]

The radial comoving distance $R$ to an event with redshift $z$ is given by solving $ds^{2}=0$, thus $dr=dt/a(t)$, leading to

[TABLE]

where $H(z)=H_{0}E(z)$ is the Hubble parameter expressed in terms of redshift, and from the Friedmann equation

[TABLE]

where $\Omega_{r,m}=\frac{8\pi G\rho_{r,m}}{3H_{0}^{2}}$, $\Omega_{\Lambda}=\frac{\Lambda}{3H_{0}^{2}}$ and $\Omega_{r}+\Omega_{m}+\Omega_{\Lambda}=1$.

A crucial quantity is the luminosity distance $d_{L}(z)$. This relates the EM luminosity of the source and the luminosity measured by the observer. In the flat FRWL metric (266) it is given by

[TABLE]

As we will see, this same distance scale determines the GW amplitude in an expanding universe.

VIII.3 FLRW metric: gravitational waveforms

Consider a GW propagating radially outwards from the source at $r=0$ and redshift $z_{s}$ with $k_{\mu}=\omega(1,-1,0,0)$. From Eq. (259) it is possible to determine how the GW amplitude decreases along the null GW geodesic. In a FLRW metric (in conformal time) Eq. (259) becomes

[TABLE]

Thus $A(\eta,r)a(\eta)r$ remains conserved during the propagation, and

[TABLE]

The constant is fixed by the known amplitude of the wave in the wave-zone approximation, close to the source, where the Minkowski results are valid. Then the remainder of the solution $\bar{h}_{\mu\nu}$ is obtained by parallel transporting this solution from the source to the observer. We now carry out these steps.

Before doing so we note that in a flat FLRW universe and focusing on the spatial TT components only, then in fact (250) reduces to

[TABLE]

since the spatial components of the Riemann and Ricci tensors vanish identically. From here the scaling of the amplitude of GWs as $1/a(\eta)$ is also immediate.

We now consider a compact binary system on circular orbits, as discussed in section VII.3, but this time in a FLRW background. In the wave-zone approximation and at a physical distance $R=a(t_{s})r$ from the source as measured by time $t_{s}$ of the source clock, the plus and cross polarisations of the GW are given in (216) and (217). Focusing on the cross polarisation,

[TABLE]

where $t^{\rm ret}_{s}=t_{s}-t_{c}$ is the time to coalescence at $t_{c}$ and ${\cal M}$ is the chirp mass (6). The time dependence of the frequency is given in Eq. (211) namely

[TABLE]

leading to

[TABLE]

so that the phase dependence is

[TABLE]

We now parallel transport this solution (278) along a null geodesic to the observer. Along the geodesic the GW phase remains constant because the time dilation effects cancel the redshifting of the frequency. Thus at the observer whose clock measures $dt=dt_{s}(1+z)$, the observed GW frequency $f=f_{s}/(1+z)$ leading to $\Phi(t)=\Phi_{s}(t_{s})$. However, at the observer, the GW amplitude is changed. From Eq. (276), and using (278)

[TABLE]

where we have included the redshifting of frequency. Let us now define the redshifted chirp mass

[TABLE]

Then (282) becomes

[TABLE]

where in the second line we have used Eq. (274) defining the luminosity distance to the source (today $a(t_{o})=1$). The dependence of the observed frequency on time $t$ is obtained by $f_{s}=(1+z)f$ into Eq. (279):

[TABLE]

Assuming that changes in $z$ are negligible during the observation time, then $z$ can be taken as constant323232See Bonvin:2016qxr for a discussion of where this assumption may lead to biases leading to

[TABLE]

namely the GW phase depends on the redshifted chirp mass,

[TABLE]

To summarize, the GW frequency depends on the redshifted chirp mass ${\cal M}_{z}$ which is therefore determined by measurements of the phase of an inspiral signal. The GW amplitude depends on both ${\cal M}_{z}$ and $d_{L}(z)$. Given that the former is determined from the phase, measurements of the amplitude of the signal determine $d_{L}(z)$. Generally speaking therefore, GW observations from individual CBC events determine the luminosity distance $d_{L}(z)$ and the so-called ‘redshifted’ masses,

[TABLE]

which are related to the ‘source’ masses $m_{1,2}$ by the same factor of $1+z$ as in Eq. (283).333333To determine each redshifted mass individually, rather than in the combination of the chirp mass, requires the waveform beyond the lowest order quadrupolar form discussed here

It is important to note that while the redshifted masses and the luminosity distance can be deduced from the waveform, the redshift $z$ of the individual CBC remains undetermined. To deduce this redshift, one possibility is to assume a cosmological model — such as $\Lambda$CDM, with given values of $H_{0}$, $\Omega_{m}$ etc, say from the Planck observations — so that $z$ can be read off the luminosity distance, see Eq. (274). Another possibility is to find, and measure, effects in the waveform that also depend on the source-frame masses as well as the redshifted masses. For binary black holes there are no such effects (even to higher orders in the PN expansion). For binaries including neutron stars, tidal effects with these properties enter, however at higher PN order (see e.g. DelPozzo:2015bna ; Messenger:2011gi and references within). Furthermore, the features depend on the — uncertain — equation of state of the nuclear matter making up the neutron star. In conclusion, to lowest order in the PN expansion for CBCs, it is not possible to determine the redshift $z$ of the source from GW observations alone: there is a perfect degeneracy between source masses, redshift, as well as spins.

Note that if $z$ is determined, then from the detected masses $m^{\rm detected}_{1,2}$ one can obtain the value of the source masses $m_{1,2}$ via Eq. (289). This was done in LIGOScientific:2016aoc , for example, assuming $\Lambda$CDM cosmological model with Planck values of $H_{0}$, $\Omega_{m}$ to find $z$: this is how the source-frame values of the two black hole masses was determined. However, the values of the cosmological parameters $H_{0}$ etc are in fact not known precisely, and a source of tension in cosmology today, see e.g. DiValentino:2021izs ; Verde:2023lmm . For these reasons, it can be interesting to use GW observations in a different way, namely as a new observable with which to measure cosmological parameters.

VIII.4 Measuring cosmological parameters with GWs: outline of GW cosmology

In this last brief subsection we outline how GW observations can be used to measure cosmological parameters.

The luminosity distance $d_{L}(z)$ in $\Lambda$CDM is given in Eq. (274), and is a function of cosmological parameters such as $H_{0}$ and $\Omega_{m}$. At low redshifts $z\ll 1$, the domain of the O3 measurements of the LVK collaboration LIGOScientific:2021aug , only $H_{0}$ enters since Eq. (274) reduces to

[TABLE]

Clearly, in order to measure $H_{0}$, not only is $d_{L}$ required (and obtained from GW observations, as we have discussed), but the redshift $z$ of the source is also needed. However, as mentioned above, this cannot be determined from GW observations: extra non-gravitational information is necessary to determine $z$. Such information could, for example, be electro-magnetic (EM).

Indeed, the most straightforward way to determine $z$ is use EM observations to uniquely identify the “host galaxy” of the GW signal, namely the galaxy in which GW event occurred. This was possible for GW event GW170817 which occurred on August 17th 2017 and which corresponded to the the merger of two neutron stars GW170817 . This GW signal was observed the two LIGO and Virgo GWdetectors, and 1.7s following the GW merger, EM observers around the globe observed a subsequent gamma-ray burst as well as multiple EM signals in different frequency bands. Using this EM data it was possible to determine the host galaxy, namely NGC 4993, a galaxy in the Hydra constellation. This constellation is receding from us with a velocity $cz=3327\pm 72$ km/s, due to the expansion of the universe. Combining this with the distance $d_{L}=43.8^{+2.9}_{-6.9}$ Mpc inferred from the GW signal led, using Eq. (290), to an estimated value for the Hubble constant of $H_{0}=70^{+12}_{-8}$ km/s/Mpc LIGOScientific:2017adf . This result, using one GW event only, is consistent with other measurements but is of course less accurate because of its larger error bars. Its interest is that it shows that the idea works. The errors would be reduced (with a $\sim 1/\sqrt{N}$ scaling) if $N$ other measurements of this kind existed, but unfortunately GW170817 was an extremely rare event as since then no further GW events with associated EM counterparts (known as standard sirens) have been detected.

However, LVK has detected GWs from hundreds of BBHs and a few NS-BH, for each of which there is a measured $d_{L}$ and $m_{1,2}^{\rm detected}$ — but no EM counterpart. Even for these dark sirens, it is possible to obtain redshift information, and therefore measure $H_{0}$. Today, two pieces of information are used together to get a statistical redshift for GW events: (i) galaxy catalogues and (ii) astrophysical modelling of the formation channels of BBHs. We refer the reader to Mastrogiovanni:2023emh for a review of these methods and results.

IX Acknowledgements

We are grateful to Abhay Ashtekar, Tom Bertheas, Cédric Deffayet, Nathalie Deruelle, Eanna Flanagan, Alejandro Perez, Pierre Piovesan, Eric Poisson, Syksy Räsänen, Carlo Rovelli, Robert Wald for discussions on the topics of this review, as well as multiple Master 2 and doctoral students whom we have had the pleasure to teach. We also thank ChatGPT for technical support (not for the calculations though). This work was supported in part by Perimeter Institute for Theoretical Physics. Research at Perimeter Institute is supported by the Government of Canada through the Department of Innovation, Science and Economic Development and by the Province of Ontario through the Ministry of Research, Innovation and Science. This work was also supported in part by CERN.

Appendix A Second order action for perturbations around any background solution

This appendix gives the details of the steps required to calculate the Einstein-Hilbert action of GR to second order in perturbations about a general background metric $\bar{g}_{\mu\nu}$:

[TABLE]

The starting point is

[TABLE]

where $S^{(i)}$ is of ${\cal{O}}(h^{i})$. The extremisation of the second order action $S^{(2)}$ with respect to $h_{\mu\nu}$ will give the linearized equations of motion for $h_{\mu\nu}$, namely those we want to calculate. The background metric $\bar{g}_{\mu\nu}$ satisfies the background Einstein equations which follow from $S^{(1)}$.

To find $S^{(2)}$ there are two main steps: calculating the Ricci scalar $R$ to second order, and then the determinant of the metric and hence $\sqrt{-g}$ to second order.

Perturbed Riemann tensor, Ricci tensor and scalar, Einstein tensor

$\bullet$ The inverse metric, or contravariant metric tensor corresponding to Eq. (291) is given at second order by

[TABLE]

where indices of $h_{\mu\nu}$ are raised and lowered with $\bar{g}_{\mu\nu}$.

$\bullet$ The perturbed Christoffel symbols are given by

[TABLE]

In the following we denote the covariant derivative with respect to $\bar{g}$ by $\bar{\nabla}$ (with of course $\bar{\nabla}_{\mu}\bar{g}_{\nu\alpha}=0$). Substitution of Eqs. (291) and (293) gives

[TABLE]

$\bullet$ Next we calculate the Riemann tensor. Let

[TABLE]

where, from (294),

[TABLE]

The definition of the Riemann tensor together with (297) gives

[TABLE]

On writing

[TABLE]

and using (298) one can read off the different orders of the Riemann tensor. To first order

[TABLE]

where we have used Eq. (295). The second order term follows from (300) and (298) and reads

[TABLE]

$\bullet$ The Ricci tensor is then obtained by contraction:

[TABLE]

From (302) it follows that

[TABLE]

whereas from (303)

[TABLE]

Its explicit form is not required below, but for completeness we give it here:

[TABLE]

$\bullet$ The Ricci scalar is obtained from Eqs. (293) and (304) and is given by

[TABLE]

where

[TABLE]

(in the last line we have used (305)), and

[TABLE]

$\bullet$ For completeness we also give the first order perturbed Einstein tensor

[TABLE]

Perturbed metric determinant

To expand $\sqrt{-g}$ to second order, we write Eq. (291) $g_{\alpha\beta}=\bar{g}_{\alpha\mu}(\delta^{\mu}_{\beta}+M^{\mu}_{\beta})$ where $M^{\mu}_{\beta}=\bar{g}^{\mu\lambda}h_{\lambda\beta}$. Thus

[TABLE]

where the matrix $\mathbf{M}$ has components $M^{\mu}_{\beta}$. To quadratic order

[TABLE]

Replacing ${\rm tr}{\mathbf{M}}=h$ and ${\rm tr}({\mathbf{M}}^{2})=h_{\mu\nu}h^{\mu\nu}$ in Eq. (311), and then taking the square root, gives

[TABLE]

EH action to second order

Substituting (307) and (312) into the perturbed Einstein Hilbert action Eq. (292) gives

[TABLE]

$\bullet$ The first order action is

[TABLE]

The first two terms are total derivatives. After integration by parts and dropping the boundary terms

[TABLE]

On including matter through the stress tensor, the variation of this gives the background Einstein equation.

$\bullet$ The second order part becomes, on substituting (309),

[TABLE]

where $R^{(1)}$ and $R^{(1)}_{\nu\sigma}$ are given in (308) and (305) respectively.

The first term $\int d^{4}x\sqrt{-\bar{g}}\bar{g}^{\nu\sigma}R^{(2)}_{\nu\sigma}$ splits into four parts on using (306). The first two parts are total derivatives and do not contribute. The last two parts give

[TABLE]

Collecting the expressions together into (314) gives

[TABLE]

Finally, after integration by parts,

[TABLE]

In flat space $\bar{g}_{\mu\nu}=\eta_{\mu\nu}$, the terms in square brackets all vanish and this reduces to the usual perturbed equations around Minkowski space.

In terms of the trace reversed perturbation

[TABLE]

this becomes

[TABLE]

Linearised equations of motion

$\bullet$ Equations of motion. Variation of action (314) with respect to $h^{{\alpha\beta}}$ gives

[TABLE]

leading to the equations of motion

[TABLE]

In terms of the trace reversed perturbation, after commuting covariant derivatives, for example,

[TABLE]

and then imposing the Lorenz gauge $\bar{\nabla}_{\mu}\bar{h}^{{\mu\nu}}=0$ these read

[TABLE]

which can be rewritten as

[TABLE]

where

[TABLE]

Including matter, and splitting its energy momentum tensor $T_{\mu\nu}=\bar{T}_{\mu\nu}+T^{\scriptscriptstyle\rm(1)}_{\mu\nu}$ into a background (determining $\bar{g}_{\mu\nu}$) plus first-order perturbation, the first-order field equations are

[TABLE]

Appendix B Green’s functions

Given a linear differential operator $\Delta$ acting on functions on $\mathbb{R}^{n}$, its corresponding Green’s function is a function on $\mathbb{R}^{n}\times\mathbb{R}^{n}$ satisfying

[TABLE]

Green’s functions are useful to solve differential equations in the presence of sources, since they allow one to write the solutions of $\Delta\Phi=J$ as

[TABLE]

where $\Phi^{\circ}$ is a solution of the homogeneous problem with $J=0$. For a given $\Delta$ the Green function is typically not unique, but a unique one can be selected via boundary conditions or other physical requirements.

For the Laplace equation, $n=3$ and $\Delta=\vec{\partial}^{2}$, there is a unique solution of (322) with vanishing boundary conditions at infinity, given by

[TABLE]

Therefore

[TABLE]

For the d’Alembert equation, $n=4$ and $\Delta=\square$, fixing vanishing boundary conditions at spatial infinity is not enough to have a single solution: there are infinitely many solutions of the homogeneous equation that can be added to any given $G$ and still satisfy the defining equation. Two notable examples are the retarded and advanced ones, which are uniquely characterised by vanishing for $x$ respectively in the past or the future of $x^{\prime}$, and given by

[TABLE]

where $\Theta$ is Heaviside’s step function. The retarded solution imposes no-incoming radiation boundary conditions, and it is the one relevant to study the emission of waves from a source.

Appendix C Landau-Lifshitz approach

The Landau-Lifshitz formulation of Einstein’s equations is a convenient approach to perturbative theory around Minkowski, and it is widely used by the community working on the PN and PM expansions. We provide here a brief description of the approach from the perspective of the main text, and refer the reader to Blanchet:2006jqj ; poisson for more details. In the Landau-Lifshitz approach one uses a density-weighted inverse metric as fundamental variable,

[TABLE]

The interest in doing so is that the quantity

[TABLE]

is related to the Einstein tensor via ordinary derivatives. More precisely,

[TABLE]

where $t_{\scriptscriptstyle\rm LL}^{\alpha\beta}$ is a pseudo-tensor of density-weight two, given explicitly by some lengthy expression in terms of first derivatives of the metric. Crucially, it is conserved on-shell $\partial_{\alpha}(-gt_{\scriptscriptstyle\rm LL}^{\alpha\beta})\,\hat{=}\,0$. The dimension-full numerical factor in front of it is included in its definition for convenience when going on-shell in the last equality above.

This arrangement of Einstein’s equations manifestly breaks covariance.343434It is like rearranging the covariant geodesic equation (112) as in (113), where neither side of the equation is covariant by itself. Not only we have partial derivatives as opposed to covariant derivatives, but also tensor densities, aka pseudo-tensors, appearing. Its usefulness is limited to situations in which there are regions of spacetimes that are approximately flat, and where one can choose a Cartesian coordinate system, so that partial derivatives can be interpreted. This is precisely the case when studying a perturbative approximation around Minkowski. Within that context, the reformulation has two useful advantages.

The first is that it provides a prescription for a gravitational energy-momentum pseudo-tensor valid to all orders, given by $t_{\scriptscriptstyle\rm LL}^{\alpha\beta}$. By analogy with the matter counterpart, one has a prescription to further split this quantity into contributions to energy, momentum and angular momentum. All these expression are gauge-dependent; in particular, we have the usual problem that the pseudo-tensor can be made to vanish at any given point, using a local inertial frame. But the logic in this approach is that we assume to have a preferred coordinate system, the Cartesian ones of the fiducial flat metric, and that is the gauge we stick to. Furthermore, since the left-hand side is a total derivative, the total energy, momentum and angular momentum can be expressed as surface integrals. This provides a prescription for these quantities that can be evaluated in a region far from the sources, where one can safely assume that spacetime is approximately flat and use Cartesian coordinates, and computed to all order in perturbation theory.

The second advantage of the formulation is that it allows to set up an iteration scheme for the perturbative resolution of the field equations in a very practical way. To see that, we change variables to $\hat{h}^{\mu\nu}:=\eta^{\mu\nu}-{\mathfrak{g}}^{\mu\nu}$, where $\eta^{\mu\nu}$ is a fiducial background metric for which the coordinates are Cartesian. Then (329) are equivalent to

[TABLE]

where $\square$ is the flat spacetime d’Alembertian, and $t_{\scriptscriptstyle\rm H}^{\mu\nu}$ are the remaining terms in the LHS of (329). Specifically, $t_{\scriptscriptstyle\rm H}^{\mu\nu}$ is quadratic in $h^{\mu\nu}$ and satisfies $\partial_{\mu}(-gt_{\scriptscriptstyle\rm H}^{\mu\nu})\equiv 0$, whereas $t_{\scriptscriptstyle\rm NH}^{\mu\nu}$ does not but contains only terms that vanish in harmonic gauge. It follows that in the harmonic gauge,

[TABLE]

and this equation includes the conservation of the matter energy-momentum tensor, namely the matter dynamics. We stress that (330) are the exact Einstein’s equations, no approximation has been done yet. It is the introduction of a fiducial flat background and a fixed choice of Cartesian coordinates on it, that allows us to rewrite the exact equations in the form of a flat spacetime wave equation with a complicated source satisfying the matter dynamics through (331). In particular, all non-linearities are recasted on the RHS of the equations. To solve the equations perturbatively, we assume that $h$ is small, and proceed iteratively as explained in (37). The idea that makes this approach particularly convenient is to first solve the ‘relaxed field equations’ (330) alone, in harmonic gauge, and afterwards impose the gauge consistency condition (331) on the matter dynamics. Notice that in this scheme, the $n$-th iteration of the gravitational potentials $h_{\mu\nu}$ are sourced by matter fields satisfying the equations of motion that use the $(n-1)$-th iteration of the gravitational field. See Blanchet:2006jqj ; poisson for further details.