Propagation of hydraulic fracture under the joint impact of bedding planes and natural fractures in shale reservoirs

The generation of the complex fracture network through hydraulic fracturing is a prerequisite for economic exploitation of shale gas. Bedding planes and natural fractures in shale reservoirs exert a direct impact on the propagation of hydraulic fracture and the formation of the complex fracture network. Based on theoretical analysis and flow‐stress‐damage coupled approach, the propagation of hydraulic fracture under the joint impact of the bedding planes and the natural fractures is studied. The results show the main factors that influence the fracture propagation are the horizontal principal stress difference, the minimum horizontal principal stress, the approaching angle of bedding plane (AABP), the intersection angle of natural fracture (IANF), and the cohesion of the reservoir matrix and the bedding plane. The bigger the horizontal principal stress difference is and the smaller AABP and the cohesion are, the more likely hydraulic fracture will propagate toward the direction of the maximum horizontal principle stress; the smaller AABP and the cohesion are, the more easily it will propagate along the bedding planes; AABP and IANF are the key factors that affect the fracture network formation and both of them bear a negative correlation on the principal stress difference for the formation of the fracture network; the smaller the stress difference and the cohesion difference between the bedding plane and the reservoir matrix are, and AABP is close to IANF and they are in the middle of its value range, the more easily the fracture network forms. The research results are of great significance to reveal the propagation law of hydraulic fracture and the formation of the complex fracture network in shale reservoirs and the efficient exploitation of shale gas.

