An overview of particle-based numerical manifold method and its application to dynamic rock fracturing

Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

Contents lists available at ScienceDirect

Journal of Rock Mechanics and Geotechnical Engineering journal homepage: www.rockgeotech.org

Review

An overview of particle-based numerical manifold method and its application to dynamic rock fracturing Xing Li a, *, Jian Zhao b, ** a b

School of Civil Engineering, Southeast University, Nanjing, 211189, China Department of Civil Engineering, Monash University, Clayton, VIC, 3800, Australia

a r t i c l e i n f o

a b s t r a c t

Article history: Received 21 August 2018 Received in revised form 7 January 2019 Accepted 18 February 2019 Available online 13 April 2019

This review summarizes the development of particle-based numerical manifold method (PNMM) and its applications to rock dynamics. The fundamental principle of numerical manifold method (NMM) is first briefly introduced. Then, the history of the newly developed PNMM is given. Basic idea of PNMM and its simulation procedure are presented. Considering that PNMM could be regarded as an NMM-based model, a comparison of PNMM and NMM is discussed from several points of view in this paper. Besides, accomplished applications of PNMM to the dynamic rock fracturing are also reviewed. Finally, some recommendations are provided for the future work of PNMM.
2019 Institute of Rock and Soil Mechanics, Chinese Academy of Sciences. Production and hosting by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/ licenses/by-nc-nd/4.0/).

Keywords: Dynamic rock fracturing Numerical modeling Numerical manifold method (NMM) Particle-based numerical manifold method (PNMM)

1. Introduction Dynamic rock fracturing involves the topics of crack nucleation, fracture propagation, rock fragmentation and rock post-failure behavior. Results of a series of experimental tests (Masuda et al., 1987; Lajtai et al., 1991) have revealed that the fracture pattern and mechanical properties of rock materials are affected by the dynamic strain rate. However, the underlying mechanism of the rate-dependent behavior is still unclear. The microstructure of rock material is recently considered as one of the influence factors leading to this phenomenon. When the mechanical behavior on microscopic scale is concerned, analytical methods are likely invalid due to the complex microstructure of rock materials, and experimental approaches will be limited as existing facilities are not sensitive enough to detect the failure process under high loading rates, especially those occurred at the inner part of rock samples. Under such a circumstance, numerical methods could be an alternative tool to study the mechanism of dynamic effect on rock materials at the microscopic scale.

* Corresponding author. ** Corresponding author. E-mail addresses: [email protected] (X. Li), [email protected] (J. Zhao). Peer review under responsibility of Institute of Rock and Soil Mechanics, Chinese Academy of Sciences.

Due to the complexity of both rock materials and rock engineering problems, rock is difficult to be modeled, compared with other solid materials. To obtain better results, a large number of numerical models have been developed, extended, and applied in this field, including in-house software, commercial software, and open source code. Existing numerical methods for rock dynamics are usually classified into three categories (Jing, 2003): continuous methods, discontinuous methods, and coupled methods. Continuous methods are suitable for those problems whose system is a continuum with infinite degree of freedom (DOF). The behavior of such system is dominated by the governing differential equation of the problem and the continuity conditions at the interfaces between adjacent elements. The well-developed continuous methods include the finite element method (FEM) (Clough, 1960), finite difference method (FDM) (Narasimhan and Witherspoon, 1976), realistic failure process analysis (RFPA) (Tang, 1997; Liang et al., 2004; Zhu, 2008; Wang et al., 2014), element-free Galerkin (EFG) method (Belytschko et al., 1994, 2000) and the related cracking particle method (Rabczuk and Belytschko, 2004, 2007; Rabczuk and Areias, 2006; Rabczuk et al., 2010a), peridynamics (PD) and dual-horizon peridynamics (DH-PD) (Silling, 2000; Ren et al., 2016, 2017; Rabczuk and Ren, 2017), partition of unity (PU) method (Babuska and Melenk, 1997; Moës et al., 1999; Sukumar et al., 2000; Rabczuk and Zi, 2007; Rabczuk et al., 2010b; Amiri et al., 2014), smoothed particle hydrodynamics (SPH) method (Gingold and Monaghan, 1977; Lucy, 1977; Fakhimi and Lanari, 2014; Zhou

https://doi.org/10.1016/j.jrmge.2019.02.003 1674-7755
2019 Institute of Rock and Soil Mechanics, Chinese Academy of Sciences. Production and hosting by Elsevier B.V. This is an open access article under the CC BYNC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

et al., 2015), and other meshless methods (Liszka et al., 1996; Liu et al., 1997; Atluri and Zhu, 1998; Ma et al., 2015). Some reviews on continuous methods can be found in Belytschko et al. (1996, 2009), Li and Liu (2002), Nguyen et al. (2008), Liu and Liu (2010), and Monaghan (2012). This type of method is mostly used for rock masses without fractures. It can also be used for the rock mass with a few or many fractures, the behavior of which is established through equivalent properties or other techniques, e.g. inserting cohesive elements along the path of fractures (Vocialta et al., 2017) and enriching approximation functions in the element containing fractures (Moës et al., 1999). Discontinuous methods are suitable for those problems whose system is a combination of a finite number of well-defined components. Usually there is no need to discretize such a system, as it has been automatically done. The behavior of such system is dominated by the well-defined inter-relations between adjacent components. This type of method is most suitable for moderately fractured rock masses where the significant fractures are for continuous methods, and/or where large displacements of individual blocks are possible. The term ‘discontinuous method’ here indicates the discrete element methods and other discontinuous methods, e.g. the lattice models (Hrennikoff, 1941; Holecek and Moravec, 2006; Zhao and Zhao, 2012; Zhao et al., 2012a). The discrete element methods are a large family of numerical methods. The distinct element method (DEM) is a class of discrete element methods that use an explicit time-domain integration scheme to solve the equations of motion for rigid or deformable discrete bodies with deformable contacts. Since DEM is almost the most well-known explicit discrete element method, researchers tend to use the terms ‘DEM’ and ‘discrete element method’ interchangeably. DEM treats the simulated material as an assembly of separate particles or blocks. The particle-based DEM (Cundall and Strack, 1979; Potyondy and Cundall, 2004) is mainly adopted to simulate the granular microstructure of the material through particles with varying diameters. Adjacent particles can contact each other through bonds. The contact is typically assigned with a normal and shear stiffness as well as a friction coefficient. Crack nucleation and fracture propagation are simulated by breaking of bonds. Consequently, blocks of arbitrary shapes can be formed as a result of the simulated fracturing process. There are typically two types of bonds in the particlebased DEM: the contact bond and the parallel bond. The contact bond is an elastic spring with a constant normal and shear stiffness function between the bonded particles, allowing only normal and shear forces to be transmitted. Whereas on the parallel bond, the moment induced by particle rotation is resisted by a set of elastic springs distributed between the bonded particles. One of the major drawbacks of the bonded-particle model is that the straightforward adoption of circular/spherical particles cannot fully capture the behaviors of complex-shaped and highly interlocked grain structures, which is common in hard rocks (Lisjak and Grasselli, 2014). To overcome this limitation, a clumped particle model (Cho et al., 2007) was later proposed. In this model, a group of particles are glued together to behave as a single rigid body. Clumped particles can act like a single particle that has an irregular shape but moves as a deformable and non-breakable body. Apart from sphere, particles in other shapes have also been developed to improve the geometry description of this method, including ellipsoids (Vu-Quoc et al., 2000; Yan et al., 2010; Zheng et al., 2013), superquadrics (Wellmann and Wriggers, 2012), polyhedral (Feng et al., 2012; Smeets et al., 2015), combination of simple shape primitives (Lu and McDowell, 2007; Ferellec and McDowell, 2010; Fang et al.,

685

2015), and mineral-grain-shaped (Li et al., 2018a,b). On the other hand, the block-based DEM discretizes the computational domain into blocks using a finite number of intersecting discontinuities (Cundall, 1971; Zhao et al., 2008; Gui et al., 2016). Each block is internally subdivided for calculation of displacement, strain and stress. Adjacent blocks interact in the normal direction with a finite stiffness together with a tensile strength criterion and in the tangential direction with a tangential stiffness together with a shear strength. In the standard block-based DEM, rock failure can only be captured either in terms of plastic yield or displacements of pre-existing discontinuities. Therefore, new discontinuities could not be driven within the continuum portion of the model, so that discrete fracturing through intact rock could not be simulated (Lisjak and Grasselli, 2014). In contrast to the explicit DEM, the discontinuous deformation analysis (DDA) method is an implicit discrete element method developed for the modeling of static and dynamic behaviors of discrete block systems, originally proposed by Shi and Goodman (1985, 1989). The blocks in DDA is deformable and can be arbitrarily shaped. The kinematic constraints of no tension and no penetration between blocks are imposed by the penalty method, the Lagrange multiplier method, or the augmented Lagrangian method. Different from the explicit DEM, the governing equations in DDA are represented by a global system of linear equations obtained by minimizing the total potential energy of the system, and a global system of equations in matrix form needs to be formed and solved. By differentiating several energy contributions, including block strain energy, contacts between blocks, displacement constraints and external loads, the stiffness matrix of the model could be assembled. After years of development, DDA has been applied to the simulation of various topics, including the rock sliding (Hatzor et al., 2004; Wu et al., 2009; Ning and Zhao, 2013; Jiao et al., 2014; Yang et al., 2014a; Zhang et al., 2015), rock fracturing (Kong and Liu, 2002; Bao and Zhao, 2013; Chen et al., 2013; Tian et al., 2014), wave propagation in rock masses (Jiao et al., 2007; Gu and Zhao, 2009), behavior of jointed rock masses (Lin et al., 1996; Tsesarsky and Talesnick, 2007; He and Zhang, 2015; He et al., 2018), and failure of masonry structures (Kamai and Hatzor, 2008; Jiang et al., 2014). A coupled method, in most cases, is a combination of one continuous method and one discrete method implicitly or explicitly, aiming to combine the characteristics of both methods. The numerical manifold method (NMM) is a typical coupled method for rock mechanics (Shi, 1991, 1995; Ma et al., 2010). It provides a framework to unify continuous and discrete methods by implicitly combining FEM and DDA. Recently, a novel numerical model, named particle-based numerical manifold method (PNMM), has been proposed as an extension of NMM and applied to several problems in the field of rock dynamics. The purpose of this review is to provide a summary of selected coupling numerical method for rock dynamics. Specifically, the NMM and PNMM are considered in this paper. As for other coupling numerical methods, e.g. the finitediscrete element method (FDEM), the reader is suggested to refer to the works (Mahabadi, 2012; Lisjak, 2013; Lisjak and Grasselli, 2014), and its applications to dynamic rock fracturing can be found in Rougier et al. (2014) and An et al. (2017). In this context, the fundamental principle of NMM is first briefly introduced. Recent improvement and applications in the field of rock dynamics are summarized and reviewed. Then, the history of the newly developed PNMM model is given. The basic idea of PNMM and its simulation procedure are also presented. Considering that PNMM could be regarded as an NMM-based model, a comparison of PNMM and NMM is discussed in this paper. Besides, accomplished

