Capability of discrete element method to investigate the macro-micro mechanical behaviours of granular soils considering different stress conditions and morphological gene mutation

Wei Xiong, Jinfeng Wng,b,*, Zhung Cheng

a Department of Architecture and Civil Engineering, City University of Hong Kong, Hong Kong, China

b Shenzhen Research Institute of City University of Hong Kong, Shenzhen, China

c School of Civil Engineering and Architecture, Wuhan University of Technology, Wuhan, China

d Sanya Science and Education Innovation Park, Wuhan University of Technology, Sanya, China

Keywords:

ABSTRACT Discrete element method (DEM) has been widely utilised to model the mechanical behaviours of granular materials.However, with simplified particle morphology or rheology-based rolling resistance models,DEM failed to describe some responses,such as the particle kinematics at the grain-scale and the principal stress ratio against axial strain at the macro-scale.This paper adopts a computed tomography(CT)-based DEM technique, including particle morphology data acquisition from micro-CT (μCT),spherical harmonic-based principal component analysis (SH-PCA)-based particle morphology reconstruction and DEM simulations, to investigate the capability of DEM with realistic particle morphology for modelling granular soils’ micro-macro mechanical responses with a consideration of the initial packing state, the morphological gene mutation degree, and the confining stress condition.It is found that DEM with realistic particle morphology can reasonably reproduce granular materials’ micro-macro mechanical behaviours, including the deviatoric stress-volumetric strain-axial strain response, critical state behaviour,particle kinematics,and shear band evolution.Meanwhile,the role of multiscale particle morphology in granular soils depends on the initial packing state and the confining stress condition.For the same granular soils, rougher particle surfaces with a denser initial packing state and a higher confining stress condition result in a higher degree of shear strain localisation.

1.Introduction

Granular materials are ubiquitous in nature, from the sand castles by little builders on the beach to the computer chip that is indispensable in our daily lives.In the domain of soil mechanics and powder technology,the study of granular materials has a pedigree going back to the last century, trying to understand and predict how the granular materials behave upon various loading conditions(Terzaghi et al., 1996; Maeda et al., 2010; Wang and Yan, 2013;Radjai et al.,2017;Zhao et al.,2017,2020;Cheng and Wang,2018a;Nie et al., 2020; Xiong et al., 2022a).Many macro-mechanical responses of granular materials were studied, for instance, the small strain stiffness (Payan et al., 2016; He et al., 2022), stressdilatancy response (Gutierrez and Wang, 2009; Yin and Chang,2013; Xiao et al., 2019), critical state behaviour (Been and Jefferies, 1985; Chu, 1995), and shear band evolution (Desrues et al., 2018; Cheng and Wang, 2019).However, the constitutive modelling of granular materials is still mysterious until today.The complicated mechanical responses of the granular material stem from its intrinsic characteristics,such as discontinuity and heterogeneity.In this context,the granular materials’responses cannot be fully described by most current constitutive theories proposed to explain continuous medium (Xing et al., 2021).Therefore, accurately modelling the granular materials’ mechanical responses is still significant.

Fig.1.Schematic diagram of the CT-based DEM modelling technique.

The discrete element method(DEM)has provided direct access to the particle kinematics and the internal fabric evolution for granular materials (Cundall and Strack, 1979; Bathurst and Rothenburg, 1988; Mueth et al., 1998).It has become a powerful technique to model granular materials’ mechanical behaviours.Early DEM studies on granular materials employed spheres with a rheology-type rolling resistance model (Iwashita and Oda, 2000;Zhou et al., 2013; Qu et al., 2022) or utilised simplified nonspherical particles to model granular materials under loading(Rothenburg and Bathurst,1991;Azéma et al.,2013;Nassauer et al.,2013; Kuhn et al., 2015; de Bono and McDowell, 2020).Although yielding some good macro-mechanical responses, the early DEM failed to accurately depict the particle kinematics and fabric evolution at the grain-scale (Zhao and Zhao, 2019) as well as some macro-scale responses like principal stress ratio against axial strain(Qu et al.,2022).It was then believed that the accurate modelling of granular particles should consider both the micro-and macro-scale characteristics.

