# Correlation Analysis of Macroscopic and Microscopic Parameters of Coal Measure Soil Based on Discrete Element Method

Authors: Hong Zhang; Wei Liao; Junbo Lin; Bang Zhang; Haoran Feng
Journal: Advances in Civil Engineering
 (2019)

---

## Abstract

Numerical simulation of the triaxial test of coal measure expansive soil distributed along the highways in Pingxiang District, Jiangxi, was carried out by means of discrete particle flow, during which the macromechanical properties and the formation and developmental patterns of shear displacement field of the coal measure expansive soil were studied from a mesoscopic perspective. The result showed that the macroscopic stress and strain of test specimens can be significantly influenced by the interparticle friction coefficient of the coal measure expansive soil. Peak value of the deviatoric stress of test specimens increased with increasing friction coefficient, and before reaching the deviatoric stress peak value, the stress-strain relationship of the soil body basically presented a linear variation trend; the soil interparticle contact stiffness varied hyperbolically with the deviatoric stress peak value of test specimens, and the increasing contact stiffness ratio led to a gradual decrease of the deviatoric stress peak value but had only a small impact on the residual strength of test specimens; confining pressure was found to have remarkable influence on both the deviatoric stress peak value and the residual strength of test specimens; when the experimental confining pressure increased from 0.2 MPa to 1.2 MPa, the deviatoric stress peak value and the residual strength of test specimens increased by 2.14 times and 5.11 times, respectively. This paper reveals the macroinstability and failure mechanism of coal measure expansive soil from a microperspective.


---

## Body

## 1. Introduction

Coal measure soil mostly appears ash black in appearance, while the soil mass exposed is greyish with weak interlayer bonding and loose structure. Since coal measure soil is easily weathered and disintegrated with poor water stability after excavation, large gullies are likely to be caused by water and soil loss during heavy rain. Therefore, an urgent problem which should be solved currently is how to handle and utilize engineering with coal measure soil. Some researchers at home and abroad have researched problems and reinforcement measures for engineering with coal measure soil [1–3]. Nevertheless, characteristics of engineering with coal measure soil and causes of the abovementioned problems are rarely documented. Focusing on the coal measure soil along the Guangzhou-Wuzhou expressway, Hu et al. [4] determined the shear strength of coal measure soil under different water contents and its influence law with the routine straight shear test. Zhu et al. [5] established a mathematical model regarding the stability of the slope infiltrated by rainfall and analyzed the influence of cracks and weathered layers on slope stability using the limit equilibrium method together with the shear strength theory of unsaturated soil. Zhang et al. [6] simulated and analyzed the slope stability of coal measure soil with different cracks and weathering degrees under rainfall infiltration using the transient unsaturated seepage software based on some secondary slope of Daqing-Guangzhou expressway in Jiangxi. Li and Liu [7] took the red clay and powdered coal measure soil in hilly and mountainous areas of southwestern China as the objects of study to investigate the influencing rule of different water contents on parameters such as the cohesion, internal friction angle, and shear strength using the common straight shear test. Yang and Zheng [8] determined the soil-water characteristic curve for samples of remolded unsaturated coal measure soil with different dry densities using a soil moisture characteristic curve meter, and they fit the data measured with the VG model. Most of the abovementioned documents tend to research from the perspective of macromechanics. However, macromechanical property is a symbol of parameter changes in the microstructure of soil mass. In that case, the slope of coal measure soil can be considered as a system comprised of soil particles in different sizes. Those soil particles in the system are discrete, and the interaction between the particles is characterized by nonlinear dissipation. To be specific, it is the micromechanical behavior of particles inside the system that determines many macromechanical properties of the slope. Hence, microscopic analysis must be performed on the composition of particle structure of the slope soil and distribution characteristics of contact force chain and movement characteristics of particles, so as to disclose the mechanism of macroscopic static and dynamic features shown by the granular system of the slope of coal measure soil [9–14].In 1979, Cundall and Strack [15, 16] proposed the distinct element method (DEM), which has gradually become an effective tool to study the microstructures of granular systems. Many scholars modeled typical soil mechanics experiments (e.g., biaxial compression and direct shear tests) using the discrete element method [17–23]. The unloading of silos is a typical bulk flow process. After the appearance of the discrete element method, Langston et al. [24–29] used DEM to study hopper unloading. In the field of complex soil dynamics or machine-soil dynamics, some scholars use the discrete element method to simulate the excavation, drilling, shoveling, penetration, and pipeline laying processes [30–36].A large amount of coal measure strata containing montmorillonite and illite are distributed along the Pingxiang-Lianhua expressway in Jiangxi. As a result, expansive soil can be found in coal measures. Nonetheless, few researches on mechanical properties of the expansive coal measure soil have been reported currently. In this paper, the discrete element method (DEM) is employed to have a 3D numerical simulation of the triaxial test on the expansive coal measure soil, to probe into the influencing rule of parameter changes of the micromechanic structure of the soil mass particles on macromechanical characteristics of the expansive coal measure soil, and to discuss the formation and development law of the shear displacement field of the expansive coal measure soil. The research findings can provide theoretical references for the protective slope design and construction of the expansive coal measure soil.


## 2. Principles of Discrete Element Method

### 2.1. Physical Equations

The discrete element method (DEM) is to consider objects of research as an integral whole consisting of a finite number of discrete particles, in which each particle is idealized as a rigid body and deemed as a discrete calculating unit. The particles are independent of each other with contact and friction. Assuming that the normal contact forcesFin between two contact solid elements are directly proportional to the normal relative displacement Un between them, the following algebraic relation can be obtained [37]:(1)Fin=KnUnni,where Kn is the normal stiffness of the distinct element; Un is the normal displacement of the distinct element; and ni is the unit normal vector over the contact surface.The tangential contact force between granular particles is defined as the shear force whose size is closely related to the movement track of the particles and their loading history. In this paper, the shear force between granular particles is expressed in the form of increment. When a contact is formed between granular particles, the total tangential contact forceFis between the particles is assumed to be initialized to zero. And the increment of the tangential contact force resulting from the relative displacement between the particles is accumulated in each cycle of calculating time steps. The calculating relationship can be expressed as follows:(2)ΔFis=−KsΔUs,where ΔFis is the shear stress increment between distinct elements in each calculation step; Ks is the tangential contact stiffness of the distinct element; and ΔUs is the increment of the tangential relative displacement of the distinct element.


### 2.2. Equations of Motion

Based on equation (2), the resultant force and the resultant moment imposed on granular particles after the accumulation of each calculation of time steps are calculated. After that, the acceleration and the rotational angular velocity of granular particles are calculated according to Newton’s second law. Thus, physical quantities such as velocity, linear displacement, and angular displacement of the granular particles in the calculation period can be known. If the resultant shearing force acting on the granular particles in the direction of x at the moment of t0 is assumed as Fx, the motion of granular particles can be expressed by the following two equations:(3)u¨xt0=Fxm,ω¨xt0=MxIx,where u¨xt0 is the average acceleration of granular particles at the moment of t0 in the direction of x; ω¨xt0 is the rotational acceleration of granular particles at the moment of t0 in the direction of x; m is the mass of granular particles; Mx is the bending moment of the resultant force acting on the granular particles; and Ix is the rotational inertia of granular particles.


