Discrete modeling of rock joints with a smooth-joint contact model

C. Lamert, C. Coll

aDepartmentofCivilandNaturalResourcesEngineering,UniversityofCanterbury,PrivateBag4800,Christchurch8140,NewZealand

bGolderAssociates,Christchurch,NewZealand

Discrete modeling of rock joints with a smooth-joint contact model

C. Lamberta,*, C. Collb

aDepartmentofCivilandNaturalResourcesEngineering,UniversityofCanterbury,PrivateBag4800,Christchurch8140,NewZealand

bGolderAssociates,Christchurch,NewZealand

A R T I C L E I N F O

Articlehistory:

Received 14 May 2013

Received in revised form

20 November 2013

Accepted 8 December 2013

Rock joint

Structural defects such as joints or faults are inherent to almost any rock mass. In many situations those defects have a major impact on slope stability as they can control the possible failure mechanisms. Having a good estimate of their strength then becomes crucial. The roughness of a structure is a major contributor to its strength through two different aspects, i.e. the morphology of the surface (or the shape) and the strength of the asperities (related to the strength of the rock). In the current state of practice, roughness is assessed through idealized descriptions (Patton strength criterion) or through empirical parameters (Barton JRC). In both cases, the multi-dimensionality of the roughness is ignored. In this study, we propose to take advantage of the latest developments in numerical techniques. With 3D photogrammetry and/or laser mapping, practitioners have access to the real morphology of an exposed structure. The derived triangulated surface was introduced into the DEM (discrete element method) code PFC3D to create a synthetic rock joint. The interaction between particles on either side of the discontinuity was described by a smooth-joint model (SJM), hence suppressing the artif i cial roughness introduced by the particle discretization. Shear tests were then performed on the synthetic rock joint. A good correspondence between strengths predicted by the model and strengths derived from well-established techniques was obtained for the fi rst time. Amongst the benef i ts of the methodology is the possibility offered by the model to be used in a quantitative way for shear strength estimates, to reproduce the progressive degradation of the asperities upon shearing and to analyze structures of different scales without introducing any empirical relation.

? 2013 Institute of Rock and Soil Mechanics, Chinese Academy of Sciences. Production and hosting by Elsevier B.V. All rights reserved.

1. Introduction

The presence of discontinuities is inherent to almost any rock mass and is a major contributor to strength and deformation of rock structures (natural or engineering). The characteristics of those discontinuities not only control structurally controlled failures but also greatly inf l uence the shear strength of the rock mass. Being able to describe the structure of a rock mass is critical to an understanding of its potential behavior. The development of various mapping techniques leads to a higher level of conf i dence on crucial characteristics such as location, orientation and persistence from which stochastic discrete fracture network (DFN) representations of the rock fabric are developed (Dershowitz, 1995; Rogers et al., 2007). Based on numerical methods, equivalent rock mass can be created and tested in order to characterize its constitutive behavior (Pierce et al., 2007; Pine et al., 2007; Deisman et al., 2010). These approaches are now able to model the engineering responses of rock and rock masses using some basic measured properties of the rock and the rock mass geometry as inputs. Offering a wider spectrum of predictions than the classical empirically-based classif i cation schemes (anisotropy, heterogeneous, etc.), the synthetic rock mass approach and equivalents (Pierce et al., 2007; Pine et al., 2007; Deisman et al., 2010) are turning to be a step forward for rock mechanics practitioners. However, the question of the shear strength of the discontinuities is in many cases poorly addressed in engineering practice despite having a signif i cant impact on the rock mass strength (Lambert, 2008).

The shear behavior of discontinuities is a combination of various complex phenomena and interactions, such as dilation, asperity failure, deformation and interaction. Direct shear tests on natural rock discontinuities quickly enhanced the inf l uence of roughness on the mechanical behavior of discontinuities. Barton (1973) proposed to assess roughness with an empirical parameter, jointroughness parameter (JRC), from which the shear strength of the discontinuity can be established. Initially estimated by visual comparison with standard roughness prof i les, correlations between JRC and various statistical parameters or fractal dimension were established (Tse and Cruden, 1979; Carr and Warriner, 1989). More recently, laser scanner and photogrammetry were used to def i ne the surface topography and estimate its roughness (Grasselli, 2001; Hans and Boulon, 2003; Haneberg et al., 2007). The dependence of shearing on the location and distribution of the three-dimensional (3D) contact area was demonstrated (Gentier et al., 2000) and new constitutive relations were developed based on a general description of roughness (Grasselli and Egger, 2003). Laser scanning and 3D photogrammetry techniques were applied in the fi eld (Fardin et al., 2004) for large-scale surface measurements. Asperity shape and distribution on a discontinuity can now be measured with a great detail and potentially incorporated in any analysis. However with the complexity of the interaction between the two walls, a complete analytical formulation remains a hard task. Since the fi rst idealized “saw-tooth” description proposed by Patton (1966), various constitutive models were developed that accommodate effect of asperities (Barton and Choubey, 1977; Saeb and Amadei, 1992) and their progressive degradation during shearing (Plesha, 1987; Hutson and Dowding, 1990; Lee et al., 2001; Misra, 2002) to name a few. Despite being each time more advanced, these models still rely on empirical relations or simplif i ed descriptions of the surface asperities.