The rapid development of high-resolution imaging techniques(Santamarina and Cho, 2004; Altuhafi and Coop, 2011; Viggiani et al., 2015; Zhao and Wang, 2016; Zheng and Hryciw, 2017)paves the way to extract and represent the multiscale and highresolution particle morphology.Extensive research on the connections of macro-and micro-scale responses for granular soils has been carried out, for instance, particle morphology effects (Zhou et al., 2015; Zhao and Wang, 2016; Fei and Narsilio, 2020; Nie et al., 2021; Thakur and Penumadu, 2021; Xiong et al., 2022b),fabric evolution (Maeda et al., 2010; Fonseca et al., 2012, 2013;Cheng and Wang, 2018a), particle translation and rotation (Hall et al., 2010; Hasan and Alshibli, 2010; Andrade et al., 2011; Chengand Wang, 2018b), and particle breakage (Alikarami et al., 2015;Zhou et al., 2015; Karatza et al., 2018; Cheng and Wang, 2021).Many particle reconstruction techniques were proposed based on the extracted particle morphology (Grigoriu et al., 2006; Liu et al.,2011; Mollon and Zhao, 2014; Zhou and Wang, 2017; Wei et al.,2018; Sun and Zheng, 2021; Xiong and Wang, 2021; Wang et al.,2022).For reconstructing the three-dimensional (3D) particle morphology, the spherical harmonic analysis (SHA) was usually adopted.By treating the spherical harmonic(SH)coefficients as the“morphological gene” of the granular particles, Xiong and Wang(2021) successfully implemented the “morphological gene mutation”by varying the particle morphological features at any specific length scale through SH-based PCA (SH-PCA).

Table 1 Morphology descriptors for gene-mutated particles.

By considering the sophisticated and realistic particle morphology,modern DEM models established the bridge between micro- and macro-scale mechanical behaviours of granular soils(Kawamoto et al., 2018; Wu et al., 2021a).This paper further investigated the capability of DEM to model the micro-macro mechanical behaviours of granular soils, considering the initial packing state, the confining stress condition, and the morphological gene mutation degree.The remainder of this paper is structured as follows.Section 2 introduces the particle reconstruction framework to generate morphological gene-mutated particles with different degrees.Thirty three DEM samples in three groups are reconstructed, and the utilised DEM modelling technique is well discussed.Section 3 discusses the DEM results of the macro-micro mechanical responses for granular soils considering the initial packing state, the degree of morphological gene mutation and the confining stress condition.Section 4 gives the conclusions.

2.Virtual sample preparation framework

The combination of mechanical loading and in situ X-ray microcomputed tomography (μCT) has become a powerful technique to study the micro-macro mechanical behaviours of granular soils considering realistic particle morphology and different stress conditions.However, acquiring an adequate understanding of the soil micro-mechanical behaviours generally requires a large number of computed tomography (CT)-based physical tests, which are expensive and time-consuming.Thus, this paper adopted a CTbased DEM modelling technique to generate morphological genemutated samples to perform numerical triaxial testing.Fig.1 presents the schematic diagram of the CT-based DEM technique,including particle morphology data acquisition from μCT, SH-PCAbased particle morphology reconstruction and DEM simulations.

2.1.Particle morphology data through X-ray μCT

The entire in situ μCT testing of the Leighton Buzzard Sand(LBS)particles was conducted at the BL13W beamline of the Shanghai Synchrotron Radiation Facility (SSRF), China.This scanner consists of an X-ray parallel beam, a rotation stage to place the triaxial apparatus, and a detector to receive the transmitted X-ray energy.This study set the input X-ray energy and the detector resolution to 25 Kev and 6.5 μm, respectively,to reach a good contrast between different phases (i.e.sand grains and void spaces) (Cheng and Wang, 2018a).The rotation stage allowed the triaxial apparatus and a small sand sample to rotate constantly during the CT scanning.The testing sample(8 mm in diameter and 16 mm in height)contained around 1500 LBS grains(diameter of 0.6-1.18 mm).With its irregularly shaped particle surface,LBS has been investigated in the single particle compression test (Zhao et al., 2015) and the triaxial test (Cheng and Wang, 2018b).During the triaxial test, the testing sample was compressed axially at 0.2% per minute and confined laterally with 1.5 MPa.The axial compression loading was paused for CT scanning at different stages of axial strain,i.e.0%,5%,10% and 15%.

