Black Hole Simulation

Developed a physics-based numerical simulation modeling light propagation in the curved spacetime surrounding a Schwarzschild black hole, visualizing null geodesic trajectories and demonstrating gravitational lensing effects through relativistic ray tracing.

Code: https://github.com/ItzSmudge/modsim-coursework

Project Overview:

Created a computational framework to simulate photon trajectories around black holes using General Relativity principles, implementing the Schwarzschild metric and geodesic equations to explore how extreme gravitational fields warp spacetime and deflect light paths. The simulation reveals critical phenomena including photon sphere dynamics, event horizon boundaries, and the sensitivity of escape/capture outcomes to initial conditions.

Theoretical Foundation:

General Relativity Framework: Applied Einstein’s field equations through the Schwarzschild metric for non-rotating black holes, describing spacetime geometry via the metric factor f(r) = 1 – r_s/r
Null Geodesic Derivation: Calculated Christoffel symbols (Γ^μ_αβ) from the metric tensor to derive coupled second-order differential equations governing photon paths in curved spacetime
State-Space Formulation: Converted geodesic equations into first-order system with state vector X = [r, φ, ṙ, φ̇] for numerical integration, expressing radial and angular accelerations through metric-dependent terms

Mathematical Implementation:

Geodesic Equation Derivation:

Computed four independent Christoffel symbols (Γ^r_tt, Γ^r_rr, Γ^r_φφ, Γ^φ_rφ) representing gravity, radial curvature, centrifugal, and angular coupling terms respectively
Derived radial acceleration: d²r/dλ² incorporating three competing effects—gravitational attraction, metric correction, and angular momentum
Derived angular acceleration: d²φ/dλ² = -2(ṙφ̇)/r maintaining angular momentum conservation

Coordinate Transformation:

Developed Cartesian-to-polar coordinate mapping for initial conditions, converting position (x, y) and velocity (v_x, v_y) into geodesic parameters (r₀, φ₀, ṙ₀, φ̇₀)
Calculated conserved quantities: angular momentum L = r²φ̇ and energy-like parameter E from null geodesic constraint ds² = 0

Numerical Solution:

Implemented MATLAB’s ode45 (adaptive Runge-Kutta) solver for robust integration of stiff differential equations near the photon sphere (r = 1.5r_s)
Integrated Gaussian noise (σ ~ 10^-6 to 10^-2) into metric components and accelerations to model environmental perturbations from nearby masses and unresolved spacetime fluctuations
Applied geometrized units (G = c = 1) to eliminate extreme numerical magnitudes (10^±22) and prevent overflow/underflow errors while preserving physical behavior

Experimental Analysis:

1. Initial Condition Sensitivity Study:

Demonstrated extreme sensitivity to starting position—0.5 unit change in y-coordinate transformed outcomes from escape to capture or multi-orbit trajectories
Mapped critical boundary between stable escape and black hole capture through systematic parameter sweeps

2. Monte Carlo Escape Probability:

Executed 1000-trial simulations per initial condition to quantify stochastic escape rates
Identified sharp transition zone: 100% escape at y=13.5, 50% at y=13.0, 0% at y=12.6 (constant x=-30), revealing narrow instability boundary of ±0.3 units
Validated that even with stochastic noise, most initial conditions deterministically lead to capture or escape