# Numerical Study of Particle Morphology Effect on the Angle of Repose for Coarse Assemblies Using DEM

Authors: Jing Chen; Rui Gao; Yangzepeng Liu
Journal: Advances in Materials Science and Engineering
 (2019)

---

## Abstract

The morphologies of coarse particles are usually irregular and play a dominant role in the mechanical behaviors of the particle assemblies. This paper quantitatively studies the effect of particle shape on the angle of repose, which is an important macroscopic parameter for ballast materials, via laboratory tests and numerical simulations by means of the discrete element method (DEM). Forty ballast particle templates and four simply created clump templates are reconstructed using an image-based method and quantified with two shape factors, sphericity and convexity. A series of simulations are conducted with the coefficient of sliding friction between particles changing from 0.2 to 0.6 at an interval of 0.1 to study its influence on various shapes of particles, and an appropriate value of sliding friction coefficient is chosen for the comparison of particle shape effect. The results show that increasing sphericity and convexity can significantly decrease the angle of repose, and the real ballast model gives a more realistic angle of repose behaviors as that of laboratory tests compared to simply created models. By analyzing the characteristics of particle motions and contacts from a microscopic perspective, the mechanism of particle shape attributed to the formation of granular aggregation is also discussed and revealed in this research.


---

## Body

## 1. Introduction

Granular materials such as powders, grains, soils, and gravels are extremely common in industrial, agricultural, and engineering fields. Therefore, thorough understandings of their mechanical properties are required for better application of these granular materials in practical. Angle of repose (both dynamic and static) is one of the important mechanical parameters for granular packings. It is usually used to indicate the ability of materials to stay self-stable under certain external conditions. In different application environments, the use of repose angle as an indicator for mechanical behaviors also varies. For example, the angle of repose is adopted as an important index in designing silos or tanks for corn stockpiles [1]; in geotechnical engineering, it is usually associated with the internal friction angle when soils or crushed rocks are under their loosest state [2]. Studies on the angle of repose for granular materials have attracted much interest among engineers and researchers over the last few decades [1, 3–10]. Despite the knowledge and understandings about the angle of repose that have been gained so far, there are still some insufficiencies that need further investigations.One major limitation about the angle of repose in the current literatures is that particle size is usually limited to less than 10 millimetre [5, 6, 9, 11–14], and research on coarse particle assemblies with large sizes is still insufficient. Compared to soils or sands, coarse large-sized matters such as crushed rocks in dam engineering or ballast used as bearing materials in the railway system tend to have irregular morphologies, and the interlocking effects between particles are much more pronounced than that of small-sized substances, leading to different behaviors and properties in typical angle of repose tests. It has been investigated that sizes of particles have an unnegligible impact on the final formation of piles for granular aggregates [3, 4]. Therefore, it is necessary to carry out studies on the angle of repose for coarse assemblies, for both theoretical and practical use.To numerically study the mechanical properties of granular materials, a majority of scholars choose the discrete element method (DEM) which was proposed by Cundall and Strack [15, 16]. In DEM, the relative movements and rotations of individual particles inside the assembly are taken into account based on Newton’s second law and specific contact models. Supplemented with the commercial programs such as particle flow code (PFC), not only the macroscopic behaviors of particle aggregates can be investigated with DEM but also some microscopic properties including particle movement and fabric anisotropy can easily be obtained. Using DEM, researchers have found that the angle of repose for granular assemblies is influenced by many factors. Zhou et al. [10, 11] stated that the angle of repose for sand pile is significantly affected by the sliding and rolling friction and particle size, while it is not sensitive to density, Poisson’s ratio, damping coefficient, and Young’s modulus. Nakashima et al. [9] have proved that the effects of gravity on the angle of repose is negligible. Li et al. [8, 17] investigated the impacts of lifting velocity and container shape on the angle of repose. However, research on the particle morphology effect on the angle of repose is limited. Generally, in the prevailing research, some simple-shape [6, 9, 11] or artificially created templates [12–14, 18–20] are usually used in order to save computational time and to improve simulation efficiency. This simplification can generate satisfying results especially for some fine particles such as sands or soils where mechanical behaviors depend more on the particle size rather than the morphologies of individual grains. For large-sized coarse aggregates on the contrary, the angularity on particle surface can result in greater interlocking effects, which is difficult to simulate with simple particle shapes. Thus, it is necessary to study the morphology effect for coarse grains in DEM simulations.In this paper, a set of irregular-shaped particle templates that resemble the real gravel shapes were created. A series of DEM simulations for the angle of repose tests were carried out for both simple-shaped and irregular-shaped particles to numerically study the particle morphology effect on the formation of piles for coarse aggregates. In addition, the underlying mechanism that contributes to the difference of macroscopic behaviors for particles with various morphologies in the angle of repose test was also explored with DEM models. Meanwhile, a laboratory experiment using the cylinder lifting test (CLT) was conducted in parallel, and the experimental results were then compared to that obtained with DEM models to verify the important effects of particle shape on the mechanical behaviors for coarse aggregates.


## 2. Lab and DEM Specification

### 2.1. Cylinder Lifting Test

To determine the angle of repose for granular assemblies, there are generally four kinds of methods, among which the cylinder lifting test [21], the slump test [12], or the hopper flow test [13] are usually used for the static angle of repose determination, while the tumbling test [22] is used for the dynamic angle of repose determination. During the whole service life of ballast, it is subjected to cyclic and impact loadings generated by fast-moving hauls. However, the static coefficient friction of ballast, which is believed to have a close relationship with its static angle of repose, is usually paid more attention by geo-researchers and engineers in the design and construction of railway tracks. So one of the three typical static repose angle tests, the cylinder lifting test (CLT), was adopted to measure the static angle of repose for ballast in this research.The ballast used in this research was made of granite and was collected from a quarry near Wuhan, China. The size distribution of the ballast is shown in Figure1. It can be seen that the sample used in the experiment meets the Chinese National Industry Standard of Railway Ballast (TB/T 2140-2008). The maximum particle size (Dmax) of ballast is 64 mm, and the coefficient of uniformity, Cu, is 1.15. D10, D30, D50, and D60 are the particle sizes below which the amounts of particles by mass are 10%, 30%, 50%, and 60%, respectively.Figure 1
Particle-size distributions in simulations and experiments.In the CLT, about 100 kilograms of ballast with sizes ranging from 10.8 mm to 64 mm was poured into an acrylic cylinder with no bottom. The size of the cylinder is 400 mm in diameter and 850 mm in height to avoid the boundary effect [21]. The cylinder was then lifted slowly to a height of 0.4 m at a velocity of 5 mm/s until the ballast was settled under gravity to form a pile. According to the research of Li et al. [17], the effect of lifting velocity cannot be ignored on the formation of piles, but it is not within the scope of this study and is thus directly drawn from the references [21]. During the lifting process, the velocity and the displacement of the cylinder is controlled by a motor shaft with strings connected to it in order not to cause fierce collision to the particle system.The lifting process and the final states of the experiment are presented in Figure2. After lifting, the pile angle was measured by using an electronic inclinometer every 45 degrees as shown in Figure 3, and the angle of repose was calculated as the average of the eight measured angles. This test was repeated 3 times, and the results are listed in Table 1. The angle of repose for ballast is determined as the average, that is, 32.5°.Figure 2
The process of the cylinder lifting test. (a) Before lifting. (b) Middle of lifting. (c) After lifting.
(a)
(b)
(c)Figure 3
Schematic diagram of the measurement of angle of repose.Table 1
The measured angles of repose in CLT and the average value for each test.
θ
1
°
θ
2
°
θ
3
°
θ
4
°
θ
5
°
θ
6
°
θ
7
°
θ
8
°
AOR
°
Test1
32.5
33.4
31.7
32.7
32.9
31.4
33.5
32.3
32.6
Test2
31.2
33.5
31.3
32.6
33.4
32.8
31.7
32.5
32.4
Test3
33.1
33.6
31.9
32.4
31.8
32.4
32.7
32.8
32.6
AOR
j
=
1
/
8
∑
i
=
1
8
θ
i
j, where θij is the ith measured angle of repose in Testj.