### 2.3. Particle Contact Constitutive Model

Since many coal particles are contained in expansive coal measure soil, especially that in the weathering layer of the slope is loose, it results in weak cohesive force between the soils. The shear strength is formed primarily by the internal friction between soil particles. Therefore, a contact sliding model can be utilized to describe the constitutive relationship between the normal force and the tangential shear force interacting between particles. The normal tension between particles is neglected in the contact sliding model. The sliding motion between particles only occurs under set conditions, that is, the shearing force actually generated between particles should be greater than the preset threshold; also, the sliding between particles must be performed within a certain range of shear strength. Assuming that the displacement overlapping between two particles is less than or equal to zero, it can be considered that the two particles have no contact or their contact force is zero. In this way, their constitutive behavior can be described as(4)Fmaxs=μFin,where μ is the friction coefficient between particles. If so, then the particles will slide. Meanwhile, in the next cycle step, Fis will be(5)Fis⟵FisFmaxsFis.


## 2.1. Physical Equations

The discrete element method (DEM) is to consider objects of research as an integral whole consisting of a finite number of discrete particles, in which each particle is idealized as a rigid body and deemed as a discrete calculating unit. The particles are independent of each other with contact and friction. Assuming that the normal contact forcesFin between two contact solid elements are directly proportional to the normal relative displacement Un between them, the following algebraic relation can be obtained [37]:(1)Fin=KnUnni,where Kn is the normal stiffness of the distinct element; Un is the normal displacement of the distinct element; and ni is the unit normal vector over the contact surface.The tangential contact force between granular particles is defined as the shear force whose size is closely related to the movement track of the particles and their loading history. In this paper, the shear force between granular particles is expressed in the form of increment. When a contact is formed between granular particles, the total tangential contact forceFis between the particles is assumed to be initialized to zero. And the increment of the tangential contact force resulting from the relative displacement between the particles is accumulated in each cycle of calculating time steps. The calculating relationship can be expressed as follows:(2)ΔFis=−KsΔUs,where ΔFis is the shear stress increment between distinct elements in each calculation step; Ks is the tangential contact stiffness of the distinct element; and ΔUs is the increment of the tangential relative displacement of the distinct element.


## 2.2. Equations of Motion

Based on equation (2), the resultant force and the resultant moment imposed on granular particles after the accumulation of each calculation of time steps are calculated. After that, the acceleration and the rotational angular velocity of granular particles are calculated according to Newton’s second law. Thus, physical quantities such as velocity, linear displacement, and angular displacement of the granular particles in the calculation period can be known. If the resultant shearing force acting on the granular particles in the direction of x at the moment of t0 is assumed as Fx, the motion of granular particles can be expressed by the following two equations:(3)u¨xt0=Fxm,ω¨xt0=MxIx,where u¨xt0 is the average acceleration of granular particles at the moment of t0 in the direction of x; ω¨xt0 is the rotational acceleration of granular particles at the moment of t0 in the direction of x; m is the mass of granular particles; Mx is the bending moment of the resultant force acting on the granular particles; and Ix is the rotational inertia of granular particles.


## 2.3. Particle Contact Constitutive Model

Since many coal particles are contained in expansive coal measure soil, especially that in the weathering layer of the slope is loose, it results in weak cohesive force between the soils. The shear strength is formed primarily by the internal friction between soil particles. Therefore, a contact sliding model can be utilized to describe the constitutive relationship between the normal force and the tangential shear force interacting between particles. The normal tension between particles is neglected in the contact sliding model. The sliding motion between particles only occurs under set conditions, that is, the shearing force actually generated between particles should be greater than the preset threshold; also, the sliding between particles must be performed within a certain range of shear strength. Assuming that the displacement overlapping between two particles is less than or equal to zero, it can be considered that the two particles have no contact or their contact force is zero. In this way, their constitutive behavior can be described as(4)Fmaxs=μFin,where μ is the friction coefficient between particles. If so, then the particles will slide. Meanwhile, in the next cycle step, Fis will be(5)Fis⟵FisFmaxsFis.


## 3. Numerical Simulation of Particle Flow of Expansive Coal Measure Soil

A particle flow calculation model is firstly established in this paper according to the actual particle size range of the soil mass by assuming that particles of the soil mass are rigid materials. Meanwhile, the above contact sliding model is used to indicate the particle interaction. After simulating the stress-strain relationship of test pieces in the triaxial test by changing microscomic parameters of particles, the results obtained are compared with those from the laboratory triaxial test. Finally, microscomic parameters of the soil particles that can reflect macromechanical properties of the soil mass are determined through numerous trial calculations and analyses. Furthermore, the calculation model is taken as a basic model for analyzing the influencing law of changes of particle microscopic parameters on macroscopic characteristics of the expansive coal measure soil.


### 3.1. Principles of Numerical Calculation of the Triaxial Test

The key to the calculation principle of numerical simulation in the triaxial test is to calculate the stress state of the test piece in the loading process along with stress control. While the test piece is loaded by controlling the movement velocity of the bottom plate and the top plate of the model constraint boundary in the triaxial test that is simulated with the discrete particle flow program [38]. The stress of the test piece can be obtained by dividing the resultant force of all the forces acting on the test piece by the area of the corresponding test piece. The strain of the test piece in the directions of x and y can be obtained as follows:(6)ε=L−L01/2L0+L,where L is the current length of the test piece in the corresponding direction and L0 is the original length of the test piece in that direction. Concerning the cylindrical test piece, the axial strain εH and the radial strain εR can be expressed as follows:(7)εH=H−H01/2H0+H,εR=R−R01/2R0+R.Servo functions can be called in the cycle of calculating the time steps to reduce the difference between the monitoring stress and the preset stress in the calculation of numerical simulation. The numerical servo functions are adopted to adjust the velocity of conducting lateral restraint against the wall of the test piece, so that the restraint stress can be inclined to a certain constant. The numerical servo functions are achieved by the following algorithms; the velocity of boundary restraint can be expressed as follows:(8)u˙ω=Gσmeasured−σrequired=GΔσ,where u˙ω is the velocity of laterally restraining the wall and G is an “incremental” parameter obtained through the following method.Assuming that the increment of the binding force generated in the movement of lateral wall restraining within a calculated step can be expressed as follows:(9)ΔFω=KnωNcu˙ωΔt,where Nc is the number of contact surfaces on the constraint boundary and Knω is the average stiffness of these contact surfaces. Therefore, the average increment of restraint stress is(10)Δσω=KnωNcu˙ωΔtA,where A is the area of regions involved in the constraint boundary. In order to achieve stable test loading in the numerical simulation of the triaxial test, the absolute value of the variable of restraint stress must be smaller than the difference between the monitored value and the preset value. Amplification of coefficients is adopted in the simulation of actual numerical simulation. α is assumed as an amplification coefficient. To guarantee the loading stability, we may have(11)Δσω<αΔσ,which leads to(12)KnωNcu˙ωΔtA<αΔσ,where Knω>0, Nc>0, Δt>0, and A>0; thus,(13)KnωNcΔtAGΔσ<αΔσ.EliminateΔσ by taking G>0 and Δσ>0, we can obtain(14)KnωNcΔtAG<α.The following equation can be obtained by transposing the inequation:(15)G<αAKnωNcΔt.Therefore, the increment can be determined by the following equation:(16)G=αAKnωNcΔt.In general,α=0.5 is taken in actual simulated calculations.


