Part V: Gravitational Waves

Gravitational waves are ripples in the fabric of spacetime, predicted by Einstein's general relativity in 1916 and directly detected a century later by the LIGO/Virgo collaboration. These waves are generated by the acceleration of massive objects—most powerfully by the inspiral and merger of compact binaries such as black holes and neutron stars. The detection of gravitational waves has opened an entirely new observational window on the universe, enabling tests of strong-field gravity, measurements of the Hubble constant, constraints on the neutron star equation of state, and the birth of multi-messenger astronomy.

1. Linearized Gravity and Gravitational Wave Generation

Gravitational waves emerge naturally from Einstein's field equations when we consider small perturbations about a flat (Minkowski) background spacetime. This linearized approach, valid in the weak-field regime, reveals that gravity propagates as a wave at the speed of light.

1.1 Weak-Field Metric Perturbation

We decompose the spacetime metric into the flat Minkowski background $\eta_{\mu\nu}$ plus a small perturbation $h_{\mu\nu}$:

$$g_{\mu\nu} = \eta_{\mu\nu} + h_{\mu\nu}, \quad |h_{\mu\nu}| \ll 1$$

Here $\eta_{\mu\nu} = \text{diag}(-1, +1, +1, +1)$ is the Minkowski metric in Cartesian coordinates. The perturbation $h_{\mu\nu}$ encodes the gravitational wave content. To first order in $h_{\mu\nu}$, the inverse metric is:

$$g^{\mu\nu} = \eta^{\mu\nu} - h^{\mu\nu} + \mathcal{O}(h^2)$$

The linearized Ricci tensor, computed to first order in $h_{\mu\nu}$, yields the linearized Einstein equations. Defining the trace-reversed perturbation:

$$\bar{h}_{\mu\nu} = h_{\mu\nu} - \frac{1}{2}\eta_{\mu\nu}h, \quad h = \eta^{\mu\nu}h_{\mu\nu}$$

and imposing the Lorenz gauge condition $\partial^\nu \bar{h}_{\mu\nu} = 0$, the linearized Einstein equations reduce to a wave equation:

$$\Box \bar{h}_{\mu\nu} = -\frac{16\pi G}{c^4} T_{\mu\nu}$$

where $\Box = -\frac{1}{c^2}\frac{\partial^2}{\partial t^2} + \nabla^2$ is the d'Alembertian operator. In vacuum ($T_{\mu\nu} = 0$), this gives a free wave equation: gravitational waves propagate at the speed of light $c$.

1.2 Transverse-Traceless (TT) Gauge and Polarizations

The residual gauge freedom in the Lorenz gauge allows us to further specialize to the transverse-traceless (TT) gauge. For a wave propagating in the $z$-direction, the TT gauge conditions impose:

Transverse: $h_{0\mu}^{TT} = 0$ and $h_{z\mu}^{TT} = 0$ (only spatial components perpendicular to propagation survive)
Traceless: $h^{TT} = \eta^{ij}h_{ij}^{TT} = 0$

This leaves exactly two independent degrees of freedom—the two polarization states. For a wave propagating along $\hat{z}$, the TT metric perturbation takes the form:

$$h_{\mu\nu}^{TT} = \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & h_+ & h_\times & 0 \\ 0 & h_\times & -h_+ & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} \cos\left[\omega\left(t - \frac{z}{c}\right)\right]$$

The two polarizations have distinct physical effects on a ring of freely falling test particles:

Plus polarization (h+): stretches along x while squeezing along y, then reverses
Cross polarization (h×): same pattern rotated by 45 degrees

A general gravitational wave is a superposition of these two polarizations. The geodesic deviation equation for nearby test particles separated by $\xi^j$ gives:

$$\frac{d^2\xi^j}{dt^2} = \frac{1}{2}\ddot{h}_{jk}^{TT}\xi^k$$

For a detector with arm length $L$, the fractional change in proper distance is $\delta L / L \sim h/2$, which directly motivates the design of interferometric detectors.

1.3 Quadrupole Formula for GW Emission

The leading-order gravitational wave emission comes from the time-varying mass quadrupole moment. Unlike electromagnetic radiation, which is dipolar at lowest order, gravitational radiation begins at the quadrupole level because mass conservation eliminates the monopole and momentum conservation eliminates the dipole contribution.

The mass quadrupole moment tensor is defined as:

$$I_{ij} = \int \rho(\mathbf{x}) x_i x_j \, d^3x$$

The reduced (trace-free) quadrupole moment is:

$$\mathcal{I}_{ij} = I_{ij} - \frac{1}{3}\delta_{ij} I_{kk}$$

The gravitational wave strain in the TT gauge at distance $r$ from the source is given by the quadrupole formula:

$$h_{ij}^{TT} = \frac{2G}{c^4 r}\ddot{\mathcal{I}}_{ij}^{TT}(t_{\rm ret})$$

where $t_{\rm ret} = t - r/c$ is the retarded time. The TT projection extracts the transverse-traceless part using the projection operator $\Lambda_{ij,kl}$ constructed from the unit vector $\hat{n}$ pointing from source to observer.

The key insight is that GW amplitude scales as the second time derivative of the quadrupole moment, so sources with rapidly changing mass distributions are the strongest emitters. The prefactor $G/c^4 \approx 8.26 \times 10^{-45}\;\text{s}^2\,\text{kg}^{-1}\,\text{m}^{-1}$ is extraordinarily small, which is why gravitational waves are so weak and why only the most violent astrophysical events produce detectable signals.

1.4 Gravitational Wave Luminosity

The total power radiated in gravitational waves is obtained by integrating the GW energy flux over a sphere surrounding the source. The result is the quadrupole luminosity formula:

$$L_{\rm GW} = \frac{G}{5c^5}\left\langle\dddot{\mathcal{I}}_{ij}\dddot{\mathcal{I}}^{ij}\right\rangle$$

where the angle brackets denote time averaging over several orbital periods. The factor $G/c^5 \approx 2.76 \times 10^{-53}\;\text{W}^{-1}$ sets the scale: to produce $1$ watt of gravitational wave power requires an enormous rate of change of the quadrupole moment.

For dimensional analysis, the characteristic GW luminosity scale is:

$$L_0 = \frac{c^5}{G} \approx 3.63 \times 10^{52}\;\text{W} \approx 10^{26}\;L_\odot$$

This is the maximum possible luminosity in general relativity—the Planck luminosity. During the peak of a binary black hole merger, the GW luminosity can briefly reach a few percent of this value, exceeding the combined electromagnetic luminosity of all stars in the observable universe.

1.5 Strain Amplitude Estimates

For a compact binary system with component masses $m_1$ and $m_2$ in a circular orbit at separation $a$, the characteristic GW strain amplitude at distance $r$ is:

