#### Abstract

The anchoring mechanism of a bolted joint subjected to a shear load was investigated using a bilinear constitutive model via the inner-embedded FISH language of particle flow code based on the discrete element method. The influences of the anchoring system on the macro-/micromechanical response were studied by varying the inclination angle of the bolt. The results indicate a clear relationship between the mechanical response of a bolted rock joint and the mechanical properties of the anchoring angle. By optimizing the anchorage angle, the peak strength can be increased by nearly 50% relative to that at an anchorage angle of 90°. The optimal anchorage angle ranges from 45° to 75°. The damage mechanism at the optimal anchorage angle joint is revealed from a macroscopic mechanical perspective. The concentration of the contact force between disks will appear in the joint and around the bolt, resulting in crack initiation. These cracks are mainly tensile cracks, which are consistent with the formation mechanism for compression-induced tensile cracks. Therefore, the macroscopic peak shear stress in the joint and the microscopic damage to the anchoring system should be considered when determining the optimal anchoring angle to reinforce a jointed rock mass.

#### 1. Introduction

Joints are extremely common in rock masses, and the stability of a jointed rock mass mostly depends on the mechanical properties of the joint. Since the 1960s, the shear behavior of rock joints has been the focus of many studies, and many achievements have been made [1–9]. In recent years, rock bolts have been widely used in engineering due to their convenience, affordability, and reliability. However, the anchoring mechanism for a jointed rock mass is still not clearly understood. As is known, the anchoring mechanism for a bolted rock mass is highly complex due to the anisotropy of the rock mass and the complexity of the joint pattern. Many scholars have conducted experimental and theoretical studies to reveal the anchoring mechanism in a bolted rock joint.

Sten Bjurström was the first to report systematic research on fully grouted rock bolts. His shear tests were conducted on full cement/grout bonded rock bolts embedded in blocks of granite. According to Bjurström [10], inclining the bolt resulted in stiffening of the shearing surface by increasing the shear strength at a small displacement. Ge and Liu [11] discussed the influence of the bolt on the joint and the mechanism behind the dowel effect via laboratory simulation tests and theoretical analyses in bolted joints and proposed a modified shear stress formulation for a bolted joint. Spang and Egger [12] conducted a series of shear tests in a bolted joint to study the deformation characteristics of a bolt and found two critical points in a bolt: one at the bolt-joint intersection and the other at the hinge point. Egger and Zabuski [13] found that the anchoring mechanism of a bolt in a reinforced joint contributes additional shear resistance to prevent the shear failure of a joint. Bolt deformation typically occurs near the joint. Ferrero [14] proposed a shear strength model for a reinforced joint that considered both the dowel effect and the incremental increase in axial force due to the bar deformation. Pellet and Egger [15] related theoretical and experimental analyses of the rock bolt shear strength and found that bolts installed perpendicular to a joint plane allowed the greatest displacement along a joint prior to failure but that the displacement at failure decreased rapidly as the angle between the bolt and joint plane decreased. They also found that harder rock led to bolt failures at smaller displacements. Jalalifar and Aziz [16] performed double shear tests and numerical simulations on five types of bolts to study the bolt force and failure mechanism. Wang et al. [17] used ANSYS to calculate the stress intensity factors on the tip of cracks for different anchoring conditions and obtained a relationship between the anchoring spacing, anchoring angle, and stress intensity factor using FLAC3D software to simulate the coalescence modes for main-control cracks at different anchor spacings. Li et al. [18] analyzed a numerical model of a fully grouted rock bolt installed in concrete with FLAC3D and concluded that the bolt’s resistance to shear was influenced by the rock strength, inclination angle, and diameter of the rock bolt. In recent years, with the development of high power computers, the numerical software PFC based on the discrete element method has become increasingly mature in studying rock and soil failure behavior [19]. Many achievements have been made in determining the mechanical properties of rock joints [20–24]. Rock joints effectively compensate for the lack of repeatable rock mass. Some studies have examined the roughness of jointed rock mass by using a DEM model. For example, Xia et al. [25] carried out a numerical simulation of shear experiments on nonbolted rock using the discrete element software PFC and analyzed the direct shear characteristics of rough joints. Using the experimental results of rock joints, Zhou et al. [26] analyzed the fracture evolution characteristics of rock mass during shearing via discrete element PFC software. Cao et al. [27] used PFC to generate shear models with different roughness values and discussed the influence of roughness on the surface morphology of nonanchor joints from a micro point of view.

