AbstractA method is presented for developing multiple particle breakage distribution functions that include the shape factor, as well as the particle size. The technique is to first develop a joint probability distribution that accounts for the size and shape factor of each child particle. This is then reduced by successive integration to the marginal probability used in the breakage equation. This technique guarantees mass conservation and exchangeability of child particles while allowing the user to choose the number of child particles formed per breakage. It also allows the user to choose the functional form for the size distribution of particles as well as independently choosing the functional form for the shape factor distribution of particles. To solve the breakage equation, a discretization method is used that conserves particle mass during breakage and correctly predicts changes in the total number of particles. Simulation results show that the shape distribution functions may yield distinguishable shape factor profiles after crushing. Comparison with experimental data indicates that these theoretical functions can be used to model real systems. © 2004 American Institute of Chemical Engineers AIChE J, 50: 937–952, 2004