technique for shale gas production in North America. [8][9][10] So far, the generation of complex fracture network through hydraulic fracturing on shale reservoir is a prerequisite for economic exploitation of shale gas. 11 However, shale reservoir that contains lots of bedding planes and natural fractures is different from other common rock (Figure 1). 12 The bedding planes and natural fractures in the shale reservoir exerts a direct impact on the propagation of hydraulic fracture as well as the formation of complex fracture network, affecting the economic and high-efficient shale gas production. Therefore, a deep understanding of influence factors and rules of the propagation of hydraulic fracture under the joint impact of the bedding planes and natural fractures in shale reservoir is crucial to the economic and high-efficient shale gas production.
So far, some researches have been made on the influence of bedding planes and natural fractures on hydraulic fracture propagation. Regarding the factors affecting hydraulic fracture propagation, Li et al 13 concluded through numerical simulation that the fracture propagation is primarily influenced by original crustal stress, bedding plane quality as well as the lithologic character. Lu et al 14 also drew a conclusion through the same means that horizontal principal stress difference, minimum horizontal principal stress and intersection angle between hydraulic fracture and stratigraphic interface which are the dominant factors that influence fracture propagation and the forms. With regard to the impact of bedding planes on propagation, Grasselli et al 15 conducted a research on the problem of crack initiation. Zhao et al 16 studied the effect of the thickness and strength of bedding planes on the fracture propagation. Men et al 17 found out that, with the increase in the bedding plane angles, the maximum principal stress and bedding planes dominantly control the initiation and propagation of the fracture crack. Zou et al 18 conducted a research by means of numerical simulation on how the bedding planes affect hydraulic fracture network propagation. Zou et al 19 analyzed the fracture network after the hydraulic fracture through CT scanning technique. Through research, Heng et al 20 found that the hydraulic fracture will ramify and sheer off while propagating vertically to bedding planes at weak bedding planes, and more importantly, it will further link up the weak bedding planes or natural fractures in the process of continuous spread, thus forming complex fracture networks. Hou et al 21 found in his physical simulation experiment of large-size triaxial hydraulic fracture that while the crack arrest, ramification, passing through and sheering off that happened to hydraulic fracture meeting weak bedding planes in the propagation process are the major factors that form complex fracture networks, the lots of weak planes themselves are the precondition for the fracture networks. As for the impact of natural fractures on the propagation, Zangeneh et al 22 carried out a study with the discrete element method. Blanton 23,24 concluded through experiment approaching angle between hydraulic fracture and natural fractures, and horizontal principal stress difference are the major factors that affect the propagation. Cheng et al 25 put forward a rule on judging whether the hydraulic fracture has penetrated the natural fractures based on his theoretic analysis and experiment study. Hou et al 26 noticed in his hydraulic fracture experiment based on large-size triaxial that a large intersection angle between the direction of the maximum horizontal principal stress and natural fractures and bedding planes, a small stress difference and a large shale elastic modulus will facilitate the forming of the fracture network. Taleghani 27 found in his simulation study by expanding finite element method that it is difficult to form the fracture network when intersection angle between natural factures and the direction of the maximum horizontal principal stress remains small, but it can be easy to grow into a complex one when the angle is large.
The work done above has made contributions to the study of the propagation of the hydraulic facture. However, there are few studies on the propagation of hydraulic fracture under the joint impact of bedding planes and natural fractures in shale reservoirs. Therefore, to better understand the factors and rules of the propagation of hydraulic fracture under the joint impact of the bedding planes and natural fractures in shale reservoirs, as well as the formation of the complex fracture network, based on elastic theory and fracture mechanics, a theoretical model of hydraulic fracture propagation in shale reservoirs with bedding planes and natural fractures is built. With the help of RFPA2D-Flow, a flow-stress-damage coupled numerical simulation software, the propagation of hydraulic fracture in shale reservoirs that contains bedding planes and natural fractures is Considering the vertical fracture along the maximum horizontal principal stress direction formed, and when hydraulic fracture is imposed on the deep reservoir, as shown in Figure 2A, the fracture height will not be considered when studying the horizontal propagation of the fracture. Instead, studying the horizontal propagation on purpose by sectioning the model into a 2D model, as is shown in Figure 2B. The hydraulic fracture under the impact of far-field stress propagates along the direction of the maximum horizontal principal stress. Because of bedding planes of shale, three situations may arise when the propagating fracture meets the bedding plane: a. the fracture will directly cross the bedding plane, forming a single fracture; b. the fracture will propagate along the bedding plane, forming a single fracture as well; c. the fracture will ramify when propagating, forming a complex fracture that crosses and moves along the bedding plane at the same time.
In the far-field, the expanding hydraulic crack is affected by the maximum horizontal principal stress, the minimum horizontal principal stress and the internal water pressure inside the fracture ( Figure 2B). The normal stress and shear stress on the wall of the hydraulic crack are: where n is the normal stress (MPa), p is the water pressure (MPa), h is the minimum horizontal principal stress (MPa), and is the shear stress (MPa).
As can be inferred from the theory of stress intensity factor, the stress intensity factor of crack tip can be put as: where K I is the stress intensity factor, and a is the fracture halflength (m).
For the linear elastic material, Irwin proposed in 1957 a fracture criterion of stress intensity factor, which is, for fracture type I, the fracture will propagate when stress intensity factor (K I ) reaches its critical value (K Ic ), namely, K I = K Ic . According to fracture mechanics, K Ic is: where K Ic is the critical value of stress intensity factor, E is the Young's modulus (MPa), is the surface energy (J/m 2 ), and v is the Poisson's ratio.
Based on fracture propagation theory, other conditions alike, Griffith fracture needs the least fluid pressure. Therefore, the research on fracture propagation based on the fracture mechanics usually assumes the form of the fracture as Griffith fracture, thus according to Equations (1)-(3), the required critical fluid pressure (p 1 net ) for propagate of hydraulic fracture in shale reservoir before encountering bedding plane or natural fracture is: When the hydraulic fracture meets bedding plane in its propagation in shale reservoir, as shown in Figure 2C, in which is the approaching angle of bedding plane (AABP, (1) namely, the approaching angle between the bedding plane and the direction of the maximum horizontal principal stress), then transform the fracture tip stress field into the polar coordinate system where rr is the radial stress (MPa), is the normal stress (MPa), r is the tangential stress (MPa), and θ is the approaching angle of bedding plane (°).
According to 2D elasticity theory, with the impact of farfield stress, the normal stress and the shear stress on the bedding plane can be calculated by where ∞ n is the far-field normal stress (MPa), ∞ is the far-field shear stress (MPa), and H is the maximum horizontal principal stress (MPa). Therefore, with the joint impact of the far-field stress and the internal water pressure inside the hydraulic fracture, the normal stress and shear stress on the bedding plane can be calculated by In theory, when the bedding plane takes tensional damage at n < 0, the needed net pressure (p 2 net ) inside the hydraulic fracture is: To intuitively and quantitatively analyze the effect of the horizontal principal stress and AABP on the fracture propagating in the bedding planes, assuming H = 20 MPa, √ a∕2r = 10 base on the field experience of hydraulic fracturing tests. When the fracture propagates to the bedding planes and Equation (8) achieves the same value on both sides, the net fluid pressure value is equal to that of the hook face in Figure 3. It can be seen that with the increase in the horizontal principal stress difference, or, to be more specific, the decrease in the minimum horizontal principal stress, the required net pressure for the tensional damage on the bedding plane almost takes a linear decrease, but with limited degree. Even AABP reaches zero and the minimum horizontal principal stress is 10 MPa, the change caused by the minimum horizontal principal stress to the net pressure only gets to 1.1 MPa. As AABP increases, contribution of decrease in the minimum horizontal principal stress is gradually approaching zero. In fact, the maximum impact factor on the fluid net pressure in Equation (8) is AABP instead of the horizontal principal stress difference. Decrease in AABP and rapid reduction of the net pressure makes it easier for tensional damage to happen to the bedding planes. This explains why opening damage can easily occur along the bedding plane when AABP is small.
Theoretically, when shear damage to the bedding planes happens at > c, the required fluid net pressure (p 3 net ) is: where c is the cohesion of the bedding planes (MPa). Similar to Equation (8), in Equation (9), the fluid net pressure at the time when both sides are equal under different cohesion of bedding planes are the pressure values of the curved surface in Figure 4. It shows that the horizontal principal stress difference contributes almost nothing to the fracture propagation in such a damage mode. Propagation now is mainly controlled by AABP and the cohesion of the bedding plane. With the strength of the bedding planes dropping, the needed fluid net pressure for fracture propagation is reduced as well. As AABP increases, the net pressure will fall rapidly, but with the continuous increase in the angle, the pace of the decrease in the net pressure will gradually slow down. When the angle passes 30°, its change nearly does not make any contribution to the fall of the net pressure needed for fracture propagation. In other words, when AABP remains small, the needed net pressure for propagation is primarily in the control of AABP. The smaller the angle is, the harder it is for the shear failure to take place, and the bigger the angle is, the easier it is for the shear failure to take place. Otherwise, it is controlled by the cohesion of the bedding plane. This explains why shear failure hardly happens along the bedding plane when AABP is small, and hydraulic fracture directly cross the bedding plane is prone to shear failure when the angle is big. At the same time, this also explains the impact brought by the cohesion to complexity of the fracture penetrating the bedding plane.

