- The ability to design a size reduction system prior to full scale experiments and to optimize existing systems has long been a goal of designers. Such a design and optimization could be achieved by correctly simulating any system under any operating condition. In this paper we present a new and innovative procedure to implement empirical comminution functions into DEM–CFD simulations. The paper is focused on the implementation procedures and not the DEM/CFD simulations, which deserve full attention. Therefore, this paper is not aimed to study any specific mill. The comminution functions include: initial strength distribution, selection function, breakage function and fatigue function. First, the traditional comminution functions (strength distribution, selection and breakage functions) and the recently investigated fatigue function are briefly described and modified. Then a procedure for implementing the functions into a DEM–CFD model or any other source to provide impact velocities and number of impacts, is described in detail. The implementation involves converting the probability comminution functions into individual particle properties by a random method and then converting the velocity dependent comminution functions into strength dependent ones. In this way, and mainly owing to the use of the fatigue function (which defines the weakening of those particles that are not breaking), a real size reduction system, in which each particle is subjected to multiple impacts at various velocities can be simulated. Three case studies for multiple impact conditions at the same average velocity (several impacts at the same velocity, various velocities at each impact and randomly selected velocities) are presented and analyzed in order to confirm qualitatively the procedure, although the comminution functions need to be further quantitatively modified. It should be emphasized that although the new procedure presents a step towards the final goal, some limitations do exist and some questions remain open.