Discrete linear population balance models have been widely used for comminution processes and particle breakage tests that involve a series of breakage events. The fundamental tenet of these models is the independence of particle breakage from the surrounding particle population. On the other hand, ample experimental data have shown that multi-particle interactions among particles of different sizes at the particle ensemble scale lead to a population-dependent breakage probability. In an attempt to analyze and model the multi-particle interactions, we formulate a time- or space-discrete and particle size-discrete non-linear population balance model using a recently introduced non-linear kinetics framework for rate-based comminution processes. The discrete model is developed as a limiting case and presented in equivalent matrix and index notation. The practical applicability, accuracy, and stability aspects of the model are discussed. With an empirical, flexible, non-linear functional accounting for multi-particle interactions explicitly, the proposed model can serve as a framework for quantitatively analyzing multi-particle interactions. Numerical simulations with a simplified discrete model have provided significant insight into the retardation effect of fines on the breakage probability of the coarser particles experimentally observed in dense-phase particulate systems such as particle bed compression tests. (C) 2011 Elsevier B.V. All rights reserved.