The raw projection slices(Fig.1a)obtained from the experiment went through a series of image processing steps to finally get the high-resolution particle morphology (Fig.1c) (Cheng and Wang,2018b).Firstly, the Gridrec algorithm (Dowd et al.,1999) and the‘ring artefact’ removal algorithm (Rivers et al., 2004; Chen et al.,2012) were used to transform the raw projection slices into greyscale slices.Secondly, a low-pass Gaussian filter and the Otsu(1979)’s thresholding method were adopted to obtain the binary images further.Thirdly, a marker-based morphological watershed algorithm can extract each particle within the given sample(Beucher and Lantuéjoul,1979).Finally, particle surface data were extracted by MATLAB function‘bwprim’and expressed by 3D point clouds.

2.2.Gene mutation through SH-PCA

Through the X-ray μCT reconstruction, each particle was expressed by a set of 3D point clouds.The morphological gene mutation on those μCT-reconstructed particles was implemented through SH-PCA.The main procedure was briefly introduced as follows (Xiong and Wang, 2021).

First,any individual particle can be expressed by SH coefficients through the following equations:

Table 2 Morphology descriptors.

Table 3 DEM model parameters after calibration.

where r（θ,φ）is the polar radius;Cmnis the SH coefficient;Ymn（θ,φ）is the SH function;Pmn（x）is the associated Legendre function;θ and φ are the zenith angleand azimuth angle of the spherical coordinates,respectively; and n and m are the frequency and the order,respectively.With given Cmn, Eq.(1) can be utilised to reconstruct the target particle surface to different length scales.In addition,the multiscale particle morphology, i.e.general form (GF), local roundness (LR) and the surface texture (ST), can be accurately reconstructed when the SH frequency is 15 or greater(Zhou et al.,2015, 2018; Xiong et al., 2020).

Since quantities in SH coefficients are complex numbers,the real and imaginary parts should be separated to make the analysis easier.Then,the whole sample B and any individual particle can be described as

Fig.2.Model validation: (a) Stress-strain curves; and (b) 3D rose map.

Fig.3.The relationship of deviatoric stress vs.axial strain for morphological gene-mutated samples under different confining stresses: (a) RDPS; and (b) RLPS.

where Cmnjdenotes the SH coefficients of particle j; real(·) and imag(·) extract the real and imaginary parts of complex numbers,respectively; uBand σ2Bare the average particle surface and the variance matrix of B,respectively;and y is the variation coefficient.

Due to the orthogonality of Ymnwith different frequencies,σ2Bin Eq.(5) was modified by setting the covariances at different frequencies to zero.After that,the principal component(PC)analysis(PCA) is implemented on the modified σ2B.The obtained eigenvalues were organised in a decreasing trend, and the eigenvalues were organised accordingly to form a matrix VPC.With VPC, any individual particle can then be generated as

where yt, vPCtand λtare the tth variation coefficient, eigenvector and eigenvalue,repectively;and Rddenotes the length scale-based PC ranges.Xiong and Wang(2021)conducted a parametric study to determine the PC division for each length scale.It was found that 1-20 PCs mainly contribute to the GF while 31-73 PCs to LR.According to this conclusion, Eq.(6)can be rewritten as follows:

2.3.DEM modelling of gene-mutated particles

Fig.4.The volumetric strain evolution of morphological gene-mutated RDPS upon different stress conditions: (a) σ3 = 0.3 MPa; (b) σ3 = 0.5 MPa; (c) σ3 = 1 MPa; and (d)σ3 = 1.5 MPa.

