Open Events in the Horizon: A Black Hole Visualization Project

Ever since the first actual image of a black hole was released by the Event Horizon Telescope team in 2019, I've been fascinated by these cosmic objects. As a programmer with a curiosity for astrophysics, I challenged myself to create a physically accurate visualization of a black hole that would run in real-time.

This project, which I've named "Open Events in the Horizon" (OEH), combines Python, CUDA, and OpenGL to create an interactive black hole simulation that anyone can run on their computer.

The Science Behind the Visualization

While I'm not an astrophysicist myself, I based my simulation on peer-reviewed research. Specifically, I implemented the magnetically dominated accretion disk model from a 2003 paper by Pariev, Blackman & Boldyrev:

"Extending the Shakura-Sunyaev approach to a strongly magnetized accretion disc model" (Astronomy & Astrophysics, 407, 403-421)

This model accounts for how strong magnetic fields affect the structure and appearance of the disk of matter spiraling into the black hole. The result is a more realistic temperature profile and visual appearance than simpler models.

Technical Components

Creating this simulation involved several technical challenges:

1. Ray Tracing in Curved Spacetime

To visualize a black hole correctly, I had to trace the paths of light rays as they bend around the intense gravitational field. I implemented a modified ray tracer in CUDA that:

Approximates light paths in Schwarzschild geometry
Computes gravitational lensing effects
Calculates relativistic effects like Doppler shifting and redshift

@cuda.jit(device=True)
def trace_ray_equatorial(
    cam_x: float, cam_y: float,
    dir_x: float, dir_y: float,
    mass_bh_cgs: float,
    max_steps: int,
    dt: float,
    horizon_radius: float,
    integrator_choice: int
) -> float:
    """
    Trace a photon in the equatorial plane around a BH with improved physics.
    """
    # Position
    x, y = cam_x, cam_y
    # Velocity (normalized to c)
    vx, dir_x * C
    vy = dir_y * C

    # Physics implementation...

2. Physically Accurate Accretion Disk

I implemented the magnetically dominated disk model, which calculates:

Proper temperature profile based on radius and magnetic field strength
Realistic emission spectrum using blackbody radiation
Electron scattering effects that modify the spectrum
Turbulence and plasma physics at different disk radii

@cuda.jit(device=True)
def disk_radiance(r_cm: float, m_bh: float, nu: float, b_field_exp: float) -> float:
    """
    Returns total specific intensity I_nu from the disk at radius r_cm for frequency nu.

    Includes both blackbody and modified blackbody effects depending on the
    optical depth regime, following the Pariev+2003 model.
    """
    # Get surface temperature at this radius
    T = disk_surface_temperature(r_cm, m_bh, b_field_exp)

    # Calculate effective optical depth to determine emission regime
    r_g = G * m_bh / (C*C)
    ratio = r_cm / (10.0 * r_g)

    # Implementation of physics model...

3. Real-time Rendering Pipeline

The OpenGL rendering pipeline includes:

GPU-accelerated simulation using CUDA and Numba
Real-time post-processing effects
Interactive camera controls
Shader-based effects for bloom, exposure, contrast, etc.

def render_frame(self):
    """Renders a single frame by running the raytracer and displaying the result."""
    # Apply auto-rotation if enabled
    self._update_camera_for_rotation()

    # Run the simulation only if not paused
    if not self.paused or self.last_image is None:
        t_start = time.time()

        # Run the simulation with the current parameters
        image = run_simulation(
            custom_camera_position=self.camera_position,
            custom_black_hole_mass=self.black_hole_mass,
            custom_fov=self.fov,
            b_field_exponent=self.b_field_exponent,
            integrator_choice=self.integrator_choice
        )

        # Update image cache and timing
        self.last_image = image
        self.render_time_ms = (time.time() - t_start) * 1000

Integration Challenges

One of the biggest challenges was integrating the physics, computation, and visualization components. I ran into issues with:

Performance optimization: Finding the balance between physical accuracy and real-time rendering
Memory management: Efficiently transferring data between the CPU, CUDA cores, and OpenGL
Numerical stability: Ensuring the integration methods worked properly in extreme conditions near the event horizon

The Development Process

This project stretched my programming skills and taught me a lot about GPU programming, physics simulation, and scientific visualization. Here's how I approached it:

Research phase: Understanding the physics papers and translating equations into code
Core simulation: Implementing the basic ray tracing and physics calculations
GPU acceleration: Porting the code to CUDA for massive parallelization
Visualization: Creating the OpenGL-based visualization pipeline
Refinement: Continuous improvement of visuals and physical accuracy

I also leveraged modern AI tools (Claude, Gemini Pro, and GPT) to help with debugging and optimization. These tools were particularly helpful in identifying numerical issues and suggesting code optimizations.

Try It Yourself

The code is open source under the GNU license. You can try it yourself:

# Clone the repository
git clone https://github.com/EricsonWillians/OEH.git
cd OEH

# Install with UV (preferred package manager)
uv pip install -e .

# Run the simulation
python -m oeh.main

You'll need a CUDA-capable GPU to run the simulation at full speed.

Conclusion

Creating this black hole visualization was both a technical challenge and a fascinating learning experience. It demonstrates how consumer-grade hardware can now simulate complex physics that previously required supercomputers.

Even if you're not an astrophysicist or specialized graphics programmer, modern tools and libraries make it possible to create impressive scientific visualizations with relatively accessible technology.

If you try out the project or have suggestions for improvement, I'd love to hear from you!

Building a Physically Accurate Black Hole Visualization with Python, CUDA, and OpenGL