Analytical Model of Wellbore Stability of Fractured Coal Seam Considering the E ﬀ ect of Cleat Filler and Analysis of Inﬂuencing Factors

: Currently, coal borehole collapses frequently occur during drilling. Considering that the coal near to the wellbore is cut into blocks, and the cleat ﬁller of the coal inﬂuences the stress distribution near the wellbore, a new theoretical solution of a near-wellbore Stress Field in coal bed wells is established. In addition, according to the limit equilibrium theory and the E.MG-C criterion, the limit sliding formula of the quadrilateral and triangular block is deduced, and the slipping direction of the blocks is further judged. Finally, the wellbore stability model of the coal seam is established. The accuracy of the theoretical model is veriﬁed through a numerical method by using the PFC software. Based upon this wellbore stability theoretical model of coal, many cleat a ﬀ ecting factors such as cleat spacing, cleat length, cleat angle and the cleat geometric position, are studied, and the results show that a quadrilateral block slides o ﬀ more easily than a triangular block under the same boundary condition; the bigger the cleat spacing and cleat length are, the lower is the risk that blocks slide o ﬀ , and increasing the cleat angle could cause blocks to slide o ﬀ easily. Under the same boundary condition, whether blocks slide o ﬀ or not is closely related to the well round angle.


Introduction
It is an essential part of the unconventional development of energy to exploit coalbed methane in China [1]. Under the condition of natural strata, the coalbed presents its discontinuous and fragile properties, due to the limitation of the structural characteristics, such as natural fracture extension, abundant cleat system, low elastic modulus and strong anisotropy [2,3]. When the drilling contacts with this kind of broken coal seam, it tends to cause the collapse of the coal wellbore by the disturbance of the drilling string and the intrusion of the drilling fluid [4][5][6]. Until now, many specialists have carried out a number of studies about the stability of the sandstone wall through the continuum theory, which is useful for quasi-continuous reservoirs, such as sandstone, but is not useful for a broken coal seam. So, it is urgent to study the mechanism of the wellbore instability of the broken coal seam.
Specialists have achieved a number of results [7][8][9][10] on the study of the wellbore stability of oil and gas wells, which can give us ideas on solving the problems of a broken coal seam. Zhang et al. [11] proposed the coupling system of both flow and solid in rock failure process analysis, in order to

Stress Analysis of the Wellbore Surrounding Rocks as a Continuous Body
In order to obtain the stress distribution of broken coal seams, firstly we assume that the coal rock is a continuous body. Then, the stress field around the wall of the well can be obtained under non-uniform, in-situ stress and wellbore pressure P0:  Figure 1. Model of coal wellbore stability.

Stress Analysis of the Wellbore Surrounding Rocks as a Continuous Body
In order to obtain the stress distribution of broken coal seams, firstly we assume that the coal rock is a continuous body. Then, the stress field around the wall of the well can be obtained under non-uniform, in-situ stress and wellbore pressure P 0 : where σ r , σ θ , τ rθ , are radial stress, cyclic stres, and shear stress under polar coordinates, respectively. 2σ m = σ H + σ h , 2σ n = σ H − σ h ; R 0 is the wellbore radius; r is the length of the calculation point to the wellbore central axis; θ is the wellbore angle; α Biot is Biot coefficient. Then wellbore rock displacement [23] can be expressed as follows, and the specific derivation process of the formula can be referred to paper 23: (12 − 6µ)K 1 r 3 − (1 + µ)K 2 r + (6 + 4µ)K 3 /r + 9(1 + µ)K 4 /r 3 u θ = 6 sin 2θ where E, µ are the elastic modulus and Poisson's ratio of coal seam under plane strain conditions, respectively; u r , u θ are the radial displacement and the cyclic displacement, respectively; R 2 is the radius of the circular outer boundary set for the solution, R 2 R 0 ; K is a constant coefficient determined by the geometric size of the model. Its specific calculation formula can be written as Equation (3): The displacement field under the polar coordinates obtained by Equation (2) can be changed to the displacement field under the orthogonal Cartesian coordinate system by the coordinate transformation method. The transformed displacement component can be expressed as u x and u y , and the transformation process is shown in Equation (4).