To simplify and highlight the effect of realistic particle morphology on DEM simulations in representing the constitutive behaviours of granular soils, one particle was chosen from the physical experiment at the initial loading state to generate the reference DEM sample.The morphological gene-mutated particles were reconstructed through SH-PCA and further utilised to form the corresponding morphological gene-mutated samples.Note that the morphological gene mutation by SH-PCA used morphology data of all particles within the physical sample.With varying χ in Eq.(7) from 0 to 1.5, different degrees of morphological gene mutation at either GF- or LR-level can be achieved (Fig.1d).Table 1 adopts four traditional morphology descriptors targeting different length scales to show the morphological variance of different morphological gene-mutated samples.The definition of each morphology descriptor can be found in Table 2.From Table 1, the morphological gene mutation at different length scales shows a different variation trend with the corresponding scaling factor.Specifically,the morphological gene mutation at GF-level(GM-GF)experiences an augmented aspect ratio with increasing χa, while the morphological gene mutation at LR-level(GM-LR) experiences a decreased roundness with increasing χc.

To avoid the inter-scale effect of particle morphology on the DEM simulations, only GM-LR was considered in this study.The resulted one reference particle and the corresponding three morphological gene-mutated particle geometries were then imported into particle flow code in 3 dimensions(PFC3D)to generate DEM particles with the bubble-packing algorithm (Fig.1d).A bonded ball membrane system was utilised to model the latex membrane from the physical triaxial test(Fig.1e).The bonded ball membrane deformed with the axial loading, and an external equivalent force Fewas employed and updated every 100 timesteps.The Feon each membrane ball was calculated as (Fazekas et al., 2006; de Bono et al., 2012; Qu et al., 2019):

where niand SAiare the unit normal vector and the surface area of the ith neighbouring triangle, respectively; and σcis the confining stress.All samples were then reiterated to an equilibrium state for further test.

This study reconstructed 33 samples in three groups, i.e.16 relatively dense packing samples(RDPS),16 relatively loose packing samples(RLPS),and a one-to-one mapping sample to calibrate the DEM model parameters.Note that the one-to-one mapping DEM sample contained particles that shared the same volume distribution, particle orientation and spatial location as the physical experiment.The DEM model parameters after calibration are shown in Table 3 and had good agreement with those in de Bono and McDowell (2014) and Wu et al.(2021a).Wu et al.(2021a)utilised the chosen DEM model parameters to conduct a parametric study on the one-to-one mapping sample and concluded that the DEM model responses match well with the X-ray μCT experiments.The stress-strain curve,volumetric strain curve and 3D rose map of the one-to-one mapping DEM sample during the triaxial test were compared with those from the physical experiment (as shown in Fig.2).The excellent agreements in all these three aspects indicate the reasonability and capability of the DEM model parameters chosen in this paper.Additionally,three degrees of morphological gene mutation at LR(i.e.χc=0,0.5,and 1.5),four different confining stresses(i.e.σ3=0.3 MPa,0.5 MPa,1 MPa,and 1.5 MPa), and two types of initial packing states were considered.Note that the two types of initial packing states were reconstructed by setting different friction coefficients during the sample compaction process,i.e.μ=0 for RDPS and μ=0.5 for RLPS.During triaxial shear tests, all the obtained DEM samples were loaded axially at a constant velocity of 0.005 m/s to an ultimate axial strain of 30%.

Fig.5.The volumetric strain evolution of morphological gene-mutated RLPS upon different stress conditions: (a) σ3 = 0.3 MPa; (b) σ3 = 0.5 MPa; (c) σ3 = 1 MPa; and (d)σ3 = 1.5 MPa.

3.Numerical results and discussion

3.1.Deviatoric stress against axial strain

Fig.3 shows the relationship of macroscopic deviatoric stress q versus axial strainlfor both RDPS and RLPS.In this paper,q andlare shown as