686

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

applications of PNMM to the dynamic rock fracturing are also reviewed. Finally, some recommendations are provided for the future work of PNMM. 2. Numerical manifold method 2.1. Fundamental principles The NMM was originally proposed by Shi (1991, 1992, 1995, 1997, 2012) at a series of conferences. NMM is derived based on the finite cover approximation theory and is named after the mathematical notion of manifolds. NMM feasibly simulates fractures by truncating elements and their shape functions, and treats continuum bodies, fractured bodies and assemblage of discrete blocks in a unified form. This characteristic of NMM is mainly due to the adoption of a dualcover-system, namely the mathematical cover and the physical cover. The mathematical cover is independent of the shape of modeling domain but covers all the spaces that the modeling domain may occupy. It is usually generated from a uniform FEM mesh, from which the finite elements that share the same node form a mathematical cover. Mathematical covers are common to be overlapping. The mathematical covers with different geometric patterns could be generated from several kinds of typical FEM meshes, as illustrated in Fig. 1, which are supposed to have a slight influence on simulation results depending on the analytical solution of concerned problems (Zhang et al., 2010). On the other hand, the physical covers are the intersection of mathematical covers and the modeling domain (Ma et al., 2010). Specifically, a physical patch is formed by cutting a mathematical cover with discontinuities, such as physical boundaries, material interfaces, and fractures. Then, the union of all physical patches constitutes a physical cover. The physical cover system is to connect the uniform mathematical cover system with the arbitrary modeling domain, providing local approximation function in NMM. The overlap of neighboring physical covers is called a manifold element. A manifold element is the basic computation unit in NMM. The advantage of adopting two cover systems for NMM is that arbitrary boundaries and internal physical features in the physical domain are allowed without meshes conforming to them. Therefore, the meshing task in NMM is greatly simplified and fracturing process can be modeled without remeshing. At the same time, NMM model fractures straightforwardly by splitting each physical cover that is completely cut by the fractures into several separate covers, assigning each cover an independent local function. By such a way, complex cracks with arbitrary number of branches can be modeled in an exactly way as the modeling of a single crack (see Fig. 2). In this figure, a branched fracture separates four physical covers into twelve physical covers, considering that the four mathematical covers in the figure will lead to four physical covers without the branched fracture. Besides, NMM implements a novel simplex integration method proposed by Shi (1996) for both NMM and DDA. This integration method is able to conveniently evaluate the weak form integration

Fig. 1. Typical patterns of the mathematical cover (the polygon indicated by red lines).

over elements intersected by internal discontinuities and/or external boundaries. For the manifold element of a general shape, the simple integration scheme will first partition the element into several sub-triangles, the same process as the element partitioning in extended FEM (XFEM), and then conduct the quadrature on each sub-triangle. Another distinct feature of NMM is that the frictional contact boundary conditions between the two sides of a crack or two discrete blocks can be accurately satisfied, due to the contact logic in NMM. It has been proven that FEM is a special case of NMM when the following conditions are satisfied (Ma et al., 2010): (1) Mathematical covers in NMM are generated from a finite element mesh; (2) Weight functions defined on mathematical covers are finite element shape functions; (3) Cover functions defined on physical covers are constants; (4) Physical features including internal discontinuities (e.g. cracks and material interfaces) and external boundaries do not intersect manifold elements. Although NMM can be converted to FEM, its distinct features in dealing with discontinuous problems cannot be covered in FEM. XFEM/GFEM (generalized FEM) is another numerical method that is capable of dealing with discontinuities. Recently, the concept of enriched functions in XFEM/GFEM has been integrated into NMM to simulate complex crack problems (Ma et al., 2009; An, 2010). It should be noticed that NMM is able to tackle cracks without enriched functions in most cases, and NMM furnished with enriched functions can tackle more complex crack problems as compared with XFEM/GFEM. Meanwhile, since NMM is developed on the basis of DDA, it preserves all the characteristics of discrete element modeling such as the kinematics constraints, contact detection, and modeling from DDA (Ma et al., 2010). However, the DOF in NMM is usually much higher than DDA as there is more than one manifold element in each block in most cases. The benefit at this cost is that NMM provides more accurate displacement and stress fields in blocks than DDA. It can be concluded that if every discrete block is a manifold element with linear displacement field, then NMM will degenerate exactly into DDA. 2.2. Recent development and applications In the original version of NMM by Shi (1991, 1992, 1995, 1997, 2012), the simplest triangular manifold element with constant cover function was adopted for two-dimensional (2D) issues, lacking criteria for crack initiation and propagation. Since then, various developments as well as applications have been made in the following decades. Shyu and Salami (1995) implemented quadrilateral isoparametric elements in NMM. Chiou et al. (2002) studied the mixed-mode fracture propagation by combining NMM with the virtual crack extension method. Cheng et al. (2002) incorporated Wilson nonconforming elements in NMM. Chen et al. (1998) developed a higher-order NMM with high-order cover functions, and then Su et al. (2003) proposed a subroutine in the commercial software Mathematica to automatically produce expressions for high-order NMM. Recently, Ghasemzadeh et al. (2014) proposed a high-order NMM for dynamic problems. Terada et al. (2003) introduced the finite cover method (FCM) as an alias of NMM. Lin et al. (2005) developed the formulations of threedimensional (3D) NMM with high-order cover functions and proposed a fast simplex integration method based on special matrix operations, without considering the linear dependence problem. Later, Cheng and Zhang (2008) proposed the 3D NMM with tetrahedron and hexahedron elements and derived the basic matrices

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

687

Fig. 2. The mathematical cover (MC) and physical cover (PC) cut by a branched crack (after Li et al., 2017). Patches are filled with different colors to indicate the relation between PCs and MCs.

Fig. 3. Road protection from rock slope failure simulated by 3D NMM (He et al., 2013).

for equilibrium equations. He and Ma (2010) also proposed a 3D NMM based on tetrahedron elements (Fig. 3). An et al. (2011) investigated the linear dependence problem of NMM approximation space using finite element covers and polynomial local

functions at both elemental and global levels. Yang and Zheng (2016) and Xu et al. (2017) recently proposed a high-order local approximation function and eliminated the linear dependence on the triangular and quadrilateral mathematical meshes, respectively.

688

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

Their model has a higher order of global approximations, better accuracy and continuous nodal stress, and was then applied to the fracture analysis as well as free and forced vibration analysis of solids (Yang et al., 2016a,b,c; Yang and Zheng, 2016). Cai et al. (2013) proposed a generalized and efficient cover generation procedure, which is applicable for dealing with interfaces, inclusions and discontinuities with complex geometry. Zheng and Xu (2014) proposed strategies for several specific issues that NMM may encounter in simulation of crack propagation, including the rank deficiency induced by the high-order cover functions, the integrals with singularity of 1/r (r indicates the distance to the crack tip), and the kinked cracks. A flat-top partition of unity (PU)-based NMM was proposed by He et al. (2015) to alleviate the linear dependence difficulty of traditional NMM. Recently, Yang et al. (2017) and Zheng and Yang (2017) developed a mass lumping scheme for NMM and other PU methods, which shows a more excellent behavior in dynamic analysis compared with the conventional consistent mass matrix since solving large-scale simultaneous algebraic equations can be avoided. For simulating fractures, Li et al. (2005a,b) and Gao and Cheng (2010) developed an enriched meshless manifold method for 2D crack modeling. Meanwhile, Ma et al. (2009) and An (2010) incorporated the enriched functions of XFEM/GFEM into standard NMM to simulate complex cracks. Zhao et al. (2012a) coupled NMM with the distinct lattice spring model to simulate the dynamic failure of rock masses. Wu and Wong (2013a,b, 2014), Wu et al. (2013) and Wong and Wu (2014) implemented NMM and studied rock fracturing in many different circumstances (Fig. 4). Yang et al. (2014b)

