Comparing Different Numerical Integrators for the N-Body Gravity Problem

Comparing Different Integrators for the N-Body Gravity Problem

In previous blogs, I have shown how to apply Parker-Sochacki to the 2-body, 3-body, and n-body problem. In this blog, we will show different results of different integrators on the n-body problem. We start with showing the system of ODEs:

n_body

 

What you will see next are 5 different test cases where I show plots that update the position in 2 dimensions. I wanted to get a visual idea of habits that the following different methods have:

Forward Euler, Semi-Implicit Euler, Velocity Verlet, Runge-Kutta 2 Midpoint Method, Runge-Kutta 4, and Parker-Sochacki. When I ran tests, I only modified the time step except for Parker-Sochacki. It was a fixed time step for the entire simulation. To be consistent with Unity, I started with a time step of 1/30. If no method could achieve that large of a time step, then I lowered it. There are a lot of images coming up so if you want a quick summary, I have it here in General Notes. In every plot I will show how many bodies there were(2 or 3), the method, the time step, and how many frames were needed to run the simulation. If a time step needed to be cut in half, then the total frames were doubled.

The order of images will go as such. It will start with Parker-Sochacki where only accuracy is adjusted to get the desired image. Then I will show images of Forward Euler, Semi-Implicit Euler, Velocity Verlet, Runge-Kutta 2 Midpoint, and Runge-Kutta 4 where the next image has a smaller time step achieving a more accurate result. For the non-Parker-Sochacki methods, a smaller time step yielded a more accurate result.

General Notes:

————————————————————————–

This test does not check for flexibility of constants(we are hoping that the gas constant will have more flexibility)

This test does not measure how long each simulation took. A larger time step does not mean that the simulation was done in less time.

This test does not check for conservation of energy. I spent all of Friday learning how to calculate the total energy and did not get results that made sense. But I think I have an explanation. To get choreographic patterns for different bodies, I set G to 1 from initial conditions I saw online. G should not be 1 if this was “real life.” I created a separate report that checks for Conservation of Energy for the Simple Pendulum model.

There is a limit to how large a time step can be for a simulation. In practice, there were simulations where at a time step of 1/30, Parker-Sochacki order 30 still blew up.

Parker-Sochacki is the only method where I can adjust the “order” or accuracy of the simulation without adjusting the time step.

With Parker-Sochacki, generally a higher order improves the accuracy but not always. Look at the third example. Parker-Sochacki 2nd Order performs better than Parker-Sochacki 3rd Order(both are still wrong). The hope is that when Parker-Sochacki is at an order higher than 4(I say 4 because that is where Runge-Kutta 4 goes) we will achieve results in SPH simulations that other methods cannot.

Semi-Implicit Euler and Velocity Verlet are more stable than Forward Euler, Runge-Kutta 2 Midpoint, and Runge-Kutta 4 but that does not make them “better” or more accurate. There were cases you will see where the simulation did not blow up but the movement was “wrong”

Parker-Sochacki at first order performs the same accuracy as Forward Euler(this is expected). Parker-Sochacki at fourth order performs about the same accuracy as Runge-Kutta 4.

I did not take time steps larger than 1/30. While there were cases where more than one method could achieve this, Parker-Sochacki(at a higher order) could take the largest time steps of any method.

If I were to base the success of SPH from this experiment I would say the following. Forward Euler needs extremely tiny time steps and it still does not keep up. Semi-Implicit Euler is a much better choice(although this is what you have already). Runge-Kutta 2 has 3 versions: Midpoint(what I have here), Heun, and Ralston. I first thought about putting those in there as well but I don’t see any of the 5 examples where Runge-Kutta 2 Midpoint Method does well when I have Velocity Verlet in there. It won’t hurt to put in SPH just to have it(to make Runge-Kutta 4 I need to go through Runge-Kutta 2 Midpoint Method anyway). Velocity Verlet is more stable than Runge-Kutta 4 so I think both of those methods should be options for SPH. Parker-Sochacki is meant to be a higher order method so I hope that at higher order methods it will outperform everything else.

 


2-BODY PROBLEM WHERE BOTH PARTICLES ORBIT EACH OTHERPS.PNGForward_Euler1.PNGForward_Euler3.PNGForward_Euler4.PNGForward_Euler5.PNGForward_Euler6.PNG

Semi_Euler1.PNGSemi_Euler2.PNGSemi_Euler3.PNGVelocity_Verlet1.PNGVelocity_Verlet2.PNGVelocity_Verlet3.PNGRK2_Midpoint.PNGRK2_Midpoint2.PNGRK2_Midpoint3.PNGRK2_Midpoint4.PNGRK2_Midpoint5.PNGRK2_Midpoint6.PNGRK2_Midpoint7.PNGRK4.PNGRK4_2.PNG