### 2.2. Particle Shape Reconstruction and Quantification

In DEM, scholars have tried many ways to recreate clumps that resemble the real particle shapes [23–26]. In this paper, an image-based method proposed by Gao et al. [24] was adopted for ballast reconstruction. 40 representative particles chosen randomly from the CLT were scanned by a X-ray CT system. During the process of scanning, the vertical image of ballast was scanned every 0.6 mm, and approximately 100–200 slice images were obtained for each particle. Figure 4(a) shows a typical cross-sectional image of ballast. The slice images of sample particles were then binarized (Figure 4(b)), and the contours of the grains were extracted (Figure 4(c)) with a built-in routine in MATLAB, after which the contour points in different slices were added up layer by layer to form the surface point clouds (Figure 4(d)). Then, the surface points were prepared for Delaunay triangulation, and the triangular mesh of the objects (Figure 4(e)) was stored in STL format which can be read directly in PFC (version 5.00.35). Finally, the 40 particle surface meshes were input into PFC to generate clump templates (Figure 4(f)) with a group of circumscribed balls (3D) combined together to approximate the morphology of the irregular particles, based on the method proposed by Taghavi [27]. In Taghavi’s method, the accuracy of the shape modeling is controlled by two dimensionless parameters, the maximum intersection angle ϕ which determines the smoothness of the surface and the diameter ratio of the smallest to the largest ball ρ which decides the basic form of the clump. A combination of ϕ=120° and ρ=0.3 was used for the generation of particle bank in this study.Figure 4
The process of particle shape reconstruction. (a) Sliced image. (b) Binarized image. (c) Grain contour. (d) Surface clouds. (e) Triangular mesh. (f) Clump template.
(a)
(b)
(c)
(d)
(e)
(f)In order to quantify the particle shape, there are generally three parameters, i.e., form, roundness, and roughness, and each represents one aspect of particle morphology. It has been pointed out by Barrett [28] that form or sphericity is the first-order property which controls the overall form of particle and roundness or convexity is the second-order property that reflects variations at the corners as shown in Figure 5. Roughness, which is the third-order property that relates directly to the surface texture, cannot be reproduced by Taghavi’s method; therefore, it is achieved by the calibration of friction coefficient between particle templates in this study.Figure 5
A grain contour (heavy solid line) with its three orders of component elements of form/sphericity, roundness, and surface roughness identified.In this paper, two shape factors, sphericity and convexity [29–32], were selected. It should be mentioned that sphericity and convexity are two individual parameters with sphericity mainly distinguished between a particle and a sphere, and a larger sphericity means that the overall shape of the particle is more approximate to an absolute sphere, while convexity is associated with a local convex area along the grain surface. Figure 6 gives the schematic diagram for the calculation of these two shape factors. Referring to Figure 6, sphericity [29] and convexity [32] can be expressed as equations (1) and (2), respectively:(1)sphericity=VpVcs3,(2)convexity=ApAp+Ach,where Vp and Vcs represent the volume of particle and its circumscribing sphere and Ap and Ach represent the area of particle projection and its convex hull, respectively.Figure 6
Schematic diagram for the definition of two shape factors. (a) Calculation of sphericity. (b) Calculation of convexity.
(a)
(b)From equations (1) and (2), it is obvious that sphericity is a three-dimensional parameter and convexity is a two-dimensional parameter. So the projection with the maximum area of every particle was selected as the target section for convexity evaluation.In this research, a total of five different types of models were developed, among which the real ballast model was created with 40 real ballast templates as shown by RB-1 to RB-40 in Figure7. The other four custom models were created using four simple templates as shown by S-1, S-2O, S-4, and S-2, respectively in Figure 7. It can be seen from Figure 7 that the S-1 template has the largest sphericity and convexity values of 1, and as the values become smaller, the particle shape becomes more angular and irregular.Figure 7
Particle bank (RB-X represents the real ballast clump template with the ID ofX, S-X represents the simply created template with the pebble number of X, and S-2O is the 2-ball with overlap clump template) with corresponding sphericity and convexity; S, C, and NP refer to sphericity, convexity, and the number of pebbles, respectively.


### 2.3. DEM Models

The rainfall filling method [19, 33] was used to place the particles into the cylinder of CLT DEM models so that the particles in models were in a natural state. A horizontal wall and an open-bottom cylinder wall with 400 mm in diameter and 8500 mm in height were first created. Then, a total of 3000 particles duplicated from Figure 7 were generated inside the cylinder with particle size determined via a distribution FISH program for five types of models mentioned above. Notably, 40 real ballast templates were chosen evenly for the creation of the real ballast model with a FISH subroutine developed by the authors. The size distribution of clumps in the DEM model is also plotted in Figure 1, which shows great coincidence with that of the real ballast assemblies. The particles settled down under gravity until the unbalanced force ratio was less than 1.0e − 5, and the container was leveled at a desired height to be prepared for lifting, as shown in Figure 8.Figure 8
DEM model of 4-ball particles in a bottom-opened cylinder prepared to be lifted.In a typical DEM model, granular materials such as ballast interact with its neighboring particles via various contact constitutive models, among which the linear model equipped with or without rolling resistance is one of the simplest and the most commonly used models. It has to be clear that rolling resistance coefficient is usually introduced into the contact model as a major means of energy dissipation caused by shape irregularity or surface texture roughness. Since the irregularity in particle morphology was considered by physical geometries of the five templates, and the surface roughness was achieved via the friction efficient between particles as mentioned above in Section2.2, the linear model without rolling resistance [16] was adopted in this research. Particularly, two damping-associated parameters, i.e., dp_nratio and damp, were used in this research. As explained by Cundall and Strack [16], the inclusion of dp_nratio in the DEM model provides viscous behavior for normal directions of each contact. It was assumed that ballast material has a low coefficient of restitution, and therefore, the dp_nratio was set to a value of 0.8 [19]. The utilization of damp, which is also known as local damping, helps DEM models to accelerate convergence towards stable configurations when solving the Newton equations of motion [34, 35]. According to [35], the value of local damping should be tuned down for dynamic analyses, and it was set to be 0.2 in this research.To measure the coefficient of friction between particle and walls, five particles with different shapes and sizes (Figure9(a)) were chosen to conduct a simple sliding test [13] as depicted in Figure 9(b). The tilting angle when the particle began to slide was recorded, and the tangent of this angle was regarded as the frictional coefficient between particles and wall (i.e., μw−p). The sliding test was conducted three times for each particle at three different placing directions. Since No. 4 ballast is in a flat shape, only two values can be obtained. The results are given in Table 2, and μw−p was calculated to be 0.5.Figure 9
Measurement of coefficient of friction between wall geometry and particles. (a) Five ballast particles. (b) The sliding test.
(a)
(b)Table 2
Tilting angles for five ballast particles.
Particle numbers
Tilting angles
No. 1
27.75
°
25.35
°
29.38
°
No. 2
25.74
°
23.77
°
28.97
°
No. 3
28.56
°
24.95
°
27.35
°
No. 4
26.14
°
28.56
°
No. 5
25.35
°
28.16
°
24.76
°A series of DEM simulations with different coefficients of sliding friction were conducted firstly to choose their proper value, and its influence on the final formation of piles for granular assemblies was also investigated. The impacts of other parameters can be negligible according to Nakashima et al. [9], so they were inherited directly from an existing paper [36]. The DEM parameters are listed in Table 3.Table 3
DEM parameter list.
Parameters
Value
Particle densitykg/m3
2700
Contact normal stiffnesskn (N/m)
1.0
×
10
7
Contact shear stiffnessks (N/m)
1.0
×
10
7
Contact normal stiffness of wall particleknp−w (N/m)
5.2
×
10
6
Contact shear stiffness of wall particleksp−w (N/m)
5.2
×
10
6
Wall-particle coefficient of frictionμw−p
0.5
Particle-particle coefficient of sliding frictionμs,p−p
0.2, 0.3, 0.4, 0.5, 0.6
Damping ratiodp_nratio
0.8
Damp
0.2


## 2.1. Cylinder Lifting Test