| Propagation of the hydraulic fracture under the influence of natural fracture
When propagating along the maximum horizontal principal stress direction and meeting the closed natural fractures in the far-field, the hydraulic fracture may intersect with the natural fracture by extending along the bedding plane or cross the plane. As indicated in Figure 5, it is assumed that is the intersection angle of natural fracture (IANF, namely, the intersection angle between the natural fracture and the direction of the maximum horizontal principal stress). Based on the above theoretical analysis, when the natural fractures remain shut, the hydraulic fracture will directly cross them and then continue to extend. As shown in Equation (7), the natural fracture will suffer tension damage at the time of the normal stress of the natural facture n < 0, namely, ∞ n < , and the required fluid net pressure (p 4 net ) for hydraulic fractures to propagate along the natural fracture is: where is the intersection angle of natural fracture (°).
When the natural fractures open, the fluid is filled with the hydraulic fracture and makes the fracture expand with the fluid pressure. At the time of a shear damage, according to Equation (9), the needed fluid net pressure (p 5 net ) is: where c 0 is the cohesion of shale reservoir matrix (MPa). As can known from Equations (10) and (11), when = and c = c 0 in Equations (8) and (9), Equations (10) and (11) can be obtained. Therefore, the parameters in Equations (10) and (11) have the same rules regarding their impact on the hydraulic fracture propagation as those in Equations (8) and (9). Only at the time of = will the shear failure occur to the bedding plane first.

