The Challenge | The Approach | Publications | Investigators | Research Sponsors |

Trajectories of 10 bodies interacting through a gravitational force over 5000 time steps. Initial positions and velocities follow a 3D Plummer distribution. |

A central problem in many physical simulations is the evolution of a system of many elements or bodies interacting over large distances. This problem arises in cosmological simulations, where the long-range force is gravity, in molecular dynamics simulations, where the long range force is primarily electrostatics, and in many other areas including vortex formation in computational fluid dynamics and radiosity computation in computer graphics.

The simplest approach to determine the total force experienced by
each body due to all other bodies in the system requires O(*N*^{2})
work per integration timestep and quickly becomes intractable for large systems
that must be evolved over many timesteps.
The celebrated fast multipole methods reduce the cost to O(*N* log *N*)
or even O(*N*) by the hierarchical approximation of the force due to a
group of bodies at a distance. These
methods have rigorous error bounds but require a spatial decomposition of the
bodies to be created and maintained.
Furthermore, the efficient implementation of the fast multipole methods
for 3D long-range forces is extraordinarily complex, particularly when the
bodies are non-uniformly distributed.
In general, several thousands of bodies are required before the methods
become competitive with alternatives with larger asymptotic complexity.

In our case, we are particularly interested in molecular dynamics
simulations, where N may be relatively small (on the order of a few thousand
atoms), the number of integration time steps is very large (on the order of 10^{6}
or more), and the accuracy required is very high.
For this class of problems in particular, as well as N-body
problems in general, we are investigating a new class of algorithms that are
adaptive in time, and have work complexity that grows as O(*N*^{4/3})
or even O(*N*) per integration step when amortized over
(*N*)
steps.

Instead of using a fixed time step for all interactions (which is
too conservative for most interactions), an approach that equalizes impulse is
followed. Mathematically, the constant
time step *t* is traded for constant impulse *I* defined as *Fij*
×*Tij*,
where *Fij* is the force between particles *i* and *j*, and *Tij*
becomes the time step used to re-evaluate *Fij*.
Using a constant impulse instead of a constant time step leads to
an amortized execution complexity of O(*N*^{4/3}) per simulation
step. Algorithmic improvements that
rely on the first and second derivatives of force further reduce the amortized
work per step to O(*N* log *N*) and O(*N*), respectively.

Distribution of the 45 pairwise distances between bodies after 100 time steps in the system shown previously. |

Schedule for the next force evaluation according to the equal impulse method. Over half of all interactions need not be re-evaluated within the next 10 timesteps, and 20% of the interactions need not be re-evaluated within the next 100 timesteps. |

- Lars Nyland (PI), research associate professor, Computer Science
- Susan Paulsen, postdoctoral researcher, Computer Science
- Jan F. Prins, professor, Computer Science
- William V. Wright, research professor emeritus, Computer Science
- Jan Hermans, professor, Biochemistry and Biophysics