  Standard Eulerian--Lagrangian (EL) methods generally employ drag force models
that only represent the mean hydrodynamic force acting upon a particle-laden
suspension. Consequently, higher-order drag force statistics, arising from
neighbor-induced flow perturbations, are not accounted for; with implications
on predictions for particle velocity variance and dispersion. We develop a
force Langevin (FL) model that treats neighbor-induced drag fluctuations as a
stochastic force within an EL framework. The stochastic drag force follows an
Ornstein-Uhlenbeck process and requires closure of the integral time scale for
the fluctuating hydrodynamic force and the standard deviation in drag. The
former is closed using the mean-free time between successive collisions,
derived from the kinetic theory of non-uniform gases. For the latter,
particle-resolved direct numerical simulation (PR--DNS) of fixed particle
assemblies is utilized to develop a correlation. The stochastic EL framework
specifies unresolved drag force statistics, leading to the correct evolution
and sustainment of particle velocity variance over a wide range of Reynolds
numbers and solids volume fractions when compared to PR--DNS of freely-evolving
homogeneous suspensions. By contrast, standard EL infers drag statistics from
variations in the resolved flow and thus under-predicts the growth and steady
particle velocity variance in homogeneous suspensions. Velocity statistics from
standard EL approaches are found to depend on the bandwidth of the projection
function used for two-way momentum coupling, while results obtained from the
stochastic EL approach are insensitive to the projection bandwidth.
