Modeling of the Particle Abrasion Process and a Discrete Element Method Study of Its Shape Effect

This study introduces a novel method for particle abrasion derived from fundamental natural phenomena and mechanical principles, allowing precise control over the degree of abrasion and more accurately mimicking natural processes. The method’s validity is confirmed using a specific shape index. Through conventional triaxial tests, the mechanical behavior of granular aggregates with varying degrees of abrasion was analyzed. The findings indicate that increased particle abrasion leads to a decrease in the average coordination number and sliding amount, while the rotation amount increases. This suggests an inverse relationship between the degree of abrasion and the structural stability and interlocking of the particle aggregate. The fabric anisotropy of the system is mainly attributed to the anisotropy of the contact normal force, which decreases as particle abrasion increases. The partial stress ratio of the particle system is influenced by fabric anisotropy and remains independent of particle shape. Additionally, the internal friction angle may be overestimated in conventional triaxial tests.


Introduction
Rock and soil particles formed in upstream rivers or by blasting have irregular shapes and sharp angles.Most pebbles in downstream rivers, after long-term scouring and abrasion by water flow, become ellipsoidal and have smooth surfaces.Since differences in particle morphology can lead to variations in the mechanical properties of their aggregates, accurate particle construction in numerical simulations is crucial to ensuring the reliability of the results.
Researchers have confirmed that aggregate shape significantly affects the mechanical properties of granular materials, such as shear strength, deformability, and breakage effects, through both indoor experiments [1][2][3][4][5][6] and numerical methods [7][8][9][10].Since spherical particles do not reflect interlocking and occlusion between particles, many researchers have developed various methods to generate non-spherical granular aggregate models for DEM simulations.Some researchers used polygons and polyhedra to represent granular materials [11][12][13][14], while others have built two-dimensional ellipsoidal particles [15,16] and expanded these to three-dimensional ellipsoidal shapes [17,18].Compared with spherical particles, these shapes can partially reflect the interlocking effect between particles.Additionally, granular aggregates with various aspect ratios have been simulated using a spherocylindrical model.However, these models do not fully capture the surface characteristics of the particles.Therefore, some scholars have proposed more complex particle models.In some cases, the surface properties of particles can be reflected by two-dimensional superquadratic curves based on implicit function expressions [19], and three-dimensional superquadratic surfaces have been used to simulate particles with random forms by improving the expressions of the superquadratic surface functions in each quadrant [20][21][22].These methods, however, do not adequately capture particle asymmetry.Additionally, by using Fourier spherical harmonic functions, some researchers [23][24][25] have created complex particles with asymmetric shapes, but the resulting concave surfaces complicate particle contact identification.Wang [26] proposed a random generation technique for 3D convex aggregates to investigate the impact of particle shape on DEM simulations.This approach involves constructing a curve skeleton using Bezier curve fitting and generating a solid surface using the Gregory surface interpolation algorithm.Nonetheless, the practical implementation of this method can be intricate, and achieving local continuity control and refinement can be challenging.With the help of laser or CT scanning equipment, researchers can extract real particle shapes [27][28][29].However, these irregular particles require a large number of overlapping spheres to simulate surface angles and depressions, making them inefficient for large-scale calculations.In general, algorithms for generating non-spherical aggregates with random shapes remain limited in DEM simulations.
This study introduces a particle generation method for convex polyhedrons based on spheres and develops a particle abrasion function utilizing mesh segmentation.Using this function, particles with varying degrees of abrasion were generated, effectively simulating the shape changes in particles in river channels due to continuous tumbling and collisions under the influence of water flow, thereby replicating the natural process of abrasion.Additionally, triaxial tests were conducted on particles with different abrasion levels to elucidate the impact of particle morphology on mechanical properties from a microscopic perspective.

