Help designing an n-body gravity simulator in Python

In summary, the author wrote a gravity simulator in Python and found that it is unstable over long time scales and often results in bodies being propelled outward in straight lines. He has read up on Runge-Kutta integration and other methods of calculating orbits, but does not know how to implement them himself. He is looking for help from others who may be more experienced.
  • #1
Wer900
2
0
Hi everyone:

I wanted to design a 2-dimensional gravity simulator as part of a larger research project. I wrote it in Python using numpy to store positions, velocities, and accelerations, and matplotlib to display the trajectories; code is available upon request.

In order to solve the problem I've used Euler's method, which I know is numerically unstable over long timescales. In particular, I've been observing that over time many bodies are propelled away, in straight lines, from their original positions. By imposing an artificial minimum limit of 0.3 meters on effective gravitational radius (the masses and gravitational constant I'm using for testing purposes are very different from what I'll actually use) I'm able to reduce this mostly, although this reduces the accuracy of predictions at small radii.

I've read up on Runge-Kutta integration and other, more exotic methods of calculating orbits, but I have no idea of how to implement them in an orbital simulator. Would anyone here be able to help? Thanks in advance to everyone.
 
Technology news on Phys.org
  • #2
If you want a fairly easy improvement over explicit Euler, you may want to try the Leapfrog method [1], which for small enough time steps do not add false energy to the system you are simulating.

[1] https://en.wikipedia.org/wiki/Leapfrog_integration
 
  • #3
How do you know that the bodies that are "propelled away in straight lines" are incorrect? Close two body interactions should result in some bodies gaining energy and being ejected. This happens in real clusters and results in the clusters losing energy and becoming more tightly bound over time.
 
  • #4
Wer900 said:
I've read up on Runge-Kutta integration and other, more exotic methods of calculating orbits, but I have no idea of how to implement them

The magic phrase to Google for is "runge kutta for system of differential equations". You'll find pages like this one:

http://www.physics.buffalo.edu/phy410-505-2009/topic3/lec-3-2.pdf

Start with Newton's Second Law for each pair of objects. Convert the second-order differential equations into pairs of first-order differential equations. Put the components of all positions and velocities into a single array or vector. The intermediate values k1, k2, etc. in the Runge-Kutta algorithm are now also vectors of the same dimension. You basically re-interpret the scalar R-K equations as vector equations. Write functions that carry out the steps of the calculation, taking vectors as input and returning vectors as output.

I did something like this about 25 years ago, for a solar-system simulator, but it was in Fortran, not Python, so I can't help with the actual code. Instead of standard 4th-order Runge-Kutta, I used a variation called Runge-Kutta-Fehlberg which varies the step-size dynamically.

I had the problem that orbits tended to drift outwards because the R-K algorithm tends to not conserve energy over long time spans. Loosely speaking, the errors introduced by the discrete steps tend to reinforce each other instead of cancelling. I understand there are algorithms in which the errors do tend to cancel and energy is at least approximately conserved, but I didn't go any further with my simulator.

[added] Ah, now I see that Filip also mentioned this "false energy" problem. Maybe a good start would be to vectorize the leapfrog algorithm that he mentioned, and look for higher-order versions if necessary.
 
  • #5
Filip Larsen said:
If you want a fairly easy improvement over explicit Euler, you may want to try the Leapfrog method [1], which for small enough time steps do not add false energy to the system you are simulating.

[1] https://en.wikipedia.org/wiki/Leapfrog_integration

Thanks for the link. I think I tried this before, but maybe reversed the order in which I set velocity and acceleration (causing eventual drifting outward of the system).

phyzguy said:
How do you know that the bodies that are "propelled away in straight lines" are incorrect? Close two body interactions should result in some bodies gaining energy and being ejected. This happens in real clusters and results in the clusters losing energy and becoming more tightly bound over time.

Yes, I understand, but I sometimes saw bodies (representing the majority of the mass of the overall system) propelled outward even in 2 or 3 body systems. This doesn't make a whole lot of sense.

Also: Thanks jtbell for the info. I'm reading up on it right now, hopefully I can get something working soon.
 
Last edited:

Related to Help designing an n-body gravity simulator in Python

1. How can I create a simulation of n-body gravity in Python?

To create a simulation of n-body gravity in Python, you will need to use the laws of gravity and motion to calculate the forces and accelerations between each body. You can then use this information to update the positions and velocities of each body over time.

2. What libraries or modules are available for creating a gravity simulator in Python?

There are several libraries and modules available for creating a gravity simulator in Python, such as PyGame, PyOpenGL, and PyODE. You can also use the built-in math and physics modules to perform the necessary calculations.

3. How can I add realistic features to my gravity simulator, such as collisions and orbital decay?

To add realistic features to your gravity simulator, you will need to incorporate more advanced algorithms and methods, such as the Barnes-Hut algorithm for efficient calculation of forces between bodies and the Runge-Kutta method for more accurate integration of equations of motion. You can also add in factors such as drag and friction to simulate collisions and orbital decay.

4. Can I customize the visual representation of my gravity simulator?

Yes, you can customize the visual representation of your gravity simulator by using graphics libraries such as PyGame or PyOpenGL. You can also import custom images or 3D models to represent the bodies in your simulation.

5. Is there a way to speed up the calculations in my gravity simulator?

Yes, there are several ways to speed up the calculations in your gravity simulator. One way is to use multi-threading or parallel processing to distribute the calculations across multiple CPU cores. You can also simplify your calculations by using approximations or simplifying assumptions, although this may reduce the accuracy of your simulation.

Similar threads

  • Programming and Computer Science
Replies
19
Views
2K
  • Programming and Computer Science
Replies
4
Views
1K
  • Programming and Computer Science
Replies
1
Views
2K
  • Programming and Computer Science
Replies
15
Views
2K
  • Programming and Computer Science
Replies
3
Views
4K
  • Programming and Computer Science
Replies
6
Views
1K
  • Classical Physics
Replies
2
Views
1K
  • Programming and Computer Science
Replies
7
Views
2K
  • Classical Physics
Replies
7
Views
914
Replies
12
Views
2K
Back
Top