Mechanical Behaviour of Anchored Rock Containing Weak Interlayer under Uniaxial Compression: Laboratory Test and Coupled DEM–FEM Simulation

Abstract

The reason for instability in a rock mass with a weak interlayer is not only the sliding failure of the rock interlayer structural plane but also the tip crack propagation caused by the stress concentration at the tip of the interlayer. In this study, a uniaxial compression test of an anchored rock mass with a weak interlayer was carried out to determine the influence of the anchorage on the failure pattern and the strength of samples with different interlayer dip angles. In addition, the coupled DEM–FEM numerical simulation method was used to study the internal stress evolution of the sample and the stress distribution of the anchor under the anchorage effect. The results showed that the anchorage effect on reinforcement and strength enhancement was greatest for the sample with an interlayer dip angle of 30°. Under the anchorage effect, crack initiation was limited and there was more shear failure in the samples. The reinforcement range of the anchorage effect for anchors with restrained ends was larger than for anchors with free ends. When the rock–anchor interface was unbonded, the effect of the free-ends anchor reflected the residual friction, but the restrained-ends anchor still worked by limiting the lateral expansion of the rocks. The stress values and deformation of the anchors decreased gradually with an increase in dip angle.

1. Introduction

A weak interlayer is a special and widely distributed weak structure plane in a rock mass, whose characteristics are low modulus, low strength and large deformation [1,2]. When subjected to stress interference, nonlinear stress transfer will occur at the interface between the rock and the weak interlayer, which further leads to irregular stress distribution, internal instability, failure and crack propagation in the rock mass [3]. At present, anchoring technology for rock masses is an important part of geotechnical engineering. An anchor has the advantages of low materials consumption and construction safety. For foundation pit and open-pit mine slope problems especially, anchoring technology can make full use of the characteristics of the rock mass itself, reinforcing the rock mass with the help of the anchor in order to achieve rock mass anti-sliding, shear and crack prevention effects [4,5].

A joint is a special unfilled and occlusive weak plane. The mechanical behaviour of an anchored fractured rock mass is currently an important research topic. Studies on the mechanical behaviour of anchored jointed rock masses have been carried out under uniaxial compression and shear, the two most common and simplest stress states [6,7,8,9,10,11,12,13,14,15]. Li et al. [6] and Yang et al. [7] conducted a laboratory test on the effect of a prestressed anchor in a fractured rock mass. The results illustrated the optimal anchorage angle and position, and also showed that the role of the anchorage was to inhibit crack initiation and expansion and increase the integrity of the fractured rock mass. For an anchored rock mass with a perforated weak interlayer, Zong et al. [8] found that the anchoring effect could change the failure mode of the sample from brittle to ductile and that tensile failure was replaced by more shear failure. A similar conclusion can be drawn from the research in [9,10]. Under shear conditions, many scholars focused on the effect of the anchorage on the shear force, displacement and deformation, and achieved good results. Ma et al. [11] proposed a full-length grouting anchorage jointed rock mass model, which can fully verify the relationship between global shear load and shear displacement and can predict the process of global shear evolution. Li et al. [12] conducted an anchorage shear test on a jointed rock mass. In this study, it was found that the residual resistance contribution value of the full-length bonded bolt was highest for an anchorage angle of 55–70 degrees. With an increase in the joint expansion angle, the maximum contribution value of the anchor also increased linearly. In addition, the effect of factors such as different anchor types [13], anchorage modes [14], joint roughness values [15], loading angles and displacement angles [16] were studied in anchorage shear tests. The above studies were based on a perforated weak interlayer or an internally closed jointed rock mass. However, the influence of geological environment, seepage and weakening in the inner rock mass could create a non-interpenetrated weak interlayer. A non-interpenetrated weak interlayer combines the characteristics of the above two types and has a binary medium and an anisotropic structural plane, and should therefore form part of anchored rock mass studies.

In recent years, numerical simulations have been increasingly widespread in geotechnical engineering. There are usually two methods for the numerical simulation of anchoring samples, namely, the finite element method [7,17,18,19] and the discrete element method [20,21,22]. Using the finite element method, the internal monitoring value of the sample can easily be obtained [7,17,18,19], as well as the deformation and stress of the anchor. However, the internal microscopic changes in the rock mass are hard to capture. Most researchers use bonded particles to simulate the anchor and rock mass in the discrete element method [20,21,22]. Although micro-crack initiation and propagation in the sample can be better simulated using DEM, the elastic–plastic behaviour of a continuous medium such as an anchor cannot be reflected. Hence, coupled DEM–FEM combines the advantages of the two methods and provides a solution to these problems [23].

In summary, current research mainly takes a perforated weak interlayer or an internally closed jointed rock mass as the test object and adopts a single method for numerical simulation. In view of this, this study adopted laboratory testing and coupled DEM–FEM numerical simulation to study the strength, deformation, failure and anchor stress in anchored sandstone samples with non-interpenetrated weak interlayers under uniaxial compression, in order to provide a useful reference for anchorage design for rock engineering with interlayers.

2. Experiment and Simulation Preparation

2.1. Experimental Samples Preparation

In the experiment, the main rock mass was white sandstone, and the weak interlayer was made of cement slurry. White sandstone samples were obtained from the Sichuan region, Southwest China. Quartz, feldspar, illite, kaolinite, etc., are the main mineral components of the sandstone in this region.