| Analysis of the fracture propagation direction and the fracture network formation
1. If min p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net = p 1 net and it differs greatly from min p 2 net ,p 3 net ,p 4 net ,p 5 net , the hydraulic fracture will primarily move along the maximum horizontal principal stress direction and form a single fracture by passing through the bedding plane and the natural fractures, in the shale matrix with relatively homogeneous mechanical properties. However, this propagation pattern only happens when the shale matrix mechanical properties are weak and the minimum horizontal principal stress small. 2. When min p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net = p 2 net and it has a huge difference from min p 1 net ,p 3 net ,p 4 net ,p 5 net , the hydraulic fracture will suffer tension damage along the bedding plane, forming a single fracture. And this propagation pattern mainly occurs when the minimum horizontal principal stress is small, the horizontal principal stress difference big, and AABP small, as shown in Figure 3. 3. If min p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net = p 3 net and it differs greatly from min p 1 net ,p 2 net ,p 4 net ,p 5 net , the hydraulic fracture will suffer shear damage when it meets the bedding planes, going through them and forming a single fracture. This propagation pattern mainly occurs at the time of big AABP and small the cohesion of bedding planes, as can be seen in Figure 4. 4. When min p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net = p 4 net and it has a huge difference from min p 1 net ,p 2 net ,p 3 net ,p 5 net , the fracturing fluid will fill the natural fractures when the hydraulic fracture meets the latter, which expands the natural fractures with the fluid pressure. At this time the hydraulic fracture propagates along the direction of the least needed water pressure. When it again intersects the bedding plane, it may, depending on AABP and the cohesion, form a single fracture or a crossover fracture. And this propagation pattern mainly occurs when the minimum horizontal principal stress is small, horizontal principal stress difference big, and IANF small. 5. If min p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net = p 5 net and it differs greatly from min p 1 net ,p 2 net ,p 3 net ,p 4 net , the hydraulic fracture will propagate beyond the natural fracture and suffer shear failure on the latter's surface, thus forming a crossover fracture. This fracture pattern usually happens when horizontal principal stress difference is big, shale matrix cohesion small and IANF big. 6. Propagation toward all directions is possible for the fractures when p 1 net ≈ p 2 net ≈ p 3 net ≈ p 4 net ≈ p 5 net . Particularly, the propagation directions are random, and when they reach bedding planes, they may deflexion, making it easy for them to form a complex fracture network. The equation of p 1 net = p 2 net = p 3 net = p 4 net = p 5 net lead to = and c = c 0 . When H = 20 MPa, c = 12 MPa and √ a∕2r = 10 are taken for an example, AABP and IANF will almost show a negative linear correlation on the horizontal principal stress difference. As shown in Figure 6, with the increase in the horizontal principal stress difference, the needed and for the formation of the complex fracture network will gradually decrease and increase with the former's decrease, but still in the middle of the value range of and . If max p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net − min p 1 net ,p 2 net ,p 3 net ,p 4 net ,p 5 net ≤ 1 MPa , then ≈ and c ≈ c 0 . This is when the condition for easily forming a complex fracture network falls in the area between the two approximately paralleled straight lines up and down in Figure 6. It can be seen that it is possible to form a complex fracture network under any horizontal principal stress differences, which leads to the conclusion that the stress difference is not a precondition for a complex fracture network, but it will affect the influence of other factors on the fracture propagation. The forming of the fracture network is very fastidious to and . As for what happens in Figure 6, with a certain horizontal principal stress difference, the value range for and are very limited regarding the formation of a complex fracture network. No network will emerge if the range is exceeded. Therefore, and are the key factors for forming a complex fracture network.

