DEM investigation of the effect of particle breakage on compact properties

Particle breakage during compaction affects compaction behavior and the quality of the formed compact. This work conducted a numerical study based on the discrete element method (DEM) to investigate the effect of particle breakage on compaction dynamics and compact properties, including particle size and density distributions, and pore properties. A force-based breakage criterion and Apollonian sphere packing algorithm were employed to characterize particle breakage behavior. The pore structures of the compacts were characterized by the watershed pore segmentation method. Calibrated with experimental data, the model was able to simulate the stress-strain relation comparable with experimental observation. During compaction, the particles were gradually broken from top to bottom with increasing pressure. Both density and pore size of the compacts had relatively uniform distribution at larger stress, while the pore size decreased sharply when the particles started to break, indicating that the smaller fragments in the compact system have a significant effect on the pore size distribution.


Introduction
Compaction of powders is widely used in fields such as pharmaceutical tableting, soil consolidation, mineral processing, and powder metallurgy [1]. At large pressure, particles may experience breakage, which affects not only compaction behavior but also the structural and mechanical properties of formed compacts [2]. A better understanding of particle breakage during compaction is critical to process design and product quality control. Previous experimental studies on particle breakage in compaction were mainly focused on individual particle breakage patterns and the effects of initial density and particle size distribution (PSD) on particle breakage [3,4]. Their results showed that the compaction curve and breakage pattern depended on the loading and initial PSD. However, particle breakage mechanisms are difficult to identify in experiments due to the lack of information at the particle scale.
On the other hand, the discrete element method (DEM) which considers particle breakage through fragment replacement method has been proposed to study the effects of failure criterion and materials on particle breakage and compaction behavior have been investigated [5][6][7]. But studies on the effect of particle breakage on compact structures during the compaction process are still limited.
In this study, we conducted a DEM study of particle compaction in which particles were subject to breakage upon large force. A force-based breakage model was implemented in the DEM model, and an image-based algorithm was adopted to analyze pore structures in a compact. In particular, the compact response and pore structure at different stages were analyzed.

DEM model
The details of the current DEM model can be found in our previous study [8]. Only the particle breakage model is described here. The breakage model in this work consisted of three stages: determination of particle breakage, fragment replacement, and relaxation. Failure criteria based on stress [5,9] and force [6,10] have been proposed to determine particle breakage. However, the stress-based criteria are unable to simulate the contact related phenomena -such as stress concentrations [11]. Thus, a force-based breakage criterion was adopted to determine if a particle of size d and strength lim breaks under a compressive force F. The particle experiences plastic deformation and is subject to breakage when the induced normal stress is larger than its strength, given by [12], where m is a parameter related to particle material properties, which is often obtained based on singleparticle breakage experiments. lim,0 is the mean limiting strength of the particles at the reference size d0, the mean limiting strengths of particles follow the normal distributions. AF is the contact area of two contacting particles. Note that Eq. (1) is derived from the Weibull theory [13], and it has been validated by other works [6,7].
*Corresponding author: r.yang@unsw.edu.au In this work, we set the minimum size of breakable particles based on the work conducted by [6] which observed that when the minimum size of breakable particles was set to 0.25 0 (soft materials) or 0.55 0 (hard materials such as silica sand), the mechanical response under compaction was less affected. Thus the same minimum breakable particle sizes were adopted in the work to reduce computational cost.
In the simulations, when a particle broke, it was replaced by smaller fragments (spheres in the work) of different sizes. The Apollonian sphere packing algorithm was used in the work as: i) it is the only structure which can effectively fill up the whole space with various sizes of spheres; ii) it has a fractal dimension of 2.47, similar to a fully fragmented particle; and iii) the positions and sizes of individual particles are pre-determined, which can significantly reduce the computational cost to replace the original particles with new particles.
In the relaxation stage, the sizes of the new particles gradually increased and their positions were adjusted based on force balance until the total mass of new particles is about 90% of the mass of the original particles. The relaxation process, normally took about 20 time steps in the simulations, significantly maintained the mass and compact properties during breakage. The 10% mass loss was treated as generation of fines which has less impact on the compact response.

Simulation conditions
A simulation began with the random generation of spheres inside a cylindrical container. After a stable packing was formed, a loading plate moved from the top at a speed of 0.05 m/s until the prescribed pressure was reached. In compaction, particle breakage was checked every 2 ms, and information, such as the pressure on the plate, particle size distribution, compact density, and inter-particle forces were monitored and recorded for post-simulation analysis. The silica sand particles were selected in the simulations and their parameters were determined based on previous experiments [4]. The initial PSDs and strengths in the simulations were similar to those in the experiments so that the results could be compared. The material parameter m and mean limit strength σlim,0 in Eq. (1) were calibrated based on the experiment data.
The tests showed that with m = 10, σlim,0 = 3.5 GPa (reference diameter, d0 = 2 mm) the simulation results agreed well with the experiments. Table 1 lists the key particle properties and the simulation conditions.