To date, only a few studies have considered the constitutive model of the bolt, particularly under bolted shear conditions. Most studies have focused on the shear performance of a bolted joint and its deformation characteristics while neglecting the influence of the anchorage angle on the mechanical properties of bolted rock joints. These studies can reveal the anchoring mechanism in a bolted joint to a certain extent. The micromechanism is the fundamental basis of the macroscopic properties. Therefore, it is of practical significance to research the macro-/microanchoring mechanism for a bolted rock joint at different anchorage angles to provide further references to guide designs of anchors in a bolted rock joint at varying anchorage angles.

Based on the previous results, in this study, we modified the bolt constitutive model and analyzed the macro-/micromechanical response of a bolted jointed rock mass using the discrete element software PFC at different anchorage angles. The disk rotation and microcracking in the thin layer on the face of the anchor joint can be analyzed at different inclinations based on the stress state of the anchor bolts and the disk contact force during the shearing process.

#### 2. PFC Software and BPM Model

PFC is a discrete element commercial software developed by the Itasca Consulting Group and has been already widely used in the rock mechanics field [19]. PFC represents a rock mass as an assemblage of circular disks with a finite thickness and connected via cohesive and frictional bonds.

A basic linear contact model describes the elastic relationship between the relative displacements and forces of disks at the point contact, as shown in Figure 1.

**(a)**

**(b)**

**(c)**

This model involves the contact normal force component, , contact overlap, , shear force increment, , and shear-displacement increment, , and is given bywhere is the normal stiffness (force/displacement) at the contact. The value of is determined by the current contact-stiffness model.

Note that the normal stiffness, , is a secant modulus in that it relates total displacement and force. The shear stiffness , on the other hand, is a tangent modulus in that it relates incremental displacement and force. An uppercase will be used to denote a secant modulus, and a lowercase will be used to denote a tangent modulus. The computation of the normal contact force from the geometry alone makes the code less prone to numerical drift and able to handle arbitrary placement of balls and changes in ball radii after a simulation has begun.

The shear contact force is computed in an incremental fashion. When the contact is formed, the total shear contact force is initialized to zero. Each subsequent relative shear-displacement increment results in an increment of elastic shear force that is added to the current value. The motion of the contact must be considered during this procedure.

The contact velocity can be resolved into normal and shear components with respect to the contact plane. Denoting these components by and for the normal and shear components, respectively, the shear component of the contact velocity can be written as

The shear component of the contact displacement-increment vector, occurring over a time step of , is calculated byand is used to calculate the shear elastic force-increment vector

The frictional resistance of the contact is given bywhere is the friction coefficient between the disks.

To simulate a relatively brittle rock-like material, it is necessary to cement these disks with a bonded model. This study uses the parallel bond model, which resists not only the contact forces but also the moments between the disks at a cemented contact (Figure 1). The function mechanism of the parallel bond model is described bywhere and are the force components about the center of the cemented-contact zone, and are the moments about the center of the cemented-contact zone, and are the normal and shear bond stiffness per unit area, respectively, and are the rotation angle components, and , , and are the area, polar moment of inertia, and moment of inertia of the bond contact cross section, respectively. The strength of the cemented contact is then given bywhere is the radius of the bonded zone between the disks, is the length of the bonded zone between the disks (Figure 1), and and are the tensile and shear strength of the bond contact, respectively.

Tensile cracks occur when the applied tensile stress exceeds the specified tensile strength of the parallel bond, . Shear cracks occur when the applied shear stress exceeds the specified shear strength of the parallel bond, , either by rotation or by the shearing of the disks. The tensile strength at the contact immediately drops to zero once the crack occurs, and the shear strength reduces to the residual friction strength (see (5)) [19, 22], as illustrated in Figure 2.

#### 3. Direct Shear Simulation of the Bolted Joint Using PFC

##### 3.1. Microparameters of the Rock Sample