where σ1and σ3are the axial stress and the confining stress,respectively; and h0and h are the initial sample height and the corresponding deformed sample height,respectively.From Fig.3a,RDPS show a different variation trend with various stress conditions.Precisely, RDPS with higher confining stresses (for instance,σ3=1.5 MPa)displays a more robust strain-hardening response to an earlier and greater peak-state deviatoric stress, followed by a stronger strain-softening till the advent of the critical condition than RDPS with lower confining stresses (for instance,σ3=0.3 MPa).In the same stress condition,the curves in RDPS have a consistent variation trend for both the pre-peak as well as postpeak states, regardless of the degrees of morphological gene mutation involved.This finding implies that particle morphology plays an insignificant role in granular soils’ deviatoric stress/axial strain response at a relatively dense packing state with minor confining stress.

Compared to RDPS, the corresponding curves in RLPS (Fig.3b)experience a more minor initial stiffness,a much smaller and more delayed peak-state deviatoric stress, and a milder post-peak state strain-softening response.GM-LR is found to have notable effects on the macroscopic deviatoric stress/axial strain response of granular soils,especially for the pre-peak state.Specifically,a more extensive scaling factor(like χc=1.5)witnesses a more prominent initial stiffness,and a more significant peak deviatoric stress than a smaller scaling factor(like χc=0)for all different stress conditions.This is because the more prominent scaling factor yields a smaller roundness,leading to a rougher particle surface that contains more sharp edges and corners.The obtained particle surface is much rougher with more minor particle rotation in grain-scale,and thus a more considerable initial stiffness and peak-state deviatoric stress result.Besides,the more significant confining stress leads to a more notable difference in the macroscopic deviatoric stress/axial strain curves for different morphological gene mutation degrees.

All the macroscopic deviatoric stress/axial strain characteristics of granular soils agree well with the basic features of physical experiments (Nie et al., 2020, 2021; Sharma et al., 2021; Qu et al.,2022), which reflects the feasibility and capability of DEM with realistic particle morphology to describe the macro-scale deviatoric stress/axial strain response of granular soils under different stress conditions.

Fig.6.Evolution of the critical state in e-log10p space with different packing states: (a) RDPS; and (b) RLPS.

Fig.7.Evolution of the critical state in e-log10p space with different stress conditions: (a) σ3 = 0.3 MPa; (b) σ3 = 0.5 MPa; (c) σ3 = 1 MPa; and (d) σ3 = 1.5 MPa.

3.2.Volumetric strain evolution

Figs.4 and 5 show the variations of volumetric strainVagainst axial strainlunder various stress conditions for both RDPS and RLPS.In this paper,Vis defined as

whereR= ΔR/R0is the radial strain; V0and R0are the sample volume and average membrane radius at the initial state, respectively; and ΔV and ΔR are the corresponding increments at the deformed state.

Fig.8.The relationship of principal stress ratio vs.axial strain for morphological gene-mutated RDPS under different confining stresses:(a)χc =0;(b)χc =0.5;(c)χc =1;and(d)χc = 1.5.

RDPS in Fig.4 exhibits a similar overall variation of volumetric strain/axial strain curves.The morphological gene mutation degrees and the confining stress conditions seem to show marginal effects.In general,samples in RDPS show a minor contraction at the initial stage, a significant dilation after that, and a nearly constant volumetric strain at the critical state.GM-LR shows marginal effects on the evolution of volumetric strain/axial strain curves.The reason is that the relatively denser packing state of RDPS leads to high interlocking forces, which impede the particle translation and particle rotation within the samples.Besides, the maximum volumetric strain shows a decreasing trend with augmenting confining stresses.The higher confining stresses restrict the granular assemblage from expanding upon triaxial shearing; hence, the maximum volumetric strain is lower.Compared to RDPS, RLPS in Fig.5 show a more significant contraction at the initial stage,and a notable difference in maximum volumetric strain for different morphological gene mutation degrees.Specifically, a more extensive scaling factor (like χc= 1.5) leads to a more prominent maximum volumetric strain than a smaller scaling factor (like χc= 0).The reason is that particles in RLPS have more space to translate and rotate than in RDPS.It is concluded that realistic particle morphology shows more notable effects on the evolution of volumetric strain/axial strain in RLPS than in RDPS.

