Multi-peak deformation behavior of jointed rock mass under uniaxial compression: Insight from particle flow modeling

Engineering Geology 213 (2016) 25–45

Contents lists available at ScienceDirect

Engineering Geology journal homepage: www.elsevier.com/locate/enggeo

Multi-peak deformation behavior of jointed rock mass under uniaxial compression: Insight from particle flow modeling Cheng Cheng a, Xin Chen b,⁎, Shifei Zhang b a b

Key Laboratory of Shale Gas and Geoengineering, Institute of Geology and Geophysics, Chinese Academy of Sciences, Beijing 100029, China State Key Laboratory for Geomechanics & Deep Underground Engineering, China University of Mining and Technology (Beijing), Beijing 100083, China

a r t i c l e

i n f o

Article history: Received 28 October 2015 Received in revised form 28 July 2016 Accepted 21 August 2016 Available online 25 August 2016 Keywords: Jointed rock mass Multi-peak deformation behavior Joint geometry Joint strength mobilization Smooth joint contact model

a b s t r a c t Evaluation of mechanical properties of jointed rock mass is very important for geological engineering problems related to underground opening, slope stability, dam foundation and hydraulic fracturing, etc. It is generally known that both strength and deformability of rock mass may change with joint geometry and arrangement. By PFC (Particle Flow Code) modeling, this paper aims at understanding the damage mechanism of the multipeak deformation behavior of jointed rock mass with different joint orientations and continuity factors under uniaxial compression observed in the physical model tests. With the smooth joint contact model, not only has the inherent roughness in conventional joint model been avoided, but also we can monitor the joint aperture, number of closed joints, and forces acting on the joints during the damage process. It is shown that mobilization of joint strength plays an important role in the multi-peak deformation behavior of rock mass. If the joint strength is slightly mobilized or immobilized (i.e., most joints keep open) during the test, strain softening occurs. However, if the joint strength is significantly or fully mobilized, yielding platform or strain hardening may happen. © 2016 Elsevier B.V. All rights reserved.

1. Introduction It is generally known that understanding the mechanical behavior of jointed rock mass is very important for the geological engineering such as rock excavation in underground or open pit mining, tunneling and dam foundation construction of hydropower station, wellbore maintaining and hydraulic fracturing in oil and gas exploration, etc. In the past decades, lots of studies have been focused on the strength estimation, deformation characteristics, and failure mechanism of rock mass with persistent joints. Nevertheless, according to the ISRM suggested method (ISRM, 1978), joint persistence implying the areal extent or size of a joint within a plane, is a very important parameter when considering the safety of rock tunnels, dam foundations, and slopes, especially if the joints are unfavorably oriented for stability. Much more work is required to be done on the behaviors of rock mass containing non-persistent joints, depicted as the co-planar joints separating volumes of intact rock (Prudencio and Van Sint Jan, 2007). In practical engineering, some empirical methods have been developed to estimate the strength and deformability of the jointed rock mass. For example, deformation modulus of rock mass (EM) can be assessed by its empirical relations with the values of Rock Mass Rating (RMR) (Bieniawski, 1978; Serafim and Pereira, 1983), Q (Barton, 2002), or Geological Strength Index (GSI) system (Hoek et al., 2002; Cai et al., 2004). These indirect methods are widely used for practical ⁎ Corresponding author. E-mail address: [email protected] (X. Chen).

http://dx.doi.org/10.1016/j.enggeo.2016.08.010 0013-7952/© 2016 Elsevier B.V. All rights reserved.

engineering design. Nevertheless, they cannot be expected to always provide accurate values of deformation modulus owing to the complicated geometrical and mechanical properties of the joints. In order to study the strength and deformation behavior of the jointed rock mass, direct test method is required. Compared with the in-situ test, physical model test is an economical and practical method to investigate mechanical behavior of rock masses. There are two ways to simulate the rock mass with non-persistent joints: (1) existing fissures are placed when the rock-like material is molded; (2) a large amount of modeling blocks are assembled to simulate rock mass with different sets of persistent or non-persistent joints. Different loading conditions were applied on these two types of samples including uniaxial, poly-axial compression (Bobet and Einstein, 1998; Brown, 1970; Chen et al., 2013, 2012; Horii and Nemat-Nasser, 1985; Ladanyi and Archambault, 1972; Prudencio and Van Sint Jan, 2007; Singh and Seshagiri Rao, 2005; Singh and Singh, 2008; Tiwari and Rao, 2006; Wong and Einstein, 2009) and direct shear (Gehle and Kutter, 2003; Ghazvinian et al., 2012; Lajtai, 1969; Li et al., 2014). In addition, large size of samples containing non-persistent joints were also applied in physical model tests to study the responses of tunnels (Gong et al., 2015, 2013; He, 2011; He et al., 2010; Singh et al., 2009) and slopes (Alejano et al., 2011). Most physical model tests focused on the influence of the joints arrangement on the failure modes and strength characteristics of rock mass, by observation on the macro failure processes and mechanical behaviors. Furthermore, some studies attempted to establish the strength criterion for rock mass containing non-persistent joints based on the observed failure mechanisms (Cording and Jamil,

26

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

1997; Jennings, 1970; Lajtai, 1969). For example, Lajtai (1969) established a combined failure envelope with three parts, i.e., tension, shear failure of rock bridges, and shear failure of the crushed material matrix, under negative or low, intermediate, and high normal stresses, respectively. In his study, a mobilization factor of the joint, C1 (0 ≤ C1 ≤ 1), was proposed to consider the load carrying capacity of the non-persistent joints with different apertures, whose friction resistance may be fully or partially mobilized when the joints are closed, or even immobilized when the joints keep open during the tests. Only a few physical model tests investigated the dependence of the stress-strain behavior on the joints arrangement. Tiwari and Rao (2006) studied the post-peak deformation of jointed blocks under true triaxial compression. The samples used in their tests were assemblies of cubes with one set of non-persistent joints and two sets of persistent joints, and their studies show the effect of both joint geometry and σ2/σ3 on the slope of post-peak stress-strain curves of the jointed rock mass. Prudencio and Van Sint Jan (2007) researched the failure modes and strength criterion of rock models with a set of non-persistent joints under uniaxial and biaxial compressive loads. It can be observed that the stress-strain curves have different features for the tested samples with different joint configurations under various confining pressures. However, the features of deformation behaviors and their relationship with the failure processes were not analyzed. What is more, the range of joint continuity factor k in their models was quite limited, k = 0.71 for most series of the samples and 0.63 for one series of the samples. Chen et al. (2011, 2012, 2013) set a wide range of the two geometrical parameters, joint inclination angle (0° to 90°) and joint continuity factor (0.2 to 0.8) in their physical models under uniaxial compression tests, and the stress-strain curves of the specimens show different types of deformation behaviors, depending on the combined variation of the two geometrical parameters. Compared with the conventional strain softening behavior after the single peak stress for the uniaxial compression test on the intact model, it was found that there may be one or several peaks for the jointed specimens. With different joint configurations, yielding platform or strain hardening may occur before the last peak followed by the final strain softening, i.e., the final strain softening may be delayed by different types of ductile behavior. In some other cases, there may also be one or several small peaks or oscillations during the general trend of strain softening. In this study, the term “multi-peak deformation behavior” is used to refer to this kind of deformation behaviors influenced by different joint configurations, and more detailed classification of the deformation behaviors will be described in Section 2.2. Chen et al. (2011, 2012, 2013) attempted to explain this kind of behaviors based on the damage processes, but it was mainly according to observations on the failure patterns of the samples. In order to learn more about underlying mechanism of the deformation characteristics and failure phenomena observed in the physical model tests, numerical studies have been applied in many researches and some novel insights have been provided. With various numerical methods such as FEM (Wasantha et al., 2014; Wong et al., 2001; Zhang et al., 2006), BEM (Gehle and Kutter, 2003), DEM (Halakatevakis and Sofianos, 2010; Jiang et al., 2009; Kulatilake et al., 1992; Park et al., 2004), NMM (An et al., 2014) and hybrid FEM/DEM approach (Karami and Stead, 2008; Morris et al., 2006), the numerical models can be built to simulate and analyze the stress redistribution and fracturing mechanism, which cannot be easily measured or characterized in the physical model tests. However, it should be noted that, the crack initiation and propagation direction around the joint tips are very complex and are always simplified in these numerical simulations. In recent years, PFC (Particle Flow Code), a DEM method developed by Itasca Co. Ltd., has shown many advantages in damage mechanism analyses of rock mass. A PFC model is generated by a group of particles with varied size, deformation modulus, normal or shear bonding strength, etc., which can be bonded together or broken up dependent on the strength of each bond, therefore no macro constitutive law or strength criterion is required while yielding, strain softening, fracture

initiation and propagation can all be simulated. It has been widely employed to build models of rock mass with non-persistent joints (Bahaaddini et al., 2013b; Fan et al., 2015; Ghazvinian et al., 2012; Park et al., 2004; Scholtès and Donzé, 2012). Especially with the development of smooth joint contact model (Bahaaddini et al., 2014, 2013a, 2013b; Chiu et al., 2013; Esmaieli et al., 2010; Hadjigeorgiou et al., 2009; Lambert and Coll, 2014; Mas Ivars et al., 2011), PFC modeling has been proved to be a promising method to make a better understanding on the mechanism of mechanical behavior of jointed rock mass. However, most of the previous PFC studies on rock mass with non-persistent joints focused on the mechanical parameters such as peak strength, Young's modulus and Poisson's ratio, as well as the failure patterns, while the stress-strain behaviors and the corresponding micromechanical mechanism, especially the opening or closing of joints, were not considered effectively. By PFC3D, Fan et al. (2015) and Yang et al. (2016) gave a successful simulation on the strength behavior and failure modes of the jointed specimens researched by Chen et al. (2011, 2012, 2013). They also presented the calculated stress-strain relationships, however, the damage mechanism was not studied. In this study, intending to understand the damage mechanism on multi-peak deformation behavior of these jointed specimens under uniaxial compression, smooth joint contact model in PFC is found to be suitable for some quantitative analyses. It is noted that these jointed specimens have geometries of plate, hence PFC2D is applied in this work for saving computer memorize. With numerical analysis, the multi-peak deformation behaviors of rock masses dependent on joint orientation, continuity factor, and aperture are studied.

