Effect of Multiple Hole Distribution and Shape Based on Particle FlowonRocklikeFailureCharacteristics andMechanicalBehavior

Based on the particle flow code, numerical models of vertical and horizontal orientations of holes with different shapes were established, and the effects of preexisting holes with different shapes and arrangement patterns on the mechanical behaviors and failure characteristics of rocklike materials were studied.(e evolution trend of the stress field is discussed by taking a circular hole as an example. (e results show that the existence of holes reduces the peak stress, peak strain, and elastic modulus of the sample, and different shapes of holes and different arrangement patterns have different effects on the mechanical properties and damage degree of the sample and significantly affect the horizontal orientation model. Before crack formation, the compressive stress and tensile stress concentration areas of each sample are located at the left and right ends and the upper and lower ends of the hole, respectively. After model failure, the compressive stress and tensile stress concentration areas of each sample are relatively scattered. In the vertical orientation model, the middle area of vertical holes is the main compressive stress concentration area, which is approximately “columnar” distribution. In the horizontal orientation model, the compressive stress concentration area between the holes is cross distribution and approximately “X” type distribution. (e vertical orientation model sample forms a “columnar” distribution to bear the applied load with a more favorable bearing orientation.


Introduction
Most natural rock masses have various defect structures such as fractures and holes, affecting the mechanical properties of rocks [1][2][3][4][5]. High tectonic stress exists in the deep rock mass. During the construction of tunnels, reservoirs, mines, and other underground projects, the original stability of rock mass is destroyed, easily leading to instability of rock mass and rock burst, as shown in Figure 1 [5][6][7][8][9][10][11][12][13][14][15][16]. At the same time, the distribution of internal stress in the rock mass is more complex due to the effect of fracture hole defects [17][18][19][20][21][22]. erefore, it is important to study the failure characteristics and mechanical behavior of defective rock mass to understand the instability failure trends of defective rock mass and predict the engineering geological hazards of rock mass.
Many studies have been conducted on the failure characteristics of circular holes and the coalescence mode of cracks, and some meaningful conclusions have been obtained [24][25][26][27][28][29]. Wong et al. [30] conducted a series of uniaxial compression physical and numerical tests on singlehole specimens with different diameters and widths and studied the splitting failure, failure mode, and strength characteristics caused by crack propagation. Huang et al. [31] studied the strength failure behavior of three noncoplanar circular hole specimens under uniaxial loading using experimental and numerical simulation methods, determined four typical crack coalescence modes, and elucidated the crack evolution mechanism around preexisting holes in granite specimens. Wong and Lin [32] studied crack initiation, coalescence mechanism, and failure behavior of granite materials with multiple circular holes under uniaxial compression by physical and numerical experiments, respectively, and proposed a coalescence mechanism criterion related to circular hole distribution. Zhang et al. [33] studied the crack propagation characteristics of rocks around pressure relief drilling with different numbers and different arrangement patterns using experimental methods, studied the unloading effect of three holes with different geometric shapes using FRACOD program, and proposed a new optimization method for pressure relief drilling arrangement parameters.
For hole defects with different shapes, some interesting results have been obtained [34][35][36][37][38][39][40]. Li et al. [41] used twodimensional particle flow code (PFC2D) to simulate the failure process of marble containing preexisting holes under uniaxial and biaxial compression conditions and analyzed the effects of preexisting hole shape, confining pressure, and rock heterogeneity on the mechanical properties and crack propagation of marble. Feng et al. [42] studied the failure characteristics and crack propagation of typical hard rock specimens during the unloading of central holes with different shapes using the finite element method (FEM/DEM) and discussed the effect of rock heterogeneity on the strength and range of failure around the central holes of specimens. Raymond and Herbert [43] conducted uniaxial compression tests on prismatic gypsum samples containing one or two inclusions and found that the shear crack propagation degree of samples containing two inclusions at the inclusion boundary significantly increased more than single inclusions. Gui et al. [44] used mixed continuous discrete element method to conduct numerical studies on the effect of circular, rectangular, and triangular holes on rock mechanical behavior, proposed a simulation method that could reasonably reproduce the crack initiation and propagation, and better predicted the mechanical property changes caused by holes.
e above studies were mainly focused on circular holes or single holes with different shapes, while studies on multiple holes with different shapes and different arrangement patterns of multiple holes are scarce. In addition, previous studies mainly focused on the crack propagation and failure modes of rocks containing defects, and relatively few studies were conducted on stress distribution and stress field evolution before crack formation and after crack failure of the sample. However, various types of defective holes exist in the natural rock mass. ese defects pose a potential threat to the stability of rock mass. erefore, PFC2D program was used in this study to establish vertical and horizontal orientation models containing four types of holes with different shapes, studied the effect of hole shape and arrangement on the failure characteristics and mechanical properties of samples, and analyzed the evolution of the stress field in a circular hole model.

Particle Flow Code.
Particle flow theory was proposed by Cundall and Strack using the discrete element method [45]. Rigid particles are usually used to characterize the rock materials in particle flow code (PFC). e force and displacement between particles are achieved through contact. Among them, the contact bond (CB) and parallel bond (PB) can be used to simulate the connection between rock particles; the PB not only transfers force but also transfers moment [46][47][48]. In this study, a uniaxial compression model of a specimen with hole defects was established using the parallel bonding method.

Determination of Microparameters.
When PFC2D is used for a numerical simulation test, the microscopic parameters of particles are important [49]. In this paper, the microscopic parameters were repeatedly adjusted through "trial and error" until the microscopic parameters satisfied the requirements of simulation analysis. e "trial and error" parameter check process of PFC model is shown in Figure 2. e parameter of mi is the Hoek-Brown strength parameter. A Shimadzu AG-X250 precision universal testing machine was used as the loading system. e universal testing machine possesses good reliability and high precision, and it can be used to perform regular compression and tensile tests. e displacement-based loading control method was used to perform the uniaxial compression test until the sample was  e peak strength and elastic modulus of the test sample and numerical sample are basically the same. Although the peak strains of the two samples are different, they are both within the acceptable range. Table 1 shows the micromechanics parameters based on PFC.

Establishment of the Numerical Model.
To study the effect of multiple holes with different shapes on the failure and mechanical behavior of samples, two types of hole models in vertical orientation and horizontal orientation were established. Each type of model had four types of samples with holes of different shapes, and each sample contained seven holes. ese two types of models were called the vertical orientation model and horizontal orientation model. e hole shapes of vertical and horizontal orientation models were circular, elliptical, square, and equilateral triangle; the areas of the holes of four shapes were equal; and eight models were established. A schematic diagram of the model is shown in Figure 4(a). e model is 100 mm high and 100 mm wide. e area of each hole is about 126.5 mm 2 . e distance between the midpoint of the hole at the center of the model and the midpoint of the surrounding hole is 25.4 mm. e diameter of the circular hole is 12.7 mm. e major axis of the elliptical hole is 18 mm, and the minor axis is 9 mm. e side length of the equilateral triangle is 17.1 mm. e side length of the square is 12.5 mm. Specific dimensions of the four hole shapes are shown in Figure 4(b). According to the size and location of holes in the model, the particles of holes were deleted, and a numerical model was established as shown in Figure 5.

Mechanical Properties of Samples.
e stress-strain curves and crack number-strain curves of intact rock and defective sample are shown in Figure 6. As shown in Figure 6(a), the peak strength, peak strain, and elastic modulus of the complete sample are 9.51 MPa, 0.0082, and 1.32 GPa, respectively. e crack development had three stages, namely, no crack formation, slow crack growth, and fast crack growth, and the stress-strain curve and crack number-strain curve are relatively smooth. Compared with the complete sample, the stress-strain curve of the defective sample fluctuates before and after the peak value. e stressstrain of each sample is different, and crack development is relatively smooth, as shown in Figure 6(b). is is mainly because of the existence of hole defects and different arrangements of holes, varying the original geometric structure of the sample and reducing its bearing capacity. e strain and peak strain at the initial occurrence of a crack in the defective sample are less than those of the intact sample, indicating that the initial crack development and final failure of the defective sample are earlier than those of the complete sample.
A comparison chart of the mechanical properties of the vertical orientation model, horizontal orientation model, and complete sample is shown in Figure 7. e peak stress, peak strain, and elastic modulus of the vertical orientation model and horizontal orientation model are all smaller than those of the complete sample, while the peak stress, peak strain, and elastic modulus of the vertical orientation model are all larger than those of the horizontal orientation model. In the vertical orientation model and horizontal orientation model, there is no obvious trend of peak stress and peak strain from the elliptical sample, triangular sample, and square sample to circular sample, but the elastic modulus gradually increases. Compared with the complete sample, the peak stress and peak strain of the elliptical sample and triangular sample in the vertical orientation model and horizontal orientation model decrease the most, by 63.4%, 35.7% and 65.5%, 41.5%, 49.4%, 26.3% and 65.9%, 45.1%, respectively, while those of the circular sample decrease the least, by 43.3%, 21.6% and 56.2%, 35.4%, respectively. e difference in peak stress, peak strain, and elastic modulus between the vertical orientation model and horizontal orientation model is at most 12.72%, 18.66%, and 5.29%, respectively, and at least 2.1%, 5.7%, and 0.6%, respectively.    and 0°in the vertical orientation and horizontal orientation, respectively. In particular, it should be noted that the partial cracks generated in the triangular defect model coalesce along a direction of about 60°. is is mainly due to the fact that the different shapes and arrangement patterns of defects affect the direction and degree of crack growth. When axial load acts on the sample, the stress distribution state inside the sample is changed, thus changing the propagation path of crack. is indicates that the initial development position and propagation trend of the crack are influenced by the arrangement pattern of holes and end position.
e number of cracks in the final failure of each sample of the vertical orientation model and horizontal orientation model is less than that of the complete sample, and the number of cracks in each sample of the vertical orientation model is greater than that of the horizontal orientation model, as shown in Figure 9. erefore, the existence of hole defects reduces the damage degree of the model, and the damage degree of the horizontal orientation model is larger than that of the vertical orientation model. According to the final development degree of cracks when the two models in Figure 8 are damaged, it can be found that elliptical holes are

Contact Force Distribution before Crack Formation and after Model
Failure. e contact force distribution before crack formation and after crack failure of the vertical orientation model and horizontal orientation model is shown in Figure 10. Before crack formation, the compressive stress of each sample is mainly concentrated at the left and right ends of the hole, and the tensile stress concentration area is mainly distributed at the upper and lower sides of the hole. Influenced by the vertical arrangement of holes, the two sides of middle vertical holes of each sample in the vertical orientation model are the main compressive stress concentration areas, approximately distributed in a "columnar" shape as shown by the yellow dashed line in Figure 10(a). is area is the main stress area and plays a major supporting role in the model. Because of the difference between the arrangement of holes in the horizontal orientation model and that in the vertical orientation model, the compressive stress concentration areas between the holes are distributed cross, which is approximately "X" type as shown in Figure 10(b). e compressive stress concentration area at the end of the hole near the left and right boundary of the sample is larger than that at the end of the other holes, indicating that the area near the left and right boundary of the sample also plays a good supporting role for the model. e stress concentration area of each sample is distributed around the hole before crack formation. When the model is damaged, because of the influence of model failure mode, although still compressive stress concentration exists around some holes, the main stress concentration area diffuses around the sample, mainly caused by the influence of model failure mode.  Advances in Civil Engineering

Stress Field Evolution of Circular Holes.
Because of a large number of established model samples and limited space, more common circular hole samples (circular hole samples of vertical and horizontal orientation models) were selected in this study as the research object to analyze the stress field evolution. Figure 11 shows the stress field evolution nephogram of circular holes (where the blue area is the main compressive stress concentration area, the red area is the main tensile stress concentration area, the tensile stress is positive, and the compressive stress is negative).
In the vertical orientation model, before the crack is initially generated, that is, when the strain is 0.368 × 10 −3 , the compressive stress and tensile stress concentration at the left and right ends and the upper and lower ends of the hole are obvious, respectively, and the stress decreases significantly near the tensile stress concentration area. With continuous loading, the compressive stress concentration range, tensile stress concentration range, and stress decreasing zone at both ends of the hole increase. When the strain is 5.55 × 10 −3 , affected by crack propagation, the compressive stress concentration area shifts to the left and right boundaries of the model. When the strain is 6.66 × 10 −3 , the model failure is further intensified, the range of compressive stress concentration is relatively dispersed, and the compressive stress concentration areas near the left and right ends of two circular holes on the right side disappear, mainly due to the release of stress caused by the propagation of cracks. When the model runs to 0.6 times of the peak stress (strain is 6.91 × 10 −3 ), although the model failure is relatively severe, still a compressive stress concentration exists between the two rows of elliptical holes, and the maximum compressive stress reaches 18 MPa at this time, still playing a supporting role on the model consistent with the distribution of contact force shown in Figure 10(a). In the horizontal orientation model, the distribution positions of compressive stress concentration area and tensile stress concentration area are the same as those in the vertical orientation model, located at the left and right ends and the upper and lower ends of holes, respectively, but the stress concentration range is relatively scattered, mainly affected by the hole arrangement patterns and vertical and horizontal spacing between the two holes. With continuous loading, the compressive stress concentration range first increases and then decreases, while the  Advances in Civil Engineering tensile stress concentration range increases all the time. However, still a high compressive stress concentration exists in the area between the holes in the vertical orientation, supporting the model. From the stress distribution nephogram of the vertical orientation model and horizontal orientation model, it can be observed that the arrangement pattern of holes significantly affects the stress field distribution of the model.
To further understand the failure mechanism of rock with holes, the stress distribution around circular holes under uniaxial compressive stress was calculated according to the elastic mechanics theory and the famous Kirsch analytical solution [50]. Kirsch's analytical solution is shown in equation (1), where σ r , σ θ , and τ rθ are radial, tangential, and shear stresses, respectively. p and q are vertical and horizontal stresses, respectively. θ is the angle between the horizontal axis and radius:  Figure 10: Distribution of contact forces before crack formation and after crack failure in vertical and horizontal orientation models (the black lines represent compressive stress. e red lines represent tensile stress. e dense and sparse lines represent the magnitude of stress. e yellow ellipse and yellow curve marked areas are compressive stress concentration areas, and the green ellipse marked areas are tensile stress concentration areas). (a) Contact force before crack formation and after crack failure of the vertical orientation model. (b) Contact force before crack formation and after crack failure of the horizontal orientation model. 8 Advances in Civil Engineering In this study, the radius of the circular hole under uniaxial compressive stress r � L, and the horizontal force q � 0. us, equation (1) can be converted into (2) e tangential stress distribution at the boundary of the circular hole under uniaxial compressive stress [51] is shown in Figure 12(a). Tensile stress concentration areas are formed at the top and bottom of the circular hole, and compressive stress concentration areas are formed at the left and right ends of the circular hole. is is consistent with the studies of compressive stress and tensile stress distribution in Section 3.3 of this paper. σ θ and p are the tangential stress and vertical stress, respectively. θ is the polar angle between the horizontal axis and polar radius, r. e tangential stress along the contour of According to the distribution of stress field in the circular hole model, the existence of holes causes compressive stress, tensile stress concentration area, and stress decreasing zone to appear around the holes, and the existence of regular holes can make the distribution of these stress concentration areas and decreasing zones more regular. is provides ideas for predicting and changing the stress magnitude and location of a certain area in engineering construction. e vertical and horizontal orientations of circular holes are taken as examples for simple analysis. To illustrate that stress distribution might be affected by the stress direction of the model, the spacing between holes, and the arrangement pattern of holes, the circular hole combination is decomposed, as shown in Figure 13, and the marked A and B areas in the figure are symmetrical areas. Figures 13(a) and 13(d), 13(b), and 13(c) have the same geometric figures in areas A, B and C, D, E, F and G, H, I, J and L, M, which are isosceles trapezoids and parallelograms, respectively, with equal areas. Because the loading method used in this study is axial loading, the geometric edges of external forces in Figures 13(a) and 13(d), 13(b), and 13(c) are different, resulting in different distributions of internal stress. is indicates that stress direction affects the stress distribution of the hole model. In Figure 13(e), the area between two circular holes encircled by the rectangle is mainly tensile stress and stress decreasing zone. In Figure 13(f ), the upper and lower ends of two holes encircled by the rectangle also have tensile stress and stress decreasing zone, but the tensile stress area between the two holes does not coalesce and still is dominated by compressive stress. is indicates that hole spacing affects the distribution of stress (although other holes will affect the stress of two holes, they have less influence on the stress in the vertical orientation, and no major analysis will be made). Although the I, J and L, M areas in Figures 13(g) and 13(h) are both geometric triangular areas, the different arrangement patterns of I and J areas in Figure 13(g) and the L and M areas in Figure 13(h) result in a large difference of stress distribution between them, indicating that the hole arrangement pattern affects the stress distribution of the model. In addition, the stress distribution of vertical and horizontal orientation models of holes with the same shape also illustrates this point.
In engineering construction, the concentration of high stress may cause rock burst, which might pose a threat to the safety of personnel in engineering construction. is analysis shows that to reduce or transfer stress, according to the general situation of engineering geology, based on the analysis of actual structural stress and in combination with engineering needs, the arrangement pattern, size, and spacing of holes can be designed to satisfy the needs of pressure relief or stress transfer.

Conclusions
(1) e existence of holes affected the mechanical properties of the sample and reduced the peak stress, peak strain, and elastic modulus of the sample. e effect of holes on the mechanical properties of the horizontal orientation model was greater than that of the vertical orientation model; that is, the different arrangement patterns of holes had different influences on the mechanical properties of the sample.
(2) e development and propagation trend of cracks were influenced by the arrangement pattern of holes and the position of their ends. e existence of holes reduced the degree of model failure, and the damage degree of the model in the horizontal orientation decreased more. e damage modes of both models were mainly tensile failure. (3) e compressive stress and tensile stress concentration areas were located at the left and right ends and the upper and lower ends of the hole, respectively. In the vertical orientation model, the area between the holes in the vertical orientation was the main compressive stress concentration area, which was approximately "columnar" distribution. In the horizontal orientation model, the compressive stress concentration area between the holes is cross distribution and approximately "X" type distribution. (4) e arrangement pattern, size, and spacing of holes were important factors that affected the stress distribution of the model, which can provide ideas for predicting and changing the stress magnitude and location of a certain area in engineering construction.

Data Availability
All data are available within the article and also from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.