Three-Body Problem Simulation (MATLAB)

This project numerically simulates the classical three-body problem using a fifth-order Runge–Kutta (RK4(5)) method. It models the gravitational interactions between three point masses under Newtonian mechanics, allowing users to modify the initial positions, velocities, and masses to observe different orbital behaviors.

---

Overview

The three-body problem studies the motion of three bodies moving under their mutual gravitational attraction. Unlike the two-body problem, there is no general analytical solution, so numerical methods like Runge–Kutta are required to approximate the trajectories over time.

This MATLAB program:

• Implements a custom Runge–Kutta 4(5) integrator with a specified Butcher Tableau.
• Computes the positions and velocities of all three bodies at each time step.
• Ensures the center of mass velocity is zero (for physical consistency).
• Outputs arrays of trajectories that can be plotted to visualize orbital motion.

---

Features

• Fifth-order accuracy (RK4(5))
• Customizable initial conditions
• Conservation of center-of-mass motion
• Configurable Butcher Tableau (alternative version included)
• Long-term trajectory simulation for chaotic systems

---

File Structure

three_body_rk45.m        # Main MATLAB script (this program)
README.md

How It Works

1. Constants and Parameters

• Gravitational constant G is set to 1 (normalized units).
• Time step dt and total iterations n determine the simulation resolution.

2. Initial conditions User can modify

M = [1 1 1];         % Masses of the three bodies
x_1(1) = -1; x_2(1) = 1; x_3(1) = 0;
y_1(1) = 0;  y_2(1) = 0;  y_3(1) = 0;
vx_1 = 0.347113; vx_2 = 0.447113; vx_3 = -0.694226;
vy_1 = 0.532727; vy_2 = 0.532727; vy_3 = -1.065454;

These define the starting configuration and momenta.

3. Runge–Kutta Integration The script uses the RK4(5) method to iteratively compute the the instantaneous position and velocity of each body.

4. Output The script produces arrays

x_1, x_2, x_3, y_1, y_2, y_3 — positions
vx_1, vx_2, vx_3, vy_1, vy_2, vy_3 — velocities

You can visualize the animation in the MATLAB script.

---

Example Visualization

---

Notes

• The simulation assumes dimensionless units (G = 1).
• The total number of steps n = 1e6 may take time, adjust for shorter tests.
• Numerical instability may appear for very small separations or long durations due to chaotic nature of the system.

---

Mathematical references

• Hairer, E., Nørsett, S. P., & Wanner, G. (1993). Solving Ordinary Differential Equations I: Nonstiff Problems.
• Butcher, J. C. (2003). Numerical Methods for Ordinary Differential Equations.

Author

Hngjesse

Undergraduate Project : Three-Body Problem Numerical Simulation using RK4(5) Developed in MATLAB