Polyhedron Particle Generation Algorithm
A particle generation method with a particular shape index is suggested in this study to simulate the randomness and angularity of the gravel particles.Figure 1 shows the modeling method of convex polyhedron particles with random shape.First, a sphere with radius R is created in space at the origin.The sphere's equation is: ls 2024, 17, x FOR PEER REVIEW Convex polyhedrons take on different shapes as a result of t of the points on the sphere, better reproducing the variety of part the degree of angularity of the particle surface gradually decreases of the convex polyhedron increases.Roundness can be used to r particles from a standard sphere; the smaller the roundness, the gre the closer it is to 1, the more similar it is to a standard sphere.Th the spherical surface can also be infinite; when the convex hull's and the particle's theoretical roundness is 1, the particle is roughly dicated that 0.7 is the most likely roundness for an angular partic  Assuming that any random point on the surface of the sphere has the coordinates (x i , y i , z i ), the equation of the plane that is tangent to the surface of the sphere and passes through this random point is: As shown in Equation ( 3), the N-hedron formed by the N tangent points {p i | i = 1, . . ., N} can be expressed as a set of inequality sets if there are N random points on the sphere.
Convex polyhedrons take on different shapes as a result of the random positioning of the points on the sphere, better reproducing the variety of particle shapes.In general, the degree of angularity of the particle surface gradually decreases as the number of faces of the convex polyhedron increases.Roundness can be used to reflect the deviation in particles from a standard sphere; the smaller the roundness, the greater the deviation, and the closer it is to 1, the more similar it is to a standard sphere.The number of points on the spherical surface can also be infinite; when the convex hull's face count is sufficient and the particle's theoretical roundness is 1, the particle is roughly spherical.Ma [30] indicated that 0.7 is the most likely roundness for an angular particle.As a result, the random convex polyhedron examined in this study has 6~14 faces, as shown in Figure 2. Convex polyhedrons take on different shapes as a result of the random p of the points on the sphere, better reproducing the variety of particle shapes.I the degree of angularity of the particle surface gradually decreases as the numb of the convex polyhedron increases.Roundness can be used to reflect the de particles from a standard sphere; the smaller the roundness, the greater the devi the closer it is to 1, the more similar it is to a standard sphere.The number of the spherical surface can also be infinite; when the convex hull's face count is and the particle's theoretical roundness is 1, the particle is roughly spherical.M dicated that 0.7 is the most likely roundness for an angular particle.As a resul dom convex polyhedron examined in this study has 6~14 faces, as shown in Fig Figure 2. Polyhedron particle template.

