Realistic planetary orbits

An easy way would be to use the 2D parametric representation of an ellipse, and then transform into 3D space. Given any time t, where t starts at 0 and ends at 2ϖ, the you can get the 2D location (x, y) using these two equations:

x = a·cos(t)
y = b·sin(t)

where a and b are the major and minor radii. So plug in a value for t, and you get (x, y).

From there, you can translate and rotate the ellipse any way you like. The easiest way would be to treat the y coordinate that you obtained above as though it were the z coordinate and treat your 3D y coordinate as 0 (so that your ellipse is lying flat in 3D space). So you would have the 3D vector <a·cos(t), 0, b·sin(t)>. Then apply the relevant 3D matrix transformations to that 3D vector to translate it so that one of its two foci are at the center, tilt it a bit, rotate it some. (Here are a couple of links that I quickly found on the subject: GameDev: 3D Matrix Math Demystified, Vector Math for 3D Computer Graphics.)

I would use the polar form of the conic section. This way, a focus of the conic section is automatically at the origin. And in case you didn't know, the Sun is located at a focus of a planet's orbit. (Well, technically not necessarily, its only the CG of the system that is located there, but the Sun is so much more massive that for all practical purposes the Sun is located there.)

The polar form in question being:

r=e*k/(1 + e*cos theta)

Where e is the eccentricity of the orbit, and k is the directrix of the ellipse. k is related to the major semiaxis, a, by:

k=a*(1/e - 1)

Or more relevantly, the perihelion, r_perihelion, is related to k by:

k=r_perihelion * (1/e + 1)

Well, forgot to mention, that obviously only gives you the shape of the orbit, but not the dynamics of the orbit. But the dynamics are also fairly easy to handle. Just calculate the total energy of the planet at perihelion and conserve that quantity over the rest of the orbit. The magnitude of the velocity can then be calculated from the kinetic energy portion of the total energy.

You would like to use 3D coordinate system with right angle axis, or you'd be SEVERELY shooting yourslef into the feet.

Planetery systems don't orbit around a star on elipses. It could be mathematically proven that these curves are not elipses, because if they would, you would enter null energy state. Google for Literka and planetary orbits, if you'd like more detailed explanation. Mathworld would give you correct equation for two body orbit. I hope you are able to implement iterative computation of required variable. (piece of cake) Google for VIVO equation to get correct speed.
Of course this isn't enough if you'd like to use two stars, or if you'd have a giant planet close. Correct precession of orbit in nonhomogenous gravity field is possible to do by integration by some smart integrator. There is need to correct integrated values (note that for 100 body system pure conservation of energy isn't enough). I'd recommend a sympletic integrator of higher order. Though for a simple simulation, v = v0 + a foloved by pos = pos0 + v would be enough, or not.
If you are doing 3D simulation, you'd be mostly in your local coordinate system. The biggest problem might be moving from planet to other.

@Agony: Your solution looks easy and straightforward, thanks a lot. Although I have been playing around a little with matrixes already, the "3D Matrix Math Demystified" link looks like the perfect opportunity to finally improve my understanding of the subject.

@Mastaba: Yes I was aware of the sun being in the focus. But to be honest, by the time of writing this reply I haven't quite grasped your math yet. But thanks none the less, I'm having a closer look at it right now.

@Raghar: I'm trying to model the Sol star system (our own). I know the orbits aren't perfect ellipses, but shouldn't they be sufficiently accurate for a game?

Quote:

Original post by transwarp @Mastaba: Yes I was aware of the sun being in the focus. But to be honest, by the time of writing this reply I haven't quite grasped your math yet. But thanks none the less, I'm having a closer look at it right now.

My method is very similar to that of Agony's. It is simply in polar form so that the Sun automatically gets placed at the origin. The conversion from polar to rectangular coordinates is straight forward.

x=r*cos theta
y=r*sin theta

Then I would rotate the XY plane into position like Agony said.

Many planetary tables list the eccentricity of the orbits of the planets, as well as the perihelion of each.

@Mastaba: Thanks to your additional description, I now understand it.
Your formula seems to be the perfect all-round solution. Thanks a lot! :)
One last stupid question: Theta is the angle from the focus pointing towards the planet, right?

Quote:

Original post by transwarp One last stupid question: Theta is the angle from the focus pointing towards the planet, right?

Correct.

I was interrupted, be phone call for help against spam virus, so I didn't add.

Note if you are doing simulation, planets would be highly static.
The second problem would be scale. It's quite hard to put all planets of this system on an one 10 m sized screen. Difference between Mercury, and Pluto orbits is quite high. 100:1

Re your question. It might, however the slightly more correct version is no brainer to implement, and if you need to have correct planetary speeds, you might need it anyway. Of course if it would look nice anyway and you will not have time to correct it... However you said you'd like to implement planets on different planes thus it would be somewhat more fancy.