In an attempt to address this problem, many authors used numerical tools to assess the shear strength of discontinuities. Twodimensional DEM (discrete element method) simulations were fi rst presented as they offer a provision for asperity degradation (Cundall, 2000; Lambert et al., 2004). They have been successfully used to investigate gouge formation and evolution upon shearing (Zhao et al., 2012; Zhao, 2013). However these simulations were at this stage limited to qualitative observations. Hybrid FEM/DEM (Karami and Stead, 2008) and FEM (Giacomini et al., 2008) methods proved their ability to reproduce typical behavior of rock joints including dilation and asperity degradation. Using 3D DEM, Kulatilake et al. (2001) showed that realistic macroscopic friction (i.e. at the joint level) could be obtained combining very small particles at the joint interface and extremely low contact friction. However this approach appears to be not very practical for engineering purposes. No formulation is available to calibrate the micro-properties of the joint model material against a given macroscopic behavior and the macroscopic friction targeted was quite high (friction coeff i cient of 0.7). In the fi eld, discontinuities often exhibit a much lower strength. The particle size required may hence increase the computational cost to unpractical levels. Park and Song (2009) performed numerical shear tests on standard roughness prof i les using the DEM code, PFC3D. This work once again highlighted the current limitations of particulate description as the discrete nature of the medium can introduce an artif i cial roughness to the discontinuity. The apparent roughness of the numerical specimen is higher than the introduced roughness (i.e. the initial roughness of the introduced surface or prof i le). The consequence is a slight overestimation of the strength and most importantly unrealistic predictions of dilation. The later point can be of major importance as joint aperture controls fl uid fl ow in the discontinuities (Hans and Boulon, 2003; Buzzi et al., 2008). The recent development of a new contact model named “smooth-joint model”(SJM) (Pierce et al., 2007) in PFC3D where particles are allowed to slide past one another without over-riding one another was a major breakthrough to represent discontinuities as planar surfaces associated to a realistic behavior for structural defects. In this study, we propose to develop in PFC3D a synthetic rock joint where a digital representation of a surface is introduced and described as a series of SJMs. The mechanical behavior of the synthetic rock joint is then analyzed performing numerical direct shear tests.

Table1Target (laboratory) and calibrated (calculation) bulk properties of the granite.

2. DEM simulations of constant normal stress shear tests

2.1.Thediscreteelementmethod

The commercially available PFC3D (Itasca, 2008) software package was used for the 3D DEM simulations presented here. Unlike continuum codes, materials are described in PFC3D as a discontinuous medium as a collection of spherical rigid particles. The particles displace independently of one another following Newton’s second law and interact with each other through contact forces that are generated at each contact point. Rock and more generally cohesive materials are represented as a bonded particle assembly, adding parallel bonds to create a synthetic material. A parallel bond acts like a conceptual cementitious material between particles. It has a fi nite dimension def i ned as a fraction of the particle diameter, a tensile and shear strength and a normal and tangential stiffness. When the contact force exceeds either tensile or shear strength, the parallel bond breaks and a micro-crack forms between the particles. Micro-cracks can eventually coalesce as external loading is applied and form fractures that can split the material into clusters. The location and the failure mode of the cracks are recorded. A detailed description of contact and bond models is provided in the user manual (Itasca, 2008).

The mechanical response of such assemblies, observed at a macroscopic level, is an emergent property of the complex interactions between the particles. Input parameters of the bonded particle model are micro-properties, contact properties and bond strength, and are not measurable with conventional laboratory apparatus. They are calibrated through an iterative process. Once a particle size distribution has been selected, cylindrical particle assemblies are generated and unconf i ned compression tests are simulated varying micro-properties until the mechanical response of the synthetic material conforms to the mechanical properties (i.e. uniaxial compressive strength, UCS; Young’s modulus; Poisson’s ratio) of the physical material (measured in the lab). A detailed description of the calibration procedure can be found in Potyondy and Cundall (2004). Once properly calibrated, such bonded assemblies proved their ability to reproduce typical behavior of rock-like materials (Kulatilake et al., 2001; Potyondy and Cundall, 2004).

Properties of the granite considered for the scope of this study are given in Table 1. The micro-properties were calibrated accordingly. Normal and shear stiffnesses for contact and parallel bonds have impact on elastic properties of the particle assembly whereas bond shear and normal strengths mainly control UCS values. Various studies by Cundall (2000), Kulatilake et al. (2001) and Park and Song (2009) illustrated the necessity to introduce low particle friction to reproduce the shearing behavior of fracture planes in cohesive materials. In this study, bond strengths were calibrated considering zero friction between particles (φp= 0°).

Besides Potyondy and Cundall (2004) showed that particle friction impacts mainly on the post peak behavior of bulk material with little effect on peak strength. The inf l uence of φpwill be discussed with more detail in Section 3.3. The result of the calibration is given in Table 2 and the emergent bulk properties of the synthetic material are listed in Table 1.

Table2

Material micro-properties for the granite sample.

Particle properties Particle mean radius (mm) Particle radius ratio Particle contact modulus (GPa) Particle normal to shear stiffness ratio Particle friction coeff i cient 1.5 1.66 56.1 2.5 0.0 Parallel-bond properties Parallel-bond radius multiplier Parallel-bond modulus (GPa) Normal to shear stiffness ratio Bond normal strength (MPa) Bond shear strength (MPa) Mean value Standard deviation Mean value Standard deviation 1.0 56.1 2.5 191 19.1 191 19.1

Fig.1.Morphology of the granite surface. All dimensions are in mm.

2.2.Descriptionoftheinterface

The interface morphology used in the simulations is based on a natural discontinuity in granite studied by Grasselli (2001). The surface is 140 mm × 140 mm and the maximum amplitude of the asperities is around 9 mm. Fig. 1 shows a general view of the surface. The 3D surface was triangulated using a Kriging gridding method with a horizontal spacing of 1.4 mm between the grid points (inX-andY-directions). Ninety nine prof i les along the sliding direction (X-direction) were extracted for which the coeff i cientZ2(root mean square of the fi rst derivative of the prof i le) was estimated:

whereziandzi+1are the elevation of two consecutive grid points on a prof i le,Nis the total number of grid points on a prof i le, and Δxis the horizontal spacing.

For each prof i le, a value of JRC was derived using the empirical relation proposed by Yang et al. (2001):

The prof i les of the triangulated surface exhibited an average JRC of 10.4, varying from 4.9 to 13.9. Grasselli et al. (2002) estimated a JRC of 12.5 using a quantitative 3D surface description. The value obtained by Grasselli et al. (2002) is in the range of variation of JRC for individual prof i les but appears to be slightly higher than the average value. However, the value ofZ2and hence the derived JRC are sensitive to sampling intervals (Yu and Vayssade, 1991). Applying the same approach varying the horizontal spacing between the grid points of the triangulated surface, the average JRC increases to a value of 11.6 for a horizontal spacing of 0.56 mm (Lambert and Coll, 2009). The average JRC value of 10.4 was considered to be a reasonable estimate of the surface roughness used for this study.

Fig.2.Smooth joint contact model between ball 1 and ball 2. Surface 1 and surface 2 denote either side of the joint lying at a dip angle of θp(after Itasca, 2008).

2.3.Thesyntheticrockjointmodel

The numerical rock joint consists of a 140 mm × 140 mm × 50 mm (respectivelyX-,Y-, andZ-directions) parallelepiped particle assembly. The specimen genesis procedure is described in detail in Potyondy and Cundall (2004).

For the fi rst series of simulations, three particle assemblies were generated, each one containing around 98,000 particles having a radius ranging from 0.5 mm (in the vicinity of the interface) to 2.4 mm. These specimens differ from one another only in their packing. A discontinuity is usually represented in PFC3D by debonding contacts along a surface. However, the particle geometry is still present and the discrete nature of the medium generates an artif i cial roughness that is added to one of the introduced surfaces, thus creating a particle size dependent joint behavior. For example, the DEM model presented by Park and Song (2009) using a standard prof i le with a JRC of 11.49 exhibited an apparent JRC of 17.55. To overcome the problem, an alternate scheme, termed as the “smooth-joint model” or SJM, initially proposed by Pierce et al. (2007), was implemented in PFC3D (Itasca, 2008). A smooth-joint model is a contact model that simulates the behavior of an interface, regardless of the local particle contact orientation along the interface. A typical smooth-joint is shown in Fig. 2. It allows particles to slide past one another without over-riding one another. A smoothjoint is created by assigning this new contact model to all the contacts between particles that lie upon opposite sides of the surface. The SJM def i nes the tangential and normal directions according to the local orientation of the surface (by opposition to the initial normal and tangential directions of the contact, see Fig. 2). The joint normal and tangential force increments (ΔFnjand ΔFtj, respectively) are derived from normal and tangential displacement increments (ΔUnjand ΔUtj), multiplying by the joint stiffnesses (ΔFnj=knjΔUnjand ΔFtj=ktjΔUtj). The joint force is then adjusted to satisfy the force-displacement relationship and mapped back into the global system. This new formulation accommodates the standard shear behavior of a joint (friction, cohesion and dilation) independently of particle induced roughness. A complete description of the formulation can be found in the manual (Itasca, 2008). An initial study by Lambert et al. (2010) on the behavior of a rockconcrete interface suggested that realistic shear behavior, shear strength and dilation, could be obtained associating the SJM with a true morphology.

An algorithm was developed for the importation of the triangulated surface presented in Section 2.2 into a bonded particle assembly. The same surface was used for each wall of the joint. To be assigned a smooth-joint model, a contact must satisfy two conditions: (1) the two contacting balls must lie on opposite side of the plane containing the triangle; and (2) the projection of the contact location onto this plane must lie within the bounds of the triangle. The orientation of the smooth joint, def i ned by a dip angle and a dip direction, corresponds to the orientation of the triangle. The process is repeated for every triangle of the surface. The joint surface is hence modeled as a collection of smooth-joint contacts. The discontinuity was considered to be purely frictional (i.e. no bond is introduced) with a friction angle set to 20°. No dilation was introduced as macroscopic dilation (i.e. at the joint level) is expected to be an emergent property of the surface topology. The SJM parameters are given in Table 3. The output is a 140 mm × 140 mm synthetic rock joint sample (SRJ) whose morphology corresponds to the natural rock joint that is being analyzed.

In this paper, “synthetic rock joint” or SRJ will refer to the discrete element model of a rock joint. Its properties such as strength or stiffness will be macro-properties (i.e. computed at the scale of sample) and are denoted using an uppercase letter (e.g. Φpeak, Ψpeak). The SJM on the other hand refers to a local contact on the joint surface. Lowercase letters will be used for micro-properties (e.g. φj,knj). A full 3D view of the numerical sample can be seen in Fig. 3.

During the direct shear tests, specimens are fi rstly subjected to a compression along axisZ(Fig. 3) and then to a shearing along axisXat constant normal stress. During the compression stage, the normal load is applied to the upper wall of the specimen while displacements of the lower wall are restraint. The required load is applied in fi ve incremental stages. For each stage, the incremental normal force is equally shared between the particles of the top layer of the specimen (i.e. between the particles whose centers are within one average diameter from the top of the specimen) and progressively applied in 100 time steps. The system is then dynamically set to equilibrium before proceeding to the next stage. During shearing, displacements alongY-axis are restrained whereas shear displacements alongX-axis are applied to the lower wall. The sum of contact forces on the periphery of the upper wall is used to compute the average normal stress and shear stress on the interface whereas relative normal and tangential displacements are monitored, averaging particle displacements on the periphery of the lower wall (Z-displacements andX-displacements, respectively). Joint aperture is def i ned as the relative normal displacement. A particle is def i ned as belonging to the periphery if the distance from its center to the closest specimen boundary is lower than average diameter. Micro-cracks due to bond breakage, contact force distribution and stress–strain path are monitored during the shear tests. The direct shear tests are run in a large strain mode. As shear displacement increases, new contacts are created along the discontinuity. These contacts are assigned a smooth-joint contact model and the orientation of the smooth joint depends on the location of the contact. A special algorithm was developed to determine which triangle of the surface morphology is intersected by the newly created contact. As shearing occurs, the mirror surface associated with the upper wall does not match the lower surface. Each contact intersects two triangles with possibly different orientations. The orientation of each new contact could be associated in reality with any of the two surfaces or be a combination of the two surfaces. In this study, the assumption was made to consider only the surface morphology associated with the upper wall.