To determine the angle of repose for granular assemblies, there are generally four kinds of methods, among which the cylinder lifting test [21], the slump test [12], or the hopper flow test [13] are usually used for the static angle of repose determination, while the tumbling test [22] is used for the dynamic angle of repose determination. During the whole service life of ballast, it is subjected to cyclic and impact loadings generated by fast-moving hauls. However, the static coefficient friction of ballast, which is believed to have a close relationship with its static angle of repose, is usually paid more attention by geo-researchers and engineers in the design and construction of railway tracks. So one of the three typical static repose angle tests, the cylinder lifting test (CLT), was adopted to measure the static angle of repose for ballast in this research.The ballast used in this research was made of granite and was collected from a quarry near Wuhan, China. The size distribution of the ballast is shown in Figure1. It can be seen that the sample used in the experiment meets the Chinese National Industry Standard of Railway Ballast (TB/T 2140-2008). The maximum particle size (Dmax) of ballast is 64 mm, and the coefficient of uniformity, Cu, is 1.15. D10, D30, D50, and D60 are the particle sizes below which the amounts of particles by mass are 10%, 30%, 50%, and 60%, respectively.Figure 1
Particle-size distributions in simulations and experiments.In the CLT, about 100 kilograms of ballast with sizes ranging from 10.8 mm to 64 mm was poured into an acrylic cylinder with no bottom. The size of the cylinder is 400 mm in diameter and 850 mm in height to avoid the boundary effect [21]. The cylinder was then lifted slowly to a height of 0.4 m at a velocity of 5 mm/s until the ballast was settled under gravity to form a pile. According to the research of Li et al. [17], the effect of lifting velocity cannot be ignored on the formation of piles, but it is not within the scope of this study and is thus directly drawn from the references [21]. During the lifting process, the velocity and the displacement of the cylinder is controlled by a motor shaft with strings connected to it in order not to cause fierce collision to the particle system.The lifting process and the final states of the experiment are presented in Figure2. After lifting, the pile angle was measured by using an electronic inclinometer every 45 degrees as shown in Figure 3, and the angle of repose was calculated as the average of the eight measured angles. This test was repeated 3 times, and the results are listed in Table 1. The angle of repose for ballast is determined as the average, that is, 32.5°.Figure 2
The process of the cylinder lifting test. (a) Before lifting. (b) Middle of lifting. (c) After lifting.
(a)
(b)
(c)Figure 3
Schematic diagram of the measurement of angle of repose.Table 1
The measured angles of repose in CLT and the average value for each test.
θ
1
°
θ
2
°
θ
3
°
θ
4
°
θ
5
°
θ
6
°
θ
7
°
θ
8
°
AOR
°
Test1
32.5
33.4
31.7
32.7
32.9
31.4
33.5
32.3
32.6
Test2
31.2
33.5
31.3
32.6
33.4
32.8
31.7
32.5
32.4
Test3
33.1
33.6
31.9
32.4
31.8
32.4
32.7
32.8
32.6
AOR
j
=
1
/
8
∑
i
=
1
8
θ
i
j, where θij is the ith measured angle of repose in Testj.


## 2.2. Particle Shape Reconstruction and Quantification

In DEM, scholars have tried many ways to recreate clumps that resemble the real particle shapes [23–26]. In this paper, an image-based method proposed by Gao et al. [24] was adopted for ballast reconstruction. 40 representative particles chosen randomly from the CLT were scanned by a X-ray CT system. During the process of scanning, the vertical image of ballast was scanned every 0.6 mm, and approximately 100–200 slice images were obtained for each particle. Figure 4(a) shows a typical cross-sectional image of ballast. The slice images of sample particles were then binarized (Figure 4(b)), and the contours of the grains were extracted (Figure 4(c)) with a built-in routine in MATLAB, after which the contour points in different slices were added up layer by layer to form the surface point clouds (Figure 4(d)). Then, the surface points were prepared for Delaunay triangulation, and the triangular mesh of the objects (Figure 4(e)) was stored in STL format which can be read directly in PFC (version 5.00.35). Finally, the 40 particle surface meshes were input into PFC to generate clump templates (Figure 4(f)) with a group of circumscribed balls (3D) combined together to approximate the morphology of the irregular particles, based on the method proposed by Taghavi [27]. In Taghavi’s method, the accuracy of the shape modeling is controlled by two dimensionless parameters, the maximum intersection angle ϕ which determines the smoothness of the surface and the diameter ratio of the smallest to the largest ball ρ which decides the basic form of the clump. A combination of ϕ=120° and ρ=0.3 was used for the generation of particle bank in this study.Figure 4
The process of particle shape reconstruction. (a) Sliced image. (b) Binarized image. (c) Grain contour. (d) Surface clouds. (e) Triangular mesh. (f) Clump template.
(a)
(b)
(c)
(d)
(e)
(f)In order to quantify the particle shape, there are generally three parameters, i.e., form, roundness, and roughness, and each represents one aspect of particle morphology. It has been pointed out by Barrett [28] that form or sphericity is the first-order property which controls the overall form of particle and roundness or convexity is the second-order property that reflects variations at the corners as shown in Figure 5. Roughness, which is the third-order property that relates directly to the surface texture, cannot be reproduced by Taghavi’s method; therefore, it is achieved by the calibration of friction coefficient between particle templates in this study.Figure 5
A grain contour (heavy solid line) with its three orders of component elements of form/sphericity, roundness, and surface roughness identified.In this paper, two shape factors, sphericity and convexity [29–32], were selected. It should be mentioned that sphericity and convexity are two individual parameters with sphericity mainly distinguished between a particle and a sphere, and a larger sphericity means that the overall shape of the particle is more approximate to an absolute sphere, while convexity is associated with a local convex area along the grain surface. Figure 6 gives the schematic diagram for the calculation of these two shape factors. Referring to Figure 6, sphericity [29] and convexity [32] can be expressed as equations (1) and (2), respectively:(1)sphericity=VpVcs3,(2)convexity=ApAp+Ach,where Vp and Vcs represent the volume of particle and its circumscribing sphere and Ap and Ach represent the area of particle projection and its convex hull, respectively.Figure 6
Schematic diagram for the definition of two shape factors. (a) Calculation of sphericity. (b) Calculation of convexity.
(a)
(b)From equations (1) and (2), it is obvious that sphericity is a three-dimensional parameter and convexity is a two-dimensional parameter. So the projection with the maximum area of every particle was selected as the target section for convexity evaluation.In this research, a total of five different types of models were developed, among which the real ballast model was created with 40 real ballast templates as shown by RB-1 to RB-40 in Figure7. The other four custom models were created using four simple templates as shown by S-1, S-2O, S-4, and S-2, respectively in Figure 7. It can be seen from Figure 7 that the S-1 template has the largest sphericity and convexity values of 1, and as the values become smaller, the particle shape becomes more angular and irregular.Figure 7
Particle bank (RB-X represents the real ballast clump template with the ID ofX, S-X represents the simply created template with the pebble number of X, and S-2O is the 2-ball with overlap clump template) with corresponding sphericity and convexity; S, C, and NP refer to sphericity, convexity, and the number of pebbles, respectively.


## 2.3. DEM Models