Due to the influence of the geological environment, structure and other factors, the compositions and strengths of white sandstone in this region are quite different [24,25,26]. The main mineral components of sandstone in Chongqing city are quartz, feldspar, clay, mica, cuttings, etc., and the mechanical parameters of this sandstone were found to be as follows: the uniaxial compression strength, Poisson’s ratio and tensile strength were 54.6 MPa, 0.24 and 5.86 MPa, respectively [24]. Sandstone sourced from Sichuan Basin in China contained 41.7% quartz, 4.7% feldspar, 2.5% calcite, 21.8% illite, 24.7% plagioclase, 3.6% zeolite and 1% other clay minerals, with mechanical properties as follows: a density of 2428 kg/m³, modulus of 18.58 GPa, uniaxial compression strength of 88.52 MPa and Poisson’s ratio of 0.16 [25]. The mechanical properties of sandstone from Zigong, Sichuan, were as follows: a density of 2329 kg/m³, modulus of 8.6 GPa, uniaxial compression strength of 55.6 MPa and tensile strength of 3.01 MPa. The components of this sandstone were 79.1% quartz, 4.4% K-feldspar, 13.3% kaolinite, 2.0% illite and 1.2% other minerals [26].

Complete white sandstone samples of size 140 mm (height) × 70 mm (length) × 40 mm (width) were cut by water knife to produce an opening flaw in the middle. The angle between the flaw and the horizontal plane was denoted α, and this was set to 30°, 45°, 60°, 75° and 90°. The size of the opening flaw was 35 mm × 5 mm. There was also a transversely perforated borehole with a diameter of 4 mm, as shown in Figure 1a.

The anchor material was 304 steel, which is a common stainless steel material with proportions of cadmium and nickel of more than 18% and 8%, respectively. In addition, 304 stainless steel is a brand of stainless steel produced according to the standard of the American Society for Testing and Materials (ASTM). The anchorage mode at each end was divided into two types: free-ends anchorage mode and restrained-ends anchorage mode. The restrained-ends anchorage fixed both ends of the anchor with a nut and gasket, and there was nothing at either end of the free-ends anchorage type, as shown in Figure 1b. The rock–anchor interface was bonded by resin. After the anchor was installed, the weak layer was added to fill the opening flaw. The cement slurry was obtained by mixing water and cement in a weight ratio of 1:2. The opening flaw was filled with the cement slurry and subjected to repeated impacts to ensure that the cement slurry fully filled the entire flaw, as shown in Figure 1c. After the cement was set, the cement surface was watered and cured for 7 days. The final anchored sample with the interlayer is shown in Figure 1d. Testing started after the cement reached a certain strength.

In order to facilitate the description, three sample name abbreviations are used. The basic abbreviation format is “anchorage method-interlayer dip angle”. Three letters were used to denote the anchorage method: N is no anchor, F is free-ends anchorage and R is restrained-ends anchorage. For example, “R-30” represents the sample with an interlayer dip angle of 30° under the restrained-ends anchorage condition.

2.2. Test Systems

Experimental tests were carried out on a multifunctional testing machine, as shown in Figure 2. Basic mechanical tests for rock and rock-like specimens can be carried out on this machine, including uniaxial compression tests and direct shear tests. The maximum normal load of the testing machine was 500 kN. When testing, the axial load was applied with a displacement-controlled loading speed of 0.12 mm/min. According to the test criteria published by the International Society for Rock Mechanics, this loading rate reaches the standard for a quasi-static test [27].

In order to better identify the failure modes of the specimens, digital image correlation (DIC) technology was used in this experiment. DIC technology has been widely used in laboratory tests in rock mechanics in recent years [28,29,30,31], and it has become an important monitoring tool. A charge-coupled device (CCD) camera was used for photography, and the collected photos were imported into GOM Correlate software (GOM Correlate 2021, ZEISS Company, Oberkochen, Germany) for analysis of the digital images of the sample surface [32]. After recording the digital images of the sample surface before and after deformation, the DIC system computes the motion of each image point by comparing the digital images of the test object surface in different states [33].

2.3. Numerical Model Setup

At present, numerical simulation in geotechnical engineering is commonly performed using PFC software or FLAC software, both from the Itasca Corporation, Minneapolis, MN, USA. The main principle of the PFC software is the discrete element method [34]. In this software, particles are used to build models, and they are connected by contact. Adjacent particles interact with each other through internal forces and moments. The contact mechanics is reflected in the laws of particle interaction and internal force and moment updates. These laws include Newton’s law of motion and the force–displacement law, as shown in Figure 3a.

The basic principle of the FLAC software is the Lagrangian fast finite difference method [35]. Two built-in basic elements exist in this software for modelling: solid element zones and structural elements. Solid element zones can be given different constitutive models in order to better simulate different materials such as soil, rock, etc. Structural elements include shell, beam, cable, pile, etc., and each structural element has its own characteristics which can be used to simulate special structures in engineering (shield tunnel segments, piles, anchors, etc.). The Lagrangian fast finite difference method is used to simulate the plastic flow and yield of materials by the mixed discrete method. As shown in Figure 3b, in this method the model is divided into element difference grids and the motion equations are solved for grid nodes. The differential equations are approximated by difference equations, and the explicit time step method is used to solve the equations. In the calculation process, the new velocity and displacement are calculated using the motion equation, initial stress and boundary conditions. The strain rate is calculated from the velocity. All grid nodes and elements in the model are updated in multiple time steps.

When coupling the calculations of PFC and FLAC (coupled PFC–FLAC), the stress and deformation of a discrete body are simulated from the mesoscopic point of view through the PFC software, while the mechanical behaviour of the continuous medium is calculated from the macroscopic point of view through the FLAC software. The calculation data for the two types of software are not of the same type, so the coupling interface Socket I/O is used for the data exchange calculation [36], as shown in Figure 4a. In PFC, particles can only contact particles and walls.