Fig.3.Visualization of the synthetic rock joint sample (upper wall in orange and lower wall in brown): (a) general 3D view; (b) horizontal cross-section through mean fracture plane; (c) vertical cross-section.

Table3Smooth-joint model contact parameters.

2.4.Mechanicalbehaviorofthediscreteinterface

Numerical shear tests under constant normal stress were performed on the SRJ for three values of normal conf i nement. Normal stress values of 0.5 MPa, 1 MPa and 1.5 MPa were applied to the sample which correspond approximately to 0.35%, 0.7% and 1.05% of the intact rock UCS, respectively. Those low values of normal stress corresponding to the order of magnitude of normal stress practitioners usually have to be dealt with for slope stability problems. Fig. 4 shows the evolution of shear stress and normal displacement with shear displacement for one particular packing. It can be seen that the classical elasto-plastic response of rock joints is well captured, thus conf i rming a good agreement with typical behaviors that can be observed experimentally. The mobilized shear stress increases to a peak value as roughness is mobilized and then decreases due to asperity degradation. The peak value def i nes the shear strength of the SRJ (the higher the normal stress, the higher the shear strength). It can be noted that the peak is reached after 1.5 mm of tangential displacement which is on the upper limit of what is usually observed. The stiffness of the smooth joint was probably underestimated and this question will be discussed in Section 3.2. The peak shear strengths were 0.51 MPa, 0.9 MPa and 1.21 MPa. The friction calculated from the ratio of peak shear stress to applied normal stress was higher at lower normal stresses (1.01, 0.90 and 0.81 at 0.5 MPa, 1 MPa and 1.5 MPa, respectively), enhancing a nonlinear strength envelope.

Normal displacement versus tangential displacement curves in

Fig. 4 show that overall dilation of the rock joint is reduced as normal stress increases. As shearing takes place and roughness is mobilized, the dilation angle def i ned as δun/δusincreases to a maximum value Ψpeak(peak dilation angle) at the peak of the shear stress.

Three different specimens were generated varying the packing of the particle assembly. The same shear tests were simulated on each specimen. Peak shear strengths τpeakand peak dilation angles Ψpeakare reported in Fig. 5. Peak shear strength shows a limited sensitivity to particle packing whereas signif i cant variation of peak dilation angle is observed.

The numerical shear tests performed under increasing normal stress def i ne the strength envelope of the model from which a Barton failure criterion (Barton and Choubey, 1977) can be expressed. In Barton’s formulation, the shear strength is expressed as a function of the JRC, joint compressive strength (JCS) and base friction Φb:

where τpis the peak shear stress and σnthe normal stress.

Barton’s failure criterion was applied to predict the strength of the SRJ. Φbrefers to the base friction of the joint which corresponds to the friction a perfectly smooth joint would have. The base friction of the SRJ corresponds to the friction angle of the SJM (φj= 20°). The JCS corresponds to the UCS of the synthetic material, 143 MPa (Table 1), and the JRC of the triangulated surface was estimated to10.4. In Barton’s formulation, the dilation angle can be estimated using the following empirical relation:

Fig.4.Stress and displacement curves of direct shear tests under constant normal stress (ranging from 0.5 MPa to 1.5 MPa) on a 140 mm × 140 mm surface. (a) Shear stress versus tangential displacement. (b) Normal displacement versus tangential displacement.

Fig.5.Comparison of peak shear strength and dilation angles of the synthetic rock joint (diamond) with their respective predicted values from Eqs. (3) and (4) (dash line).

The predictions of shear strength and peak dilation angle of the SRJ with Eqs. (3) and (4) were compared with the numerical results obtained from the simulations and can be seen in Fig. 5. The measured strength appears to be in very good agreement with the prediction obtained with a widely used relationship such as Barton’s failure criterion. Fig. 5 shows some differences between the measured dilation angles and Barton’s predictions. If for the range of normal stress applied in this study, the measured and predicted values are of the same order, the general trend in the decrease of dilation with normal stress is signif i cantly different. An overestimation of the dilation can be expected at high normal stress. However, the SRJ seems to well capture the mechanical behavior of a natural rock joint.

3. Parametric study

The SRJ, as described in Section 2, does not rely on any empirical scheme or any particular assumption on surface roughness. It is generated using a 3D measurement of the surface morphology and intact rock properties.

Such discrete model seems to well capture the effect of surface roughness on the mechanical behavior of rock joints. Results of shear test simulations show a very good agreement with Barton’s prediction, based on JRC. This suggests that predictive estimations of shear strength should be possible combining 3D surface measurements with a smooth-joint contact model. A number of SRJ samples were generated varying the properties of the SJM and the particle friction angle to enhance the relation between some contact properties and the emergent macroscopic behavior. The same surface morphology (140 mm × 140 mm) was introduced and the same particle size distribution was used to generate the particle assembly representing the rock. Scale dependency and particle discretization will be discussed in Section 4. A SJM is de fi ned through fi ve parameters, i.e. friction, cohesion, dilation and stiffness (normal and tangential). Only purely frictional joint was considered at this stage. No local dilation (i.e. at the contact level) was introduced through the SJM as macroscopic dilation (i.e. at the joint level) is expected to be an emergent property of the surface’s roughness. Effects of joint friction angle φj, joint stiffnessesknjandktj, and particle friction φpwere analyzed.

