Meshless method simulation and experimental investigation of crack propagation of CBM hydraulic fracturing

For simulating CoalBed Methane (CBM) hydraulic fracturing using 3-D meshless method, this paper analyzed the hydraulic fracturing mechanism and cracking form for coal rock and established the geometric and mathematical models of hydraulic fracturing propagation in coal rock in terms of the Hillerborg model on crack opening displacement theory. With the theoretical basis of hydromechanics, the formulas for calculating hydraulic pressure inside the fracture by numerical simulation were deduced from the analysis on this fluid-structure interaction problem. The geometric and mathematical models established above were described by 3-D meshless Galerkin (EFG, Element-Free Galerkin) method and compiled into the numerical simulation program using VB and FORTRAN programming language to simulate the fracture propagation for an actual coal rock sample with a drilling hole as an example. Then the physical simulation experiment of hydraulic fracturing propagation of coal seam was conducted on the same coal rock sample. Through the direct observation with naked eyes and detection by advanced instruments of ESEM and Micro-CT, the shape and parameters of cracks on the surface of and inside the coal rock sample were achieved, which indicated that experimental results are reasonably consistent with numerical simulation results.


Introduction
In recent decades CoalBed Methane (CBM) has become an important source of energy all over the world. Hydraulic fracturing is a popular and powerful technology mainly used in the CBM industry to stimulate reservoirs for improving the productivity. Although the productive benefit from hydraulically fractured wells is potentially large, it is hard to predict the fracture propagation and perform an optimal hydraulic fracturing process. The special structures and properties of coal seam (such as natural cleat, fracture growing, lower elastic modulus, higher Poisson ratio, frangibility, and strong heterogeneity), together with drilling disturbance and fluid flow make it difficult to figure out the complex and diverse 3-D (three dimensional) fracture propagation using mathematical analysis method without the help of computer.
Aimed to the principle and regulation of fracture propagation of oil and gas recovery for an optimal hydraulic fracturing process, a number of studies are carried out recently by computer numerical simulation or physical experimental simulation. Finite element methods are usually employed to the numerical investigation of fracture shape and propagation laws on hydraulic fracturing crack of various kinds of rock mass. According to the relevant literatures, FRANC-3D, HYFRANC-3D, or ABAQUS [1][2][3][4] are popular numerical simulation software based on the finite-element methods for simulating the propagation of hydraulic fractures.
These studies are instrumental to understand the complex fracture growth to some extent. Results from Rahman show that two or more multiple transverse fractures may propagate but they involve different treatment pressure behaviors compared to a single fracture. Mofazzal Hossain and M.K. Rahman developed insights to explain how the complex fracture geometry grows under different stress and well conditions, and their effects on the injection pressure. Zuorong Chen developed the cohesive zone model to avoid the singularity in the crack tip stress field that is present in classic fracture mechanics. Yao Yao discussed the fundamental mechanical difference of line elastic fracture mechanics and cohesive fracture mechanics-based methods in according to the size of the fracture process zone and its effect on crack extension in ductile rock.
As a matter of fact, however, it was just an expedient way to simulate the fracture propagation of rock cracking, especially resulting from hydraulic fracture inside CBM wellbore, using finite element methods. As we well known, finite element methods were originally defined on meshes of data points. In such a mesh, each point has a fixed number of predefined neighbors, and this connectivity between neighbors can be used to define mathematical operators like the derivative for constructing the equations to simulate. But in simulations where the material being simulated can move around (as in computational fluid dynamics in hydraulic fracturing) or where large deformations of the material can occur (as in simulations of plastic or linear elasto materials of coal rock), the connectivity of the mesh can be difficult to maintain without introducing error into the simulation. If the mesh becomes tangled or degenerate during simulation, the operators defined on it may no longer give correct values. The mesh may be recreated during simulation (called remeshing), but this can also introduce error, since all the existing data points must be mapped onto a new and different set of data points.
As a much better alternative solution for simulating fracture propagation resulting from hydraulic fracturing, meshfree methods (also called element free methods), may remedy above mentioned shortfalls using finite element methods [5][6][7][8] since they are very suitable for: (1) simulations where creating a useful mesh from the geometry of a complex 3D object may be especially difficult or require human assistance, (2) simulations where nodes may be created or destroyed, such as in cracking simulations, (3) simulations where the problem geometry may move out of alignment with a fixed mesh, such as in bending simulations, and (4) simulations containing nonlinear material behavior, discontinuities or singularities.
Based on the use of moving least-squares interpolants with a Galerkin method, the Element Free Galerkin (EFG) method is particularly effective in line elastic fracture problems with the advantages of the circumvention of the need for element data, no any volumetric locking, the rate of convergence much faster than that of finite elements significantly, and a high resolution of localized steep gradients [9]. The implementation of the EFG method for problems of fracture and static crack growth based on an infinite plate with a central circular hole subjected to a unidirectional tensile load in the x direction shows that accurate stress intensity factors can be obtained without any enrichment of the displacement field by a near-crack-tip singularity and that crack propagation can be easily modeled since it requires hardly any remeshing [10,11]. More simulating examples for arbitrary three-dimensional propagating cracks dynamically in elastic bodies also further verify that EFG method is particularly attractive in dynamic crack propagation under various tensile loading cases, such as quasi-static loading and step load, and even the combination of tension and torsion [12][13][14][15].
An expansive application on a crack within an infinite square plate using element free method attempted to simulate the hydraulic fracturing of rock mass under the external unidirectional compression loading without or with water pressure inside the rock mass [16]. Although this kind of condition seems like but totally not equal to the hydraulic fracturing in the wellbore for oil and gas recovery, the relevant results still verify the meshless method can commendably deal with the problems of rock fracturing under hydraulic pressure. In addition, a coupled hydro-mechanical analysis for prediction of hydraulic fracture propagation in saturated porous media using EFG meshless method successfully simulates the leak-off of fluid from the fracture into the surrounding homogeneous and heterogeneous saturated soils [17].
From the analysis above, all of these significant researches indicate that all properties of EFG method seem to be customized for the requirements for simulating dynamic fracture propagation of hydraulic fracturing in coal seam for CBM recovery. However, these researches are hardly involved in directly applying EFG method to analyze fracture initiation and propagation of CBM hydraulic fracturing. Furthermore, up to now, many research scholars are still accustomed to numerical modeling of hydraulic fracturing in rock mass by finite element method [18,19].
The main objective of this paper is to investigate the feasibility and validity of the EFG method on the crack propagation of CBM hydraulic fracturing by geometric and mathematic modeling after theoretic analysis on the cracks of an actual coal rock sample with a drilling hole. Then the physical simulation experiments are conducted on the same coal rock sample. Through the direct observation with naked eyes or detection by advanced instruments on the shape and parameters of cracks, the actual physical process of fracture propagation is scientifically described to verify and improve theoretical model and numerical model. The research results will be beneficial to the design of hydraulic fracturing in coal seam.
2 Geometric and mathematic modeling 2.1 Judging criteria for coal rock cracking Although coal rock belongs to fragile material, the tip of existing crack inside the wellbore under hydraulic pressure will not be in the form of simple single brittle fracture, but rather in the form of the combination of brittle-elastic-plastic due to the crustal stress from coal reservoirs, especially buried deeply [20,21]. This process can be considered as three procedures: (1) the fracture tip produces brittle fracture while it is subjected to the maximum fracture strength resulting from hydraulic pressure, (2) due to the constraint of in-situ stress, the fracture will propagate under the continuous action of hydraulic pressure till the pressure is unloaded or fracture propagation is complete, and (3) the fracture occurred will be partially closed after pressure unloading, thus forming the plastic fracture. This also can be vividly called a ''spring-splint'' model illustrated in Figure 1. K 1 and K 2 are the equivalent spring coefficients of the coal rock between roof and floor, w is the crack width, and a is the opening angle.
Within this model, the mechanic effect for a fracture propagation of coal rock can be stated as that the extended fractures are subjected to a closing force r from a spring drive fractures. In this case, when the stability coefficient of tip f = r/w exceeds the recovery coefficient K of spring, plastic fractures will be completed finally. Hence the fracture propagation laws based on plastic criteria can be used as the judging criteria in theory for hydraulic fracture propagation laws.
The Hillerborg model based on Crack Opening Displacement (COD) theory has been proved that it is consistent with the mechanism of fracture propagation and quite suitable to the numerical simulation analysis for it [22][23][24]. Therefore the mechanic model based on the Hillerborg model is built to simulate the hydraulic fracture propagation of coal rock using meshless method. According to the spring-splint model, for a certain kind of coal rock, the critical opening displacement resulting from the brittle fracture or elasto plastic deformation can be derived from Equations (1) and (2) respectively.
Since the brittle fracture is dominated by the plane stress, the critical opening displacement is Similarly, since the elasto plastic deformation is dominated by the plane strain, the critical opening displacement is where K IC is the fracture toughness, f t is the tensile strength, E is the elastic modulus, and m is the Poisson ratio of coal petrography. Then with the COD criterion of w ! w c , it can be judged that whether the fracture propagation directions follow the original fracture surfaces or not.

