eggrobin

Which is why I'm trying to find a way to avoid doing that. I'd like my inaccuracies to come from something I programmed (I could add radiation pressure if you want the Pioneer anomaly), rather than cancellations. At least, I'd like the errors from cancellations not to impede decent trajectory planning (especially since things can be quite unstable in the N-body problem).

I didn't do the math, Mattasmack did the simulation, and although his integrator was not symplectic, his energy drift was low enough and his results consistent enough with mine for me to conclude that everything except the Jool system is stable for at least a century. I'll try that and a few other things. I'm not sure what you mean by 'cancel each other out', but I'm not sure to what extent these gee forces actually exist either. They come from KSP (I'm just tracing perturbation and a few other things), so maybe bad calculations on KSP's side contribute to that. Bear in mind though that Unity uses single-precision floats, so with cancellations this can explain a lot here. Small enough accelerations can be dwarfed by cancellations. More testing is in order.

So will N-body physics Minmus is stable (the wobbling of the orbit as seen in-game is just the difference between an orbit around the barycenter and a an orbit around Kerbin). The big three moons of the Jool system are in a resonance, it's just that I don't use barycentric coordinates for the initial conditions. A mess ensues. I'm not sure if the system would be stable even if the resonance were properly set up though. I used to think that too, but here's a hint: a division and a square root are ~10 clock cycles each. A trigonometric function is ~100 clock cycles. In real life, the constants matter. Yes, I don't feel like explaining basic calculus. I can always point to Wikipedia for that. Anyway, unless you're really into topology, you don't want a math undergrad to teach you calculus (I'd probably go off a tangent onto algebraic topology anyway ). I'm rather surprised by the number of people in the 'I don't know what a derivative is' group. This is a good thing, KSP gets people into math. There are definitely some unacceptable fluctuations, just look at the peaks of .gee_immediate. You make some really good points here (drag is not an issue in space, but collisions and decoupling are); If I could grab the calls to AddForce and the like this would be easy, but yes, otherwise I'm basically talking about making my own physics engine here (which would be a good thing, but I'd have to be able to dig into SQUAD's code. I'm not sure C++ would be a good choice though, the .NET framework has its advantages. Maybe C#, interfaced with unmanaged bits for the performance-critical stuff? Of course if I could pick a language I'd say Ada, but that comes with its own set of problems.). It would probably be trivial to just insert calls to my methods in RigidBody, but an installation process that requires the user to do disassembly/reassembly is cumbersome, not to mention legal issues (not even with SQUAD, that's Unity code). Until I can come up with something smarter, I'll just experiment to see how far upstream I can grab the acceleration, and whether I can clamp it to kill numerical drift. Low accelerations (long-duration burns) will be specially handled anyway. What accelerations should be considered negligible?

The orbit on the second image is the same as the one on the first image. It's just a Kerbin-Synchronous Orbit (KSO). The weird shape comes from the choice of reference frame: as seen in a nonrotating Kerbin-centric reference frame (say, somebody standing on a pole of Kerbin rotating as the day goes by to always look at the same distant star) it's just an ellipse (the cyan line, as predicted by stock), but as seen in the rotating reference frame (by an observer on the equator of Kerbin, looking up), it's this sort of figure-eight (it's an analemma).

It's all very unstable, look at the drift of the apoapsis/periapsis when out of warp (including the data from .perturbation, which wobbles around 1e-5 when it should be 0, and is actually an average. .geeForce_immediate wobbles a lot more, with peaks above 1 g0 when it should be 0). This is probably all coming from Unity physics. The question is not really efficiency (this is completely negligible compared to the other things I have to do) but compatibility: for every module that applies an external force in space, I need to have code that grabs the force (in stock, I think this is ModuleRCS, ModuleEngines, ModuleEnginesFX; can anyone think of what I'm missing in this list?). Physics engines are not made to handle physics: no-one cares if your bullets are slightly too energetic or pick up angular momentum in a first-person shooter.

I suspect this might just be numerical instability on Unity's side. Conservation of angular momentum is on the roadmap for my mod (admittedly rather far away), and since I'll have to do most of the integration myself, my approach should fix the issue. In the meantime, perhaps you can try applying a few corrections in the FixedUpdates after coming out of warp to get the angular momentum to its expected value. Good luck!

To some extent, yes. I'll be managing the movement of the center of mass of every vessel in space (maybe in the atmosphere too, the ground looks hard however). It looks like this mod is not going to be too friendly on compatibility though (grabbing precise enough accelerations is hard, I'll have to get the forces from the engines/RCS/&c. directly), so compatibility with any module that exerts thrust will have to be implemented through a plugin.