| Method for hydraulic fracturing simulation
Rocks that contain a large number of geological structures like joints and cracks are a complex heterogeneous material, before reaching the maximum carrying capacity and losing its stability under the water pressure, the rock will go through a series of damage processes such as original fracture development, microfissure conceiving, growing and interconnecting. Meanwhile, the conceiving, extending and linking of the internal microfissures inside the rock increase the water diversion paths as well as the rock's permeability that can be compromised by the water proof area caused by the partial compressive deformation zone of the rock. It is because of the above mentioned complex characteristics contained in the process of coupling and deformation of the water and rocks that makes it difficult to build a flow-stress-damage coupling model. The classical rock constitutive theory cannot express the whole process of the rock deformation and damage, nor can it be simply used for the research on the failure mechanism. As the research on this mechanism, the idea that degradation of mechanical properties of the rock is caused by internal structure damage and crack growth under the impact of external forces has been gradually accepted, from the point of rock microstructure. Tang, 28-30 a professor from rock failure and unstability center of Northeastern University of  China, adheres to this idea and leads a team in developing a software called RFPA2D-Flow for analyzing the damage and fracture process of the rocks. It is used to effectively analyze the discontinuous deformation and flow-stress-damage coupling mechanism and simulate the flow-stress coupling as well as the permeability evolvement rules in the process of rock failure and crack extension.
Considering the rock mineral composition, cement mechanical property differences and other flaws in the rock itself like original microcracks, the obvious inhomogeneity of the rock in the microscopic physical properties is inevitable. RFPA2D-Flow uses Weibull power function of threshold value to describe the space distribution characteristics of the mechanical properties of the material, and elasticity finite element method to analyze the stress so as to have an analytical calculation of the stress field and displacement field of the object. Usually, the maximum tensile stress (or strain) criterion and Mohr-Coulomb criterion are seen as the threshold value of the damage and destruction. When stress (or strain) of the crack unit reaches the maximum tension stress (or tension strain) and meets the failure conditions of Mohr-Coulomb, the unit will start to suffer initial damage and destruction of tension and shear. In the uniform deformation field, besides the stress concentration caused by such structure factors as load inequality, the intensity unevenness of the microscopic elements is the more important factor that leads the microcrack to develop. Thus, for the problem of inhomogeneity, the numerical simulation analysis method based on statistical theory is considered an effective method. The governing equations and Weibull power function that describes space distribution characteristics of mechanic properties of materials used in RFPA2D-Flow are as follows: Function for space distribution characteristics of mechanic properties of materials: where ( ) is the statistical distribution density of , and 0 are the mechanical property parameters of material medium element and its average value, and m is the property parameter of Weibull power function, which reflects the uniformity of the material medium.
Equilibrium equation: where ij,j is the partial derivative of stress tensor on j, is the density (kg/m 3 ), and X j is the component of body force in the j direction. Geometric equation: where ij is the strain tensor, u i,j is the partial derivative of the component of displacement in i direction on j, u j,i is the partial derivative of the component of displacement in j direction on i, v is the volumetric strain, and xx , yy , and zz are the normal strain in the x, y and z directions, respectively.
Constitutive equation: where ′ ij is the derivative of stress tensor, ij is the Kronecker constant, is the pore water pressure coefficient, is the Lame constant, G is the shear modulus (MPa).
Seepage equation: where K is the coefficient of permeability (m/d), ∇ 2 is the Laplace operator, and Q is the Biot's constant.

Seepage-stress coupling equation:
where is the mutation coefficient, K 0 is the initial coefficient of permeability, and is the coupling parameter.