3.1.Effectofjointfrictionangle

Four 140 mm × 140 mm joints were generated with a joint friction angle φjranging from 15°to 30°. Direct shear tests were performed at a constant normal stress of 1.5 MPa. The evolution of shear stress and normal displacement with tangential displacement is shown in Fig. 6. As joint friction increases, peak shear strength and residual shear strength increase from 1.06 MPa to 1.60 MPa and from 0.82 MPa to 1.16 MPa, respectively. The residual shear strength can be characterized by relatively stable shear and normal stresses with degradation on joint surfaces still occurring (Gentier et al., 2000). Direct shear tests by Grasselli (2001) show that residual strength is reached slightly before 3 mm of tangential displacement. In this study, however, residual shear strength was def i ned as the shear stress after 3 mm of tangential displacement.

As φjincreases, the mechanical behavior becomes more brittle. As expected, dilation remains unchanged and emerges as independent of the smooth-joint friction angle. In comparison simulations by Park and Song (2009) exhibited a signif i cant increase in dilation when friction coeff i cient increased from 0 to 0.3 (from 0°to 16.7°).

Mobilized peak friction angles Φpeak(ratio between peak shear stress and normal stress) and mobilized residual friction angles Φres(ratio between residual shear stress and normal stress) can be seen in Fig. 7. Φpeakvaries from 35.4°to 46.8°whereas Φresvaries from 28.8°to 37.8°. A very good linear relation can be drawn between the mobilized friction angles and the smooth-joint friction angle (coeff i cient of determination of 0.999 for Φpeakand 0.987 for Φres). This result is consistent with the idealized decomposition of rock joint strength as the addition of a frictional component and an asperity component as suggested by Barton and Bandis (1982) and conf i rms that the smooth-joint friction angle φjshould be calibrated according to base friction (i.e. friction angle of a planar surface).

Fig.6.Stress and displacement curves of direct shear tests under constant normal stress (σn= 1.5 MPa) on a 140 mm × 140 mm surface with a smooth-joint friction angle ranging from 15°to 30°. (a) Shear stress versus tangential displacement. (b) Normal displacement versus tangential displacement.

3.2.Effectofjointnormalandtangentialstiffness

Fig.7.Peak and residual mobilized friction angle vs. smooth-joint friction angle φjat normal stress of 1.5 MPa. Diamonds and disks represent model values and dash lines best linear fi t.

Table4Normal and tangential stiffnesses for smooth-joint contact model.

Fig.8.Stress and displacement curves of direct shear tests under constant normal stress (σn= 1.5 MPa) on a 140 mm × 140 mm surface with different smooth-joint normal and tangential stiffnesses. (a) Shear stress versus tangential displacement. (b) Normal displacement versus tangential displacement.

3.3.Effectofparticlefriction

Fig.9.Contact shear force distribution across the joint surface at a tangential displacement of 1.4 mm.

Park and Song (2009) suggested that for planar surfaces particle–particle friction coeff i cient (μ = tan φp) should be determined between 0.0 and 0.15 as signif i cant increase in the shear strength is obtained and very little variation of the dilation angle isobserved. In this study, direct shear tests on rough synthetic rock joints were performed under a constant normal stress of 1.5 MPa. Particle friction angles of 0°, 5°and 10°were used, corresponding to friction coeff i cients of 0.0, 0.087 and 0.176, respectively. For each particle friction angle, the micro-properties were calibrated to generate particle assemblies exhibiting the same mechanical behavior (deformation and strength) on unconf i ned compression test. The evolutions of the shear stress and normal displacement upon shearing are shown in Fig. 10.

Peak and residual shear strengths increase with particle friction (φjkept constant). No signif i cant difference materializes on the shear stress curves before a tangential displacement of 1 mm. After 1 mm, micro-cracking (i.e. bond breakage) in the particle assembly becomes signi fi cant as can be seen in Fig. 10. The increase of microcracks leads to the formation of a gouge (materialized as single particles or clusters of several particles) between the two walls of the joint. Differences in the number of clusters forming the gouge appear after a tangential displacement of 1 mm (see Fig. 10d). It can be noted that signi fi cant degradation occurs long after the peak shear stress. The rate of degradation however is maximal at or immediately after the peak. With the formation of a gouge, forces across the interface are transmitted not only through smooth-joint contacts but also through clusters of particles. Shear strength of fi lled joints is highly in fl uenced by the strength of the in fi ll material. As the particle friction angle is increased, the strength of the newly formed gouge is increased. The dilation shows little variationfor the range of particle friction angles studied here. However signif i cant increase can be expected for higher values of φpas stated in Park and Song (2009).

Fig.10.Stress and displacement curves of direct shear tests under constant normal stress (σn= 1.5 MPa) on a 140 mm × 140 mm surface for different particle friction angles (σpranging from 0°to 10°). (a) Shear stress versus tangential displacement. (b) Normal displacement versus tangential displacement. (c) Number of micro-cracks versus tangential displacement. (d) Number of clusters in the gouge vs. tangential displacement.

Fig.11.Scale effects in the shear strength components of non-planar defects (after Bandis et al. (1981)).

4. Large-scale discontinuities

Bandis et al. (1981) identif i ed two contributors to rock joint strength: a basic frictional component (base friction) and a roughness component (Fig. 11).