Volumetric strain curves of both RDPS and RLPS agree with the primary characteristic of real granular materials upon triaxial shearing (Mitchell and Soga, 2005).The effects of different stress conditions and morphological gene mutation degrees on the evolution of volumetric strain/axial strain match well with the previous experimental and numerical works(Xiao et al.,2019;Nie et al.,2020; Sharma et al., 2021; Wu et al., 2021b).Therefore, DEM with realistic particle morphology can feasibly depict the macroscopic volumetric strain evolution.

3.3.Critical state behaviour

Under a constant confining stress, the granular material will ultimately reach a critical state after a large shear-induced volume change.Fig.6 shows the relationship of void ratio e against mean effective stress p.The void ratio experiences a similar variation trend regardless of the stress conditions, the initial packing state,and the morphological gene mutation degrees.In general, RDPS(Fig.6a)shows a slight decrease in e with the increasing p at early stages, followed by a sharp increase in the critical state value.Compared to RDPS, RLPS (Fig.6b) witnesses a similar overall variation trend but a more substantial decrease in void ratio at early stages.Besides, the critical state void ratio exhibits a decreasing trend with augmenting mean effective stress.Particle morphology is found to play a significant role in the critical state void ratio for both RDPS and RLPS.The larger the morphological gene mutation degree, the more excellent the critical void ratio.For all different stress conditions in Fig.7, the void ratio of RDPS and the corresponding RLPS will ultimately arrive at a similar critical state value,which implies that the critical state void ratio is independent of the initial packing state.Additionally, the difference among the void ratio curves of RDPS as well as RLPS is getting narrower with the increase of morphological gene mutation degree.The more considerable morphological gene mutation degree with a rougher particle surface makes the granular assemblage harder to compact for both RDPS and RLPS, primarily upon high confining stresses.

Fig.9.The relationship of principal stress ratio vs.axial strain for morphological gene-mutated RLPS under different confining stresses:(a)χc =0;(b)χc =0.5;(c)χc =1;and(d)χc = 1.5.

Fig.10.Particle rotation and translation upon triaxial shearing.

The evolution curves of the void ratio/mean effective stress and the variation trend of critical state void ratio are consistent with other experimental or numerical results (Been and Jefferies,1985;Mitchell and Soga, 2005; Huang et al., 2014; Qu et al., 2022), indicating that DEM with realistic particle morphology is feasible to depict the critical state behaviour of granular soils.

3.4.Principal stress ratio against axial strain

Figs.8 and 9 show the variations of the principal stress ratio(i.e.σ1/σ3) against axial strain for both morphological gene-mutated RDPS and RLPS under different confining stresses.It is found that RDPS exhibits a similar pattern for the principal stress ratio varying with the axial strain, regardless of the morphological gene mutation degrees and the confining stress conditions.In addition,compared to RDPS,RLPS shows a similar variation trend of principal stress ratio at minor axial strains but a different trend at more significant axial strains.Under the minor confining stress condition,both RDPS and RLPS witness a more prominent principal stress ratio at all stages of axial strain.By adopting realistic particle morphology of granular soils, the simulation results of principal stress ratio/axial strain show a good agreement with previous experimental works on real granular soils (Lade, 1977; Kolymbas and Wu,1990; Xiao et al., 2014; Alshibli and Cil, 2018; Sun et al.,2019), which proves the feasibility and capability of DEM with realistic particle morphology to illustrate the evolution of principal stress ratio/axial strain.

Recently,Qu et al.(2022)reported that the DEM simulations had limited capabilities in reproducing some important behaviours of real granular soils,such as the relationship of principal stress ratio/axial strain.They found significant intersection points in those principal stress ratio/axial strain curves under different confining stress conditions.However,such an abnormality is not found in our results.This abnormal phenomenon might stem from the adoption of a simplified particle shape(i.e.simple two-ball clumps)or rolling resistance model (sphere with Wensrich and Katterfeld (2012)’s model,or even sphere without any rolling resistance)in their study.By comparing the DEM results from simulations adopting spheres and several kinds of clumps with different simplified shapes,Zhou et al.(2013) found that samples with different particle shapes exhibited different kinds of particle anti-rotation mechanisms,leading to notable difference in the localisation of particle rotation and shear strain behaviours.They concluded that the particle shape effects for modelling the realistic mechanical behaviour of irregularly shaped granular materials could not be fully replaced by the rheology-type rolling resistance model.In the current study, the realistic particle morphology of granular soils has been approximated to the best possible degree.It is believed that any such inconsistency between the DEM results and real granular soil behaviour could be minimised.

