# Is there a way to reduce aberration in computations of planets' trajectories?

I don't think the title is very accurate , sorry for that.

I simulate bodies in space using two timestep:

the TIMESTEP is the Δt wich I use to make the calculation and XTIME is the number of times I make the calculation

I want the each visible step to be 7 days , the Δt will be 86400 (seconds) and XTIME 7 , so it'll be calculate seven days.

I made that because if $Δt = 86400*7$ there is a lot of aberration (especially on moon orbit).

Since I'll use less than 20 planets each time , calculate 7 times Δt is not really a problem , but more could be.

And this is the problem , with $Δt = 86400$ I still have aberation when planet come too close to the massive mass (the sun), the result is the planet going away at the speed of light ! It's because the Δt is too large to calculate the "counter force" .

I heard of leapfrog but I must be stupid , I didn't succeed to adapt it to my code (the result was worst than anything, like the moon going away from earth) and i'm don't sure this will solve my problem.

I precise that I want something that look believable but not realistic , so if there is not "not-eating-cpu" method , I'll certainly going to cheat a little with those cases, but I prefer not to.

I ask this question here because i think it's more computing than physics since I already got the physics

edit: maybe a lead ? something like if (distanceToOldPosition > FIXED_VALUE) calculate more times

edit2: precision on calculations asked on comment: \begin{align} F = \frac{GM_{1} M_{2}}{r} \end{align}

the acceleration \begin{align} a = \frac{F}{M} \end{align}

then the linear vector \begin{align} v_{t+1} = v_{t} + aΔt \end{align}

and finally the new position

\begin{align} x_{t+1} = x_{t} + v_{t+1}Δt \end{align}

• welcome to SciComp. Could you be more specific in what algorithm you use to integrate the movement of the planets? Jan 31, 2013 at 9:36
• sorry for the delay , I added the formulas Jan 31, 2013 at 9:50
• I am going to guess you are having problems due to the instability of Euler's method with your large time step. If possible, could you switch to backwards Euler to see if that improves your results? Jan 31, 2013 at 16:18
• After reading docs on backward euler , I don't really get how to adapt it to my formulas Jan 31, 2013 at 16:50
• Feb 8, 2013 at 13:34

To get something that looks realistic for planetary orbits, you shouldn't use the forward or backward Euler methods. These will cause your planets to spiral outward or inward. You should use a symplectic method.

You may also need to adjust the timestep to be smaller when two bodies are very close to each other.

Read Chambers (1999) A hybrid symplectic integrator that permits close encounters between massive bodies to get started with such methods.

• If you want a simple one, go for leapfrog integration. May 15, 2013 at 15:47

In general I would suggest not writing your own ODE solver. This is a well solved problem, read about how popular ones work and pick one suitable to your problem.

If you let us know what language you are working in I am sure someone can suggest their favorite package. :)

Remember, these are never exact solutions, but some are less wrong that others.

If don't need physics, you can just cheat and make your planets orbit in a perfect circumference or even better in a ellipse.

$$pos_{planet} = circumference(angle(time),radius,center)$$

The moon and other satellites would be a bit different:

$$pos_{satellite} = pos_{planet} + circumference(angle(time),radius,center)$$