The intact synthetic model material selected in this study is represented by compacted disks cemented with the parallel bond model. Therefore, the micromechanical parameters consist of two categories: one for disks and the other for parallel bonds. The PFC code allows one to simulate the macromechanical behavior of the selected synthetic model material by selecting appropriate values for the micromechanical parameters, such as the disk size distribution and packing, disk and parallel bond stiffness, disk friction coefficient, and bond strengths [19, 22]. To select the appropriate micromechanical parameter values for the synthetic intact material, one needs to basically go through a trial and error procedure by iteratively varying the micromechanical parameter values to match the required macromechanical behaviors of the selected synthetic material. This procedure is known as the calibration of the intact synthetic material based on the micromechanical parameters. Herein, the required macromechanical behaviors are the strength response, represented by the uniaxial compressive strength (UCS), and the deformability response, represented by Young’s modulus and Poisson’s ratio of a cylindrical specimen with a 50 mm diameter and a height/diameter ratio of 2 : 1 under unconfined compression.

In this study, a special calibration sequence was followed to rationalize the micromechanical parameter calibration and to minimize the number of iterations. First, the disk and parallel bond moduli and the ratios of normal to shear stiffness were set equal between the disks and parallel bonds to reduce the number of independent parameters [22]. Then, Young’s modulus was calibrated by setting the material strengths to a large value and varying disk Young’s modulus and parallel bond Young’s modulus to match Young’s modulus between the numerical and laboratory specimens. Next, by changing the ratio of disk normal stiffness/disk shear stiffness and parallel bond normal stiffness/parallel bond shear stiffness, Poisson’s ratio of the numerically simulated intact synthetic cylindrical specimen was matched to that of the laboratory specimen.

After calibrating the aforementioned deformation micromechanical parameters, the peak strength between the numerical and laboratory specimens was matched by gradually reducing the normal and shear bond strengths of the parallel bonds. During this procedure, it is important to fix the ratio of normal to shear bond strength , because these parameters affect the failure mode of the specimen. Therefore, a series of numerical simulations were conducted to match the failure mode of the cylindrical specimen between the numerical and laboratory specimens by varying , while keeping the other parameters unchanged. The obtained failure modes with varying for the cylindrical specimens used in the numerical simulations are displayed in Figure 3.

**(a) Sample from PFC**

**(b) Samples of rock-like materials**

As seen in Figure 3, as the specimen is compressed under uniaxial stress, numerous cracks (tensile and shear cracks) are produced through the breakage of parallel bonds. The failure mode of the numerically simulated synthetic cylindrical specimen with agrees well with the results from the laboratory specimen. Thus, setting the ratio of normal to shear bond strength equal to one seems to increase the confidence of the numerical model in simulating the appropriate failure behavior of the synthetic material.

The numerically obtained stress-strain curves were compared with the stress-strain curve obtained from the laboratory sample as shown in Figure 4. In this comparison, in addition to comparing the prepeak behavior, it is necessary to match the macro UCS and Young’s modulus values obtained for different numerical simulations of the synthetic material with those obtained for the laboratory specimen. Therefore, the disk contact modulus values need to be modified simultaneously with the value. Figure 4 shows that the prepeak and postpeak trend of the laboratory stress-strain curve is approximately parallel to that of all the stress-strain curves of numerical specimens. Therefore, the prepeak and postpeak behavior of numerical specimens did not provide any guidance for the selection of disk friction value. However, the coefficient of friction obtained from the laboratory triaxial tests on the synthetic intact rock was approximately 0.50. Therefore, the friction coefficient value of 0.50 was selected to represent the disk friction. Table 1 lists the determined micromechanical parameter values to simulate the intact synthetic rock with the parallel bond model. The failure mode, the stress-strain curves, and the comparison of the macromechanical parameters indicate that the results are in good agreement with the laboratory test results.

##### 3.2. Bilinear Constitutive Model of the Bolt

The bolt is represented by a cluster of balls with parallel bonds that meet the requirements of the force-displacement law. Parallel bonds can transmit both forces and moments between disks, which is consistent with the stress state in a bolt. The stress state in a bolt is highly complex when it is reinforced in a jointed rock mass. The microparameters of the bolt are determined via the bolt pull-out test, as shown in Figure 5(a). The microparameters of the bolt are shown in Table 2. Figure 5(b) describes the bolt constitutive model, which is not modified. One shortcoming of this model is that it can represent the characteristics of a metallic bar only during elastic deformation; the model cannot represent the large plastic deformations during yielding. Therefore, this constitutive model has been modified for use in the internal FISH language in PFC to make the modified bolt constitutive model meet the stiffness and deformation requirements. As shown in Figure 5(c), the yield strength of the bolt is set such that it remains in the linear elastic stage prior to the axial stress reaching the yielding strength. However, when the axial stress exceeds the yield strength, the bonds will not break, meaning that the bolt maintains a certain capacity. If the magnitude of the axial stress in the bolt exceeds the maximum normal stress, the bonds are broken, and the bolt will lose its bearing capacity.