$$h \sim \frac{4G^2 m_1 m_2}{c^4 a r} = \frac{4G}{c^4 r}\frac{\mu M}{a}$$

where $\mu = m_1 m_2 / M$ is the reduced mass and $M = m_1 + m_2$ is the total mass. Using Kepler's third law $\omega^2 = GM/a^3$ to eliminate the separation in favor of the orbital angular frequency $\omega = \pi f_{\rm GW}$ (where $f_{\rm GW} = 2f_{\rm orb}$ for circular orbits):

$$h = \frac{4}{r}\left(\frac{G\mathcal{M}}{c^2}\right)^{5/3}\left(\frac{\pi f_{\rm GW}}{c}\right)^{2/3}$$

where $\mathcal{M}$ is the chirp mass (defined below). For a binary black hole system with $\mathcal{M} \sim 30\;M_\odot$ at $r \sim 400$ Mpc emitting at $f_{\rm GW} \sim 100$ Hz:

$$h \sim 10^{-21}$$

This corresponds to a length change of about $10^{-18}$ meters across a 4 km arm—a thousand times smaller than a proton radius. Detecting such minuscule distortions is the extraordinary engineering feat achieved by LIGO.

2. Binary Inspiral

The inspiral of a compact binary system is the most well-understood source of gravitational waves, thanks to the separation of scales between the orbital dynamics and the radiation reaction timescale. The system slowly loses energy and angular momentum to gravitational wave emission, causing the orbit to shrink, the frequency to increase, and the amplitude to grow—producing the characteristic "chirp" signal.

2.1 Keplerian Binary: Energy and Frequency

For a Newtonian binary in a circular orbit with separation $a$, the orbital energy and angular frequency are:

$$E_{\rm orb} = -\frac{Gm_1 m_2}{2a} = -\frac{G\mu M}{2a}$$

$$\omega_{\rm orb} = \left(\frac{GM}{a^3}\right)^{1/2}, \quad f_{\rm GW} = 2f_{\rm orb} = \frac{\omega_{\rm orb}}{\pi}$$

The GW frequency is twice the orbital frequency because the quadrupole moment has a period of half the orbital period (the mass distribution looks the same after a half-orbit rotation). The angular momentum is:

$$J = \mu\sqrt{GMa} = \mu(GMa)^{1/2}$$

2.2 Peters Formula: Orbital Decay

Peters (1964) derived the rate of orbital decay due to gravitational wave emission for circular orbits. Equating the GW luminosity to the rate of energy loss:

$$\frac{dE}{dt} = -L_{\rm GW} = -\frac{32G^4}{5c^5}\frac{m_1^2 m_2^2(m_1+m_2)}{a^5}$$

Since $E = -G\mu M/(2a)$, we have $dE/dt = (G\mu M/2a^2)\dot{a}$, giving Peters' formula for the semi-major axis decay:

$$\frac{da}{dt} = -\frac{64G^3}{5c^5}\frac{m_1 m_2(m_1+m_2)}{a^3}$$

This can be integrated to give the time remaining until coalescence from an initial separation $a_0$:

$$t_{\rm merge} = \frac{5c^5 a_0^4}{256 G^3 m_1 m_2(m_1+m_2)}$$

For eccentric orbits, Peters showed that gravitational radiation circularizes the orbit—eccentricity decreases faster than the semi-major axis, so by the time the binary enters the LIGO band, the orbit is nearly circular. The generalized Peters formula includes eccentricity-dependent enhancement factors$f(e)$ and $g(e)$ that increase the energy and angular momentum loss rates for eccentric orbits.

2.3 Chirp Mass

The chirp mass is the particular combination of the component masses that determines the leading-order gravitational waveform during the inspiral. It is defined as:

$$\mathcal{M} = \frac{(m_1 m_2)^{3/5}}{(m_1+m_2)^{1/5}} = \mu^{3/5}M^{2/5} = M\eta^{3/5}$$

where $\eta = \mu/M = m_1 m_2/(m_1+m_2)^2$ is the symmetric mass ratio ($0 < \eta \leq 1/4$). The chirp mass is the best-determined parameter from the GW signal because it enters the phase evolution at lowest post-Newtonian order. From a detected waveform, $\mathcal{M}$ can typically be measured to better than 0.1% accuracy.

The individual masses require higher-order post-Newtonian corrections and are therefore measured with larger uncertainties. For equal masses ($m_1 = m_2 = m$), $\mathcal{M} = 2^{-1/5}m \approx 0.871\,m$.

2.4 Frequency Evolution

The GW frequency sweeps upward during inspiral. From Kepler's law and Peters' formula, the frequency evolution is:

$$\dot{f}_{\rm GW} = \frac{96}{5}\pi^{8/3}\left(\frac{G\mathcal{M}}{c^3}\right)^{5/3}f_{\rm GW}^{11/3}$$

This equation can be inverted to give the time to coalescence as a function of GW frequency:

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

For a $30+30\;M_\odot$ binary entering the LIGO band at $f = 20$ Hz, the time to merger is approximately 0.2 seconds. For a $1.4+1.4\;M_\odot$ neutron star binary, the signal lasts about 100 seconds in band from 20 Hz to merger. This frequency evolution directly encodes the chirp mass, enabling its precise measurement.

The GW phase evolution to leading (Newtonian) order is:

$$\Phi(t) = -2\left(\frac{\tau(t)}{5G\mathcal{M}/c^3}\right)^{5/8} + \Phi_0$$

Higher-order post-Newtonian corrections introduce dependence on the symmetric mass ratio $\eta$, enabling extraction of individual component masses from the waveform. The PN expansion is currently known to 3.5PN (fourth post-Newtonian) order in phase and 3PN order in amplitude.

2.5 Time to Coalescence

The Newtonian approximation gives the time remaining until coalescence from a given GW frequency:

$$\tau = 2.18\;\text{s}\;\left(\frac{1.21\;M_\odot}{\mathcal{M}}\right)^{5/3}\left(\frac{100\;\text{Hz}}{f_{\rm GW}}\right)^{8/3}$$

Representative values for systems entering the LIGO band at 10 Hz:

System	Chirp Mass (M&odot;)	Time in band (10 Hz)	Cycles in band
NS-NS (1.4+1.4)	1.22	~17 minutes	~15,000
BH-NS (10+1.4)	3.06	~100 seconds	~3,000
BH-BH (30+30)	26.1	~0.5 seconds	~10
BH-BH (50+50)	43.5	~0.1 seconds	~5

2.6 Hulse-Taylor Binary PSR B1913+16

