[WIP][1.8.1, 1.9.1, 1.10.1, 1.11.0–2, 1.12.2–5] Principia—version ‎‎Kronecker, released 2024-11-01—n-Body and Extended Body Gravitation

[eggrobin]

That might not work for performance reasons. The best I've managed to get 100k steps of RK5(4) executed is in about 1 second if my memory serves me - that it's on powerful i7 3930K CPU with significant portions of code written manually on SSE2/3 intrinsics (=assembly).

And in your case it will have to be done in addition to regular updates calculations. You could probably run it for t and (t+dt) in parallel on separate cores if you write the code carefully not to use any shared state, but experience with my mod (which is heavily parallel and extensively uses all available cores) indicates that there are enough users here with relatively weak 2-core or even 1-core CPUs.

To clarify:

The main integration currently uses a pretty dumb 5th order SPRK (it has 6 RHS evaluations per step); as it is written in C# it's about an order of magnitude slower than yours (do you use Dormand-Prince, Cash-Karp, Fehlberg, or something else? RK5(4) is a rather ambiguous name. Also, what is the dimension of your system? Most of the time is spent computing the RHS anyway, so your mileage may vary quite a bit depending on the specifics of the problem). I intend to switch to an integrator designed for N-body gravitation, which should allow me to run 1_000_000x timewarp easily.

The trajectories ferram4 was talking about are those from the dynamic integration (I haven't found a better name, see the 'terminology' section). Those will only be computed while the active vessel is under thrust (not counting the ion engine stuff, which I'll have to handle differently), so out of timewarp. When out of timewarp, the main integration has negligible performance cost (recall that this is the same integration that we have to be able to run at 1_000_000x; even with my current SPRK, this is virtually free outside timewarp).

Rough estimates for the performance cost of the dynamic integration over 10 years yield something on the order of 50 ms: use 10 days as a timestep, which may sound crazy, but look at the timesteps in (Chambers 1999), assume only 10_000 steps per second, which is probably pessimistic. Even if I have to double that to show two trajectories, I can just update the predictions every couple of physics update, or make them shorter, or, as ferram4 said, make one with extrapolated initial conditions less precise.

[eggrobin]

Well, its an official statement from MS regarding cli/c++ assemblies on vs2010+ that target <net4, so you might want to check out your assembly in ILSpy or sth similar. You will most likely find a "Runtime: .NET 4.0" and reference of mscorlib and System.dll to both version 2.0 and 4.0, like this:

http://i.imgur.com/NlnyuD7.png?1

Ofc that might not be a problem and Unity-Mono simply ignores those references to non-existent librarys... until you accidentally use any feature that requires .NET 4.0 and you just don't get whats wrong, especially since it compiles just fine. Thats the same reason we frequently tell newer modders to target a framework lower than 4.0 that is more compatible with Ksp's Unity-Mono.

I found about about those CLI/C++ limitation when i tried to use my library in a C# project that properly targeted a pure .NET 2 runtime. But hey, if it works for you, no one will force you to change anything. And in case you experience "crazy" issues you now know an additional possible source.

KSP's Assembly-CSharp already targets the .NET framework version 4.5 (and so do I), so I'm not sure why this should be a problem...

[thorfinn]

The thing is, once you get used to how they work, 2-body orbital mechanics are pretty simple. This will hopefully provide an interesting challenge and new possibilities, while not fundamentally changing most of the game.

I guess we'll know how hard it becomes when we get there (I'm still a bit scared, but after all satellites moving around would be a really big problem only for Remotetech, true.)

But if you ever got to this:

Gravitational moment (quadrupole) for planets, enabling such amusing things as Sun-synchronous orbits

I would be very curious to see what orbits would look like around a realistic Haumea

[asmi]

The main integration currently uses a pretty dumb 5th order SPRK (it has 6 RHS evaluations per step); as it is written in C# it's about an order of magnitude slower than yours (do you use Dormand-Prince, Cash-Karp, Fehlberg, or something else? RK5(4) is a rather ambiguous name. Also, what is the dimension of your system? Most of the time is spent computing the RHS anyway, so your mileage may vary quite a bit depending on the specifics of the problem). I intend to switch to an integrator designed for N-body gravitation, which should allow me to run 1_000_000x timewarp easily.

