<p>Recently, I designed an n-body engine and tested it with a binary star system. After a minute of running, the orbit began noticeably expanding:</p>
<p><a href="https://i.sstatic.net/pejhJ.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/pejhJ.png" alt="enter image description here" /></a></p>
<p>Why is it happening, and is this natural? The source code is based on <a href="https://github.com/lemmingapex/N-bodyGravitySimulation" rel="nofollow noreferrer">https://github.com/lemmingapex/N-bodyGravitySimulation</a>.</p>

CC BY-SA 4.0
2021-05-05T16:55:14.027Z
Did you use an Euler integration scheme? If so, you probably want to use a Leapfrog. The latter will conserve energy and angular momentum on average, the former doesn't. Welcome to the wondrous world of numerical artifacts...
AtmosphericPrisonEscape
CC BY-SA 4.0
2021-05-05T17:13:48.463Z
That definitely looks like energy conservation failure due to a non-symplectic integrator. If you try a smaller time step, and less eccentricity, it will improve things, but eventually stuff will start spiralling. I also recommend Leapfrog, especially the variation by Yoshida (briefly described on the Wikipedia page). You could also use Verlet, which is also a symplectic integrator, like Leapfrog. That readme of that GitHub repo says it has Verlet.
PM 2Ring
CC BY-SA 4.0
2021-05-05T17:54:01.800Z
@AtmosphericPrisonEscape I am using Runge-Kutta 4 as my integrator.
WarpPrime
CC BY-SA 4.0
2021-05-05T21:30:52.940Z
Take into account what we said. RK4 can still easily fail with large timesteps. Don't start with a binary, simplify even more, do a extreme mass-ratio binary (i.e. testmass and star). Then check that a circular and an elliptic orbit are conserved. Once you have established that, you can go and complicate things again. Also see how things change when using the verlet integrator. If things look odd with verlet, then the program you are using is inherently flawed.
AtmosphericPrisonEscape
CC BY-SA 4.0
2021-05-05T21:51:21.070Z
@AtmosphericPrisonEscape Turns out, objects in highly eccentric orbits tend to have massive changes after periapsis (both increase and decrease). I think the timestep is too high, which is causing the issues. Any workaround?
WarpPrime
CC BY-SA 4.0
2022-02-12T00:29:27.260Z
If you're using RK4 you are doing something wrong. RK4 doesn't behave anywhere close to that badly if used correctly. If you're using double precision, you want about 4000 to 6000 steps per orbit for a circular orbit, more for an eccentric orbit. It's always a good idea to vary the step size to see where the optimal step size lies.
David Hammen
<orbit><binary-star><n-body-simulations>
2-body orbits expanding in nbody simulation?
WarpPrime
# 2-body orbits expanding in nbody simulation?

(Asked by: username_0 on 2021-05-05)

Recently, I designed an n-body engine and tested it with a binary star system. After a minute of running, the orbit began noticeably expanding:

[![enter image description here](https://i.sstatic.net/pejhJ.png)](https://i.sstatic.net/pejhJ.png)

Why is it happening, and is this natural? The source code is based on <https://github.com/lemmingapex/N-bodyGravitySimulation>.

**Tags:** orbit, binary-star, n-body-simulations

> **Comments:**
> **username_1**: Did you use an Euler integration scheme? If so, you probably want to use a Leapfrog. The latter will conserve energy and angular momentum on average, the former doesn't. Welcome to the wondrous world of numerical artifacts...  
> **username_2**: That definitely looks like energy conservation failure due to a non-symplectic integrator. If you try a smaller time step, and less eccentricity, it will improve things, but eventually stuff will start spiralling. I also recommend Leapfrog, especially the variation by Yoshida (briefly described on the Wikipedia page). You could also use Verlet, which is also a symplectic integrator, like Leapfrog. That readme of that GitHub repo says it has Verlet.  
> **username_0**: @username_1 I am using Runge-Kutta 4 as my integrator.  
> **username_1**: Take into account what we said. RK4 can still easily fail with large timesteps. Don't start with a binary, simplify even more, do a extreme mass-ratio binary (i.e. testmass and star). Then check that a circular and an elliptic orbit are conserved. Once you have established that, you can go and complicate things again. Also see how things change when using the verlet integrator. If things look odd with verlet, then the program you are using is inherently flawed.  
> **username_0**: @username_1 Turns out, objects in highly eccentric orbits tend to have massive changes after periapsis (both increase and decrease). I think the timestep is too high, which is causing the issues. Any workaround?  
> **username_3**: If you're using RK4 you are doing something wrong. RK4 doesn't behave anywhere close to that badly if used correctly. If you're using double precision, you want about 4000 to 6000 steps per orbit for a circular orbit, more for an eccentric orbit. It's always a good idea to vary the step size to see where the optimal step size lies.  

---

## 0 answers

