Damage and mass transport behavior of particulate materials depend on the spatial distribution and orientation of the phase fractions and the properties of their interfaces. Discretization of this three-dimensional setting is the first critical step in the modeling of such behavior. In this paper, a general procedure is developed for discretizing particulate materials, such as particles embedded within a cementing matrix, in which the particle inclusions are assumed to be randomly generated polyhedra. Polyhedron discretization is based on the Voronoi tessellation of a set of points, which is obtained through a multistep process that relies on the concept of point saturation. The matrix-particle and/or particle-particle interfaces are explicitly constructed. To demonstrate the utility of the approach, uniaxial tensile tests (of a single particle, embedded in a matrix) are simulated using a lattice model. It is seen that interface debonding is strongly affected by inclusion shape and orientation with respect to the tension field. A procedure is developed for the discretization and analysis of multi-particle systems. (C) 2011 Elsevier B.V. All rights reserved.