The binary pulsar PSR B1913+16, discovered by Hulse and Taylor in 1974, provided the first indirect evidence for gravitational wave emission, earning them the 1993 Nobel Prize in Physics. This system consists of two neutron stars ($m_1 \approx 1.44\;M_\odot$, $m_2 \approx 1.39\;M_\odot$) in a highly eccentric orbit ($e \approx 0.617$) with period $P_b \approx 7.75$ hours.

General relativity predicts the orbital period decay rate due to gravitational wave emission:

$$\dot{P}_b^{\rm GR} = -\frac{192\pi}{5}\left(\frac{2\pi G\mathcal{M}}{c^3 P_b}\right)^{5/3}\frac{1}{(1-e^2)^{7/2}}\left(1 + \frac{73}{24}e^2 + \frac{37}{96}e^4\right)$$

The predicted value is $\dot{P}_b^{\rm GR} = -2.402531 \times 10^{-12}\;\text{s/s}$. After more than 40 years of precise timing observations, the measured value agrees with the GR prediction to better than 0.2%. The cumulative orbital phase shift exceeds 40 seconds. This system will merge in approximately 300 million years.

Additional binary pulsars (PSR J0737-3039, the double pulsar) have further confirmed GR predictions with even greater precision, testing multiple post-Keplerian parameters simultaneously.

Simulation: Binary Inspiral Chirp Waveform

This simulation computes the gravitational wave strain for an inspiraling compact binary using the leading-order (Newtonian) quadrupole approximation. The waveform shows the characteristic chirp: increasing frequency and amplitude as the binary spirals inward toward merger.

Binary Inspiral Chirp Waveform

Python

Compute GW strain h(t) for an inspiraling compact binary showing the chirp signal

Parameters

m₁ (M☉)

Primary mass

m₂ (M☉)

Secondary mass

Distance (Mpc)

Luminosity distance

f₀ (Hz)

Initial GW frequency

binary_inspiral_chirp.py129 lines

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# Physical constants
G = 6.67430e-11       # m^3 kg^-1 s^-2
c = 2.99792458e8      # m/s
M_sun = 1.98892e30    # kg
pc = 3.08567758e16    # m

# Binary parameters
m1_solar = M1_VAL     # solar masses
m2_solar = M2_VAL     # solar masses
d_Mpc = DIST_VAL      # distance in Mpc
f0 = F0_VAL           # initial GW frequency in Hz

m1 = m1_solar * M_sun
m2 = m2_solar * M_sun
M_total = m1 + m2
mu = m1 * m2 / M_total           # reduced mass
eta = mu / M_total                # symmetric mass ratio
M_chirp = M_total * eta**(3.0/5.0)  # chirp mass
d = d_Mpc * 1e6 * pc             # distance in meters

print(f"=== Binary Inspiral Chirp Waveform ===")
print(f"m1 = {m1_solar:.1f} M_sun, m2 = {m2_solar:.1f} M_sun")
print(f"Total mass M = {(m1_solar + m2_solar):.1f} M_sun")
print(f"Chirp mass Mc = {M_chirp/M_sun:.2f} M_sun")
print(f"Symmetric mass ratio eta = {eta:.4f}")
print(f"Distance = {d_Mpc:.1f} Mpc")
print(f"Initial GW frequency = {f0:.1f} Hz")

# Time to coalescence from initial frequency (Newtonian)
# tau = (5/256) * (c^3 / (G*Mc))^(5/3) * (pi*f0)^(-8/3)
Mc_geo = G * M_chirp / c**3  # chirp mass in seconds
tau_merge = (5.0 / 256.0) * (np.pi * f0)**(-8.0/3.0) * Mc_geo**(-5.0/3.0)
print(f"Time to coalescence = {tau_merge:.4f} s")

# Compute waveform using post-Newtonian inspiral approximation
# Frequency evolution: f(t) = (1/pi) * (5/(256*(t_c - t)))^(3/8) * Mc_geo^(-5/8)
# We integrate backwards from coalescence

N_points = 80000
t_end = tau_merge * 0.999  # stop just before merger

# Create time array
t = np.linspace(0, t_end, N_points)
t_to_merger = tau_merge - t

# Avoid division by zero
t_to_merger = np.maximum(t_to_merger, 1e-10)

# GW frequency as function of time
f_gw = (1.0 / np.pi) * (5.0 / (256.0 * t_to_merger))**(3.0/8.0) * Mc_geo**(-5.0/8.0)

# GW strain amplitude (leading order quadrupole)
# h = (4/d) * (G*Mc/c^2)^(5/3) * (pi*f/c)^(2/3)
h_amp = (4.0 / d) * (G * M_chirp / c**2)**(5.0/3.0) * (np.pi * f_gw / c)**(2.0/3.0)

# GW phase (integrated from frequency evolution)
# Phi(t) = -2 * (t_to_merger / (5 * Mc_geo))^(5/8) + Phi_0
Phi = -2.0 * (t_to_merger / (5.0 * Mc_geo))**(5.0/8.0)

# Plus polarization (face-on, iota=0)
h_plus = h_amp * np.cos(Phi)

# Cut off when frequency exceeds ISCO frequency
f_isco = c**3 / (6.0**(3.0/2.0) * np.pi * G * M_total)
print(f"ISCO frequency = {f_isco:.1f} Hz")
mask = f_gw <= f_isco
t_plot = t[mask]
h_plot = h_plus[mask]
f_plot = f_gw[mask]
h_amp_plot = h_amp[mask]

print(f"Peak strain amplitude = {np.max(h_amp_plot):.2e}")
print(f"Final GW frequency = {f_plot[-1]:.1f} Hz")

# Create figure
fig, axes = plt.subplots(3, 1, figsize=(12, 10), facecolor='#0f172a')

# Panel 1: GW strain waveform
ax1 = axes[0]
ax1.set_facecolor('#1e293b')
ax1.plot(t_plot, h_plot, color='#38bdf8', linewidth=0.3, alpha=0.9)
ax1.fill_between(t_plot, -h_amp_plot, h_amp_plot, alpha=0.15, color='#38bdf8')
ax1.set_xlabel('Time (s)', color='#cbd5e1', fontsize=11)
ax1.set_ylabel('Strain h(t)', color='#cbd5e1', fontsize=11)
ax1.set_title(f'Gravitational Wave Chirp: {m1_solar:.0f}+{m2_solar:.0f} M☉ at {d_Mpc:.0f} Mpc',
              color='white', fontsize=13, fontweight='bold')
ax1.tick_params(colors='#cbd5e1')
for spine in ax1.spines.values():
    spine.set_color('#334155')
ax1.grid(True, alpha=0.15, color='#64748b')
ax1.ticklabel_format(axis='y', style='scientific', scilimits=(-2,2))