**(a)**

**(b)**

**(c)**

##### 3.3. Direct Shear Test Model of the Bolted Joint

The direct shear test model is built using PFC. The model is composed of a top and bottom rock sample, grout, joint, and a bolt with a size of 100 mm × 50 mm, with a total of 2,926 disks. The microparameters of the rock sample used are shown in Table 1. The grout microparameters are based on the rock parameters with the parallel bond strength changed to 40 MPa. The joint is created using its own command in PFC, and the disks in the joint are unbounded. Therefore, the parallel bond stiffness and the strength of the disks in the joint are set to zero. A shear box is generated by the wall element. Among these walls, the sixth wall is the loading wall, with the first wall, fifth wall, and the sixth wall together forming the down active shear box, which can move in the horizontal direction. The upper shear box is fixed in the horizontal direction, and only the second wall can move up or down, as it is the servo wall, which keeps the normal stress constant during shearing. As shown in Figure 6, the rock mass is represented by light gray disks, and the bolt is composed of a cluster of light gray disks, which are tangential to each other. Around the bolt, the dark gray disks represent the grout. Along the horizontal direction in the middle of the model, the joint is composed of black disks.

**(a) Numerical model of a bolted joint at 45°**

**(b) Numerical model of a bolted joint at 75°**

##### 3.4. Direct Shear Test of the Anchoring Model

During shearing, the second wall is controlled by the servo mechanism to apply a constant normal stress, with the vertical displacement monitored as the normal displacement. The fifth wall is subjected to a constant velocity to apply horizontal shear loading, with the horizontal displacement monitored as the shear displacement. The shear stress is calculated by dividing the average force on the fifth and sixth walls. The test is ended when the shear displacement reaches the present value. Microcracking due to bond breakage, shear stress, normal displacement, stress in the bolt, and displacement variations in the bolt disks are monitored during the shear test.

#### 4. Analysis of the Results

##### 4.1. Analysis of the Shear Performance of the Joint with Different Anchorage Angles

Direct shear simulations of the bolted joint were conducted using the bolt at different anchorage angles under a constant normal stress to examine the macro-/micromechanical responses of the bolted joint during shearing. Figure 7 shows the evolution of the shear stress with the shear displacement at different anchorage angles. The shear simulation at a normal stress of 5 MPa is achieved after a shear displacement of 2.5 mm is reached. The peak shear strength for different anchorage angles is shown in Table 3.

Our data show that () the peak shear strength increased to a maximum and then decreased with further increase in the anchorage angle and () the peak shear strength is the largest at an anchorage angle of 75°. Compared to the 90° anchorage angle, changing the angle of the anchor can greatly improve the bolt’s antishear abilities. Therefore, it is important to determine the optimum design for a bolted rock joint in practice.

##### 4.2. Analysis of the Bolt Stress in a Compression-Shear Stress State

###### 4.2.1. Stress Analysis of the Anchor Joint Plane in a Compression-Shear Stress State

Under compressive shear stress, the rock bolts experience not only shear deformation through the thickness of the weak plane but also notable shear deformation over a considerable section (the effective length). To simplify the analysis, we set the effective length of the rock bolt to be and the effective length of the bolt to be , as shown in Figure 8.

The average shear stress distribution on the plane of the section is related to the contact conditions of the rod body and hole wall. Three different distribution patterns for the simplified calculation are shown in Figure 9(a). The distribution pattern is a curve but can be measured in three basic forms. Similar to the section, the axial stress distribution can be used, as shown in Figure 9(b).

**(a) Shear force**

**(b) Axial force**

From Figures 8 and 9 and from the joint surface equilibrium conditions, the stress tensor of the rock mass can be obtained from the normal and tangential components of the stress tensor along the joint plane, as shown in