The rainfall filling method [19, 33] was used to place the particles into the cylinder of CLT DEM models so that the particles in models were in a natural state. A horizontal wall and an open-bottom cylinder wall with 400 mm in diameter and 8500 mm in height were first created. Then, a total of 3000 particles duplicated from Figure 7 were generated inside the cylinder with particle size determined via a distribution FISH program for five types of models mentioned above. Notably, 40 real ballast templates were chosen evenly for the creation of the real ballast model with a FISH subroutine developed by the authors. The size distribution of clumps in the DEM model is also plotted in Figure 1, which shows great coincidence with that of the real ballast assemblies. The particles settled down under gravity until the unbalanced force ratio was less than 1.0e − 5, and the container was leveled at a desired height to be prepared for lifting, as shown in Figure 8.Figure 8
DEM model of 4-ball particles in a bottom-opened cylinder prepared to be lifted.In a typical DEM model, granular materials such as ballast interact with its neighboring particles via various contact constitutive models, among which the linear model equipped with or without rolling resistance is one of the simplest and the most commonly used models. It has to be clear that rolling resistance coefficient is usually introduced into the contact model as a major means of energy dissipation caused by shape irregularity or surface texture roughness. Since the irregularity in particle morphology was considered by physical geometries of the five templates, and the surface roughness was achieved via the friction efficient between particles as mentioned above in Section2.2, the linear model without rolling resistance [16] was adopted in this research. Particularly, two damping-associated parameters, i.e., dp_nratio and damp, were used in this research. As explained by Cundall and Strack [16], the inclusion of dp_nratio in the DEM model provides viscous behavior for normal directions of each contact. It was assumed that ballast material has a low coefficient of restitution, and therefore, the dp_nratio was set to a value of 0.8 [19]. The utilization of damp, which is also known as local damping, helps DEM models to accelerate convergence towards stable configurations when solving the Newton equations of motion [34, 35]. According to [35], the value of local damping should be tuned down for dynamic analyses, and it was set to be 0.2 in this research.To measure the coefficient of friction between particle and walls, five particles with different shapes and sizes (Figure9(a)) were chosen to conduct a simple sliding test [13] as depicted in Figure 9(b). The tilting angle when the particle began to slide was recorded, and the tangent of this angle was regarded as the frictional coefficient between particles and wall (i.e., μw−p). The sliding test was conducted three times for each particle at three different placing directions. Since No. 4 ballast is in a flat shape, only two values can be obtained. The results are given in Table 2, and μw−p was calculated to be 0.5.Figure 9
Measurement of coefficient of friction between wall geometry and particles. (a) Five ballast particles. (b) The sliding test.
(a)
(b)Table 2
Tilting angles for five ballast particles.
Particle numbers
Tilting angles
No. 1
27.75
°
25.35
°
29.38
°
No. 2
25.74
°
23.77
°
28.97
°
No. 3
28.56
°
24.95
°
27.35
°
No. 4
26.14
°
28.56
°
No. 5
25.35
°
28.16
°
24.76
°A series of DEM simulations with different coefficients of sliding friction were conducted firstly to choose their proper value, and its influence on the final formation of piles for granular assemblies was also investigated. The impacts of other parameters can be negligible according to Nakashima et al. [9], so they were inherited directly from an existing paper [36]. The DEM parameters are listed in Table 3.Table 3
DEM parameter list.
Parameters
Value
Particle densitykg/m3
2700
Contact normal stiffnesskn (N/m)
1.0
×
10
7
Contact shear stiffnessks (N/m)
1.0
×
10
7
Contact normal stiffness of wall particleknp−w (N/m)
5.2
×
10
6
Contact shear stiffness of wall particleksp−w (N/m)
5.2
×
10
6
Wall-particle coefficient of frictionμw−p
0.5
Particle-particle coefficient of sliding frictionμs,p−p
0.2, 0.3, 0.4, 0.5, 0.6
Damping ratiodp_nratio
0.8
Damp
0.2


## 3. Results and Discussion

The impacts of the coefficient of sliding friction between particles and particle morphology on the angle of repose are analyzed in this section. By comparing the DEM modeling results for different particle shapes and that of the experimental tests, the mechanism that contributes to the difference of repose angle between various granular assemblies concerning the particle shapes is also discussed in this section.


### 3.1. Influence of Friction Coefficient between Particles on the Angle of Repose

To qualify the influence of sliding friction coefficientμs,pp on the angle of repose, simulations were conducted for five types of material models with different μs,pp ranging from 0.2 to 0.6 at an interval of 0.1. The final states for 2-ball with overlap models under different μs,pp were picked and depicted in Figure 10 as examples. As shown in Figure 10, increasing sliding friction coefficient μs,pp can increase the angle of repose for particle assemblies. Figure 11 presents the results of angle of repose for different models. The smallest and the greatest angle of repose were observed in the 1-ball model with μs,pp at 0.4 and the 2-ball without overlap model with μs,pp at 0.6, respectively. This indicates that the angle of repose for ballast not only associates with the sliding friction coefficient between particles but also depends on the shapes of grains. Here, the angles of repose for 1-ball models with μs,pp at 0.2 and 0.3 were too small to measure as all the particles were extremely scattered on the base plane and no distinct pile was formed, so the results for the two simulations were not plotted in Figure 11.Figure 10
Particle piles for 2-ball with overlap model under differentμs,pp. (a)μs,pp=0.2; (b)μs,pp=0.3; (c)μs,pp=0.4; (d)μs,pp=0.5; (e)μs,pp=0.6.
(a)
(b)
(c)
(d)
(e)Figure 11
The effect of sliding coefficient of friction on the angle of repose.In DEM, scholars tend to use spheres or clumps with fewer amount of pebbles in order to shorten the simulation time and to improve the efficiency and lots of satisfying results have been obtained in the existing researches [9, 10, 13, 19, 37]. However, this simplicity in particle shape requires stricter parameter calibration for the reliability of simulation results and might lead to inaccurate results in some situations. As shown in Figure 11, it was impossible for 1-ball or 2-ball with overlap models to reach the expected value with only one parameter, i.e., μs,pp, being used in this study. This is due to the reason that the effect of interlocking is much more dominant on particle movements (both translational and rotational) for coarse grains such as ballast or crushed rocks. Therefore, it can be indicated that rolling friction has to be equipped as a way to account for asperities or nonspherical shape factors of particles in DEM, when simplified particle shapes are utilized. For the other three models, although a certain value of μs,pp can always be determined to obtain the same angle of repose as CLT, they reacted differently to the variation of μs,pp. The 2-ball without overlap models have a wider range of repose angle when μs,pp changes from 0.2 to 0.6 compared to the real ballast models and the 4-ball models, which means it is much easily affected by μs,pp variations. Notably, even if the clumps in models have same numbers of pebbles, such as 2-ball with overlap models and 2-ball without overlap models, the relationships between μs,pp and the angles of repose for them are different. It can be found that 2-ball without overlap model produces much larger repose angles than 2-ball with overlap models due to aspect ratio provides geometric resistance to particle rotation. In general, particle shape plays an important role in the angle of repose simulations, and their specific impacts on it have to be fully studied.Whenμs,pp is at a relatively small value, the input value of μs,pp can have a significant impact on the final formation of grain piles and larger μs,pp produces larger angles of repose for all five models. These parameter-dependent characteristics of the angle of repose test make it a prevailing way for parameter calibration in current numerical studies. However, when μs,pp exceeds or equals 0.5, the growth in the angle of repose becomes flat with μs,pp increasing. The phenomenon indicates that the angle of repose is mainly determined by the particle shape at this stage. So 0.5 is adopted as the value of μs,pp to study the influence of particle shape in the following research.


### 3.2. Influence of Particle Shape on the Angle of Repose

Five DEM models withμs,pp of 0.5 and the other parameters kept to be the same as listed in Table 2 were prepared in PFC, and the final states of the granular piles for five particle templates are shown in Figure 12. It is observed from Figure 12 that sphere-shaped particles, i.e., 1-ball particles, are more likely to roll down the pile until the friction between particles and the base plane dissipated all the abundant energy stored in the system, and this hinders their forward movement. By contrast, the number of particles rolling down the piles decreases in the granular assemblies of the other four types of models, i.e., the 2-ball with overlap model, the real ballast model, the 4-ball model, and the 2-ball without overlap model, which demonstrates that more particles would keep stable on the surface of piles, resulting in larger slope angles than that of spherical particles. The largest angle of repose is obtained in the model with a particle shape of 2-ball without overlap, followed by a slight decrease in 4-ball and real ballast models. The angles of repose for 1-ball and 2-ball with overlap models are relatively smaller, and the 1-ball granular assembly has the smallest value.Figure 12
The final states of the granular piles for five models. (a) 1-ball model; (b) 2-ball with overlap model; (c) real ballast model; (d) 4-ball model; (e) 2-ball without overlap model.
(a)
(b)
(c)
(d)
(e)As presented in Figure7, the sphericity and convexity values for 40 ballast particles are evaluated. The means at 0.77 and 0.936 were calculated as the representative values of real ballast sphericity and convexity, respectively, in the following research. Figure 13 shows the relationship between the angle of repose and sphericity of five different particle templates. Obviously, as particle sphericity increases, the angle of repose decreases. This trend is in agreement with what was observed by Jo et al. [38]. When it comes to particle convexity, it can be seen from Figure 14 that increasing convexity can significantly decrease the angle of repose for granular assemblies. Particularly, it is worth noting that the negative correlations between particle shape parameters and the angle of repose begins to disappear and stable values of repose angle can be attained when particle sphericity and convexity are lower than 0.77 and 0.96, respectively, as shown by the fitting curves in Figures 13 and 14.Figure 13
The effect of particle sphericity on the angle of repose.Figure 14
The effect of particle convexity on the angle of repose.In this research, the results of angle of repose obtained by various DEM models and laboratory tests are also compared to investigate the importance of particle morphology in DEM simulations. It can be seen from Figures13 and 14 that models using the real ballast and the 4-ball templates are in best accordance with the experiments. Moreover, the results are more deviated with that obtained in CLT if the shape factors of particle templates are more discrepant with that of the real ballast templates. The results shown in Figures 13 and 14 imply that particle sphericity and convexity could also be taken into consideration when conducting parameter calibration in DEM simulations for coarse grains, besides those regular parameters such as friction coefficient and damping.