2. Laboratory physical model experiments 2.1. Geometry and properties of the specimens The geometry of specimen, geometrical parameters of joints and the photograph of a sample are presented in Fig. 1. The specimen has the geometry of plate with a size of 150 mm × 150 mm × 50 mm, and a single set of non-persistent joints penetrates through the thickness of the specimen. A global coordinate system ox1x2x3 is related to the specimen dimension, where x2 is direction of thickness. In addition, a local coordinate system oζ1ζ2ζ3 is related to the arrangement of non-persistent joints, where ζ2 shows the direction of specimen thickness and cannot be plotted in Fig. 1(b), while ζ3 presents the direction of joint plane. For the non-persistent joints, s is joint spacing; c is the joint center distance, which is the distance between the centers of two neighboring joints in a joint plane; Lj is the length of each single joint; β is the joint inclination angle, defined by the angle between the joint plane (oζ3 direction) and the horizontal plane (ox3 direction); and k is the continuity factor, defined by the ratio of the length of a single joint to the joint center distance, i.e., k = Lj/c. In this study, joints are regularly arranged with constant joint spacing s = 30 mm and constant joint center distance c = 30 mm. The combined variation of the continuity factor k (k = 0.2, 0.4, 0.6, and 0.8) and the joint inclination angle β (β = 0°, 15°, 30°, 45°, 60°, 75°, and 90°) is considered, and the corresponding joint length Lj is 6 mm, 12 mm, 18 mm, and 24 mm. For comparison, the intact specimen is also studied, which is presented as k = 0. The mixture of plaster and water with a ratio of 1:0.6 by weight was casted into the mold made of PMMA (polymethyl methacrylate) to prepare the specimens. To form the pre-existing non-persistence joints, a group of 0.3 mm-thick nickel steel sheets were inserted into the mixture through the precut slots in the upper and bottom PMMA plates, and then were removed after the setting of the liquid mixture started. The specimens were cured for 21 days in room temperature before being subjected to mechanical testing. In order to observe repeatable results, the mixing, molding and curing of the material were carefully controlled, and at least three samples with the same joint configuration

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

27

Fig. 1. Sketches of (a) specimen geometry, and (b) joint geometry. (c) photograph of a sample (k = 0.8 and β = 45°). (Chen et al., 2011)

were prepared for the tests. The physical and mechanical properties of the intact specimen are listed in Table 1. 2.2. Test setup and results Uniaxial compressive tests were carried out in a servo controlled INSTRON 8506 test machine. Displacement control was applied with a constant loading rate of 0.0025 mm/s. Fig. 2 shows the complete axial stress-strain curves of the specimens with different joint geometries. It can be seen that all the curves contain the early three stages before the first peak, i.e., the upward concave curve induced by micro-crack closure of the matrix, linear elastic stage and plastic yielding stage. Thereafter, they exhibit different stages and four types of brittle or ductile behaviors after the first peak have been classified (see Fig. 3 for their typical curves with descriptions of key points in Table 2): Table 1 Physical and mechanical properties of the intact specimen (Chen et al., 2011). Physical and mechanical properties 3

Unit weight γ (kN/m ) UCS (MPa) Young's Modulus E (GPa) Cohesion c (MPa) Internal friction angle φ (°) Joint friction angle φ (°)

Values 11.58 8.51 2.56 2.2 38 37

Type I: strain softening. After the first peak stress (which is also the maximum peak stress in this case), the axial stress decreases until the end of the test with the increasing strain. This type of curve shows the typical brittle behavior that similar to the complete stress-strain curve of intact rock. Type II: general strain softening with oscillations. It is similar to Type I, and the only difference is that there are one or several peak stress points during the general process of strain softening, but they are always lower than the first peak stress. This type of curve shows the ductile behavior to some extent. Type III: yielding platform – strain softening. With the increasing strain, the stress keeps almost constant for a period before decreasing until the end of the test. The maximum peak stress may occur during the yielding platform. The ductile behavior is very obvious in this type of curve. Type IV: yielding platform – strain hardening – strain softening. With the increasing strain, the stress keeps almost constant for a period, and then increases again to the maximum peak value followed by decreasing until the end of the test. This type of curve also shows obvious ductile behavior. 3. Numerical models 3.1. A brief introduction of PFC2D code and smooth joint contact model PFC2D is a DEM code by which the macro plastic behavior and failure process of rock can be simulated based on the micro responses of the particles such as deformation, breaking up, rotation, slip etc. under the

28

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

10

10

=0°

=15° 8

k=0 k=0.2 k =0.4 k=0.6 k=0.8

6 4

(MPa)

(MPa)

8

4

2 0 0.0

k=0 k=0.2 k=0.4 k=0.6 k=0.8

6

2

0.5

1.0

1.5

2.0

2.5

3.0

0 0.0

3.5

0.5

1.0

1.5

(%)

2.0

2.5

3.0

3.5

(%) 10

10

=45°

=30° 8 k=0 k=0.2 k=0.4 k=0.6 k=0.8

6 4

(MPa)

(MPa)

8

4 2

2 0 0.0

k=0 k=0.2 k=0.4 k=0.6 k=0.8

6

0.5

1.0

1.5

2.0

2.5

3.0

0 0.0

3.5

0.5

1.0

1.5

(%)

2.0

2.5

3.0

10

10

=75°

=60° 8

k=0 k=0.2 k=0.4 k=0.6 k=0.8

6 4 2

(MPa)

(MPa)

8

0 0.0

3.5

(%)

k=0 k=0.2 k=0.4 k=0.6 k=0.8

6 4 2

0.5

1.0

1.5

2.0

2.5

3.0

3.5

(%)

0 0.0

0.5

1.0

1.5

2.0

2.5

3.0

3.5

(%)

10

=90°

(MPa)

8

k=0 k=0.2 k=0.4 k=0.6 k=0.8

6 4 2 0 0.0

0.5

1.0

1.5

2.0

2.5

3.0

3.5

(%) Fig. 2. Axial stress-strain curves of the specimens obtained from the laboratory tests (Chen et al., 2011).

loading or unloading conditions (Potyondy and Cundall, 2004). In PFC2D models, an assembly of particles can be bonded together, and the bond will break if the tensile or shear strength is reached, forming a tensile or shear crack. The movement and deformation of particles follow the basic Newton's Law and force-displacement law, respectively.

Consequently, the strain hardening or softening can be simulated without a prescribed constitutive model or strength criterion. During the calculation, cracking process and chains of force distribution can both be monitored to analyze the mechanism of the macro deformation or failure behavior.

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

Type I 3

B

B

(MPa)

4 2

Type II

F1( F)

R

F2

2

1

(MPa)

1

6

4

F1( F)

8

29

S

R

S

1 A

A

0 O 0.0

0.5

1.0

1.5

2.0

0 O 0.0

2.5

0.5

1.0

(%) 1

1.5

2.0

2.5

(%) 1

6 4

F2( D)

3 R S

2

F1

2

R

D

B

S

1

1 A 0 O 0.0

Type IV

F2 ( F)

3

(MPa)

B

1

4

Type III

F1( F)

1

(MPa)

5

0.5

1.0

1.5 1

2.0

2.5

A 0 O 0.0

3.0

0.5

1.0

(%)

1.5 1

2.0

2.5

3.0

(%)

Fig. 3. Four types of axial stress-strain curves (Chen et al., 2012).

There are usually three ways to simulate the joints in PFC2D models: (1) creating joint planes with different geometries by modifying the properties of the contacts existing between the balls that fall on opposite sides of the planes (Fig. 4(a)); (2) creating a band of at least 3-particles width with weak properties (Fig. 4(b)); and (3) deleting one or several thin blocks of particles to form fissures or non-persistent open joints (Fig. 4(c)). The main problem for these methods is the inherent joint roughness induced by the circular particles with various sizes and arrangements, in addition, for method (1) the open joint with different apertures cannot be considered; for method (2) joint closing is limited owing to a certain width of the joint band filled with particles; for method (3) the joint aperture is always too large to be closed. Fortunately, smooth-joint contact model provides a novel way for joint simulation (Fig. 4(d)), by which the shortcomings of the former methods can be avoided. On one hand, the particle pairs joined by a smooth joint contact can slide through each other along the joint orientation regardless of the local contact orientation, instead of being forced to move around one another (Han et al., 2012; Mas Ivars et al., 2011). On the other hand, an aperture a (a N 0 indicates that the joint is open) can be given to the smooth joint contact, and joint opening/closing as well as the forces on the joint surfaces can both be monitored during the loading process. Table 2 Description of the key points in the stress-strain curves. Key points

Description

O A B D F F1 F2 R S

Start of test End of micro-crack closure stage and start of linear elastic stage End of linear elastic stage End of yielding platform Maximum peak stress First peak Last peak A point near the end of test End of test

The smooth joint contact behaves like a set of elastic springs existing uniformly over a cross section between a pair of particles, and the cross section has the center at the contact point, the orientation parallel with the joint plane, as well as the area A defined as:
 A ¼ 2Rt; R ¼ λ min RðAÞ ; RðBÞ

ð1Þ

where, R and t are the smooth joint radius and thickness (t = 1.0 in PFC2D models), respectively, and λ is the radius multiplier (λ = 1.0, usually); R(A) and R(B) are the radii of the pair of particles. During the calculation, the force F acting on the contact and the relative translational displacement U between the two particles can be ^ j and tangential unit vector decomposed along the normal unit vector n ^t j (presented in Fig. 4(d)) as follows: ^ j þ F s^t j F¼F n n

ð2Þ

^ j þ U s^t j U¼U n n

ð3Þ

where, Fn and Un are the magnitudes of force and displacement at the normal direction of the joint plane (Fn N 0 means compressive force, and Un N 0 indicates particles overlapping), respectively; Fs and Us are the magnitudes of force and displacement at the shear direction, respectively. For the unbounded smooth joint, the contact force-displacement law follows Coulomb sliding with dilation. In each cycle, the forces update according to the displacements increment, and the relationship is elastic if Un N 0 (otherwise the displacement change cannot induce any force variation because the joint remains open) as follows: F n :¼ F n þ knj AΔU en 0