# Panel 2: GW frequency evolution
ax2 = axes[1]
ax2.set_facecolor('#1e293b')
ax2.plot(t_plot, f_plot, color='#f472b6', linewidth=1.2)
ax2.axhline(y=f_isco, color='#ef4444', linestyle='--', alpha=0.7, label=f'f_ISCO = {f_isco:.0f} Hz')
ax2.set_xlabel('Time (s)', color='#cbd5e1', fontsize=11)
ax2.set_ylabel('GW Frequency (Hz)', color='#cbd5e1', fontsize=11)
ax2.set_title('Frequency Chirp', color='white', fontsize=13, fontweight='bold')
ax2.tick_params(colors='#cbd5e1')
ax2.legend(fontsize=10, facecolor='#1e293b', edgecolor='#334155', labelcolor='white')
for spine in ax2.spines.values():
    spine.set_color('#334155')
ax2.grid(True, alpha=0.15, color='#64748b')
ax2.set_yscale('log')

# Panel 3: Strain amplitude envelope
ax3 = axes[2]
ax3.set_facecolor('#1e293b')
ax3.semilogy(t_plot, h_amp_plot, color='#4ade80', linewidth=1.5)
ax3.set_xlabel('Time (s)', color='#cbd5e1', fontsize=11)
ax3.set_ylabel('Strain Amplitude |h|', color='#cbd5e1', fontsize=11)
ax3.set_title('Amplitude Growth', color='white', fontsize=13, fontweight='bold')
ax3.tick_params(colors='#cbd5e1')
for spine in ax3.spines.values():
    spine.set_color('#334155')
ax3.grid(True, alpha=0.15, color='#64748b')

plt.tight_layout()
plt.savefig('output.png', dpi=150, bbox_inches='tight', facecolor='#0f172a')
plt.close()
print("\nPlot saved successfully.")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

3. Compact Binary Mergers

The complete gravitational waveform from a compact binary coalescence consists of three phases: inspiral, merger, and ringdown. While the inspiral is well-described by post-Newtonian theory and the ringdown by black hole perturbation theory, the merger phase requires full numerical relativity simulations.

3.1 Types of Compact Binary Mergers

Binary Black Hole (BBH)

The cleanest GW source: no matter effects, purely gravitational dynamics. The waveform is determined entirely by 8 intrinsic parameters (two mass values, two spin vectors). BBH mergers produce the strongest signals and were the first detected (GW150914). Mass range spans from stellar-mass ($\sim 5$–$100\;M_\odot$) to supermassive ($\sim 10^6$–$10^9\;M_\odot$).

Binary Neutron Star (BNS)

Involves tidal interactions during late inspiral that encode the neutron star equation of state through the tidal deformability parameter $\Lambda$. The merger produces a hypermassive or supramassive neutron star remnant (or prompt collapse to BH) and is accompanied by electromagnetic counterparts: short gamma-ray bursts, kilonovae, and multi-wavelength afterglows.

Neutron Star–Black Hole (NSBH)

The neutron star may be tidally disrupted outside the BH's innermost stable circular orbit (ISCO) if the BH mass is sufficiently low and/or its spin is aligned with the orbital angular momentum. Disruption produces an accretion disk and potential EM counterpart; otherwise the NS plunges directly into the BH with minimal EM signature.

3.2 Inspiral-Merger-Ringdown Waveform

The three phases of a compact binary coalescence waveform are:

Inspiral: The binary slowly spirals inward due to GW emission. The waveform is accurately described by post-Newtonian (PN) theory or the effective-one-body (EOB) formalism. The frequency sweeps upward as $f \propto \tau^{-3/8}$ and the amplitude grows as$h \propto f^{2/3}$. This phase dominates the signal duration.
Merger: The binary reaches the innermost stable orbit and plunges together. Post-Newtonian theory breaks down and full numerical relativity is required. The waveform reaches peak amplitude and frequency. For BBH, the merger produces a highly distorted, rapidly spinning BH.
Ringdown: The remnant BH settles to a stationary Kerr state by radiating quasi-normal modes (QNMs). The waveform is a superposition of exponentially damped sinusoids.

Modern waveform models (IMRPhenom, SEOBNRv4, NRSur7dq4) combine all three phases into complete inspiral-merger-ringdown (IMR) templates used for matched filtering searches and parameter estimation.

3.3 Numerical Relativity

The binary merger problem requires solving the full nonlinear Einstein equations on a computer. After decades of effort, the breakthrough came in 2005 when three groups (Pretorius; Campanelli et al.; Baker et al.) independently achieved stable long-term BBH merger simulations using different techniques:

Generalized Harmonic formulation (Pretorius): coordinate conditions as wave equations with constraint damping
Moving puncture method (Campanelli et al.; Baker et al.): BSSN formulation with 1+log slicing and Gamma-driver shift conditions

Modern NR catalogs (SXS, RIT, Georgia Tech) contain thousands of simulations spanning the BBH parameter space. These simulations are essential for calibrating semi-analytical waveform models and for building NR surrogate models that interpolate across parameters at a fraction of the computational cost.

3.4 Ringdown: Quasi-Normal Modes

After merger, the remnant black hole rings down by emitting quasi-normal modes (QNMs)—the characteristic vibrational frequencies of a Kerr black hole. Each QNM is labeled by indices $(l, m, n)$ and has a complex frequency:

$$\omega_{lmn} = 2\pi f_{lmn} + i/\tau_{lmn}$$

The ringdown waveform is a superposition of damped sinusoids:

$$h(t) = \sum_{lmn} A_{lmn}\, e^{-t/\tau_{lmn}}\cos(2\pi f_{lmn}\, t + \phi_{lmn})$$

For the dominant $(l,m,n) = (2,2,0)$ mode, the frequency and damping time depend only on the remnant mass $M_f$ and dimensionless spin $\chi_f = a_f c/(GM_f)$:

$$f_{220} \approx \frac{c^3}{2\pi G M_f}\left[1.5251 - 1.1568(1-\chi_f)^{0.1292}\right]$$

$$Q_{220} = \pi f_{220}\tau_{220} \approx 0.7000 + 1.4187(1-\chi_f)^{-0.4990}$$

For GW150914 ($M_f \approx 62\;M_\odot$, $\chi_f \approx 0.67$), $f_{220} \approx 251$ Hz and $\tau_{220} \approx 4$ ms. The detection of multiple QNM frequencies enables "black hole spectroscopy"—testing the no-hair theorem by verifying that all mode frequencies are consistent with a single Kerr BH.

3.5 Final Mass and Spin

A significant fraction of the total mass-energy is radiated as gravitational waves during the merger. For non-spinning equal-mass binaries, numerical relativity gives:

$$M_f \approx (1 - E_{\rm rad}/Mc^2)\,M, \quad E_{\rm rad}/Mc^2 \approx 0.05\;\text{(non-spinning equal mass)}$$