### 3.3. Microscopic Behaviors of Coarse Assemblies

To further study the mechanism of particle morphology influence on the angle of repose for coarse assemblies, some microscopic behaviors, including characteristics of particle movement and contacts for all five models, are explored in this section.


#### 3.3.1. Characteristics of Particle Movement

In this research, the translational and rotational velocities of particles located at three different heights (i.e., the upper, the middle, and the bottom parts of the cylinder) over the first eight seconds were monitored. At each height, 4 particles located in the circumference of the assemblies were picked as they are more likely to roll along the surface of the pile compared to those located in the middle part of the assemblies. Their velocities were averaged and fitted as presented by the movement curves for five different material shapes in Figures15–17, and their characteristics can be summarized as follows:(1)
For translational movement depicted in Figure15, all five types of particles initially increase to its peak velocity and then gradually decrease nearly to zero. However, the stabilization period is much shorter for nonspherical models with only about 2 seconds after the first 8 seconds than that of 1-ball model with more than 40 seconds which has been discarded in the graph for a better comparison. Notably, for particles located at the middle and top parts, the time of peak velocity is relatively earlier for 1-ball models than that for other shapes of models especially those with small particle shape factors (i.e., real ballast, 4-ball, and 2-ball particles). The difference is due to the lack of sufficient interlocking and hence the resistance between particles with large sphericity and convexity. These features are more pronounced on rotational movement for five models as shown in Figure 16.(2)
The average translational velocity curves in Figure17(a) show that the translational velocity has a close relationship with particle shape. Particles located at the upper part of cylinder have a greater average translational velocity, and this value decreases as the height decreased for all five types of models. For the four simply created clump templates, the average translational velocity becomes smaller with particle shape factors increasing. At the bottom of cylinder, the translational velocity of 1-ball templates is significantly higher than that of other three custom clump templates, which leads to the severe scatter of particles on the base plane as observed in the final simulation states. Figure 17(a) also shows that the average translational velocities for real ballast are relatively lower compared with those artificial ones in all three heights, which is in accordance with what was observed in the research by Zeng et al. [21].(3)
When it comes to the average rotational velocity as shown in Figure17(b), it is obvious that, at all three levels, the spherical or the nearly spherical models exhibit a much more distinctive rotation. They have a relatively larger average rotational velocity than those less spherical models. With a decrease in particle shape factors, the average rotational velocity drops significantly. More importantly, for 1-ball models and 2-ball with overlap models, the rotational velocity for particles located at the bottom of the cylinder is much greater than that located at the middle and the upper parts. By contrast, for particles with less shape factors, the rotational velocity gradually decreases from the upper parts to the bottom. This observation explains why particles with small sphericity values can form excellent grain piles on the base plane while it was extremely hard for spherical or the 2-ball with overlap clump templates to be stable when released from the cylinder. During the lifting process, spherical or nearly spherical particles from the bottom are much easier to roll away with a higher rotational velocity when suffered from horizontal push by upper particles, resulting in the local collapse of particle assemblies. However, in the cases of nonspherical particles, the bottom particles are difficult to rotate and they can build up a strong base for those particles rolling from upper parts, leading to a larger angle of repose compared to 1-ball and 2-ball with overlap models.Figure 15
The translational velocities for particles located at different parts of the cylinder. (a) Bottom part. (b) Middle part. (c) Upper part.
(a)
(b)
(c)Figure 16
The rotational velocities for particles located at different parts of the cylinder. (a) Bottom part. (b) Middle part. (c) Upper part.
(a)
(b)
(c)Figure 17
The average velocities for particles located at different heights of the cylinder. (a) Translational movement. (b) Rotational movement.
(a)
(b)The histories of kinetic energy in the system against simulation time for five models are presented in Figure18. It is obvious from Figure 18(a) that the kinetic energy for all five types of particle templates went through an initial fluctuation with a local maximum at the very beginning of the simulation. This is due to the onset of the movement and the fact that the lifting cylinder breaks the balance that has been maintained in the preparation phase. During this initial stage, the kinetic energy of the 1-ball model is relatively higher than other models and no significant difference of the kinetic energy is observed between the remaining 4 types of models. This shows that the 1-ball model is much more susceptible to external disturbances for the lack of interlocking and rolling resistance between particles. The kinetic energy climbs up to its peak value when all the particles in the cylinder are released out 3 to 4 seconds after the lifting starts, followed by a long period of rebalance and a decrease nearly close to zero. For four custom clump templates, the peak kinetic energy becomes greater with larger shape factors. It should be mentioned that the 1-ball model needs a much longer time, up to about 50 seconds, to reachieve the same balance state as other models, which indicates that the abundant energy stored in the system is hard to dissipate for the lack of rolling resistance, and the last 40 seconds were abandoned for a better comparison. Compared to the other four manually created template models, the real ballast model has a relatively smaller kinetic energy as shown in Figure 18(b) . This can be explained as the result of shape diversity in the real ballast systems. The similarity in particle shapes tends to result in larger void ratios for artificial assemblies even if having the same size distributions as that of real ballast systems. By contrast, the voids of real ballast can be better filled with clumps with different shapes and hence lead to a denser assemblies, in which particles movement can be restricted.Figure 18
The kinetic energy variation of particle assemblies for five models. (a) 0–0.062 s; (b) 0–10 s.In this research, we set the particle densities of five templates (i.e. 2700kg/m3) same rather than the bulk densities of five specimens. However, the variation in bulk densities may change the dynamics of particle flow and hence the angle of repose. The impacts of bulk density on the angle of repose for coarse assemblies will be further investigated in the future.


#### 3.3.2. Characteristics of Particle Contacts

The concave area along the outlines of the clump templates increases the contact points between particles, allowing for an additional generation of resistance torque; thus the movement of particles is hindered, and a much larger angle of repose is generated. For spherical templates, i.e., 1-ball clump, there is only one contact point between two particles. The contact normal is collinear with the branch vector, which is the vector that connects the centroids of two particles. Thus, the normal contact will not impart any additional resistance torque on two bodies. In the cases of nonspherical bodies, more than one contact point will be observed and these extra contact forces will resist the free rotation of particles to some extent. Additionally, in a polydisperse system, small size particles tend to rest in the concave area of large size particles, which leads to the compactness of assemblies and hence a large angle of repose will be observed.To prove the idea, the variation of coordinate number, which is defined as the average number of contacts around a single particle, of the five models were monitored. The coordinate number (CN) is a parameter that reflects the contact density of assemblies at the scale of particles. According to the definition used by O’Sullivan [34], the CN can be calculated as follows:(3)CN=2×NcNp,where Nc is the total number of contacts and Np is the number of particles of the assemblies. Here, the divisor is multiplied by 2 as each contact is shared between two bodies. During the lifting process, the CN was recorded every second for five models, and the results are shown in Figure 19. For all five models, the CN decreases typically during first 4 seconds of simulation and then gradually increases until a stable value is reached. This is because that the number of contacts between particles drops as the collapse of the initial balance state, and more contacts are recreated inside the particle assemblies in order to reach another balance state at the final stage. It is particularly noteworthy that as particle convexity increases, the CN significantly decreases, which means that the convex area along the surface of particles increases the number of contacts of the assemblies and thus the angle of repose is affected as a result. However, the variation in CN for the real ballast models, the 4-ball models, and the 2-ball without overlap models is insignificant during the whole process. This is due to the reason that the CN is influenced not only by particle shapes but also by their sizes. There is an ultimate limit of contact intensity for individual particle when its size is determined as explained by O’Sullivan [34], and the 3 nonspherical models in this research just reach this critical state.Figure 19
The effect of particle convexity on the coordinate number.


