Large scale simulations of particulate systems using the PME method

This work describes the development of a novel O (N_p ln N_p) Particle-Mesh-Ewald method applied to concentrated systems of hydrodynamically interacting rigid particles in Stokes flow. The algorithm is based on the P3M approach (particle-particle, particle-mesh) commonly employed in plasma physics and applications concerned with electrostatic particle interactions.; The method calculates interactions between particles in a fashion similar to the popular Stokesian dynamics algorithm, but in contrast to Stokesian dynamics, the many-body interactions do not require costly matrix inversion to compute. More specifically, near-field interactions are included through the addition of the exact pair-wise resistance functions, while the many-body interactions are included using the method of reflections, combined with a particle-mesh solution of the Stokes equations. The far-field interactions are computed by translating the particle force moments to a fixed grid and computing the induced velocity field using FFT solution of the Stokes equations. The velocities at the particle positions are evaluated by interpolating from the known values at the cubic grid. Rigorous error estimates are presented, demonstrating the equivalence of this procedure with the explicit evaluation of Ewald sums for the exact multipole distribution.; The procedure may be used to compute individual particle forces on a collection of particles with known velocities. To avoid solution of the linear system required to determine particle velocities from the known forces, the method uses the Langevin equation, retaining the acceleration term. Given this choice, the time stepping procedure is subject to a severe stability constraint based on the magnitude of the lubrication forces.; The computational load imposed by this constraint is eased by emplying a hierarchical time stepping procedure with different size time steps for the near field and far field particle interactions. The overall computational effort of the algorithm scales as O(N_p ln N_p) where N_p is the number of particles.; Comprehensive tests of the algorithm have been conducted demonstrating its accuracy and speed relative to Stokesian dynamics and spectral boundary integral computations. The speed of the algorithm is such that simulations of 4096 particles may be run on an IBM RS6000 workstation. Memory requirements are modest, and extension of the method to multiple-processors architectures is relatively straightforward. Results are discussed for both Monte Carlo and dynamic simulations for the sedimentation and rheology of concentrated hard-sphere suspensions.