A constant-number direct simulation Monte Carlo (DSMC) model was developed for the analysis of nanoparticle (NP) agglomeration in aqueous suspensions. The modeling approach, based on the "particles in a box" simulation method, considered both particle agglomeration and gravitational settling. Particle-particle agglomeration probability was determined based on the classical Derjaguin-Landau-Verwey-Overbeek (DLVO) theory and considerations of the collision frequency as impacted by Brownian motion. Model predictions were in reasonable agreement with respect to the particle size distribution and average agglomerate size when compared with dynamic light scattering (DLS) measurements for aqueous TiO(2), CeO(2), and C(60) nanoparticle suspensions over a wide range of pH (3-10) and ionic strength (0.01-156 mM). Simulations also demonstrated, in quantitative agreement with DLS measurements, that nanoparticle agglomerate size increased both with ionic strength and as the solution pH approached the isoelectric point (IEP). The present work suggests that the DSMC modeling approach, along with future use of an extended DLVO theory, has the potential for becoming a practical environmental analysis tool for predicting the agglomeration behavior of aqueous nanoparticle suspensions.