Induced Stress Analysis
According to the research, we find that the internal filling material is not completely filled, and it appears to be intermittent along the orientation of the cleat direction. Therefore, it could be assumed that the filler is composed of multiple rod structures, which is like a spring, as shown in Figure 2. There is no connection between these rods. The strong effect of the rod structure is only connected with the upper side Γ 1 and the lower side Γ 2 on the cleat. When the normal displacement occurs on both sides of the cleat, the filler can be considered as a rod structure. When the cleat surface undergoes a tangent displacement, the filler can be considered as a beam structure. it appears to be intermittent along the orientation of the cleat direction. Therefore, it could be assumed that the filler is composed of multiple rod structures, which is like a spring, as shown in Figure 2. There is no connection between these rods. The strong effect of the rod structure is only connected with the upper side 1  and the lower side 2  on the cleat. When the normal displacement occurs on both sides of the cleat, the filler can be considered as a rod structure. When the cleat surface undergoes a tangent displacement, the filler can be considered as a beam structure. Therefore, under the external loading, the relationship between the filler and the normal load P and the tangent load Q generated on the cleat surface, and the relative displacement on both sides of the cleat surface, can be written as [24]: Therefore, under the external loading, the relationship between the filler and the normal load P and the tangent load Q generated on the cleat surface, and the relative displacement on both sides of the cleat surface, can be written as [24]: where, E F is the elastic modulus of the cleat filler; x , y are directions along the cleat and the vertical cleat, respectively; u x and u y are displacements in the x and y directions, respectively; d b is the width of the filler, which can be set as the width of the unit; due to the length of the cleat, along the cleat direction, the normal load P and the tangent load Q will change with different positions, so the cleat load and the cleat distribution load are the function of positions, seen in Figures 3 and 4.
Appl. Sci. 2020, 10, x FOR PEER REVIEW  5 of 18 where, EF is the elastic modulus of the cleat filler; x , y are directions along the cleat and the vertical cleat, respectively; x u  and y u  are displacements in the x and y directions, respectively; db is the width of the filler, which can be set as the width of the unit; due to the length of the cleat, along the cleat direction, the normal load P and the tangent load Q will change with different   According to the strength factor under a single normal load and a tangent load in the stress field intensity factor manual, and using the superposition principle, the strength factor under an asymmetric distribution load is obtained by integrating [25], seen in Equations (6)and (7):   According to the strength factor under a single normal load and a tangent load in the stress field intensity factor manual, and using the superposition principle, the strength factor under an asymmetric distribution load is obtained by integrating [25], seen in Equations (6)and (7): Figure 4. A tangential load model with asymmetric distribution for cleats. According to the strength factor under a single normal load and a tangent load in the stress field intensity factor manual, and using the superposition principle, the strength factor under an asymmetric distribution load is obtained by integrating [25], seen in Equations (6)and (7): Then the stress field expression near the crack tip can be written [25] in Equation (9): where, σ x , σ y and τ xy are the normal stress and the shearing stress in the additional stress field, respectively. The stress state is shown in Figure 5, where ρ, ϕ are the radius and angle of the point to the cleat tip, respectively.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 6 of 18 where, x   , y   and xy   are the normal stress and the shearing stress in the additional stress field, respectively. The stress state is shown in Figure 5, where ρ, φ are the radius and angle of the point to the cleat tip, respectively.
where G represents the shear elastic modulus; Thus, it can be obtained that the total displacement field of coal seam considering the cleats can be expressed as: Through the change of coordinates, the relationship between the displacement field and the total displacement field of coal rocks can be obtained: Cleat tip The induced displacement field distribution near the cutting tip can be written as [23]: where G represents the shear elastic modulus; κ = 3−µ 1+µ . Thus, it can be obtained that the total displacement field of coal seam considering the cleats can be expressed as: Through the change of coordinates, the relationship between the displacement field and the total displacement field of coal rocks can be obtained: Appl. Sci. 2020, 10, 1169 where φ is the cutting inclination. Through the simultaneous Equations (2), (5), (8), (10) and (11), the specific solution can be obtained. In this paper, we propose a program to gain the result of the above formula, because Equations (5)- (8) and (10) keep positive correlation with displacement u x , u y . So, we firstly set the values of u x and u y , and next through the above formula, we can get other values of u x and u y . Then taking the result of the previous iteration as the initial value of the continued iteration, we stop iterating until the error between the former and the latter is within the allowable range.

The Total Stress Field of the Broken Coal Seam Wellbore Surrounding Rocks
The Equations (1) and (9) can be obtained by the stress field distribution of a non-continuous coal seam containing a single shearing. Equations (1) and (9) are in different coordinate systems. For this purpose, the calculation results need to be unified into the same coordinate system for superposition. The coordinate system of the axis X and Y with the σ H and σ h directions is the final coordinate system. The stress in Equations (1) and (9) are all converted to this coordinate system, and the stress expression of the non-continuous coal seam is obtained by superposition.
By extension, if n cleats exist in the coal seam, the stress States of σ x , σ y and τ xy at any point in the coal seam can be expressed as: where σ i x , σ i y , and τ i x y are the induced stress fields caused by the i cutting, respectively; φ i is the direction of section i.

