We develop a Fluctuating Immersed Boundary (FIB) method for performing Brownian dynamics simulations of confined particle suspensions. Unlike traditional methods which employ analytical Green's functions for Stokes flow in the confined geometry, the FIB method uses a fluctuating finite-volume Stokes solver to generate the action of the response functions "on the fly." Importantly, we demonstrate that both the deterministic terms necessary to capture the hydrodynamic interactions among the suspended particles, as well as the stochastic terms necessary to generate the hydrodynamically correlated Brownian motion, can be generated by solving the steady Stokes equations numerically only once per time step. This is accomplished by including a stochastic contribution to the stress tensor in the fluid equations consistent with fluctuating hydrodynamics. We develop novel temporal integrators that account for the multiplicative nature of the noise in the equations of Brownian dynamics and the strong dependence of the mobility on the configuration for confined systems. Notably, we propose a random finite difference approach to approximating the stochastic drift proportional to the divergence of the configuration-dependent mobility matrix. Through comparisons with analytical and existing computational results, we numerically demonstrate the ability of the FIB method to accurately capture both the static (equilibrium) and dynamic properties of interacting particles in flow.

Including the effect of thermal fluctuations in traditional computational fluid dynamics requires developing numerical techniques for solving the stochastic partial differential equations of fluctuating hydrodynamics. These Langevin equations possess a special fluctuation-dissipation structure that needs to be preserved by spatio-temporal discretizations in order for the computed solution to reproduce the correct long-time behavior. In particular, numerical solutions should approximate the Gibbs-Boltzmann equilibrium distribution, and ideally this will hold even for large time step sizes. We describe finite-volume spatial discretizations for the fluctuating Burgers and fluctuating incompressible Navier-Stokes equations that obey a discrete fluctuation-dissipation balance principle just like the continuum equations. We develop implicit-explicit predictor-corrector temporal integrators for the resulting stochastic method-of-lines discretization. These stochastic Runge-Kutta schemes treat diffusion implicitly and advection explicitly, are weakly second-order accurate for additive noise for small time steps, and give a good approximation to the equilibrium distribution even for very strong fluctuations. Numerical results demonstrate that a midpoint predictor-corrector scheme is very robust over a broad range of time step sizes.

We introduce numerical methods for simulating the diffusive motion of rigid bodies of arbitrary shape immersed in a viscous fluid. We parameterize the orientation of the bodies using normalized quaternions, which are numerically robust, space efficient, and easy to accumulate. We construct a system of overdamped Langevin equations in the quaternion representation that accounts for hydrodynamic effects, preserves the unit-norm constraint on the quaternion, and is time reversible with respect to the Gibbs-Boltzmann distribution at equilibrium. We introduce two schemes for temporal integration of the overdamped Langevin equations of motion, one based on the Fixman midpoint method and the other based on a random finite difference approach, both of which ensure that the correct stochastic drift term is captured in a computationally efficient way. We study several examples of rigid colloidal particles diffusing near a no-slip boundary and demonstrate the importance of the choice of tracking point on the measured translational mean square displacement (MSD). We examine the average short-time as well as the long-time quasi-two-dimensional diffusion coefficient of a rigid particle sedimented near a bottom wall due to gravity. For several particle shapes, we find a choice of tracking point that makes the MSD essentially linear with time, allowing us to estimate the long-time diffusion coefficient efficiently using a Monte Carlo method. However, in general, such a special choice of tracking point does not exist, and numerical techniques for simulating long trajectories, such as the ones we introduce here, are necessary to study diffusion on long time scales.