The radiated energy fraction can reach up to $\sim 10$% for maximally spinning, aligned BH binaries. For GW150914: $M_1 + M_2 \approx 65\;M_\odot$, $M_f \approx 62\;M_\odot$, so $E_{\rm rad} \approx 3\;M_\odot c^2 \approx 5.4 \times 10^{47}$ J was radiated in about 0.2 seconds.

The final spin for non-spinning equal-mass mergers is $\chi_f \approx 0.69$. Fitting formulas calibrated to NR simulations (e.g., Hofmann, Barausse & Rezzolla 2016) give $M_f$ and $\chi_f$ as functions of the mass ratio and component spins. Gravitational wave recoil ("kicks") from asymmetric emission can reach up to $\sim 5000$ km/s for certain spin configurations, potentially ejecting the remnant BH from its host galaxy.

4. LIGO/Virgo Detection

The Laser Interferometer Gravitational-Wave Observatory (LIGO) and its European partner Virgo are kilometer-scale Michelson interferometers designed to measure spacetime strain at the level of one part in a thousand billion billion.

4.1 Michelson Interferometer Principle

A Michelson interferometer splits a laser beam into two perpendicular arms, reflects the beams off mirrors (test masses) at the ends of the arms, and recombines them at the beamsplitter. A passing gravitational wave differentially modulates the arm lengths:

$$\Delta L = \frac{1}{2}h\,L$$

where $L$ is the arm length and $h$ is the strain amplitude. For LIGO ($L = 4$ km) and $h \sim 10^{-21}$:

$$\Delta L \sim \frac{1}{2}\times 10^{-21}\times 4000\;\text{m} \approx 2\times 10^{-18}\;\text{m}$$

This is about $1/1000$ the diameter of a proton! LIGO achieves this sensitivity through several key technologies: Fabry-Perot arm cavities (effective arm length $\sim 1200$ km from$\sim 300$ bounces), power recycling (circulating power $\sim 750$ kW), signal recycling, and the world's most stable laser and mirror suspension systems.

4.2 Noise Sources

The sensitivity of a GW detector is characterized by its noise power spectral density (PSD), which varies with frequency. The main noise sources are:

Seismic Noise (f < 10 Hz)

Ground vibrations from microseisms, traffic, wind, earthquakes. Mitigated by multi-stage pendulum suspensions providing $f^{-2}$ isolation per stage. The seismic wall at$\sim 10$ Hz defines the low-frequency limit of ground-based detectors. Active seismic isolation systems push this below 10 Hz.

Thermal Noise (10–100 Hz)

Brownian motion of mirror surfaces and suspension fibers. Minimized using high-Q materials (fused silica fibers, crystalline mirror coatings) and cryogenic cooling (planned for KAGRA, Einstein Telescope). Thermal noise from mirror coatings is the dominant noise source in the most sensitive frequency band.

Quantum Noise (f > 100 Hz)

Shot noise: Poisson fluctuations in photon arrival times. Scales as $\propto 1/\sqrt{P}$ where$P$ is the circulating laser power. Dominates at high frequencies.
Radiation pressure noise: Random momentum kicks from photon fluctuations on the mirrors. Scales as $\propto \sqrt{P}$. Dominates at low frequencies. The trade-off between these two defines the Standard Quantum Limit (SQL). Frequency-dependent squeezing can beat the SQL and is deployed in Advanced LIGO since O4.

Other Noise Sources

Gravity gradient (Newtonian) noise from local mass fluctuations, residual gas scattering, laser frequency and intensity noise, scattered light, control system noise, and transient instrumental artifacts ("glitches"). Each requires dedicated mitigation strategies.

4.3 Matched Filtering and Signal-to-Noise Ratio

When the signal waveform is known (as for compact binary inspirals), matched filtering is the optimal detection strategy. The detector output $s(t) = n(t) + h(t)$ is cross-correlated with template waveforms. The optimal signal-to-noise ratio (SNR) is:

$$\rho^2 = 4\int_0^\infty \frac{|\tilde{h}(f)|^2}{S_n(f)}\,df$$

where $\tilde{h}(f)$ is the Fourier transform of the GW signal and $S_n(f)$ is the one-sided noise power spectral density. A detection typically requires $\rho \geq 8$ in each detector, with a network SNR $\rho_{\rm net} = \sqrt{\sum_i \rho_i^2} \geq 12$.

The template bank for CBC searches contains $\sim 10^5$–$10^6$ templates covering the mass and spin parameter space. The false alarm rate is assessed by time-sliding detector outputs relative to each other. For GW150914, the false alarm probability was less than $2 \times 10^{-7}$, corresponding to a $5.1\sigma$ detection significance.

For a compact binary inspiral, the SNR can be expressed in terms of the chirp mass and distance:

$$\rho \propto \frac{\mathcal{M}^{5/6}}{d_L}\left[\int_{f_{\rm low}}^{f_{\rm ISCO}} \frac{f^{-7/3}}{S_n(f)}\,df\right]^{1/2}$$

4.4 GW150914: The First Detection

On September 14, 2015, LIGO made the first direct detection of gravitational waves from the merger of two black holes. Key parameters of this landmark event:

Source masses	$m_1 = 36^{+5}_{-4}\;M_\odot$, $m_2 = 29^{+4}_{-4}\;M_\odot$
Chirp mass	$\mathcal{M} = 28.3^{+1.8}_{-1.5}\;M_\odot$
Final mass	$M_f = 62^{+4}_{-4}\;M_\odot$
Final spin	$\chi_f = 0.67^{+0.05}_{-0.07}$
Radiated energy	$E_{\rm rad} \approx 3.0\;M_\odot c^2 \approx 5.4\times10^{47}\;\text{J}$
Peak luminosity	$L_{\rm peak} \approx 3.6\times10^{49}\;\text{W} \approx 200\;M_\odot c^2/\text{s}$
Luminosity distance	$d_L = 410^{+160}_{-180}\;\text{Mpc}\;(z \approx 0.09)$
Peak strain	$h_{\rm peak} \approx 1.0 \times 10^{-21}$
Network SNR	24

The signal lasted about 0.2 seconds in the LIGO band, sweeping from 35 Hz to 250 Hz. The 7 ms time delay between Hanford and Livingston was consistent with a source in the southern sky. This detection confirmed the existence of stellar-mass binary black holes and validated general relativity in the strong-field, highly dynamical regime.

4.5 Observing Run Catalog Highlights

Through four observing runs, the LIGO-Virgo-KAGRA collaboration has assembled a rich catalog of gravitational wave events:

O1 (2015–2016): 3 BBH detections including GW150914, GW151226, and candidate LVT151012
O2 (2016–2017): 8 events including the first BNS merger GW170817 and the first triple-detector BBH GW170814
O3 (2019–2020): ~80 candidates in GWTC-3, including GW190521 (the most massive BBH with remnant in the intermediate-mass range), GW190814 (highly asymmetric with a 2.6 M&odot; companion), and GW200105/GW200115 (first confirmed NSBH mergers)
O4 (2023–2025): Greatly expanded catalog with improved sensitivity using frequency-dependent squeezing, detection rate exceeding 1 event per week

The observed BBH mass distribution reveals features including a peak near $\sim 35\;M_\odot$, a possible gap between $\sim 50$–$120\;M_\odot$ from pair-instability supernovae, and a population of low-spin systems suggestive of dynamical formation channels.

Simulation: LIGO Sensitivity and Detection Horizon

This simulation plots the Advanced LIGO noise amplitude spectral density alongside the characteristic strain from compact binary inspirals at various distances. It computes the optimal matched-filter SNR and detection horizon distance.

LIGO Sensitivity Curve & Detection Horizon

Python

Plot LIGO noise curve with compact binary signals and SNR vs distance

Parameters

Chirp Mass (M☉)

Chirp mass of binary

Distance (Mpc)

Luminosity distance

SNR Threshold

Detection SNR threshold

ligo_sensitivity.py141 lines

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# Physical constants
G = 6.67430e-11
c = 2.99792458e8
M_sun = 1.98892e30
pc = 3.08567758e16

# Parameters
Mc_solar = MC_VAL     # chirp mass in solar masses
d_Mpc = DIST2_VAL     # distance in Mpc
snr_threshold = SNR_VAL  # detection threshold

Mc = Mc_solar * M_sun
d = d_Mpc * 1e6 * pc

print(f"=== LIGO Sensitivity & Detection Horizon ===")
print(f"Chirp mass = {Mc_solar:.1f} M_sun")
print(f"Source distance = {d_Mpc:.1f} Mpc")
print(f"SNR threshold = {snr_threshold:.1f}")

# Frequency array (LIGO band: 10 Hz to 5000 Hz)
f = np.logspace(np.log10(10), np.log10(5000), 5000)

# aLIGO design sensitivity (analytic approximation)
# Based on LIGO-T1800044 noise budget
f0 = 215.0  # Hz
S0 = 1.0e-49  # Hz^-1

# Phenomenological fit to aLIGO design
x = f / f0
Sn = S0 * (x**(-4.14) - 5.0*x**(-2) + 111.0*(1.0 - x**2 + 0.5*x**4) / (1.0 + 0.5*x**2))
# Ensure physical values
Sn = np.maximum(Sn, 1e-50)

# Add low-frequency seismic wall
seismic_wall = S0 * 1e8 * (f / 10.0)**(-20)
Sn = Sn + seismic_wall

# Characteristic strain noise: h_n(f) = sqrt(f * S_n(f))
h_n = np.sqrt(f * Sn)

# GW signal: characteristic strain from compact binary inspiral
# h_c(f) = (G*Mc)^(5/6) / (pi^(2/3) * c^(3/2) * d) * f^(-1/6)
# valid up to ISCO
f_isco_arr = []
h_c_signals = []
distances = [50, 200, 500, 1000, 3000]
colors_sig = ['#4ade80', '#38bdf8', '#f472b6', '#fbbf24', '#a78bfa']

for dist_mpc in distances:
    d_i = dist_mpc * 1e6 * pc
    h_c_i = (1.0 / (np.pi**(2.0/3.0))) * (G * Mc)**(5.0/6.0) / (c**(3.0/2.0) * d_i) * f**(-1.0/6.0)
    # Convert to characteristic strain: h_c = 2*f*|h_tilde(f)|
    # For inspiral: h_c(f) = sqrt(4*f^2 * |h(f)|^2) ~ A * f^(-1/6) effectively
    h_c_signals.append(h_c_i)

# Also compute for user-specified distance
h_c_user = (1.0 / (np.pi**(2.0/3.0))) * (G * Mc)**(5.0/6.0) / (c**(3.0/2.0) * d) * f**(-1.0/6.0)

# Compute optimal matched-filter SNR for user distance
# SNR^2 = 4 * integral |h(f)|^2 / S_n(f) df
# For inspiral: |h(f)|^2 = (1/(16*pi^(4/3))) * (G*Mc)^(10/3) / (c^5 * d^2) * f^(-7/3)
integrand = (1.0 / (16.0 * np.pi**(4.0/3.0))) * (G * Mc)**(10.0/3.0) / (c**5 * d**2) * f**(-7.0/3.0) / Sn
# Trapezoidal integration
snr_sq = 4.0 * np.trapezoid(integrand, f)
snr = np.sqrt(np.maximum(snr_sq, 0))
print(f"Optimal SNR at {d_Mpc:.0f} Mpc = {snr:.1f}")

# Detection horizon: distance at which SNR = threshold
# SNR scales as 1/d, so d_horizon = d * (snr / snr_threshold)
if snr > 0:
    d_horizon_Mpc = d_Mpc * (snr / snr_threshold)
    print(f"Detection horizon (SNR={snr_threshold:.0f}) = {d_horizon_Mpc:.0f} Mpc = {d_horizon_Mpc/1000:.2f} Gpc")
    z_approx = d_horizon_Mpc * 70.0 / 3e5  # rough Hubble law
    print(f"Approximate redshift z ~ {z_approx:.2f}")

# Create figure
fig, axes = plt.subplots(2, 1, figsize=(12, 10), facecolor='#0f172a')

# Panel 1: Noise ASD and signal characteristic strain
ax1 = axes[0]
ax1.set_facecolor('#1e293b')
ax1.loglog(f, np.sqrt(Sn), color='#ef4444', linewidth=2.5, label='aLIGO Design ASD', zorder=5)

for i, dist_mpc in enumerate(distances):
    alpha = 0.9 if dist_mpc == d_Mpc else 0.5
    lw = 2.0 if dist_mpc == d_Mpc else 1.2
    ax1.loglog(f, h_c_signals[i], color=colors_sig[i], linewidth=lw, alpha=alpha,
               label=f'{Mc_solar:.0f} M☉ @ {dist_mpc} Mpc', linestyle='--')

ax1.loglog(f, h_c_user, color='#ffffff', linewidth=2.0, alpha=0.9,
           label=f'{Mc_solar:.0f} M☉ @ {d_Mpc:.0f} Mpc (selected)', zorder=6)