| Numerical model and numerical simulation scheme
The numerical simulation of the propagation of hydraulic fracture under the joint impact of the bedding planes and the natural fractures is conducted by RFPA2D-Flow. The model size has been set as 10 m × 10 m and divided into 500 × 500 = 250 000 units. The horizontal principal stress has been simplified as uniform load imposed on the model boundaries. The inclined strip belt in the model stand for the bedding planes of the shale reservoir, whose physical and mechanical parameters are given according to the shale outcrop in the shale specimen area. In the middle of the model was dug an oval whose long axis is 1 m and minor axis is 0.03 m, which represent the propagating hydraulic fracture whose internal water pressure is inflicted based on boundary conditions. In addition, a closed fracture is set in front of the extending direction of the hydraulic fracture and intersects with the bedding planes, which is regarded as a natural fracture, as shown in Figure 7.
The above analysis shows that the horizontal principal stress difference, the minimum horizontal principal stress, AABP, IANF, and the cohesion of shale matrix and bedding planes all have an important impact on the propagation pattern and morphological characteristics of the fractures. Based on the field conditions of the present shale hydraulic fracturing and the literature search, a numerical simulation scheme is designed as shown in Table 1. The physical and mechanical parameters of shale reservoir in the numerical simulation are shown in Table 2. 31

| Impact of the horizontal principal stress difference
The impact of different horizontal principal stress differences on the hydraulic fracture propagation is shown in Figure 8, and it is notable. When the horizontal principal stress difference is small, the hydraulic fracture upon meeting the bedding plane will directly cross the bedding plane because of the big AABP. Meanwhile, some of the fractures still extend along the bedding plane. With the gradual increase in the horizontal principal stress difference, the hydraulic fracture gradually turn to the direction of the maximum horizontal principal stress in its propagation, and finally travel totally along it in a way similar to a straight line. Due to = = 55 • in the four groups in Figure 8, which is within the scope of the angle where a complex fracture network can easily take shape under the low stress difference. Thus, cross fractures arise in 1#, 2# of low stress difference, which can easily lead to forming a complex fracture network. However, with the increase in the horizontal principal stress difference, for the same value of and a single fracture will build up, which is no longer good for a complex fracture network formation, as shown in 4#.