### 3.2. Test Materials and Laboratory Triaxial Test

The test soil samples source from the slope of original coal measures at K12 + 110 in A2 section of Wanzai-Yichun expressway in Jiangxi. Natural physical and mechanical indexes of the samples are as follows: the average density of soil particles is 2230 kg/m3; the average natural density is 1680 kg/m3; the natural porosity is 0.32; and the moisture content W of the test pieces is 18.2%. The distribution of particle contents in different particle size ranges of the expansive coal measure soil can be obtained by conducting particle grading test on the samples with the sieving method and the mixed method, as shown in Table 1. Soil particles mainly range from 0.25 to 10 mm, accounting for 85.41%; with d10 = 0.129, d30 = 0.58, and d60 = 1.27; the nonuniform coefficient Cu is 9.84; and the curvature coefficient Cc is 2.08. The SLB-1 triaxial shear permeameter manufactured by the Nanjing Soil Instrument Factory is used in the test which takes unconsolidated-undrained (UU) shearing. The test pieces are 40 mm in diameter and 80 mm in height. The triaxial shear test is conducted under the ambient pressure of 0.1 MPa, 0.3 MPa, and 0.6 MPa, respectively.Table 1
Test results of particle size distribution of coal measure expansive soils.
Test material
Particle size (mm)
Content (%)
Coal measure expansive soil
>20
2.26
10∼20
7.73
5∼10
14.63
2∼5
36.31
0.5∼2
24.3
0.25∼0.5
10.17
0.075∼0.25
3.02
>0.075
1.58


### 3.3. Numerical Modeling of Triaxial Test and Determination of Microscopic Parameters

Sizes of numerical modeling and particle simulation are magnified to save time for calculation and improve the efficiency of analysis. The numerical model of the triaxial test is 2 meters in height and 0.8 meters in diameter. As can be known from the results of the laboratory particle grading test for expansive coal measure soil in Table1, the particle size ranges from 0.25 to 10 mm. The particle diameter should be set within this range in modeling. And a total of 3,058 solid particles can be generated in the form of uniform distribution. The calculation model is shown in Figure 1.Figure 1
Computational model for triaxial numerical experiments. (a) Boundary wall. (b) Particle set.
(a)
(b)The walls on the top and at the bottom of the particle flow soil samples are driven to conduct relative movement via the numerical servo system under consolidated and ambient pressure of 0.1 MPa, 0.3 MPa, and 0.6 MPa, so as to unload the test pieces. At the same time, the ambient pressure of the test pieces can be guaranteed unchanged through adjustment of the side wall displacement before the end of loading. Numerical calculation results are approximated to the results of the laboratory triaxial test via repeated trial calculation by changing microscopic parameters of the particles. The final calibrated results of microscopic parameters of the particles are obtained, as shown in Table2. Three-axis simulation results calculated by DEM are compared with the laboratory test results, as shown in Figure 2. It can be seen from the figure that the peak value of the deviatoric stress calculated by numerical simulation of the triaxial test is relatively close to the laboratory triaxial test results with consistent change trends of the stress-strain curve. Thus, it indicates that the microscopic parameters calibrated finally can relatively approach the real model.Table 2
Microscopic parameters used in the discrete element method simulations.
Particle density
2230 kg/m3
Porosity
0.32
Model particle diameter
20 mm∼200 mm
Friction coefficient
0.3
Normal stiffness of particle contact
1 × 109 N/m
Tangential stiffness of particle contact
1 × 107 N/mFigure 2
The results of deviation stress and strain of the triaxial test are compared with those of numerical simulation.


## 3.1. Principles of Numerical Calculation of the Triaxial Test

The key to the calculation principle of numerical simulation in the triaxial test is to calculate the stress state of the test piece in the loading process along with stress control. While the test piece is loaded by controlling the movement velocity of the bottom plate and the top plate of the model constraint boundary in the triaxial test that is simulated with the discrete particle flow program [38]. The stress of the test piece can be obtained by dividing the resultant force of all the forces acting on the test piece by the area of the corresponding test piece. The strain of the test piece in the directions of x and y can be obtained as follows:(6)ε=L−L01/2L0+L,where L is the current length of the test piece in the corresponding direction and L0 is the original length of the test piece in that direction. Concerning the cylindrical test piece, the axial strain εH and the radial strain εR can be expressed as follows:(7)εH=H−H01/2H0+H,εR=R−R01/2R0+R.Servo functions can be called in the cycle of calculating the time steps to reduce the difference between the monitoring stress and the preset stress in the calculation of numerical simulation. The numerical servo functions are adopted to adjust the velocity of conducting lateral restraint against the wall of the test piece, so that the restraint stress can be inclined to a certain constant. The numerical servo functions are achieved by the following algorithms; the velocity of boundary restraint can be expressed as follows:(8)u˙ω=Gσmeasured−σrequired=GΔσ,where u˙ω is the velocity of laterally restraining the wall and G is an “incremental” parameter obtained through the following method.Assuming that the increment of the binding force generated in the movement of lateral wall restraining within a calculated step can be expressed as follows:(9)ΔFω=KnωNcu˙ωΔt,where Nc is the number of contact surfaces on the constraint boundary and Knω is the average stiffness of these contact surfaces. Therefore, the average increment of restraint stress is(10)Δσω=KnωNcu˙ωΔtA,where A is the area of regions involved in the constraint boundary. In order to achieve stable test loading in the numerical simulation of the triaxial test, the absolute value of the variable of restraint stress must be smaller than the difference between the monitored value and the preset value. Amplification of coefficients is adopted in the simulation of actual numerical simulation. α is assumed as an amplification coefficient. To guarantee the loading stability, we may have(11)Δσω<αΔσ,which leads to(12)KnωNcu˙ωΔtA<αΔσ,where Knω>0, Nc>0, Δt>0, and A>0; thus,(13)KnωNcΔtAGΔσ<αΔσ.EliminateΔσ by taking G>0 and Δσ>0, we can obtain(14)KnωNcΔtAG<α.The following equation can be obtained by transposing the inequation:(15)G<αAKnωNcΔt.Therefore, the increment can be determined by the following equation:(16)G=αAKnωNcΔt.In general,α=0.5 is taken in actual simulated calculations.


## 3.2. Test Materials and Laboratory Triaxial Test

