n-bodies: a parallel TBB solution: realizing addAcc(i,j)

By Robert Reed (Intel) (25 posts) on September 25, 2009 at 3:42 pm

Having settled the question of whether I should accumulate forces or accelerations last time, now it’s time to build the accumulation function.

void
addAcc(int i, int j) {

i and j are indices selecting elements of the body array. First task is to compute the distance between them.

        double dx = body[i].pos[0]-body[j].pos[0];
        double dy = body[i].pos[1]-body[j].pos[1];
        double dz = body[i].pos[2]-body[j].pos[2];
        double distsq = dx*dx + dy*dy + dz*dz;

Pythagorean Theorem in three dimensions gets me the square of the hypotenuse, but before doing the square root, I’ll avoid the singularity:

if (distsq < MINDIST) distsq = MINDIST;
double dist = sqrt(distsq);

That is, if the point masses get too close together, act like they’re not. But wait! Why do I even need the square root, if I’m working with gravitation, an inverse-squared law? Well, because acceleration is a vector so I need the next step.

        double ud[3];
        ud[0] = dx/dist;
        ud[1] = dy/dist;
        ud[2] = dz/dist;

Array ud represents the unit vector (length 1 direction vector) pointing from body j to body i. I need just one more thing, the magnitude of those accelerations.

        double Gdivd = GFORCE/distsq;
        double ai = Gdivd*body[j].mass;
        double aj = Gdivd*body[i].mass;

All that’s left is to compute the acceleration vector components and apply them to the bodies.