Modeling hydraulic fracturing shape of coal rock
The single fracture propagation model is the basis of analyzing and calculating complex fractures. There are three basic modes of crack extension, i.e., opening mode (Mode I), shearing mode (Mode II) and tearing mode (Mode III). Good agreement is found between the analytically and numerically estimated Mode I for the single fracture scenario [1]. That is the case the crack of coal rock by hydraulic fracturing is mainly Mode I.
Considering the differences on the natural properties of between the coal roof and coal floor resulting in high aspect ratio hydraulic fracture in coal seam, the constant-height rectangular section (KGD mode) can be employed to study fracture propagation based on fracture Mode I of coal rock illustrated in Figure 2. The so-called KGD model (plane strain) was independently developed Geertsma and de Klerk at the same time the PKN model was formulated by Nordgren based on the PK model developed by Perkins and Kern. Variations of the KGD and PKN models were used routinely for treatment designs as recently as the 1990s, and are still used today for computer simulation of hydraulic fractures. The PKN model is applicable to long fractures of limited height and elliptical vertical crosssection, whereas the KGD model for width calculation is height independent, and is used for short fractures where plane strain assumptions are applicable to horizontal sections [25,26].

Computing hydraulic pressure inside the fracture
For a specific coal seam, the fracture initiation and propagation law can be expressed by the functional relationship of the surrounding ground stress and the water pressure in coal seam since it belongs to typical fluid-structure interaction problem. As far as KGD hydraulic fracture model, the computing formula of pressure distribution of the water flow inside the fracture can be deduced using mass conservation and momentum conservation theorems.