refined the mathematical covers near crack tips. Ning et al. (2011), An et al. (2014) and Zheng et al. (2014a) analyzed the stability of rock slopes with 2D NMM, while He et al. (2013) and Liu et al. (2017) extended the analysis to 3D cases. Zheng et al. (2014b) adopted the influencing domains of nodes in the moving least squares (MLS) interpolation when constructing the mathematical cover and proposed the MLS-based NMM to better address the problems with large movement. Then, they adopted the MLS-based NMM and proposed a so-called projection-contraction algorithm to tackle the complementarity problem in the static growth of multiple fractures (Zheng et al., 2015). Recently, Liu et al. (2018a,b) studied the micro/macro failure of rock by developing a Voronoi element based-NMM. For stress wave propagation issues, Zhao et al. (2014) improved the performance of NMM in the simulation of wave propagation in rock masses by importing Newmark system equations, edge-to-edge contact scheme and non-reflection boundary condition. Wu and Fan (2014) developed a timedependent absorbing boundary condition for wave propagation problems. Wei et al. (2018) also proposed some new boundary conditions for NMM to study the seismic response of geomechanics problems. Yang et al. (2018) developed boundary settings for the seismic dynamic response analysis of rock masses in NMM. Moreover, Zhang et al. (2014, 2017, 2018) studied the thermo-mechanical fractures and transient heat conduction problems using NMM in different solids. NMM has been a promising numerical method for rock mechanics in recent years due to its great advantages in theoretical basement for continua-discontinua analysis. However, as a quasi-

Fig. 4. Simulation of the failure of rock samples using NMM: (a) Failure of Brazilian tensile disc test in pre-fractured sample (Wong and Wu, 2014); and (b) Cracking of rock containing two inclusions (Wu and Wong, 2013b).

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

block-based model, its geometrical and topological operations, and then its contact detection and operation, can be complicated, especially in 3D cases. Furthermore, NMM is generally suitable for the analysis in macroscopic scale, but not in microscopic scale. 3. Particle-based numerical manifold method The PNMM is a newly proposed extension of NMM for rock dynamics. The basic idea of PNMM is to incorporate the particle concept into NMM. The motivation of this development is to simplify the geometrical Boolean operations and contact operations in NMM, as well as to better understand the heterogeneity of rock materials and rock fractures. 3.1. Particle manifold method The model originally proposed by Zhao et al. (2012a,b) and Sun et al. (2013, 2014) was entitled particle manifold method (PMM) during its early stage. PMM replaced the polyhedron representation of physical cover with a particle representation. The basic idea of this improvement is to reduce the complexity of the topological operations between polygons in NMM. Considering that the physical cover could be in general shapes due to the intersection between the mathematical cover and discontinuities (e.g. boundaries, material interfaces, and fractures), topological operations in NMM are usually complicated and time-consuming in practice, especially in 3D cases. For example, in Fig. 5, Zone I in polygon representation is the result of a Boolean operation ‘AND’ between the triangular simulation domain and the rectangular mathematical cover M. Using the particle representation, Zone I is described by a group of red particles located in M. Therefore, in generation of the physical cover I (same as Zone I in geometry), the topological operation between two polygons is simplified to the operation between several points and a polygon. This polygon-to-polygon topological operation, sometimes called block cutting procedure in NMM, could be more complicated when fractures are taken into consideration. The particle manifold element, as the basic element in PMM, was defined as the combination of the mathematical cover and particles. A relationship termed ‘link’ was defined to describe the continuous state between neighboring particles. PMM adopted a particle simplex integration scheme for numerical calculation. Although the integration scheme is analytical, there could still exist slight accuracy loss due to the particle representation of polygonal domains (Sun, 2012). A virtual integration model and a realistic integration model were developed for homogeneous and inhomogeneous materials, respectively. The most noticeable improvement of PMM compared to NMM was that the polygon-to-polygon contact was simplified to particle-to-particle contact, since the contact between particles was easier to be detected and computed. The contact forces between particles as well as between particle

and rigid wall were developed. Only normal contact force was considered in the model. The maximum tensile stress criterion and Mohr-Coulomb criterion were adopted in PMM to study the brittle failure of rock. The link was taken as the basic failure component in the model, and the directivity of the failure was discussed. Application of PMM includes a preliminary study of the stress wave propagation across joints, the spalling of rock bars (Fig. 6), and some other fracturing problems. Furthermore, a graphics processing unit (GPU) parallelization and a coupling of PMM and NMM were also performed, attempting to further improve the efficiency of computation.

3.2. Fundamental principles An improvement to PMM was later performed and consequently the numerical model PNMM was proposed. This section briefly presents some fundamental principles of PNMM. More details of the components, formulations and implementation of the model can be found in Li (2017) and Li et al. (2017). PNMM inherits the dual-layer-cover system from NMM, i.e. the mathematical cover and physical cover. On each mathematical cover, a weight function 4i is defined to combine the interpolation functions from different physical covers, which satisfies

X 4i ðxÞ ¼ 1 i

9 > = (1)

4i ðxÞ  0 ðcx˛Pi Þ > ; 4i ðxÞ ¼ 0 ðcx;Pi Þ

where Pi is the ith mathematical cover. The physical cover is to provide the local approximation function by defining the cover function ui:

ui ¼

n X

Tij ðxÞdij ¼ Tdi

(2)

j¼1

where n is the number of DOFs, T is the basis of the cover function, and di is the vector of DOFs on the ith physical cover. Since a manifold element is the overlap of neighboring physical covers, its displacement field is generated by combing the cover functions of related physical covers using weight functions as

ue ðxÞ ¼

m X

4i ðxÞui ðxÞ ¼

i¼1

m X

N i ðxÞdi ¼ N e de

(3)

i¼1

where ue is the displacement field of the manifold element e; m indicates the number of physical covers; Ni and Ne are the shape functions on the ith physical cover and manifold element e, respectively; and de is the vector of DOFs on the manifold element e. The generation of manifold elements, as the result of the duallayer-cover system, forms the first discretization in PNMM. With this discretization, the global equation of PNMM for dynamic analysis can be obtained in matrix form as

Kd þ M d€ ¼ F

Fig. 5. Physical cover in PMM: Replacing the polygon representation with a particle representation (Sun et al., 2013).

689

(4)

where K is the global stiffness matrix, d is the global vector of unknowns, M is the global mass matrix, and F is the global vector of equivalent loads. This governing equation is the same as that in NMM, since particles have not yet been introduced by far. An extra level of discretization is incorporated by dividing each manifold element into a group of particles, as illustrated in Fig. 7. As the governing equation and DOFs are defined on manifold elements, introduction of the particles does not bring in additional

690

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

Fig. 6. The spalling of rock bar under impulse stress wave predicted by PMM (Sun et al., 2013).

Fig. 7. The dual-level discretization in PNMM (Li et al., 2018c).

DOFs and does not affect the form of the governing equation. However, the material properties, boundary conditions and body forces are all assigned to individual particles independently. Particles in a manifold element are possible to have varying material properties and boundary conditions. For calculation of the manifold element matrices, a particle integration scheme was proposed by taking the particle centroids and areas as integration points and weights, respectively, as

Z f ðxÞdA ¼ A

p p h X X   i   c  f xci r 2i f xi Ai ¼ p

i¼1

(5)

i¼1

sL ¼ ðs1 þ s2 Þ=2



εL ¼ ðε1 þ ε2 Þ=2

where f(x) is the integrand, p is the number of particles, xci is the coordinates of the centroid of the ith particle, Ai is the area of the ith particle, and ri is the radius of the ith particle. The particle integration scheme is a numerical integration technique, and it is suitable for both convex and concave integration domains and any integrands, including polynomial, exponential, and trigonometric functions (Li, 2017). Using the particle integration scheme, the global matrices in the governing equation could be assembled. Then, the mechanical fields of manifold elements are obtained by solving the equation, and the mechanical results of each particle can be simply derived as

      ui ¼ ue xci ; εi ¼ εe xci ; si ¼ se xci

conditions are satisfied at the same time: (1) the two particles are next to each other geometrically; (2) the two particles belong to the same block/object; and (3) there are no micro/macro fractures between the two particles. In implementation, the Munjiza-NBS contact detection algorithm (Munjiza and Andrews, 1998) is utilized for the generation of links. A link is the basic failure component in PNMM, as that in PMM. The stress and strain of a link are taken as the average stress and strain of the two particles it connects:

(6)

where ui, εi and si are the displacement, strain and stress of the ith particle in the manifold element e, respectively; and ue, εe and se are respectively the displacement, strain and stress of the manifold element e. Among particles in a continuum, links are defined as the continuous status between neighboring particles. A link is generated between any two particles as long as the following three

(7)

where sL and εL are the stress and strain tensors of the link, respectively; s1 and ε1 are the stress and strain tensors of the first particle, respectively; and s2 and ε2 are the stress and strain tensors of the second particle, respectively. For rock dynamics issues, the Johnson-Holmquist-Beissel (JHB) model (Johnson et al., 2003) is chosen as the strength criterion to reproduce the rate-dependent behaviors of rock:

sc ¼ s0 ð1 þ C ln_ε*Þ

(8)

where sc is the dynamic strength at the strain rate of ε_ *; s0 is the strength at ε_ * ¼ 1; and C is the dimensionless strain rate constant. More details of the implementation of JHB model in PNMM and the determination process of JHB parameters for rock materials can be found in Ma and An (2008) and Li et al. (2017). The strength criterion is applied on each link, and once the stress of a link reaches the dynamic strength from the criterion, the link is supposed to be failed, representing that a microcrack is generated. In a continuous model, it could be proven that any two particles within the model could be connected by either one or more links. When there are failed links generated in a manifold element at a time step in

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