The test soil samples source from the slope of original coal measures at K12 + 110 in A2 section of Wanzai-Yichun expressway in Jiangxi. Natural physical and mechanical indexes of the samples are as follows: the average density of soil particles is 2230 kg/m3; the average natural density is 1680 kg/m3; the natural porosity is 0.32; and the moisture content W of the test pieces is 18.2%. The distribution of particle contents in different particle size ranges of the expansive coal measure soil can be obtained by conducting particle grading test on the samples with the sieving method and the mixed method, as shown in Table 1. Soil particles mainly range from 0.25 to 10 mm, accounting for 85.41%; with d10 = 0.129, d30 = 0.58, and d60 = 1.27; the nonuniform coefficient Cu is 9.84; and the curvature coefficient Cc is 2.08. The SLB-1 triaxial shear permeameter manufactured by the Nanjing Soil Instrument Factory is used in the test which takes unconsolidated-undrained (UU) shearing. The test pieces are 40 mm in diameter and 80 mm in height. The triaxial shear test is conducted under the ambient pressure of 0.1 MPa, 0.3 MPa, and 0.6 MPa, respectively.Table 1
Test results of particle size distribution of coal measure expansive soils.
Test material
Particle size (mm)
Content (%)
Coal measure expansive soil
>20
2.26
10∼20
7.73
5∼10
14.63
2∼5
36.31
0.5∼2
24.3
0.25∼0.5
10.17
0.075∼0.25
3.02
>0.075
1.58


## 3.3. Numerical Modeling of Triaxial Test and Determination of Microscopic Parameters

Sizes of numerical modeling and particle simulation are magnified to save time for calculation and improve the efficiency of analysis. The numerical model of the triaxial test is 2 meters in height and 0.8 meters in diameter. As can be known from the results of the laboratory particle grading test for expansive coal measure soil in Table1, the particle size ranges from 0.25 to 10 mm. The particle diameter should be set within this range in modeling. And a total of 3,058 solid particles can be generated in the form of uniform distribution. The calculation model is shown in Figure 1.Figure 1
Computational model for triaxial numerical experiments. (a) Boundary wall. (b) Particle set.
(a)
(b)The walls on the top and at the bottom of the particle flow soil samples are driven to conduct relative movement via the numerical servo system under consolidated and ambient pressure of 0.1 MPa, 0.3 MPa, and 0.6 MPa, so as to unload the test pieces. At the same time, the ambient pressure of the test pieces can be guaranteed unchanged through adjustment of the side wall displacement before the end of loading. Numerical calculation results are approximated to the results of the laboratory triaxial test via repeated trial calculation by changing microscopic parameters of the particles. The final calibrated results of microscopic parameters of the particles are obtained, as shown in Table2. Three-axis simulation results calculated by DEM are compared with the laboratory test results, as shown in Figure 2. It can be seen from the figure that the peak value of the deviatoric stress calculated by numerical simulation of the triaxial test is relatively close to the laboratory triaxial test results with consistent change trends of the stress-strain curve. Thus, it indicates that the microscopic parameters calibrated finally can relatively approach the real model.Table 2
Microscopic parameters used in the discrete element method simulations.
Particle density
2230 kg/m3
Porosity
0.32
Model particle diameter
20 mm∼200 mm
Friction coefficient
0.3
Normal stiffness of particle contact
1 × 109 N/m
Tangential stiffness of particle contact
1 × 107 N/mFigure 2
The results of deviation stress and strain of the triaxial test are compared with those of numerical simulation.


## 4. Correlation Analysis of Microscopic Parameters and Macroscopic Mechanics of Coal Measure Soil

### 4.1. Relationship between Friction Coefficients and Stress-Strain

The friction coefficientμ between particles is set as 0.1, 0.3, 0.5, 0.7, and 0.9 for numerical simulation in order to discuss the influence of friction coefficients on the macroscopic stress and strain of the expansive coal measure soil. Microscopic parameters of other particles in the calculation model can have assignment according to Table 1. The influencing relationship curve of different friction coefficients on stress-strain of expansive coal measure soil is shown in Figure 3.Figure 3
Graph showing the relationship between bias stress and axial strain with different friction coefficients: (a) bias stress and axial strain; (b) peak deviating stress and friction coefficient.
(a)
(b)According to the deviatoric stress-axial strain curve in Figure3(b), as the friction coefficient μ between the particles grows, the peak value of the deviatoric stress of the test piece will also increase continuously. The stress-strain relationship of the soil mass tends to be in linear change before the test piece reaches the peak. The soil mass enters the strain softening stage right after that. As can be seen from Figure 3(a), strain softening of the soil is apparent. As the friction coefficient increases, a growth of the interaction force between soil particles can be observed; and the force for overcoming the sliding between particles is also improved. Thus, there is an increase of the modulus of deformation of the soil. Certain enhancement of the residual strength of the soil can be seen with the friction coefficient after the test piece is loaded and unloaded, which, however, is insignificant. Moreover, the residual strength of the soil increases insignificantly when the friction coefficient is large.Figure4 presents the relationship curve between the volumetric strain and the axial strain of the test piece under different friction coefficients, indicating that the volumetric strain decreases with the increase of axial strain at the initial stage of loading, whereas particle deformation changes into tangential sliding from normal compression and overlapping between particles. Consequently, the volume of the test piece becomes smaller, macroscopically resulting in shrinkage. The volumetric strain increases with the axial strain at the later stage of loading. The test piece is slid tangentially with the increase in volume, macroscopically leading to dilatancy. The volumetric strain of the test piece is enlarged with the increase of friction coefficients.Figure 4
Graph showing the relationship between volumetric strain and axial strain with different friction coefficients.


### 4.2. Relationship between Contact Stiffness and Stress-Strain of Particles

Based on the original calculation model, the ambient pressure is 0.3 MPa; and the contact stiffnessKs is 1 GPa, 3 GPa, 5 GPa, 7 GPa, 9 GPa, and 12 GPa, respectively. Other parameters are constant. Numerical simulation of the triaxial test specific to the test piece is conducted in 6 working conditions. The stress-strain relationships of the test piece with different contact stiffness are shown in Figures 5 and 6.Figure 5
Graph showing the relationship between different tangential contact stiffness and stress-strain: (a) deviation stress and axial strain; (b) peak deviating stress and tangential contact stiffness.
(a)
(b)Figure 6
Graph showing the relationship between volumetric strain and axial strain with different tangential contact stiffness.As can be seen from Figure5, the deviatoric stress peak of the test piece decreases with the increase of the contact stiffness of the particles. When the contact stiffness Ks increases from 1 GPa to 3 GPa, the decline of the deviatoric stress peak of the test piece is significant; when the contact stiffness Ks increases from 5 GPa to 12 GPa, the deviatoric stress peak of the test piece is slow; when the contact stiffness Ks is large, the residual strength of the test piece tends to be the same value (as shown in Figure 5(a)). According to the test results of numerical simulation, with the increase of the contact stiffness Ks between particles of the test piece, the macroscopic initial tangential modulus of the test piece will also experience an increase, while the deviatoric stress peak and the axial strain corresponding to the strength reaching the peak value become smaller, suggesting that the contact stiffness has little influence on the residual strength of the test piece.It can be observed from Figure6 that the deformation of the test piece is constituted primarily by normal particle compression and overlapping as the volumetric strain decreases with the increase of the axial strain at the initial stage of loading. As a result, the test piece is shrunk. The influence of the contact stiffness of particles on the shrinkage of the test piece is not obvious at the early stage. As the test piece is further loaded, its volumetric strain becomes positive and increases with the axial strain, whereas particle deformation changes into tangential sliding from normal compression and overlapping between particles and causes volume expansion, which is macroscopically manifested as dilatancy. The volumetric strain of the test piece rises with the contact stiffness of the test piece. Such increase, however, is unapparent at the later stage.As can be seen from Figure7 that presents calculation results of the stiffness ratio Kn/Ks between different particles, the initial tangential modulus of the test piece at the early stage is decreased and the strength of the peak value of the deviatoric stress gradually decreases with the increase of the particle stiffness ratio, exerting little influence on the residual strength of the test piece. Furthermore, the stiffness ratio Kn/Ks has an invisible influence on characteristics such as shrinkage and dilatancy of the test piece.Figure 7
Graph showing the relationship between stress and strain with different contact stiffness ratios of particles: (a) deviation stress and axial strain; (b) volumetric strain and axial strain.
(a)
(b)


