Simple orbit simulation problem

[fluoZor]

Hello, i have been trying to program REALLY simple simulation of an orbit, i currently just use a simple kinectic energy vector that i give originally some speed and then i just proceed to add the force of gravity to it on every step.

Only problem i have with this system is that orbit is always centered on the domain body and there seem to be 2 apoapses and 2 periapses, like this

What force or variable have i failed to account for when this always keeps happening?

[metaphor]

Are you using Cartesian or radial coordinates? If you just have velocity and you add the acceleration from gravity to it at every step you should have an actual orbit. Kinetic energy is not a vector, it doesn't have a direction.

[tavert]

I think you can get a centered ellipse like that (instead of the center of the planet being at a focus of the ellipse, like it should be) when you assume a constant gravity force instead of an inverse-square force magnitude. Or maybe it's when you leave out the Coriolis term in polar coordinates, so the angular velocity is constant instead of the angular momentum.

Edited June 22, 2013 by tavert

[K^2]

Ellipse like that corresponds to a body in central harmonic potential. I suspect it's due to the way you are dealing with force as a vector, but you'll have to paste a portion of your code to say with certainty where the problem is.

[fluoZor]

Okay, i shall now explain the system i have used more clearly. Because of the simple nature of this i use two components to represent to motion on the 2D coordinate system, these are hspeed for horizontal speed and vspeed for vertical speed. In the beginning of the program i have hspeed on the value of 0.7 so that the craft will have its motion to keep it on the orbit.

Every step of the program i calculate the gravitational force and calculate its components and just simply add the components of the hspeed and vspeed variable. Using these hspeed and vspeed variables i move the spaceship every step.

This is the code i use to calculate the gravitational force:

grav=((Gcon*Mse*(5.9*10^24))/(distance_to_point(earth.x,earth.y)*50)^2)

Gcon is the gravitational constant and Mse is the mass of the craft.

[softweir]

Are you remembering to calculate A=F/m where A is the acceleration due to gravity, F is the gravitational force and m is the mass of the spaceship? You don't add force to velocity, you add acceleration to velocity!

Another way is to calculate A directly as:

Acc = Gcon*5.9*10^24 / (distance_to_point(earth.x,earth.y)*50)^2

.... in other words, just leave out the mass of the spaceship.

[fluoZor]

I did the things above and this still doesn't seem right

at the moment i have no idea what force or variable i have forgot to account for

[Epthelyn]

You should be able to keep the mass of the spacecraft. After all, although the mass of the spacecraft is negligible relative to the planet, if you want to replace it with larger object later you're better off leaving it in and getting it to work. Are you sure the problem you're having isn't due to some other part of the simulation? The formula is correct.

[fluoZor]

The problematic thing is, that this is everything that there is to it at this state, there is just a static planet and an orbiter that only does get its speed at the start and counts the gravity to itself and applies its acceleration.

[maltesh]

Without seeing the rest of your code, it's going to be pretty much impossible to diagnose.

That said, you aren't going to get an absolutely perfect ellipse using the method you detailed, even if you implement everything correctly. The smaller your timestep, the closer the resulting path will be to a approximating a Keplerian orbit, but the step-and-calculate method is always going to result in shifting of the orbit over time.

[fluoZor]

I actually think i have found the problem here, which is that the program was cutting corners and variables got overflown and the program counted errors.

Now i either need a better way of simulating this or some good way to cut down the numbers and have it still work rather same.

[maltesh]

The question is: What are you trying to simulate?

Are you ultimately headed for multi-body simulation, where every body gravitationally interacts with every other one?

Or do you want to just put non-gravitating satellites in the gravitational field of the primary body, generate the Keplerian orbit that is defined by their initial position and velocity,and work out their positions as a function of time based on that orbit?

The former gets somewhat complicated for long-term accuracy. The latter is mostly straightforward.

[fluoZor]

I am just trying to have an simple satellite orbiting one domain object, so the latter

[K^2]

What are the exact formulae you evaluate on each step. This isn't just numerical mess, you definitely aren't computing relevant quantities correctly.

[fluoZor]

I am sorry everyone, the problem i had here was simply that i had too big numbers in these calculations, i compensated every number down and now have correct orbits.

Thanks for the help.

[MBobrik]

I don't claim to know what is wrong with your simulation ( you would have to show your code to be able to), but the "naive" method of numerical integration ( simply summing up acceleration x timestep into velocity and velocity x timestep into position ) has very poor stability when applied to orbital dynamics. And as a result, the orbit will always drift (though most probably not much as what you've shown ), even for small timesteps. Google some numerical methods, the leapfrog method is a bare minimum.

[PakledHostage]

the problem i had here was simply that i had too big numbers in these calculations, i compensated every number down and now have correct orbits.

Congratulations on making the effort! I'm glad to hear that you figured it the source of your bug. An article you might be interested in is "Lies my Calculator and Computer Told Me". It is a bit old, but it still relevant to the issue you experienced.

[fluoZor]

I don't claim to know what is wrong with your simulation ( you would have to show your code to be able to), but the "naive" method of numerical integration ( simply summing up acceleration x timestep into velocity and velocity x timestep into position ) has very poor stability when applied to orbital dynamics. And as a result, the orbit will always drift (though most probably not much as what you've shown ), even for small timesteps. Google some numerical methods, the leapfrog method is a bare minimum.

I know it is not totally accurate but for my use, i do well enough using this method and forcing the craft into same position and speed on the end of each orbit to keep it from going too far. I am making space RTS with semi realistic orbit mechanics and try to avoid having too complex calculations for the units.

[K^2]

the leapfrog method is a bare minimum.

Leapfrog/Verlet only marginally improve stability of the orbital simulation, unfortunately. If you want to improve precision significantly, you have to go for higher order implicit RK integration, and that's a royal mess. And even then, you are going to have energy drifts. As far as I know, a general conservative method for n-body 1/r potential is not known.

Of course, for a 2-body problem, you can always just fall back on analytic solutions.

[grawl]

Here is a thread that might help you out.

http://forum.kerbalspaceprogram.com/showthread.php/32590-Maths-Physics-Equation-of-an-orbit

[metaphor]

If you just have a single massive body generating gravity, you can calculate where a satellite will be at any time if you know its position and velocity at any other time, using the orbit equations. I think this is what KSP uses, which is much faster than numerical integration and will give a stable unchanging orbit over any amount of time.