calculation, the topological relation between particles would be redetected. And if disconnected particles are found to separate particles into more than one group, the manifold element would be subdivided accordingly, representing that the fracture has propagated across the manifold element. A cover function enrichment is implemented around the fracture tip to consider the singular fields. The contact operation in PMM is kept, i.e. the particle-toparticle and particle-to-plane contacts are defined respectively. Particles on the surface of different continuums/blocks may form a pair of particle-to-particle contacts when they are adjacent to each other. The particle-to-plane contacts are to simulate the interaction between deformable bodies and rigid static/moving walls or infinite planes. The difference between a contact and a link is that there is no force defined on the link whereas the contact force is applied. In the current version of PNMM, only normal contact force is defined and conducted. Although the direction of a pair of contact is determined by the position of two particles at particle scale, PNMM accumulates pairs of contacts to simulate the complex block-to-block contact at macro-scale, and the contact force between blocks is the sum of contact forces between particles, as illustrated in Fig. 8. PNMM implements the same open-close iteration procedure as NMM (Liu et al., 2018c) to prevent the penetration between blocks. It should be noted that although only normal contact is considered in PNMM, the model still gains the ability to simulate shear failure, since the failure is defined on the link instead of the contact. The flowchart of PNMM is illustrated in Fig. 9 and the calculation steps are concluded as follows:

(4)

(5)

(6) (7)

(8) (9) (10)

(11) (12)

(1) Input the model and generate mathematical covers. The model is created in other software and read by PNMM from a text file. The text file imported into PNMM includes the model geometry information, material properties, and some other setting variables. (2) Generate manifold elements. The generation of manifold elements is based on the topological operation between the imported model and mathematical covers. (3) Generate particles. Each manifold element is first divided into several inner triangular elements. Then, one particle is generated in each inner triangular element. The center of the particle is the centroid of the inner triangular element. The area of the particle is equal to the area of the inner triangular element. In such a way, the sum of particle areas remains equal to the area of the polygonal element. It is common to find adjacent particles overlapped due to the limitation of geometry. However, particles in a same manifold element are mainly used to carry material parameters and serve as the numerical integration points. There is no contact force

(13) (14)

(15)

(16)

(17)

691

defined between adjacent particles, which are distinctly different from other particle-based methods. Therefore, particles in a same element could be regarded as point collocation, and the overlap between particles should have no physical meaning (Li et al., 2017). Generate links. The Munjiza-NBS contact detection algorithm (Munjiza and Andrews, 1998) is implemented to detect the pairs of neighboring particles in PNMM and the links are generated accordingly. Details of this procedure are available in Li (2017). Generate blocks. This is to detect how many blocks are involved in the model. A pair of neighboring particles is connected by a link. A block is constituted of a group of particles, where any two of them could be connected through numbers of links. A so-called seed filling method is adopted to accomplish this task (Sun, 2012; Li, 2017). Detect pairs of contacts. Apply body forces, external loads, and boundary conditions. These parameters are defined on particles and applied to the element matrices. Generate element matrices using the particle integration scheme. Assemble global matrices. Solve global equation. The preconditioned conjugate gradient (PCG) method is implemented to solve the global equation, which is a system of linear equations. Derive DOF results on manifold elements from the global vector of unknowns, which are solved in the previous step. Obtain particle mechanical results from the mechanical fields of a manifold element straightforwardly. Obtain link mechanical results from the two particles it connects as the average value. Apply the JHB criterion. The JHB model is implemented in PNMM and applied on each link at every time step. Once the JHB criterion is met, the status of the link will be set to be failed. A group of material constants in the criterion need to be determined before the simulation. Remove failed links, and redetect connectivity between particles. This task is accomplished by the seed filling method as well. Generate new manifold element and re-group particles. Additional manifold elements need to be generated when fractures have propagated through an element. Post-process. PNMM writes the simulation results into an external text file, which will be later processed by other postprocess software.

3.3. Applications

Fig. 8. From particle-to-particle contact to block-to-block contact.

The application of PNMM mainly focuses on the dynamic fracturing of rock. The proposed PNMM was first validated by modeling some fundamental problems, including the cantilever beam, Brazilian disc test, and rock blasting under dynamic loading (Li et al., 2017). Moreover, stress wave propagation is an important issue in studying the dynamic fracturing of rock. The example illustrated in Fig. 10 presents the wave propagation in a rock cavern which is semi-circle in shape. The propagated wave was induced by an explosion load uniformly applied on the surface of the rock cavern. Simulation results showed that the wave arrives at the upper bound of the model at the time of 9.72 ms, being in good agreement with the theoretical result of 9.61 ms. Three monitoring points were set in the model. The time-displacement results of the three monitoring points obtained by PNMM agreed well with those simulated by ABAQUS (Dassault Systèmes, 2017). Later, the spalling

692

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

Fig. 9. The flowchart of PNMM.

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

693

Fig. 10. The stress wave propagation in a rock cavern: (a) Numerical model; (b) Simulated vertical displacement results at t ¼ 10 ms (unit: m); and (c) Displacement histories at monitoring points e comparison with FEM results.

of rock bars under different loading rates and different types of triangular compression waves were simulated, as shown in Fig. 11 (Li et al., 2018d). It was found that the spalling pattern is closely related to the type of the wave applied on the bar. When a sawtooth or reverse sawtooth wave was imposed, the second spall fracture would possibly be created on the first spall piece or the remaining part of the bar; whereas in the case where a symmetrical triangular

wave is applied, a spalling zone rather than a single spalling fracture will be generated. The spalling process in plate impact test with different impact velocities and different flyer sizes were also studied (Li, 2017; Li et al., 2018d). Two flyers were adopted in the test, a long flyer with the same length as the target and a short flyer with half the length of the target (Fig. 12). In the long flyer test, spalling fractures were generated by the reflected stress in a narrow

694

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

Fig. 11. The spalling of rock bar under different types of triangular compression waves: (a) Reverse sawtooth compression wave; (b) Sawtooth compression wave; and (c) Symmetrical triangular compression wave (Li et al., 2018d).

Fig. 12. The plate impact test with (a) long flyer (Li et al., 2018d) and (b) short flyer.

area near the center of the target. The length of the fractured area was almost the same as the length of the target. Then, spalling fractures were extended to the left and right boundaries of the target, and the target was consequently separated into two parts, where a marked difference between the flying speeds existed. In the short flyer test, spalling fractures were generated by the reflected stress in a narrow area near the center of the target as well; however, the fracturing zone was shorter than the flyer. Then, some longitudinal fractures were generated from the spalling zone to the lower bound of the target, creating several spalling pieces from the target. The speed of these spalling pieces was non-uniform at their creation but became uniform with the increase of the simulation time. There were also some fractures originating from the upper bound of the target, due to the mismatch between the length of the flyer and target. Based on the spalling phenomenon in rock, the rockburst of tunnel induced by static in situ stress and far-field dynamic disturbance was further simulated in Li et al. (2018a). This simulation was divided into two stages. In the first stage, static in situ stresses were imposed on the model to obtain the initial stress state. Then, a dynamic disturbance was uniformly imposed on the left boundary of the model to trigger the spalling and

rockburst. The applied dynamic disturbance was in the form of an impulse function. Numerical results revealed that if the vertical and horizontal in situ stresses are close to each other, the surface of tunnel will be under the state of compression. A compressive in situ stress field helped to offset the reflected tensile stress wave; however, spalling could still occur when the peak value of the dynamic disturbance was high enough. Spalling fractures that propagated to the tunnel formed the rockburst with the rock ejection from the surface. Other fractures remaining in the rock mass did not enhance the rockburst but could be a potential hazard to the safety of the tunnel in the following use. A parametric study was conducted to discuss the effect of in situ stresses on the rockburst of tunnel, as shown in Fig. 13. In the last two cases, the lateral pressure coefficient was increased to 10, which had exceeded the typical value range in practice, but it helped to reveal the mechanism of rockburst under coupled in situ stress and dynamic disturbance and made our results comparable to those of Zhu et al. (2010). Simulation results showed that in these cases, the horizontal surfaces and surrounding regions of the tunnel were initially under the state of tension instead. The superposition between the initial tensile stress and the reflected tensile stress weakened the ability

Fig. 13. The rockburst of tunnels subjected to a same dynamic disturbance and different in situ horizontal stresses (PH) and vertical stresses (PV) (enlarged view around the tunnel): (a) PH ¼ 20 MPa, PV ¼ 10 MPa; (b) PH ¼ 10 MPa, PV ¼ 1 MPa; and (c) Case 5 in Li et al. (2018d), PH ¼ 20 MPa, PV ¼ 2 MPa.

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

to withstand the dynamic disturbance, therefore severer rockbursts occurred. It is concluded that the most dangerous situation for deep tunnels is that there is a large difference between the horizontal and vertical in situ stresses and the dynamic disturbance comes from the direction of higher in situ stress. Of course, it would also be possible to have a very large difference between the horizontal and vertical in situ stresses caused by other activities instead of sole geological conditions. Li et al. (2018e) used PNMM to simulate the rock scratch tests. Three different failure patterns, i.e. brittle failure, ductile failure, and brittle-ductile combined failure, under decreasing cutting depths were successfully modeled (Fig. 14). It is found that there is a

695

range of cutting depth at which the mode of rock scratching is in the transition from ductile to brittle. The transitional cutting depth range could be estimated from the result of mechanical specific energy (MSE). Besides, a detailed parametric study was performed by a series of PNMM simulations to investigate the effect of operational parameters of the cutter, including the cutting depth, cutting speed, and rake angle. The effect of cutting depth was divided into three phases: the value of MSE decreases rapidly when the cutting depth is shallow, remains approximately constant when the cutting depth is intermediate, and decreases nonlinearly at an intermediate rate in deep cuts. The effect of cutter rake angle was significant on cutting force but moderate on MSE, and an actual

