Speeding up particle simulation

[fargoth]

Hi, I've got a list of particles (where each particle has it's coordinates in space and velocity)
these particles can interact with each-other ($$f=\frac{a}{r^4}-\frac{b}{r^8}$$)
(i'll add 6 - 12 soon too)
anyway, as it is now, the program has a complexity of $$O(\frac{n(n-1)}{2})$$
because every pair has to be calculated to get the force.

the problem is - its too slow for realtime simulation - and i want it to be interactive (changing the temperature, zooming in and out, etc.)

so, i started thinking about approximations i could do to make it work with O(n), and still get good results...
I think calculating the force for near neighbors is enough - but then i need a quick way to know who those neighbors are for a given particle (in less then O(n)).

I thought i'd use an octree (for 3D case) or quadtree (for 2D case), but i can't think of a way to maintain it - because the particles are always moving from place to place... i need to rebuild it each time...

I hope i was clear enough.. my question is - do you think my approach is the best one?
do you have any suggestion as to how to maintain the tree, or alternatively use a better approximation?

 

[las3rjock]

http://math.nyu.edu/faculty/greengar/shortcourse_fmm.pdf

 

[Jheriko]

don't forget to optimise the code as well as the algorithm, this way you can squeeze the most out of which ever algorithm you choose (algorithmic optimisation is more important in general though)

use a vector math library if possible as it might optimise to use threads and simd instruction

also try to reduce divisions and other expensive math operations, e.g. to calculate f:

var = 1 / r*r;
var = var * var;
f = a * var;
var = var * var;
f = f - b * var;

3 multiplications is faster than 2 calls to power functions and there is only one division, the other has been exchanged for 2 multiplications in essence...

you can optimise further by accumulating the force from all of the particles before moving each one, you can remove one of the multiplies if you use 1/r^4 - (b/a)/r^8 as well, then multiply the total force by a afterwards (calculate b/a once in a variable)

in the inner loop, don't take the square root of r if possible, since you probably use its square a lot, you might be able to get away with it.

another simple trick is to not count distant particles. you are calculating r in the inner loop, if you do this first you can early out if it is above some value, gives the same effect as rounding small values to zero.

if precision doesn't matter much then use only 32-bit floats for a minor performance gain.

hope this helps.

 

[fargoth]

don't forget to optimise the code as well as the algorithm, this way you can squeeze the most out of which ever algorithm you choose (algorithmic optimisation is more important in general though)

use a vector math library if possible as it might optimise to use threads and simd instruction

also try to reduce divisions and other expensive math operations, e.g. to calculate f:

var = 1 / r*r;
var = var * var;
f = a * var;
var = var * var;
f = f - b * var;

3 multiplications is faster than 2 calls to power functions and there is only one division, the other has been exchanged for 2 multiplications in essence...

you can optimise further by accumulating the force from all of the particles before moving each one, you can remove one of the multiplies if you use 1/r^4 - (b/a)/r^8 as well, then multiply the total force by a afterwards (calculate b/a once in a variable)

in the inner loop, don't take the square root of r if possible, since you probably use its square a lot, you might be able to get away with it.

another simple trick is to not count distant particles. you are calculating r in the inner loop, if you do this first you can early out if it is above some value, gives the same effect as rounding small values to zero.

if precision doesn't matter much then use only 32-bit floats for a minor performance gain.

hope this helps.

thank you very much.
I will spare one division in the inner loop thanks to you..

(EDIT: after testing to see what impact it has on the speed i can tell you that using multiplications instead of divisions doesn't help - on the contrary - on my system it gave me 5% penalty... and sqrtf is VERY efficient, if i use it in the inner loop i only lose about 5% performance)

i'm already using the distance limit you've suggested - but when i cool the system down, and everything is squeezed to liquid or solid - the distance limit includes them all - i could make my distance limit change according to temperature, but i think using an efficient algorithm for finding nearest neighbors will get me better accuracy, and it might even be a better...

i'm looking into K-D trees now... see if i can use it to find nearest neighbors efficiently.

 

[christianjb]

This is probably an order n problem- because there are no long range forces.

You only have to sum over interactions which are less than rmax apart. For instance- it might be that the simulation converges by summing over a sphere of ~100 particles around each target particle.

The problem is how to work out which particles are within the sphere. This can be done using neighbor lists.

 

[fargoth]

http://math.nyu.edu/faculty/greengar/shortcourse_fmm.pdf

thanks, I've gone over the text.
I knew I'm in the right direction.
(using the adaptive nearest neighbor approach to get $$O(n log(n))$$)
but I'm having trouble implementing it... and the paper doesn't talk about implementation...
oh well, i guess I'll have to tackle it some more =)