Particle Abrasion Method
As shown in Figure 3, the Catmull-Clark subdivision algorithm is widel computer graphics to create smooth surfaces from coarse polygonal meshes [ primary function is to iteratively refine a mesh, transforming it into a smoot visually appealing surface.This is achieved by calculating new points for faces, e vertices in each iteration and then subdividing the original faces using these p algorithm is particularly valued for its ability to handle complex topologies and surfaces that are smooth and continuous, making it a standard tool in 3D mod animation for generating high-quality, realistic shapes.

Particle Abrasion Method
As shown in Figure 3, the Catmull-Clark subdivision algorithm is widely used in computer graphics to create smooth surfaces from coarse polygonal meshes [31,32].Its primary function is to iteratively refine a mesh, transforming it into a smoother, more visually appealing surface.This is achieved by calculating new points for faces, edges, and vertices in each iteration and then subdividing the original faces using these points.The algorithm is particularly valued for its ability to handle complex topologies and produce surfaces that are smooth and continuous, making it a standard tool in 3D modeling and animation for generating high-quality, realistic shapes.The algorithm operates iteratively, refining the mesh at each step to produce a smoother result.The process begins with the initial mesh, where each face, edge, and vertex is considered.In each iteration, the following steps are performed: 1.
Calculate the new face point, noted as f -point, which is the average of all vertex positions of the defined face, according to Equation (4).
Calculate the new edge point, denoted as e-point, defined as the average of the two endpoints of the edge and the two f -points of the face sharing the edge, according to Equation ( 5). 3.
Calculate the new vertex point, denoted as v 1 point, defined as the vertex generated for each interior vertex of a given mesh; the new vertex is a linear combination of adjacent face points, adjacent edge midpoints, and vertices, calculated according to the equation: where Q is the average of the new face points sharing the initial vertices.R is the midpoint average of all polyhedron edge positions calculated last time that incident on the polyhedron vertex position calculated last time; S indicates the old vertex; and N is the order of the vertex.Taking the ortho-hexahedral mesh shown in Figure 4 as an example, after generating new vertices, the boundaries are generated in the following order: The grid edge is represented by green lines, and the orthohexahedral is approximated to a sphere after 3 abrasions.It can be found that the method in this study corrodes the mesh vertices the most, the edges the second, and the surface the least, which can better reflect the abrasive process of particles in nature.Figure 5 shows the particle models of the nine convex polyhedron particles in Figure 2 formed after 1 to 3 abrasions, respectively.(In the following study, the initial polyhedron particles will be represented by I_P, and the particles after 1 to 3 abrasions will be represented by Ab_1, Ab_2, and Ab_3, respectively.)The randomly convex polyhedron particles go through abrasion, which is akin to the abrasion process of the hexahedral mesh, to reduce the angularity to some extent and make the particles more rounded overall.The rounding of the particles after abrasion is faster for particles with few faces.Furthermore, it should be noted that this abrasion process is infinite.However, a grid that is too dense will reduce the computational efficiency of numerical simulation; therefore, only 1 to 3 abrasions are examined in this study.Figure 5 shows the particle models of the nine convex polyhedron particles in Figure 2 formed after 1 to 3 abrasions, respectively.(In the following study, the initial polyhedron particles will be represented by I_P, and the particles after 1 to 3 abrasions will be represented by Ab_1, Ab_2, and Ab_3, respectively.)The randomly convex polyhedron particles go through abrasion, which is akin to the abrasion process of the hexahedral mesh, to reduce the angularity to some extent and make the particles more rounded overall.The rounding of the particles after abrasion is faster for particles with few faces.Furthermore, it should be noted that this abrasion process is infinite.However, a grid that is too dense will reduce the computational efficiency of numerical simulation; therefore, only 1 to 3 abrasions are examined in this study.
reduce the angularity to some extent and make the particles more rounded overall.The rounding of the particles after abrasion is faster for particles with few faces.Furthermore, it should be noted that this abrasion process is infinite.However, a grid that is too dense will reduce the computational efficiency of numerical simulation; therefore, only 1 to 3 abrasions are examined in this study.

Particle Shape Quantification
According to Barrett's [33] definition, Equation ( 7) can be used to calculate the roundness evaluation index of particle shape.This evaluation method is straightforward and can more accurately reflect the degree of deviation between particles and standard spheres.As represents the surface area of a sphere equal in volume to the particle; Ap and Vp stand for the particle's surface area and volume, respectively.Figure 6 illustrates how the number of surfaces after abrasion affects the roundness of the polyhedron particles.After abrasion, the roundness of the initial polyhedron particles with varying numbers of surfaces increases.The roundness statistics of convex polyhedron particles with various numbers of surfaces after four regenerations are shown in Figure 7 in order to make the findings more universal.When there are fewer than ten surfaces, the polyhedron's roundness increases quickly with the number of surfaces; however, when there are more than eleven surfaces, the roundness of the particles almost no longer increases with the number of surfaces.The dispersion of the distribution of roundness increases as the surface number of polyhedrons decreases.Following abrasion, the dispersion gradually decreases while the roundness of polyhedron with fewer surfaces grows more quickly.Although the initial polyhedron particles' roundness varies greatly, after three abrasions, all of the particles' roundness stabilize at around 0.91, which amply proves that the method described in this study can accurately simulate the abrasion process of polygonal particles.

Particle Shape Quantification
According to Barrett's [33] definition, Equation ( 7) can be used to calculate the roundness evaluation index of particle shape.This evaluation method is straightforward and can more accurately reflect the degree of deviation between particles and standard spheres.A s represents the surface area of a sphere equal in volume to the particle; A p and V p stand for the particle's surface area and volume, respectively.Figure 6 illustrates how the number of surfaces after abrasion affects the roundness of the polyhedron particles.After abrasion, the roundness of the initial polyhedron particles with varying numbers of surfaces increases.The roundness statistics of convex polyhedron particles with various numbers of surfaces after four regenerations are shown in Figure 7 in order to make the findings more universal.When there are fewer than ten surfaces, the polyhedron's roundness increases quickly with the number of surfaces; however, when there are more than eleven surfaces, the roundness of the particles almost no longer increases with the number of surfaces.The dispersion of the distribution of roundness increases as the surface number of polyhedrons decreases.Following abrasion, the dispersion gradually decreases while the roundness of polyhedron with fewer surfaces grows more quickly.Although the initial polyhedron particles' roundness varies greatly, after three abrasions, all of the particles' roundness stabilize at around 0.91, which amply proves that the method described in this study can accurately simulate the abrasion process of polygonal particles.
s 2024, 17, x FOR PEER REVIEW

Sample Preparation
It is assumed that the particles in the four gr Figures 2 and 5, with the same amount of particles in of particle size between 10 mm and 40 mm.As frag into account in this study, the rigid block model in P particles.The linear model in PFC6.0 is used for the since the bonding force between particles is not tak different specimens varies because of the wide rang make the results of each group of specimens compar of sample preparation, ensuring the greater relative terial.Although there is contact sliding friction be [35], in the discrete element method (DEM), to av stress of the sample and to produce specimens with tion, the sliding friction coefficient between the part ies have shown that this approach is feasible unde The DEM parameters for simulations obtained from shown in Table 1.

Sample Preparation
It is assumed that the particles in the four groups of specimens are nine forms in Figures 2 and 5, with the same amount of particles in each form and a uniform distribution of particle size between 10 mm and 40 mm.As fragmentation of the particles is not taken into account in this study, the rigid block model in PFC6.0 [34] can be used to simulate the particles.The linear model in PFC6.0 is used for the inter-particle contact intrinsic model since the bonding force between particles is not taken into account.The initial porosity of different specimens varies because of the wide range of particle morphology.In order to make the results of each group of specimens comparable, a coefficient was used at the time of sample preparation, ensuring the greater relative density of the deposited granular material.Although there is contact sliding friction between the particles and the sidewalls [35], in the discrete element method (DEM), to avoid gradient changes in the internal stress of the sample and to produce specimens with a relatively uniform particle distribution, the sliding friction coefficient between the particles and the walls is set to zero.Studies have shown that this approach is feasible under low micropressure conditions [36].The DEM parameters for simulations obtained from the previous calibration work [37] are shown in Table 1.

Conventional Triaxial Test
As shown in Figure 8, four groups of granular aggregates are generated in a frustumshaped area of a cone with a base diameter of 300 mm and a height of 600 mm.The boundary is a shell-type element with a thickness of 2 mm and an elastic modulus of 7.8 MPa to simulate the rubber mold.The top and bottom loading plates are rigid wall elements.The interaction between the discrete particles and the lateral boundary can be realized by the discrete-continuous coupling method built into PFC6.0.A constant 100 kPa surrounding pressure σ 3 is applied to the sidewalls, and top and bottom loading plates are applied to all four groups of specimens after they reach a relative density of 0.8 in order to simulate the consolidation of the specimens.After the specimens reached stability, the confining pressure was kept constant, and the specimens were loaded by applying a relative velocity on the upper and lower loading plates, and the loading was stopped when the axial strain reached 15%.At the end of loading, an obvious uneven bulging deformation in the middle of the specimen was observed.
discrete-continuous coupling method built into PFC6.pressure 3 σ is applied to the sidewalls, and top and bo all four groups of specimens after they reach a relative the consolidation of the specimens.After the specime pressure was kept constant, and the specimens were loa on the upper and lower loading plates, and the loading reached 15%.At the end of loading, an obvious uneven of the specimen was observed.

Deviatoric Stress and Strain
In this study, a stress tensor is used to quantify the of granular materials in the shear process.The calculat

Deviatoric Stress and Strain
In this study, a stress tensor is used to quantify the macroscopic mechanical response of granular materials in the shear process.The calculation method is as follows: where σ ′ ij is the deviator part of the stress tensor σ ij , V is the volume of the specimen, f c i is the component of the contact force vector f , and d c j is the component of the contact branch vector d.Mean stress p and deviatoric stress q can be defined as: The axial strain of the specimen is calculated as ε a = ln ∆H/H 0 , where ∆H is the loading displacement of the load plate and H 0 is the initial height of the specimen.According to Ma's [31] method, the volumetric strain ε v is determined by the relative placement of the loading plate and the lateral shell nodes during the loading process.
At a confining pressure of 100 kPa, Figures 9 and 10 depict the macroscopic mechanical response of the four groups of materials.The simulation results demonstrate that the deviatoric stress-axial strain patterns are similar across all specimens and that after a brief linear increase at the start of loading, the growth rate quickly declines until the axial strain peaks at about 5%.The peak strength of the corroded granular specimens was marginally less than that of the initial polyhedron specimens, and it decreased as the degree of abrasion increased.At the same time, the dilatancy effect of particle specimens with a higher degree of angularity is more obvious, and increasing the abrasion degree of particles will inhibit the dilatancy of the material.The model test results show that this law is convincing for coarse sand materials composed of dense non-cohesive particles.
peaks at about 5%.The peak strength of the corroded g less than that of the initial polyhedron specimens, and sion increased.At the same time, the dilatancy effect o degree of angularity is more obvious, and increasing th inhibit the dilatancy of the material.The model test res ing for coarse sand materials composed of dense non-c

Strength and Deformation
Figures 11 and 12 demonstrate the effect of the d secant modulus and strength of the granular material that the mechanical properties of the granular materia own particle shape.Increasing the degree of particle er granular material and increases its Poisson's ratio.The less than that of the initial polyhedron specimens, and sion increased.At the same time, the dilatancy effect o degree of angularity is more obvious, and increasing th inhibit the dilatancy of the material.The model test res ing for coarse sand materials composed of dense non-c

Strength and Deformation
Figures 11 and 12 demonstrate the effect of the d secant modulus and strength of the granular material that the mechanical properties of the granular material own particle shape.Increasing the degree of particle er granular material and increases its Poisson's ratio.The angles of the granular material decrease with increasin

Strength and Deformation
Figures 11 and 12 demonstrate the effect of the degree of particle abrasion on the secant modulus and strength of the granular material.The results of both figures show that the mechanical properties of the granular material are significantly influenced by its own particle shape.Increasing the degree of particle erosion weakens the modulus of the granular material and increases its Poisson's ratio.The peak friction and shear expansion angles of the granular material decrease with increasing abrasion.This is due to the fact that the internal structure of the corroded granular material is more fragile and less able to resist external loads, a phenomenon that can be explained by the stress-dilatancy theory [38,39].

Coordination Number
The coordination number is an important indicato properties and stability of granular materials; it is defi contact with a given particle [40,41].Usually, researc average coordination number of the particle system.T nation number of the four groups of specimens during Figure 13.The initial coordination number of the four

Coordination Number
The coordination number is an important indicato properties and stability of granular materials; it is defi contact with a given particle [40,41].Usually, research average coordination number of the particle system.T nation number of the four groups of specimens during Figure 13.The initial coordination number of the four and the specimens showed a slight increase in the coor

Coordination Number
The coordination number is an important indicator in evaluating the microstructural properties and stability of granular materials; it is defined as the number of particles in contact with a given particle [40,41].Usually, researchers are more concerned with the average coordination number of the particle system.The changes in the average coordination number of the four groups of specimens during the loading process are shown in Figure 13.The initial coordination number of the four sets of specimens was around 6.2, and the specimens showed a slight increase in the coordination number at the beginning of loading, indicating a slight shear contraction of the specimens.Subsequently, the average coordination number of the four groups of specimens decreased rapidly, but the rate of decrease in coordination number was proportional to the degree of particle abrasion.This indicates that the ability of the particle system to resist external loading decreases as the degree of abrasion increases.When the axial strain reaches about 11%, the coordination numbers of the four groups of tests remain stable.The critical coordination number, which is thought to be the minimum necessary required to maintain the specimen's mechanical stability during the shear process, decreases slightly as the degree of abrasion increases but not noticeably.The critical coordination number for the four test groups remains constant between 3.8 and 4.1.

Internal Sliding and Rotation
The macroscopic deformation of the specimen is man relative sliding and rotation between particles.Whether rel ticles occurs depends on the degree of friction play., relative sliding of the p ft are the normal and tangential contact forces, respectively, tion coefficient.The sliding rate of the particle system is de where Nsc is the number of particles in which sliding occurs contacts.Figure 14 shows the variation in the sliding rate d of specimens.Sp increases rapidly with increasing ε a for the reaches a peak at a 1.5% ε = .It then decreases slowly.Over for the initial polyhedron group, and the average Sp during sion increases.

Internal Sliding and Rotation
The macroscopic deformation of the specimen is manifested microscopically as the relative sliding and rotation between particles.Whether relative sliding between the particles occurs depends on the degree of friction play.The friction mobility index I m = | f t |/µ f n can be used to judge whether sliding occurs between particles in contact with each other.When I m > 0.9999, relative sliding of the particles occurs, where fn and f t are the normal and tangential contact forces, respectively, and µ is the inter-particle friction coefficient.The sliding rate of the particle system is defined as S p = N sc /N c × 100%, where N sc is the number of particles in which sliding occurs and N c is the total number of contacts.Figure 14 shows the variation in the sliding rate during loading for four groups of specimens.S p increases rapidly with increasing ε a for the four groups of specimens and reaches a peak at ε a = 1.5%.It then decreases slowly.Overall, S p increases more rapidly for the initial polyhedron group, and the average S p during loading is smaller as the abrasion increases.

Internal Sliding and Rotation
The macroscopic deformation of the specimen is man relative sliding and rotation between particles.Whether rel ticles occurs depends on the degree of friction play.where Nsc is the number of particles in which sliding occurs contacts.Figure 14 shows the variation in the sliding rate d of specimens.Sp increases rapidly with increasing ε a for the reaches a peak at a 1.5% ε = . It then decreases slowly.Ove for the initial polyhedron group, and the average Sp during sion increases.Another microscopic cause of the deformation of the specimen is the rotation of the particles, and the average particle rotation w p is used here to describe the degree of particle rotation during the loading process.It is calculated as follows: where w a is the rotation angle of a single particle during loading and N r is the number of particles in the material.Figure 15 shows the average rotation angle of particles during the loading process in the numerical simulations described in Section 3. When ε a was less than 1%, the initial polyhedron specimens produced almost no increment in rotation angle, while the group with the greatest abrasion produced significant rotation of the particles from the beginning of loading.Thereafter, the incremental rate of w p in the four groups of specimens gradually accelerated during loading, and after ε a reached 8%, w p increased linearly.The rate of increase in w p increased significantly with the increase in abrasion degree, while the rate of increase in the rotation of the specimens in the initial polyhedron group was significantly lower than that of the other three groups.
, x FOR PEER REVIEW groups of specimens gradually accelerated increased linearly.The rate of increase in w abrasion degree, while the rate of increase polyhedron group was significantly lower The results from the two experimenta and rotation between particles and indicat at low shear strain is primarily driven by rel increases, the percentage of material deform increases.Specifically, in the low-roundne mation resulting from sliding is higher dur The results from the two experimental sets reveal the occurrence of mutual sliding and rotation between particles and indicate that the deformation of granular aggregates at low shear strain is primarily driven by relative sliding between particles.As shear strain increases, the percentage of material deformation attributed to particle rotation gradually increases.Specifically, in the low-roundness granular system, the percentage of deformation resulting from sliding is higher during shear deformation.This is attributed to the stronger interlocking effect among particles with lower roundness, which compels particles to slide in order to adjust the deformation of the specimen.Moreover, the particle system with lower roundness exhibits a larger average coordination number during loading, which enhances the particle resistance to rotation.These findings are consistent with the conclusions drawn by Estrada [42] and Azéma [43].

Fabric Anisotropy
(1) Anisotropy evaluation method Fabric, which refers to the macroscopic statistics of particle arrangement or contact, is the term used to describe the microstructure of granular materials [44][45][46].The fundamental factor in determining the macro mechanical properties of granular materials is fabric anisotropy, a significant micromechanical and structural index of granular materials.The anisotropy of granular materials during loading has always been a research focus.Figure 16 shows the spatial distribution of the normal contact force in the granular system of four groups of samples before loading, at the peak state, and at the end of loading.The average contact force of each group of samples is marked on the upper right corner of the corresponding three-dimensional histogram.The normal contact force distribution of each group of samples is approximately isotropic under the same initial pressure in all directions.In the process of loading, the normal contact force of each group of specimens in the vertical direction increases fastest, indicating that shear causes the anisotropy of the normal contact force.Additionally, the rate of anisotropy is influenced by the degree of particle abrasion.The histogram becomes more slender as a result of the vertical normal contact force of particles with high angularity increasing more quickly during loading.The results demonstrates that the roundness of particles has a significant impact on the anisotropy brought on by shear, and the initial polyhedron sample's average contact force is noticeably higher than that of the other samples.Based on the results in Figure 13, particles with higher angularity form tighter contact structures between particles.Furthermore, the results in Figure 15 indicate that more angular particles, because of their irregular shapes, can interlock better with each other.This interlocking effect increases the mechanical interlock force between particles, making it more difficult for the particles to slide relative to each other under external forces.Quantitative indicators can be used to describe the anisotropy of granular materials.Anisotropy is typically divided into geometric and mechanical types.Geometric anisotropy refers to the fabric's vector distribution of particle spacing, while mechanical anisotropy refers to the directivity of contact force under external load.Five fabric tensors are Quantitative indicators can be used to describe the anisotropy of granular materials.Anisotropy is typically divided into geometric and mechanical types.Geometric anisotropy refers to the fabric's vector distribution of particle spacing, while mechanical anisotropy refers to the directivity of contact force under external load.Five fabric tensors are used to describe these indicators in order to analyze them quantitatively, including the fabric tensor of contact normal Φ c ij : the fabric tensioners of normal and tangent branch vectors Ψ b n ij and Ψ b t ij due to irregular particle shape: and the fabric tensors of normal and shear contact forces Ω f n ij and Ω f t ij due to external forces: In Equations ( 11)-( 13), n i and t i are the unit vectors in the contact normal and tangential directions, and b n , b t , f n , and f t are the projections of the branch vectors and contact force vectors in the contact normal and contact shear directions.The degree of contact anisotropy can be described by the deviatoric tensor A c ij in Equation (19).
The anisotropy evaluation method of normal and tangential branch vectors can be described by the deviatoric part of the tensor A b n ij and A b t ij in Formula (15): Similarly, the anisotropy evaluation method of normal and tangential contact force anisotropy can be described by the deviatoric part of the tensor A f n ij and A f t ij in Equation ( 16): where are the deviatoric parts of the fabric tensors in Equations ( 15) and ( 16), respectively.
The scalar value A * is introduced to describe the degree of fabric anisotropy.The calculation method is as follows: where the superscript "*" indicates c, b n , b t , f n , and f t in Equations ( 11)- (17), sign() is a symbolic function, and σ ′ ij is the deviatoric part of Equation (9).When the principal direction of A * ij is consistent with the direction of the stress tensor, sign(A * ij σ ij ′ ) is a positive number, and vice versa. (

2) Evolution of anisotropy
The evolution patterns of A c , A b n , A f n , and A f t during loading are shown in Figure 17.Since there is not a distinct long and short axis for the particles in this study, the variation pattern in the A b t term is not given here because of the minimal contribution of the branch vector in the contact shear direction (about 1%).The results show that the individual anisotropy exhibits softening characteristics as loading proceeds.In each group of results, the initial polyhedron particles have the highest degree of anisotropy, and as the abrasion level increases, the anisotropies all show varying degrees of decline.It is interesting to note that the anisotropy of A c and A f n during loading closely resembles the granular system's deviatoric stress.A b n grows slowly in the loading process and does not show a significant peak, and the degree of anisotropy is far less than the other three groups.The shear force anisotropy A f t first appears as a peak value in the loading process and rapidly decreases after the peak value.(3) Stress-Force-Fabric relationship According to some studies, the shear process of cohesionless granular mater their fabric tensors have the following correspondence: Equation ( 13) establishes the stress-force-fabric relationship, which is know relationship between the granular system's fabric anisotropy and stress tensor.Co ing that the anisotropy tensor product in square brackets is a second-order micro nent, it can be omitted.Therefore, Equation ( 18) can be combined with Equatio derive:  (3) Stress-Force-Fabric relationship According to some studies, the shear process of cohesionless granular materials and their fabric tensors have the following correspondence: Equation ( 13) establishes the stress-force-fabric relationship, which is known as the relationship between the granular system's fabric anisotropy and stress tensor.Considering that the anisotropy tensor product in square brackets is a second-order micro component, it can be omitted.Therefore, Equation ( 18) can be combined with Equation ( 9) to derive: Equation (19) demonstrates that the degree of anisotropy A * can be used to describe the deviator stress ratio.
According to the measurement principle of the macro triaxial experiment, the confining pressure σ 3 is fixed at 100kPa and the normal stress σ 1 is defined as the total axial contact force of the particles acting on the loading plate divided by the area of the loading plate.The generalized average principal stress p g and deviatoric stress q g can be expressed as follows: In Equation ( 20), the deviatoric stress ratio of the sample in the conventional triaxial test is 3(σ 1 − σ 3 )/(σ 1 + 2σ 3 ).It should be noted that there are differences between the deviatoric stress ratios obtained using Equations ( 9) and (20).In contrast to Equation (20), which defines the shear stress ratio from a macro perspective, Equation ( 9) defines it from a micro contact perspective.The deviatoric stress ratio obtained by Equation ( 9) is denoted as b, while the one defined by Equation ( 20) is denoted as c for the purpose of distinction.In the three computational models based on Equations ( 9), (19), and (20), Figure 18 compares the variation in the deviatoric stress ratio throughout the loading process.The theoretical deviatoric stress ratio curve is represented in the figure as a red straight line.The deviatoric stress ratio scatter plot and the theoretical stress ratio straight line projection in the horizontal plane during the loading process of each group of tests are very similar, showing that the approximate deviatoric stress ratio calculated by Equation ( 19) and the theoretical deviatoric stress ratio of Equation ( 9) are in good agreement.The degree of anisotropy can be used to define the deviatoric stress ratio of the specimen.The large deviatoric stress ratio found in Equation ( 20) may be attributed to the specimen's additional constraint by the rubber mold.
Materials 2024, 17, x FOR PEER REVIEW 15 In Equation ( 20), the deviatoric stress ratio of the sample in the conventional tria test is σ σ σ σ − + .It should be noted that there are differences between the viatoric stress ratios obtained using Equations ( 9) and (20).In contrast to Equation which defines the shear stress ratio from a macro perspective, Equation ( 9) defines it f a micro contact perspective.The deviatoric stress ratio obtained by Equation ( 9) is den as b, while the one defined by Equation ( 20) is denoted as c for the purpose of distinc In the three computational models based on Equations ( 19), (9), and (20), Figure 18 c pares the variation in the deviatoric stress ratio throughout the loading process.The oretical deviatoric stress ratio curve is represented in the figure as a red straight line.deviatoric stress ratio scatter plot and the theoretical stress ratio straight line projectio the horizontal plane during the loading process of each group of tests are very sim showing that the approximate deviatoric stress ratio calculated by Equation ( 19) and theoretical deviatoric stress ratio of Equation ( 9) are in good agreement.The degre anisotropy can be used to define the deviatoric stress ratio of the specimen.The large viatoric stress ratio found in Equation ( 20) may be attributed to the specimen's additi constraint by the rubber mold.

Results
In this study, a modeling method for random convex polyhedron particles is firs veloped, and subsequently, a particle abrasion method is proposed, which can better re

Results
In this study, a modeling method for random convex polyhedron particles is first developed, and subsequently, a particle abrasion method is proposed, which can better reflect the abrasion process of particles in their natural state.Based on the discrete element method, conventional triaxial tests are conducted on four granular aggregates with different abrasion degrees, and the macroscopic-microscopic connection between the granular materials is better captured.The main findings of this study are summarized as follows: • A modeling method of randomly convex polyhedron particles was established.In addition, a method of generating different degrees of abrasion of the particles based on the corners and edges of the particles and the relative positions of the surfaces was established, which better reflects the abrasion process of the particles in their natural state.According to the results of the roundness index S3d, the particles have higher sphericity after abrasion.The particle abrasion function is derived from graphical and fundamental mechanical principles.However, it is necessary to couple it with the flow field to obtain a more detailed understanding of the particle abrasion mechanics.• The macroscopic mechanical behavior of cohesionless granular materials is reproduced using traditional triaxial tests, including the nonlinear stress-strain relationship and dilatancy.With an increase in the degree of particle abrasion, the stiffness and macro shear strength of the material decrease, and the bulging deformation of the material is also restrained.• The granular aggregate with a higher degree of abrasion has lower structural stability, a quicker decline in coordination number, and a lower average coordination number when achieving stability, according to microscopic analysis.At the initial loading stage, the granular material with a higher abrasion degree is more susceptible to relative rotation, while the granular aggregate with a lower abrasion degree is more susceptible to relative sliding.These micro mechanical phenomena can, to a certain extent, reveal the macro mechanical characteristics of granular materials.• The major contribution to fabric anisotropy comes from normal contact force anisotropy, which is followed by contact normal anisotropy.The results of the S-F-F relationship for the four materials demonstrate that there is a good correspondence between the macro mechanical properties and the micro fabric of the materials and that this correspondence is independent of the particle abrasion level.The abrasion level has an inverse relationship with each fabric's parameters.The additional constraint of the rubber mold will cause the internal friction angle of granular materials obtained from conventional triaxial tests to be overestimated even if the friction constraint between the particles and boundary is not taken into consideration.

Figure 1 .
Figure 1.Modeling principle of a polyhedron particle.
Random location points on a sphere Tangent plane across the point The numerical model of gravel Intersection line of planes

Figure 1 .
Figure 1.Modeling principle of a polyhedron particle.

Figure 1 .
Figure 1.Modeling principle of a polyhedron particle.

Figure 3 .
Figure 3.The application of the Catmull-Clark subdivision method in computer modeling.

Figure 6 .
Figure 6.The roundness of initial particles with different numbers o tions.

Figure 6 .
Figure 6.The roundness of initial particles with different numbers of faces after abrasion calculations.

Figure 6 .
Figure 6.The roundness of initial particles with different tions.

Figure 7 .
Figure 7. Influence of the surface number and abrasion d

Figure 7 .
Figure 7. Influence of the surface number and abrasion degree on the roundness of particles.

Figure 8 .
Figure 8. Conventional triaxial test of the DEM.

11 .
Variation in the secant modulus E 50 and Poisson's ratio v with the abrasion degree.

12 .
Variation in the peak friction angle φ p and dilatancy angle Ψ with the abrasion degree.
to judge whether sliding occurs with each other.When 0.9999 m I > , relative sliding of the ft are the normal and tangential contact forces, respectively, tion coefficient.The sliding rate of the particle system is de

Figure 15 .
Figure 15.Average rotation angle of particles du

Figure 15 .
Figure 15.Average rotation angle of particles during loading.
Materials 2024, 17, x FOR PEER REVIEW 12 of 18 ular shapes, can interlock better with each other.This interlocking effect increases the mechanical interlock force between particles, making it more difficult for the particles to slide relative to each other under external forces.

Figure 16 .
Figure 16.Column coordinate histogram of contact normal force.

Figure 16 .
Figure 16.Column coordinate histogram of contact normal force.

Figure 17 .
Figure 17.Evolution of anisotropy.(a) Contact normal vectors c A ; (b) normal branch vect

Figure 17 .
Figure 17.Evolution of anisotropy.(a) Contact normal vectors A c ; (b) normal branch vectors A b n ; (c) contact normal force vectors A f n ; and (d) contact shear force vectors A f t ;.

Figure 18 .
Figure 18.Validation in the stress-force-fabric relationship in terms of the deviatoric stress rati the 4 types of granular materials.

Figure 18 .
Figure 18.Validation in the stress-force-fabric relationship in terms of the deviatoric stress ratio for the 4 types of granular materials.