Fig. 14. Typical failure patterns in rock scratch test: (a) Ductile failure at shallow cut; (b) Brittle-ductile failure at intermediate cut; and (c) Brittle failure at deep cut (after Li et al., 2018e).

696

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

cutting depth was found to be noticeably greater than that when the cutter has a large rake angle, leading to a significant difference between the cutting force and MSE. The shortcoming of this numerical study was that the cutting speed was set to be 4 m/s, apparently higher than that in an experimental test, which is usually in the order of mm/s. The reason is to shorten the scratching process and reduce the computation time to a practical level. The same problem had been encountered by Jaime et al. (2015) in their FEM analysis and we adopted the same cutting speed as them. The cutting speed was found to greatly affect the cutting efficiency since rock had a much higher strength at such high strain rates. He and Xu (2016) observed experimentally that the effects of cutting speed can be regarded to be negligible for realistic cutting speeds. We suggest that it is necessary to apply a more efficient numerical model to this simulation.

PNMM. Compared with the block-DEM, e.g. universal distinct element code (UDEC) and 3D distinct element code (3DEC), PNMM has an advantage that it is easier to have the fracture propagation driven into blocks due to the existence of internal particles. This advantage also applies when comparing PNMM with FDEM, especially the widely used code Y-Geo (Lisjak and Grasselli, 2014). However, compared with previously mentioned numerical models, PNMM still lacks improvement in the following aspects. First, PNMM is currently still a 2D model, which limits its application to many problems. Second, application of PNMM is still limited, since it fails to simulate heats, fluids, and coupled problems. Third, the efficiency of PNMM should be further improved by incorporating parallel computation. Based on the development of PNMM and simulation results obtained over the past years, the following recommendations for future work are outlined:

3.4. Discussion Since PNMM is developed by incorporating the particle concept into NMM, it is reasonable to be compared it with other particlebased models, e.g. the DEM (Potyondy and Cundall, 2004), SPH (Lucy, 1977), PD (Rabczuk and Ren, 2017), LM (Zhao et al., 2011) and reproducing kernel particle method (RKPM) (Liu et al., 1995). However, particles in PNMM are in fact based on a quite different concept. Serving as an extension and improvement to NMM, PNMM has the same theoretical foundation and global governing equation as NMM. In PNMM, the DOFs are defined on manifold elements rather than particles. The element matrices are to be assembled and the mechanical fields of elements are to be solved from the global equation. Recalling that the link between adjacent particles only presents the continuous status and carries no additional force, particles in a continuum do not interact with each other, which can distinguish PNMM from other particle methods. Particles in PNMM can be regarded as a point collocation under most circumstances. The force between particles is only defined at the contact between different continuums or blocks on their surfaces. Due to this characteristic, introduction of particles as well as links does not bring in many micro-parameters as other methods. Macro-mechanical parameters are assigned to particles as those assigned to elements in NMM. Parameters at the contact are similar to NMM as well, since the only difference is that contact is defined between particles instead of blocks and the same open-close iteration is utilized. For fracturing analysis, JHB parameters are introduced on the links in order to take the rate-dependent behavior of rock into consideration. The determination of JHB parameters could be found in the literature (Li et al., 2017). In summary, particles in PNMM are functioned in the following aspects. First, particles carrying parameters assist in simulating the heterogeneity of rock materials. Second, the numerical integration scheme is implemented by particles. Third, fractures could be initiated and propagate on the particle level. Finally, the contacts between different blocks are conducted by particles. Compared with NMM, PNMM simplifies the contact operation between blocks, improves the flexibility in initiating and determining the propagation path of fractures, and gains the ability to simulate the heterogeneity of rock materials. Compared with the particle-DEM, e.g. particle flow code (PFC) (Potyondy and Cundall, 2004), PNMM is a continumm-discontinumm coupled method, which is supposed to provide a better solution to continuous problems. Moreover, PNMM does not have the micro-parameters on particles that need to be calibrated prior to the simulation as PFC does. However, PNMM would probably have a lower accuracy than PFC if their models have the same number of particles, since the DOFs in PNMM are defined on elements instead of particles. Besides, PFC is apparently more suitable for granular materials than

(1) Investigating the heterogeneity of rock materials. The heterogeneity could be defined by either a random distribution of material properties among particles or a meso-scale structure. In the second scenario, the geometry information of aggregates needs to be provided first. Then, particles located in aggregates are assigned with different material parameters. This work could be an important application of PNMM in the future. In the next phase, a PNMM model should be able to be generated from a computed tomography (CT) scanning image. The effect of the varying particle size and particle allocation should also be studied in the scope of this work. (2) Enhancing the computational efficiency of PNMM by developing a parallelized version. Both central processing unit (CPU)-based and GPU-based technologies are possibly adopted. The first attempt could be made on solving the global equation. (3) Developing an explicitly coupled method using PNMM and NMM. A major difference between PNMM and NMM is the introduction of particles. It is possible to introduce particles as the second level of discretization in a part of the model, and conducts only the first discretization in the rest portion of the model with less interest, e.g. in far field, in the domain with small deformation, in the domain where no fracture and contact occurs. (4) Extending the PNMM code to 3D cases. Some formulae proposed for 3D PMM in Sun (2012) could be introduced. The extension of the proposed particle integration scheme to 3D cases is straightforward. The 3D PNMM could be applied to the problems that do not obey the plane strain and plane stress assumptions, e.g. the projectile penetration in rock plates. (5) Incorporating thermal and fluid components into PNMM to investigate the coupled problems, which have been increasingly important in rock engineering. (6) Application of PNMM could be extended to the stress wave propagation in jointed rock masses, which is also an important topic in rock dynamics. 4. Concluding remarks Over the last decade, PNMM has emerged as a promising numerical tool for rock dynamics. PNMM inherits the cover system from NMM and introduces a dual-level discretization using particles. Compared with NMM, PNMM simplifies the contact operation between blocks, improves the flexibility in initiating and determining the propagation path of fractures, and gains the ability to simulate the heterogeneity of rock materials. Compared with PMM,

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

several aspects have been re-clarified and improved in PNMM, including adoption of high-order interpolation functions, development of a particle integration scheme, an enrichment function around fracture tips, an improved algorithm for the generation of links, unambiguous failure of links, incorporation of a ratedependent strength model, expanded applications in rock dynamics, and improved calculation performance. Due to its performances, PNMM could be easily utilized to study the heterogeneity of rock materials, the initiation, coalescence and propagation of fractures, the detachment and post-failure behavior of rock fragments, the contact between rock blocks, and the rate-dependent behavior of rock, making it a promising tool for modeling of rock dynamic fracturing problems. Application of PNMM mainly is now focused on rock fracturing under different dynamic loads, including rock spalling, plate impact, rockburst in tunnel, and rock scratching. Furthermore, the advantages and limitations of PNMM, compared with other widely used models, are summarized. Advices to the future development of PNMM are also outlined. Conflicts of interest The authors wish to confirm that there are no known conflicts of interest associated with this publication and there has been no significant financial support for this work that could have influenced its outcome. Acknowledgements The authors would like to acknowledge the financial support to the development of PNMM from the Laboratory of Rock Mechanics at École Polytechnique Fédérale de Lausanne (LMR-EPFL), Monash University, and the National Natural Science Foundation of China (Grant No. 11802058) in the past few years. References Amiri F, Anitescu C, Arroyo M, Bordas SPA, Rabczuk T. XLME interpolants, a seamless bridge between XFEM and enriched meshless methods. Computational Mechanics 2014;53(1):45e57. An HM, Liu HY, Han H, Zheng X, Wang XG. Hybrid finite-discrete element modelling of dynamic fracture and resultant fragment casting and muck- piling by rock blast. Computers and Geotechnics 2017;81:322e45. An X. Extended numerical manifold method for engineering failure analysis. PhD Thesis. Singapore: Nanyang Technological University; 2010. An X, Ning Y, Ma G, He L. Modeling progressive failures in rock slopes with nonpersistent joints using the numerical manifold method. International Journal for Numerical and Analytical Methods in Geomechanics 2014;38(7):679e701. An XM, Li LX, Ma GW, Zhang HH. Prediction of rank deficiency in partition of unitybased methods with plane triangular or quadrilateral meshes. Computer Methods in Applied Mechanics and Engineering 2011;200(5e8):665e74. Atluri SN, Zhu T. A new meshless local Petrov-Galerkin (MLPG) approach in computational mechanics. Computational Mechanics 1998;22(2):117e27. Babuska I, Melenk JM. The partition of unity method. International Journal for Numerical Methods in Engineering 1997;40(4):727e58. Bao H, Zhao Z. Modeling brittle fracture with the nodal-based discontinuous deformation analysis. International Journal of Computational Methods 2013;10(6):1350040. https://doi.org/10.1142/S0219876213500400. Belytschko T, Gracie R, Ventura G. A review of extended/generalized finite element methods for material modeling. Modelling and Simulation in Materials Science and Engineering 2009;17(4):043001. https://doi.org/10.1088/0965-0393/17/4/ 043001. Belytschko T, Krongauz Y, Organ D, Fleming M, Krysl P. Meshless methods: an overview and recent developments. Computer Methods in Applied Mechanics and Engineering 1996;139(1e4):3e47. Belytschko T, Lu YY, Gu L. Element-free Galerkin methods. International Journal for Numerical Methods in Engineering 1994;37(2):229e56. Belytschko T, Organ D, Gerlach C. Element-free Galerkin methods for dynamic fracture in concrete. Computer Methods in Applied Mechanics and Engineering 2000;187(3e4):385e99. Cai Y, Zhuang X, Zhu H. A generalized and efficient method for finite cover generation in the numerical manifold method. International Journal of Computational Methods 2013;10(5):1350028. https://doi.org/10.1142/ S021987621350028X.