Hence, in order to establish communication between these two software types, a wall zone (wall structure) for coupling was generated at the surface of the FLAC element, and particles could contact the wall zone, as seen in Figure 4b. PFC particles could transfer forces and moments to the wall through particle–wall contacts, and the wall exerted a corresponding mechanical effect on the FLAC element through nodes, and vice versa [37]. It is worth noting that when the coupling is carried out, the large-deformation switch must be turned on. In order to establish the same model as the experimental model, the same geometry and size was used in the numerical model. The rock and interlayer consisted of bonded particles in PFC, and the zone element in FLAC was used to simulate the anchor. Particles were generated in the specified sample domain and divided into two particle groups according to the location. Figure 5 shows the numerical sample model and contacts in the coupled PFC–FLAC program.

2.4. Parameter Determination and Calibration

2.4.1. Parameter Determination

Two types of cylindrical samples were produced from white sandstone and cement. The sizes of the two types of cylindrical samples were 50 mm (diameter) × 100 mm (height) and 50 mm (diameter) × 25 mm (height), respectively. The ends of the specimens were flat to ±0.01 mm [27]. Since all white sandstone samples were obtained by cutting a large white sandstone cube, two of each type of cylindrical white sandstone sample were required, considering the low discreteness of rock samples. However, considering errors in fabrication, five of each type of cylinder cement sample were required, to ensure statistically more reliable results. A uniaxial compression test was carried out to obtain the uniaxial compression strength (UCS), elastic modulus and Poisson’s ratio. The Brazilian splitting test was used in the tensile strength determination. In order to obtain the strain data for the rock and cement cylindrical samples for calculating the Poisson’s ratio, a lateral and a vertical strain gauge were attached at the middle position of each cylinder sample, as shown in Figure 6b. These were model BFH120-10AA-D-D150 strain gauges, with other characteristics as follows: an electric resistance of 120 Ω and a coefficient of sensitivity of 2.0 ± 1%. This type of strain gauge can be seen in Figure 6a. The Poisson’s ratio was obtained from the ratio of the lateral strain increment and the vertical strain increment in the elastic stage. The strain–time curves of the cylindrical samples are presented in Figure 6c.

One test result for each type of cement sample showed a large error compared with the other results, so only four of the five test results were identified as valid. The results of the tests can be seen in

Table 1. Parameters of white sandstone and cement in this study.

Parameter	Density kg/m3	UCS /MPa	Elastic Modulus /GPa	Poisson’s Ratio	Tensile Strength /MPa
Sandstone1	2331	32.9	5.32	0.160	2.93
Sandstone2	2373	34.6	5.14	0.166	2.99
Sandstone (mean values)	2352	33.75	5.23	0.163	2.96
Sandstone (standard deviation)	21	0.85	0.09	0.003	0.03
Sandstone (coefficient of variation)	0.89%	2.52%	1.72%	1.84%	1.01%
Cement1	1982	10.1	1.13	0.23	0.94
Cement2	1918	9.03	0.96	0.21	1.02
Cement3	1936	9.96	1.11	---	1.06
Cement4	1964	9.21	1.056	0.24	0.99
Cement (mean values)	1950	9.58	1.064	0.227	1.01
Cement (standard deviation)	24.7	0.462	0.066	0.012	0.044
Cement (coefficient of variation)	1.27%	4.82%	6.20%	5.28%	4.36%

. Two test results for the white sandstone rock presented good consistency, with the two specimens having similar uniaxial compressive strength, modulus and tensile strength, together with coefficients of variation of the rock parameters not in excess of 3%. The uniaxial compressive strengths of the cement samples showed some differences (the values were in the range of 9.03–10.1 MPa). The moduli of the cement samples were similar (the values were in the range of 0.96–1.13 GPa), but there was no stress–strain curve coincidence, owing to the influence of the compaction stage (caused by fabrication errors in the preparation of the cement samples). In general, although there were some small discrepancies in the cement parameters, the coefficients of variation of the parameters were still not in excess of 7%, which was within an acceptable range. Therefore, the averaged values of each parameter were taken as the parameter value and the parameter calibration standard in this study, as shown in Table 1. The anchor material was directly obtained from the steel factory, so the parameters of the anchor material could be obtained from the correlation standards, as shown in

Table 2. Parameters of anchor.

Parameter	Density kg/m3	Tensile Yield Strength /MPa	Elastic Modulus /GPa	Poisson’s Ratio	Ultimate Tensile Strength /MPa
Value	7900	450	198	0.3	620

.

2.4.2. Contact Model Selection and Parameter Calibration

In the PFC software, two contact models are usually used to simulate rock or rock-like materials: the parallel bond contact model and the flat-joint contact model [34]. A comparison between the parallel bond contact model and the 3D flat-joint model is shown in Figure 7. There are two parts in the parallel bond contact model: the linear part and the parallel bond part. The parallel bond is equivalent to a set of springs with constant stiffness, which can transmit forces and moments to particles or walls at both ends of the contact. The linear part is non-bonded, and this is manifested as an infinitely small, linear elastic and frictional interface. After parallel bond breakage, proper rotational resistance cannot be reproduced [38], as seen in Figure 7a. The contact plane of the flat-joint contact model produces multiple discrete elements, where each element can be bonded or unbonded. According to this feature, the states of the flat-joint contact plane are divided into three types (completely bonded, locally bonded and unbonded). After the fracture of the contact plane, these elements still provide interlocking and rotational resistance of particles, as shown in Figure 7b.

The tensile-to-compressive strength ratio is a key parameter for rock and rock-like materials. If the parallel bond contact model of particles is chosen in the numerical modelling, the tensile-to-compressive strength ratio of the sample obtained in the numerical simulation will not match the real situation [39,40]. However, every contact in the flat-joint model consists of some independent element, and the bond state of every element is also independent. Therefore, this improves the interlocking action of circular particles, which means that the tensile-to-compressive strength ratio of rock and rock-like materials using the flat-joint contact model is satisfactory [41].