Pore structure analysis
The pore structure was characterized by the watershed pore segmentation method which is a region-based method for pore analysis [14]. The procedure of the watershed method was as follows: Firstly, a compact structure was converted to the binary volume data to separate pores and particles (the grid size is 1/20 of the minimum particle diameter). Then the distance between the void pixels and the solid pixel were calculated. Finally, the watershed algorithm was applied to the distance data to obtain the ridge points (defined as the local minimum values between two nearly deepest values (Catchment basins). The pixels of a segmented pore (a pore body) were labeled with the same numbers and the pore volume is the sum of pixels in pores body, its center position x, y, z corresponding to the center of mass of the horizontal, vertical, and planar coordinates of pore body.
The watershed pore segmentation method was tested by analyzing the structures of an ordered hexagonal close-packed (HCP) structure. For the HCP structure, two main pore structures, which have sizes of 0.36d (d is particle diameter) and 0.63d, were identified. The two types of pore were similar to the theoretically obtained analytic values [15] (0.37d and 0.63d respectively), indicating the validity of the method.

Sensitivity test and model validation
Sensitivity tests were conducted to investigate the effect of the number of fragments replacing the original particles on the compaction curve and final PSD. 8, 20, and 57 spheres, corresponding to generation 1, 2, and 3 of the Apollonian packing (see inset of Fig. 1(a)), were adopted in the simulations. As shown in Fig. 1, the compaction curves are less affected by the number of fragments, while the final PSDs under at 5 MPa pressure are clearly different, particularly for particles replaced with 8 spheres. On the other hand, the PSDs for 20 and 57 spheres are similar, consistent with the previous findings [6]. To balance model accuracy and computational cost, 20 fragments are adopted in the simulations. Fig. 2(a) shows the variations of pressure with porosity for silica sands when particle breakage is or is not considered. The compaction curves with and without particle breakage initially follow the similar trends but diverge at 10 MPa. The compact response with particle breakage is much softer than when no particle breakage is considered. The decrease in porosity with pressure is much faster due to the generation of small fragments that fill up pores. The figure also shows the simulation results have a good agreement with the experiment data [4]. Fig. 2(b) shows the comparisons of the simulated and experimentally measured PSDs compressed under different pressures. With larger pressures, more particles are broken, leading to finer size distributions. The trends in the simulations show a good agreement with the experiment [4], indicating the current model is able to predict the breakage behavior of particles during compactions.

Evolutions of compacts
After the model is calibrated and validated with the experimental data, it is used to simulate particle compaction under different pressures, as shown in Fig.  3. At the initial state (compact A) before compaction, the particles are evenly distributed in the compact. When the compaction starts (compact B), the particles near the compact plane start to break, and the smaller particles are generated and interacted with the other particles in the compact. With compaction proceeds (compact C), the breakage gradually spreads to the bottom of the compact. Finally, at large pressure (compact D), most of the large particles in the compact are broken.
The densities for compacts A-D are 0.572, 0.596, 0.667, 0.743, respectively (the mass loss during particle breakage was included in the calculations), and the mean particle sizes are 1.595 mm, 1.295 mm, 0.83 mm, and 0.56 mm respectively. It shows that the mean particle size decreases with pressure and the density increases with pressure except for the final compact.   4 shows the evolution of particle size distribution during the compaction. The particle sizes in the compact gradually decrease and the percentage of smaller particles increases with the pressure. This figure is also consistent with the results in Fig. 3. Fig. 5 shows the density distribution of the compacts at different heights (not include the fine particles). The compact density increases with the height for the initial compact and compact B. Compact B has a larger variation between the top and bottom of the compact because the particles near the compaction plane are broken. With increasing pressure (compact C), the compact has a relatively evenly density distribution while the density at bottom of the compact is still smaller due to the large particles at the bottom of the compact (see Fig. 3(c)). At larger compact stress, the compact has an evenly density distribution with compact height because all the weak particles are broken in the compact.   6 shows the pore size distribution with the compact height using watershed pore segmentation. The pore distribution shows some similarities with the density distribution. For example, compacts C and D have uniform distribution for both compact density and pore size. For compact B, the pore size and compact density have larger variations between the top and bottom of the compact. Also, the pore size decreases with the compact height while the density increases with the compact height. Besides, the pore size sharply decreases when the particles start to break, indicating that the smaller particles in the compacts have significant effect on the pore size distribution.

Conclusions
A DEM model coupled with a particle breakage model was developed to simulate particle fragmentation under compaction. The model was calibrated based on the experimental results in terms of compaction curve and PSD. Subsequently, this model was applied to compacts with different pressures. During the compaction, the particles are gradually broken from the top (near the compaction plane) to the bottom of the compacts with increasing pressure. In addition, the compacts had relatively uniform distributions for both density and pore size at larger stress. The pore size decreases with compact height and the density increased with compact height, while the pore size sharply decreases when the particles start to break, indicating that the smaller particles have a great effect on pore size distribution. Pore distribution with the compact height, the inset is the pore segmentation for compact A, the black color represents the particles.