The Theory of Block Slip and the Determination of Movable Block
The rock in the engineering structure tends to be divided into different shapes by faults, joints, or weak interlayers. Because the strength of the joints and weak interlayers is weaker than the strength of the rock body, the damage of the engineering structure can be attributed to the sliding damage of rocks along the joint surface and the weak interlayer plane, and both the damage and the deformation of rocks are basically ignored.
The block slip theory assumes that the through structure surface is a spatial plane, the mosaic block body is regarded as a convex body, and the various loads of the structure surface are regarded as a spatial vector. Under the conditions of each spatial plane, the geometric method is used to pass simple static force calculation. Then, using the vector operation algorithm, the force of each plane of the mosaic block is summed up by vectors, and the sliding force of the mosaic block instability is obtained. Due to natural cutting, the coal rock is cut into finite rock bodies and infinite rock bodies. The finite rock blocks are divided into movable rock bodies and immovable rock bodies. The key to block theory is to determine the movable blocks.
Through the research, it is found that the parameters that determine the motility of the blocks mainly include the shear plane parameters and the geometric parameters of the near-empty plane. By moving the half-space interface of the constituent block to the coordinate zero point, a series of cones can be formed: the crack cone (CP) is a cone formed by the semi-space of the thermal fracture plane, and the excavation cone (EP) is a cone composed of barely the semi-space of the near-empty plane and the block cone (BP). Combining the block finiteness theorem and the mobility theorem, we can determine the movable block [26] as in Equation (14): Appl. Sci. 2020, 10, 1169 7 of 17 where ∅ represents an empty set.

Coal and Rock Cleat Failure Criterion
A large number of rock strength tests show that the strength of rock materials in high stress areas is hyperbolic. Therefore, we use the Extend Mogi-Coulomb (E.MG-C) criterion proposed by A.M. Ai-Ajm [27] as a criterion for the shear strength failure of dry, hot rock cracks. It can be expressed as: where τ oct is octahedral shear stress; σ 1 , σ 2 , σ 3 are the first, second and third main stresses respectively; σ m is the average main stress; a, b are parameters related to the cohesive force and internal friction angle of the rock, respectively, and c is the parameter of the nonlinear nature of the rock strength curve in high stress areas.
The criterion for tensile damage to cracks is shown as: where σ t is the positive tensile stress and [σ t ] is an allowable positive stress.

Block Sliding Direction
The block sliding direction is mainly determined by the structure of the block itself and the force of the block. Assuming that α is the angle between the side of the block and the direction → k , which is from the center of the block to the center of the wellbore, and α is marked as positive in a clockwise direction from the direction → k , α i bc is the angle between the boundary i and → k , and the sliding angle range determined by the block's own structure can be expressed as: where α sl is the sliding angle of the block; min α i bc , max α i bc are the minimum and maximum values of the α i bc composition set when it can slide alone along the boundary i. The equation for calculating the force on the cleat surface and the weak sandwich plane of the movable block is: where F ix bc and F iy bc are the combined forces along the x and y directions on the boundary I, respectively. x direction is along the boundary i of the block, and y direction is perpendicular to the boundary i of the block; n is the number of segments of the boundary i of the block. Since the weak surfaces such as cracks, cleat and weak interlayers, are generally very long, the stress is unevenly distributed along the thin and weak surfaces, and the calculation accuracy of the result by summation after segmentation calculation is higher in solid segments; σ j , τ j , l j are the positive stresses, shear stresses, and lengths of segment j above the boundary i respectively; d is for the thickness of weak surfaces such as cracks, joints and the weak mezzanine. The equation for calculating the force on the thermal fracture of a movable block can be shown as: where f is the coefficient of the friction on the surface of the hot fracture.
Combining the above two equations, the combined force on the boundary i of the block is: Assuming the x direction is along the direction → k , and the y direction is perpendicular to the direction → k , the two components of force → F i bc are calculated as: The resultant of forces → F bc on the entire block is: where N is the number of boundaries on the block. The size and direction of the joint force → F bc on the block can be expressed as: where F bc is the size of the resultant of forces; γ f is the angle between the resultant of forces and the direction → k . If γ f ∈ min α i bc , max α i bc , the block slides along the angle γ f in the direction