In the previous formula, and are the weak levels of tangential and normal displacement, respectively; and are the bolt axial and transverse deformations, respectively; is the angle between the bolt and joint surface; and and are Young’s modulus and shear modulus for the material, respectively. The shear coefficient is associated with the rod body section ( for a circular solid rod). and are the inner pole sections’ average shear stress and axial stress distribution of the shape coefficient, respectively, which are 1, 3/4, and 1/2. and are the axial and shear forces on the bolt, respectively. and are the normal and tangential stiffness coefficients for the weak level. is the anchor section reinforcement ratio. Therefore,where

Therefore, the relationship between the stress and the stress component of the element can be determined as shown in

###### 4.2.2. Analysis of the Bolt Stress for a Particular Anchorage Angle

A bolt installed at a particular anchorage angle will have a better anchoring effect than a bolt installed perpendicular to the joint plane. Figures 10 and 11 show illustrations of the stresses in a bolt at anchorage angles of 75° and 90°, respectively, during shearing.

**(a) Axial stress diagram**

**(b) Bending moment diagram**

**(c) Shear stress diagram**

**(a) Axial stress diagram**

**(b) Bending moment diagram**

**(c) Shear stress diagram**

The diagram showing the two different anchorage angle bolts indicates that the anchor force and bending moment are relatively large in the vertically installed bolt, leading to bolt axial and shear stress that are larger than those in the tilted bolt, whereas, for the tilted bolt, the axial and shear stress are small and do not reach the yield strength. Therefore, the tilted anchor bolts not only improve the shear capacity of the joints but also maintain an appropriate stress state in the anchor. The rock bolt axial stress diagram indicates that the axial stress increases as the shear displacement increases for different anchorage angles. The axial stress is the highest when the anchor is vertically installed on the joint surface. The results of these two cases reveal that the optimal anchoring angle is approximately 75°. This result corresponds to the previous conclusions that the optimal anchoring angle should be between 45° and 75°.

##### 4.3. Evolution of the Contact Force and Microscopic Cracking of Disks

The anchorage angle of 75° was chosen for a more detailed analysis on the micromechanical response and damage mechanism during shearing. Figure 12 shows the evolution of the shear stress, number of cracks, and normal displacement for a normal stress of 5 MPa.

To accurately determine the evolution features, six different monitoring points (a to f) were set on the shear stress curve, which correspond to shear displacements of 0, 0.198, 0.217, 0.513, 1.500, and 2.500 mm, respectively. The variations in the contact force and cracking can be observed from the six monitoring points. As shown in Figure 12, “a” indicates the initial point of the shear displacement and “c” indicates the shear displacement at which the shear stress reaches the peak shear strength. As shown in Figure 12, no shear crack is generated during the initial stage of the shear test. The number of shear cracks increases gradually as the shear displacement increases. When the shear stress reaches the peak shear strength, the number of shear cracks increases significantly, and the shear crack growth rate tends to stabilize as the shear stress increases beyond the peak shear strength. The trends in tension cracking and normal stress are highly similar and increase monotonically.

During the shearing of a bolted joint, the loaded side of the bolt is in a state of compression, whereas the other side of the bolt is in a state of tension due to the bending deflection of the bolt. The adjacent grout is in a state of compression, and the rock mass is in a state of tension because of the bonds between the bolt and grout, as shown in Figure 13. Bending deflection will occur in the bolt as the shear displacement increases. Simultaneously, cracks will occur continuously due to extrusion and slippage between disks. In the compression zone, there are a large number of cracks, which are mainly tension cracks generated by the extrusion between disks. This result agrees with the formation mechanism for compression-induced tension cracks shown in Figure 13.

**(a)**

**(b)**

Figure 14 shows the illustration of the contact force and crack propagation between disks during shearing. The contact force is initially concentrated on the joint and the loading end, and cracks are generated in the joint, whereas there are only a few cracks around the bolt. As the shear displacement increases, the disk’s contact force around the bolt becomes gradually larger. Cracks occur on the joint, and those around the bolt extend further. The cracks around the bolt propagate most rapidly. The cracks are initially concentrated at the intersection of the bolt and joint and then propagate to the ends of the bolt. Finally, the majority of the cracks are distributed at the joint and around the bolt after a larger displacement occurs. Finally, more cracks occur in the compression zone around the bolt than in the tension zone. In addition, the distribution of cracks corresponds to the contact force in the stress concentration area. The detailed maximum contact force and crack number in the specimen at the different monitoring points during shearing are provided in Tables 4 and 5, respectively.

