n-bodies: a parallel TBB solution: computing accelerations

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

Last time I fleshed out a very simple simulator loop that alternately computes body interactions and then applies them. I coded a version of the application part, using a loop to update the vector components of velocity and position.

This time we’ll shift to the interactions part. The gravitational force equation mentioned a couple posts ago acts on pairs of bodies: the magnitude of the force between body i and body j is just the same as that between body j and body i, just in the opposite direction. My first optimization is to avoid that duplication: computations for each pair of bodies will only be done once: thinking of the bodies i and j, I can eliminate the duplications by respecting a relationship between i and j. Just restrict it so that either one is always larger or it’s always smaller to ensure each pair is touched only once. I’ll choose j smaller, as indicated by the green triangle to the right. Future posts may refer to this ij-space as interaction space and call the green triangle the interaction triangle. I don’t care about the case when i and j are equal: the gravitational problems there are beyond the scope of this point model. :-)

Turning that into code should produce something close to this:

    for (double s = 0.; s < STEPLIMIT; s += TIMESTEP) {
        int i, j;

        // Compute the accelerations of the bodies
        for (i = 0; i < n – 1; ++i)
            for (j = i + 1; j < n; ++j)
                addAcc(i, j);

        // Apply the accelerations and update body positions
        for (i = 0; i < n; ++i) {
            for (j = 0; j < 3; ++j) {
                body[i].vel[j] += body[i].acc[j] * TIMESTEP;
                body[i].pos[j] += body[i].vel[j] * TIMESTEP;
                body[i].acc[j] = 0.;
            }
        }
    }

Note the addition of a function taking i and j as arguments, in which I’ll encapsulate the body interaction code. I’ll defer until next time whether this should be called addAcc() or addForce(). ;-)

Next time: computing accelerations? or forces?