BTW don't forget the Sun is rotating in different plane than the Earth.

Please try out the "Planetary Gravitation Movement Demo" at

http://www.game-physics-engine.info/4.html (source code included)

Quote:

The second problem would be scale. It's quite hard to put all planets of this system on an one 10 m sized screen. Difference between Mercury, and Pluto orbits is quite high. 100:1

Yes, I've been thinking about that problem as well (this program gives you a very nice overview to the topic: http://www.shatters.net/celestia/). However I think I have come up with a solution:
Much of the game-play is planned to be done in those large empty spaces between the planets. The normal player view will only show the closer surroundings and thus never more than one planet at a time.
Additionally there will be an overview map which will use real distances as well but will scale up the planets by huge factor to make them visible.

Quote:

Re your question. It might, however the slightly more correct version is no brainer to implement, and if you need to have correct planetary speeds, you might need it anyway.

I don’t quite get that. Why do I need to calculate the sun rotation to get the correct planetary speeds?

Quote:

Original post by AlexanderRudolph Please try out the "Planetary Gravitation Movement Demo" at http://www.game-physics-engine.info/4.html (source code included)

I just had a look at the physics engine and found it very impressive. However I am afraid it might be a bit of an overkill for my game, which only requires very basic physics calculations.

@Mastaba:
I have now implemented your formula in my simulation and it seems to be working fine. However I have two questions:

1. As far as I know, a planet's velocity increases as it gets closer to the perihelion and decreases as it moves back to the aphelion.
Using your formula and a constant angular velocity the behaviour is the opposite. So I guess I need to continuously adapt/change the angular velocity. Is that correct and if so, how can I calculate the angular velocity of an orbiting body?

2. When I insert zero for the eccentricity (which would be a perfect circle) the radius simply stays zero because of the multiplication. This is no real problem, since I wont use any perfect circles and even if I did I could use a separate, simpler formula for that. I'm just curious why a formula for ellipses doesn't cover circles as well.

[Edited by - transwarp on August 8, 2006 10:55:11 AM]

With regards to question 1, you are correct. Like I wrote earlier, this gives you the shape of the path, but not the dynamics of the path. To get the speed right you will need to conserve the total energy of the planet as it orbits. Or since you really only need it to look correct, the linear velocity would need to be inversely proportional to the square root of r. Just choose some appropriate scale factor and divide it by the square root of r. Now, to get the angular velocity from the linear velocity, just divide it by r (not exactly correct, but since most of your orbits will be nearly circular it should be good enough), and then the change in theta would be your time delta * angular velocity.

With regards to question 2, r should not go to 0 when e = 0 (well, it does if you didn't replace k with one of the formulas for it). But if you solve for k in terms of the perihelion, and substitute that into the formula for r, and then set e to 0, r should equal the perihelion.

r = r_perihelion * (1 + e) / (1 + e * cos theta)

1. Thanks for the hint. I think I will go for the "inversely proportional to the square root of r" method.

2. Sorry for wasting your time on that one.
I simply didn't turn
r = e*p*(1/e+1)/(1+e*Cos(Theta))
into
r = p*(1+e)/(1+e*Cos(Theta))
like you did on your last post. Thanks once more for the support.

I'd suggest checking those links:

Computing planet positions using mean orbital elements

Positions of the Sun and the planets

Wikipedia page on orbital elements

The first link in particular should be enough to implement "realistic" planetary orbits, and it even provides code for doing so.

Y.

@Ysaneya: Thanks for those links, I'll have a close look at them soon.

In the mean time I'm using this, not quite perfect, but functioning code:

// Universal gravitational constantconst float G = 6.672e-20f;// Calculate current position within 2D-planeradius = periapsis * (1 + eccentricity) / (1 + eccentricity * Math.Cos(Theta));x = radius * Math.Cos(Theta);y = radius * Math.Sin(Theta);// Update angular velocityv = Math.Sqrt(G * orbitCenter.Mass / radius);theta += ElapsedTime * v / radius;

Notes:
ElapsedTime = DeltaTime
Radius = Distance from Sun
orbitCenter.Mass = The sun's mass

Look at this site first. It has nice Java applets that would show you what you can do. It has also descriptions of correct integrations with a limits on them.

http://www.burtleburtle.net/bob/

At least one of links down on this page, will get you to place where you can look at speed equations.

http://en.wikipedia.org/wiki/Orbit

Majority of people did

angleO = angleO + step


E = angleO + e*sin(E) .... kepler equation


Then they converted E to real coordinates, and used VIVO equation for real speed. It might be bad in reality, but it doesn't drift away as real orbits do, and it's quicker than integration.


re Ysaneya

Considering how he was uneasy when confronted with kepler equation, he might be perfectly happy with correctly spaced kinda elliptic orbits.