3.5.Particle kinematics

Particles within a given sample can translate and rotate during loading, as shown in Fig.10.In this study, the particle displacements and particle rotations for any load increment were quantified based on the particle locations and orientations of the samples at the beginning and the final states of the shear increment.The particle kinematics were calculated for each shear strain increment of 1%according to a methodology documented in Cheng and Wang(2018a).

Figs.11 and 12 show the average incremental particle displacement of morphological gene-mutated RDPS and RLPS,respectively.It is found that RDPS shows a similar variation of the average incremental particle displacement, regardless of the morphological gene mutation degrees and the stress conditions.In contrast, RLPS show a notable difference of average incremental particle displacement for samples with different morphological gene mutation degrees.The difference is getting narrower with the increase of confining stress.Generally, a more extensive scaling factor with a rougher particle surface leads to a higher average incremental particle displacement than a smaller scaling factor with a smoother particle surface.The average incremental particle displacement evolution matches well with the volumetric strain evolution (Figs.4 and 5).

Fig.13 shows the average incremental particle rotation of morphological gene-mutated RDPS under different stress conditions.For all the samples, the average particle rotation increases sharply to the peak at an axial strain of about 10% and then fluctuates until approaching the critical state.Besides, the average incremental particle rotation at the axial strain of 1%decreases with increasing confining stress.Compared to RDPS, RLPS experiences notably different average incremental particle rotations at the early shear stage for the samples with different morphological gene mutation degrees, as shown in Fig.14.The rougher the particle surface,the smaller the average incremental rotation.This particle morphology-induced difference is getting narrower with the increase of confining stress.Similar to RDPS, RLPS experiences a nearly unified variation of average incremental particle rotations at large axial strains, regardless of the morphological gene mutation degrees and the confining stress conditions.Thus, particle morphology at LR-level shows more important effects on the particle rotation at small axial strains than at large axial strains.

Fig.12.Average incremental particle displacement of morphological gene-mutated RLPS under different stress conditions:(a)σ3 =0.3 MPa;(b)σ3 =0.5 MPa;(c)σ3 =1 MPa;and(d) σ3 = 1.5 MPa.

The particle translation and rotation features agree with the basic features of real granular soils and match well with other macroscopic behaviours, such as the volumetric strain evolution and the deviatoric stress responses.Therefore, DEM with realistic particle morphology is considered suitable to depict the particle kinematics of granular materials upon triaxial shearing.

3.6.Shear band evolution

The shear band evolution,an intrinsic characteristic of granular materials,can be treated as an instability problem from an initially homogeneous mode to a localised mode (Rudnicki and Rice,1975;Rice,1976; Wang et al., 2007a, b).According to Cheng and Wang(2019), the shear strainsof any point within the granular sample is defined as

whereij（i,j=x,y,z）is the strain tensor of the tetrahedron element formed by four neighbouring particles.Fig.15 shows the shear band evolution of RDPS with χc= 1 under different confining stress conditions.Note that the incremental shear strain field of different axial strain ranges,i.e.0%-2%,2%-5%,5%-10%,10%-15%and 15%-20%,is determined.No clear evidence of strain localisation is found in all samples at 0%-2% axial strain range, regardless of the confining stress condition adopted.At the axial strain range of 2%-5%, strain localisation is observed in all samples, and the higher confining stress condition witnesses a more concentrated shear strain localisation.With the ongoing deformation of the sample,the strain localisation is further developed, and samples are failed through a single shear band.The shear band evolution results match the experimental works of real granular soils (Desrues and Andò, 2015; Cheng and Wang, 2019).

