Optimizing gravity calculations

Question

I've got a bunch of objects of varying size and velocity which gravitate towards each other. On every update, I have to go over every object and add up the forces due to gravity of every other object. It doesn't scale very well, is one of two big bottlenecks I've found in my game, and I'm not sure what to do to improve the performance.

It feels like I should be able to improve the performance. At any given time, probably 99% of the objects in the system will only have a negligible influence on an object. I of course can't sort the objects by mass and only consider the top 10 largest objects or something, because the force varies with the distance more than with the mass (the equation is along the lines of force = mass1 * mass2 / distance^2). I think a good approximation would be to consider the largest objects and the closest objects, ignoring the hundreds of tiny fragments of rock on the other side of the world that can't possibly affect anything -- but in order to find out which objects are closest I have to iterate over all the objects, and their positions are changing constantly, so it's not like I can just do it once.

Currently I'm doing something like this:

private void UpdateBodies(List<GravitatingObject> bodies, GameTime gameTime)
{
    for (int i = 0; i < bodies.Count; i++)
    {
        bodies[i].Update(i);
    }
}

//...

public virtual void Update(int systemIndex)
{
    for (int i = systemIndex + 1; i < system.MassiveBodies.Count; i++)
    {
        GravitatingObject body = system.MassiveBodies[i];

        Vector2 force = Gravity.ForceUnderGravity(body, this);
        ForceOfGravity += force;
        body.ForceOfGravity += -force;
    }

    Vector2 acceleration = Motion.Acceleration(ForceOfGravity, Mass);
    ForceOfGravity = Vector2.Zero;

    Velocity += Motion.Velocity(acceleration, elapsedTime);
    Position += Motion.Position(Velocity, elapsedTime);
}