## 3.1. Influence of Friction Coefficient between Particles on the Angle of Repose

To qualify the influence of sliding friction coefficientμs,pp on the angle of repose, simulations were conducted for five types of material models with different μs,pp ranging from 0.2 to 0.6 at an interval of 0.1. The final states for 2-ball with overlap models under different μs,pp were picked and depicted in Figure 10 as examples. As shown in Figure 10, increasing sliding friction coefficient μs,pp can increase the angle of repose for particle assemblies. Figure 11 presents the results of angle of repose for different models. The smallest and the greatest angle of repose were observed in the 1-ball model with μs,pp at 0.4 and the 2-ball without overlap model with μs,pp at 0.6, respectively. This indicates that the angle of repose for ballast not only associates with the sliding friction coefficient between particles but also depends on the shapes of grains. Here, the angles of repose for 1-ball models with μs,pp at 0.2 and 0.3 were too small to measure as all the particles were extremely scattered on the base plane and no distinct pile was formed, so the results for the two simulations were not plotted in Figure 11.Figure 10
Particle piles for 2-ball with overlap model under differentμs,pp. (a)μs,pp=0.2; (b)μs,pp=0.3; (c)μs,pp=0.4; (d)μs,pp=0.5; (e)μs,pp=0.6.
(a)
(b)
(c)
(d)
(e)Figure 11
The effect of sliding coefficient of friction on the angle of repose.In DEM, scholars tend to use spheres or clumps with fewer amount of pebbles in order to shorten the simulation time and to improve the efficiency and lots of satisfying results have been obtained in the existing researches [9, 10, 13, 19, 37]. However, this simplicity in particle shape requires stricter parameter calibration for the reliability of simulation results and might lead to inaccurate results in some situations. As shown in Figure 11, it was impossible for 1-ball or 2-ball with overlap models to reach the expected value with only one parameter, i.e., μs,pp, being used in this study. This is due to the reason that the effect of interlocking is much more dominant on particle movements (both translational and rotational) for coarse grains such as ballast or crushed rocks. Therefore, it can be indicated that rolling friction has to be equipped as a way to account for asperities or nonspherical shape factors of particles in DEM, when simplified particle shapes are utilized. For the other three models, although a certain value of μs,pp can always be determined to obtain the same angle of repose as CLT, they reacted differently to the variation of μs,pp. The 2-ball without overlap models have a wider range of repose angle when μs,pp changes from 0.2 to 0.6 compared to the real ballast models and the 4-ball models, which means it is much easily affected by μs,pp variations. Notably, even if the clumps in models have same numbers of pebbles, such as 2-ball with overlap models and 2-ball without overlap models, the relationships between μs,pp and the angles of repose for them are different. It can be found that 2-ball without overlap model produces much larger repose angles than 2-ball with overlap models due to aspect ratio provides geometric resistance to particle rotation. In general, particle shape plays an important role in the angle of repose simulations, and their specific impacts on it have to be fully studied.Whenμs,pp is at a relatively small value, the input value of μs,pp can have a significant impact on the final formation of grain piles and larger μs,pp produces larger angles of repose for all five models. These parameter-dependent characteristics of the angle of repose test make it a prevailing way for parameter calibration in current numerical studies. However, when μs,pp exceeds or equals 0.5, the growth in the angle of repose becomes flat with μs,pp increasing. The phenomenon indicates that the angle of repose is mainly determined by the particle shape at this stage. So 0.5 is adopted as the value of μs,pp to study the influence of particle shape in the following research.


## 3.2. Influence of Particle Shape on the Angle of Repose

Five DEM models withμs,pp of 0.5 and the other parameters kept to be the same as listed in Table 2 were prepared in PFC, and the final states of the granular piles for five particle templates are shown in Figure 12. It is observed from Figure 12 that sphere-shaped particles, i.e., 1-ball particles, are more likely to roll down the pile until the friction between particles and the base plane dissipated all the abundant energy stored in the system, and this hinders their forward movement. By contrast, the number of particles rolling down the piles decreases in the granular assemblies of the other four types of models, i.e., the 2-ball with overlap model, the real ballast model, the 4-ball model, and the 2-ball without overlap model, which demonstrates that more particles would keep stable on the surface of piles, resulting in larger slope angles than that of spherical particles. The largest angle of repose is obtained in the model with a particle shape of 2-ball without overlap, followed by a slight decrease in 4-ball and real ballast models. The angles of repose for 1-ball and 2-ball with overlap models are relatively smaller, and the 1-ball granular assembly has the smallest value.Figure 12
The final states of the granular piles for five models. (a) 1-ball model; (b) 2-ball with overlap model; (c) real ballast model; (d) 4-ball model; (e) 2-ball without overlap model.
(a)
(b)
(c)
(d)
(e)As presented in Figure7, the sphericity and convexity values for 40 ballast particles are evaluated. The means at 0.77 and 0.936 were calculated as the representative values of real ballast sphericity and convexity, respectively, in the following research. Figure 13 shows the relationship between the angle of repose and sphericity of five different particle templates. Obviously, as particle sphericity increases, the angle of repose decreases. This trend is in agreement with what was observed by Jo et al. [38]. When it comes to particle convexity, it can be seen from Figure 14 that increasing convexity can significantly decrease the angle of repose for granular assemblies. Particularly, it is worth noting that the negative correlations between particle shape parameters and the angle of repose begins to disappear and stable values of repose angle can be attained when particle sphericity and convexity are lower than 0.77 and 0.96, respectively, as shown by the fitting curves in Figures 13 and 14.Figure 13
The effect of particle sphericity on the angle of repose.Figure 14
The effect of particle convexity on the angle of repose.In this research, the results of angle of repose obtained by various DEM models and laboratory tests are also compared to investigate the importance of particle morphology in DEM simulations. It can be seen from Figures13 and 14 that models using the real ballast and the 4-ball templates are in best accordance with the experiments. Moreover, the results are more deviated with that obtained in CLT if the shape factors of particle templates are more discrepant with that of the real ballast templates. The results shown in Figures 13 and 14 imply that particle sphericity and convexity could also be taken into consideration when conducting parameter calibration in DEM simulations for coarse grains, besides those regular parameters such as friction coefficient and damping.


## 3.3. Microscopic Behaviors of Coarse Assemblies

To further study the mechanism of particle morphology influence on the angle of repose for coarse assemblies, some microscopic behaviors, including characteristics of particle movement and contacts for all five models, are explored in this section.


### 3.3.1. Characteristics of Particle Movement