Geometry (or morphology) of the discontinuity (shape of the asperities) and asperity failure (the strength of the asperities) are the basis of the roughness component. When base friction appears to be scale-independent and can be estimated on laboratory-scale experiments, the roughness component is highly scale-dependent. Roughness decreases as scale increases (Bandis et al., 1981). Numerous studies were carried out trying to quantify the scale dependence of joint strength from which empirical relations were proposed (Barton and Bandis, 1982):

Because of the scale dependency observed in the mechanical behavior of discontinuities, their properties should be assessed at the relevant scale. In a rock mass, the scale of the discontinuities ranges from meters to hundreds of meters (and more). Laboratory methods, where scale is usually restricted to meter and below, cannot be directly extended for fi eld estimates.

4.1.ScaledependencyofSRJbehavior

In this study, scale dependency of the model was investigated performing numerical shear tests on samples of various sizes. Two smaller scales were tested, 70 mm × 70 mm and 46.7 mm × 46.7 mm, splitting the initial surface into respectively four and nine sub-surfaces which were imported into a bonded-particle assembly. Same micro-properties and particle size distribution were used to represent the intact rock. Direct shear tests under a constant normal stress of 1.5 MPa were performed on each of the four + nine newly created synthetic rock joints. Fig. 12 shows peak shear stress (mean value and range of distribution) and peak dilation angle (mean value and range of distribution) versus sample size for all the tests.

Peak shear strength and peak dilation decrease signif i cantly for a sample size increasing from 46.7 mm to 140 mm. Mean peak shear stress dropped from 1.67 MPa to 1.37 MPa, corresponding to a 17.8% decrease. Peak dilation angle dropped from 25.4°to 19.9°, corresponding to a 21.8% decrease. Variability in peak and dilation angles is reduced as sample size increases. Combining the empirical relations Eqs. (3)–(6), predictions on the scale dependency of peak friction angle and peak dilation angle are shown in Fig. 12. The SRJ exhibits a scale dependency of its mechanical properties in good agreement with predictions based on empirical relations. However, peak dilation angles of the SRJ appear slightly lower than those predicted with Barton and Bandis’ relations.

Fig.12.Variation of peak shear stress and peak dilation angle with joint size (σn= 1.5 MPa): mean value (diamond) and variability (plain line). Comparison with combined empirical relations (Eqs. (3)–(6)) (after Barton and Choubey, 1977; Barton and Bandis, 1982).

4.2.Effectofparticlesize

With the development of laser measurement systems and 3D photogrammetry techniques, practitioners can now have access to topological descriptions of large discontinuities (meter and above). These large-scale morphologies could be used to generate large synthetic rock joints. Estimates of their shear strengths could hencebe derived without the need for any empirical relations. However with the current computer limitations, a limited number of particles can be used, and testing these large-scale discontinuities with DEM would require the use of larger particles. As particles become larger than the smallest asperities, surface roughness can be artifi cially reduced. Direct shear tests under constant normal stress (1.5 MPa) were performed on 70 mm × 70 mm samples using different particle size distributions, with an average diameter in the vicinity of the interface ranging from 1.29 mm to 2.56 mm. The shape of the particle size distribution (ratio between minimum and maximum radius of 1.66) and micro-properties were kept unchanged.

Fig. 13 shows the peak shear stress and the peak dilation angle versus the average particle diameter. Dilation is reduced as particle size increases which tends to con fi rm that the surface roughness of the joint is reduced as particles becomes larger. The horizontal spacing between two grid points of the surface is 1.4 mm whereas the lowest average diameter used for these simulations is 1.29 mm, resulting in approximately one particle per interval. No asymptotic value appears when decreasing the particle diameter, thus suggesting that a minimum of two (or more) particles per interval should be introduced to capture the full roughness of the surface. Interestingly, shear strength exhibits no such trend. In comparison, shear test simulations by Park and Song (2009) show no conclusive effect on peak friction angle and peak dilation angle. A full understanding of the effect of particle size with the SJM requires additional analyses.

Fig.13.Evolution of peak shear stress τpand peak dilation angle ψpeakwith average particle diameter. 70 mm × 70 mm sample and σn= 1.5 MPa.

4.3.Signif i canceforlargediscontinuities

The surface shape or topology of a discontinuity is seen at the micro-scale level as a series of asperities and results from two different components of surface texture, roughness and waviness (Belem et al., 2009). The roughness component is termed“secondary” or second order surface roughness, and the waviness component is termed “primary” or fi rst order surface roughness. Both orders of asperities have to be taken into consideration when considering joint roughness and thus joint strength (Plesha, 1987; Yang et al., 2001; Haneberg et al., 2007). The second order asperities exhibit high angles and narrow base lengths (or wave length) in opposition to the fi rst order asperities that have lower angles and longer base lengths. The behavior of rock joints is controlled primarily by the second order asperities during small displacements, and fi rst order asperities govern the shearing behavior for large displacements. Barton and Choubey (1977) fi rst stated that at low normal stress levels, the second order asperity controls the shearing process. With increasing normal stress, the second order asperity is sheared off and the fi rst order asperity takes over as the controlling factor. Fardin et al. (2001) suggested that a resolution of 0.2 mm in the roughness measurement was required to correctly capture the second order asperities whereas a resolution of 20 mm seems suff i cient to capture the fi rst order asperities. Yang et al. (2001) obtained similar conclusions using analytical decompositions.

With the current computer capacities, capturing the effect of the second order asperities for large joints is currently not achievable. However, the fi rst order asperities could be accurately described for joint surfaces of 1 m2and above. Strength characterizations would then be restrained to situations where primary asperities are the controlling factor.