ax1.set_xlabel('Frequency (Hz)', color='#cbd5e1', fontsize=12)
ax1.set_ylabel('Strain (1/√Hz) or h_c', color='#cbd5e1', fontsize=12)
ax1.set_title('LIGO Sensitivity & Compact Binary Signals', color='white', fontsize=14, fontweight='bold')
ax1.set_xlim(10, 5000)
ax1.set_ylim(1e-25, 1e-18)
ax1.legend(fontsize=8, facecolor='#1e293b', edgecolor='#334155', labelcolor='white', loc='upper right')
ax1.tick_params(colors='#cbd5e1')
for spine in ax1.spines.values():
    spine.set_color('#334155')
ax1.grid(True, alpha=0.15, color='#64748b', which='both')

# Panel 2: SNR vs distance
ax2 = axes[1]
ax2.set_facecolor('#1e293b')
d_range = np.logspace(1, 4, 500)  # 10 to 10000 Mpc
snr_range = snr * (d_Mpc / d_range)  # SNR scales as 1/d

ax2.loglog(d_range, snr_range, color='#38bdf8', linewidth=2.5, label=f'Mc = {Mc_solar:.0f} M☉')
ax2.axhline(y=snr_threshold, color='#ef4444', linestyle='--', linewidth=1.5, alpha=0.8,
            label=f'Detection threshold (SNR={snr_threshold:.0f})')
ax2.axvline(x=d_Mpc, color='#fbbf24', linestyle=':', linewidth=1.5, alpha=0.7,
            label=f'Source distance = {d_Mpc:.0f} Mpc')
if snr > 0:
    ax2.axvline(x=d_horizon_Mpc, color='#4ade80', linestyle=':', linewidth=1.5, alpha=0.7,
                label=f'Horizon = {d_horizon_Mpc:.0f} Mpc')
    ax2.plot(d_Mpc, snr, 'o', color='#fbbf24', markersize=10, zorder=5)
    ax2.plot(d_horizon_Mpc, snr_threshold, 's', color='#4ade80', markersize=10, zorder=5)

ax2.set_xlabel('Distance (Mpc)', color='#cbd5e1', fontsize=12)
ax2.set_ylabel('Optimal SNR', color='#cbd5e1', fontsize=12)
ax2.set_title('Detection Signal-to-Noise Ratio vs Distance', color='white', fontsize=14, fontweight='bold')
ax2.set_xlim(10, 10000)
ax2.set_ylim(0.5, 1e4)
ax2.legend(fontsize=10, facecolor='#1e293b', edgecolor='#334155', labelcolor='white')
ax2.tick_params(colors='#cbd5e1')
for spine in ax2.spines.values():
    spine.set_color('#334155')
ax2.grid(True, alpha=0.15, color='#64748b', which='both')

plt.tight_layout()
plt.savefig('output.png', dpi=150, bbox_inches='tight', facecolor='#0f172a')
plt.close()
print("\nPlot saved successfully.")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

5. Multi-Messenger Astronomy

The simultaneous observation of gravitational waves and electromagnetic radiation from the same astrophysical event represents the dawn of multi-messenger astronomy—a transformative capability that combines the strengths of different observational channels to extract physics inaccessible to either alone.

5.1 GW170817: The First Multi-Messenger Event

On August 17, 2017, LIGO and Virgo detected gravitational waves from a binary neutron star merger (GW170817) with component masses $m_1 \in [1.36, 1.60]\;M_\odot$ and $m_2 \in [1.17, 1.36]\;M_\odot$ at a distance of $\sim 40$ Mpc. The signal lasted approximately 100 seconds in the LIGO band, accumulating over 3000 cycles. This was followed by a cascade of electromagnetic counterparts:

t + 1.7 s: Short gamma-ray burst GRB 170817A detected by Fermi-GBM and INTEGRAL, confirming the long-hypothesized connection between BNS mergers and sGRBs
t + 11 hours: Optical/IR counterpart AT2017gfo (kilonova/macronova) discovered in NGC 4993 by the Swope telescope, showing thermal emission from r-process nucleosynthesis ejecta
t + 9 days: X-ray afterglow detected by Chandra, consistent with an off-axis structured jet
t + 16 days: Radio afterglow detected by VLA

Over 70 observatories across the electromagnetic spectrum and in neutrinos observed the aftermath, producing a wealth of scientific results from a single event.

5.2 Electromagnetic Counterparts

Short Gamma-Ray Burst

The sGRB from GW170817 was underluminous by a factor of $\sim 1000$ compared to cosmological sGRBs, consistent with an off-axis viewing angle of $\theta_{\rm obs} \approx 20$–$30°$. This confirmed that sGRBs are produced by relativistic jets from BNS mergers and that the jet has angular structure (a narrow core with a wider cocoon).

Kilonova

The kilonova AT2017gfo showed a blue component (lanthanide-poor, $T \sim 10^4$ K, lasting days) and a red component (lanthanide-rich, $T \sim 2500$ K, lasting weeks). The blue emission arises from shock-heated polar ejecta; the red from tidally-stripped equatorial material rich in heavy r-process elements. The total ejecta mass was estimated at $\sim 0.05\;M_\odot$, producing$\sim 10$ Earth masses of gold and platinum-group elements.

Afterglow

The delayed X-ray and radio afterglow, rising over months, was consistent with a structured relativistic jet viewed off-axis. Superluminal motion of the radio source confirmed the jet nature. The afterglow modeling provided constraints on the jet energy, opening angle, and the surrounding interstellar medium density.

5.3 Standard Sirens and the Hubble Constant

Gravitational waves provide a direct measurement of the luminosity distance to the source, independent of the cosmic distance ladder. The GW strain amplitude is inversely proportional to the luminosity distance:

$$h \propto \frac{\mathcal{M}^{5/3}f^{2/3}}{d_L}$$

Since the chirp mass is independently determined from the frequency evolution, the distance is directly measured from the signal amplitude. Combined with a redshift measurement from the electromagnetic counterpart (host galaxy NGC 4993, $z = 0.0099$), this gives the Hubble constant:

$$H_0 = \frac{cz}{d_L} = 70.0^{+12.0}_{-8.0}\;\text{km/s/Mpc}$$

While this single measurement has large uncertainties (primarily from the GW distance-inclination degeneracy), it is consistent with both the Planck CMB value ($H_0 = 67.4$) and the local distance ladder value ($H_0 = 73.0$). With $\sim 50$ BNS events with counterparts, the Hubble tension can be resolved at the 2% level. "Dark sirens" (BBH events without counterparts) can also contribute using statistical host galaxy associations.

5.4 Neutron Star Equation of State from Tidal Deformability

During the late inspiral, tidal interactions between the neutron stars modify the GW phase. The key parameter is the dimensionless tidal deformability:

$$\Lambda = \frac{2}{3}k_2\left(\frac{c^2 R}{Gm}\right)^5$$

where $k_2$ is the Love number and $R$ is the neutron star radius. The tidal deformability measures how easily the star is deformed by the companion's tidal field: larger $\Lambda$ corresponds to a stiffer EOS and larger radius. The tidal contribution to the GW phase enters at 5PN order:

$$\delta\Phi_{\rm tidal} = -\frac{39}{2}\tilde{\Lambda}\left(\frac{G\mathcal{M}\pi f}{c^3}\right)^{10/3}$$

where $\tilde{\Lambda}$ is the mass-weighted tidal deformability. From GW170817, the constraint $\tilde{\Lambda} = 300^{+420}_{-230}$ (at 90% credibility) was obtained, implying a neutron star radius of $R_{1.4} \approx 11.0$–$13.0$ km for a 1.4 M&odot; star. This disfavors very stiff equations of state and rules out some exotic matter models.

Combined with the NICER X-ray timing measurements of pulsars PSR J0030+0451 and PSR J0740+6620, the NS mass-radius relation is being precisely constrained, probing the behavior of matter at supranuclear densities.

6. Future Gravitational Wave Detectors

The gravitational wave spectrum spans many decades in frequency, from nanohertz to kilohertz, each band requiring different detection technologies and probing different astrophysical sources.

6.1 LISA: Space-Based Millihertz Observatory

The Laser Interferometer Space Antenna (LISA), adopted by ESA with a planned launch in the mid-2030s, will consist of three spacecraft in a triangular formation with $2.5 \times 10^6$ km arm lengths, trailing Earth in a heliocentric orbit. LISA will observe in the millihertz band ($10^{-4}$–$10^{-1}$ Hz), inaccessible from the ground due to seismic and gravity gradient noise.

Key LISA science targets include:

Massive black hole binaries ($10^4$–$10^7\;M_\odot$): mergers at cosmological distances with SNR up to thousands, enabling precision tests of GR and mapping the cosmic history of SMBH assembly
Extreme mass-ratio inspirals (EMRIs): stellar-mass objects spiraling into SMBHs, producing $\sim 10^5$ orbit cycles that map the Kerr spacetime geometry to exquisite precision
Galactic binaries: thousands of known ultra-compact WD-WD binaries serving as verification sources, plus discovery of the Galactic double-compact-object population
Stochastic backgrounds: from unresolved sources, early-universe phase transitions, and cosmic strings

The LISA Pathfinder mission (2015–2017) demonstrated the required drag-free control and acceleration noise performance, achieving results exceeding LISA requirements by a factor of 5.

6.2 Third-Generation Ground-Based Detectors

Next-generation ground-based detectors will achieve a factor of $\sim 10$ improvement in strain sensitivity over Advanced LIGO, corresponding to a factor of $\sim 1000$ increase in detection volume:

Einstein Telescope (ET)

European project: triangular underground detector with 10 km arms. Uses a xylophone configuration with separate high-frequency (room temperature) and low-frequency (cryogenic, 10 K silicon mirrors) interferometers. Extends sensitivity down to $\sim 3$ Hz. Expected to detect BBH mergers throughout the observable universe.

Cosmic Explorer (CE)

US project: L-shaped surface detector with 40 km arms. Longer arms provide proportionally better strain sensitivity. Will detect stellar-mass BBH mergers to $z \sim 100$ and BNS mergers to $z \sim 5$. Two detectors planned for a network with ET.

With 3G detectors, the detection rate will increase to $\sim 10^5$–$10^6$ events per year, enabling precision population studies, cosmography with standard sirens to percent-level accuracy, and detection of intermediate-mass black holes and primordial gravitational wave backgrounds.

6.3 Pulsar Timing Arrays: Nanohertz Band

Pulsar Timing Arrays (PTAs) use networks of millisecond pulsars as a galaxy-scale gravitational wave detector. A gravitational wave passing between Earth and a pulsar modulates the pulse arrival times with a characteristic pattern:

$$\frac{\delta t}{t} \sim h \quad \text{at frequencies } f \sim \frac{1}{T_{\rm obs}} \sim 1\text{--}100\;\text{nHz}$$

The hallmark of a GW background is the Hellings-Downs angular correlation: correlations between timing residuals of pulsar pairs depend on their angular separation $\zeta$ as:

$$\Gamma(\zeta) = \frac{1}{2} - \frac{1}{4}\left(\frac{1 - \cos\zeta}{2}\right) + \frac{3}{2}\left(\frac{1 - \cos\zeta}{2}\right)\ln\left(\frac{1 - \cos\zeta}{2}\right)$$

The primary PTA source is the stochastic background from the cosmic population of inspiraling supermassive black hole binaries formed after galaxy mergers, with a characteristic strain spectrum:

$$h_c(f) = A\left(\frac{f}{f_{\rm yr}}\right)^{-2/3}$$

where $f_{\rm yr} = 1/(1\;\text{year})$ is the reference frequency.

6.4 NANOGrav and IPTA: Evidence for the GW Background

In 2023, the NANOGrav collaboration (North American Nanohertz Observatory for Gravitational Waves), along with the European PTA, the Parkes PTA, and the Chinese PTA, independently announced evidence for a stochastic gravitational wave background in the nanohertz frequency band.

The NANOGrav 15-year dataset, analyzing timing residuals from 67 millisecond pulsars, found a common-spectrum red noise process with amplitude $A \approx 2.4 \times 10^{-15}$ at$f_{\rm yr}$ and spectral index consistent with $-2/3$ (the expected SMBH binary value). Crucially, the inter-pulsar correlations are consistent with the Hellings-Downs pattern expected for a gravitational wave origin, with evidence at the $\sim 3$–$4\sigma$ level.

While the most natural interpretation is a background from SMBH binaries, alternative sources include cosmological phase transitions, cosmic strings, domain walls, primordial gravitational waves from inflation, and scalar-induced GWs from enhanced primordial perturbations. Continued observations with more pulsars and longer baselines will improve the significance and enable spectral characterization to distinguish among these sources.

6.5 The Gravitational Wave Spectrum: A Complete Picture

Band	Frequency	Detector	Sources
Ultra-low	$10^{-18}$–$10^{-16}$ Hz	CMB B-modes	Inflationary GWs
Nanohertz	$10^{-9}$–$10^{-7}$ Hz	PTAs	SMBH binaries, cosmological
Millihertz	$10^{-4}$–$10^{-1}$ Hz	LISA	SMBH mergers, EMRIs, galactic binaries
Decihertz	$10^{-1}$–$10$ Hz	DECIGO, BBO	Intermediate-mass BHs, early inspiral
Audio	$10$–$10^4$ Hz	LIGO/Virgo, ET, CE	Stellar-mass CBC, supernovae, pulsars

Together, these detectors will cover more than 20 decades in frequency, providing a complete view of the gravitational wave universe from the Big Bang to the present day.

← Part IV: High-Energy Astrophysics Part VI: Extragalactic Astrophysics →