Getting back to Orbital Mechanics, I've started taking the hard route and doing it the KSP way - writing a separate physics system for them. The problem comes in the interactions between the Orbital Physics and the Rigidbody Physics.
1. How can I make Rigidbody interactions have proper effects on the Orbit? Difficulty: Adjusting the orbit to small collisions and thrust. I have a theory that I should just track the changes in velocity from one frame to the next and consider that kinetic energy, but I'd also have to subtract the orbital motion from that, and oyoyoy...
2. How can I let Orbital changes and Gravity have the right effects on the Rigidbodies? Difficulty: I somehow suspect that rigidbodies don't take kindly to being constantly pushed around by non-unity-physics forces, since I've already seen some strange behaviour when a rigidibody in motion was translated between two frames and it ended up somewhere it didn't expect to be.
I'm not sure about how this should work with Unity in particular, but essentially what you're doing by putting the object on a on-rails orbital track is you're splitting its total motion/physics into two pieces. One piece comes from the solution of the equations of motion for an object under the effect of gravity; the other piece comes from whatever detailed motions the object performs due to collisions/rotation/etc within the context of that orbiting reference frame.
One way to think about the orbiting reference frame is that its a choice of reference frame such that the fictitious force due to the frame's acceleration exactly counter-balances the force of gravity, so the motion of a single rigid object within the orbiting reference frame is in free-fall (e.g. you don't need to worry about gravity because its taken care of by choosing the orbiting reference frame).
Once you have two objects, the problem is that they are in slightly different orbits, so if you go to a shared reference frame containing both objects then the result is that gravity once again matters and needs to be calculated carefully. As far as I'm aware, KSP uses 'on-rails' orbits only when objects are not interacting, and when they are interacting it switches to something in which gravity/etc are fully taken into account - that's why objects sometimes jump around when you enter non-physics-timewarp.
What you can in principle do is to create a 'center of mass orbit frame' when two objects come near to eachother. This is basically the orbit of the average of the objects' motions. However, there isn't that much mathematical benefit in doing so because you still are going to have to take into account differential gravity (e.g. tidal forces) due to different distances within the frame from the primary. Of course if the interaction range is much smaller than the distance to the primary, you may choose to just assume that these tidal forces are zero.
When objects left collision distances of eachother, you'd then recalculate their individual orbits and update their orbital parameters.
3. Orbits require large, massive bodies, which greatly increases the distances in play. Even if it's just an earthlike 15,000 km or so, I'm guessing that the physics simulation will get pretty wonky that far away from the origin...
Double precision means that roughly speaking you've got 12-14 orders of magnitude of smooth resolution and then a few orders of magnitude that may have small, visible glitching. So with double-precision floating points you can resolve motions on the order of millimeters while still retaining distances of the order of 15000km. Its generally better however to separate things into two parts if there's a huge gap in the numbers like that - store the 'orbital-scale' numbers in one variable while storing the local position of sub-objects and nearby objects relative to their shared center-of-mass or something like that. That way, you don't have the issue that a 15000km orbit screws up the millimeter-scale rigid body stuff. Separating out the orbital physics like you're intending to should help you do this, though you'll have to be careful at the point where the two scales intersect (e.g. asking whether or not two orbiting objects are within 1km of eachother to initialize the detailed positional information).
A separate issue is that I have no idea on how to actually calculate orbital mechanics. If anyone knows a good page to read - this non-mathematician would be quite grateful.
I don't want to open the can of worms that is writing an all-new physics engine, unless you really think that's the best way ;P
And obviously none of this is anyone's speciality here, so miscellaneous input is just as appreciated as qualified knowledge
The idea with orbital mechanics is that you can separate out the shape of the orbit (a particular ellipse or hyperbola) from your position along the orbit. Calculating the specific position that corresponds to a specific time still requires numerical integration or solving a transcendental equation, but the benefit is that basically you're always just solving versions of the same problem and so the error doesn't grow with time.
It looks like 'Orbital Elements' is one of the keywords to search for in finding information about this kind of thing. What you're going to want is a way to calculate orbital elements given 'orbital state vectors'. There's a stackexchange for this particular question which might help:
http://space.stackexchange.com/questions/1904/how-to-programmatically-calculate-orbital-elements-using-position-velocity-vectoThe other part you need is to be able to calculate the time evolution of the angle around the orbit. To do this you need to solve the equation for the 'true anomaly'. It seems like this is generally done by starting with the 'mean anomaly' M, which is the fraction of the orbital period that has passed since perigee (so sort of like an interpolating parameter). This is then related by a pair of transcendental equations to the true anomaly:
M = E - e sin(E) (the 'e' here is eccentricity; you want to solve this for E)
cos(E) = (e + cos(nu))/(1+ e cos(nu)) (the 'nu' here is the true anomaly, which is basically the angle around the orbit; you want to solve this for nu)
This page discusses it in more detail and gives a handy approximation, depending on how accurate a solution you need:
http://www.braeunig.us/space/orbmech.htm