In this research, the translational and rotational velocities of particles located at three different heights (i.e., the upper, the middle, and the bottom parts of the cylinder) over the first eight seconds were monitored. At each height, 4 particles located in the circumference of the assemblies were picked as they are more likely to roll along the surface of the pile compared to those located in the middle part of the assemblies. Their velocities were averaged and fitted as presented by the movement curves for five different material shapes in Figures15–17, and their characteristics can be summarized as follows:(1)
For translational movement depicted in Figure15, all five types of particles initially increase to its peak velocity and then gradually decrease nearly to zero. However, the stabilization period is much shorter for nonspherical models with only about 2 seconds after the first 8 seconds than that of 1-ball model with more than 40 seconds which has been discarded in the graph for a better comparison. Notably, for particles located at the middle and top parts, the time of peak velocity is relatively earlier for 1-ball models than that for other shapes of models especially those with small particle shape factors (i.e., real ballast, 4-ball, and 2-ball particles). The difference is due to the lack of sufficient interlocking and hence the resistance between particles with large sphericity and convexity. These features are more pronounced on rotational movement for five models as shown in Figure 16.(2)
The average translational velocity curves in Figure17(a) show that the translational velocity has a close relationship with particle shape. Particles located at the upper part of cylinder have a greater average translational velocity, and this value decreases as the height decreased for all five types of models. For the four simply created clump templates, the average translational velocity becomes smaller with particle shape factors increasing. At the bottom of cylinder, the translational velocity of 1-ball templates is significantly higher than that of other three custom clump templates, which leads to the severe scatter of particles on the base plane as observed in the final simulation states. Figure 17(a) also shows that the average translational velocities for real ballast are relatively lower compared with those artificial ones in all three heights, which is in accordance with what was observed in the research by Zeng et al. [21].(3)
When it comes to the average rotational velocity as shown in Figure17(b), it is obvious that, at all three levels, the spherical or the nearly spherical models exhibit a much more distinctive rotation. They have a relatively larger average rotational velocity than those less spherical models. With a decrease in particle shape factors, the average rotational velocity drops significantly. More importantly, for 1-ball models and 2-ball with overlap models, the rotational velocity for particles located at the bottom of the cylinder is much greater than that located at the middle and the upper parts. By contrast, for particles with less shape factors, the rotational velocity gradually decreases from the upper parts to the bottom. This observation explains why particles with small sphericity values can form excellent grain piles on the base plane while it was extremely hard for spherical or the 2-ball with overlap clump templates to be stable when released from the cylinder. During the lifting process, spherical or nearly spherical particles from the bottom are much easier to roll away with a higher rotational velocity when suffered from horizontal push by upper particles, resulting in the local collapse of particle assemblies. However, in the cases of nonspherical particles, the bottom particles are difficult to rotate and they can build up a strong base for those particles rolling from upper parts, leading to a larger angle of repose compared to 1-ball and 2-ball with overlap models.Figure 15
The translational velocities for particles located at different parts of the cylinder. (a) Bottom part. (b) Middle part. (c) Upper part.
(a)
(b)
(c)Figure 16
The rotational velocities for particles located at different parts of the cylinder. (a) Bottom part. (b) Middle part. (c) Upper part.
(a)
(b)
(c)Figure 17
The average velocities for particles located at different heights of the cylinder. (a) Translational movement. (b) Rotational movement.
(a)
(b)The histories of kinetic energy in the system against simulation time for five models are presented in Figure18. It is obvious from Figure 18(a) that the kinetic energy for all five types of particle templates went through an initial fluctuation with a local maximum at the very beginning of the simulation. This is due to the onset of the movement and the fact that the lifting cylinder breaks the balance that has been maintained in the preparation phase. During this initial stage, the kinetic energy of the 1-ball model is relatively higher than other models and no significant difference of the kinetic energy is observed between the remaining 4 types of models. This shows that the 1-ball model is much more susceptible to external disturbances for the lack of interlocking and rolling resistance between particles. The kinetic energy climbs up to its peak value when all the particles in the cylinder are released out 3 to 4 seconds after the lifting starts, followed by a long period of rebalance and a decrease nearly close to zero. For four custom clump templates, the peak kinetic energy becomes greater with larger shape factors. It should be mentioned that the 1-ball model needs a much longer time, up to about 50 seconds, to reachieve the same balance state as other models, which indicates that the abundant energy stored in the system is hard to dissipate for the lack of rolling resistance, and the last 40 seconds were abandoned for a better comparison. Compared to the other four manually created template models, the real ballast model has a relatively smaller kinetic energy as shown in Figure 18(b) . This can be explained as the result of shape diversity in the real ballast systems. The similarity in particle shapes tends to result in larger void ratios for artificial assemblies even if having the same size distributions as that of real ballast systems. By contrast, the voids of real ballast can be better filled with clumps with different shapes and hence lead to a denser assemblies, in which particles movement can be restricted.Figure 18
The kinetic energy variation of particle assemblies for five models. (a) 0–0.062 s; (b) 0–10 s.In this research, we set the particle densities of five templates (i.e. 2700kg/m3) same rather than the bulk densities of five specimens. However, the variation in bulk densities may change the dynamics of particle flow and hence the angle of repose. The impacts of bulk density on the angle of repose for coarse assemblies will be further investigated in the future.


### 3.3.2. Characteristics of Particle Contacts

The concave area along the outlines of the clump templates increases the contact points between particles, allowing for an additional generation of resistance torque; thus the movement of particles is hindered, and a much larger angle of repose is generated. For spherical templates, i.e., 1-ball clump, there is only one contact point between two particles. The contact normal is collinear with the branch vector, which is the vector that connects the centroids of two particles. Thus, the normal contact will not impart any additional resistance torque on two bodies. In the cases of nonspherical bodies, more than one contact point will be observed and these extra contact forces will resist the free rotation of particles to some extent. Additionally, in a polydisperse system, small size particles tend to rest in the concave area of large size particles, which leads to the compactness of assemblies and hence a large angle of repose will be observed.To prove the idea, the variation of coordinate number, which is defined as the average number of contacts around a single particle, of the five models were monitored. The coordinate number (CN) is a parameter that reflects the contact density of assemblies at the scale of particles. According to the definition used by O’Sullivan [34], the CN can be calculated as follows:(3)CN=2×NcNp,where Nc is the total number of contacts and Np is the number of particles of the assemblies. Here, the divisor is multiplied by 2 as each contact is shared between two bodies. During the lifting process, the CN was recorded every second for five models, and the results are shown in Figure 19. For all five models, the CN decreases typically during first 4 seconds of simulation and then gradually increases until a stable value is reached. This is because that the number of contacts between particles drops as the collapse of the initial balance state, and more contacts are recreated inside the particle assemblies in order to reach another balance state at the final stage. It is particularly noteworthy that as particle convexity increases, the CN significantly decreases, which means that the convex area along the surface of particles increases the number of contacts of the assemblies and thus the angle of repose is affected as a result. However, the variation in CN for the real ballast models, the 4-ball models, and the 2-ball without overlap models is insignificant during the whole process. This is due to the reason that the CN is influenced not only by particle shapes but also by their sizes. There is an ultimate limit of contact intensity for individual particle when its size is determined as explained by O’Sullivan [34], and the 3 nonspherical models in this research just reach this critical state.Figure 19
The effect of particle convexity on the coordinate number.


## 3.3.1. Characteristics of Particle Movement

