# Simulations of Planetary Motions

I wrote a spreadsheet that simulates the trajectory of 3 planets in 2D space. The method is simple: for each moment in space, calculate the force felt and velocity of each planet, then for the next $\delta$t period, assume that the planet continues to move with the same acceleration. Then recalculate the position, velocity, and acceleration of each object. Proceeds until desired time frame is simulated.

The computation works roughly correct, but there is a problem with it, that is: I can not produce closed curves typical seen in a planet's orbit. Regardless of what velocity and mass I chose, the orbit's radius always gradually increases.

This is probably caused by the approximation method I chose. I am wondering, what are some improvements I can make to the numerical method?

-
See related scicomp question: scicomp.stackexchange.com/q/3566/3224 –  Chris White Apr 4 '13 at 1:55

Your method is known as the forward Euler method, the simplest DEQ-solver but a really bad one. I suppose it's the approximation everyone tries first, at least I did, too, when I started with such simulations... and soon was troubled by similar "numerical explosions", which can't really be avoided with this algorithm which is why it's virtually never used in real applications.

The most standard, "proper" alternative is the 4th-order Runge-Kutta solver, I suggest you read up on that. It's going to become a little unwieldy in a spreadsheet (you may want to switch to a proper programming language; I'd recommend Python or Haskell for such stuff).

-
Be careful - while Runge-Kutta methods have their uses (like non-periodic N-body simulations), they are not entirely stable and will gain or lose energy at least linearly with time in 2-body orbital problems. Still, they are better than forward Euler, and you recommended switching to a proper language, so +1. –  Chris White Apr 4 '13 at 2:21