F s :¼ F s −ksj AΔU es

ð4Þ ð5Þ

30

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

a)

c)

b)

y x

d) Fig. 4. Four methods to create joint in the rock model. Note that the joint plane is presented as the boundary between the yellow and green particles in (a), the band of green particles with weak properties in (b) and a series of non-persistent slots in (c). The blue and red lines in (d) show the joint plane and smooth joint contacts, respectively. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

dilation angle of the joint (ψj):

where, knj and ksj are the normal and shear stiffness of the smooth joint contact, respectively; ΔUenand ΔUes are the elastic portion of normal and shear displacement increments, respectively. It should be noted that the updated shear force F's in each cycle should be checked by a maximum value F⁎s given by: F s ¼ F n μ j

ð6Þ

j F s j ¼ F s

ð7Þ


 F n :¼ F n þ ΔU s tanψ j knj A

ð8Þ

where, ΔU⁎s is the shear displacement increment during the sliding process. The smooth joint can also be bonded, and the bonding strength follows the coulomb criterion with a tension cut-off. The force-displacement relationship is also elastic before the bond breaks.

where, μj is the friction coefficient of the smooth joint. If |F's| ≤F⁎s , then F's = Fs; otherwise, sliding happens and the shear force keeps constant while the normal force may be updated because of the

B

B

vx

B

vy

y

y

y x

A

x

A

b) step 2000

x

c) step 4000

0.18

0.24 Horizontal moving

Vertical moving

0.12 Normal force Shear force Joint aperture

0.06

0.12

0.00 0.00 -0.06 Sliding

-0.12 Joint open

Joint closed (particle overlapping)

-0.12

Joint aperture (mm)

Normal and shear forces (N)

a) step 0

A

-0.18 0

1000

2000

3000

4000

Timestep

d) Fig. 5. “Direct shear test” on double particles with a smooth joint. In (a)–(c), the red line indicates the smooth joint contact, and the black bar shows the resultant contact force. (d) the evolution of normal and shear force components on the joint as well as the joint aperture during the whole test. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

31

Table 3 Micro-properties of particles and bonds. Ball parameters

Parallel bond parameters

Ball density (kg/m3)

1158

Minimum ball radius (mm) Ball radius ratio Rmax/Rmin Contact modulus Ec(GPa) Coefficient of friction Normal to shearing stiffness ratio kn/ks

0.6 1.66 2.0 0.5 2.5

Bond modulus Ec (GPa) Normal bond strength (MPa) S.D.⁎ normal bond strength (MPa) Shearing bond strength (MPa) S.D. shearing bond strength (MPa) Normal to shearing bond stiffness ratio kn /ks

2.0 5.5 2.1 8.0 3.3 2.5

⁎ S.D.: standard deviation.

To demonstrate the behavior of the smooth joint, a “direct shear test” on two particles (with the same radius of 1.0 m) with a smooth joint contact was simulated, and the evolution of forces on the joint surface was monitored. Fig. 5 shows the model, where the smooth joint is parallel with the x axis with an original aperture of 0.2 mm. The normal and shear stiffness of the joint are both 1.0 N/m3, and the joint friction coefficient is 1.0. During the test, particle A is fixed, while particle B moves downward at a constant velocity (vy = − 8.0 × 10−5 m/s) for 2000 steps, thereafter moves horizontally to the left at a constant velocity (vx = −8.0 × 10−5 m/s) for 2000 steps, see Fig. 5(a)–(c). It should be noted that the model is calculated with a time step algorithm, so the velocities here are not exactly the same as those in the physical world. The change of the joint aperture, the normal and shear forces along with steps were monitored and drawn in Fig. 5(d). In the first 2000 steps, the normal force remains zero until the joint is closed, thereafter it starts to increase in an linear elastic way; during the following 2000 steps, the normal force keeps the constant while the shear force increases elastically to the maximum value until sliding begins. This simple example illustrates that with the smooth joint contact model, joint opening/ closing and the forces on the joint can be monitored during the whole test process, and this supplies a method to analyze the micro mechanism of the rock mass with non-persistent joints in the following study. 3.2. Model calibration and validation In order to build a valid PFC2D model for the further research compared with the laboratory tests, it is required to calibrate the microproperties of particles, parallel bonds, and smooth joints so as to obtain a numerical model with comparative macro-properties with the physical models. The PFC2D model is built with a size of 150 mm × 150 mm, and a series of uniaxial and biaxial numerical tests are carried out on the intact specimens with an axial strain rate of 1.0 s− 1. Note that this strain rate in time stepping algorithm of this PFC2D model is equal to the axial displacement rate of about 2.864 × 10−8 m per time step, so it is low enough to guarantee that the specimen is always in quasi-static equilibrium. Finally,a calibrated model with 10,347 particles is obtained, and the micro-properties of particles and parallel bonds for the matrix material are listed in Table 3. A comparison of macro-properties of the intact specimen by the numerical test and physical test is presented in Table 4. It should be noted here that the calibrated internal friction angle (30°) is a little lower than the measured value (38°) of the physical model. As pointed in references (Itasca Co. Ltd., 2004; Potyondy and

Table 4 Calibrated mechanical properties of the intact specimen by numerical models.

UCS (MPa) Young's Modulus E (GPa) Cohesion c (MPa) Internal friction angle φ (°)

Experimental

Numerical

8.51 2.56 2.2 38

8.05 2.57 2.0 30

Cundall, 2004), conventional PFC models give unrealistic lower friction angle owing to the rotation of circular particles. Although irregular grains modeled by cluster, clump or Grain Based Model (GBM) may increase the macro friction angle, they will also bring difficulties in observing the damage process owing to the complex parameters and complicated distributions of inter-grain and inner-grain cracks. Fortunately, the conventional PFC model can capture the key features of failure modes and stress-strain behaviors of the physical models. In addition, in uniaxial compressive tests dominated by tensile failure, the influence of internal friction angle is limited. Consequently, it is reasonable to use this numerical model to analyze the main characteristics and mechanism of the multi-peak deformation behaviors in this study. According to the results of the physical model test presented in Section 2, the joint friction angle is available and the joint dilation angle can be considered as zero because no apparent joint asperities have been found in the physical model test, but the normal and shear stiffness are not known as no direct shear testing studies on the joints have been reported. However, try and error method can be applied during the validation of the numerical models with non-persistent joints to solve this problem. All cases of the non-persistent joint arrangements in the physical model tests are calculated in the numerical model validation. For illustration of the numerical model, the case of β = 15°, k = 0.8 is shown in Fig. 6. In the models, non-persistent joints consist of several co-planar joint segments separated by the solid bridges. For each joint segment, it is composed by a series of smooth joint contacts, which are defined as joint elements in this paper. The specimens with different joint inclination angles and joint continuity factors have different numbers of joint segments and joint elements, which are listed in Table 5. With the validated parameters of smooth joint presented in Table 6, a series of axial stress-strain curves are obtained and presented in Fig. 7. It is observed that the results for the jointed specimens obtained by numerical tests have covered the basic characteristics of those by the physical tests. There are also four types of brittle or ductile behaviors, similar to those presented in the physical tests. Fig. 8 shows the corresponding typical axial stress-strain curves found in the numerical tests. The key points marked here (Fig. 8) are in accordance with those in the physical model tests (Fig. 3), and it should be noted that point O and A coincide in the numerical tests because there is no micro-crack closure in the densely packed and well-connected assembly of particles. It is observed that these curves are not exactly the same as those from the physical model tests, e.g., the physical models have larger axial strains, more stress vibrations and lower peak strengths in some cases. The above differences between the results of the physical and numerical tests may be attributed to the following two factors: (1) the physical model is a porous material while the PFC2D model is built with an assembly of densely packed particles; (2) 2D numerical models cannot behave exactly like the 3D physical model tests. These are also the reasons why the validated joint aperture (0.14 mm) in the PFC2D models is lower than that in the physical models (b0.3 mm, because of the shrinkage of the material after removing the inserted sheet). However, generally speaking, it is believed that the numerical models are reasonable and valid enough to conduct some further researches compared with the physical model tests.

32

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

Joint segment

Joint segment

b)

Joint element

a)

c)

Fig. 6. Numerical model of the specimen with non-persistent joints (β = 15°, k = 0.8). The blue lines indicate the joint planes, and the red lines show the smooth joint contacts, which are defined as the joint elements in this model. (b) shows a magnified joint segment from (a), and some of the joint elements are zoomed in and presented in (c). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

4. Mechanism of multi-peak deformation behavior influenced by joint persistence and orientation In the physical model test, it is not easy to measure the forces distribution in the matrix and on the joint surfaces, as well as the joint opening or closing process. The mechanism of multi-peak deformation behavior can only be speculated according to relationship between the observed failure patterns and the corresponding types of stress-strain curves. Fortunately, by using PFC2D model with smooth joint, four types of stress-strain curves and the multi-peak deformation behavior due to the variation of joint inclination angle and joint continuity factor have been reproduced, as shown in Figs. 7 and 8. In PFC2D models, force chains distribution and the joint aperture evolution can both be monitored at each time step, which are very useful to understand the mechanism of the multi-peak deformation behavior of the jointed specimens with different joint orientations and persistence. As we know, with the increase of the axial strain in the uniaxial compression test, the aperture of each joint element may change, and some of the joint elements may close. Since the closed joint elements can bear loads while the open ones cannot, the mobilization of the strength of a specific joint segment depends on the close or open status of all joint elements in the joint segment. Namely, if none of

the joint elements are closed, the strength of the joint segment is immobilized; if all the joint elements are closed, the strength of the joint segment is fully mobilized; if only some of the joint elements are closed, the strength of the joint segment is partially mobilized. As an example, Fig. 9 presents a partially closed joint segment (from the case of β = 0° and k = 0.6) where closed joint elements are located in the center to right part of the joint segment, and the force can only be transmitted through these closed joint elements, which leads to partially mobilization of the strength of the joint segment. Therefore, the evolution of force chains distribution during the test is helpful for understanding the role of the joints played in the whole process of ductile behavior. In addition to the force chains distribution, some other parameters are proposed and monitored to describe the behavior of the joint system in the specimens. The parameters are given as follows: (1) Average joint aperture (aave). It is defined as the sum of the joint element apertures (ai) divided by the total number of joint elements (N):