Water flow motion equation inside the fracture
With regard to the water flow problem of high aspect ratio hydraulic fracture in coal seam, the control volume method is usually an efficient way to figure out the solution [27].
In terms of the mass conservation theorem, there is where m is the mass of control volume, q is the mass density of control volume, v is the volume of control volume, X is the surface of control volume, u is the velocity vector of control volume, and n is outer normal vector of control volume. Also, in terms of the momentum conservation theorem, there is:  where M is the momentum of control volume and P F is the vector sum of forces on the control volume.

Modeling theoretic hydraulic pressure distribution inside hydraulic fracture
To simplify the computing model rationally, some recognized rules inside coal petrography fractures should be adopted on the premise of no excessive effect on computational accuracy as follows: (1) The fluid volume should be incompressible; (2) The fluid belongs to Newtonian fluid; (3) The fluid loss should be ignored; and (4) The water flow should obey the rules of onedimensional (1-D) laminar flow.
The velocity distribution of 1-D laminar flow conforms to Poiseuille flow at time t in any section illustrated in Figure 3. There is: where u x,y is the flow velocity at any point inside the section x, u 0 is the maximum fluid velocity in the section x, and w x is the fracture width of the section x.
As the fracture propagation shape of coal rock is chosen to be the KGD model, the longitudinal profile of the fractures section should be the constant-height rectangular section while the transverse profile of the fractures section should be prolate ellipse expressed by expression (6) at time t is where a and b are the long axis and minor axis of fracture ellipse section at time t respectively. By taking an infinitesimal segment dx for the control volume x into the mass conservation equation and momentum conservation equation of fluid, the differential equation of hydraulic pressure inside the fracture can be obtained after ignoring higher-order differentials illustrated in Equation (7): where w x is the fracture width of section x, l is the kinetic viscosity coefficient of fluid, q x and p x is the quantity of flow and the average head pressure of section x respectively. Substituting Equation (6) into Equation (7) gives where w 0 is the fracture width, w 0 0 is the change rate of fracture width with time, and a 0 is the change rate of fracture length with time.
For the brittle-elastic-plastic fracture of coal rock, the change rate of fracture length with time can be set as zero, i.e., a 0 = 0. Then by introducing the boundary condition, i.e., i x | x=0 = p 0 (p 0 is the water injection pressure), Equation (8) can be simplified to Equation (9) that statements the distribution of hydraulic pressure inside the fracture of coal rock: where m is the kinematic viscosity coefficient of fluid and 1/m = q/l.