Case Analysis
Taking vertical well HG61-4-X in the Qinshui Basin (Shanxi) as an example, it is affected by the fault zone and the tectonic stress around. Therefore, the cleats grow extremely. During the drilling, there is a collapse in the range of 834.90-845.12 m in the coal seam. The numerical model of the coal seam blocks is established based on the slide theoretical model. As shown in Figure 6a, quadrilateral ABCD is the quadrilateral block. As shown in Figure 6b, triangle AED is the triquetrous block. Coal and rock mechanics parameters are shown in Table 1. The basic parameters are listed as follows: the depth is 841 m, the maximum horizontal principal stress σ H = 17.87 MPa, the minimum horizontal principal stress σ h = 16.68 MPa, the drilling fluid pressure P 0 = 8.3 MPa, the formation pore pressure is 8.15 MPa and the borehole diameter R 0 = 153.1 mm. seam blocks is established based on the slide theoretical model. As shown in Figure 6a, quadrilateral ABCD is the quadrilateral block. As shown in Figure 6b, triangle AED is the triquetrous block. Coal and rock mechanics parameters are shown in Table 1. The basic parameters are listed as follows: the depth is 841 m, the maximum horizontal principal stress σH = 17.87 MPa, the minimum horizontal principal stress σh = 16.68 MPa, the drilling fluid pressure P0 = 8.3 MPa, the formation pore pressure is 8.15 MPa and the borehole diameter R0 = 153.1 mm.

The Stress Field Comparative Analysis of Coal Seams Containing Multiple Cleats
In order to verify the accuracy of the coal seam stress field containing the multiple cleats analytical solution, the identical mechanical parameters are used, and the finite element model is established. Meanwhile, the principal stress map is obtained, as shown in Figure 7.

The Stress Field Comparative Analysis of Coal Seams Containing Multiple Cleats
In order to verify the accuracy of the coal seam stress field containing the multiple cleats analytical solution, the identical mechanical parameters are used, and the finite element model is established. Meanwhile, the principal stress map is obtained, as shown in Figure 7. As shown in Figure 7, the existence of cleats makes the stress yielded around the coal bed in the wellbore extraordinarily complicated. According to the calculation consequence, the stress value of any point in the coal seam is vitally impacted by the cleats around. The following assumption is made, that is, the stress values are completely different among two points which are closed in the coal seam. To verify the accuracy of the analytical solution, the finite element results in the wellbore should be As shown in Figure 7, the existence of cleats makes the stress yielded around the coal bed in the wellbore extraordinarily complicated. According to the calculation consequence, the stress value of any point in the coal seam is vitally impacted by the cleats around. The following assumption is made, that is, the stress values are completely different among two points which are closed in the coal seam. To verify the accuracy of the analytical solution, the finite element results in the wellbore should be taken into account. Additionally, after calculation of the consequence at the counterpart position, the curve comparison is shown in Figure 8. According to Figure 8, the maximum and minimum principal stresses calculated by the theory method are slightly larger than those calculated by the finite element method, and the relative error of those methods ranges from 0.64% to 9.16%. When the round angle ranges from 120° to 300°, the error is minimum. Based on these analyses, the proposed analytic model of the coal bed stress yield is reliable.

The Wellborn Stability Analysis of the Coal Seam
If SH 1   , the cleats present strength failure, but on the contrary, the cleats do not present strength failure. If EQ 1   , the load balance of the block along the slide is broken, which means that the block begins sliding. Otherwise, the load balance of the block along its sliding is maintained. For the quadrilateral block, the only way to slide is by overcoming the strength in both sides along the sliding, and breaking the limit equilibrium in the side which is perpendicular to the sliding direction. So, there are three conditions required to judge the falling of the quadrilateral block. For the triquetrous block, the only way to slide is in overcoming the strength fracture on one side and breaking the limit equilibrium on the other side. Therefore, there are two conditions required to judge the falling of the triquetrous block. According to Figure 8, the maximum and minimum principal stresses calculated by the theory method are slightly larger than those calculated by the finite element method, and the relative error of those methods ranges from 0.64% to 9.16%. When the round angle ranges from 120 • to 300 • , the error is minimum. Based on these analyses, the proposed analytic model of the coal bed stress yield is reliable.