In this research, the translational and rotational velocities of particles located at three different heights (i.e., the upper, the middle, and the bottom parts of the cylinder) over the first eight seconds were monitored. At each height, 4 particles located in the circumference of the assemblies were picked as they are more likely to roll along the surface of the pile compared to those located in the middle part of the assemblies. Their velocities were averaged and fitted as presented by the movement curves for five different material shapes in Figures15–17, and their characteristics can be summarized as follows:(1)
For translational movement depicted in Figure15, all five types of particles initially increase to its peak velocity and then gradually decrease nearly to zero. However, the stabilization period is much shorter for nonspherical models with only about 2 seconds after the first 8 seconds than that of 1-ball model with more than 40 seconds which has been discarded in the graph for a better comparison. Notably, for particles located at the middle and top parts, the time of peak velocity is relatively earlier for 1-ball models than that for other shapes of models especially those with small particle shape factors (i.e., real ballast, 4-ball, and 2-ball particles). The difference is due to the lack of sufficient interlocking and hence the resistance between particles with large sphericity and convexity. These features are more pronounced on rotational movement for five models as shown in Figure 16.(2)
The average translational velocity curves in Figure17(a) show that the translational velocity has a close relationship with particle shape. Particles located at the upper part of cylinder have a greater average translational velocity, and this value decreases as the height decreased for all five types of models. For the four simply created clump templates, the average translational velocity becomes smaller with particle shape factors increasing. At the bottom of cylinder, the translational velocity of 1-ball templates is significantly higher than that of other three custom clump templates, which leads to the severe scatter of particles on the base plane as observed in the final simulation states. Figure 17(a) also shows that the average translational velocities for real ballast are relatively lower compared with those artificial ones in all three heights, which is in accordance with what was observed in the research by Zeng et al. [21].(3)
When it comes to the average rotational velocity as shown in Figure17(b), it is obvious that, at all three levels, the spherical or the nearly spherical models exhibit a much more distinctive rotation. They have a relatively larger average rotational velocity than those less spherical models. With a decrease in particle shape factors, the average rotational velocity drops significantly. More importantly, for 1-ball models and 2-ball with overlap models, the rotational velocity for particles located at the bottom of the cylinder is much greater than that located at the middle and the upper parts. By contrast, for particles with less shape factors, the rotational velocity gradually decreases from the upper parts to the bottom. This observation explains why particles with small sphericity values can form excellent grain piles on the base plane while it was extremely hard for spherical or the 2-ball with overlap clump templates to be stable when released from the cylinder. During the lifting process, spherical or nearly spherical particles from the bottom are much easier to roll away with a higher rotational velocity when suffered from horizontal push by upper particles, resulting in the local collapse of particle assemblies. However, in the cases of nonspherical particles, the bottom particles are difficult to rotate and they can build up a strong base for those particles rolling from upper parts, leading to a larger angle of repose compared to 1-ball and 2-ball with overlap models.Figure 15
The translational velocities for particles located at different parts of the cylinder. (a) Bottom part. (b) Middle part. (c) Upper part.
(a)
(b)
(c)Figure 16
The rotational velocities for particles located at different parts of the cylinder. (a) Bottom part. (b) Middle part. (c) Upper part.
(a)
(b)
(c)Figure 17
The average velocities for particles located at different heights of the cylinder. (a) Translational movement. (b) Rotational movement.
(a)
(b)The histories of kinetic energy in the system against simulation time for five models are presented in Figure18. It is obvious from Figure 18(a) that the kinetic energy for all five types of particle templates went through an initial fluctuation with a local maximum at the very beginning of the simulation. This is due to the onset of the movement and the fact that the lifting cylinder breaks the balance that has been maintained in the preparation phase. During this initial stage, the kinetic energy of the 1-ball model is relatively higher than other models and no significant difference of the kinetic energy is observed between the remaining 4 types of models. This shows that the 1-ball model is much more susceptible to external disturbances for the lack of interlocking and rolling resistance between particles. The kinetic energy climbs up to its peak value when all the particles in the cylinder are released out 3 to 4 seconds after the lifting starts, followed by a long period of rebalance and a decrease nearly close to zero. For four custom clump templates, the peak kinetic energy becomes greater with larger shape factors. It should be mentioned that the 1-ball model needs a much longer time, up to about 50 seconds, to reachieve the same balance state as other models, which indicates that the abundant energy stored in the system is hard to dissipate for the lack of rolling resistance, and the last 40 seconds were abandoned for a better comparison. Compared to the other four manually created template models, the real ballast model has a relatively smaller kinetic energy as shown in Figure 18(b) . This can be explained as the result of shape diversity in the real ballast systems. The similarity in particle shapes tends to result in larger void ratios for artificial assemblies even if having the same size distributions as that of real ballast systems. By contrast, the voids of real ballast can be better filled with clumps with different shapes and hence lead to a denser assemblies, in which particles movement can be restricted.Figure 18
The kinetic energy variation of particle assemblies for five models. (a) 0–0.062 s; (b) 0–10 s.In this research, we set the particle densities of five templates (i.e. 2700kg/m3) same rather than the bulk densities of five specimens. However, the variation in bulk densities may change the dynamics of particle flow and hence the angle of repose. The impacts of bulk density on the angle of repose for coarse assemblies will be further investigated in the future.


## 3.3.2. Characteristics of Particle Contacts

The concave area along the outlines of the clump templates increases the contact points between particles, allowing for an additional generation of resistance torque; thus the movement of particles is hindered, and a much larger angle of repose is generated. For spherical templates, i.e., 1-ball clump, there is only one contact point between two particles. The contact normal is collinear with the branch vector, which is the vector that connects the centroids of two particles. Thus, the normal contact will not impart any additional resistance torque on two bodies. In the cases of nonspherical bodies, more than one contact point will be observed and these extra contact forces will resist the free rotation of particles to some extent. Additionally, in a polydisperse system, small size particles tend to rest in the concave area of large size particles, which leads to the compactness of assemblies and hence a large angle of repose will be observed.To prove the idea, the variation of coordinate number, which is defined as the average number of contacts around a single particle, of the five models were monitored. The coordinate number (CN) is a parameter that reflects the contact density of assemblies at the scale of particles. According to the definition used by O’Sullivan [34], the CN can be calculated as follows:(3)CN=2×NcNp,where Nc is the total number of contacts and Np is the number of particles of the assemblies. Here, the divisor is multiplied by 2 as each contact is shared between two bodies. During the lifting process, the CN was recorded every second for five models, and the results are shown in Figure 19. For all five models, the CN decreases typically during first 4 seconds of simulation and then gradually increases until a stable value is reached. This is because that the number of contacts between particles drops as the collapse of the initial balance state, and more contacts are recreated inside the particle assemblies in order to reach another balance state at the final stage. It is particularly noteworthy that as particle convexity increases, the CN significantly decreases, which means that the convex area along the surface of particles increases the number of contacts of the assemblies and thus the angle of repose is affected as a result. However, the variation in CN for the real ballast models, the 4-ball models, and the 2-ball without overlap models is insignificant during the whole process. This is due to the reason that the CN is influenced not only by particle shapes but also by their sizes. There is an ultimate limit of contact intensity for individual particle when its size is determined as explained by O’Sullivan [34], and the 3 nonspherical models in this research just reach this critical state.Figure 19
The effect of particle convexity on the coordinate number.


## 4. Conclusion

Based on an image-processing method, four types of simple clumps together with 40 real ballast templates have been reconstructed and quantified with sphericity and convexity to study the impact of particle shape and the interparticle friction coefficient on the angle of repose on a ballast-size scale in DEM simulations. By comparing the results of laboratory cylinder lifting tests and that of the DEM simulations, it can be concluded that particle shape has a profound effect on the angle of repose for coarse aggregates. Despite that rolling friction was excluded in the current research, it is still a promising method to model angular particles and will be studied in further research. The results can be summarized as follows:(1)
Generally, the angle of repose increases as the coefficient of friction between particles increases for all models no matter whether the simple clump templates or the real ballast templates are used. However, all five types of models react variously to the variation ofμs,pp. Models using the real ballast template and the 4-ball template perform better than the other three types of models in the parameter calibration, without rolling friction coefficient being used in the models.(2)
The angle of repose decreases with the increase of particle sphericity and convexity. Models using real ballast clumps produce the most satisfying results compared to the other four simply created ones. The more deviated of particle shape factors from that of real ballast, the more discrepant of simulation results with that of the experimental tests. This implies that particle shape factors ought to be included into the parameter calibration for coarse particles in DEM simulations.(3)
Spherical particles tend to have greater translational and rotational velocities than nonspherical ones do at three different heights. Particularly, the rotational velocity for 1-ball model and 2-ball with overlap model at the bottom of the cylinder is much larger than that at the middle and upper parts, while the situation is just opposite for the other three types of models. The higher rotation at the bottom instead of middle or upper part is more likely to result in a local collapse of grain assemblies, and this explains why an excellent pile cannot be formed when clump templates with larger particle shape factors are used. The kinetic energy of real ballast assemblies is the least during the lifting process, which is due to particle shape diversity and clumps with various real ballast shapes can better fill the voids between particles, making the assemblies denser.The coordination number during the lifting process for five models gradually drops to its minimum and then rises to a constant value in the remaining period of the simulation. Moreover, the CN in models with spherical or nearly spherical particle templates is apparently lower than that with small convexity ones, which means that a larger concave area along the surface of particles would increase the possibility of contact and hence the angle of repose.(4)
Particle shape factors are supposed to be taken into consideration for parameter calibration in future research. By using 4-ball or 2-ball without overlap particles which have similar sphericities and convexities to that of real ballast, adequate results could be provided at a shorter simulation time and a higher computational efficiency. However, calibration of other necessary parameters such as friction coefficient has to be carefully carried out if such simplified particles are used.



---

*Source: 8095267-2019-09-18.xml*