The particle contact model was chosen in the flat-joint model. The parallel bond contact model was used between particles and the wall zone for simulating resin, and the parameters of this bond model are given in [42]. A linear model is widely used for the contact between particles and the wall for loading. Therefore, in addition to the contact between the wall zone and the particles, the contact between the wall and the particles of the loading plate was selected for the linear model.

Among the structure elements in FLAC, the pile and cable structural elements can simulate an anchor, but the structural element has no width or thickness, only length. Although there is a parameter assignment for a cross-sectional area in the model, the cross-sectional area is regarded as a point without size, so this structure element is only suitable for large-scale engineering simulations. The ratio of the sample width to the anchor diameter in this test was a finite value, so the influence of this ratio cannot be ignored. In summary, the elastic model was selected for the anchor, and the specific parameters were consistent with the actual parameters, as shown in Table 2.

The “trial and error” method has been applied to the parameter calibration by many researchers [20,42,43,44,45]. The microscopic parameter values are continuously adjusted to match the experimental and numerical macroscopic parameter values and failure modes. When the macroscopic strength, deformation parameters and failure modes obtained in the simulation are consistent with those in laboratory test, these microscopic parameters are regarded as the parameters that can be used in the simulation. The matching results for strength and failure are shown in Figure 8. The parameters of this test are shown in

Table 3. Microscopic parameter values used in simulation.

Microscopic Parameters	White Sandstone	Cement
Particle minimum radius Rmin/mm	0.9	0.9
Particle radius ratio Rmax/Rmin	1.5	1.5
Effective modulus E*/GPa	4.5	0.8
Stiffness ratio kn/ks	1.8	2.5
Tensile strength fj_ten/MPa	6.5	1.5
Bond strength fj_coh/MPa	8.5	3.2
Friction angle fj_fa/°	38	35
Friction coefficient fj_fric	0.15	0.1

, and the macro parameters and errors obtained from the test and simulation results are shown in

Table 4. Parameter value comparison between experiment and simulation.

	White Sandstone			Cement
Parameters	Experiment (Mean Value)	Simulation	Error	Experiment (Mean Value)	Simulation	Error
UCS/MPa	33.75	32.2	4.59%	9.58	9.62	0.42%
Tensile strength/MPa	2.96	2.93	1.01%	1.01	0.96	4.95%
Elastic modulus/GPa	5.23	5.67	8.41%	1.064	1.01	5.08%
Poisson‘s ratio	0.163	0.171	4.91%	0.227	0.205	9.69%

.

3. Experimental Results Analysis

3.1. Strength Analysis

The force–displacement curves and peak load relationships of the three groups are shown in Figure 9. White sandstone is a brittle rock, and the strength and stiffness of the anchor were much larger than those of the rock and the weak interlayer, so the load curve dropped suddenly after the peak. In the restrained-ends anchored sample, the pre-peak plasticity was larger, which meant that the anchoring effect was stronger. Under the anchoring effect, the elastic modulus of the anchored sample was similar, which indicated that the elastic modulus of the sample with a small dip angle increased more than that of samples with a larger dip angle. It can be seen from Figure 9d that the peak strength of unanchored and anchored samples improved with an increase in interlayer dip angle. In the case of free-ends anchorage, the maximum rate of increase of the peak strength was 21.05% when the dip angle was 30°. However, with α ≥ 45°, the improvement rates in strength were both flat and low and had a value of about 8%. This indicated that, under the condition of free-ends anchorage, strength improvement was significant only for the sample with an interlayer dip angle of 30°. In the case of restrained-ends anchorage, the rate of increase of the peak strength decreased with an increase in α, indicating that the anchoring mode had a greater effect on the interlayer with a small dip angle. When the interlayer dip angle was large, the strength of the non-anchored sample was closer to that of the intact sample, so the anchorage effect was relatively small. Regardless of the anchorage method, an intuitive effect was reflected in the improvement in the peak strength.

3.2. Failure Pattern

3.2.1. Failure Mode

The failure mode represents the failure characteristics shown by a class of samples. The summary represented by the failure mode can often help in efficiently studying samples with the same failure mode, so the effect of other factors on the sample failure can be quickly and comprehensively obtained. In the unanchored samples, three failure modes were found, with specific failure details as follows:

(1)

Failure mode I: Combined with the instantaneous failure image captured by the high-speed camera, it can be seen that a shear crack occurred at the tip of the interlayer, and then the tensile stress dominated the failure. The blocks at the left and right tip moved to the outside under the action of the tensile stress to give the final failure. There were not many powders or particles after friction and sliding. This failure mode occurred in the samples with interlayer dip angles of 30° and 45°, as shown in Figure 10a.

(2)

Failure mode II: A tensile crack was initiated and propagated at the tip of interlayer, leading to upper crack coalescence and tensile failure; from the lower left tip to the lower left corner, a whole shear friction slip occurred. Shear failure occurred at the right tip, and the tensile crack propagated to the lower right corner. The failure of the right part was a mixed failure. This failure mode occurred in samples with interlayer dip angles of 60° and 75°, as shown in Figure 10b.

(3)

Failure mode III: The tensile crack was hard to propagate from the tip of the interlayer. Only one mixed crack extended upward from the bottom interlayer tip. Surface spalling and block collapse occurred in the right part of the sample. The failure degree of samples showing failure mode III was closer to that of intact samples. This failure mode occurred in the samples with interlayer dip angles of 90°, as shown in Figure 10c.

3.2.2. Crack Type

Researchers have classified the crack types in fractured rock under uniaxial compression [46,47,48]. The crack types include direct tensile cracks at the tip, secondary shear cracks after a shear crack at the tip, mixed tensile–shear cracks at the tip and so on. The most important method of distinguishing between tensile and shear cracks involves judging whether friction sliding has occurred at both sides of the crack. According to the experimental results, friction sliding was also used to judge tensile and shear crack types in this study. Crack types in all failure patterns of the samples were identified and classified, as shown in Figure 11. These crack types can be described as follows:

(1)

Tensile crack: a type 1 tensile crack was a wing crack, which can be observed in the samples with dip angles of 30°, 45° and 60°; a type 2 tensile crack was an anti-wing crack; type 3 tensile cracks mainly appeared in samples with dip angles of 75° and 90°, where the initiation angle was small and the crack was initiated directly from the tip along the outer normal direction of the interlayer.

(2)

Shear cracks: type 1 shear cracks propagated along the dip direction of the interlayer, mostly caused by the overall sliding shear of the upper and lower parts; type 2 shear cracks and type 2 tensile cracks occurred at the same location, but S₂ cracks were driven by shear stress.

(3)

Mixed cracks: the initiation of a type 1 mixed crack was driven by shear stress, and propagation was caused by tensile stress; a type 2 mixed crack consisted of an S₂ shear crack and an anti-wing crack; the tip of an M₃ mixed crack was driven by tensile stress, and the shear stress was dominant in the lateral shear propagation after crack initiation.

3.2.3. Anchorage Effect on the Final Failure Pattern

In order to show the crack types more clearly, the completely destroyed specimens were observed, and the failure diagram was drawn and marked with the crack types identified in the previous section. Each sample failure pattern consisted of three parts: the experimental pattern (top), the DIC pattern (middle) and the drawn failure diagram (bottom). The results for all the samples’ failure patterns are shown in Figure 12.

Figure 12a shows the schematic failure diagram for unanchored samples. When the interlayer dip angle was <90°, with increasing dip angle the crack propagating downward at the lower left tip gradually changed from a tensile crack (30°) to a mixed crack (45°) and a shear crack (60° and 75°). At the same time, the crack propagating upward at the right tip gradually changed from a mixed crack (30° and 45°) to a tensile crack (60° and 75°). Anti-wing cracks propagated well, and the anti-wing crack propagating downward on the right tip also gradually changed from a tensile crack to a shear crack with an increase in dip angle. Anti-wing crack coalescence occurred at both sides in sample N-75. When the interlayer dip angle was 90°, the sample integrity was best, and it was most difficult to destroy the sample.

Figure 12b shows the schematic failure diagram of the free-ends anchored samples. For failure mode I (α = 30°, 45°), compared with the unanchored samples, the T₂ tensile crack at the left tip was suppressed, and the shear sliding coalescence failure was more evident at the lower left side. The T₁ tensile crack at the tip was initiated but not propagated. The shear failure of the free-ends anchored samples was significantly increased. For failure mode II (α = 60°, 75°), the T₃ tensile crack at the upper tip was suppressed, and it was unable to propagate upward. The left anti-wing crack type changed from tensile to shear, and the number of cracks was significantly reduced. Crack initiation, propagation and coalescence at the tip were difficult. For failure mode III (α = 90°), the tip tensile crack could not be initiated. Only one oblique shear crack and its secondary crack appeared, and the crack could not be clearly observed on the surface of the sample. It can be seen that the effect of the free-ends anchorage was reflected in the inhibition of crack initiation at the tip, preventing the tensile effect, increasing the normal stress of the failure surface and turning tensile cracks into tensile–shear mixed cracks.

As seen in Figure 12c, for R-30–R-75, the shear stress effect dominated, and the tensile crack at the tip was completely suppressed. The S₁ shear cracks at the upper and lower tips completely coalesced, and a large amount of powder and surface spalling occurred during the shear process. Due to the restrained-ends anchorage effect, it was hard for the rock mass on each side of the interlayer to expand outward. A rock–anchor co-deformation area was formed in the middle part of the samples, and a shear anti-wing crack replaced the tensile anti-wing crack. For R-90, the anchorage force acted directly in the normal direction of the interlayer, which meant that the tip of the interlayer was unable to crack. At the same time, the internal strain energy of the samples accumulated and could not be released. With increasing load, a sudden collapse occurred, when the strain energy reached the limit of the sample capacity.

Clearly, with the help of DIC technology, more micro-cracks (S₂ shear cracks in samples F-75, R-45 and R-60) could be observed, making the failure pattern analysis more comprehensive. Similarly, the DIC technique also played a key role in the observation of micro-cracks, surface point displacement and strain fields [28,29,31]. This technique provides great help in the study of fracture damage and crack propagation in rocks and rock-like materials.

3.3. Anchorage Effect on Crack Prevention

Figure 13a shows the crack initiation load for all the samples under the three anchorage conditions. It can be seen that when α is 30°, there were obvious differences in the crack initiation load under the three anchorage conditions. The crack initiation load was smallest for sample N-30 and largest for sample R-30. When 45° ≤ α ≤ 75°, the crack initiation load of the anchored samples was larger than that of the unanchored samples, but there was no great difference in the crack initiation load for the two anchorage conditions. This indicated that the effects of the two anchorage types on the inhibition of crack initiation were similar. The crack initiation load of sample N-90 was almost equal to the sample’s peak value, while crack initiation was not observed in samples F-90 and R-90. When α was 90°, the additional force generated by the anchorage was applied in the direction normal to the interlayer, resulting in the complete closure of the rock’s weak layer structural plane.