5. Conclusions

In this paper, a new DEM representation of rock joints was presented. Numerical joints were generated combining a real 3D surface morphology and the smooth-joint contact model. Particles lying on opposite side of the joint surface were assigned a smoothjoint contact. At macro-scale level, the behavior of a natural joint is a combination of surface roughness, intact rock properties and frictional contact behavior. The behavior of the synthetic rock joint is an emergent property of the surface morphology, micro-properties of the particle assembly and micro-properties of the smooth-joint contact model. The surface morphology introduced was measured by Grasselli (2001) using a laser scanner. Micro-properties of the particle assembly were calibrated on the basis of measured intact rock properties (UCS, Young’s modulus).

Direct shear tests under constant normal stress were simulated and the mechanical response of the discrete model was analyzed. The shear behavior was compared to the expected behavior of a joint with the same morphology, the latter being assessed with conventional JRC based estimation methods. A relatively good agreement could be established. The effect of roughness was consistently captured throughout the simulations, for various normal stresses, though the model appears to slightly overestimate dilation under high normal stresses.

Surface morphology was responsible for surface roughness and micro-properties of the particle assembly def i ned intact rock behavior; sensitivity analyses were presented to assist with the calibration of the smooth-joint micro-properties. As expected, the smooth-joint friction angle exhibited a direct inf l uence only on peak friction angle and residual friction angle only and should be calibrated according to base friction angle. Dilation was an emergent property of the surface roughness only whereas aperture was inf l uenced by the normal stiffness of the smooth joint. The macroscopic shear stiffness was related to the shear stiffness of the smooth joint.

The main benef i t of this approach is the possibility to assess the shear strength of a discontinuity on the basis of directly measurable properties such as rock mechanical properties (UCS, Young’s modulus), planar surface friction and surface morphology. The model offers the perspective of considering the 3D effect of roughness for strength assessment. Such approach could provide a means of establishing a constitutive behavior that is not reliant on any empirical formulation or classif i cation scheme and def i ning strength parameters for rock joints that can be used in more conventional analyses.

Conf l ict of interest

We wish to conf i rm that there are no known conf l icts of interest associated with this publication.

Acknowledgements

The authors would like to thank Itasca Consulting Group Inc. for their technical advice and Dr. Giovanni Grasselli for providing the natural morphology of a discontinuity in granite whose work was performed with funding provided by the Swiss Federal Off i ce for Water and Geology.

Bandis S, Lumsden A, Barton N. Experimental study of scale effects on the shear behavior of rock joints. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 1981;18(1):1–21.

Barton N, Choubey V. The shear strength of rock joints in theory and practice. Rock Mechanics 1977;10(1–2):1–54.

Barton N. Review of a new shear strength criterion for rock joints. Engineering Geology 1973;7(4):287–332.

Barton N, Bandis SC. Effect of block size on the shear behavior of jointed rocks. In: The 23rd US Symposium on Rock Mechanics. Berkeley: American Rock Mechanics Association; 1982. p. 732–60.

Belem T, Souley M, Homand F. Method for quantif i cation of wear of sheared joint walls based on surface morphology. Rock Mechanics and Rock Engineering 2009;42(6):883–910.

Buzzi O, Hans J, Boulon M, Deleruyelle F, Besnus F. Hydromechanical study of rock-mortar interfaces. Physics and Chemistry of the Earth 2008;32(8–14): 820–31.

Carr JR, Warriner JB. Relationship between the fractal dimension and joint roughness coeff i cient. Bulletin of the Association of Engineering Geologists 1989;26(2):253–63.

Cundall P. Numerical experiments on rough joints in shear using a bonded particle model. In: Lehner FK, Urai JL, editors. Aspects of Tectonic Faulting. Berlin: Springer-Verlag; 2000. p. 1–9.

Deisman N, Ivars MD, Darcel C, Chalaturnyk RJ. Empirical and numerical approaches for geomechanical characterization of coal seam reservoirs. International Journal of Coal Geology 2010;82(3–4):204–12.

Dershowitz W. Interpretation and synthesis of discrete fracture orientation, size, shape, spatial structure and hydrologic data by forward modeling. In: Cook NGW, Goodman RE, Myer LR, Tsang CF, editors. Fractured and Jointed Rock Masses. Netherlands: A.A. Balkema; 1995. p. 579–86.

Fardin N, Stephansson O, Jing L. The scale dependence of rock joint surface roughness. International Journal of Rock Mechanics & Mining Sciences 2001;38(5):659–69.

Fardin N, Feng Q, Stephansson O. Application of a new in situ 3D laser scanner to study the scale effect on the rock joint surface roughness. International Journal of Rock Mechanics & Mining Sciences 2004;41(2):329–35.

Gentier S, Riss J, Archambault G, Flamand R, Hopkins D. Inf l uence of fracture geometry on shear behavior. International Journal of Rock Mechanics & Mining Sciences 2000;37(1–2):161–74.

Giacomini A, Buzzi O, Krabbenhoft K. Modelling the asperity degradation of shear rock joint using FEM. In: The 8th World Congress on Computational Mechanics (WCCM8) and the 5th European Congress on Computational Methods in Applied Sciences and Engineering (ECCOMAS 2008); 2008.

Grasselli G. Shear strength of rock joints based on quantif i ed surface description. PhD thesis. Switzerland: Ecole Polytechnique Federale de Lausanne; 2001.

Grasselli G, Wirth J, Egger P. Quantitative three-dimensional description of a rough surface and parameter evolution with shearing. International Journal of Rock Mechanics & Mining Sciences 2002;39(6):789–800.

Grasselli G, Egger P. Constitutive law for the shear strength of rock joints based on three-dimensional surface parameters. International Journal of Rock Mechanics & Mining Sciences 2003;40(1):25–40.