(note that I've removed a lot of code -- for example the collision tests, I do not iterate over the objects a second time to detect collisions).

So I'm not always iterating over the whole list -- I only do that for the first object, and every time the object finds the force it feels toward another object, that other object feels the same force, so it just updates both of them -- and then that first object doesn't have to be considered again for the rest of the update.

The Gravity.ForceUnderGravity(...) and Motion.Velocity(...), etc. functions just use a bit of XNA's built in vector math.

When two objects collide, they create massless debris. It's kept in a separate list and the massive objects do not iterate over the debris as part of their velocity calculation, but each piece of debris must iterate over the massive particles.

This doesn't have to scale to incredible limits. The world isn't unlimited, it contains a border which destroys objects that cross it -- I'd like to be able to handle maybe a thousand or so objects, currently the game starts to choke around 200.

Any thoughts on how I might improve this? Some heuristic I can use to shave the loop length from hundreds down to just a few? Some code I can execute less often than every update? Should I just multithread it until it's fast enough to allow for a decent sized world? Should I try to offload the velocity calculations to the GPU? If so, how would I architect that? Can I keep static, shared data on the GPU? Can I create HLSL functions on the GPU and call them arbitrarily (using XNA) or do they have to be part of the draw process?

Answer 1

This sounds like a job for a grid. Divide your game space into a grid and for each grid cell keep a list of the objects currently in it. When objects move across a cell boundary, update which list they're in. When updating an object and searching for others to interact with, you can look at just the current grid cell and a few neighboring ones. You can tweak the size of the grid for the best performance (balancing the cost of updating the grid cells - which is higher when the grid cells are too small - with the cost of doing the searches, which is higher when the grid cells are too large).

This will, of course, cause objects that are farther apart than a couple of grid cells not to interact at all, which is probably an issue because a large accumulation of mass (either a big object, or a cluster of many small objects) should, as you mentioned, have a larger region of influence.

One thing you could do is keep track of the total mass within each grid cell, and treat the whole cell as a single object for the purposes of farther-away interactions. That is: when you calculate the force on an object, calculate the direct object-to-object acceleration for the objects in a few nearby grid cells, then add in a cell-to-cell acceleration for each farther-away grid cell (or maybe just the ones with a non-negligible amount of mass in them). By a cell-to-cell acceleration, I mean a vector computed using the two cells' total masses and the distance between their centers. That should give a reasonable approximation of the summed gravity from all the objects in that grid cell, but much more cheaply.

If the game world is very large, you could even use a hierarchical grid, like a quadtree (2D) or octree (3D), and apply similiar principles. Longer-distance interactions would correspond to higher levels of the hierarchy.

Answer 2

The Barnes-Hut algorithm is the way to go with this one. It's been used in supercomputer simulations to solve your exact problem. It's not too hard to code, and it's very efficient. I actually wrote a Java applet not too long ago to solve this problem.

Visit http://mathandcode.com/programs/javagrav/ and press "start" and "show quadtree".

In the options tab, you can see that the particle count can go all the way up to 200,000. On my computer, the calculation finishes in about 2 seconds (the drawing of 200,000 dots takes about 1 second, but the calculation runs on a separate thread).

Here's how my applet works:

1. Create a list of random particles, with random masses, positions, and starting velocities.
2. Construct a quadtree from these particles. Each quadtree node contains the center of mass of each subnode. Basically, for each node you have three values: massx, massy, and mass. Every time you add a particle to a given node, you increase massx and massy by particle.x*particle.mass and particle.y*particle.mass respectively. The position (massx/mass, massy/mass) will end up as the node's center of mass.
3. For each particle, calculate the forces (fully described here). This is done by starting at the top node and recursing through each subnode of the quadtree until the given subnode is small enough. Once you stop recursing, you can calculate the distance from the particle to the node's center of mass, and then you can calculate the force using the node's mass and the particle's mass.

Your game should easily be able to handle a thousand mutually attracting objects. If each object is "dumb" (like the bare-bones particles in my applet), you should be able to get 8000 to 9000 particles, maybe more. And this is assuming single-threading. With multi-threaded or parallel computing applications, you can get many more particles than that updating in real time.

See also: http://www.youtube.com/watch?v=XAlzniN6L94 for a large rendering of this

Answer 3

Nathan Reed has an excellent answer. The short version of it is to use a broadphase technique that fits your simulation's topology, and to only run the gravity calculations on pairs of objects that are going to have a noticeable effect on each other. It's really no different than what you'd do for regular collision detection broadphase.

Carrying on from that, though, another possibility is to only update objects intermittently. Basically, each time step (frame) only update a fraction of all the objects, and leave the velocity (or acceleration, depending on your preference) the same for the other objects. The user is unlikely to notice any delay in the updates from this so long as the intervals aren't too long. This will give you a linear speed up of the algorithm, so definitely look into the broadphase techniques like Nathan suggested as well, which can give much more significant speedups if you have a ton of objects. While not at all modeled the same in the slightest, this is kind of like having "gravity waves." :)

Also also, you could generate a gravity field in one pass, then update the objects in a second pass. The first pass you're basically filling a grid (or a more complex spatial data structure) with the gravity influences of each object. The result now is a gravity field that you can even render (looks pretty cool) to see what acceleration will be applied to an object at any given location. Then you iterate over the objects and simply apply the gravity field's effects to that object. Even cooler, you can do this on a GPU by rendering the objects as circles/spheres to a texture, then reading the texture (or using another transform-feedback pass on the GPU) to modify the object's velocities. So long as the GPU has a decent fillrate and your objects have a reasonable (if technically unrealistic) cut-off distance for their gravitational effect, this will likely be the fastest approach (but profile and test, as always).

Answer 4

I would recommend using a Quad Tree. They allow you quickly and efficiently look up all objects in an arbitrary rectangular area. Here's the wiki article on them: http://en.wikipedia.org/wiki/Quadtree

And a shameless link to my own XNA Quad Tree project on SourceForge: http://sourceforge.net/projects/quadtree/

I would also maintain a list of all the large objects so that they can interact with everything regardless of distance.

Answer 5

Just a small bit of (possibly naive) input. I don't do game programming, but what I am feeling is that your fundamental bottleneck is the gravity-due-to-gravity calculation. Instead of iterating over each object X and then finding the gravitational effect from each object Y and adding it, you can take each pair X,Y and find the force between them. That should cut the number of gravity calculations from O(n^2). Then you will be doing a lot of addition (O(n^2)), but this is normally less expensive.

Also at this point you can implement rules such as "if the gravitational force will be less than \epsilon because these bodies are too small, set the force to zero". It may be advantageous to have this structure for other purposes too (including collision detection).

Answer 6

In extending seanmiddleditch's answer, I thought I might shed some light (irony?) on the gravity field idea.

Firstly, don't think of it as a texture, but a discrete field of values that can modified (a two-dimensional array, as it were); and the subsequent accuracy of the simulation could be the resolution of that field.

When you introduce an object into the field, its gravitation potential can be calculated for all surrounding values; thereby creating a gravitation sink in the field.

But how many of these points should you calculate before it becomes more or as ineffective as before? Probably not many, even 32x32 is a substantial field to iterate for each object. Therefore break the entire process into multiple passes; each with varying resolutions (or accuracy).

Ie, the first pass may calculate the objects gravity represented in a 4x4 grid, with each cell value representing a 2D coordinate in space. Giving an O(n * 4 * 4) sub-total complexity.

The second pass may more accurate, with a 64x64 resolution gravity field, with each cell value representing a 2D coordinate in space. However, as the complexity is very high, you can restrict the radius of surrounding cells affected (perhaps, only the surrounding 5x5 cells are updated).

An additional third pass could be used for high accuracy calculations, with maybe a resolution of 1024x1024. Remembering at no time are you actually performing 1024x1024 separate calculations, but operating only on portions of this field (perhaps 6x6 sub-sections).

In that way, your overall complexity for the update is O(n * (4*4 + 5*5 + 6*6)).

To then calculate the velocity changes to each of your objects, for each gravity field (4x4, 64x64, 1024x1024) you just map the point masses position to a grid cell, apply that grid cells overall gravitational potential vector to a new vector; repeat for each "layer" or "pass"; then add them together. This should give you a good resultant gravitational force vector.

Therefore, the overall complexity is: O(n * (4*4 + 5*5 + 6*6) + n). What really counts (for complexity) is how many surrounding cells you update when calculating the gravitational potential in the passes, not the overall resolution of the gravity fields.

The reason for low resolution fields (first passes) is to obviously encompass the universe as a whole, and ensure outlying masses are attracted to more dense areas despite distance. Then use higher resolution fields as separate layers to increase accuracy for neighboring planets.

I hope this made sense.

Answer 7

How about another approach:

Assign an area of influence to objects based on their mass--they're simply too small to have a measurable effect beyond that range.

Now divide your world up into a grid and put each object in a list of all the cells it has an influence on.

Do your gravity calculations only on objects in the list attached to the cell an object is in.

You only have to update the lists when an object moves into a new grid cell.

The smaller the grid cells the less calculation you'll do per update but the more work you'll do updating lists.