Figure 13b shows the crack initiation-to-peak load ratio for all samples. When α was 30°, the load ratio amplification was a maximum under anchorage, illustrating that anchorage worked best for crack prevention in samples with an interlayer dip angle of 30°. When α was 45°, the crack initiation-to-peak load ratio of F-45 was roughly equal to that of N-45, and the ratio of R-45 was even lower than that of the former two. This shows that the effect of anchorage on crack initiation was not straightforward, and additional stress caused by the restrained-ends anchor occurred under a higher load situation when α was 45°, only preventing further crack propagation after crack initiation. When α ≥ 75°, the initiation load was closer to the peak load. The tangential stress at the interlayer surface decreased, and the normal stress increased. The rock–interlayer surface tended to be more closed, which also meant that the crack could rapidly expand after the peak. Therefore, this showed that for the samples with large dip angles, the role of the anchorage was mainly reflected in preventing early crack initiation, and the anchorage effect was not evident after crack initiation.

4. Numerical Results Analysis

4.1. Experimental and Numerical Results Comparison

Laboratory testing and numerical simulation are two methods of studying the same problem. When the variables and parameters are the same, the results should be roughly the same. Therefore, in order to verify the accuracy of the results of the two methods and seek more realistic results, the force–displacement curves and failure diagrams of the test and simulation results are compared in Figure 14.

Some samples with typical failure characteristics were selected for comparison. There were no micro-fissure closure stages in the force–displacement curves in the simulation results. When establishing the model, suspended particles often exist, as seen in Figure 15a. Suspended particles may lead to a large and unbalanced force in the initial model. At the same time, the suspended particles contact the surrounding particles in the compression process, resulting in an increase in modulus and strength. In fact, this increase is more in line with the compaction process in real rock compression. However, the increase in modulus and strength in the compression process is uncontrollable, and its specific value cannot be calculated. Therefore, it is impossible to determine the elastic modulus value of the model when it enters the elastic stage, which affects the parameter calibration. Therefore, in the modelling, suspended particles must be removed, and the method involves making suspended particles contact the surrounding particles as much as possible. When the particles are not suspended, the surrounding pores are incompressible, as seen in Figure 15b. Whatever contact model is used, compressible pores are completely eliminated in the modelling process, so the simulation curve enters the elastic stage at the beginning [49,50]. In fact, there were still some errors in the curves and failure patterns between the numerical and test results, but they were within the acceptable range. Due to the large difference in the modulus between the anchors and the samples, splitting failures occurred at the post-peak stage, and this could be observed in most of the anchored samples. This is shown in the comparison of the lateral failures of the F-75, R-60 and R-90 samples. Overall, most of the simulation results matched the experimental results, which verifies the correctness of all the results.

4.2. Analysis of Anchorage Effect on Internal Stress Evolution

“Measure circle” is a monitoring tool in the PFC software. In the measure circle tool, the stress, strain rate, porosity and other measurement quantities can be computed and read off for further study [34]. Three stress states can be computed: stress in the X direction (horizontal stress), stress in the Z direction (axial stress) and stress in the Y direction. In addition, the shear stress in the XOZ plane can be computed. In order to monitor the internal stress evolution of the sample, taking the centre of the sample as the origin, 400 measurement circles with radii of 0.007 m were arranged in a square area of side 0.07 m in the XOZ plane. The measurement circles had a certain amount of overlap to ensure that particles in the range were monitored completely. The arrangement of measurement circles is shown in Figure 16.

It can be seen from Section 3 that the anchorage effect was more obvious for the samples with small interlayer dip angles. When the dip angle was large, the interlayer gradually lost the trend of slip, so the difference in stress evolution was not clear. Therefore, a sample with an interlayer dip angle of 30° was selected to study the evolution process of horizontal stress, axial stress and shear stress at different key time points in the failure process. When the force or stress of a contact exceeded its capacity limit, then the contact was defined as a breakage contact. At the same time, a crack (blue area in the failure process sample) was initiated at the position of the breakage contact in order to better observe the failure process.

4.2.1. N-30 Sample

Figure 17 shows the failure process and stress evolution diagram of the N-30 sample. Wing cracks initiated and propagated sufficiently; when the axial load reached the peak value (point C), an anti-wing crack propagated at the right lower tip of the interlayer, and shear coalescence began at the left lower tip. The sample showed complete failure with a small axial displacement after the peak.

At point A, there was a stress concentration phenomenon at the interlayer tips. The horizontal stress at the wing crack initiation position was a tensile stress (tensile stress is positive, compression stress is negative), which also revealed that the wing crack was a tensile crack. An axial compressive stress difference occurred at each end of the initiated wing crack, and there was a large compressive stress concentration on both sides of the interlayer. Once the stress exceeded the capacity limit at this point, the sample would undergo crack coalescence failure. There were shear stresses of opposite sign at the sides of the crack tip, indicating that there was a friction trend here.

The wing crack propagated completely at point B, and failure occurred at the tips. This led to the horizontal tensile stress concentration gradually moving away from the tips, and the compressive stresses on both sides of the interlayer were further aggregated. Crack propagation concentrated the shear stress from a large area to a small area, and the value of the shear stress showed no significant improvement.

At peak point C, the horizontal stress concentration disappeared at the right tip but remained at the left tip, because crack propagation had occurred in the right part of the sample, and stress was released. In contrast, a tensile anti-wing crack was initiated at the right tip. The axial stress changed considerably. The upper and lower cracks coalesced on the right side, so the stress concentration area disappeared. It can be seen from the shear stress on the left tip that there were two opposite directions of shear stress concentration, and shear crack propagation was driven by these two stress concentrations. A small part of the shear trend existed at the tips, and this was consistent with the test results.

At point D, the cracks completely coalesced, and the stress concentration area on both sides of the interlayer tips disappeared. Both sides of the sample completely lost bearing capacity, and the upper and lower sides of the interlayer–rock surface closed, so the horizontal stress was concentrated in the interlayer area. Similarly, the minimum axial compressive stress value occurred on the extended crack path. At this point, there was still a stress concentration at the lower left tip, which may have been caused by the interlocking of particles. After sample failure, the distribution of the shear stress showed disorder, and the stress value was also lower than the values at the three points A, B and C, which represented the completion of the shear crack propagation.