The Wellborn Stability Analysis of the Coal Seam
η SH , η EQ are defined as block strength coefficient and equilibrium slide coefficient, respectively, in this paper.
If η SH ≥ 1, the cleats present strength failure, but on the contrary, the cleats do not present strength failure. If η EQ ≥ 1, the load balance of the block along the slide is broken, which means that the block begins sliding. Otherwise, the load balance of the block along its sliding is maintained. For the quadrilateral block, the only way to slide is by overcoming the strength in both sides along the sliding, and breaking the limit equilibrium in the side which is perpendicular to the sliding direction. So, there are three conditions required to judge the falling of the quadrilateral block. For the triquetrous block, the only way to slide is in overcoming the strength fracture on one side and breaking the limit equilibrium on the other side. Therefore, there are two conditions required to judge the falling of the triquetrous block.
In order to verify the accuracy of the analytical model, the numerical model is established by PFC software. The boundary condition (B.C.) of the numerical model is the same as that of the analytical model. The PFC numerical model is presented in Figure 9, and the red lines represent the face cleats, while the blue lines represent the butt cleats. As seen from Figure 10, the PFC has calculated the sliding result of the block around the wellbore. The calculated result shows that the sliding coefficients of the analytical model prediction are listed as follows: the quadrilateral block, . Since all sliding coefficients are greater than 1, it indicates that the block will fall off. The prediction made by the discrete element method is in accord with the theoretical prediction. To further verify the accuracy of the analytical model, the minimum drilling fluid pressure of maintaining the stability of the wellbore will be calculated spontaneously by the analytical model and the discrete element model. From the results, for the quadrilateral block, the minimum drilling fluid pressure by the analytical mode is 9.4 MPa, and that by the discrete element model is 9.7 MPa and the gap is 3.19%; for the triquetrous block, the minimum drilling fluid pressure by the analytical model is 11.1 MPa, and that by the discrete element is 11.5 MPa, and the gap is 2.68%. Those results indicate that the analytical model and the discrete element model keep highly identical.

The Effect Factor Analysis
The stability of the coal seam wellbore is impacted greatly by the geometry and the position of the block. Therefore, this paper focuses on discussing how the face cleat spacing, the butt cleat spacing, the cleat length, the cleat inclination angle and the block geometric position have an effect upon the strength sliding coefficient and the equilibrium slip coefficient of the quadrilateral block and the triquetrous block. The calculation parameters and the analysis parameters can be referred to Figure 6. As seen from Figure 10, the PFC has calculated the sliding result of the block around the wellbore. The calculated result shows that the sliding coefficients of the analytical model prediction are listed as follows: the quadrilateral block, η SH 1 = 1.58, η SH 2 = 3.31, η EQ = 69.79; the triquetrous block, η SH = 2.49, η EQ = 1.43. Since all sliding coefficients are greater than 1, it indicates that the block will fall off. The prediction made by the discrete element method is in accord with the theoretical prediction. To further verify the accuracy of the analytical model, the minimum drilling fluid pressure of maintaining the stability of the wellbore will be calculated spontaneously by the analytical model and the discrete element model. From the results, for the quadrilateral block, the minimum drilling fluid pressure by the analytical mode is 9.4 MPa, and that by the discrete element model is 9.7 MPa and the gap is 3.19%; for the triquetrous block, the minimum drilling fluid pressure by the analytical model is 11.1 MPa, and that by the discrete element is 11.5 MPa, and the gap is 2.68%. Those results indicate that the analytical model and the discrete element model keep highly identical. As seen from Figure 10, the PFC has calculated the sliding result of the block around the wellbore. The calculated result shows that the sliding coefficients of the analytical model prediction are listed as follows: the quadrilateral block, . Since all sliding coefficients are greater than 1, it indicates that the block will fall off. The prediction made by the discrete element method is in accord with the theoretical prediction. To further verify the accuracy of the analytical model, the minimum drilling fluid pressure of maintaining the stability of the wellbore will be calculated spontaneously by the analytical model and the discrete element model. From the results, for the quadrilateral block, the minimum drilling fluid pressure by the analytical mode is 9.4 MPa, and that by the discrete element model is 9.7 MPa and the gap is 3.19%; for the triquetrous block, the minimum drilling fluid pressure by the analytical model is 11.1 MPa, and that by the discrete element is 11.5 MPa, and the gap is 2.68%. Those results indicate that the analytical model and the discrete element model keep highly identical.

The Effect Factor Analysis
The stability of the coal seam wellbore is impacted greatly by the geometry and the position of the block. Therefore, this paper focuses on discussing how the face cleat spacing, the butt cleat spacing, the cleat length, the cleat inclination angle and the block geometric position have an effect upon the strength sliding coefficient and the equilibrium slip coefficient of the quadrilateral block and the triquetrous block. The calculation parameters and the analysis parameters can be referred to Figure 6.

The Effect Factor Analysis
The stability of the coal seam wellbore is impacted greatly by the geometry and the position of the block. Therefore, this paper focuses on discussing how the face cleat spacing, the butt cleat spacing, the cleat length, the cleat inclination angle and the block geometric position have an effect upon the strength sliding coefficient and the equilibrium slip coefficient of the quadrilateral block and the triquetrous block. The calculation parameters and the analysis parameters can be referred to Figure 6.