### 4.3. Stress-Strain Relationship under Different Ambient Pressures

The stress state around the test piece in the triaxial test is changed to explore the changing law of deviatoric stress and axial strain under different ambient pressures. Ambient pressures calculated are valued as 0.2 MPa, 0.4 MPa, 0.6 MPa, 0.8 MPa, 1.0 MPa, and 1.2 MPa, while other microscopic parameters remain constant. The calculation results are shown in Figures8 and 9.Figure 8
Graph showing the relationship between stress and strain with different confining pressures: (a) deviation stress and axial strain; (b) peak deviating stress and confining pressure.
(a)
(b)Figure 9
Graph showing the relationship between volumetric strain and axial strain with different confining pressures.From Figure8, the initial tangential modulus of the test piece is growing (the slope of the curve increases gradually) before the softening of the curve strain with the increase of the ambient pressure. In the meantime, the deviatoric stress peak and the residual strength of the test piece are also improved constantly to realize linear change. As a result, ambient pressure has appreciable influence on the residual strength of the test piece.The relationship curve of the volumetric strain and the axial strain under different ambient pressures of the test piece is shown in Figure9. It can be observed that the volumetric strain of the test piece is negative at the early loading stage. As the axial strain increases, the test piece is shrunk gradually, which becomes increasingly obvious with the increase of the ambient pressure. When the ambient pressure is small, the test piece is diverted to dilatancy from shrinkage. It can also be observed that the volumetric strain of the test piece constantly decreases with the increase of the ambient pressure. When the ambient pressure is small, the declining range is significant; when the ambient pressure is large, the declining range of the volumetric strain is insignificant.


### 4.4. Analysis of the Influence of Porosity on Shear Strength Parameters

Based on the original calculation model, the ambient pressure is set as 0.8 MPa; and the porosity of the test piece is set as 0.3, 0.35, 0.4, 0.45, and, 0.5, respectively, while other parameters maintain constant, thus to analyze the stress-strain relationships of the test piece under different porosity conditions. The simulated calculation results are shown in Figure10. As can be seen from Figure 10(a), the porosity has certain influence on the deviatoric stress peak of the test piece. The deviatoric stress peak of the test piece is large when the porosity is small; that is, the larger the porosity is, the smaller the deviatoric stress peak of the test piece will be. Nonetheless, the porosity has an insignificant influence on the residual strength of the test piece.Figure 10
Graph showing the relationship between stress and strain with different porosity: (a) deviation stress and axial strain; (b) volumetric strain and axial strain.
(a)
(b)The relationship curve of volumetric strain and axial strain under different porosity conditions of the test piece is shown in Figure10(b). It can be observed when the test piece with low porosity (samples are dense) is sheared, a shrinkage will be witnessed for the test piece with the decreased volume at the early loading stage. After that, the volume of the test piece will rapidly increase, macroscopically manifesting as dilatancy in an obvious manner, as the particles in the test piece are sheared and moved. When the test piece with high porosity (samples are loose) is sheared, its volume will decrease before the subsequent rebound, mainly manifesting shrinkage. Moreover, the volume increase of the test piece gradually decreases with the increase of porosity.


### 4.5. Discussion on the Evolution Law of Shear Displacement Field of Expansive Coal Measure Soil

The displacement and velocity change paths of particles in the test piece can be monitored and tracked in real time in the loading process of numerically simulating the triaxial test with the particle flow code. Moreover, the formation of the shear displacement field and the change law of the stress field of the test piece can be analyzed through elaboration of the displacement vector. Figure11 shows the vector diagram presenting the formation and evolution of the shear displacement field of the test piece under different axial strains (calculation of time steps) when the ambient pressure is 0.6 MPa. To be specific, Figures 11(a) and 11(b) show the formation of the shear displacement vector at the initial loading stage of the test piece. It can be seen from Figure 11(a) that particles in the test piece are moved towards the central section at this time with the insignificant shape of the shear zone, indicating that the test piece is at the early stage of dilatancy. It can be observed from Figure 11(b) that the shear zone of the test piece is formed initially as loading continues. At this time, the motion of particle displacement vector approaches the direction of 45°. Whereas Figures 11(c) and 11(d) show the evolution of the shear displacement vector after the deviatoric stress of the test piece reaches the peak. Since particle motion changes from the low potential energy to the high potential energy, the directional trend of the shear displacement field becomes significant with the increase of the axial strain of the test piece. Meanwhile, the thickness of the shear zone is getting smaller in the state of instability.Figure 11
Graph showing the formation and variation of shear displacement field of specimen particles during loading (the confining pressure is 0.6 MPa): (a)ε1=1.35%; (b) ε1=2.32%; (c) ε1=8.28%; (d) ε1=11.15%.
(a)
(b)
(c)
(d)Distribution of the shear displacement field of the test piece under different ambient pressures is presented in Figure12. It can be observed that the width of the shear displacement belt of particles is gradually narrowed as the ambient pressure increases. Analyzing the microscopic mechanism shows that it is mainly because the resistance of particle motion in the restrained test piece increases with the ambient pressure of the test piece. That is to say, the difficulty of displacing the relative sliding of the particle has also increased, which is further manifested by the enhanced deviatoric stress peak of the test piece from the macroperspective.Figure 12
Graph showing the distribution of shear displacement field of specimens with different confining pressures: (a)σ3=0.1MPa; (b) σ3=0.4MPa; (c) σ3=0.8MPa; (d) σ3=1.6MPa.
(a)
(b)
(c)
(d)


## 4.1. Relationship between Friction Coefficients and Stress-Strain

The friction coefficientμ between particles is set as 0.1, 0.3, 0.5, 0.7, and 0.9 for numerical simulation in order to discuss the influence of friction coefficients on the macroscopic stress and strain of the expansive coal measure soil. Microscopic parameters of other particles in the calculation model can have assignment according to Table 1. The influencing relationship curve of different friction coefficients on stress-strain of expansive coal measure soil is shown in Figure 3.Figure 3
Graph showing the relationship between bias stress and axial strain with different friction coefficients: (a) bias stress and axial strain; (b) peak deviating stress and friction coefficient.
(a)
(b)According to the deviatoric stress-axial strain curve in Figure3(b), as the friction coefficient μ between the particles grows, the peak value of the deviatoric stress of the test piece will also increase continuously. The stress-strain relationship of the soil mass tends to be in linear change before the test piece reaches the peak. The soil mass enters the strain softening stage right after that. As can be seen from Figure 3(a), strain softening of the soil is apparent. As the friction coefficient increases, a growth of the interaction force between soil particles can be observed; and the force for overcoming the sliding between particles is also improved. Thus, there is an increase of the modulus of deformation of the soil. Certain enhancement of the residual strength of the soil can be seen with the friction coefficient after the test piece is loaded and unloaded, which, however, is insignificant. Moreover, the residual strength of the soil increases insignificantly when the friction coefficient is large.Figure4 presents the relationship curve between the volumetric strain and the axial strain of the test piece under different friction coefficients, indicating that the volumetric strain decreases with the increase of axial strain at the initial stage of loading, whereas particle deformation changes into tangential sliding from normal compression and overlapping between particles. Consequently, the volume of the test piece becomes smaller, macroscopically resulting in shrinkage. The volumetric strain increases with the axial strain at the later stage of loading. The test piece is slid tangentially with the increase in volume, macroscopically leading to dilatancy. The volumetric strain of the test piece is enlarged with the increase of friction coefficients.Figure 4
Graph showing the relationship between volumetric strain and axial strain with different friction coefficients.


