Abstract A novel parallel simulation algorithm on the GPU, implemented in CUDA and C++, is presented for the simulation of Brownian particles that display excluded volume repulsion and interact with long and short range forces. When an explicit Euler–Maruyama integration step is performed to take into account the pairwise forces and Brownian motion, particle overlaps can appear. The excluded volume property brings up the need for correcting these overlaps as they happen, since predicting them is not feasible due to the random displacement of Brownian particles. The proposed solution handles, at each time step, a Delaunay triangulation of the particle positions because it allows us to efficiently solve overlaps between particles by checking just their neighborhood. The algorithm starts by generating a periodic Delaunay triangulation of the particle initial positions on CPU, but after that the triangulation is always kept on GPU memory. We used a parallel edge-flip implementation to keep the triangulation updated during each time step, checking previously that the triangulation was not rendered invalid due to the particle displacements. We designed and implemented an exact long range force simulation with an all-pairs N -body simulation, tiling the particle interaction computations based on the warp size of the target device architecture. The resulting implementation was validated with two models of active colloidal particles, also showing a speedup of up to two orders of magnitude when compared to a sequential implementation. A short range forces simulation using Verlet lists for neighborhood handling was also developed and validated, showing similar performance improvements.
In this paper we discuss the nonequilibrium shear viscosity of a suspension of hard spheres that is modeled by neglecting hydrodynamic interactions in a consistent way. The aim is to establish the true capabilities of this model in predicting the properties of real suspensions. A Brownian dynamics algorithm is used to simulate the movements of hard spheres immersed in a Newtonian solvent in a nonequilibrium steady shear flow. A new development is the treatment of the overlap of spheres as elastic collisions, to simulate the no-flux boundary condition on the surfaces of rigid particles. This algorithm is compared with other algorithms suggested in the literature, and is shown to be simple and accurate even for two spheres at close distance. This provides an algorithm that is very suitable for calculating the pair distribution function and especially its hard-sphere contact value, both in equilibrium and nonequilibrium simulations. The algorithm is used to study the nonequilibrium stationary shear flow in the low shear limit. The simulations correctly reproduce the exact low-density limit of the perturbation of the pair distribution function. The perturbation of the pair distribution function in shear flow can be extracted from the simulation data and used to compute the stationary shear viscosity for a system of diffusing hard spheres without hydrodynamic interactions. This yields a flow curve for this model system including the low shear limit. It is found that the model shear viscosity fails at intermediate and high shear rates as can be expected from the neglect of hydrodynamic interactions, but also in the low shear limit at small and moderate volume fractions.