Elliptical orbit modeling

I'm playing with orbits in a simple 2-d game where a ship flies around in space and is attracted to massive things. The ship's velocity is stored in a vector and acceleration is applied to it every frame as appropriate given Newton's law of universal gravitation. The point masses don't move (there's only 1 right now) so I would expect an elliptical orbit.

Instead, I see this:

I've tried with nearly circular orbits, and I've tried making the masses vastly different (a factor of a million) but I always get this rotated orbit.

Here's some (D) code, for context:

``````void accelerate(Vector delta)
{
velocity = velocity + delta; // Velocity is a member of the ship class.
}

// This function is called every frame with the fixed mass. It's a
// method of the ship's.
void fall(Well well)
{
// f=(m1 * m2)/(r**2)
// a=f/m
// Ship mass is 1, so a = f.
float mass = 1;
Vector delta = well.position - loc;
float rSquared = delta.magSquared;
float force = well.mass/rSquared;
accelerate(delta * force * mass);
}
``````
-
woo. yeaaaah. D. Unit test that math code against known results; and all will be well. – dcousens Jun 23 '12 at 6:43

2 Answers

The bug is in the `fall` function. We have

1. `delta`: a vector from the well to the ship
2. `force`: the magnitude of the gravity between these two bodies.

`|force|` is G * m1 * m2 / r^2

but `|delta|` is already r! so you are actually accelerating too fast. You need to divide by `r` again (basically normalizing the delta vector) before calling accelerate.

``````accelerate(delta * well.mass * mass / rSquared / Math.sqrt(rSquared))
``````
-

Note that even with the math bug(s) fixed, you're using Euler integration (i.e. `velocity += delta` and presumably `position += velocity`), so you're probably going to get some odd effects like rotation of the orbital ellipse over time, and perhaps the ellipse getting larger/smaller since Euler integration isn't guaranteed to conserve energy.

You might want to switch to leapfrog integration, which is energy-conserving and should work better for orbital mechanics. Also, you should include the frame time in your equations so that the speed of your simulation is framerate-independent.

-