Face Cleat Spacing
In order to analyze the influence of the face cleat spacing for block sliding coefficient, the face cleat spacing between SL 1 and SL 2 is modified, keeping the other parameters invariant. The face cleat spacing ranges from 2.09 cm to 3.19 cm.
According to Figure 11, the addition of the face cleat distance will cause the decrease of the strength coefficient and equilibrium slip coefficient of the quadrilateral block, and the decrease of η SH2 is greater than η SH1 . When the distance is above 3.13 cm, the quadrilateral block stops falling off. On the contrary, the addition of the face cleat distance will cause the increase of the strength coefficient and the equilibrium slip coefficient of the triangular block, and the increase of η SH2 is greater than that of η SH1 . When the distance is above 2.87 cm, the triangular block starts to drop. The increase of the face cleat distance decreases the interference among the cleats. Additionally, it reduces the stresses of the lines AB and CD on the face cleats, and it causes the reduction of the strength slip coefficient. The difference of the direction and the length of the face cleats SL 1 and SL 2 are both the reason why the ranges of the discount are different between η SH2 and η SH1 . The increase of the face cleat distance causes the addition of the resistance to the triquetrous block sliding, and it reduces the un-equilibrium force in blocks, so the equilibrium slip coefficient decreases. For the triquetrous block, the shape becomes blunt and wide from a cuspidal, and the face cleat distance becomes narrow. Meanwhile, the angle of cleat SL 4 decreases. Based on Equations (4)- (6), with the increase of shearing stress and tensile stress on the cleat SL 4 , the stresses on cleats AE, DE, the strength coefficient and the equilibrium slip coefficient, all increase.

Face Cleat Spacing
In order to analyze the influence of the face cleat spacing for block sliding coefficient, the face cleat spacing between SL1 and SL2 is modified, keeping the other parameters invariant. The face cleat spacing ranges from 2.09 cm to 3.19 cm. According to Figure 11, the addition of the face cleat distance will cause the decrease of the strength coefficient and equilibrium slip coefficient of the quadrilateral block, and the decrease of SH2  is greater than SH1  . When the distance is above 3.13 cm, the quadrilateral block stops falling off. On the contrary, the addition of the face cleat distance will cause the increase of the strength coefficient and the equilibrium slip coefficient of the triangular block, and the increase of

Butt Cleat Distance
In order to analyze how the butt cleat distance makes an impact on the block slip coefficient, the length of line AF is modified, keeping other factors invariant. The length of the line AF ranges from 4 cm to 14 cm.
According to Figure 12, with the increase of the butt cleat distance, the strength slip coefficient of the quadrilateral block increases, and the equilibrium slip coefficient of the quadrilateral block climbs up and then declines. But for the triquetrous block, the strength slip coefficient and the equilibrium slip coefficient decline. When the butt cleat distance is longer than 5 cm, the triquetrous block stops falling off. The reason is that under the local coordinate system of cleat SL3, the dimensionless radius among the midpoints of cleats AB and CD increases as the addition of the butt cleat distance. In addition, this causes the increase of the stress on cleats AB and CD, and at the same time, it affects the increase of the strength slip coefficient. It is hard to conclude the change law of the equilibrium slip coefficient, because the change of the butt cleat distance causes the stress changes on all the sides of the quadrilateral block. Therefore, the calculation is required to define these values.

Butt Cleat Distance
In order to analyze how the butt cleat distance makes an impact on the block slip coefficient, the length of line AF is modified, keeping other factors invariant. The length of the line AF ranges from 4 cm to 14 cm.
According to Figure 12, with the increase of the butt cleat distance, the strength slip coefficient of the quadrilateral block increases, and the equilibrium slip coefficient of the quadrilateral block climbs up and then declines. But for the triquetrous block, the strength slip coefficient and the equilibrium slip coefficient decline. When the butt cleat distance is longer than 5 cm, the triquetrous block stops falling off. The reason is that under the local coordinate system of cleat SL 3 , the dimensionless radius among the midpoints of cleats AB and CD increases as the addition of the butt cleat distance. In addition, this causes the increase of the stress on cleats AB and CD, and at the same time, it affects the increase of the strength slip coefficient. It is hard to conclude the change law of the equilibrium slip coefficient, because the change of the butt cleat distance causes the stress changes on all the sides of the quadrilateral block. Therefore, the calculation is required to define these values. Because of the addition of the butt cleat distance, the body areas of the quadrilateral block and the triquetrous block increase, as well as the angle of cleat SL 4 . By Equations (4) and (5), we can know that it causes the decline of the shearing stresses and the tensile stresses on cleat SL 4 , as well as stresses on cleats AE and DE. As a result, the strength coefficient and equilibrium slip coefficient decline.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 14 of 18 Because of the addition of the butt cleat distance, the body areas of the quadrilateral block and the triquetrous block increase, as well as the angle of cleat SL4. By Equations (4) and (5), we can know that it causes the decline of the shearing stresses and the tensile stresses on cleat SL4, as well as stresses on cleats AE and DE. As a result, the strength coefficient and equilibrium slip coefficient decline.