Furthermore,Figs.16 and 17 illustrate the shear band evolution of morphological gene-mutated RDPS and RLPS under different confining stress conditions, respectively.Note that only the shear strains of the samples at the axial strains of 10%-15%and 15%-20%are presented due to the limited space.It is found that: (1) RDPS generally show a more concentrated and developed shear band than RLPS;(2)rougher particle surface with higher confining stress conditions will lead to a more concentrated shear band; and (3)higher confining stress condition always witnesses a more inclined sample geometry, and this incline angle is alleviated with the increase of morphological gene mutation degree.

Fig.13.Average incremental particle rotation of morphological gene-mutated RDPS under different stress conditions: (a)σ3 =0.3 MPa;(b)σ3 =0.5 MPa;(c)σ3 =1 MPa;and(d)σ3 = 1.5 MPa.

4.Conclusions

This study utilised a novel CT-based DEM modelling technique to generate morphological gene-mutated samples to perform numerical triaxial testing.This technique included particle morphology data acquisition from X-ray μCT, particle morphology reconstruction and gene mutation by SH-PCA, and DEM simulations.Thirty three DEM samples in three groups were reconstructed and simulated.With different initial packing densities,morphological gene mutation degrees and confining stresses, the capability of DEM in reproducing the macro-micro mechanical behaviour of granular materials with realistic particle morphology was investigated.The main conclusions are summarised as follows:

(1) DEM with realistic particle morphology is capable of reproducing the macro-micro mechanical behaviours of granular soils,including the relationships of the deviatoric stress/axial strain and the volumetric strain/axial strain,the critical state behaviour, the particle kinematics, and the shear band evolution.

(2) The effects of particle morphology on the deviatoric stress,volumetric strain, particle translation and particle rotation depend on the initial packing state and the confining stress condition.It is found that particle morphology has a more critical effect on the looser samples than the denser samples.For the samples with different particle morphologies,higher confining stress leads to a more considerable difference in peak-state deviatoric stress, yet a more minor difference in volumetric strain and incremental particle translation and rotation.

(3) RDPS has a more concentrated and developed shear band than RLPS.It is found that samples with rougher particle surfaces under higher confining stresses have a more focused shear band.Meanwhile, samples subjected to higher confining stresses always witness a more inclined sample geometry, which alleviates this incline angle with the increase of morphological gene mutation degree.

These results highlight the critical role of particle morphology in the micro-macro mechanical behaviours of granular soils.With realistic particle morphology, DEM can model the micro-macro mechanical responses of granular soils under different stress conditions and, more importantly, reproduce many behaviours that DEM cannot realise with spherical or simplified particles.The following limitations are acknowledged:(1)only gene mutation at LR-level was considered to simplify the morphological variance;and (2) all the tested samples share the same particle size distribution.It is worth noting that gene mutation at different length scales and particle size distributions might also affect soil behaviours.Future work will be undertaken to handle the coupling effects of particle size distribution and morphological gene mutation.

Fig.14.Average incremental particle rotation of morphological gene-mutated RLPS under different stress conditions: (a)σ3 = 0.3 MPa;(b)σ3 = 0.5 MPa;(c)σ3 = 1 MPa;and (d)σ3 = 1.5 MPa.

Fig.15.Shear band evolution of RDPS with different confining stresses.

Fig.16.Shear band evolution of morphological gene-mutated RDPS at axial strains of 10%-15%.

Fig.17.Shear band evolution of morphological gene-mutated RLPS at axial strains of 15%-20%.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

This study was supported by the General Research Fund from the Research Grant Council of the Hong Kong SAR, China (Grant Nos.CityU 11201020 and CityU 11207321), the National Science Foundation of China (Grant No.42207185), the Contract Research Project from the Geotechnical Engineering Office of the Civil Engineering Development Department of Hong Kong SAR, China(Project Ref.No.CEDD STD-30-2030-1-12R), and the BL13W beamline of Shanghai Synchrotron Radiation Facility (SSRF).The first author acknowledges the financial support from Hong Kong PhD fellowship scheme (HKPFS).