## 4.2. Relationship between Contact Stiffness and Stress-Strain of Particles

Based on the original calculation model, the ambient pressure is 0.3 MPa; and the contact stiffnessKs is 1 GPa, 3 GPa, 5 GPa, 7 GPa, 9 GPa, and 12 GPa, respectively. Other parameters are constant. Numerical simulation of the triaxial test specific to the test piece is conducted in 6 working conditions. The stress-strain relationships of the test piece with different contact stiffness are shown in Figures 5 and 6.Figure 5
Graph showing the relationship between different tangential contact stiffness and stress-strain: (a) deviation stress and axial strain; (b) peak deviating stress and tangential contact stiffness.
(a)
(b)Figure 6
Graph showing the relationship between volumetric strain and axial strain with different tangential contact stiffness.As can be seen from Figure5, the deviatoric stress peak of the test piece decreases with the increase of the contact stiffness of the particles. When the contact stiffness Ks increases from 1 GPa to 3 GPa, the decline of the deviatoric stress peak of the test piece is significant; when the contact stiffness Ks increases from 5 GPa to 12 GPa, the deviatoric stress peak of the test piece is slow; when the contact stiffness Ks is large, the residual strength of the test piece tends to be the same value (as shown in Figure 5(a)). According to the test results of numerical simulation, with the increase of the contact stiffness Ks between particles of the test piece, the macroscopic initial tangential modulus of the test piece will also experience an increase, while the deviatoric stress peak and the axial strain corresponding to the strength reaching the peak value become smaller, suggesting that the contact stiffness has little influence on the residual strength of the test piece.It can be observed from Figure6 that the deformation of the test piece is constituted primarily by normal particle compression and overlapping as the volumetric strain decreases with the increase of the axial strain at the initial stage of loading. As a result, the test piece is shrunk. The influence of the contact stiffness of particles on the shrinkage of the test piece is not obvious at the early stage. As the test piece is further loaded, its volumetric strain becomes positive and increases with the axial strain, whereas particle deformation changes into tangential sliding from normal compression and overlapping between particles and causes volume expansion, which is macroscopically manifested as dilatancy. The volumetric strain of the test piece rises with the contact stiffness of the test piece. Such increase, however, is unapparent at the later stage.As can be seen from Figure7 that presents calculation results of the stiffness ratio Kn/Ks between different particles, the initial tangential modulus of the test piece at the early stage is decreased and the strength of the peak value of the deviatoric stress gradually decreases with the increase of the particle stiffness ratio, exerting little influence on the residual strength of the test piece. Furthermore, the stiffness ratio Kn/Ks has an invisible influence on characteristics such as shrinkage and dilatancy of the test piece.Figure 7
Graph showing the relationship between stress and strain with different contact stiffness ratios of particles: (a) deviation stress and axial strain; (b) volumetric strain and axial strain.
(a)
(b)


## 4.3. Stress-Strain Relationship under Different Ambient Pressures

The stress state around the test piece in the triaxial test is changed to explore the changing law of deviatoric stress and axial strain under different ambient pressures. Ambient pressures calculated are valued as 0.2 MPa, 0.4 MPa, 0.6 MPa, 0.8 MPa, 1.0 MPa, and 1.2 MPa, while other microscopic parameters remain constant. The calculation results are shown in Figures8 and 9.Figure 8
Graph showing the relationship between stress and strain with different confining pressures: (a) deviation stress and axial strain; (b) peak deviating stress and confining pressure.
(a)
(b)Figure 9
Graph showing the relationship between volumetric strain and axial strain with different confining pressures.From Figure8, the initial tangential modulus of the test piece is growing (the slope of the curve increases gradually) before the softening of the curve strain with the increase of the ambient pressure. In the meantime, the deviatoric stress peak and the residual strength of the test piece are also improved constantly to realize linear change. As a result, ambient pressure has appreciable influence on the residual strength of the test piece.The relationship curve of the volumetric strain and the axial strain under different ambient pressures of the test piece is shown in Figure9. It can be observed that the volumetric strain of the test piece is negative at the early loading stage. As the axial strain increases, the test piece is shrunk gradually, which becomes increasingly obvious with the increase of the ambient pressure. When the ambient pressure is small, the test piece is diverted to dilatancy from shrinkage. It can also be observed that the volumetric strain of the test piece constantly decreases with the increase of the ambient pressure. When the ambient pressure is small, the declining range is significant; when the ambient pressure is large, the declining range of the volumetric strain is insignificant.


## 4.4. Analysis of the Influence of Porosity on Shear Strength Parameters

Based on the original calculation model, the ambient pressure is set as 0.8 MPa; and the porosity of the test piece is set as 0.3, 0.35, 0.4, 0.45, and, 0.5, respectively, while other parameters maintain constant, thus to analyze the stress-strain relationships of the test piece under different porosity conditions. The simulated calculation results are shown in Figure10. As can be seen from Figure 10(a), the porosity has certain influence on the deviatoric stress peak of the test piece. The deviatoric stress peak of the test piece is large when the porosity is small; that is, the larger the porosity is, the smaller the deviatoric stress peak of the test piece will be. Nonetheless, the porosity has an insignificant influence on the residual strength of the test piece.Figure 10
Graph showing the relationship between stress and strain with different porosity: (a) deviation stress and axial strain; (b) volumetric strain and axial strain.
(a)
(b)The relationship curve of volumetric strain and axial strain under different porosity conditions of the test piece is shown in Figure10(b). It can be observed when the test piece with low porosity (samples are dense) is sheared, a shrinkage will be witnessed for the test piece with the decreased volume at the early loading stage. After that, the volume of the test piece will rapidly increase, macroscopically manifesting as dilatancy in an obvious manner, as the particles in the test piece are sheared and moved. When the test piece with high porosity (samples are loose) is sheared, its volume will decrease before the subsequent rebound, mainly manifesting shrinkage. Moreover, the volume increase of the test piece gradually decreases with the increase of porosity.


## 4.5. Discussion on the Evolution Law of Shear Displacement Field of Expansive Coal Measure Soil