Cleat Length
Keeping other parameters invariant, we change cleat SL1 length and analyze the effect of the cleat length on the sliding coefficient of the block. SL1 ranges from 15 cm to 24 cm.
From Figure 13, it can be seen that the increase in the length of the cleat SL1 reduces the strength and equilibrium sliding coefficient. When the length of the cleat SL1 is greater than 21 cm, the quadrilateral block stops falling off, and the triquetrous block never falls off. The increase in the length of the cleat reduces the dimensionless radius of the midpoints of the cleat, thereby reducing the stress of the midpoints of the cleat, so the strength slip coefficient decreases. However, the internal pressure of the well wall has not changed. From the limit equilibrium equation, it is known that the decrease of the block slide force is larger than that of the resistant force, so the uneven force of the block decreases and then the equilibrium slip coefficient decreases.

Cleat Length
Keeping other parameters invariant, we change cleat SL 1 length and analyze the effect of the cleat length on the sliding coefficient of the block. SL 1 ranges from 15 cm to 24 cm.
From Figure 13, it can be seen that the increase in the length of the cleat SL 1 reduces the strength and equilibrium sliding coefficient. When the length of the cleat SL 1 is greater than 21 cm, the quadrilateral block stops falling off, and the triquetrous block never falls off. The increase in the length of the cleat reduces the dimensionless radius of the midpoints of the cleat, thereby reducing the stress of the midpoints of the cleat, so the strength slip coefficient decreases. However, the internal pressure of the well wall has not changed. From the limit equilibrium equation, it is known that the decrease of the block slide force is larger than that of the resistant force, so the uneven force of the block decreases and then the equilibrium slip coefficient decreases.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 14 of 18 Because of the addition of the butt cleat distance, the body areas of the quadrilateral block and the triquetrous block increase, as well as the angle of cleat SL4. By Equations (4) and (5), we can know that it causes the decline of the shearing stresses and the tensile stresses on cleat SL4, as well as stresses on cleats AE and DE. As a result, the strength coefficient and equilibrium slip coefficient decline.

Cleat Length
Keeping other parameters invariant, we change cleat SL1 length and analyze the effect of the cleat length on the sliding coefficient of the block. SL1 ranges from 15 cm to 24 cm.
From Figure 13, it can be seen that the increase in the length of the cleat SL1 reduces the strength and equilibrium sliding coefficient. When the length of the cleat SL1 is greater than 21 cm, the quadrilateral block stops falling off, and the triquetrous block never falls off. The increase in the length of the cleat reduces the dimensionless radius of the midpoints of the cleat, thereby reducing the stress of the midpoints of the cleat, so the strength slip coefficient decreases. However, the internal pressure of the well wall has not changed. From the limit equilibrium equation, it is known that the decrease of the block slide force is larger than that of the resistant force, so the uneven force of the block decreases and then the equilibrium slip coefficient decreases.

Cleat Inclination
Keeping other parameters invariant, we change the inclination α SL 2 of the cleat SL 2 and analyze the effect of the cleat inclination on the sliding coefficient of the block. Inclination α SL 2 ranges from 15 • to 21 • . From Figure 14, it can be seen that with the increase of the inclination of the cleat, the strength slip coefficient of the quadrilateral block increases, the equilibrium slip coefficient decreases, and the strength and equilibrium slip coefficient of the triangular block increases. When α SL 2 ≥ 20 • , the triquetrous block changes start falling off. As α SL 2 increases, from Equations (4)- (6), it can be seen that the shear stress and the tensile stress on the cleat SL 2 increase, and this results in an increase in the stress on the cleat AB and CD, and an increase in the strength slip coefficient of the quadrilateral block. Similarly, the strength slip coefficient of the triquetrous block increases as α SL 2 increases. The increase of α SL 2 shortens the length of the cleat BC, and reduces the block slide force, so the equilibrium slip coefficient of the quadrilateral block is reduced; and for the triquetrous block, the effect of the α SL 2 increase on the strength and equilibrium slip coefficient is the same as the increase in the distance between face cleats (see Section 5.3.1).

Cleat Inclination
Keeping other parameters invariant, we change the inclination 2 SL  of the cleat SL2 and analyze the effect of the cleat inclination on the sliding coefficient of the block. Inclination 2 SL  ranges from 15° to 21°.
From Figure 14, it can be seen that with the increase of the inclination of the cleat, the strength slip coefficient of the quadrilateral block increases, the equilibrium slip coefficient decreases, and the strength and equilibrium slip coefficient of the triangular block increases. When