The 2-D processing method for 3-D fractures
For success solution of mechanics and mathematical models built above, coal rock is sliced along its bedding direction and then projected so that the coal seam with a certain thickness can be disposed as a plane approximately illustrated in Figure 4. This can be called 3-D to 2-D process in which complex 3-D fracture problem are transformed to easier 2-D plane fracture problem. Before the 3-D to 2-D process, the following operations should be considered: 1. In order to prevent the nodes inside the solution domain overflow the solution domain due to the large deformation caused by crack propagation during hydraulic fracturing, the circumcircle of rectangular section in bedding direction is served as the boundary for background integral mesh. Since the circumcircle is theoretically the minimum regular boundary including the rectangular solution domain, any other envelope curves will be unreflected. 2. As the thickness of selected bedding slice along bedding direction is far less than that of coal reservoir, i.e., m << c, it can be regarded as an actual plane. Furthermore, it is assumed that the surface of bedding section is rather smooth without obvious uplift and fault, thus the vertical ground stress can be set as zero. As shown in Figures 5 and 6, the model is very explicit after the 3-D to 2-D process on both mechanics and geometric sides.
According to the principle of geometric transformation, the geometric transformation of 3-D data points can be expressed as The solution domain is The boundary condition is The transformation matrix of vertical view is At the step i, the fracture propagation length will be As the fracture initiation always starts at the wall of wellbore, the initial solution domain is In Equations (12), (13), (15) and (16), x i is x-coordinate of point at the solution domain, y i is y-coordinate of point, z i is z-coordinate of point, R is the constraint boundary radius, r is the radius of water injection hole, R i is the curvature radius at any node of the solution domain, b i is the dip angle of sought point, a i is the azimuth angle of sought point, x 0 is the x-coordinate of point of fracture initiation, y 0 is the y-coordinate of point of fracture initiation, z 0 is the z-coordinate of point of fracture initiation, b 0 is the dip angle of point of fracture initiation, a 0 is the azimuth angle of point of fracture initiation, i is the step number for i infinitesimal segment dx, and l i is the fracture length at i step.
After fracture propagation is complete or the step is ended, the 2-D fracture graph in bedding slice section will be transformed back to the 3-D fracture propagation model by multiplying an inverse matrix " A expressed in Equations (17) and (18) 3 Program design for numerical simulation Since many studies have conducted on the theory of EFG method and its effectiveness in fracture problems, this section will not repeat them but only simply indicate the program flow of combined programming in FORTRAN and Visual Basic programming languages based on that.