4.2.2. F-30 Sample

The crack behaviour of the F-30 sample was similar to that of the N-30 sample, as shown in Figure 18, but the wing crack propagation of the F-30 sample was less marked than that of the N-30 sample. During the wing crack propagation, a shear crack also started to propagate. This hindered the development of wing cracks and gradually transformed the coalescence type of the samples from tension to shear.

From the horizontal stress at points A and B, it can be seen intuitively that there is a small compressive stress area on both sides of the middle point of the weak interlayer, which was caused by the friction resistance of the horizontal anchorage. There was no difference in the distribution of horizontal stress at other sample locations between the N-30 and F-30 samples at point A. The wing crack propagation at point B led to the dispersion of the horizontal tensile stress. Due to the addition of the anchor, the axial compressive stress on the line Z = 0 was the smallest, so that the middle part could withstand compressive stress without failure. This was the reason for the increase in the modulus and strength of the anchored samples. In the case of a large axial load, the shear stress of the F-30 sample was smaller than that of the N-30 sample at the two points A and B, which intuitively reflects the weakening effect of the anchorage on the shear stress.

Under the anchorage effect, at peak point C there was a still larger compressive stress concentration at the middle of the interlayer, preventing the lateral deformation of the F-30 sample. Even during the period of crack coalescence, the axial stress on both sides of the interlayer at point C remained large. It was observed that the high strength and stiffness of the anchor could continue to provide axial load-bearing capacity for samples when sample failure occurred. In fact, the stress concentration at the tip was weakened by the anchorage effect, which further enhanced the bearing capacity of the samples. The shear stress of the F-30 sample was concentrated at the peak and after the peak, which proved that the anchorage effect weakened the shear stress and delayed shear failure.

The distribution of the horizontal tensile stress at point D was not uniform, which meant that the anchor–rock interface was unbonded. Only the middle anchorage area could bear the residual axial stress when the sample was broken, so the axial stress was concentrated on the upper and lower sides of the interlayer. There were two parallel and reverse shear stresses on the upper and lower sides of the anchor, indicating that horizontal friction sliding occurred. It is worth noting that, although the anchor–rock interface was unbonded, the shear stress at point D was still larger than at the three points A, B and C.

4.2.3. R-30 Sample

The choice of points in the R-30 sample was different from that in the N-30 and F-30 samples. There was no wing crack during the crack initiation process, and when the crack was initiated the axial load had reached 68.15 kN. The main crack propagation stage existed after the peak, as seen in Figure 19.

The horizontal tensile stress value for the R-30 sample at point A was similar to that of the first two samples, but the axial load was much greater at this time. This means that the effect of the restrained-ends anchorage was more obvious in weakening the tensile stress, so the wing crack was difficult to initiate. Similarly, the axial stress concentration at the tip was not obvious in the R-30 sample. The restrained-ends anchorage also weakened the axial compressive stress concentration at the tips, making it difficult for the crack to initiate. Compared with the N-30 and F-30 samples, the shear stress distribution and values for the R-30 sample did not change much under a higher axial load. It can be seen that at the stage of crack initiation, the reduction effect of the restrained-ends anchorage on each stress was greater than that of the free-ends anchorage.

At points B and C, the horizontal compressive stress area in R-30 was wider than in the F-30 sample. The restrained-ends anchor absorbed most of the axial stress, allowing the sample to bear only a little axial stress, which indirectly improved the axial load-bearing capacity of the R-30 sample. The shear stress of the sample was distributed along the upper and lower sides of the anchor, extending to the interior of the weak interlayer at point C. An anchorage zone formed between anchor and the sample, resisting lateral expansion deformation.

The horizontal stress distribution in the R-30 sample at point D can be used to show the difference between the restrained-ends anchorage and the free-ends anchorage. At this point, most of the horizontal stress area was a compressive stress area (green and red areas), and the stress was tensile only in the four corners. However, in the F-30 sample, the horizontal compressive stress area was only distributed near the anchor. This showed that the scope of action of the restrained-ends anchorage was greater than that of the free-ends anchorage.

From point A to point D, the shear stress had a larger value in the internal weak interlayer, while the larger value of the shear stress in the F-30 sample only existed on the surface of the weak interlayer. The restrained-ends anchorage effect led to complete closure in the interface between the rock and the weak interlayer in the anchorage zone. In the case of a large modulus difference, the stress could be transmitted through the rock–interlayer interface, which further improved the integrity of the sample. However, due to the limited free-ends anchorage effect, the interface closure was not complete, and the shear stress was only distributed along the surface or in the sliding direction.

4.3. Stress Analysis for Anchor

The stress monitoring of a zone element in FLAC provides great help in obtaining the anchor stress. Owing to the existence of the interlayer, the stress distribution may be complex, so more monitoring points were arranged near the interlayer. The monitoring-point arrangement is shown in Figure 20. In order to study the stress distribution for different interlayer dip angles, three groups of samples were chosen, with dip angles of 45°, 75° and 90°. There were three axial and shear stress distribution curves for each sample, and the sample states in the three curves showed three stages: the elastic stage, the elastic–plastic stage and near the peak, respectively. In order to distinguish the three stages, the sample axial displacement was used to mark the curves, as shown in Figure 21, Figure 22 and Figure 23.

4.3.1. Stress Analysis for Anchors in Samples F-45 and R-45

