Calculations of Gravitational Waves from Binary Black Hole Systems

Simulating a binary black hole merger properly means solving the Einstein field equations numerically, which is the work of large collaborations. This project asks a narrower and more tractable question: how much of the observed signal can be recovered from classical mechanics alone?

The strategy is to give up on computing the waveform propagating through space and compute only the orbit. Two masses circle each other; they lose energy at the rate general relativity predicts for gravitational radiation; the orbit therefore shrinks and the orbital frequency rises, until the two event horizons touch. If gravitational waves are emitted more strongly the faster the bodies move, then the frequency sweep of the orbit should resemble the frequency sweep of the signal — and that is something a laptop can integrate. The resulting “chirp” is then set beside the real spectrograms of GW150914 and GW190521.

§Research question

How far can the frequency evolution of a gravitational-wave signal from a binary black hole merger be reproduced using classical orbital mechanics together with the relativistic energy-loss rate, without solving the field equations?

§Theoretical background

For two masses $m_1$ and $m_2$ in a circular orbit of separation $r$, Peters and Mathews (1963) give the rate at which the system radiates energy as gravitational waves:

$$\frac{\dd E}{\dd t} = -\frac{32}{5}\,\frac{G^4}{c^5}\,\frac{(m_1 m_2)^2 (m_1 + m_2)}{r^5}.$$

Energy loss shrinks the orbit. Relating $\dd E/\dd t$ to $\dd r/\dd t$ through the orbital energy gives a differential equation for the separation alone, of the form $\dd r/\dd t \propto -r^{-3}$, so the inspiral accelerates sharply as the bodies close.

The orbit cannot shrink indefinitely, because black holes have a size. For mass $M$ the Schwarzschild radius — the event horizon, where escape velocity reaches the speed of light — is

$$R_s = \frac{2GM}{c^2},$$

and taking merger to occur when the two horizons meet fixes the end of the integration at

$$r_{\min} = \frac{2G(m_1 + m_2)}{c^2}.$$

The rotation rate follows from Kepler’s third law: for a binary of total mass $M$ at separation $r$,

$$\omega = \sqrt{\frac{G M}{r^3}},$$

so as $r$ falls, $\omega$ rises — the chirp.

§Methodology

Assumptions. The model is explicitly a simplification, and the manuscript states its assumptions:

1. The binary’s orbital motion follows the same equations as a classical binary star system.
2. The black holes merge into one immediately on contact.
3. Orbits are circular — eccentricity zero.
4. Gravitational waves grow stronger as the black holes move faster.
5. Each black hole emits a pulse once as it approaches Earth and once as it recedes.

Apart from the energy-loss rate and the Schwarzschild radius, nothing relativistic enters.

Procedure. Given two masses and an initial separation, $\dd r/\dd t$ was integrated numerically until $r$ reached $r_{\min}$. From $r(t)$ the angular velocity followed by Kepler’s law, and from that the trajectory of each body. With Earth taken to lie far away along the $x$-axis, the $x$-component of each body’s velocity was extracted as the quantity standing in for wave amplitude. Because assumption 5 puts two pulses per body per revolution — four per revolution in total — the frequency axis was multiplied by four. A spectrogram of the summed $v_x$ was then computed with SciPy and compared against the published LIGO spectrograms.

Implementation. Python 3.10 with NumPy and SciPy for computation, and Matplotlib for visualisation.

§Results and observations

GW150914. Initial conditions matched to the first detection — 36 and 29 solar masses, initial separation 1000 km. The integration shows the separation contracting slowly and then collapsing, with a sharp acceleration late in the run. Against the observed event the manuscript reports three points of agreement: the strongest signal lasts about 0.2 s; the frequency rises from about 35 Hz to about 250 Hz; and the bodies reach 30–60% of the speed of light.

GW190521. Masses of 85 and 66 solar masses, initial separation adjusted to 1700 km. Both simulation and observation show a strong, rapidly rising signal within roughly 0.1 s. The simulated frequency rose through about 60–150 Hz against roughly 50–100 Hz observed. For the heavier pair, simulation and observation agree in direction: a shorter signal at lower frequency. Specific frequency values and spectrogram shape differ.

The honest summary is the manuscript’s own: that classical mechanics plus one relativistic energy-loss formula produces “quite plausible” frequency evolution, with the agreement being qualitative — duration, direction and rough band — rather than a reproduction of the waveform.

§Challenges and limitations

• This is not a waveform. The model computes orbital frequency, not strain. It cannot produce the merger or ringdown that dominates the end of a real signal, and the comparison is between a frequency sweep and a spectrogram, not between two waveforms.
• The last orbits are the least trustworthy. At 30–60% of the speed of light and separations of a few hundred kilometres, post-Newtonian corrections are large and a Newtonian orbit is well outside its domain of validity. The regime the model is compared against is exactly the regime where it is weakest.
• Horizons touching is not merger. Real binaries become dynamically unstable at the innermost stable circular orbit, well before the horizons meet, so the cut-off is a convenient stopping rule rather than a physical one.
• Circular orbits and no spin. Eccentricity and spin are both set aside; both matter for real systems.
• The factor of four is an assumption, not a derivation. Quadrupole radiation from a circular binary emerges at twice the orbital frequency. The four-pulses-per-revolution rule is a modelling choice, and the agreement with the observed band depends on it.