Thanks for taking the time to do that! Well, if you were to use it to play, you would find that rockets are useless when using Newtonian physics: when the craft goes off-rails it uses stock physics, and when it goes back on-rails everything that happened off-rails is ignored. In order to play you'd have to 'Switch to Keplerian' for each burn... Don't worry, I'm working on this right now. It feels like I have to handle just about everything Unity is supposed to do... A problem with game physics engines is that they are terrible at physics. I'd prefer to keep the planets under N-body physics for several reasons. This is the first reason: it's the other way around; it would actually be slower to perform the integration with planets on Keplerian orbits (I'll explain later). The second reason is that I can probably tweak the initial conditions a bit to make things stable. If Scott Manley (illectro)'s intuitions are verified, interpreting the orbits as around the barycenter would fix things, otherwise I might have to tweak the masses a bit, but it should be doable (after all, it works for Jupiter). The third reason is that I would lose symplecticity (I'll explain later). The fourth reason is that if the bodies followed the Keplerian approximation (orbiting the primary as in stock rather than the barycenter), Lagrange points and other interesting features of the CR3BP/ER3BP wouldn't work as they should. I don't think I can help, I use VS2013. Feel free to fork and maintain a branch with your modified solution and project files for the people out there who use MonoDevelop. Building from source is only going to get harder when I get to optimisation, you'll have to install clang. Okay, it seems not everyone on this thread (and no-one at SQUAD) is really familiar with numerical integration. For the sake of everybody's sanity, I think I'll make a short series of posts roughly explaining Hamiltonian mechanics, RK and SPRK integrators, as well as the Kepler problem (which is more complicated than it looks) and the Saha & Tremaine integrator. Can I get a rough idea of where to start? In particular, can I assume that people know what a derivative is? Otherwise, can an engineer or physicist tell me how I should explain basic calculus? The best approach would be to start with point-set topology, but in the interest of brevity and not boring everyone to death, 10 posts about the properties of Hausdorff spaces would not be the best way to start. EDIT: Added a poll for those purposes.

PSA: the 'Simulator' project probably won't compile, it's leftover from my Alternis simulations. It's not needed. Alas, that's a scalar (it's the norm of total acceleration - acceleration from fictitious forces---I'm counting gravity as a fictitious force, even though I'm not implementing GR; it's just more convenient terminology). I think I found the answer though. From the OP (I'm not assuming people will read the entire thread, but I am assuming they'll read the OP), So if you want to look at the current prototype you can pull the source from GitHub and compile it; don't expect it to be playable soon (the current state is not an alpha; it's a messy way of generating nice screenshots). Once this is done, you too can generate nice screenshots! (I guess you can use it to get a feeling for what N-body trajectories look like, but it's highly likely that the current plotting is a bit off.) As far as providing the DLLs directly, configuring Git to do that without risking to accidentally publish a KSP assembly is more work than I want to do at the moment. Ask around if you don't know how to build from source.

Is gForce the gravitational force, or the nongravitational acceleration (the acceleration that is seen on the g-meter)? The name g-force is just terribly confusing, as it typically refers to the difference between the acceleration and the acceleration due to gravity... I think there's some Coriolis stuff in the vessel.flightintegrator, so that counts as reference frame weirdness. Also, what about torque? Anyway, it's 6 in the morning for some reason, so I'll get some sleep.

Cross-posting from the OP: Now that I have fixed trajectories lagging, the next step is going to be handling off-rails stuff. This is going to be interesting. Does anybody know if there is a way to find the acceleration on a vessel excluding gravity and reference frame weirdness?

That worked! The specifics of the solution: Have a class T inherit MonoBehaviour, call MapView.MapCamera.gameObject.AddComponent<T>() in your addon's Start() and do the redraw in T.OnPreCull().

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, 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).

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).

It still gets slower, it's just a question of linear vs. quadratic, see my reply to Dragon01. 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.

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: 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).

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).

I tried setting them in Update(), they still lag. 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.

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

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...

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). EDIT: 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.

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...

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.

The 'exams' thing will last until 2014-02-14 (only to come back in six months, and be replaced by 'studying' in the meantime). The other two reasons are here to stay. Why would I want to target .NET Framework assemblies below 4.0? I will only be using the C++/CLI part for easier interfacing with the unmanaged C++ part (probably compiled with clang due to clang's interpretation of long double). Anyway, I'm not downgrading this project to VS2008 (though I have it somewhere), bad things would be likely to happen. I'd rather do the marshaling myself if the C++/CLI scheme doesn't work out. Rewriting stuff in C++ is not my first priority at the moment though; a tremendous speedup will already come from using a better integrator. The expected order of development can be seen in the 'further modding' section (though neat visualisations, such as drawing gravitational potentials, are fun to write and make for great screenshots, so I might implement those sooner than I should). One of my first tasks will therefore be to find out a way to grab the acceleration of the vessel (not due to gravity) and plug that in my integrator.

That's what I was thinking about. Contours should be fairly trivial, I was wondering whether a mesh would be more readable though (it wouldn't be much harder, but some care would have to be taken in order to make the reference plane readable). Regarding the choice of plane, how about the Keplerian orbital plane? That makes things usable for inclined orbits, and it will pretty much do what you want in the CR3BP (I guess there will always be some edge cases, but if things work well in both the 2-body limit and the CR3BP I'll be happy for now). The centrifugal potential is arises from the choice of a rotating reference frame, it's not a property of the system. So that's pretty easy, I'll just add the centrifugal potential for whatever reference frame I'm drawing the orbits in. If I draw the potential, you'll get something more useful than the Hill sphere (and you'll see something that looks like the Hill sphere as one of the contours). This sounds like a very good idea, it would help with visualising complicated transfers. It would only double the performance cost of dynamic integration, so that should not be too much of a problem. This is sadly a fairly good assumption in general. I won't refrain from expanding the list of open questions as new ones come up.