697

Chen G, Ohnishi Y, Ito T. Development of high-order manifold method. International Journal for Numerical Methods in Engineering 1998;43(4):685e712. Chen H, Zhao Z, Sun J. Coupled hydro-mechanical model for fractured rock masses using the discontinuous deformation analysis. Tunnelling and Underground Space Technology 2013;38:506e16. Cheng YM, Zhang YH. Formulation of a three-dimensional numerical manifold method with tetrahedron and hexahedron elements. Rock Mechanics and Rock Engineering 2008;41(4):601e28. Cheng YM, Zhang YH, Chen WS. Wilson non-conforming element in numerical manifold method. International Journal for Numerical Methods in Biomedical Engineering 2002;18(12):877e84. Chiou YJ, Lee YM, Tsay RJ. Mixed mode fracture propagation by manifold method. International Journal of Fracture 2002;114(4):327e47. Cho N, Martin CD, Sego DC. A clumped particle model for rock. International Journal of Rock Mechanics and Mining Sciences 2007;44(7):997e1010. Clough RW. The finite element method in plane stress analysis. In: Proceedings of the 2nd ASCE conference on electronic computation. American Society of Civil Engineers (ASCE); 1960. p. 345e78. Cundall PA. A computer model for simulating progressive, large scale movements in blocky rock system. In: Proceedings of the international symposium on rock fracture. International Society of Rock Mechanics (ISRM); 1971. p. 129e36. Cundall PA, Strack ODL. A discrete numerical model for granular assemblies. Géotechnique 1979;29(1):47e65. Dassault Systèmes. ABAQUS 2017. Dassault Systèmes; 2017. Fakhimi A, Lanari M. DEM-SPH simulation of rock blasting. Computers and Geotechnics 2014;55:158e64. Fang ZQ, Hu GM, Du J, Fan Z, Liu J. A contact detection algorithm for multi- sphere particles by means of two-level-grid-searching in DEM simulations. International Journal for Numerical Methods in Engineering 2015;102(13):1869e93. Feng YT, Han K, Owen DRJ. Energy-conserving contact interaction models for arbitrarily shaped discrete elements. Computer Methods in Applied Mechanics and Engineering 2012;205e208:169e77. Ferellec JF, McDowell GR. A method to model realistic particle shape and inertia in DEM. Granular Matter 2010;12(5):459e67. Gao H, Cheng Y. A complex variable meshless manifold method for fracture problems. International Journal of Computational Methods 2010;7(1):55e81. Ghasemzadeh H, Ramezanpour MA, Bodaghpour S. Dynamic high order numerical manifold method based on weighted residual method. International Journal for Numerical Methods in Engineering 2014;100(8):596e619. Gingold RA, Monaghan JJ. Smoothed particle hydrodynamics: theory and application to non-spherical stars. Monthly Notices of the Royal Astronomical Society 1977;181:375e89. Gu J, Zhao Z. Considerations of the discontinuous deformation analysis on wave propagation problems. International Journal for Numerical and Analytical Methods in Geomechanics 2009;33(12):1449e65. Gui YL, Bui HH, Kodikara J, Zhang QB, Zhao J, Rabczuk T. Modelling the dynamic failure of brittle rocks using a hybrid continuum-discrete element method with a mixed-mode cohesive fracture model. International Journal of Impact Engineering 2016;87:146e55. Hatzor YH, Arzi AA, Zaslavsky Y, Shapira A. Dynamic stability analysis of jointed rock slopes using the DDA method: king Herod’s Palace, Masada, Israel. International Journal of Rock Mechanics and Mining Sciences 2004;41(5):813e32. He L, An X, Liu X, Zhao Z, Yang S. Augmented numerical manifold method with implementation of flat-top partition of unity. Engineering Analysis with Boundary Elements 2015;61:153e71. He L, An XM, Ma GW, Zhao ZY. Development of three-dimensional numerical manifold method for jointed rock slope stability analysis. International Journal of Rock Mechanics and Mining Sciences 2013;64:22e35. He L, Ma G. Development of 3D numerical manifold method. International Journal of Computational Methods 2010;7(1):107e29. He L, Tian Q, Zhao Z, Zhao X, Zhang Q, Zhao J. Rock slope stability and stabilization analysis using the coupled DDA and FEM method: NDDA approach. International Journal of Geomechanics 2018;18(6):04018044. https://doi.org/10.1061/ (ASCE)GM.1943-5622.0001098. He L, Zhang QB. Numerical investigation of arching mechanism to underground excavation in jointed rock mass. Tunnelling and Underground Space Technology 2015;50:54e67. He X, Xu C. Specific energy as an index to identify the critical failure mode transition depth in rock cutting. Rock Mechanics and Rock Engineering 2016;49(4): 1461e78. Hole cek M, Moravec F. Hyperelastic model of a material which microstructure is formed by “balls and springs”. International Journal of Solids and Structures 2006;43(24):7393e406. Hrennikoff A. Solution of problems of elasticity by the framework method. Journal of Applied Mechanics 1941;8(4):169e75. Jaime MC, Zhou Y, Lin JS, Gamwo IK. Finite element modeling of rock cutting and its fragmentation process. International Journal of Rock Mechanics and Mining Sciences 2015;80:137e46. Jiang H, Wang L, Li L, Guo Z. Safety evaluation of an ancient masonry seawall structure with modified DDA method. Computers and Geotechnics 2014;55: 277e89. Jiao YY, Zhang HQ, Tang HM, Zhang XL, Adoko AC, Tian HN. Simulating the process of reservoir-impoundment-induced landslide using the extended DDA method. Engineering Geology 2014;182:37e48.

698

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700

Jiao YY, Zhang XL, Zhao J, Liu QS. Viscous boundary of DDA for modeling stress wave propagation in jointed rock. International Journal of Rock Mechanics and Mining Sciences 2007;44(7):1070e6. Jing L. A review of techniques, advances and outstanding issues in numerical modelling for rock mechanics and rock engineering. International Journal of Rock Mechanics and Mining Sciences 2003;40(3):283e353. Johnson GR, Holmquist TJ, Beissel SR. Response of aluminum nitride (including a phase change) to large strains, high strain rates, and high pressures. Journal of Applied Physics 2003;94(3):1639e46. Kamai R, Hatzor YH. Numerical analysis of block stone displacements in ancient masonry structures: a new method to estimate historic ground motions. International Journal for Numerical and Analytical Methods in Geomechanics 2008;32(11):1321e40. Kong X, Liu J. Dynamic failure numeric simulations of model concrete-faced rock-fill dam. Soil Dynamics and Earthquake Engineering 2002;22(9e12):1131e4. Lajtai EZ, Duncan EJS, Carter BJ. The effect of strain rate on rock strength. Rock Mechanics and Rock Engineering 1991;24(2):99e109. Li S, Liu WK. Meshfree and particle methods and their applications. Applied Mechanics Reviews 2002;55(1):1e34. Li S, Cheng Y, Wu YF. Numerical manifold method based on the method of weighted residuals. Computational Mechanics 2005a;35(6):470e80. Li SC, Li SC, Cheng YM. Enriched meshless manifold method for two-dimensional crack modeling. Theoretical and Applied Fracture Mechanics 2005b;44(3): 234e48. Li X. Development of particle-based numerical manifold method (PNMM) for dynamic rock fracturing. PhD Thesis. Monash University; 2017. Li X, Zhang QB, He L, Zhao J. Particle-based numerical manifold method to model dynamic fracture process in rock blasting. International Journal of Geomechanics 2017;17(5):E4016014. https://doi.org/10.1061/(ASCE)GM.1943-5622. 0000748. Li XF, Li X, Li HB, Zhang QB, Zhao J. Dynamic tensile behaviours of heterogeneous rocks: the grain scale fracturing characteristics on strength and fragmentation. International Journal of Impact Engineering 2018a;118:98e118. Li XF, Zhang QB, Li HB, Zhao J. Grain-based discrete element method (GB-DEM) modelling of multi-scale fracturing in rocks under dynamic loading. Rock Mechanics and Rock Engineering 2018b;51(12):3785e817. Li X, Zhang QB, Zhao J, Li JC. A numerical study of spalling under different dynamic loads using PNMM. In: Li CC, Li X, Zhang ZX, editors. Rock dynamics e experiments, theories and applications 3: proceedings of the 3rd international confrence on rock dynamics and applications (RocDyn-3). CRC Press; 2018c. p. 299e303. Li X, Li XF, Zhang QB, Zhao J. A numerical study of spalling and related rockburst under dynamic disturbance using a particle-based numerical manifold method (PNMM). Tunnelling and Underground Space Technology 2018d;81:438e49. Li X, Zhang Q, Li J, Zhao J. A numerical study of rock scratch tests using the particlebased numerical manifold method. Tunnelling and Underground Space Technology 2018e;78:106e14. Liang ZZ, Tang CA, Li HX, Xu T, Zhang YB. Numerical simulation of 3D failure process in heterogeneous rocks. International Journal of Rock Mechanics and Mining Sciences 2004;41(Supp. 1):323e8. Lin CT, Amadei B, Jung J, Dwyer J. Extensions of discontinuous deformation analysis for jointed rock masses. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 1996;33(7):671e94. Lin S, Qi Y, Su H. Formulation of high-order numerical manifold method and fast simplex integration based on special matrix operations. In: MacLaughlin MM, Sitar N, editors. Proceedings of the 7th international conference on the analysis of discontinuous deformation (ICADD-7); 2005. p. 183e90. Honolulu, USA. Lisjak A. Investigating the influence of mechanical anisotropy on the fracturing behaviour of brittle clay shales with application to deep geological repositories. University of Toronto; 2013. PhD Thesis. Lisjak A, Grasselli G. A review of discrete modeling techniques for fracturing processes in discontinuous rock masses. Journal of Rock Mechanics and Geotechnical Engineering 2014;6(4):301e14. Liszka TJ, Duarte CAM, Tworzydlo WW. hp-meshless cloud method. Computer Methods in Applied Mechanics and Engineering 1996;139(1e4):263e88. Liu GY, Zhuang XY, Cui ZQ. Three-dimensional slope stability analysis using independent cover based numerical manifold and vector method. Engineering Geology 2017;225:83e95. Liu MB, Liu GR. Smoothed particle hydrodynamics (SPH): an overview and recent developments. Archives of Computational Methods in Engineering 2010;17(1):25e76. Liu Q, Jiang Y, Wu Z, He J. A Voronoi element based-numerical manifold method (VE-NMM) for investigating micro/macro-mechanical properties of intact rocks. Engineering Fracture Mechanics 2018a;199:71e85. Liu Q, Jiang Y, Wu Z, Qian Z, Xu X. Numerical modeling of acoustic emission during rock failure process using a Voronoi element based e explicit numerical manifold method. Tunneling and Underground Space Technology 2018b;79: 175e89. Liu WK, Jun S, Li S, Adee J, Belytschko T. Reproducing kernel particle methods for structural dynamics. International Journal for Numerical Methods in Engineering 1995;38(10):1655e79. Liu WK, Li S, Belytschko T. Moving least-square reproducing kernel methods (I) methodology and convergence. Computer Methods in Applied Mechanics and Engineering 1997;143(1e2):113e54.