| Impact of the approaching angle of bedding plane
The impact of different the approaching angle of bedding plane (AABP) on the hydraulic fracture propagation is shown in Figure 9. As can be seen, different AABP have more notable impact than that of the different horizontal principal stress differences. In Figure 9, when the AABP is small, the hydraulic fracture upon meeting the bedding plane will cross it and suffer tension damage, forming a single fracture. With the increase in the AABP, the hydraulic fracture propagation will gradually turn to follow the maximum horizontal principal stress direction from its previous moving direction, until it completely moves along the turned direction (7#), which still forms a single fracture. In Figure 9, the four groups of horizontal principal stress differences are set as 4 MPa, only 2# in which AABP and horizontal principal stress difference are matched has a cross fracture, while the AABP is small (5#) or big (7#) can only form a single fracture. This demonstrates the importance of AABP to the hydraulic fracture propagation as well as the forming of a complex network.

| Impact of the intersection angle of natural fracture
It is essential to study the impact of the very existence of the natural fractures on the hydraulic fracture propagation before to study the intersection angle of natural fracture (IANF) as it complicates the underground stress field. The concentrated stress at the tip of the natural fractures induces the hydraulic fracture to propagate toward the tip. Nevertheless, in order for the hydraulic fracture to connect the tip, it must go through the bedding plane, the weak plane of the shale. Thus the hydraulic fracture usually connect the tip through the bedding plane, leading to a partial fracture networking close to the shaft, which can be seen in the result of the numerical simulation in this study.
As for the impact of the IANF on the hydraulic fracture propagation, Figure 10 gives us a good description. When the IANF is small, the hydraulic facture will basically propagate toward the maximum horizontal principal stress direction from natural fracture tip after meeting natural fracture, thus forming a single fracture. As the IANF increases, around 55°, a complex fracture network will take shape with the suitable horizontal principal stress difference and the IANF. However, when the IANF gradually magnifies and exceeds the allowed scope needed for the forming of the complex fracture network, a single fracture will still be produced along the maximum horizontal principal stress direction in its propagation. This shows the importance of the IANF to the hydraulic fracture propagation as well as the forming of a complex network.
Based on the joint impacts of AABP and IANF on the complex fracture network formation and the conclusion reached above that the fracture network can easily take form in the situation where the two angles are the same, a numerical simulation is carried out when AABP and IANF are equal. As shown in Figure 11 the horizontal principal stress difference in the four groups is 4 MPa. With this difference, when the two angles are too small (11#) or big (13#), they are not beneficial for the formation of a complex fracture network. When the two angles are relatively small, the hydraulic fractures will propagate extensionally along the bedding planes, forming a single fracture. With the gradual increase of them, the hydraulic fracture will propagate along both the bedding plane and the maximum horizontal principal stress direction from its previous travel along the plane alone, forming a cross fracture. With the continuous increase in the two angles, the cross fracture will shrink into a single fracture along the maximum horizontal principal stress direction.

| Impact of the cohesion of bedding plane
As shown in Figure 12, = = 55 • in the four groups is still within the scope of angle in which a complex fracture network can easily take shape. However, closing to the shaft in 14# and 15# there is basically no complex fracture network, except the partial networking, while in 2# and 16# obvious cross complex fractures take shape, which explains the impact of different bedding plane cohesion. Therefore, it can be concluded that the cohesion has a significant impact on the forming of the hydraulic fracture in shale reservoirs. Based on previous analysis of complex fracture network formation, the formation of the complex fracture network can be meet at ≈ and c = c 0 . As the cohesion of the bedding plane (c ) is usually smaller than that of shale matrix (c 0 ), only when c is close to c 0 , the complex fracture network can be easily  formed (2#, 16#), or it is not good for the fracture network to take shape (14#, 15#).

| DISCUSSION
According to the theoretical analysis, the main factors that influence the fracture propagation are the horizontal principal stress difference, the minimum horizontal principal stress, AABP, IANF, the cohesion of shale reservoir matrix and bedding planes. The bigger the horizontal principal stress difference is, the more likely the hydraulic fracture will propagate toward the direction of the maximum horizontal principle stress, forming a single fracture and not beneficial for a complex fracture network to take shape. Due to the certain randomness of the fracture propagation, small horizontal principal stress difference makes it easier for the formation of the complex facture network. When AABP remains small, the hydraulic fracture upon meeting bedding planes will propagate along them and suffer tension damage. With the increase in AABP, the propagation direction of the fracture totally along the bedding plane at first will gradually turn to the maximum horizontal principal stress direction. During this process, for a certain level of horizontal principal stress difference, the complex fracture network can only take place with an appropriate AABP. Therefore, the formation of the complex facture network demands highly regarding AABP. Intersection angle of natural fracture has the same impact on the hydraulic fracture propagation rules as that of AABP. The two angles, when close in value and within the middle of its value range, are of greater benefit to the complex fracture network. The horizontal principal stress difference bears a negative correlation on AABP and IANF required for the formation of the fracture network. That is, the greater the horizontal principal stress difference is, the smaller the two angles are, which is the same as the rules disclosed in Figure 6. The very existence of natural fractures and bedding planes also complicate the underground stress field and lead to the partial fracture networking close to the shaft, causing the hydraulic fracture to deflect, turn and ramify, and thus good for the forming of the fracture network. Therefore, the bedding planes and natural fractures are very important to forming complex fracture network, AABP and IANF are the key factors that affect the fracture network formation. The cohesion of shale matrix and bedding planes also have an impact on the fracture propagation, which is mainly shown that the cohesion of the bedding plane is usually smaller than that of the shale matrix. When the two differ greatly, the hydraulic fracture will easily propagate along the bedding plane, thus not beneficial for the forming of the fracture network. When they differ little, the opposite takes place. Besides, study of the joint impact of bedding planes and natural fractures of shale reservoir on hydraulic fracture propagation also have a guidance for situations of rock mass rich in joints and fractures, and it is valuable for the production of such unconventional gases as coalbed methane and tight sandstone gas as well.