aave ¼

N X

ai =N

ð9Þ

i¼1

Table 5 Number of joint segments and joint elements for each specimen in the PFC2D models. β

k

0° 0° 0° 0° 15° 15° 15° 15° 30° 30° 30° 30° 45° 45°

0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4

Number of joint segments

Number of joint elements

β

k

25 25 25 25 25 25 29 29 25 25 29 31 25 25

167 354 526 699 179 340 507 675 179 362 532 696 173 352

45° 45° 60° 60° 60° 60° 75° 75° 75° 75° 90° 90° 90° 90°

0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8

Number of joint segments

Number of joint elements

25 25 25 27 29 31 25 27 29 29 25 25 25 25

517 695 181 352 528 704 177 342 506 675 176 355 528 706

The increasing or decreasing of ai indicates the joint element is opening or closing, respectively. Hence aave is used to evaluate whether the general trend of the joints is opening or closing; (2) Number of closed joint elements (Njc). Joint element aperture ai N 0 means the joint element is open, while ai ≤ 0 means it is closed. In each time step, the number of open joint elements

Table 6 Validated parameters for the smooth contact joint. Parameters

Values

Joint friction angle φj (°) Joint dilation angle ψj (°) Joint normal stiffness (N/m3) Joint shear stiffness (N/m3) Joint aperture (mm)

38 0 1.0 × 1012 1.0 × 1012 0.14

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

10

33

10 = 0°

4

(MPa)

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

6

1

= 15°

8

4 2

2 0 0.0

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

6

1

(MPa)

8

0.3

0.6

0.9

1.2

1.5

0 0.0

1.8

0.3

0.6

(%) 1

0.9

1.2

10

= 45°

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

1

6 4

(MPa)

8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

6

1

(MPa)

8

2

4 2

0.3

0.6

0.9

1.2

1.5

0 0.0

1.8

0.3

0.6

(%) 1

0.9

1.2

10 8

4 2 0.3

0.6

0.9 1

1.2

1.5

1.8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

6

1

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

6

(MPa)

8

(MPa)

1.8

= 75°

= 60°

1

1.5

(%) 1

10

0 0.0

1.8

10 = 30°

0 0.0

1.5

(%) 1

4 2 0 0.0

0.3

(%)

0.6

0.9 1

1.2

1.5

1.8

(%)

10

= 90°

1

(MPa)

8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

6 4 2 0 0.0

0.3

0.6

0.9 1

1.2

1.5

1.8

(%)

Fig. 7. Axial stress-strain curves of the specimens in the numerical tests with the joint aperture a = 0.14 mm.

(Njo) and the number of closed joint elements (Njc) can be counted and the summation of them equals the total number of joint elements (N) in the specimen, which can be given by:

N ¼ N jc þ Njo

ð10Þ

Since N keeps constant during loading on a given specimen, changes in the number of closed joint elements Njc can be served as an index to estimate the overall close state of the joint system.

(3) Ratio of closed joint elements (Rjc). It is defined as the number of closed joint elements (Njc) divided by the total number of joint elements (N):

Rjc ¼ N jc =N

ð11Þ

Therefore Rjc is a normalized index to estimate the close status of the joint system.

34

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

10

10

a) Type I F1(F)

(MPa)

1

6

B

4

6

F2 B

0.6 1

R

2

S

2 0.4

F1(F)

4

R 0 O(A) 0.0 0.2

b) Type II

8

1

(MPa)

8

0.8

1.0

1.2

S

O(A) 0 0.0 0.2

1.4

0.4

(%)

0.6

0.8

10

10

1.4

8

(MPa)

8 6

F1

4

F2(D/F)

1

(MPa)

1.2

d) Type IV

c) Type III

1

1.0

(%)

1

6

F2(F)

4

F1

B R

2 0 O(A) 0.0 0.2

0.4

0.6 1

0.8

1.0

1.2

2

S 1.4

(%)

D

R

B

0 O(A) 0.0 0.2

S 0.4

0.6 1

0.8

1.0

1.2

1.4

(%)

Fig. 8. Four typical axial stress-strain curves in the numerical test.

(4) Normalized normal and shear forces acting on the joints (FN nj and FN sj). The normal and shear forces on the joints (Fnj and Fsj) are defined as the sum of normal or shear contact forces on all joint elements in a specimen, respectively. For comparison, the normalized values are adopted. Here, Fmax, the maximum normal force on the joints in the specimen with β = 0° and k = 0.8 (Type IV), is set as the basic value for normalization:

FN nj ¼ F nj =F max

ð12Þ

FN sj ¼ F sj =F max

ð13Þ

With the above four parameters, i.e., average joint aperture (aave), number of closed joint elements (Njc), ratio of closed joint elements (Rjc), and normalized normal and shear forces acting on the joints (FN nj and FN sj), the deformation, close status and contact forces of the joint system in a specific specimen can be characterized averagely. Namely, the strength mobilization of the joint system can be quantitatively evaluated.

Fig. 9. Partially mobilization of the strength of a joint segment. The black bars indicate the resultant contact forces, with the thickness proportional to the force magnitude.

Based on the evolution of force chains distribution and the monitored parameters mentioned above, the mechanism of each type of ductile behavior is analyzed. Moreover, a comparison analysis is made and the effect of joint orientation and persistence on ductile behavior of rock masses is discussed in the end of this section. 4.1. Mechanism of Type I behavior Based on the stress-strain curves presented in Fig. 7, both the specimens with high inclination joint angles (β = 60°–90°) and the specimens with low inclination angles (β = 0°–30°) at low continuity factors (k = 0.2–0.4) may show Type I behavior, i.e., typical strain softening behavior. However, in the numerical study, it is found that the two kinds of specimens have different damage mechanisms. Therefore, Type I stress-strain curve can be further classified into two sub-types, namely, Type Ia and Ib, according to their different mechanisms. 4.1.1. Mechanism of Type Ia behavior The jointed specimen with β = 90° and k = 0.8 is taken as an example to study the mechanism of Type Ia behavior. Three key points are selected from the stress-strain curve of Type Ia in Fig. 8(a) to show the evolution of force chains distribution during the test (Fig. 10(a)–(c)). It can be seen that the contact forces distribute almost uniformly over the matrix of the specimen when the axial stress reaches peak strength (point F1, Fig. 10(a)), and then becomes localized in the matrix during the post-peak stage (point R and S, Fig. 10(b) and (c)), implying that the matrix is gradually damaged. During the whole test, very few contact forces across the joint elements can be observed, meaning that very little joint strength has been mobilized. The photograph of the fractured specimen in Fig. 10(d) also shows that failure occurs mainly in the matrix and all the joints are open. More detailed information of the joint behavior can be obtained from the data monitored in Fig. 11. The entire test process can be divided into three stages according to the varied mobilization of joint strength: (1) Stage 1 coincides with the elastic stage of the stress-strain curve. In this stage, the curve of average joint aperture (aave) indicates that the joints open slightly, and no closed joint elements or contact forces on the joints are observed. A few tensile cracks initiated from

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

35

(b) Point R

(a) Point F1

(c) Point S

(d)

Fig. 10. (a)–(c) Evolution of force chains distribution during the test ( k = 0.8 and β = 90°), and the selected key points correspond to those of Type Ia in Fig. 8(a). (d) the failure of the physical model.

10 1

1

8

0.25

4

0.20

1

6

(MPa)

0.30

Stage 1

Stage 2

8

40

2

20 Stage 2

Stage 3

Stage 3

0 0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

6

Fnj Fsj

(d)

1 N

1.0 0.8

N

0.6 0.4

4

N

4

(MPa)

80 60

10

Stage 2

1

6

Stage 1

20

Fj

1

Rjc

4

10

100

Rjc (%)

(MPa)

8

(c)

1

1

1

30

0 0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

0.10 0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

10

6

Stage 1

Stage 3

50 40

Njc

2

0.15

2

(b)

1

Njc

1

aave

10

0.35

aave (mm)

(MPa)

8

(a)

1

0.2 2

0.0 Stage 1

Stage 2

Stage 3

0 -0.2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

