Neutron Stars in Numerical Relativity

Numerical relativity is the backbone of my research career: it is the art of teaching supercomputers to solve Einstein's equations so that we can watch neutron stars and black holes collide in virtual universes. In my PhD work on "Neutron Stars in Numerical Relativity," I learned how to turn abstract equations into concrete simulations that predict gravitational waves and matter dynamics in some of the most extreme events known in the cosmos.​

What numerical relativity is

In everyday physics, many problems can be solved with pen and paper or simple formulas, but Einstein's theory of general relativity is different: it describes gravity as the curvature of space and time, governed by highly nonlinear partial differential equations. For complex, dynamical systems like merging neutron stars, these equations are far too complicated to solve exactly, so numerical relativity uses algorithms and computers to approximate the solution step by step in time.​

At the heart of general relativity are Einstein's field equations, which in compact form read

$$G_{\mu\nu} = 8\pi T_{\mu\nu}.$$

Here \(G_{\mu\nu}\) describes how space-time is curved, and \(T_{\mu\nu}\) describes the energy and matter content; in simple terms, "matter tells space-time how to curve, and curvature tells matter how to move." In numerical relativity, these equations are rewritten in a form suitable for computers by slicing space-time into "3+1" dimensions: three dimensions of space evolving in one dimension of time.​

Slicing space-time: the 3+1 picture

In my dissertation, a key step is the 3+1 decomposition, where space-time is broken into a stack of spatial slices (like the frames of a movie) and the equations are split into evolution equations and constraint equations. This allows the computer to start from an initial configuration of two neutron stars and then evolve their geometry and matter fields forward in time, enforcing consistency with Einstein's equations at every step.

Technically, the spatial geometry is described by a 3-metric and an extrinsic curvature, which encode the shape of space and how it is embedded in space-time. These quantities are evolved using formulations such as the BSSN scheme, which reorganizes the Einstein equations into a numerically stable system suitable for long-term simulations of inspiral, merger, and post-merger phases.​

From equations to algorithms

To put general relativity "on a computer," my work employed numerical methods like the method of lines, finite differences, spectral methods, and discontinuous Galerkin schemes. The method of lines treats time and space differently: space is discretized on a grid, and the resulting system of ordinary differential equations in time is integrated using time-stepping methods such as Runge–Kutta algorithms.​

Finite-difference methods approximate derivatives by comparing neighboring grid points, while spectral and spectral-element methods represent functions as sums of basis functions, often achieving very high accuracy with relatively few points. In my thesis, I also explored discontinuous Galerkin methods for general relativistic hydrodynamics, which combine the locality and shock-handling of finite-volume schemes with the high accuracy and flexibility of finite-element methods.​

Supercomputing and extreme simulations

All of this mathematics would be useless without large-scale computing power. My simulations ran on high-performance computing clusters, where the computational domain (the region of space containing the neutron stars and their surroundings) is divided among many processors that work in parallel. Parallelization and domain decomposition are essential: resolving the fine structure of neutron stars and their gravitational fields requires millions to billions of grid points, and each timestep involves solving large coupled systems of equations across those points.​

From a lay perspective, a typical simulation is like running a high-resolution 3D movie of a neutron-star merger, except each "frame" is computed from scratch using Einstein's equations rather than drawn by an artist. The code must also respect physical constraints, for example by using high-resolution shock-capturing techniques in the fluid sector to handle shocks and discontinuities in the neutron-star matter without producing unphysical artefacts.​

Connecting to my research themes

During my PhD, numerical relativity was the framework that allowed me to explore highly eccentric and precessing binary neutron star systems in full 3D with consistent initial data. Constructing initial data means solving the constraint part of Einstein's equations—using techniques like the conformal thin-sandwich decomposition—to find physically realistic starting configurations that match the desired masses, spins, orbital shapes, and internal structure of the stars.​

Once the evolution begins, numerical relativity tracks the inspiral, tidal interactions, merger, and post-merger remnant, simultaneously evolving the space-time geometry and the neutron-star matter fields. From these simulations, I could extract gravitational-wave signals, study the ejected matter, and compare the results to approximate models, thus building a bridge between fundamental theory and the actual data observed by gravitational-wave detectors.