Liu XW, Liu QS, He J, Liu B. Modified contact model with rock joint constitutive in numerical manifold method. Engineering Analysis with Boundary Elements 2018c;93:63e71. Lu M, McDowell GR. The importance of modelling ballast particle shape in the discrete element method. Granular Matter 2007;9:69e80. Lucy LB. A numerical approach to the testing of the fission hypothesis. The Astronomical Journal 1977;82:1013e24. Ma G, An X, He L. The numerical manifold method: a review. International Journal of Computational Methods 2010;7(1):1e32. Ma GW, An XM. Numerical simulation of blasting-induced rock fractures. International Journal of Rock Mechanics and Mining Sciences 2008;45(6):966e75. Ma GW, An XM, Zhang HH, Li LX. Modeling complex crack problems using the numerical manifold method. International Journal of Fracture 2009;156(1): 21e35. Ma ZH, Wang H, Pu SH. A parallel meshless dynamic cloud method on graphic processing units for unsteady compressible flows past moving boundaries. Computer Methods in Applied Mechanics and Engineering 2015;285:146e65. Mahabadi OK. Investigating the influence of micro-scale heterogeneity and microstructure on the failure and mechanical behaviour of geomaterials. University of Toronto; 2012. PhD Thesis. Masuda K, Mizutani H, Yamada I. Experimental-study of strain-rate dependence and pressure dependence of failure properties of granite. Journal of Physics of the Earth 1987;35(1):37e66. Moës N, Dolbow J, Belytschko T. A finite element method for crack growth without remeshing. International Journal for Numerical Methods in Engineering 1999;46(1):131e50. Monaghan JJ. Smoothed particle hydrodynamics and its diverse applications. Annual Review of Fluid Mechanics 2012;44:323e46. Munjiza A, Andrews KRF. NBS contact detection algorithm for bodies of similar size. International Journal for Numerical Methods in Engineering 1998;43(1):131e49. Narasimhan TN, Witherspoon PA. An integrated finite difference method for analyzing fluid flow in porous media. Water Resources Research 1976;12(1): 57e64. Nguyen VP, Rabczuk T, Bordas S, Duflot M. Meshless methods: a review and computer implementation aspects. Mathematics and Computers in Simulation 2008;79(3):763e813. Ning Y, Zhao Z. A detailed investigation of block dynamic sliding by the discontinuous deformation analysis. International Journal for Numerical and Analytical Methods in Geomechanics 2013;37(15):2373e93. Ning YJ, An XM, Ma GW. Footwall slope stability analysis with the numerical manifold method. International Journal of Rock Mechanics and Mining Sciences 2011;48(6):964e75. Potyondy DO, Cundall PA. A bonded-particle model for rock. International Journal of Rock Mechanics and Mining Sciences 2004;41(8):1329e64. Rabczuk T, Areias PMA. A new approach for modelling slip lines in geological materials with cohesive models. International Journal for Numerical and Analytical Methods in Geomechanics 2006;30(11):1159e72. Rabczuk T, Belytschko T. Cracking particles: a simplified meshfree method for arbitrary evolving cracks. International Journal for Numerical Methods in Engineering 2004;61(13):2316e43. Rabczuk T, Belytschko T. A three-dimensional large deformation meshfree method for arbitrary evolving cracks. Computer Methods in Applied Mechanics and Engineering 2007;196(29e30):2777e99. Rabczuk T, Zi G. A meshfree method based on the local partition of unity for cohesive cracks. Computational Mechanics 2007;39(6):743e60. Rabczuk T, Zi G, Bordas S, Nguyen-Xuan H. A simple and robust three- dimensional cracking-particle method without enrichment. Computer Methods in Applied Mechanics and Engineering 2010a;199(37e40):2437e55. Rabczuk T, Bordas S, Zi G. On three-dimensional modelling of crack growth using partition of unity methods. Computers and Structures 2010b;88(23e24): 1391e411. Rabczuk T, Ren H. A peridynamics formulation for quasi-static fracture and contact in rock. Engineering Geology 2017;225:42e8. Ren H, Zhuang X, Cai Y, Rabczuk T. Dual-horizon peridynamics. International Journal for Numerical Methods in Engineering 2016;108(12):1451e76. Ren H, Zhuang X, Rabczuk T. Dual-horizon peridynamics: a stable solution to varying horizons. Computer Methods in Applied Mechanics and Engineering 2017;318:762e82. Rougier E, Knight EE, Broome ST, Sussman AJ, Munjiza A. Validation of a threedimensional Finite-Discrete Element Method using experimental results of the split Hopkinson pressure bar test. International Journal of Rock Mechanics and Mining Sciences 2014;70:101e8. Shi GH, Goodman RE. Two dimensional discontinuous deformation analysis. International Journal for Numerical and Analytical Methods in Geomechanics 1985;9(6):541e56. Shi GH, Goodman RE. Generalization of two-dimensional discontinuous deformation analysis for forward modelling. International Journal for Numerical and Analytical Methods in Geomechanics 1989;13(4):359e80. Shi GH. Manifold method of material analysis. In: Transactions of the 9th army conference on applied mathematics and computing. Minneapolis, USA: University of Minnesota; 1991. p. 51e76. Shi GH. Modeling rock joints and blocks by manifold method. In: Rock mechanics: proceedings of the 33th US symposium. American Rock Mechanics Association (ARMA); 1992. p. 639e48.

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700 Shi GH. Numerical manifold method. In: Li JC, Wang CY, Sheng J, editors. Proceedings of the 1st international conference on analysis of discontinuous deformation (ICADD-1); 1995. p. 187e222. Chungli, Taiwan, China. Shi GH. Simplex integration for manifold method, FEM, DDA and analytical analysis. In: Discontinuous deformation analysis (DDA) and simulations of discontinuous media: proceedings of the 1st international forum on discontinuous deformation analysis (DDA) and simulations of discontinuous media. Albuquerque, USA: TSI Press; 1996. p. 205e62. Shi GH. Numerical manifold method. In: Ohnishi Y, editor. Proceedings of the 2nd international conference on analysis of discontinuous deformation (ICADD-2); 1997. p. 1e35. Kyoto, Japan. Shi GH. Rock block stability analysis of slopes and underground power houses. In: Zhao J, Ohnishi Y, Zhao GF, Sasaki T, editors. Advances in discontinuous numerical methods and applications in geomechanics and geoengineering. CRC Press; 2012. p. 3e16. Shyu K, Salami MR. Manifold with four-node isoparametric finite element method. In: Proceedings of the 1st working forum on the manifold method of material analysis. US Army Corps of Engineers; 1995. p. 165e82. Silling SA. Reformulation of elasticity theory for discontinuities and long-range forces. Journal of the Mechanics and Physics of Solids 2000;48(1):175e209. Smeets B, Odenthal T, Vanmaercke S, Ramon H. Polygon-based contact description for modeling arbitrary polyhedra in the discrete element method. Computer Methods in Applied Mechanics and Engineering 2015;290:277e89. Su H, Xie X, Liang Q. Automatic programming for high-order numerical manifold method. In: Lu M, editor. Development and application of discontinuous modelling for rock engineering, proceedings of the 6th international conference on analysis of discontinuous deformation (ICADD-6). CRC Press; 2003. p. 153e7. Sukumar N, Moës N, Moran B, Belytschko T. Extended finite element method for three-dimensional crack modelling. International Journal for Numerical Methods in Engineering 2000;48(11):1549e70. Sun L. Particle manifold method (PMM) for multiscale continuous- discontinuous analysis. PhD Thesis. Ecole Polytechnique Federale de Lausanne; 2012. Sun L, Zhao GF, Zhao J. Particle manifold method (PMM): a new continuumdiscontinuum numerical model for geomechanics. International Journal for Numerical and Analytical Methods in Geomechanics 2013;37(11):1711e36. Sun L, Zhao GF, Zhao J. A multiscale manifold method using particle representations of the physical domain. Geomechanics and Geoengineering 2014;9(2): 124e32. Tang C. Numerical simulation of progressive rock failure and associated seismicity. International Journal of Rock Mechanics and Mining Sciences 1997;34(2):249e61. Terada K, Asai M, Yamagishi M. Finite cover method for linear and non-linear analyses of heterogeneous solids. International Journal for Numerical Methods in Engineering 2003;58(9):1321e46. Tian Q, Zhao Z, Bao H. Block fracturing analysis using nodal-based discontinuous deformation analysis with the double minimization procedure. International Journal for Numerical and Analytical Methods in Geomechanics 2014;38(9): 881e902. Tsesarsky M, Talesnick ML. Mechanical response of a jointed rock beam e numerical study of centrifuge models. International Journal for Numerical and Analytical Methods in Geomechanics 2007;31(8):977e1006. Vocialta M, Richart N, Molinari JF. 3D dynamic fragmentation with parallel dynamic insertion of cohesive elements. International Journal for Numerical Methods in Engineering 2017;109(12):1655e78. Vu-Quoc L, Zhang X, Walton OR. A 3D discrete-element method for dry granular flows of ellipsoidal particles. Computer Methods in Applied Mechanics and Engineering 2000;187(3e4):483e528. Wang SY, Sloan SW, Tang CA. Three-dimensional numerical investigations of the failure mechanism of a rock disc with a central or eccentric hole. Rock Mechanics and Rock Engineering 2014;47(6):2117e37. Wei W, Zhao Q, Jiang Q, Grasselli G. Three new boundary conditions for the seismic response analysis of geomechanics problems using the numerical manifold method. International Journal of Rock Mechanics and Mining Sciences 2018;105:110e22. Wellmann C, Wriggers P. A two-scale model of granular materials. Computer Methods in Applied Mechanics and Engineering 2012;205e208:46e58. Wong LNY, Wu Z. Application of the numerical manifold method to model progressive failure in rock slopes. Engineering Fracture Mechanics 2014;119:1e20. Wu JH, Lin JS, Chen CS. Dynamic discrete analysis of an earthquake-induced largescale landslide. International Journal of Rock Mechanics and Mining Sciences 2009;46(2):397e407. Wu Z, Wong LNY. Elastic-plastic cracking analysis for brittle-ductile rocks using manifold method. International Journal of Fracture 2013a;180(1):71e91. Wu Z, Wong LNY. Modeling cracking behavior of rock mass containing inclusions using the enriched numerical manifold method. Engineering Geology 2013b;162:1e13. Wu Z, Wong LNY, Fan L. Dynamic study on fracture problems in viscoelastic sedimentary rocks using the numerical manifold method. Rock Mechanics and Rock Engineering 2013;46(6):1415e27. Wu Z, Fan L. The numerical manifold method for elastic wave propagation in rock with time-dependent absorbing boundary conditions. Engineering Analysis with Boundary Elements 2014;46:41e50. Wu Z, Wong LNY. Extension of numerical manifold method for coupled fluid flow and fracturing problems. International Journal for Numerical and Analytical Methods in Geomechanics 2014;38(18):1990e2008.