**(a) Monitoring point a**

**(b) Monitoring point b**

**(c) Monitoring point c**

**(d) Monitoring point d**

**(e) Monitoring point e**

**(f) Monitoring point f**

##### 4.4. Evolution of the Disk Rotation

During the shearing slip process, the behavior is determined by the rotation of the disks. Relative disk dislocations will occur due to the rotation of the disks. If the rotation angle is sufficiently large, it will lead to shear stresses and the breaking of disk bonds. Meanwhile, cracks occurring in the rock mass and shear cracking are caused directly by disk rotation. Figure 15 shows a graph of the evolution of the disk rotation in terms of radians for rock disks at different monitoring points.

**(a) Monitoring point a**

**(b) Monitoring point b**

**(c) Monitoring point c**

**(d) Monitoring point d**

**(e) Monitoring point e**

**(f) Monitoring point f**

Figure 15 shows that when the shear displacement is small at the initial stages, the disk rotation angle is relatively small, and disks with large rotational angles are distributed around the bolt. The disk rotation angles gradually increase with increasing shear displacement. The disks containing larger rotation angles are mainly distributed in the two loading ends and around the bolt, and the bolt and the disks around the bolt contact move relatively violently during shearing, leading to a larger rotation angle. Further increases in the shear displacement and larger rotation disks also appeared on the joint surfaces, which are caused by the joint surface microconvex body cutting into the free disks. The distribution of the large disk rotation angles is highly consistent with the crack distribution, demonstrating that the rotation of the disks leads to the generation of shear cracks.

Figure 16 shows the variations in the proportion of rotations exceeding 0.1 radians. As the shear displacement increases, the number of clockwise and counterclockwise disk rotations greater than 0.1 radians increases continuously. Additionally, the number of clockwise rotations is greater than the number of counterclockwise rotations, mainly due to the shear loading directions.

#### 5. Conclusions

In this paper, the anchoring mechanism for a tilted anchor bolt subjected to a shear load has been studied using a bilinear constitutive model with the inner-embedded FISH language of PFC based on DEM. The anchoring mechanism of the bolt at different anchorage angles has been discussed at the macro-/microlevels, and the following conclusions have been drawn.

() The bolt can improve the shear resistance of the joint, mainly due to cohesion. Tilted anchor bolts not only improve the shear strength of the joint surface but also result in a lower stress in the anchor bolt. The peak shear strength can be 50% greater than that of a vertically installed anchor bolt.

() A tilted bolt allows for fully mobilized anchor performance. This configuration not only restricts increment in the normal displacement but also indirectly improves the stress state of the anchorage system and increases the shear capacity of the anchor joint. The comprehensive shear displacement and shear stress data along with the bolt stress state analysis reveal that the optimal anchorage angle should be between 45° and 75°.

() During the shearing process, when the rock, bolt, and grout are commonly constrained, the stress concentration phenomenon around the bolt and joint increases the shear displacement and contact force. After reaching the residual phase, the contact force exhibits slight fluctuations due to the sliding of the joint surface.

() Based on the disk flow software PFC, crack initiation for an anchorage angle of 75° will start from the contact force concentration due to the interactions between the rock, grout, and bolt. Finally, the cracks are mainly distributed in the joint and around the bolt for larger displacements. The cracks are mainly tensile cracks, which correspond to the formation mechanism for a compression-induced tensile crack. The shear cracking was mainly caused by the rotation of the disks during the shearing process.

#### Disclosure

Gang Wang and Yongzheng Zhang should be regarded as co-first authors.

#### Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

#### Acknowledgments

This study was supported by the National Natural Science Foundation of China (nos. 51479108, 51379117, and 51509261), China Postdoctoral Science Foundation (2015M580599), the Taishan Scholar Talent Team Support Plan for Advantaged & Unique Discipline Areas, the Postdoctoral Science Foundation of Shandong Province (no. 201402014), and the Shandong Province Outstanding Young Scientist Research Award Fund (BS2014NJ006).