FUNKY 3-BODY PATTERN:PS.PNGPS_2.PNGPS_3.PNGForward_Euler_1.PNGForward_Euler_2.PNGForward_Euler_3.PNGForward_Euler_4.PNGForward_Euler_5.PNGForward_Euler_6.PNGForward_Euler_7.PNGSemi_Implicit_Euler_1.PNGSemi_Implicit_Euler_2.PNGSemi_Implicit_Euler_3.PNGSemi_Implicit_Euler_4PNG.PNGSemi_Implicit_Euler_5.PNGVelocity_Verlet.PNGVelocity_Verlet2.PNGVelocity_Verlet3.PNGVelocity_Verlet4.PNGRK2_Midpoint_1.PNGRK2_Midpoint_2PNG.PNGRK2_Midpoint_3PNG.PNGRK2_Midpoint_4.PNGRK2_Midpoint_5.PNGRK4_1.PNGRK4_2.PNG

 

 

3-BODY WHERE TWO BODIES ROTATE AROUND A HEAVY MASS(THE THIRD MASS IS A DOT IN THE CENTER THAT DOES NOT MOVE):PS_8.PNGPS_1.PNGPS_2.PNGPS_3.PNGPS_4.PNGPS_5.PNGPS_6.PNGPS_7.PNG

Forward_Euler_1.PNG

Forward_Euler_2.PNG

Forward_Euler_3.PNG

Forward_Euler_4.PNGForward_Euler_5.PNG

Forward_Euler_6.PNG

Forward_Euler_7.PNG

Forward_Euler_8.PNG

Forward_Euler_9.PNG

Forward_Euler_10.PNG

Forward_Euler_11.PNG

Forward_Euler_12.PNG

Forward_Euler_13.PNG

Forward_Euler_14.PNG

Forward_Euler_15.PNG

Semi_Implicit_Euler_1.PNG

Semi_Implicit_Euler_2.PNG

Semi_Implicit_Euler_3.PNG

Semi_Implicit_Euler_4.PNG

Semi_Implicit_Euler_5.PNG

Semi_Implicit_Euler_6.PNGSemi_Implicit_Euler_7.PNG

Semi_Implicit_Euler_8.PNG

Semi_Implicit_Euler_9.PNG

Semi_Implicit_Euler_10.PNGVelocity_Verlet_1.PNG

Velocity_Verlet_2.PNG

Velocity_Verlet_3.PNG

Velocity_Verlet_4.PNG

Velocity_Verlet_5.PNG

Velocity_Verlet_6.PNGRK2_Midpoint_1.PNG

RK2_Midpoint_2.PNG

RK2_Midpoint_3.PNG

RK2_Midpoint_4.PNG

RK2_Midpoint_5.PNGRK2_Midpoint_6.PNG

RK2_Midpoint_7.PNG

RK2_Midpoint_8.PNG

RK4_1.PNG

RK4_2.PNG

RK4_3.PNG

RK4_4.PNG

3-BODY FIGURE 8:

PS_1.PNG

PS_2.PNG

PS_3.PNG

PS_4.PNG

PS_5.PNG

PS_6.PNG

Forward_Euler_1.PNG

Forward_Euler_2.PNG

Forward_Euler_3.PNG

Forward_Euler_4.PNG

Forward_Euler_5.PNG

Forward_Euler_6.PNG

Forward_Euler_7.PNG

Forward_Euler_8.PNGForward_Euler_9.PNG

Forward_Euler_10.PNG

Forward_Euler_11.PNG

Forward_Euler_12.PNG

Forward_Euler_13.PNG

Forward_Euler_14.PNG

Semi_Implicit_Euler_1.PNG

Semi_Implicit_Euler_2.PNG

Semi_Implicit_Euler_3.PNGSemi_Implicit_Euler_4.PNG

Semi_Implicit_Euler_5.PNG

Semi_Implicit_Euler_6.PNG

Semi_Implicit_Euler_7.PNG

Semi_Implicit_Euler_8.PNGVelocity_Verlet_1.PNG

Velocity_Verlet_2.PNG

Velocity_Verlet_3.PNGVelocity_Verlet_4.PNG

Velocity_Verlet_5.PNG

RK2_Midpoint_1.PNG

RK2_Midpoint_2.PNG

RK2_Midpoint_3.PNG

RK2_Midpoint_4.PNGRK2_Midpoint_5.PNG

RK4_1.PNG

RK4_2.PNG

 

 

3-BODY ROTATING FIGURE 8

PS_1.PNG

PS_2.PNG

PS_3.PNG

PS_4.PNG

Forward_Euler_1.PNG

Forward_Euler_2.PNG

Forward_Euler_3.PNG

Forward_Euler_4.PNG

Forward_Euler_5.PNG

Forward_Euler_6.PNG

Forward_Euler_7.PNG

Forward_Euler_8.PNG

Forward_Euler_9.PNG

Forward_Euler_10.PNG

Semi_Implicit_Euler_1.PNG

Semi_Implicit_Euler_2.PNGSemi_Implicit_Euler_3.PNG

Semi_Implicit_Euler_4.PNG

Semi_Implicit_Euler_5.PNG

Velocity_Verlet_1.PNG

Velocity_Verlet_2.PNG

Velocity_Verlet_3.PNG

RK2_Midpoint.PNG

RK2_Midpoint_2.PNG

RK2_Midpoint_3.PNG

RK4_1.PNG

RK4_2.PNG