699

Xu D, Yang Y, Zheng H, Wu A. A high order local approximation free from linear dependency with quadrilateral mesh as mathematical cover and applications to linear elastic fractures. Computers and Structures 2017;178:1e16. Yan B, Regueiro RA, Sture S. Three-dimensional ellipsoidal discrete element modeling of granular materials and its coupling with finite element facets. Engineering Computations 2010;27(4):519e50. Yang Q, Cai F, Su Z, Ugai K, Xu L, Huang R, Xu Q. Numerical simulation of granular flows in a large flume using discontinuous deformation analysis. Rock Mechanics and Rock Engineering 2014a;47(6):2299e306. Yang S, Ma G, Ren X, Ren F. Cover refinement of numerical manifold method for crack propagation simulation. Engineering Analysis with Boundary Elements 2014b;43:37e49. Yang Y, Chen L, Xu D, Zheng H. Free and forced vibration analyses using the fournode quadrilateral element with continuous nodal stress. Engineering Analysis with Boundary Elements 2016a;70:1e11. Yang Y, Sun G, Zheng H, Fu X. A four-node quadrilateral element fitted to numerical manifold method with continuous nodal stress for crack analysis. Computers and Structures 2016b;177:69e82. Yang Y, Xu D, Zheng H. Application of the three-node triangular element with continuous nodal stress for free vibration analysis. Computers and Structures 2016c;169:69e80. Yang Y, Zheng H. A three-node triangular element fitted to numerical manifold method with continuous nodal stress for crack analysis. Engineering Fracture Mechanics 2016;162:51e75. Yang Y, Zheng H, Sivaselvan MV. A rigorous and unified mass lumping scheme for higher-order elements. Computer Methods in Applied Mechanics and Engineering 2017;319:491e514. Yang YT, Guo HW, Fu XD, Zheng H. Boundary settings for the seismic dynamic response analysis of rock masses using the numerical manifold method. International Journal for Numerical and Analytical Methods in Geomechanics 2018;42(6):1095e122. Zhang HH, Chen YL, Li LX, An XM, Ma GW. Accuracy comparison of rectangular and triangular mathematical elements in the numerical manifold method. In: Ma GW, Zhou YX, editors. Analysis of discontinuous deformation e new developments and applications. Research Publishing; 2010. p. 297e303. Zhang HH, Ma GW, Ren F. Implementation of the numerical manifold method for thermo-mechanical fracture of planar solids. Engineering Analysis with Boundary Elements 2014;44:45e54. Zhang HH, Ma GW, Fan LF. Thermal shock analysis of 2D cracked solids using the numerical manifold method and precise time integration. Engineering Analysis with Boundary Elements 2017;75:46e56. Zhang HH, Han SY, Fan LF, Huang D. The numerical manifold method for 2D transient heat conduction problems in functionally graded materials. Engineering Analysis with Boundary Elements 2018;88:145e55. Zhang Y, Wang J, Xu Q, Chen G, Zhao JX, Zheng L, Han Z, Yu P. DDA validation of the mobility of earthquake-induced landslides. Engineering Geology 2015;194:38e51. Zhao GF, Fang J, Zhao J. A 3D distinct lattice spring model for elasticity and dynamic failure. International Journal for Numerical and Analytical Methods in Geomechanics 2011;35(8):859e85. Zhao GF, Khalili N, Fang JN, Zhao J. A coupled distinct lattice spring model for rock failure under dynamic loads. Computers and Geotechnics 2012a;42:1e20. Zhao GF, Zhao XB, Zhu JB. Application of the numerical manifold method for stress wave propagation across rock masses. International Journal for Numerical and Analytical Methods in Geomechanics 2014;38(1):92e110. Zhao J, Sun L, Zhu JB. Modelling P-wave transmission across rock fractures by particle manifold method (PMM). Geomechanics and Geoengineering 2012b;7(3):175e81. Zhao SF, Zhao GF. Implementation of a high order lattice spring model for elasticity. International Journal of Solids and Structures 2012;49(18):2568e81. Zhao XB, Zhao J, Cai JG, Hefny AM. UDEC modelling on wave propagation across fractured rock masses. Computers and Geotechnics 2008;35(1):97e104. Zheng H, Liu F, Li C. The MLS-based numerical manifold method with applications to crack analysis. International Journal of Fracture 2014b;190(1e2):147e66. Zheng H, Xu DD. New strategies for some issues of numerical manifold method in simulation of crack propagation. International Journal for Numerical Methods in Engineering 2014;97(13):986e1010. Zheng H, Liu F, Du X. Complementarity problem arising from static growth of multiple cracks and MLS-based numerical manifold method. Computer Methods in Applied Mechanics and Engineering 2015;295:150e71. Zheng H, Yang Y. On generation of lumped mass matrices in partition of unity based methods. International Journal for Numerical Methods in Engineering 2017;112(8):1040e69. Zheng QJ, Zhou ZY, Yu AB. Contact forces between viscoelastic ellipsoidal particles. Powder Technology 2013;248:25e33. Zheng W, Zhuang X, Tannant DD, Cai Y, Nunoo S. Unified continuum/discontinuum modelling framework for slope stability assessment. Engineering Geology 2014a;179:90e101. Zhou XP, Zhao Y, Qian QH. A novel meshless numerical method for modeling progressive failure processes of slopes. Engineering Geology 2015;192:139e53. Zhu WC. Numerical modelling of the effect of rock heterogeneity on dynamic tensile strength. Rock Mechanics and Rock Engineering 2008;41(5):771e9. Zhu WC, Li ZH, Zhu L, Tang CA. Numerical simulation on rockburst of underground opening triggered by dynamic disturbance. Tunnelling and Underground Space Technology 2010;25(5):587e99.

700

X. Li, J. Zhao / Journal of Rock Mechanics and Geotechnical Engineering 11 (2019) 684e700 Dr. Xing Li is Associate Professor in the School of Civil Engineering at Southeast University, China since 2018. He received his BSc and MSc degrees in Mechanics at Xi’an Jiaotong University, China in 2011 and 2013, respectively. Xing conducted his doctoral research in Rock Mechanics at Ecole Polytechnique Fédérale de Lausanne (EPFL) of Switzerland during 2013e2015, and obtained his PhD

degree in the Department of Civil Engineering at Monash University, Australia in 2017. His research interests cover the dynamic rock fracturing, numerical manifold method, and discrete element method.