2

I am making an n-body simulation in python. There are many different methods to numerically solve the system of differential equations governing the gravitational interactions between the n particles; the challenge of n-body simulations is to find an integration scheme which is both accurate and fast.

For the moment, I am working with n small (say n ≤ 5). I began with the most naive method of Euler integration, which is first order. Even with two bodies, this method can produce pretty bad errors if the bodies come very close (since then the force between them becomes very large, so their accelerations become very large). Of course, one can always reduce the time-step to get more accuracy, but at the expense of computation time.

Next I would like to implement a higher order method which has better stability properties etc like leapfrog integration or Runge-Kutta.

I expect these methods to give much better results, but the problem is I am unsure how to measure quantitatively the accuracy of my simulation. What metric should I use?

One of the issues is that I can't tell if the simulation is giving realistic results, because I haven't been using the real value of the gravitational constant G and none of the distances/times/masses in my simulation have units.

The reason why I haven't been paying attention to units is because, say if I wanted to do a scale sun-earth-moon simulation; since the sun-earth distance is about 370 times the earth-moon distance, it's not possible for everything to fit inside the display window and also be able to distinguish the earth and the moon.

The other problem is using the real value of G and real masses etc leads to computations involving big numbers. But, keeping everything else fixed, choosing lower values of G reduces the errors (because the forces are then smaller). So I have no real way of telling if the simulation is accurate.

So how can I rigorously test the accuracy of the simulation, and what can I do concerning units?

Arya McCarthy
  • 8,554
  • 4
  • 34
  • 56
math_lover
  • 886
  • 7
  • 20
  • Sorry but that doesn't help me. The problem I'm having is that I have no method of comparison. I don't know how to evaluate the accuracy of my programs, because I don't know which metric(s) to use. Moreover, I can't tell if the results I'm getting are realistic because I'm not sure which units to use. – math_lover May 12 '17 at 21:13
  • 1
    Another standard tool you could do, then, is use the converge plot. Take increasingly small step sizes and see if the result converges to a solution—whether the error between the 'true' and computed values is continuing to decrease. – Arya McCarthy May 12 '17 at 21:16
  • download our star system's orbital data for the last 10 years. compare it to your program. There are just a bunch of planets out there. Shouldn't be hard find necessary data from astronomics web sites. Or you can just start with moon and earth. – huseyin tugrul buyukisik May 13 '17 at 18:57
  • @huseyintugrulbuyukisik LoL just Moon and Earth is so wrong ... it is Sun-Earth-Moon and computation involved behind it is not so easy even approximations are a nightmare. Also getting real data on our solar system is not that easy as it should be ... there are more reference frames and most data on the web is not accurate enough. – Spektre May 14 '17 at 07:06
  • @JoshuaBenabou I would compare your results with iterative numerical integration like this: [Is it possible to make realistic n-body solar system simulation in matter of size and mass?](http://stackoverflow.com/a/28020934/2521214) and use small enough time steps and the same constants as for your own integrations. But you need to use the same scale for all I am using SI so to see how to use bigger floating values see the last part of that QA or use higher precision variables. Of coarse reality and simulation will be off for non rigid bodies like Earth ... or near high mass (Mercury) – Spektre May 14 '17 at 07:09
  • 1
    Excuse my ignorance but would you not just monitor the total energy in the system. As you are simulating a closed system any problems with the calculations would manifest as a lose or gain in angular momentum, it may be a progressive accumulation or a cyclic variation, both of which would be a sure sign something is wrong. – Blindman67 May 14 '17 at 07:55
  • @Blindman67 : Indeed measuring the energy drift is one metric of accuracy, but it's not the only one. You can have integrators which more or less conserve energy but still give wrong results. – math_lover May 14 '17 at 08:17
  • Mmm precision v accuracy As i understand it there is no way known to get the correct result, nor can you know if a given state reflects flawed calculations or a flawed initial state. Another option would be the use of different models relying on the divergence between them to quantify an error. Unless you use Einstein's field equations you are just doing an approximation anyway, so accuracy is lost, so why not look at the precision of your integration in comparison to others. I sometime use a little noise in the initial state and see how the integration amplifies it to know if I am on track. – Blindman67 May 14 '17 at 09:29

0 Answers0