The displacement and velocity change paths of particles in the test piece can be monitored and tracked in real time in the loading process of numerically simulating the triaxial test with the particle flow code. Moreover, the formation of the shear displacement field and the change law of the stress field of the test piece can be analyzed through elaboration of the displacement vector. Figure11 shows the vector diagram presenting the formation and evolution of the shear displacement field of the test piece under different axial strains (calculation of time steps) when the ambient pressure is 0.6 MPa. To be specific, Figures 11(a) and 11(b) show the formation of the shear displacement vector at the initial loading stage of the test piece. It can be seen from Figure 11(a) that particles in the test piece are moved towards the central section at this time with the insignificant shape of the shear zone, indicating that the test piece is at the early stage of dilatancy. It can be observed from Figure 11(b) that the shear zone of the test piece is formed initially as loading continues. At this time, the motion of particle displacement vector approaches the direction of 45°. Whereas Figures 11(c) and 11(d) show the evolution of the shear displacement vector after the deviatoric stress of the test piece reaches the peak. Since particle motion changes from the low potential energy to the high potential energy, the directional trend of the shear displacement field becomes significant with the increase of the axial strain of the test piece. Meanwhile, the thickness of the shear zone is getting smaller in the state of instability.Figure 11
Graph showing the formation and variation of shear displacement field of specimen particles during loading (the confining pressure is 0.6 MPa): (a)ε1=1.35%; (b) ε1=2.32%; (c) ε1=8.28%; (d) ε1=11.15%.
(a)
(b)
(c)
(d)Distribution of the shear displacement field of the test piece under different ambient pressures is presented in Figure12. It can be observed that the width of the shear displacement belt of particles is gradually narrowed as the ambient pressure increases. Analyzing the microscopic mechanism shows that it is mainly because the resistance of particle motion in the restrained test piece increases with the ambient pressure of the test piece. That is to say, the difficulty of displacing the relative sliding of the particle has also increased, which is further manifested by the enhanced deviatoric stress peak of the test piece from the macroperspective.Figure 12
Graph showing the distribution of shear displacement field of specimens with different confining pressures: (a)σ3=0.1MPa; (b) σ3=0.4MPa; (c) σ3=0.8MPa; (d) σ3=1.6MPa.
(a)
(b)
(c)
(d)


## 5. Macroinstability of Coal Measure Soil and Micromechanism Analysis of Shear Band Localization

Cundall and Strack [39] as well as Bardet and Proubet [40] studied the localization of shear strain with the PFC discrete element analysis software. Vardoulakis and Aifantis [41] and Chambon et al. [42] theoretically explored the mechanism of production of shear zones based on the second gradient. In this paper, the author analyzed the mesoscopic mechanism of localization of shear zones with the volume fraction of particles, coordination number, and other mesoscopic parameters. The gradation of test samples in the triaxial test originated from the undisturbed coal-bearing soil samples on the side slope at mark K12 + 110 of section A2 on the Wanzai-Yichun expressway in Yichun, Jiangxi. The gradation is shown in Table 1. The numerical model in the triaxial test was 2.0 m in height and 0.8 m in diameter; the average soil particle density was determined as 2230 kg·m−3, and the natural porosity was 0.32. The values of mesoscopic parameters of the test samples were in accordance with the calibrated values in Table 2. To track the parametric variation in the particle system, a measuring sphere (MS1) outside the shear zone and one (MS2) inside the shear zone were selected to highlight the variation of volume fraction, coordination number, and other parameters within and outside the shear zone.


### 5.1. Macroscopic Instability of Coal Measure Soil

The results of deviator stress and volumetric strain of the test samples simulated by the triaxial test are shown in Figure13. Figure 13(a) shows the complete process of elastoplasticity, which the particle system underwent. When the axial strain was ε < 2.0%, the deviator stress and axial strain exhibited a linear change at the macrolevel, and the particle system exhibited elastic properties. When ε increased, the particle system entered the plastic stage. When the axial strain was 2% ≤ ε ≤ 3.5%, the deviator stress of the test samples underwent violent fluctuations and peaked at the axial strain value of 2.4%. Later, the deviator stress gradually decreased; during this process, the test samples hardened. When the axial strain was ε > 3.5%, the particle system entered the softening stage, where the deviator stress considerably decreased, which indicated the presence of internal damage to the system. According to Figure 13(b), when the particle system was under pressure, the volumetric strain εv decreased to a negative value; when the particle volume fraction n increased, the stress in the system increased. When the axial strain was ε < 2.0%, the system was in an elastic strain. When the axial strain was ε = 3.2%, the volumetric strain was εv = −1.8% while the system achieved its highest density and maximum volume fraction. When the axial strain was ε > 3.2%, the volumetric strain εv increased and the particle system entered the plastic stage. When the axial strain was ε = 5.3%, the volumetric strain became εv = 0, which indicated that the volumetric strain of the system returned to the initial state. When ε > 5.3%, the volume fraction of the particle system gradually decreased while the volumetric strain experienced slight fluctuations. At the same time, the change in slope remained constant. The system experienced a dilative shear failure.Figure 13
Simulated results of the triaxial test for specimens: the relation curve between (a) deviating stress and axial strain and (b) volumetric strain and axial strain.
(a)
(b)


### 5.2. Mesomechanism Analysis of Shear Band Localization of Coal Measure Soil

Figure14 illustrates the relationship curve between the axial strain and volume fraction, which were measured using the measuring spheres. Figure 14 shows that the volume fraction n of both MS1 and MS2 gradually increased at first, while the particle system remained at the elastic stage, which was represented as shrinkage. When the axial strain was 2% ≤ ε ≤ 3.2%, the volume fraction at the shear zone (MS2) remained the same at 65.38%, with the system shifting to the hardening stage. When ε = 3.2%, the MS2 volume fraction n exhibited a slight decrease with a small range of variation and the particle system may have internally developed small cracks. The MS2 volume fraction n exhibited only minor changes before ε< 4.5%. When ε = 4.5%, the MS2 volume fraction n progressed in a sustained decline and the system entered the critical state, which indicated the growth of the dominant shear zone. When the axial strain increased, the in-shear zone (MS2) volume fraction n exhibited an overall downward trend with occasional small upward and downward deviations, which implies that the dilatancy of the particle system was primarily caused by the loose arrangement of particles within the shear zone. The varying pattern of the outer-shear zone (MS1) volume fraction was similar to that of MS2.Figure 14
The curve of relationship between volume fraction and axial strain: (a) the relationship between volume fraction and axial strain at the measuring sphere 1 (outside shear band); (b) the relationship between volume fraction and axial strain at the measuring sphere 2 (in shear band).
(a)
(b)Figure15 shows the variation curve of the coordination number and axial strain. It is observed that when the axial strain was ε < 2%, as the loading process of the entire particle system continued, the coordination number Z gradually increased. When 2% ≤ ε ≤ 4.5%, the coordination number Z of MS1 and MS2 exhibited small oscillations, and the maximum coordination number of MS1 was 5.5 and that of MS2 was 5.84. When ε = 4.5%, the coordination number Z of MS1 and MS2 clearly decreased, which indicates the increase in porosity in the central area of the particle system. Meanwhile, the volumetric strain gradually increased from the minimum, which indicates a clear dilatancy of the entire particle system during the softening stage. When the axial strain was ε = 5.8%, the coordination number of the particle system increased again. Therefore, it was concluded that the dilative shear failure of the particles included slippage between the particles and angular displacement, which was generated by particle rotation. When the axial strain was 7.5%, the coordination number of the particle system remained practically unchanged and the coordination numbers of MS1 and MS2 were 4.7 and 5.1, respectively.Figure 15
Relation between coordination number and axial strain: (a) the relationship between coordination number and axial strain at the measuring sphere 1 (outside shear band); (b) the relationship between coordination number and axial strain at the measuring sphere 2 (in shear band).
(a)
(b)