The total program flow
The total program flowchart is illustrated in Figure 7. Visual Basic is used for the interactive window of input and output (shown in Fig. 13) while FORTRAN is for the fracture propagation program based on EFG method in the dashed box running in the background.

Program based on EFG method
The whole numerical simulation program of hydraulic fracturing propagation based on meshless method (including EFG method) is normally divided into three parts shown in Figure 8: the meshless method executive process, meshless method program module, and external program module. Both meshless method program module and external program module are developed in terms of the execution process of the program. Meshless program module is used to describe the solving domain mathematically based on meshless method to establish meshless method basic equations of hydraulic fracture propagation for coal rock sample. It consists of general submodules: InPut for the external data input, Gauss Coefficient for the setting of Gauss integral point coefficient, CellGaussPoints for the setting of all Gauss points in the background element and calculation of Jacob ratio, SuPPortDomain for the determination of single interpolation point support area, MLS_ShaPeFunc_2D for the calculation of shape function and derivative of single interpolation points, PointstiffessMatrix for the calculation of stiffness matrix of single Gauss integral point, EssentialBC for the implementation of penalty function   External program module, also called add-in program, is used to figure out the distribution of hydraulic pressure inside the existing fracture and simulate the fracture propagation for coal rock sample. It consists of several customized submodels: 3DChange2D for the transformation of three dimensional fracture processing, XChangeNode for the modification of original fracture plane node information and generation of new fracture plane, XInsertRNode for the deletion of node information in background grid cell for the existing fracture plane and node insertion in fracture area, XNODE for the acquisition and loading of fracture surface node and modification of boundary condition, FMONODE for the calculation of fracture morphology and orientation parameters, MFLOAD for the calculation of fluid pressure on fracture surface, and XLOAD for the calculation of loading force in the fracture surface based on Hillerborg model.

Program for computing hydraulic pressure inside the fracture
The flowchart of hydraulic pressure inside the fracture illustrated in Figure 9 is established based on the mathematical model of Equation (9). Two submodules of FMONODE and MFLOAD are combined together to complete the loading of fluid pressure.

Program for simulating hydraulic fracture propagation
According to EFG process, hydraulic fracture propagation based on Hillerborg model can be programmed as the flowchart illustrated in Figure 10. Five submodules of 3DChange2D, XChangeNode, XlnsertRNode XNODE, and XLOAD are involved during this process.

Numerical simulation example
After the coding and compiling of all modules, the program can be used to simulate hydraulic fracturing.

Preparing coal rock sample for simulation
For comparison with physical experimental results at laboratory later, an actual coal rock sample of approximately 28 · 28 · 33 mm [3] with /5 · 27 wellbore was designed to verify the developed mathematical model and EFG program as an example. The physical properties of this coal rock sample are as follows: elastic modulus parallel to bedding direction is E = 0.27 · 10 5 MPa, elastic modulus perpendicular to bedding direction is E = 4.3 · 10 5 MPa, Poisson ratio is c = 0.2, tensile strength parallel to bedding surface is f t = 0.3 MPa, tensile strength perpendicular to bedding plane is f t = 1.2 MPa, fracture toughness parallel to bedding direction is K IC = 0.8 MPa m 1/2 , and fracture toughness perpendicular to bedding direction is K IC = 2.4 MPa m 1/2 .

The Environmental Scanning Electron Microscope (ESEM) of Quanta 200 and Micron Level Computerized Tomography
Scanner (Micro-CT) of SkyScan1172 were used to detect the distribution and characteristics of existing fractures, including original fractures and artificial fractures resulting from drilling the hole, in the coal rock sample. After the detection, two representative primary fractures were found in the coal rock sample illustrated in Figure 11. One is transverse fracture (named No. 1) at the depth of 10 mm in the wellbore, its dip angle is 0.43 rad, azimuth angle of 0, crack length is 12 mm, and crack width is 0.07 mm; the other is longitudinal fracture (named No. 2) at the depth of 10.5 mm in the wellbore, its dip angle is 1.28 rad, azimuth angle is 0, crack length is 13 mm, and crack width is 0.05 mm. These initial parameters are shown in Table 1.