Fig. 11. Evolution of (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal and shear forces acting N on the joints (FN nj and Fsj) with increasing axial strain (k = 0.8 and β = 90°).

36

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

10

(a) Type I

1

(MPa)

8 6

F1(F)

Type I

B

Type I

a b

4 2

R S

0 O(A) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1 Fig. 12. Typical axial stress-strain curves for Type Ib behavior (k = 0.2 and β = 15°) compared with Type Ia behavior.

the joint tips as well as some disperse micro-cracks in the matrix can be observed in this stage; (2) In stage 2, the decrease of aave shows that joints have a minor closing trend, implying that a little dilation occurs in the solid matrix. Some joint elements with a low percentage are closed, and very few forces on joints are induced. In this stage, with linkage between the tensile cracks and the pre-existing joints, the axial stress reaches the maximum value, thereafter strain softening occurs with the further fracturing of the solid matrix; (3) In stage 3, aave

(a) Point F1

(c) Point S

increases significantly indicating a salient opening of the joints. With the severer damage of the matrix, more joint elements are getting closed, however, the percentage of the closed joint elements is so low and the mobilization of joint strength is so limited that the general process of strain softening cannot be changed. 4.1.2. Mechanicsm of Type Ib behavior The jointed specimen with β = 15° and k = 0.2 is taken as an example to study the mechanism of Type Ib behavior. Type Ib behavior has similar strain softening curves (Fig. 12) with Type Ia, while there is something different during the damage process of the matrix. The contact forces are also distributed almost uniformly in the material matrix at peak strength (point F1, Fig. 13(a)), and then turn to be localized in the post-peak stage (point R and S, Fig. 13(b) and (c)). However, the contact forces can be observed across some joint elements, meaning that these joint elements are closed and can carry loads, i.e. joint strength has been mobilized partially. We can also see closed joints in the physical test (Fig. 13(d)). Nevertheless, the matrix has been seriously damaged at this point. The axial stress-strain curve can also be divided into three stages based on the mobilization of joint strength (Fig. 14): (1) stage 1 also almost coincides with the elastic stage, however, the decreasing aave indicates a closing trend of the joints, but no closed joint elements have been observed and accordingly there is no forces transmitted though the joint system. A few wing cracks originated from the tips of the inclined joints and some disperse micro-cracks in the matrix can be observed in this stage; (2) in stage 2, the joints continue closing and some joint elements are closed. With the coalescence of the tensile or

(b) Point R

(d)

Fig. 13. (a)–(c) Evolution of force chains distribution during the test ( k = 0.2 and β = 15°), and the selected key points correspond to those of Type Ib in Fig. 12. (d) Photograph of the fractured physical model. Note that the red circles in (b) and (c) mark some joint segments bearing contact forces, and in (d) show some joint segments closed.

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

aave

8

-

4

0.05

2 Stage 2

Stage 3

10 1

(MPa)

20

Stage 2

Stage 1

0 0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

1

(MPa)

2

20

-

8

Fnj

6

Fsj

(d)

1

1.0

N

0.8

N

0.6 0.4

N

1

40 4

Fj

40

Rjc (%)

4

60

6

Stage 3

60

Ncj

2

0 0.00 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1 100 10 (c) 1 1 Rjc 80 8

6

(b)

1

Ncj

0.10

1

8

0.15

6

Stage 1

80

10 1

aave (mm)

(MPa)

0.20

(a)

1

(MPa)

-

1

1

10

37

4 0.2

Stage 1 Stage 2

Stage 3

0 0 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

2

0.0 Stage 1

Stage 2

Stage 3

0 -0.2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 (%) 1

Fig. 14. Evolution of (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal and shear forces acting N on the joints (FN nj and Fsj) with increasing axial strain (k = 0.2 and β = 15°).

(a) Point F1

(b) Point R

(c) Point S

(d)

Fig. 15. (a)–(c) Evolution of force chains distribution during the test (k = 0.6 and β = 45°), and the selected key points correspond to those of Type II in Fig. 8(b). (d) Photograph of the corresponding fractured physical model. Note that the red circles in (b) and (c) mark some joint segments bearing contact forces, and in (d) show some joint segments closed. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

38

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

shear cracks in the solid matrix, the axial stress reaches the maximum strength. More joint elements turn to be closed after the peak strength, but the normal and shear forces carried by them are so limited that the axial stress keeps dropping; (3) in stage 3, because of the failure of the matrix and movement of the fractured blocks, there is a general trend of joint opening. Although there are some joint elements closed, the mobilized strength of joint system cannot compensate the reduced capacity of the damaged matrix, hence strain softening continues. It should be noted that Type Ia behavior occurs if the joints have a high inclination angles (β = 60°–90°) and the joint systems of the specimens are generally open under vertical compression, while Type Ib occurs when the joints have low inclination angles (β = 0°–30°) and the joint continuity factor is not very high (k = 0.2–0.4), therefore the joints are mainly closing but only a limited number of them can be closed. In both cases, the loads are mainly carried by the material matrix, and strain softening occurs with the gradually damage of the matrix. However, in Type Ib behavior, there is a little more mobilization of joint strength that induces a slower drop of axial stress compared with Type Ia behavior (Fig. 12).

According to the mobilization of joint strength, the axial stress-strain curve (shown in Fig. 16) is divided into three stages as follows: (1) stage 1 almost covers the pre-peak stage. The joints keep closing in this stage, but no joint elements are closed. It is found that tensile and shear cracks initiated from the solid bridges in this stage; (2) in stage 2, with the formation of shear planes resulted from the linkage between the tensile or shear cracks and the pre-existing co-planar joints, an increasing number of joint elements are found to be closed after the peak strength. After a while, the joints start to open with the propagation of tensile cracks in the solid matrix. During the general process of joint opening, some joint elements are closed, and the curve of aave decreases a little as well. The increasing number of closed joint elements and forces on joints show a mobilization of joint strength, and that is why the axial stress increases to form a secondary peak during the general trend of strain softening; (3) in stage 3, with the further fracturing of the solid matrix, the joints continue to keep opening, and the forces carried by the joints also decreases, so the axial stress-strain relationship continues to show strain softening behavior. 4.3. Mechanism of Type III behavior

4.2. Mechanism of Type II behavior The jointed specimen with β = 45°and k = 0.6 is taken as an example to study the mechanism of Type II behavior. As presented in Section 2.2, the axial stress-strain curve of Type II is similar to that of Type I, and the main difference is that there are one or several secondary peak stress points during the general process of strain softening. This feature can also be explained by the force distribution and mobilization of joint strength during the loading. The evolution of contact force chains for the case of k = 0.6 and β = 45° is presented in Fig. 15. The contact forces are mainly distributed around the solid bridges and very few joint elements can carry loads at peak strength (point F1, Fig. 15(a)), thereafter the solid matrix becomes more and more damaged and the contact forces turn to be localized. During the damage process, part of the joint elements are closed and can carry some loads (Fig. 15(b) and (c)). That is why the axial stress goes up a little; however, this cannot change the general trend of strain softening. 10

0.5 0.4

4

0.2 0.1 Stage 1

0 0.0

Stage 2

0.2

0.4

10 1

0.6 0.8 (%) 1

1.0

90

4

60

2

0.0 1.2

100

6

4

40

2

20

0.2

0.6 0.8 (%) 1

1.0

0 1.2

0.4

Stage 3

0.6 0.8 (%) 1

1.0

0 1.2

1

1

-

1.0

N

0.8

N

0.6

Fnj Fsj

(d)

1

0.4 4

N

0.4

Stage 2

Fj

0.2

Stage 3

(MPa)

60

Rcj (%)

6

0 0.0

Stage 1

10 8

Stage 2

30

0 0.0

80

Stage 1

120

6

(c)

1

Rcj

8 (MPa)

-

Stage 3

(b)

1

Ncj

8 (MPa)

0.3

150

-

Ncj

6

2

1

10 1

aave (mm)

(MPa)

(a)

1

aave

8

1

-

1

1

The jointed specimen with β = 0°and k = 0.6 is taken as an example to study the mechanism of Type III behavior. Fig. 17 shows the evolution of contact forces distribution and the photograph of the fractured physical model. Apparently, the contact forces distribute mainly around the solid bridges while very few contact forces can be observed acting on the joints at point F1 (Fig. 17(a)). This condition changes totally at point F2, when the contact forces transmit through lots of joint elements (Fig. 17(b)). This evolution demonstrates that at the beginning of the yielding platform most of the joints are still open and the stresses are mainly carried by the solid bridges, thereafter along with the gradual damage of the material matrix, more and more joint elements are closed and therefore higher loads can be carried by the joints. Subsequently with the further failure of the matrix, the contact forces turn to be localized (Fig. 17(c)). The photograph in Fig. 17(d) shows the fractured physical models with closed joints. The analysis presented above is accordance with the data monitoring in Fig. 18, which is also divided into three stages according to the

0.2 2 0 0.0

0.0 Stage 1

0.2

Stage 2

0.4

Stage 3

0.6 0.8 (%) 1

1.0

-0.2 1.2

Fig. 16. Evolution of (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal and shear forces acting N on the joints (FN nj and Fsj) with increasing axial strain (k = 0.6 and β = 45°).

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

39

(a) Point F1

(b) Point F2

(c) Point S

(d)

Fig. 17. (a)–(c) Evolution of force chains distribution during the test (k = 0.6 and β = 0°), and the selected key points correspond to those of Type III in Fig. 8(b). (d) Photograph of the corresponding fractured physical model.

1

0.06

2

0.03

0.2

Stage 3

0.4

1.0

0.6 0.8 (%) 1

100

10

80

8

100

4

40

2

20

0.4

1

-

Fsj

0.6

0.8 (%) 1

1.0

1.2

0 1.4

(d)

1

1.0

N

0.8

N

0.6

Fnj

6

Stage 3

0.4 4

N

60

0.2

Stage 2

Fj

6

Rcj (%)

(MPa)

200

0 0.0

(c)

1

Rcj

8

1

-

4 2

0.00 1.4

1.2

500 400

Stage 1

10 1

Ncj

300

(MPa)

0 0.0

Stage 2

(MPa)

4

(b)

1

6

1

0.09

-

Ncj

6

Stage 1

8

0.12

aave (mm)

(MPa)

(a)

1

aave

8

1

-

1

1

600

10

0.15

10

0.2

Stage 1

0 0.0

0.2

Stage 2

0.4

0.6

0.8 (%) 1

0.0 Stage 1

Stage 3

1.0

2

1.2

0 1.4

0 0.0

0.2

Stage 2

0.4

0.6 0.8 (%) 1

Stage 3

1.0

1.2

-0.2 1.4

Fig. 18. Evolution of (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal and shear forces acting N on the joints (FN nj and Fsj) with increasing axial strain (k = 0.6 and β = 0°).

40

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

(a) Point F1

(b) Point F2

(c) Point S

(d)

Fig. 19. (a)–(c) Evolution of force chains distribution during the test (k = 0.8 and β = 0°), and the selected key points correspond to those of Type IV in Fig. 8(b). (d) Photograph of the corresponding fractured physical model.

10

Stage 2

10

0.4

-

1

Rcj

(MPa)

100

4

40

2

20

600

4

400

2

200 Stage 1

Stage 2

0.2

0.4

Stage 3

0.6 0.8 (%) 1

0 1.2

1.0

1

8 6

-

1.0

N

0.8

N

0.6

Fnj Fsj

(d)

1

0.4 4

N

60

800

10

80

6

Ncj

6

0 0.0

0.00 1.2

(c)

1

(b)

1

Fj

1

1.0

Rcj (%)

(MPa)

8

Stage 3

0.6 0.8 (%) 1

1

0.03

(MPa)

(MPa)

2

-

Ncj

0.06

aave (mm)

4

0.2

8

0.12 0.09

Stage 1

1000 1

6

0 0.0

10

0.15

(a)

1

aave

8

1

-

1

1

0.2

0 0.0

Stage 1

0.2

Stage 2

0.4

Stage 3

0.6 0.8 (%) 1

1.0

0 1.2

2 0 0.0

0.0 Stage 1

0.2

Stage 2

0.4

Stage 3

0.6 0.8 (%) 1

1.0

-0.2 1.2

Fig. 20. Evolution of (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal and shear forces acting N on the joints (FN nj and Fsj) with increasing axial strain (k = 0.8 and β = 0°).

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

0.4

Ia

II

Ib

aave (mm)

II III IV

0.2

Ia

Ib

Ib

0.8 (%) 1

III

1.0

1.2

Ia

1.4

N

II

Ib 0.2

1.0

Ia

0.8

Ib

0.6

II III IV

0.4

0.6 0.8 (%) 1

1.0

1.2

1.4

(d)

IV

0.4 0.2

II

Ia

0 0.0

Fnj

Rcj (%)

0.0 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 (%) 1 100 Ia (c) 80 Ib IV II 60 III III IV 40

0.6

IV

II III IV

100

III

20

(b)

200

IV

0.4

500

300

Ib

0.2

Ia

400

0.1

0 0.0

600

Ncj

0.3

(a)

41

III

Ia Ib

II 0.0 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 (%) 1

Fig. 21. (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal force acting on the joints (FN nj) during the tests on all the cases studied in Section 4.

mobilization of joint strength: (1) In stage 1, the joints keep closing under the uniaxial load, but no joint elements are closed. It can be found cracks emanated from both the tips and middle part of the preexisting joints; (2) In stage 2, although the solid bridges continue to be damaged by the propagation of tensile cracks, a large number of joint elements are closed and they can carry considerable normal forces (but very low shear force as the joints are perpendicular to the loading direction), which shows a significant mobilization of joint normal strength and it can support the specimen to remain its strength for a while, so a yielding platform occurs; (3) In stage 3, with the more serious damage of the material matrix, some joints turn to be open and lower forces can be carried by the joint system, therefore the axial stress drops and strain softening occurs. 4.4. Mechanism of Type IV behavior The jointed specimen with β = 0° and k = 0.8 is taken as an example to study the mechanism of Type IV behavior. The evolution of contact forces distribution and the photograph of the fractured physical model are presented in Fig. 19. Similar to the Type III, the contact forces distribute mainly around the solid bridges at point F1 (Fig. 19(a)), and then transfer to be carried by almost all the joints (Fig. 19(b)).Thereafter with the failure of the material matrix, the contact forces turn to be localized (Fig. 19(c)). The photograph in Fig. 19(d) shows the fractured physical models with closed joints. Comparing the force distribution with the monitoring of joint behavior (Fig. 20), more details on this type of axial stress-strain relationship can be obtained. Three stages are also divided according to the mobilization of joint strength: (1) in stage 1, all the joint elements are still open. Cracks initiated from the tips and middle part of the joints are observed. Due to the high persistence of joints, the lengths of solid bridges are very short and they are much easier to fail than the specimen with low persistence joints. As a result, the axial stress reaches the yielding point at a very low value; (2) in stage 2, the solid bridges are seriously damaged by the propagation of tensile cracks. However, almost all the joint elements are closed, and they can carry the loads which are high enough to increase the strength of the specimen to an even higher value,

therefore a strain hardening follows the yielding stage; (3) in stage 3, with more serious failure happening in the matrix, combined capacities of the matrix and joints decreases and the axial stress drops with the increasing axial strain, which leads to strain softening. 4.5. Discussions Based on the analyses presented above, the curves of average joint aperture (aave), number of closed joint elements (Njc), ratio of closed joint elements (Rjc), and normalized normal force acting on the joints (FN nj) during the loading for the four types of deformation behavior are put together in Fig. 21 for comparison. It is found that mobilization of joint strength plays an important role in the ductile behaviors of the specimen containing non-persistent joints under uniaxial compression. Table 7 presents the types of axial stress-strain curves observed in the numerical tests for all kinds of the specimens with different joint orientations and persistence under uniaxial compression. For the specimens with medium to high inclination angles (β = 45°– 90°), it is easy for the joints to open under uniaxial loading, so very limited joint strength can be mobilized. With the damage of the material matrix, strain softening occurs after reaching the peak strength, which is the typical Ia behavior. It should be noted that in some cases with higher continuity factor (k = 0.6 to 0.8), some joints may be closed and axial stress will increase a little during the general process of strain softening, and this is the main characteristics of typical II behavior.

Table 7 Ductile transition behaviors of the specimens with different joint arrangements. k= 0.2

k = 0.4

β = 0°

Ib

Ib

k = 0.6

β = 15°

Ib

Ib

β = 30°

Ib

Ib

Ib

β = 45°

Ia

β = 60°

Ia

Ia

Ia

β = 75°

Ia

Ia

Ia

β = 90°

Ia

Ia

Ia

k = 0.8

II Ia

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

For the specimens with lower inclination angles (β = 0°–30°), the joints are easier to be closed under uniaxial compression. If the continuity factor is low, i.e., k = 0.2 or 0.4, the number of closed joint elements and mobilization of joint strength are limited, so strain softening occurs with the damage of material matrix, showing Type Ib behavior; If the continuity factor is higher, i.e., k = 0.6,

0 0.0

0 0.0

1.8

1.5

1.8

0.3

0.6

=60°

8 6

1

4

1.5

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0.3

0.6

10

=75°

8

1

4

6

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

1.5

1.8

0.3

0.6

1.5

1.8

=90°

8

1

4

10

6 4

0 0.0

0 0.0

0 0.0

0.9 1.2 (%)

1.5

1.8

(MPa)

6

0.3

0.6

0.9 1.2 (%)

1

1.5

1.8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0 0.0

1.8

0.3

0.6

0.9 1.2 (%)

1.5

1.8

1

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

=45°

8

1.5

6

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0 0.0

1.8

0.3

0.6

0.9 1.2 (%)

1.5

1.8

1

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

1.5

10

=60°

8 6

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0 0.0

1.8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

0.9 1.2 (%)

1.5

0.3

0.6

10

0.9 1.2 (%) 1

=75°

8 6

1.5

1.8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0 0.0

1.8

0.3

0.6

=90°

0.9 1.2 (%)

1.5

1.8

1

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

=90°

8

4 2

1

0.6

6

2

0.6

0.3

8

2 0.3

1.8

2

10 k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

(MPa)

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

=30°

1

1

6

0.9 1.2 (%)

0.9 1.2 (%)

=75°

1

(MPa)

=90°

1.5

4

0 0.0

8

0.6

6

0 0.0

10

0.3

8

0 0.0

0.9 1.2 (%)

1.5

2

10

2

1

0.9 1.2 (%)

8

4

2

0.6

=60°

6

2 0.3

10

1

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

0.6

1

(MPa)

6

1.8

1

=75°

8

1.5

0.9 1.2 (%) 1

1

(MPa)

10

0.9 1.2 (%)

0.6

8

0 0.0

1.8

0.3

10

0 0.0

1

0.3

2

0 0.0

0 0.0

1.5

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

2

0.9 1.2 (%)

0.9 1.2 (%)

6

2

0.6

0.6

=45°

2 0.3

0.3

8

1.8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

2

10

1

6

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

(MPa)

=60°

8

0.9 1.2 (%) 1

1.8

1

2

1

10

0 0.0

1.8

(MPa)

0 0.0

0.9 1.2 (%)

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0 0.0

1.5

1

6

2

0.6

=45°

8

1.5

4

1

(MPa)

10

2

6

1

4

0.3

1.8

1

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

1

6

1.5

=30°

1

(MPa)

=45°

8

0.9 1.2 (%)

=15°

8

4

0 0.0

0.6

0.9 1.2 (%) 1

6

0 0.0

0.3

0.6

8

4

0.9 1.2 (%)

2

0.3

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

0.6

10

(MPa)

1

6

1.8

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

(MPa)

=30°

8

1.5

0 0.0

10

(MPa)

0.9 1.2 (%) 1

0.3

1

1

0.6

2

1.5

0 0.0

1.8

4

2 0.9 1.2 (%)

1.5

(MPa)

0.3

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

2

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

=15°

6

2 0.6

0.9 1.2 (%) 1

8

4

0 0.0

1.8

(MPa)

(MPa)

1.5

4

0.3

0.6

1

=30°

0.3

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

(MPa)

0.9 1.2 (%) 1

1

1

1.8

1

0.6

6

1

=15°

6

(MPa)

0.3

8

(MPa)

1.5

2

10

(MPa)

0.9 1.2 (%)

8

2

1

0.6

6

1

10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4

0 0.0

0.3

1

1

(MPa) 1

1.8

(MPa)

=15°

1.5

1

0.9 1.2 (%) 1

=0°

8

1

0 0.0

0.6

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

2

0 0.0

0.3

a = 0.0mm (closed) 10

4 2

6

1

6

1

4

=0°

8

2

8

(MPa)

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

2

10

1

a = 0.1mm 10

(MPa)

6

1

4

(MPa)

1

(MPa)

6

=0°

8 (MPa)

=0°

8

a = 0.14mm 10

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

(MPa)

a = 0.2mm 10

there are enough closed joint elements to carry loads and remain the axial stress for a period while the material matrix is damaged, so a yielding platform occurs before strain softening, namely, it shows Type III behavior; If the continuity factor is even higher, i.e., k = 0.8, the solid bridges with smaller sizes lead to a lower yielding axial stress, but there are more joints closed to carry the loads and

(MPa)

42

6

k=0 k = 0.2 k = 0.4 k = 0.6 k = 0.8

4 2

0.3

0.6

0.9 1.2 (%)

1

1.5

1.8

0 0.0

0.3

0.6

0.9 1.2 (%) 1

1.5

1.8

Fig. 22. Axial stress-strain curves of the numerical uniaxial compressive tests on the specimens with different joint inclination angles, continuity factors, and apertures (a = 0.2 mm, 0.14 mm, 0.1 mm, and 0.0 mm).

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

43

Table 8 Types of the axial stress-strain curves for the jointed specimens that varied with joint inclination angles, continuity factors, and apertures (a = 0.2 mm, 0.14 mm, 0.1 mm, and 0.0 mm) under uniaxial compression.

k = 0.2

k = 0.4

k = 0.6

0.2

0.14

0.1

0.0

0.2

0.14

0.1

0.0

0.2

β = 0°

Ib

Ib

Ib

Ib

Ib

Ib

Ib

Ib

III

β = 15°

Ib

Ib

Ib

Ib

Ib

Ib

Ib

Ib

β = 30°

Ib

Ib

Ib

Ib

Ib

Ib

Ib

Ib

β = 45°

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

β = 60°

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

β = 75°

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

β = 90°

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

The analyses in the previous sections show that it is necessary to investigate the influence of joint aperture systematically, and it is conducted by using the calibrated PFC models mentioned above with different original joint apertures. Axial stress-strain curves of the numerical uniaxial compressive tests on the specimens with different joint inclination angles (β = 0°, 15°, 30°, 45°, 60°, 75°, and 90°), continuity factors (k = 0.2, 0.4, 0.6, 0.8, and 1.0 (intact)), and apertures (a = 0.2 mm, 0.14 mm, 0.1 mm, and 0.0 mm) are presented in Fig. 22. Types of the axial stress-strain curves of the specimens are classified and listed in Table 8. If the joint is vertical (β = 90°), joint aperture has very little influence on the deformation characteristic of the specimens as the joint strength is difficult to be mobilized, and all the cases show Type Ia behavior. For the cases that the joint inclination angle β b 90°, if the joint continuity factors are low (k = 0.2, 0.4), the joint aperture has very limited influence on the deformation characteristic of the specimens because the mobilization of joint strength is quite limited. They almost show Type Ia or Ib behavior depending on β ≥ 45° or β b 45°, respectively. With increasing joint continuity factors (k = 0.6, 0.8), the joint aperture has more significant influence on the deformation behaviors. If 45° ≤ β ≤ 75°, lower aperture leads to Type Ia behavior, while Type II behavior may happen in some cases with higher aperture; if 0° ≤ β ≤ 30°,

0.0

0.2

0.14

0.1

0.0

Ib

Ib

Ib

Ib

Ib

Ib

Ib

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ia

Ib

Ib

Ia

Ia

Ia

Ia

ductile behavior transfers from Type Ib to IV and III with the increasing joint aperture. One typical example is the case of k = 0.8 and β = 15°, which shows Type Ib behavior for a = 0.0 mm, Type III behavior for a = 0.1 mm and 0.14 mm, and Type II behavior for a = 0.2 mm. The axial stress-strain curves are presented in Fig. 23. Curves of average joint aperture (aave), number of closed joint elements (Njc), ratio of closed joint elements (Rjc), and normalized normal force acting on the joints (FN nj) during the tests on the case of k = 0.8 and β = 15° with different apertures (a = 0.2 mm, 0.14 mm, 0.1 mm, and 0.0 mm) are shown in Fig. 24. When a = 0.0 mm, as most of the joints are closed initially, strain softening occurs after peak strength (Type Ib behavior); with increasing aperture (a = 0.1 mm, 0.14 mm), yielding occurs with the joint closing and damage of the matrix, followed by strain hardening when most of the joint elements are closed and higher loads can be carried (Type IV behavior); If the aperture is too high, say a = 0.2 mm, it takes longer time for the joints to be closed, and the material matrix may be seriously damaged before enough joint strength can be mobilized to carry higher stress, and strain hardening is suppressed, consequently, the ductile behavior transfers from Type IV to Type III. 6. Conclusion In this study, PFC2D analysis is successfully used to give an insight on the multi-peak deformation behavior of rock mass with non-persistent joints under uniaxial compression observed in the physical model test.

10 a=0.2mm a=0.14mm a=0.1mm a=0.0mm

8 (MPa)

5. Influence of joint aperture on brittle or ductile behaviors of the jointed specimen

0.1

Ib

1

strain hardening occurs after the yielding platform, which presents Type IV behavior. According to Table 7, for the specimens with medium to high joint inclination angles (β = 45°–90°), the ductile behavior transfers from Type Ia to II with the increasing continuity factor. For the specimens with lower joint inclination angles (β = 0°–30°), the behavior transfers from Type Ib to III and finally to IV with the increasing continuity factor. This conclusion is accordance with the analysis presented above. The different deformation behaviors of the jointed specimens with the variation of the joint inclination angle and joint continuity factor show that the mobilization of the joint strength plays an important role in the multi-peak deformation behavior of rock mass with non-persistent joints under uniaxial compression. In the cases that the joint strength cannot be mobilized, the rock mass will have a strain softening behavior after the peak strength point along with the damage of the matrix. However, if the joint strength is partially or fully mobilized, the strain softening stage will be changed. Joint aperture is an important factor affecting the strength mobilization of the joints. Therefore, the influence of joint aperture will be studied in the next section.

0.14

k = 0.8

6 4 2 0 0.0

0.3

0.6

0.9

1.2

1.5

1.8

(%) 1 Fig. 23. Axial stress-strain curves of the specimens with β = 15°and k = 0.8 at different joint apertures (a = 0.2 mm, 0.14 mm, 0.1 mm, and 0.0 mm).

44

C. Cheng et al. / Engineering Geology 213 (2016) 25–45

700

0.20

aave (mm)

0.15 0.10

a=0.2mm a=0.14mm a=0.1mm a=0.0mm

600 500 Ncj

a=0.2mm a=0.14mm a=0.1mm a=0.0mm

400 300

0.05

200

0.00

(a) 0.3

0.6

0.9 1.2 (%) 1

1.5

(b) 0.3

0.6

0.9 1.2 (%) 1

100

1.0

a=0.2mm a=0.14mm a=0.1mm a=0.0mm

0.8

Fnj

N

0.6 0.4

1.5

1.8

a=0.2mm a=0.14mm a=0.1mm a=0.0mm

80 60 40 20

0.2 0.0 0.0

0 0.0

1.8

Rcj (%)

-0.05 0.0

100

(c) 0.3

0.6

0.9 1.2 (%) 1

1.5

1.8

0 0.0

(d) 0.3

0.6

0.9 1.2 (%) 1

1.5

1.8

Fig. 24. Evolution of (a) average joint aperture (aave), (b) number of closed joint elements (Njc), (c) ratio of closed joint elements (Rjc), and (d) normalized normal force acting on the joints (FN nj) for the specimen with β = 15° and k = 0.8 at different apertures (a = 0.2 mm, 0.14 mm, 0.1 mm, and 0.0 mm).

With the smooth joint contact model, not only can the inherent roughness induced by the particle size in the conventional joint models of PFC2D be avoided, but also joint opening or closing and joint strength mobilization can be investigated. The data including joint average aperture, contact forces acting on the joints, number and ratio of the closed joint elements, etc. can be monitored to research the mechanism of the multi-peak deformation behavior instead of only observing the failure patterns. It proves to be a promising approach for the mechanism study on rock mass behaviors. Compared with the physical model tests, a series of numerical tests have been conducted on the samples with different joint arrangements, including the joint inclination angle and continuity factor. Four types of specimen behaviors are obtained and it is found that for the specimens with medium to high inclination angles (β = 45°–90°), the joints are prone to open and their strengths are difficult to be mobilized, so strain softening always occur after the peak strength point. If the joints have lower inclination angles (β = 0°–30°) and higher continuity factor (k = 0.6 to 0.8), the strengths of the joints are easy to be mobilized which are adequate to compensate the strength reduction of the damaged matrix and remain the total strength of the jointed specimen, hence a yielding platform or even a strain hardening after the yielding platform may happen before the strain softening. Since the numerical modeling has revealed that the strength mobilization of the joint system plays a very important role in the multi-peak deformation behavior of the specimens with non-persistent joints under uniaxial compression. Furthermore, the influence of joint aperture on the joint response is considered in the numerical tests. It is shown that: (1) For the samples with vertical joints (β = 90°), joint aperture has no influence on the ductile behavior, because it is very difficult for the joint strength to be mobilized in this case; (2) For the other cases (β b 90°), with an increasing joint continuity factor, joint aperture has more significant influence on the ductile behaviors. If the rock mass has high joint continuity factor (k = 0.8) and low joint inclination angles (β = 0°–30°), joint strength is more and more easy to be mobilized with the decreasing joint aperture. The findings in this research may help us understand the deformability of rock masses with different joint orientations, persistence, and apertures, which are important for stability analysis of rock

structures. Mobilization of joint strength is found to play an important role in the deformation behavior of jointed rock mass. This study also supplies an efficient modeling approach to analyze the multi-peak deformation behavior of rock mass. Further studies on influences of other geometrical parameters of joints, such as set number, spacing, roughness of rock joints, as well as the confining pressure, on the multi-peak deformation behaviors of rock mass, will be considered in the future work. Based on more detailed researches, new constitutive equations or yielding/failure criterion may be proposed and may be useful directly for the design of engineering geologists. Acknowledgments This work is financially supported by National Natural Science Foundation of China (no. 41402286 and no. 11572344) and the Strategic Priority Research Program of the Chinese Academy of Sciences, grant no. (XDB10030301). Dr. H. Qian from Institute of Crustal Dynamics, China Earthquake Administration is appreciated here for his support in PFC modeling. References Alejano, L., Ferrero, A.M., Ramírez-Oyanguren, P., Fernández, M.Á., 2011. Comparison of limit-equilibrium, numerical and physical models of wall slope stability. Int. J. Rock Mech. Min. Sci. 48 (1), 16–26. An, X., Ning, Y., Ma, G., He, L., 2014. Modeling progressive failures in rock slopes with nonpersistent joints using the numerical manifold method. Int. J. Numer. Anal. Methods Geomech. 38 (7), 679–701. Bahaaddini, M., Sharrock, G., Hebblewhite, B.K., 2013a. Numerical direct shear tests to model the shear behaviour of rock joints. Comput. Geotech. 51, 101–115. Bahaaddini, M., Sharrock, G., Hebblewhite, B.K., 2013b. Numerical investigation of the effect of joint geometrical parameters on the mechanical properties of a non-persistent jointed rock mass under uniaxial compression. Comput. Geotech. 49, 206–225. Bahaaddini, M., Hagan, P.C., Mitra, R., Hebblewhite, B.K., 2014. Parametric study of smooth joint parameters on the shear behaviour of rock joints. Rock Mech. Rock. Eng. 48 (3), 923–940. Barton, N., 2002. Some new Q-value correlations to assist in site characterisation and tunnel design. Int. J. Rock Mech. Min. Sci. 39 (2), 185–216. Bieniawski, Z.T., 1978. Determining rock mass deformability: experience from case histories. Int. J. Rock Mech. Min. Sci. Geomech. Abstr. 15 (5), 237–247. Bobet, A., Einstein, H.H., 1998. Fracture coalescence in rock-type materials under uniaxial and biaxial compression. Int. J. Rock Mech. Min. Sci. 35 (7), 863–888.

C. Cheng et al. / Engineering Geology 213 (2016) 25–45 Brown, E.T., 1970. Strength of models of rock with intermittent joints. J. Soil Mech. Found. Div. 96, 1935–1949. Cai, M., Kaiser, P.K., Uno, H., Tasaka, Y., Minami, M., 2004. Estimation of rock mass deformation modulus and strength of jointed hard rock masses using the GSI system: Int. Int. J. Rock Mech. Min. Sci. 41, 3–19. Chen, X., Liao, Z., Li, D., 2011. Experimental study of effects of joint inclination angle and connectivity rate on strength and deformation properties of rock masses under uniaxial compression. Chin. J. Rock Mech. Eng. 30 (4), 781–789 (in Chinese). Chen, X., Liao, Z., Peng, X., 2012. Deformability characteristics of jointed rock masses under uniaxial compression. Int. J. Min. Sci. Technol. 22 (2), 213–221. Chen, X., Liao, Z.-h., Peng, X., 2013. Cracking process of rock mass models under uniaxial compression. J. Cent. South Univ. 20, 1661–1678. Chiu, C.-C., Wang, T.-T., Weng, M.-C., Huang, T.-H., 2013. Modeling the anisotropic behavior of jointed rock mass using a modified smooth-joint model. Int. J. Rock Mech. Min. Sci. 62, 14–22. Cording, E., Jamil, M., 1997. Slide geometries on rock slopes and walls. Proceedings of the Fourth South American Congress on Rock Mechanics, Santiago, pp. 199–221. Esmaieli, K., Hadjigeorgiou, J., Grenon, M., 2010. Estimating geometrical and mechanical REV based on synthetic rock mass models at Brunswick mine. Int. J. Rock Mech. Min. Sci. 47 (6), 915–926. Fan, X., Kulatilake, P.H.S.W., Chen, X., 2015. Mechanical behavior of rock-like jointed blocks with multi-non-persistent joints under uniaxial loading: a particle mechanics approach. Eng. Geol. 190, 17–32 (0). Gehle, C., Kutter, H.K., 2003. Breakage and shear behaviour of intermittent rock joints. Int. J. Rock Mech. Min. Sci. 40 (5), 687–700. Ghazvinian, A., Sarfarazi, V., Schubert, W., Blumel, M., 2012. A study of the failure mechanism of planar non-persistent open joints using PFC2D. Rock Mech. Rock. Eng. 45 (5), 677–693. Gong, W.L., Wang, J., Gong, Y.X., Guo, P.Y., 2013. Thermography analysis of a roadway excavation experiment in 60° inclined stratified rocks. Int. J. Rock Mech. Min. Sci. 60 (48), 134–147. Gong, W., Peng, Y., Sun, X., He, M., Zhao, S., Chen, H., Xie, T., 2015. Enhancement of lowcontrast thermograms for detecting the stressed tunnel in horizontally stratified rocks. Int. J. Rock Mech. Min. Sci. 74, 69–80 (0). Hadjigeorgiou, J., Esmaieli, K., Grenon, M., 2009. Stability analysis of vertical excavations in hard rock by integrating a fracture system into a PFC model. Tunn.Undergr. Sp. Tech. 24 (3), 296–308. Halakatevakis, N., Sofianos, A.I., 2010. Strength of a blocky rock mass based on an extended plane of weakness theory. Int. J. Rock Mech. Min. Sci. 47 (4), 568–582. Han, Y., Damjanac, B., Nagel, N., 2012. A microscopic numerical system for modeling interaction between natural fractures and hydraulic fracturing. Proceedings of the 46th US Rock Mechanics/Geomechanics Symposium, Chicago, IL, ARMA, pp. 12–238. He, M., 2011. Physical modeling of an underground roadway excavation in geologically 45° inclined rock using infrared thermography. Eng. Geol. 121 (3–4), 165–176. He, M., Jia, X., Gong, W., Faramarzi, L., 2010. Physical modeling of an underground roadway excavation in vertically stratified rock using infrared thermography. Int. J. Rock Mech. Min. Sci. 47 (7), 1212–1221. Hoek, E., Carranza-Torres, C., Corkum, B., 2002. Hoek-Brown failure criterion-2002 edition. Proceedings of the 5th North American Rock Mechanics Symp., Toronto, Canada, pp. 267–273. Horii, H., Nemat-Nasser, S., 1985. Compression-induced microcrack growth in brittle solids: axial splitting and shear failure. J. Geophys. Res. 90 (B4), 3105–3125. ISRM, 1978. International society for rock mechanics commission on standardization of laboratory and field tests: suggested methods for the quantitative description of discontinuities in rock masses. Int. J. Rock Mech. Min. Sci. Geomech. Abstr. 15 (6), 319–368. Jennings, J., 1970. A mathematical theory for the calculation of the stability of slopes in open cast mines. Proceedings of the Symposium on Theoretical Background to the Planning of Open Pit Mines, Johannesburg, pp. 87–102.

45

Jiang, Y., Li, B., Yamashita, Y., 2009. Simulation of cracking near a large underground cavern in a discontinuous rock mass using the expanded distinct element method. Int. J. Rock Mech. Min. Sci. 46 (1), 97–106. Karami, A., Stead, D., 2008. Asperity degradation and damage in the direct shear test: a hybrid FEM/DEM approach. Rock Mech. Rock. Eng. 41 (2), 229–266. Kulatilake, P.H.S.W., Ucpirti, H., Wang, S., Radberg, G., Stephansson, O., 1992. Use of the distinct element method to perform stress analysis in rock with non-persistent joints and to study the effect of joint geometry parameters on the strength and deformability of rock masses. Rock Mech. Rock. Eng. 25 (4), 253–274. Ladanyi, B., Archambault, G., 1972. Evaluation de la resistance au cisaillement d'un massif rocheux fragmente. Proc. 24th Int. Geological Congress, Montreal, pp. 249–260. Lajtai, E.Z., 1969. Shear strength of weakness planes in rock. International Int. Int. J. Rock Mech. Min. Sci. Geomech. Abstr. 6 (5), 499–515. Lambert, C., Coll, C., 2014. Discrete modeling of rock joints with a smooth-joint contact model. J. Rock Mech. Geotech. Eng. 6 (1), 1–12. Li, B., Jiang, Y., Mizokami, T., Ikusada, K., Mitani, Y., 2014. Anisotropic shear behavior of closely jointed rock masses. Int. J. Rock Mech. Min. Sci. 71, 258–271 (0). Mas Ivars, D., Pierce, M.E., Darcel, C., Reyes-Montes, J., Potyondy, D.O., Paul Young, R., Cundall, P.A., 2011. The synthetic rock mass approach for jointed rock mass modelling. Int. J. Rock Mech. Min. Sci. 48 (2), 219–244. Morris, J.P., Rubin, M.B., Block, G.I., Bonner, M.P., 2006. Simulations of fracture and fragmentation of geologic materials using combined FEM/DEM analysis. Int. J. Impact Eng. 33 (1−12), 463–473. Park, E.-S., Martin, C.D., Christiansson, R., 2004. Simulation of the mechanical behavior of discontinuous rock masses using a bonded-particle model. Proceedings of the 6th North America Rock Mechanics Symposium (NARMS), Houston, Texas, p. 8. Potyondy, D.O., Cundall, P.A., 2004. Abonded-particle model for rock. Int. J. Rock Mech. Min. Sci. 41, 1329–1364. Prudencio, M., Van Sint Jan, M., 2007. Strength and failure modes of rock mass models with non-persistent joints. Int. J. Rock Mech. Min. Sci. 44 (6), 890–902. Scholtès, L., Donzé, F.-V., 2012. Modelling progressive failure in fractured rock masses using a 3D discrete element method. Int. J. Rock Mech. Min. Sci. 52, 18–30. Serafim, J.L., Pereira, J.P., 1983. Considerations of the geomechanics classification. In: Bieniawski (Ed.)Proceedings of International Symposium Engineering Geology and Underground Construction vol. 1, pp. 1133–1142. Singh, M., Seshagiri Rao, K., 2005. Empirical methods to estimate the strength of jointed rock masses. Eng. Geol. 77 (1–2), 127–137. Singh, M., Singh, B., 2008. High lateral strain ratio in jointed rock masses. Eng. Geol. 98 (3– 4), 75–85. Singh, M., Choudhari, J., Samadhiya, N.K., 2009. An experimental study to predict ground response curve for tunnels. In: Kulatilake, P.H.S.W. (Ed.), Procedings of International Conference on Rock Joints and Jointed Rock Masses. Tucson, Arizona, USA, p. 8. Tiwari, R.P., Rao, K.S., 2006. Post failure behaviour of a rock mass under the influence of triaxial and true triaxial confinement. Eng. Geol. 84, 112–129. Wasantha, P.L.P., Ranjith, P.G., Xu, T., Zhao, J., Yan, Y.L., 2014. A new parameter to describe the persistency of non-persistent joints. Eng. Geol. 181, 71–77. Wong, L.N.Y., Einstein, H.H., 2009. Crack coalescence in molded gypsum and Carrara marble: part 1. Macroscopic observations and interpretation. Rock Mech. Rock. Eng. 42, 475–511. Wong, R., Chau, K., Tang, C., Lin, P., 2001. Analysis of crack coalescence in rock-like materials containing three flaws—part I: experimental approach. Int. J. Rock Mech. Min. Sci. 38 (7), 909–924. Yang, X.X., Kulatilake, P.H.S.W., Chen, X., Jing, H.W., Yang, S.Q., 2016. Particle flow modeling of rock blocks with Nonpersistent open joints under uniaxial compression. Int. J. Geomech. http://dx.doi.org/10.1061/(ASCE)GM.1943-5622.0000649, 04016020. Zhang, H.Q., Zhao, Z.Y., Tang, C.A., Song, L., 2006. Numerical study of shear behavior of intermittent rock joints with different geometrical parameters. Int. J. Rock Mech. Min. Sci. 43 (5), 802–816.