## 5.1. Macroscopic Instability of Coal Measure Soil

The results of deviator stress and volumetric strain of the test samples simulated by the triaxial test are shown in Figure13. Figure 13(a) shows the complete process of elastoplasticity, which the particle system underwent. When the axial strain was ε < 2.0%, the deviator stress and axial strain exhibited a linear change at the macrolevel, and the particle system exhibited elastic properties. When ε increased, the particle system entered the plastic stage. When the axial strain was 2% ≤ ε ≤ 3.5%, the deviator stress of the test samples underwent violent fluctuations and peaked at the axial strain value of 2.4%. Later, the deviator stress gradually decreased; during this process, the test samples hardened. When the axial strain was ε > 3.5%, the particle system entered the softening stage, where the deviator stress considerably decreased, which indicated the presence of internal damage to the system. According to Figure 13(b), when the particle system was under pressure, the volumetric strain εv decreased to a negative value; when the particle volume fraction n increased, the stress in the system increased. When the axial strain was ε < 2.0%, the system was in an elastic strain. When the axial strain was ε = 3.2%, the volumetric strain was εv = −1.8% while the system achieved its highest density and maximum volume fraction. When the axial strain was ε > 3.2%, the volumetric strain εv increased and the particle system entered the plastic stage. When the axial strain was ε = 5.3%, the volumetric strain became εv = 0, which indicated that the volumetric strain of the system returned to the initial state. When ε > 5.3%, the volume fraction of the particle system gradually decreased while the volumetric strain experienced slight fluctuations. At the same time, the change in slope remained constant. The system experienced a dilative shear failure.Figure 13
Simulated results of the triaxial test for specimens: the relation curve between (a) deviating stress and axial strain and (b) volumetric strain and axial strain.
(a)
(b)


## 5.2. Mesomechanism Analysis of Shear Band Localization of Coal Measure Soil

Figure14 illustrates the relationship curve between the axial strain and volume fraction, which were measured using the measuring spheres. Figure 14 shows that the volume fraction n of both MS1 and MS2 gradually increased at first, while the particle system remained at the elastic stage, which was represented as shrinkage. When the axial strain was 2% ≤ ε ≤ 3.2%, the volume fraction at the shear zone (MS2) remained the same at 65.38%, with the system shifting to the hardening stage. When ε = 3.2%, the MS2 volume fraction n exhibited a slight decrease with a small range of variation and the particle system may have internally developed small cracks. The MS2 volume fraction n exhibited only minor changes before ε< 4.5%. When ε = 4.5%, the MS2 volume fraction n progressed in a sustained decline and the system entered the critical state, which indicated the growth of the dominant shear zone. When the axial strain increased, the in-shear zone (MS2) volume fraction n exhibited an overall downward trend with occasional small upward and downward deviations, which implies that the dilatancy of the particle system was primarily caused by the loose arrangement of particles within the shear zone. The varying pattern of the outer-shear zone (MS1) volume fraction was similar to that of MS2.Figure 14
The curve of relationship between volume fraction and axial strain: (a) the relationship between volume fraction and axial strain at the measuring sphere 1 (outside shear band); (b) the relationship between volume fraction and axial strain at the measuring sphere 2 (in shear band).
(a)
(b)Figure15 shows the variation curve of the coordination number and axial strain. It is observed that when the axial strain was ε < 2%, as the loading process of the entire particle system continued, the coordination number Z gradually increased. When 2% ≤ ε ≤ 4.5%, the coordination number Z of MS1 and MS2 exhibited small oscillations, and the maximum coordination number of MS1 was 5.5 and that of MS2 was 5.84. When ε = 4.5%, the coordination number Z of MS1 and MS2 clearly decreased, which indicates the increase in porosity in the central area of the particle system. Meanwhile, the volumetric strain gradually increased from the minimum, which indicates a clear dilatancy of the entire particle system during the softening stage. When the axial strain was ε = 5.8%, the coordination number of the particle system increased again. Therefore, it was concluded that the dilative shear failure of the particles included slippage between the particles and angular displacement, which was generated by particle rotation. When the axial strain was 7.5%, the coordination number of the particle system remained practically unchanged and the coordination numbers of MS1 and MS2 were 4.7 and 5.1, respectively.Figure 15
Relation between coordination number and axial strain: (a) the relationship between coordination number and axial strain at the measuring sphere 1 (outside shear band); (b) the relationship between coordination number and axial strain at the measuring sphere 2 (in shear band).
(a)
(b)


## 6. Conclusions

Based on the theory of discrete element method (DEM), this paper performed numerical calculations on the triaxial test of expansive coal measure soil as well as simulated the macromechanical properties of the coal measure soil as well as the formation and development of the shear displacement field from the perspective of particle microscopic parameters, from which, results similar to the laboratory test can be obtained. By analyzing microscopic parameters of the particles and the influencing law of the macroscopic stress-strain curve of the expansive coal measure soil, the following conclusions can be made.(1)
The friction coefficient between soil particles exerts a significant influence on the macroscopic stress and strain of the test piece. As the friction coefficientμ between particles grows, the deviatoric stress peak of the test piece also increases continuously. The stress-strain relationship of the soil mass tends to be in linear change before the test piece reaches the peak. The soil mass enters the strain softening stage immediately after reaching the peak, whereas the increase in the residual strength of soil mass is not obvious when the friction coefficient is large.(2)
The contact stiffness between particles increases, whereas the deviatoric stress peak of the test piece decreases, basically presenting a hyperbolic change. When the contact stiffness ratioKn/Ks increases, the deviatoric stress peak of the test piece will decrease gradually with nearly no effect on the residual strength of the test piece. Meanwhile, the stiffness ratio Kn/Ks has insignificant effect on characteristics such as shrinkage and dilatancy of the test piece.(3)
The deviatoric stress peak and the residual strength are significantly affected by the ambient pressure in the test. By analyzing from the microscopic mechanism, it is mainly a result of the resistance increase of particle motion in the restrained test piece with the ambient pressure of the test piece. That is to say, the difficulty of displacing the relative sliding of the particle has also increased, which is further manifested by the enhanced deviatoric stress peak of the test piece from the macroperspective.(4)
The monitoring of volume fraction and coordination number in shear band shows that the volume fraction in shear band of coal measure soil decreased with the increase of axial strain when the granular system was in the softening stage, which indicated that the dilatancy behavior of the granular system was basically caused by the loose distribution of particles in the shear band. At the same time, the coordination number began to decrease significantly, indicating that the porosity in the central region of the particle system began to increase, and the dilatancy failure behavior of the particles included not only the slip between the particles but also the angular displacement caused by the rotation of the particles.



---

*Source: 1302682-2019-09-22.xml*