Fastest way to calculate gravity-like forces experienced between large set of particles?

Question

Given a set of n particles in 2-dimensional space, where each particle experiences an attractive force towards all other particles in the set, proportional to the inverse square of the distances between them, can you imagine any clever methods of calculating the net force experienced by each particle, without performing n^2 distance calculations?

From the Particle class of my program:

    public void calculateAttractions(List<Particle> entities)
    {
        foreach (Particle p in entities)
        {
            if (!p.Equals(this))
            {
                Vector2 dx = this.Position - p.Position;
                this.ApplyForce(-dx / dx.LengthSquared());
            }
        }
    }

The brute-force method can only simulate ~150 particles in my program (using the XNA framework and Farseer physics engine) before running into performance issues.

If there is no other way, how might this code be optimized, or else how might I cheat? I've tried randomizing the call to calculateAttractions() on a particle-to-particle basis, such that at any given time only around 1/3 of the particles are actually being updated. This resulted in some improvement in performance, but at the cost of accuracy, of course.

What I'm trying to do is similar to gravitational force calculations. Someone in this thread suggested working with complex numbers and using the equation:

(z-z')/abs(z-z')**3

However, I don't know what z refers to and I cant find reference to this equation anywhere else.

EDIT: Though I've accepted an answer below, I ended up using a solution that nobody mentioned: a look-up table. Basically, I pre-calculate the force exerted on one particle by another particle at all distances (with a resolution of 1 pixel) up to 2000 pixels away in any direction. Forces can then be determined without any runtime distance calculations. Though it's not an exact solution, the distortion is not perceivable and it has allowed me to increase the number of simulation particles three-fold (to the point where the number is now limited by the physics engine and not by the n-body problem.)

Answer 1

1. First up, the force is symmetric so you only need n squared divided by 2 calculations.

2. Second up the force falls off with the distance squared, so you can approximate this by setting a threshold beyond which there is no need to calculate. You only need check this in 1 dimension. As in, if x - x' < t and y - y' < t calculate, otherwise no.

3. You can use two occupancy trees to put the particles into little squares, or buckets, and only do k squared divided by 2 distance calculations for the pairs in each bucket in one tree. The other tree is translated up and to the left, to take care of any pairs on the bounds of adjacent squares. For the second tree just don't repeat the distance calcs you did in the first tree.

4. Between updates it is possible that some particles may not change their position, therefor you do not needn't calculate again any forces that involve two particles stationary since te last update.

5. One may assume the masses also contribute. A further optimization is if the mass is too small then the threshold t beyond which to not calculate can shrink, too.

6. You may also want to quantize the translations of the particles so that any move below the minimum is not registered as a change in coordinates .

Hope this helps!

Answer 2

To calculate the forces between all N particles exactly you must determine the O(N^2) particle-particle distances. There's no way to do this better than O(N^2), although you can achieve (1/2)N^2 by noting symmetry (i.e. dist(i,j) == dist(j,i)).

To do almost the same thing by cheating you can instead calculate only the "near-by" particle interactions (these are the ones of largest magnitude anyway) and make use of spatial indexing structures.

One idea is to make use of quadtrees (octtrees in 3d) to allow you to efficiently determine a small set of M particles that are "close" to a particular particle. When calculating the particle forces you only consider this reduced set. You will need to experiment with the threshold used to determine whether the set M is large enough - obviously small sets will be very fast, but inaccurate, large sets will be accurate but slow.

The spatial index will also need updating as the particles move. One way to do this with quadtrees is to set an upper and lower particle count per leaf and refine the tree when these constraints are violated.

Overall, based on reasonably well distributed particle sets, the height of the tree will be O(log(N)) and you should expect O(N*log(N)) performance from the tree-based method.

Answer 3

Since you mention that cheating is OK, it might be possible to fake things by calculating the attraction between a single particle and the center of mass of the entire group (including the particle itself). This would be a O(2N) problem (N for the center of mass, plus another N for each particle's experienced force). This might deviate significantly from the actual laws of motion for certain arrangements of particles, but it might be a reasonable approximation for most (random) initial arrangements.