Geometric description of existing fractures
In order to intuitionally calculate and analyze the distribution of existing fractures, their shape and parameters were described as geometric model illustrated in Figure 12.

Numerical simulation using developed program
Within the solving domain of bedding plane of projection, 11 · 10 + 10 + 10 nodes were placed. Each of them After all setting, the numerical simulation was started in the developed program mentioned above. Figure 13 illustrates the simulating interface which is displaying the distribution and related parameters of hydraulic fractures with triaxial confining pressure under the action of 3 MPa water injection.

Results of numerical simulation
Various water injection pressures of 0-6 MPa at interval of 1 MPa were adjusted to simulating hydraulic fracture propagation on coal rock sample under the action of increasing pressure. The simulated parameters of the existing fractures     Tables 2-4 respectively. By cubic polynomial fitting in Matlab, data of the width, length, dip angle and azimuth angle changed with various hydraulic pressure in above tables were fitted to the curve graphs illustrated in Figure 14.

Analysis on numerical simulation results
Major findings from the analysis on the simulated data and result curves of hydraulic fracturing in coal rock can be concluded as follows.
(1) The expansion of hydraulic fractures is of step response on the water injection pressure since the existing fractures in coal rock sample wouldn't keep expanding in a certain water injection pressure, but rather stop at a certain stable state. Only when the water injection pressure increases to a certain higher situation, the expansion will continue to proceed again. (2) The expansion of hydraulic fractures is of inheritance since it is mainly based on the existing fractures and in the form of straight line without distortion normally. (3) For a given step hydraulic fracturing (this step is usually greater than that of the crack tip's fracture time of occurrence), the fracture propagation in the numerical simulation process is embedded in two stages. During the first stage, called fracture occurrence stage, a new long fracture will suddenly appear at the area without original fractures. During the second stage, called development stage, the fracture will propagate slowly in the water pressure loading process until the water pressure loading is ended or the fracture propagation is complete. It is clear that two stages were involved for the new fracture of No. 3 while only the second stage was for the existing fractures of No. 1 and No. 2.

Hydraulic fracturing experiments
For comparing with and verifying the numerical simulation results, the physical simulation experiments were conducted on the coal rock sample used in numerical simulation above in the laboratory.

The experimental procedure
The basic experimental procedure is outlined as follows and illustrated in Figure 15.  Step The first two steps have been done for the numerical simulation. The water injection is a similar process as many other studies. The following parts will mainly focus on the results and analysis.

Coal rock sample scanning after hydraulic fracturing
After each water injection, ESEM and Micro-CT were used to acquire the properties for the hydraulic fractures. Figure 16 illustrates the surface characteristics of three fractures scanned by ESEM after 3 MPa water injection experiment. For ease of human observation, all of them are magnified 28 times, 56 times and 112 times. Figure 17 illustrates the internal characteristics of three fractures scanned by Micro-CT after 3 MPa and 6 MPa water injection experiments.

Geometric description of hydraulic fractures
In the same way as in the numerical simulation, the shape and parameters of hydraulic fractures after water injection were described as geometric model illustrated in Figure 18 while the detailed parameters are listed in Tables 5 and 6.

Interpolating process for data before and after hydraulic fracturing experiments
By cubic spline interpolating in Matlab [28] for the data before and after hydraulic fracturing experiments, the curve graphs of the shape and orientation parameters of three fractures with water injection pressure in above tables were illustrated in Figure 19.