Figure 21 shows the distribution of axial stress and shear stress for the F-45 and R-45 samples. The axial stress was, essentially, large in the middle and small at both ends. The maximum axial stress value was not located at the middle position of anchor but in the central position of the connection between the interlayer’s left edge and the anchor’s middle point. The axial stress distribution of the anchor in the sample with an inclined interlayer was not symmetrical. A higher axial stress occurred at the side with shear sliding. When the sample load was near the peak, the axial stress of the anchor increased from −35 mm to −10 mm in sample R-45, and the stress at the ends was not equal to zero, as shown in Figure 21b. This is because the restrained-ends anchorage fixed the two end points of the anchor, and after the peak, cooperative deformation occurred between the restrained-ends anchor and the sample. In contrast, the free-ends anchor could move freely after the rock–anchor interface was unbonded, so the axial stress did not increase in this area, as shown in Figure 21a.

The shear stress distribution showed three peaks, in the middle and at both ends, indicating that the anchor shear stress distribution was relatively complex. The shear stress at both ends of the restrained-ends anchor was significantly larger than in the free-ends anchor, which was also the result of the fixed end points of the restrained-ends anchor.

4.3.2. Stress Analysis of Anchors in Samples F-75 and R-75

The regularity of the axial stress distribution in the anchor in the F-75 and R-75 samples was similar to that in the F-45 and R-45 samples, as shown in Figure 22. The axial stress of the anchor in the F-75 sample increased sharply when the sample axial displacement was 0.65 mm, but the stress value of the anchor in the R-75 sample increased steadily with an increase in sample axial displacement. This reflected the cooperative deformation of the restrained-ends anchorage sample.

The maximum shear stress value of the anchors was obtained near the peak, and two negative peak values occurred at both sides of the anchor’s middle point. From the positive peak value to the negative peak value, the point where the shear stress was zero was a plastic hinge point. The restrained-ends anchorage effect limited the anchor ends’ deformation and shear stress transmission, so the differential value of the positive–negative peak shear stress was larger, as seen in Figure 22b. The anchor bending moment had its extreme values at the two plastic hinge points, so the anchor had the largest bending deformation between these two points.

4.3.3. Stress Analysis of Anchors in Samples F-90 and R-90

When the interlayer dip angle was 90°, the axial tensile stress of the anchor was caused wholly by the horizontal deformation of the samples. In the elastic and elastic–plastic stages, the axial tensile stress values in the middle part of the anchors were roughly equal, as seen in Figure 23. However, due to the low stiffness of the interlayer, the stress transmission in the interlayer–anchor interface was small. Hence, two strong contact points were formed between the anchor and the rock on both sides of the interlayer, and the axial tensile stress of the anchors at these two points was also largest when the sample load was near the peak. The stress transmission at the anchor–interlayer interface was small, so the anchor axial tensile stress distribution in the middle was reduced. Owing to the restrained-ends anchorage effect, the two ends of the anchor in the R-90 sample could not move freely, but both ends of the anchor could withstand a greater tensile stress. Therefore, as seen in Figure 23b, the anchor axial stress at the middle point in the R-90 sample was lower than that in the F-90 sample.

The shear stress of the rock–interlayer interface was zero, and there was no sliding trend at this interface when the dip angle was 90°, so the shear stress values of the anchors were small. The distribution of shear stress was disordered and fluctuated greatly, which may have been caused by rotation and the stress fluctuation of the anchor–rock interface bonded friction.

4.4. Comparison of Anchor Deformation

The comparison of the anchor deformation between the simulation results and the experimental results is shown in Figure 24. In Figure 24, the anchor displacement vector, the anchor displacement in the Z direction and the anchor deformation in the experimental results are arranged from top to bottom. At the same time, different colour maps of the deformation in the numerical results are placed and marked on both sides. Overall, there was agreement between the anchor deformation and stress distribution, except for the anchor in the R-90 sample, where the maximum deformation of the anchor appears on both sides of the midpoint. The bending of the middle part was caused by the shear effect, and a plastic hinge was formed in the middle of the anchor.

5. Conclusions

In this paper, an experimental and coupled PFC–FLAC simulation study was used to conduct uniaxial compression tests on three types of anchored white sandstone samples with weak interlayers. The strength variation laws and failure modes of the samples were obtained. Moreover, analyses of internal stress variation, anchor stress distribution and deformation were performed. The specific conclusions are as follows:

(1)

The peak strength values of all samples improved with an increase in the interlayer dip angle. When the interface dip angle was small, the anchorage effect on the strength and elastic modulus increase was greater. The effect on strength of the restrained-ends anchorage was always greater than that of the free-ends anchorage. The crack initiation strength and initiation-to-peak strength ratio of the sample with a dip angle of 30° was increased by the anchorage effect, but the improvement effect on other samples was not obvious.

(2)

There were three failure modes for the unanchored samples. The effect of anchorage on modes Ⅰ and Ⅱ was to limit tensile crack initiation and to drive shear failure, and the effect of the restrained-ends anchorage was more obvious. For mode III, a crack could not initiate before the peak, and the internal strain energy could not be released under the anchorage effect, so the sample broke instantly after the peak.

(3)

The restrained-ends anchorage effect on the lateral reinforcement range was larger than that of the free-ends anchorage. Shear cracks were driven by two shear stresses in opposite directions on each side. The anchor contributed to bearing the axial load, and thus the stiffness and strength of the samples were increased. The anchorage effect resulted in the shear stress passing through the interlayer. However, when the rock–anchor interface was unbonded, the shear stress in the interlayer of the free-ends anchored samples disappeared. The restrained-ends anchor could provide a continuous anchorage force, so the shear stress in the interlayer did not disappear even if the interface was unbonded.

(4)

The larger the dip angle, the smaller the axial and shear stress of the anchor at the middle point. The restrained-ends anchor had an obvious tensile stress distribution at both ends. The axial tensile stress distribution in the middle of the anchor was reduced when α was 90°, and the degree of reduction in the R-90 sample was larger. Larger dip angles made it hard for interface sliding to occur, so the shear stress value was less.