I used code for algorithm from that post as baseline: http://forum.kerbalspaceprogram.com/threads/49121-On-Newtonian-trajectories-vs-patched-conics and that post for more detailed description: http://www.roesle.org/cms25/index.php/projects/81-general/95-on-newtonian-trajectories-in-kerbal-space-program The code uses Dormand and Prince's Runge Kutta RK5(4)7M. It uses adaptive time step and also provides coefficients for fifth order interpolation. The idea was to precalculate trajectory to set amount of steps ahead, and use these points and interpolation to actually move spacecraft. This way we would only need to recalculate the whole trajectory is there is some sort of dynamic event (like applying the thrust, or entering atmosphere). It provided pretty good accuracy, and would be relatively cheap since you could calculate additional data points in the background. Unfortunately due to spare time constraints I never got to the point of getting it into KSP, but was focusing on profiling algorithm's performance and used simple WinForms app for ease of testing. If you're interested I could try to dig out that code and send it to you. Oh, and I used whole KSP celestials (orbital and physical data is taken from wiki).

The trajectories ferram4 was talking about are those from the dynamic integration (I haven't found a better name, see the 'terminology' section). Those will only be computed while the active vessel is under thrust (not counting the ion engine stuff, which I'll have to handle differently), so out of timewarp. When out of timewarp, the main integration has negligible performance cost (recall that this is the same integration that we have to be able to run at 1_000_000x; even with my current SPRK, this is virtually free outside timewarp).

Rough estimates for the performance cost of the dynamic integration over 10 years yield something on the order of 50 ms: use 10 days as a timestep, which may sound crazy, but look at the timesteps in (Chambers 1999), assume only 10_000 steps per second, which is probably pessimistic. Even if I have to double that to show two trajectories, I can just update the predictions every couple of physics update, or make them shorter, or, as ferram4 said, make one with extrapolated initial conditions less precise.

Time step have to be adaptive, as sometimes even 10s is too much, while in other cases 10 days is OK. And why ion engines are so special? There are plenty of other low-thrust engines, besides it's not the thrust is what's matters, but acceleration they can provide.

[thorfinn]

And why ion engines are so special?

For very low thrust engines ("ion" is what we have in stock now in that class) you'd want to be able to apply thrust during high timewarp, so craft physics will be off and you'd have to have some replacement simplified model to apply thrust to, decrease mass et cetera. I guess you would have to handwave thrust offsets and things like that.

I suppose there would be an acceleration cutoff somewhere low (0.01G?) below which your ship is considered "ion class" and can thrust during high warp.

[eggrobin]

Time step have to be adaptive

I used to think that too. It's more complicated. It's always more complicated.