Block Geometric Position
Keeping other parameters invariant, we change α 1 and analyze the effect of the block geometry position on the strength and equilibrium slip coefficient. α 1 ranges from 0 • to 180 • .
From Figure 15a, we can see that with the increase of α 1 , the strength slip coefficient and the equilibrium slip coefficient of the quadrilateral block are approximately sinusoidal changes. When α 1 changes from 30 • to 120 • , the quadrilateral block starts falling off. From Figure 15b, we can also observe that with the increase of α 1 , the equilibrium slip coefficient of the triquetrous block is approximately sinusoidal, but the strength slip coefficient is undulating. When α 1 changes from −20 • to 10 • , and from 60 • to 120 • , the triquetrous block starts falling off. The change of α 1 causes the change of the cleat angle around the block, which causes the change in the induced stress field and at last causes a fluctuation in the strength and equilibrium slip coefficient.

Discussion
This paper presents a novel stress field analytical model of the wellbore coal rock in the broken coal seam. This analytical model considers the influence of the cleat filler, and by simplifying the cleat filler as the multiple rod structures (see Figure 2), we use the mechanics of materials and fracture mechanics to establish the induced stress field and the displacement field.
According to the limit equilibrium theory and the E.MG-C criterion, we further establish the stable analysis of the wellbore coal rock. The steps are: Firstly, we determine the motility of the block mainly including the shear plane parameters and the geometric parameters of the near-empty plane; secondly, we gain coal and rock cleat failure criterion; and lastly, we gain the block sliding direction. Taking vertical well HG61-4-X in the Qinshui Basin as an example, this paper uses a numerical method (see Figures 7 and 9) to judge the validity of the analytical model. This paper also analyzes the effect factors such as the face cleat spacing, the butt cleat spacing, the cleat length, the cleat inclination angle and the block geometric position. The results show that a quadrilateral block slides off easier than a triangular block under the same boundary condition; the bigger are the cleat spacing and cleat length, the lower is the risk at which blocks slide off, and the increasing cleat angle could cause blocks to slide off easily. Under the same boundary condition, it is closely related to the well round angle at which blocks slide off.
In addition, we merely consider the drilling fluid fluctuation in this paper. Thus, an improved and more general approach should be proposed in the future work to completely resolve the issue.

Conclusions
This paper establishes a stress field analytical model of the wellbore coal rock by considering the irregularity of the cleat distribution and the influence of the cleat filler. According to the limit equilibrium theory and the E.MG-C criterion, the wellbore stability model of the coal seam is established by determining, firstly, the motility of the block, secondly, the cleat failure criterion and thirdly, the block sliding direction. In order to judge the validity of the analytical model, this paper uses a numerical method to make a comparative analysis. At the same time, cleat affecting factors are studied in the paper.

Discussion
This paper presents a novel stress field analytical model of the wellbore coal rock in the broken coal seam. This analytical model considers the influence of the cleat filler, and by simplifying the cleat filler as the multiple rod structures (see Figure 2), we use the mechanics of materials and fracture mechanics to establish the induced stress field and the displacement field.
According to the limit equilibrium theory and the E.MG-C criterion, we further establish the stable analysis of the wellbore coal rock. The steps are: Firstly, we determine the motility of the block mainly including the shear plane parameters and the geometric parameters of the near-empty plane; secondly, we gain coal and rock cleat failure criterion; and lastly, we gain the block sliding direction. Taking vertical well HG61-4-X in the Qinshui Basin as an example, this paper uses a numerical method (see Figures 7 and 9) to judge the validity of the analytical model. This paper also analyzes the effect factors such as the face cleat spacing, the butt cleat spacing, the cleat length, the cleat inclination angle and the block geometric position. The results show that a quadrilateral block slides off easier than a triangular block under the same boundary condition; the bigger are the cleat spacing and cleat length, the lower is the risk at which blocks slide off, and the increasing cleat angle could cause blocks to slide off easily. Under the same boundary condition, it is closely related to the well round angle at which blocks slide off.
In addition, we merely consider the drilling fluid fluctuation in this paper. Thus, an improved and more general approach should be proposed in the future work to completely resolve the issue.

Conclusions
This paper establishes a stress field analytical model of the wellbore coal rock by considering the irregularity of the cleat distribution and the influence of the cleat filler. According to the limit equilibrium theory and the E.MG-C criterion, the wellbore stability model of the coal seam is established by determining, firstly, the motility of the block, secondly, the cleat failure criterion and thirdly, the block sliding direction. In order to judge the validity of the analytical model, this paper uses a numerical method to make a comparative analysis. At the same time, cleat affecting factors are studied in the paper.