Haneberg WC. Directional roughness prof i les from three-dimensional photogrammetric or laser scanner cloud points. In: Rock Mechanics: Meeting Society’s Challenges and Demands, Proceedings of the 1st Canada-US Rock Mechanics Symposium. Vancouver, Canada: Taylor & Francis; 2007. p. 101–6.

Hans J, Boulon M. A new device for investigating the hydromecanical properties of rock joints. International Journal for Numerical and Analytical Methods in Geomechanics 2003;27(6):513–48.

Hutson R, Dowding C. Joint asperity degradation during cyclic shear. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 1990;27(2):109–19.

Itasca Consulting Group Inc. . Particle fl ow code in 3 dimensions version 4.0: User manual. Minneapolis, MN: Itasca Consulting Group Inc; 2008.

Karami A, Stead D. Asperity degradation and damage in the direct shear test: a hybrid FEM/DEM approach. Rock Mechanics and Rock Engineering 2008;42(2):229–66.

Kulatilake P, Malama B, Wang J. Physical and particle fl ow modeling of jointed rock block behavior under uniaxial loading. International Journal of Rock Mechanics & Mining Sciences 2001;38(5):641–57.

Lambert C, Darve F, Nicot F. Rock slope stability from microscale to macroscale level. In: Pietruszczak P, editor. Numerical Models in Geomechanics. Ottawa: A.A. Balkema; 2004. p. 115–20.

Lambert C. Variability and uncertainty on rock mass strength via a synthetic rock mass approach. In: 1st Southern Hemisphere International Rock Mechanics Symposium. Perth, Australia: Australian Centre for Geomechanics; 2008. p. 355–66.

Lambert C, Coll C. Santiago A DEM approach to rock joint strength estimates Santiago. In: Rock Slope Stability in Open Pit Mining and Civil Engineering; 2009.

Lambert C, Buzzi O, Giacomini A. Inf l uence of calcium leaching on the mechanical behavior of a rock mortar interface: a DEM analysis. Computers and Geotechnics 2010;37(3):258–66.

Lee H, Park Y, Cho T, You K. Inf l uence of asperity degradation on the mechanical behavior of rough rock joints under cyclic shear loading. International Journal of Rock Mechanics & Mining Sciences 2001;38(7):967–80.

Misra A. Effect of asperity damage on shear behavior of single fracture. Engineering Fracture Mechanics 2002;69(17):1997–2014.

Park JW, Song JJ. Numerical simulation of a direct shear test on a rock joint using a bonded-particle model. International Journal of Rock Mechanics & Mining Sciences 2009;46(8):1315–28.

Patton FD. Multiple modes of shear failure in rock. In: Proc. 1st Int. Cong. Rock Mech. Lisbon: International Society for Rock Mechanics; 1966. p. 509–13.

Pierce M, Cundall P, Potyondy D, Ivars MD. A synthetic rock mass model for jointed rock. In: Rock Mechanics: Meeting Society’s Challenges and Demands, Proceedings of the 1st Canada-US Rock Mechanics Symposium. Vancouver, Canada: Taylor & Francis; 2007. p. 341–9.

Pine R, Owen D, Coggan J, Rance J. A new discrete fracture modelling approach for rock masses. Geotechnique 2007;57(9):757–66.

Plesha M. Constitutive models for rock discontinuities with dilatancy and surface degradation. International Journal for Numerical and Analytical Methods in Geomechanics 1987;11(4):345–62.

Potyondy DO, Cundall PA. A bonded-particle model for rock. International Journal of Rock Mechanics & Mining Sciences 2004;41(8):1329–64.

Rogers S, Kennard D, Dershowitz W, van As A. Characterising the in situ fragmentation of a fractured rock mass using a discrete fracture network approach. In: Rock Mechanics: Meeting Society’s Challenges and Demands, Proceedings of the 1st Canada-US Rock Mechanics Symposium. Vancouver, Canada: Taylor & Francis; 2007. p. 137–43.

Saeb S, Amadei B. Modelling rock joints under shear and normal loading. International Journal of Rock Mechanics & Mining Sciences 1992;29(3):267–78.

Tse R, Cruden C. Estimating joint roughness coeff i cients. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 1979;16(5):303–7.

Yang ZY, Di CC, Yen KC. The effect of asperity order on the roughness of rock joints. International Journal of Rock Mechanics & Mining Sciences 2001;38(5):745–52.

Yu X, Vayssade B. Joint prof i les and their roughness parameters. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts 1991;28(4):333–6.

Zhao Z, Jing L, Neretnieks I. Particle mechanics model for the effects of shear on solute retardation coeff i cient in rock fractures. International Journal of Rock Mechanics & Mining Sciences 2012;52:92–102.

Zhao Z. Gouge particle evolution in a rock fracture undergoing shear: a microscopic DEM study. Rock Mechanics and Rock Engineering 2013;46(6):1461–79.

C. Lambert received his ME in Hydraulic Engineering (2000) and PhD in Geomechanics (2005) from Institut National Polytechnique de Grenoble. In 2006, he joined the Division Exploration and Mining at the Commonwealth Scientif i c and Industrial Research Organisation (CSIRO) as a postdoctoral fellow. He joined the University of Canterbury in 2010 as a lecturer in geotechnical engineering and rock mechanics.

*Corresponding author.

E-mailaddresses:cedric.lambert@canterbury.ac.nz (C. Lambert),

ccoll@golder.co.nz (C. Coll).

Peer review under responsibility of Institute of Rock and Soil Mechanics, Chinese Academy of Sciences.

Shear strength

DEM simulations

Smooth-joint model