In this case we care about symplecticity, and adaptive timesteps kill symplecticity. See the quote from Chambers in the OP. You can however have different timesteps for different bodies, as long as each timestep divides the next one (Saha and Tremaine 1994). Chambers has a form of adaptive timestep (it's more of a dual-method thing), but it's far more complicated than the usual business either looking at convergence or using an embedded method like Dormand-Prince. All these 'smarter' integrators use the Keplerian approximation to some extent to allow for longer timesteps, so it's not your average RK or SPRK method.

And why ion engines are so special? There are plenty of other low-thrust engines, besides it's not the thrust is what's matters, but acceleration they can provide.

Ion engines are not special. I was referring to low-acceleration long duration planned burns, in contrast with the real-time burns we currently have. That's a pretty bad misnomer though, I guess I'll call these long-duration burns from now on (unless you can find a better name).

EDIT:

I used code for algorithm from that post as baseline: http://forum.kerbalspaceprogram.com/threads/49121-On-Newtonian-trajectories-vs-patched-conics and that post for more detailed description: http://www.roesle.org/cms25/index.php/projects/81-general/95-on-newtonian-trajectories-in-kerbal-space-program The code uses Dormand and Prince's Runge Kutta RK5(4)7M. It uses adaptive time step and also provides coefficients for fifth order interpolation.

Oh, that's Matt Roesle (Mattasmack)! He's done some interesting work on N-body simulations in KSP (see the acknowledgements), but his choice of integrator is suboptimal. A symplectic one is much better (and one specifically developed for the gravitational N-body problem would be leaps and bounds above that).

EDIT2: Actually, the fact that he doesn't use compensated summation is more critical than the lack of symplecticity. This would be easy to add though.

Edited February 6, 2014 by eggrobin
Close parenthesis.

[BubbaWilkins]

I used to think that too. It's more complicated. It's always more complicated.

In this case we care about symplecticity, and adaptive timesteps kill symplecticity. See the quote from Chambers in the OP. You can however have different timesteps for different bodies, as long as each timestep divides the next one (Saha and Tremaine 1994). Chambers has a form of adaptive timestep (it's more of a dual-method thing), but it's far more complicated than the usual business either looking at convergence or using an embedded method like Dormand-Prince. All these 'smarter' integrators use the Keplerian approximation to some extent to allow for longer timesteps, so it's not your average RK or SPRK method.

Ion engines are not special. I was referring to low-acceleration long duration planned burns, in contrast with the real-time burns we currently have. That's a pretty bad misnomer though, I guess I'll call these long-duration burns from now on (unless you can find a better name.

I expected you to trot out Cappa Fero, but seeing as Thibault cancels that out unless you've studied Itsa Glippa, we'll stick with Bonetti's defense.

I suppose you aren't left handed either.

In all seriousness, this is a very ambitious project. I look forward to see it's development progress.

Edited February 6, 2014 by BubbaWilkins

[Bluekommeh]

"We cant do N-Bodies, they would fry your computer :^)" -Squad, 2011

Mods will blow Squead into bankruptcy.

[eggrobin]

Mods will blow Squead into bankruptcy.

I hope not, there wouldn't be anything left to mod...

Mods only make the game more attractive if anything, so I think we're safe though.

Squad's attitude on this particular matter is understandable, as this is really hard to do right, and the Keplerian approximation is mostly pretty good. Stock aerodynamics are harder to justify...

Edited February 6, 2014 by eggrobin

[eggrobin]

Itsa Glippa

From the fencing context, I'm pretty sure this should read 'Camillo Agrippa'.

Cappa Fero

Capo Ferro?

[lincourtl]

"We cant do N-Bodies, they would fry your computer :^)" -Squad, 2011

Mods will blow Squead into bankruptcy.

To the contrary, mods are great for showing Squad how to do things they might not thought possible. All they have to do is to continue supporting a robust modding community, assimilate the mods that add the most to the game, and they're fine. I can understand Squad's objections to n-body orbital mechanics on a pure game-play basis though. I don't think most people would find that fun. For those of us who do want that added realism, there's always the mods.

But yeah, the aerodynamics issue is much less defensible. Similarly, I'd like to see realistic reentry heating, axial tilt, magnetary fields, and radiation belts in the stock game. I'm not going to cry into my Cheerios if that never happens though.

[BubbaWilkins]

From the fencing context, I'm pretty sure this should read 'Camillo Agrippa'.

Capo Ferro?

Correct on both accounts, however I was actually referencing the Princess Bride duel scene. Apparently the site I used to look it up had most of it spelled wrong.

[eggrobin]

Are you setting it in Update or FixedUpdate? Make sure you use the former, I've seen effects lagging behind when using FixedUpdate. If that's not it, try seeing if you can do things in localSpace.

I tried setting them in Update(), they still lag.

axial tilt

This would be really nice. I'm not sure I can find a way to mod it in, and if I cannot, conservation of angular momentum for planets is lost.

I'll add it to the open questions.

Edited February 6, 2014 by eggrobin

[jxi24]

I have a clarification question. When for example calculating the N-body dynamics how are you handling the effects of for example the moons around Jool, effecting the orbits of other bodies in the system?

I know that the way large scale simulations for the galaxy are done at the research level is to subdivide the universe into regions, and in places where there is a lot of matter the grids are smaller and done with more accuracy.

For example, to simplify the number of calculations you need to preform, are you going to take the sum of Jool and it's moons and treat it as a point source to calculate the effects on the orbits of things far away from Jool, and then near Jool do more detailed calculations for the moons orbits?

[eggrobin]

I have a clarification question. When for example calculating the N-body dynamics how are you handling the effects of for example the moons around Jool, effecting the orbits of other bodies in the system?

I know that the way large scale simulations for the galaxy are done at the research level is to subdivide the universe into regions, and in places where there is a lot of matter the grids are smaller and done with more accuracy.

For example, to simplify the number of calculations you need to preform, are you going to take the sum of Jool and it's moons and treat it as a point source to calculate the effects on the orbits of things far away from Jool, and then near Jool do more detailed calculations for the moons orbits?

No, for the moment this is truly N-body. The answer is a bit more complicated once I start using integrators specifically designed for this problem (see Saha & Tremaine 1994), but mostly it will remain truly N-body (it's just that the Sun's influence will get special treatment).

Bear in mind the number of massive bodies is very small (17 in stock, conceivably ~50 if somebody decides to add lots of planets and moons; if I need to manage stuff like asteroids or rings I might want to try to be smarter, but the Saha & Tremaine integrator should be pretty good even for that).

Edited February 6, 2014 by eggrobin

[Guest]

How will the amount of non-massive bodies impact performance? People have launched quite a lot ships in a single save, I imagine this could easily hit 100 for some.

[asmi]

How will the amount of non-massive bodies impact performance? People have launched quite a lot ships in a single save, I imagine this could easily hit 100 for some.

If it's going to be true N-body, than the more stuff is in space, the slower performance would be. In practice it's quite common to cut off bodies that will have negligible impact on trajectory. But that would allow you some very cool stuff like putting a satellite around your huge space station

[eggrobin]

How will the amount of non-massive bodies impact performance? People have launched quite a lot ships in a single save, I imagine this could easily hit 100 for some.

The main cost is the RHS computation. If the B massive bodies exert forces on the B + V bodies, the cost of the RHS computation is proportional to B(B + V), so asymptotically linear in V and quadratic in B.

From the OP:

Near-unit-roundoff (using IEEE 754 binary64) errors can be achieved while sustaining 100_000x timewarp with B = 17, V < 5. Handwavy asymptotic calculations indicate that for V = 100, this integrator would slow down to about 20_000x timewarp. Using better integrators (and saner tolerances), performance could easily be increased by a couple of orders of magnitude, making the simulation usable for RSS as well as stock.

I'll only know once I have implemented it, but I think even 10_000_000x with B = 20, V = 100 should easily be feasible (see Chambers 1999 for the sort of timesteps I'll be able to use with better integrators and my answer to asmi on post #51 of this thread for very rough estimates of the speed for a single vessel).

[Spartwo]

I thought Unity couldn't support N-body physics.My mod list is getting longer and longer.

[eggrobin]

If it's going to be true N-body, than the more stuff is in space, the slower performance would be.

It still gets slower, it's just a question of linear vs. quadratic, see my reply to Dragon01.

In practice it's quite common to cut off bodies that will have negligible impact on trajectory. But that would allow you some very cool stuff like putting a satellite around your huge space station

Ok, I forgot to mention that. The V I define in the OP is not supposed to be 0, it stands for 'vessel'. Maybe I'll add an option for massive vessels, but my current integrator definitely can't handle that. Maybe the Saha & Tremaine one can.

[eggrobin]

I thought Unity couldn't support N-body physics.My mod list is getting longer and longer.

You might be confused.

Unity doesn't support 2-body orbital mechanics, nor does it support large distances or velocities (for some mindbogglingly stupid reason it uses single-precision floating point), rotating ground, &c.

KSP handles the orbit propagation, the floating origin stuff, the rotating reference frames &c. so that Unity can do the small amount of physics it can handle---such as (badly) managing collisions, joints, and various other forces.

What this mod will do is handle orbit propagation in a different way (using an integrator, rather than using the 2-body approximation).

[thorfinn]

for some mindbogglingly stupid reason it uses single-precision floating point

Well, almost nobody has the requirements of KSP, and it's a budget game engine at that.

[eggrobin]

Well, almost nobody has the requirements of KSP, and it's a budget game engine at that.

Nevertheless, single precision floating point is almost always a bad idea. On x87 FPUs it's slower, and unless your code happens to be highly vectorisable the performance improvement with SSE2/3 is negligible (it's not going to be more than twice as fast anyway, even if everything can be vectorised). Memory usage concerns from physics are not an issue when you have stuff like textures to worry about.

Anyway,

The method I currently use for drawing trajectories (LineRenderers, object in layer 9, transform.parent = null, useWorldSpace=true, updated every time the GUI is drawn) has some issues: when rotating the camera in map mode, the lines lag behind. Does this have a known solution? (I guess so, RT2 doesn't have this problem).

I'm mostly clueless as far as this is concerned. So far I have tried putting the object in layer 31 as Cilph does, drawing the trajectories at the same time as the GUI, in Update, in OnPreCull (again from the RT2 code). Nothing changes. I'll try making a dedicated renderer class and adding it as a component of MapView.MapCamera.gameObject like Cilph, since his code works (and he seems to have a clue about how Unity works).

Edited February 7, 2014 by eggrobin

[J.Random]

If this thing works, if it works fast enough, and/or if you manage to spread the processing load on multicore systems, I will want to pay for such mod. Because that's going to be an incredible, astonishing improvement for KSP, and it deserves the support.

[FanaticalFighter]

Any idea when you're releasing a prototype/alpha test? Because if you want testers, I'll love to do it!