Analysis on physical experiment results
Major findings from the analysis on above images, data and curves before and after the experiment can be concluded as follows.
(1) The expansion of hydraulic fractures of the coal rock will mainly initiate from the existing fractures on the The results of CT scanning of the fractures after 6MPa-hydraulic-fracturing experiment The results of CT scanning of the fractures after 3MPa-hydraulic-fracturing experiment wall of wellbore and then intersect with the internal porosity and fractures to form the fracture network.
Less new fracture is able to be observed. (2) For coal rock with only one single existing fracture, either transverse or longitudinal, hydraulic fracturing propagation will be primarily in the form of straight line along the existing fracture without the occurrence of larger distortion.
(3) For coal rock with both transverse and longitudinal fractures coexisting but no intersection, various phenomena resulting from increasing hydraulic pressure will be appeared: (1) Neither fracture propagated for lower water injection pressure (less than 1 MPa).
(2) While water injection pressure is increased (1-3 MPa), the transverse fracture will propagate significantly but longitudinal one won't. (3) While water    injection pressure continues to be increased (3-6 MPa), both transverse and longitudinal fractures will propagate, and after they intersect with each other, the transverse propagation will be dominated and longitudinal one will cease basically. (4) When water injection pressure is up to a certain value (about 6 MPa), new longitudinal fracture will be prone to occur for lower three axial confining pressure (about 1 MPa) while new transverse fracture for higher three axial confining pressure (about 3 MPa).

Conclusion
This paper conducted the numerical simulation and physical experiments for the fracture propagation of hydraulic fracturing in the coal seam. For the numerical simulation, a kind of new geometric model of spring-splint was built for the numerical simulation using EFG method and was applied on an actual coal rock sample with a drilling hole as an example. For the physical simulation, the same one coal rock sample was used for the experiment as well as the advanced detection instruments of ESEM and Micro-CT were used for the nuanced observation on the surface of and inside the coal rock sample. Many findings are concluded from either numerical simulation or physical simulation experiments. Through comparison and analysis, the validity, feasibility and capability of spring-splint model and EFG method are illustrated by the reasonably good agreement between them. The results also show that since the development degree of transverse fractures (parallel to the bedding surface) is higher than that of longitudinal fractures (perpendicular to the bedding surface) in the coal seam, the extension of transverse fracture is dominated while the longitudinal crack extension is subsidiary during the process of hydraulic fracturing.
Based on the above results and analysis, the following major suggestions will be instrumental to some extent to the design of hydraulic fracturing process in actual production of CBM: 1. The jet hole should be placed at the direction abounding in transverse fractures in the wellbore to the greatest extent, which will facilitate full extension of transverse fractures with lower water injection pressure, thus improving the efficiency of CBM recovery as well as reducing the cost of hydraulic fracturing. 2. Correspondingly, the jet hole should avoid the direction with longitudinal fractures in the wellbore as far as possible. This will not only be against the CBM recovery efficiency with higher cost of higher water injection pressure, but also induce some accidents such as collapse of the roof resulting from excessive expansion of longitudinal fracture. 3. If the direction of jet hole in the wellbore has both transverse and longitudinal fractures in the trend of separation, water injection pressure should be reasonably adjusted in a range to facilitate the extension of transverse fractures and inhibit the extension of longitudinal fractures. 4. If the direction of jet hole in the wellbore has both transverse and longitudinal fractures in the trend of intersection, water injection pressure should be increased to facilitate the sufficient extension of transverse fractures without influencing the stability of the roof of coal reservoir. 5. For thicker coal seam, hydraulic fracturing should firstly be conducted in the upper part of the reservoir to produce new transverse fracture. Then a certain water injection pressure will be loaded in the deeper part for producing longitudinal fracture after sealing the jet hole in the upper part. In this way, transverse and longitudinal fractures will intersect to form fracture network, thus benefiting the CBM recovery efficiency. 10. In the case of the existence of larger fault resulting from tectonic movement nearby the wall of wellbore,  no jet hole should be arranged at this place regardless of any kind of fractures to avoid the extension of hydraulic fracture along the fault, which can result in relative slip of fault plane to induce roof unstability of coal reservoir.