Skip to main content
1 of 5

How to numerically integrate Kepler Problem?

I was "redirected" from Physics Stack Exchange to this place - it is my first question:

https://physics.stackexchange.com/questions/739224/how-to-numerically-integrate-kepler-problem

Again (but note: I'm not a mathematician and this is really a naïve approach):

I tried to solve a simple Kepler Problem numerically.

I have discrete time steps, a starting position $(x_0,y_0)$ and starting velocity $(u_0, v_0)$.

I used this iteration by calculating the forces and apply Newton's law:

$$F = k/r^2 = k/(x^2+y^2) $$ $$F_x = F \cdot x/\sqrt{x^2+y^2}$$ $$F_y = F \cdot y/\sqrt{x^2+y^2}$$

$$u_{i+1} = u_{i}+ F_x/m \cdot \Delta t$$ $$v_{i+1} = v_{i}+ F_y/m \cdot \Delta t$$

$$x_{i+1} = x_{i}+ (u_{i+1}+u_i)/2 \cdot \Delta t$$ $$y_{i+1} = y_{i}+ (v_{i+1}+v_i)/2 \cdot \Delta t$$

I used mean velocity between last and actual time to calculate new (x, y).

Naively, I thought, that by making the timestep small enough this gives good results. But I observed, that this method is completely bad, because the total energy tends to increase ant what I get is not what I expect: enter image description here

What method for integrating such problem is the standard? My idea is obviously complete garbage because it is numerically most unstable. When I calculate total energy it is not constant at all, even when I make the timestep incredible small, I get totally dissatisfying results.